Skip to content

Commit

Permalink
Merge pull request #207 from cta-observatory/time_array
Browse files Browse the repository at this point in the history
Store photoelectron arrival information in two 1d arrays
  • Loading branch information
watsonjj authored May 25, 2020
2 parents d8c1abc + 777f112 commit a9829ef
Show file tree
Hide file tree
Showing 4 changed files with 30 additions and 34 deletions.
3 changes: 3 additions & 0 deletions .gitignore
Original file line number Diff line number Diff line change
Expand Up @@ -129,3 +129,6 @@ dmypy.json

# Pyre type checker
.pyre/

# IDE
.idea
6 changes: 3 additions & 3 deletions eventio/iact/objects.py
Original file line number Diff line number Diff line change
Expand Up @@ -314,7 +314,7 @@ def parse(self):
data = self.read()

pe.update(PhotoElectrons.parse_1208(
data, pe['n_pixels'], pe['non_empty'], self.header.version, flags
data, pe['n_pixels'], pe['non_empty'], self.header.version, flags, pe['n_pe']
))

return pe
Expand Down Expand Up @@ -401,8 +401,8 @@ def parse(self):
n_strings = read_int(self)
input_card = bytearray()
for i in range(n_strings):
input_card.extend(read_string(self))
input_card.append(ord('\n'))
input_card += read_string(self)
input_card += b'\n'
return input_card


Expand Down
43 changes: 20 additions & 23 deletions eventio/var_int.pyx
Original file line number Diff line number Diff line change
Expand Up @@ -363,27 +363,26 @@ def parse_1208(
uint32_t n_pixels,
uint32_t nonempty,
uint32_t version,
uint32_t flags
uint32_t flags,
uint32_t total_n_pe
):

cdef uint32_t length

cdef cnp.npy_intp[1] shape = [n_pixels]
cdef cnp.ndarray[float, ndim=1] photoelectrons = cnp.PyArray_ZEROS(1, shape, cnp.NPY_FLOAT32, False)
cdef cnp.ndarray[float, ndim=1] mean_time = cnp.PyArray_SimpleNew(1, shape, cnp.NPY_FLOAT32)
cnp.PyArray_FillWithScalar(mean_time, NAN)
cdef cnp.npy_intp[1] pixel_shape = [n_pixels]
cdef cnp.ndarray[float, ndim=1] photoelectrons = cnp.PyArray_ZEROS(1, pixel_shape, cnp.NPY_FLOAT32, False)
cdef cnp.ndarray[int32_t, ndim=1] photon_counts = None

cdef list time = [[] for _ in range(n_pixels)]
cdef list amplitude
cdef cnp.npy_intp[1] pe_shape = [total_n_pe]
cdef cnp.ndarray[uint32_t, ndim=1] pixel_id = cnp.PyArray_ZEROS(1, pe_shape, cnp.NPY_UINT32, False)
cdef cnp.ndarray[float, ndim=1] time = cnp.PyArray_ZEROS(1, pe_shape, cnp.NPY_FLOAT32, False)
cdef cnp.ndarray[float, ndim=1] amplitude = cnp.PyArray_ZEROS(1, pe_shape, cnp.NPY_FLOAT32, False)

cdef bint has_amplitudes = flags & 1
cdef bint has_photon_counts = flags & 4

if has_amplitudes:
amplitude = [[] for _ in range(n_pixels)]
cdef bint has_photons = flags & 4

cdef uint64_t pos = 0
cdef uint64_t i_pe = 0
cdef uint32_t i
cdef int32_t j
cdef dict result = {}
Expand All @@ -403,31 +402,29 @@ def parse_1208(
pos += 4

photoelectrons[pix_id] = n_pe
if n_pe > 0:
mean_time[pix_id] = 0.0

for j in range(n_pe):
t = (<float*> &data[pos])[0]
time[pix_id].append(t)
mean_time[pix_id] += t
pixel_id[i_pe + j] = pix_id
time[i_pe + j] = t
pos += 4

mean_time[pix_id] /= n_pe

if has_amplitudes:
for j in range(n_pe):
amplitude[pix_id].append((<float*> &data[pos])[0])
amplitude[i_pe + j] = (<float*> &data[pos])[0]
pos += 4

i_pe += n_pe

result['photoelectrons'] = photoelectrons
result['pixel_id'] = pixel_id
result['time'] = time
result['mean_time'] = mean_time

if has_amplitudes:
result['amplitude'] = amplitude

if has_photon_counts:
photon_counts = cnp.PyArray_ZEROS(1, shape, cnp.NPY_INT32, False)
if has_photons:
photons = cnp.PyArray_ZEROS(1, pixel_shape, cnp.NPY_INT32, False)

nonempty = (<int32_t*> &data[pos])[0]
pos += 4
Expand All @@ -436,10 +433,10 @@ def parse_1208(
pix_id = (<int16_t*> &data[pos])[0]
pos += 2

photon_counts[pix_id] = (<int32_t*> &data[pos])[0]
photons[pix_id] = (<int32_t*> &data[pos])[0]
pos += 4

result['photon_counts'] = photon_counts
result['photons'] = photons

return result

Expand Down
12 changes: 4 additions & 8 deletions tests/iact/test_iact_objects.py
Original file line number Diff line number Diff line change
Expand Up @@ -24,16 +24,12 @@ def test_photo_electrons():
assert 0 <= data['non_empty'] <= data['n_pixels']
assert len(data['photoelectrons']) == data['n_pixels']
assert data['photoelectrons'].sum() == data['n_pe']
assert len(data['time']) == 2368
assert sum(len(pixel) for pixel in data['time']) == data['n_pe']
assert len(data['pixel_id']) == data['n_pe']
assert len(data['time']) == data['n_pe']

# times should be within 200 nanoseconds
assert all(0 <= t <= 200 for pixel in data['time'] for t in pixel)
mean_time = np.array([
np.mean(t) if len(t) > 0 else np.nan
for t in data['time']
])
assert np.allclose(data['mean_time'], mean_time, equal_nan=True)
assert np.all(0 <= data['time'])
assert np.all(data['time'] <= 200)

not_read = o.read()
assert len(not_read) == 0 or all(b == 0 for b in not_read)
Expand Down

0 comments on commit a9829ef

Please sign in to comment.