From 00c4956e97c6971515332e70d1554904b338a2a8 Mon Sep 17 00:00:00 2001
From: Tamino Huxohl <thuxohl@techfak.uni-bielefeld.de>
Date: Fri, 18 Nov 2022 14:05:19 +0100
Subject: [PATCH] rename recon module to osem and add postfilter support
---
mu_map/recon/{recon.py => osem.py} | 47 ++++++++++++++++++------------
1 file changed, 28 insertions(+), 19 deletions(-)
rename mu_map/recon/{recon.py => osem.py} (89%)
diff --git a/mu_map/recon/recon.py b/mu_map/recon/osem.py
similarity index 89%
rename from mu_map/recon/recon.py
rename to mu_map/recon/osem.py
index 3614817..d4abd98 100644
--- a/mu_map/recon/recon.py
+++ b/mu_map/recon/osem.py
@@ -16,6 +16,7 @@ from mu_map.file.interfile import (
KEY_SPACING_2,
TEMPLATE_HEADER_IMAGE,
)
+from mu_map.recon.filter import GaussianFilter
TEMPLATE_RECON_PARAMS = """
@@ -65,21 +66,6 @@ OSMAPOSLParameters :=
number of subiterations:= {N_SUBITERATIONS}
Save estimates at subiteration intervals:= {SAVE_INTERVALS}
- ; keywords that specify the filtering that occurs after every subiteration
- ; warning: do not normally use together with a prior
- inter-iteration filter subiteration interval := 4
- inter-iteration filter type := Separable Gaussian
- ;post-filter type := Separable Gaussian
- separable gaussian filter parameters :=
- x-dir filter fwhm (in mm) := 9.6
- y-dir filter fwhm (in mm) := 9.6
- z-dir filter fwhm (in mm) := 9.6
- x-dir maximum kernel size := 128
- y-dir maximum kernel size := 128
- z-dir maximum kernel size := 33
- Normalise filter to 1 := 1
- end separable gaussian filter parameters :=
-
END :=
"""
@@ -112,6 +98,7 @@ def reconstruct(
init: Optional[Tuple[Dict[str, str], np.ndarray]] = None,
n_subsets: Optional[int] = 4,
n_iterations: Optional[int] = 10,
+ postfilter: Optional[GaussianFilter] = None,
**kwargs,
):
# sanitize parameters
@@ -155,6 +142,9 @@ def reconstruct(
write_interfile(filename_init, *init)
params = params.replace("{INIT_FILE}", filename_init)
+ if postfilter:
+ params = postfilter.insert(params)
+
filename_params = os.path.join(dir_tmp.name, "OSEM_SPECT.par")
with open(filename_params, mode="w") as f:
f.write(params.strip())
@@ -162,7 +152,6 @@ def reconstruct(
recon = stir.OSMAPOSLReconstruction3DFloat(filename_params)
recon.reconstruct()
- print(params)
return load_interfile(filename_out)
@@ -209,6 +198,18 @@ if __name__ == "__main__":
default=10,
help="the number of iterations for OSEM reconstruction",
)
+ parser.add_argument(
+ "--postfilter",
+ choices=["gaussian", "none"],
+ default="gaussian",
+ help="apply a postfilter to the reconstruction",
+ )
+ parser.add_argument(
+ "--postfilter_width",
+ type=float,
+ default=1.0,
+ help="the filter witdth for the postfilter is based on the spacing in the projection and can be modified with this factor",
+ )
parser.add_argument(
"-v", "--verbosity", type=int, default=0, help="configure the verbosity of STIR"
)
@@ -230,11 +231,19 @@ if __name__ == "__main__":
else:
projection = load_as_interfile(args.projection)
+ postfilter = None
+ if args.postfilter == "gaussian":
+ postfilter = GaussianFilter(projection, args.postfilter_width)
+
kwargs = vars(args)
del kwargs["mu_map"]
del kwargs["projection"]
- print(kwargs)
- header, image = reconstruct(projection, mu_map=mu_map, **kwargs)
- image = image[:, :, ::-1]
+ del kwargs["postfilter"]
+ header, image = reconstruct(
+ projection, mu_map=mu_map, postfilter=postfilter, **kwargs
+ )
+ image = image[
+ :, :, ::-1
+ ] # STIR creates reconstructions flipped on the x-axes, this is reverted here
write_interfile(args.out, header, image)
--
GitLab