From 4e659d4ed54425310ecfd7a17b3a09f55b255c36 Mon Sep 17 00:00:00 2001 From: Maximilian Linhoff Date: Wed, 26 Oct 2022 12:36:26 +0200 Subject: [PATCH 1/2] Add option to not parse blocks --- corsikaio/file.py | 39 ++++++++++++++++++++++++++++++--------- 1 file changed, 30 insertions(+), 9 deletions(-) diff --git a/corsikaio/file.py b/corsikaio/file.py index 2d66d32..3ba69da 100644 --- a/corsikaio/file.py +++ b/corsikaio/file.py @@ -24,12 +24,19 @@ ['header', 'particles', 'longitudinal', 'end'] ) +def _to_floatarray(block): + return np.frombuffer(block, dtype=np.float32) + class CorsikaFile: + """ + A file to iterate over events in a CORSIKA binary output file + """ - def __init__(self, path): + def __init__(self, path, parse_blocks=True): self.EventClass = Event + self.parse_blocks = parse_blocks self._buffer_size = read_buffer_size(path) self._f = open_compressed(path) @@ -72,7 +79,10 @@ def __next__(self): if block[:4] != b'EVTH': raise IOError('EVTH block expected but found {}'.format(block[:4])) - event_header = parse_event_header(block)[0] + if self.parse_blocks: + event_header = parse_event_header(block)[0] + else: + event_header = _to_floatarray(block) data_bytes = bytearray() long_bytes = bytearray() @@ -87,9 +97,14 @@ def __next__(self): block = self.read_block() - event_end = parse_event_end(block)[0] - data = self.parse_data_blocks(data_bytes) - longitudinal = parse_longitudinal(long_bytes) + if self.parse_blocks: + event_end = parse_event_end(block)[0] + data = self.parse_data_blocks(data_bytes) + longitudinal = parse_longitudinal(long_bytes) + else: + event_end = _to_floatarray(block) + data = _to_floatarray(data_bytes).reshape(-1, 7) + longitudinal = _to_floatarray(long_bytes) return self.EventClass(event_header, data, longitudinal, event_end) @@ -150,13 +165,16 @@ def close(self): class CorsikaCherenkovFile(CorsikaFile): - def __init__(self, path, mmcs=False): - super().__init__(path) + def __init__(self, path, mmcs=False, **kwargs): + super().__init__(path, **kwargs) self.EventClass = PhotonEvent self.mmcs = mmcs def parse_data_blocks(self, data_bytes): + if not self.parse_blocks: + return _to_floatarray(data_bytes) + photons = parse_cherenkov_photons(data_bytes) if not self.mmcs: return photons @@ -175,9 +193,12 @@ def parse_data_blocks(self, data_bytes): class CorsikaParticleFile(CorsikaFile): - def __init__(self, path): - super().__init__(path) + def __init__(self, path, **kwargs): + super().__init__(path, **kwargs) self.EventClass = ParticleEvent def parse_data_blocks(self, data_bytes): + if not self.parse_blocks: + return _to_floatarray(data_bytes) + return parse_particle_data(data_bytes) From aae021f06e0524e661c8b7ec0b69b91b79d6e8e2 Mon Sep 17 00:00:00 2001 From: Maximilian Linhoff Date: Wed, 15 Feb 2023 13:56:22 +0100 Subject: [PATCH 2/2] Add test for parse_blocks=False, remove unneeded code --- corsikaio/file.py | 6 ------ tests/test_file.py | 17 +++++++++++++++++ 2 files changed, 17 insertions(+), 6 deletions(-) diff --git a/corsikaio/file.py b/corsikaio/file.py index 3ba69da..65f25a1 100644 --- a/corsikaio/file.py +++ b/corsikaio/file.py @@ -172,9 +172,6 @@ def __init__(self, path, mmcs=False, **kwargs): self.mmcs = mmcs def parse_data_blocks(self, data_bytes): - if not self.parse_blocks: - return _to_floatarray(data_bytes) - photons = parse_cherenkov_photons(data_bytes) if not self.mmcs: return photons @@ -198,7 +195,4 @@ def __init__(self, path, **kwargs): self.EventClass = ParticleEvent def parse_data_blocks(self, data_bytes): - if not self.parse_blocks: - return _to_floatarray(data_bytes) - return parse_particle_data(data_bytes) diff --git a/tests/test_file.py b/tests/test_file.py index cdfddb1..42b5c81 100644 --- a/tests/test_file.py +++ b/tests/test_file.py @@ -1,4 +1,5 @@ import pytest +import numpy as np def test_version(): @@ -71,6 +72,22 @@ def test_particle_longi(): assert e.header['event_number'] == 1 +def test_particle_no_parse(): + from corsikaio import CorsikaParticleFile + + with CorsikaParticleFile('tests/resources/corsika757_particle', parse_blocks=False) as f: + n_read = 0 + for e in f: + n_read += 1 + # second entry in header is event_number + assert e.header[1] == n_read + assert e.header.dtype == np.float32 + assert len(e.header) == 273 + assert e.particles.size % 273 == 0 + assert n_read == 10 + + + def test_truncated(tmp_path): '''Test we raise a meaningful error for a truncated file