Skip to content

Commit

Permalink
mv tof sinogram to lm generation to separate function
Browse files Browse the repository at this point in the history
  • Loading branch information
gschramm committed Nov 16, 2023
1 parent 2f07cc5 commit e21c7dd
Show file tree
Hide file tree
Showing 3 changed files with 187 additions and 127 deletions.
50 changes: 9 additions & 41 deletions python/00_simulate_lm_data.py
Original file line number Diff line number Diff line change
Expand Up @@ -6,6 +6,7 @@
import prd
import parallelproj
import parallelproj_utils
import utils
import array_api_compat.numpy as np
import matplotlib.pyplot as plt
from array_api_compat import to_device
Expand Down Expand Up @@ -84,6 +85,9 @@
)
projector.tof = True # set this to True to get a time of flight projector

# repeat number of TOF bin times here
num_tof_bins = projector.tof_parameters.num_tofbins

# forward project the image
noise_free_sinogram = projector(img)

Expand All @@ -101,50 +105,14 @@
noisy_sinogram = xp.asarray(
np.random.poisson(np.asarray(to_device(noise_free_sinogram, "cpu"))), device=dev
)
# ravel the noisy sinogram and the detector start and end "index" sinograms
noisy_sinogram = xp.reshape(noisy_sinogram, (noisy_sinogram.size,))

# get the two dimensional indices of all sinogram bins
start_mods, end_mods, start_inds, end_inds = lor_descriptor.get_lor_indices()

# generate two sinograms that contain the linearized detector start and end indices
sino_det_start_index = (
lor_descriptor.scanner.num_lor_endpoints_per_module[0] * start_mods + start_inds
)
sino_det_end_index = (
lor_descriptor.scanner.num_lor_endpoints_per_module[0] * end_mods + end_inds
)

# repeat number of TOF bin times here
num_tof_bins = projector.tof_parameters.num_tofbins
# -------------------------------------------------------------------------------------
# -------------------------------------------------------------------------------------
# -------------------------------------------------------------------------------------

sino_det_start_index = xp.reshape(
xp.stack(num_tof_bins * [sino_det_start_index], axis=-1),
sino_det_start_index.size * num_tof_bins,
event_det_id_1, event_det_id_2, event_tof_bin = utils.noisy_tof_sinogram_to_lm(
noisy_sinogram, lor_descriptor, xp, dev
)

sino_det_end_index = xp.reshape(
xp.stack(num_tof_bins * [sino_det_end_index], axis=-1),
sino_det_end_index.size * num_tof_bins,
)

# ----------------------------------------------------------------------------
# ----------------------------------------------------------------------------
# --- convert the index sinograms in to listmode data ------------------------
# ----------------------------------------------------------------------------
# ----------------------------------------------------------------------------

# generate listmode data from the noisy sinogram
event_sino_inds = np.repeat(np.arange(noisy_sinogram.shape[0]), noisy_sinogram)
# shuffle the event sinogram indices
np.random.shuffle(event_sino_inds)
# convert event sino indices to xp array
event_sino_inds = xp.asarray(event_sino_inds, device=dev)

event_det_id_1 = xp.take(sino_det_start_index, event_sino_inds)
event_det_id_2 = xp.take(sino_det_end_index, event_sino_inds)
event_tof_bin = event_sino_inds % num_tof_bins

print(f"number of simulated events: {event_det_id_1.shape[0]}")

# ----------------------------------------------------------------------------
Expand Down
Loading

0 comments on commit e21c7dd

Please sign in to comment.