update README

README.md

@@ -14,6 +14,21 @@ cd python
 python 01_analytic_petsird_lm_simulator.py
 ```
 
+The simulation script creates a binary petsird LM file, but also many other
+files (e.g. a reference sensitivity image) that are all stored in the output
+directory
+```
+tree --charset=ascii my_lm_sim
+
+my_lm_sim
+|-- reference_histogram_mlem_50_epochs.npy
+|-- reference_sensitivity_image.npy
+|-- scanner_geometry.png
+|-- sim_parameters.json
+|-- simulated_lm_file.bin
+`-- tof_profile_and_sensitivity_image.png
+```
+
 The simulation can be customized in many ways (number of counts to simulate,
 uniform or non-uniform efficiencies ...) via command line options.
 
@@ -22,14 +37,18 @@ These option can be listed via
 python 01_analytic_petsird_lm_simulator.py
 ```
 
+**Note:** The "reference" MLEM using histogrammed data is only run if a 
+value > 0 is given via `--num_epochs_mlem`. Otherwise it is skipped to save
+time.
+
 ## Run a listmode OSEM recon on the simulated
 
 ```
-python 02_recon_block_scanner_listmode_data.py
+python 02_lm_osem_recon_simulated_data.py
 ```
 
 Thes command line optione for the LM OSEM recon script can be listed via
 ```
-python 02_recon_block_scanner_listmode_data.py -h
+python 02_lm_osem_recon_simulated_data.py -h
 ```
 

python/01_analytic_petsird_lm_simulator.py

@@ -305,6 +305,9 @@ def parse_float_tuple(arg):
 else:
     emission_data = img_fwd_tof
 
+# save the ground truth image
+xp.save(output_dir / "ground_truth_image.npy", img)
+
 # %%
 if num_epochs_mlem > 0:
     recon = xp.ones(img_shape, dtype=xp.float32, device=dev)

python/02_recon_block_scanner_listmode_data.py → python/02_lm_osem_recon_simulated_data.py

@@ -18,8 +18,6 @@
 from pathlib import Path
 import argparse
 
-# %%
-
 
 def transform_to_mat44(
     transform: petsird.RigidTransformation,
@@ -87,34 +85,40 @@ def parse_float_tuple(arg):
 # %%
 
 parser = argparse.ArgumentParser()
-parser.add_argument("--fname", type=str, default="my_lm_sim/simulated_lm_file.bin")
+parser.add_argument("--lm_fname", type=str, default="my_lm_sim/simulated_lm_file.bin")
 parser.add_argument("--num_epochs", type=int, default=5)
 parser.add_argument("--num_subsets", type=int, default=20)
 parser.add_argument("--img_shape", type=parse_int_tuple, default=(100, 100, 11))
 parser.add_argument("--voxel_size", type=parse_float_tuple, default=(1.0, 1.0, 1.0))
 parser.add_argument("--fwhm_mm", type=float, default=1.5)
+parser.add_argument("--output_dir", type=str, default="my_lm_sim")
 
 args = parser.parse_args()
 
-fname = args.fname
+lm_fname = args.lm_fname
 num_epochs = args.num_epochs
 num_subsets = args.num_subsets
 img_shape = args.img_shape
 voxel_size = args.voxel_size
 fwhm_mm = args.fwhm_mm
+output_dir = Path(args.output_dir)
 
 dev = "cpu"
+
+if not output_dir.exists():
+    output_dir.mkdir(parents=True)
+
 # %%
-if not Path(fname).exists():
+if not Path(lm_fname).exists():
     raise FileNotFoundError(
-        f"{args.fname} not found. Create it first using the generator."
+        f"{args.lm_fname} not found. Create it first using the generator."
     )
 
 # %%
 # read the scanner geometry
 
 
-reader = petsird.BinaryPETSIRDReader(fname)
+reader = petsird.BinaryPETSIRDReader(lm_fname)
 header = reader.read_header()
 
 # %%
@@ -317,3 +321,7 @@ def parse_float_tuple(arg):
         recon *= tmp / sens_img
 
 print("")
+
+opath = output_dir / f"lm_osem_{num_epochs}_{num_subsets}.npy"
+xp.save(opath, recon)
+print(f"LM OSEM recon saved to {opath}")