From e30cf184a764a13d968f5c2553e424fd51b64bd4 Mon Sep 17 00:00:00 2001 From: Kris Thielemans Date: Mon, 18 Sep 2023 00:23:33 +0100 Subject: [PATCH] allow using full FOV for Parallelproj projector introduce restrict_to_cylindrical_FOV variable/keyword (defaulting to true), identical to what is used in the ray-tracing matrix This variable is a rename of the _use_truncation variable (which was essentially disabled), but now made to work. --- .../BackProjectorByBinParallelproj.h | 4 ++ .../ForwardProjectorByBinParallelproj.h | 9 ++--- .../ProjectorByBinPairUsingParallelproj.h | 5 +++ .../BackProjectorByBinParallelproj.cxx | 24 ++++++++++-- .../ForwardProjectorByBinParallelproj.cxx | 28 ++++++++++--- .../ProjectorByBinPairUsingParallelproj.cxx | 39 ++++++++++++++++--- 6 files changed, 89 insertions(+), 20 deletions(-) diff --git a/src/include/stir/recon_buildblock/Parallelproj_projector/BackProjectorByBinParallelproj.h b/src/include/stir/recon_buildblock/Parallelproj_projector/BackProjectorByBinParallelproj.h index 9038cab98..73b9148fb 100644 --- a/src/include/stir/recon_buildblock/Parallelproj_projector/BackProjectorByBinParallelproj.h +++ b/src/include/stir/recon_buildblock/Parallelproj_projector/BackProjectorByBinParallelproj.h @@ -85,6 +85,9 @@ class BackProjectorByBinParallelproj : void set_num_gpu_chunks(int num_gpu_chunks) {_num_gpu_chunks = num_gpu_chunks; } int get_num_gpu_chunks() { return _num_gpu_chunks; } + bool get_restrict_to_cylindrical_FOV() const; + void set_restrict_to_cylindrical_FOV(bool val); + protected: virtual void actual_back_project(const RelatedViewgrams&, @@ -100,6 +103,7 @@ class BackProjectorByBinParallelproj : void set_helper(shared_ptr); bool _cuda_verbosity; int _num_gpu_chunks; + bool _restrict_to_cylindrical_FOV; }; END_NAMESPACE_STIR diff --git a/src/include/stir/recon_buildblock/Parallelproj_projector/ForwardProjectorByBinParallelproj.h b/src/include/stir/recon_buildblock/Parallelproj_projector/ForwardProjectorByBinParallelproj.h index fceeb1c53..94756b468 100644 --- a/src/include/stir/recon_buildblock/Parallelproj_projector/ForwardProjectorByBinParallelproj.h +++ b/src/include/stir/recon_buildblock/Parallelproj_projector/ForwardProjectorByBinParallelproj.h @@ -80,14 +80,13 @@ virtual void set_up( /// Set verbosity void set_verbosity(const bool verbosity) { _cuda_verbosity = verbosity; } - /// Set use truncation - truncate before forward - /// projection and after back projection - void set_use_truncation(const bool use_truncation) { _use_truncation = use_truncation; } - // set/get number of gpu chunks to use void set_num_gpu_chunks(int num_gpu_chunks) {_num_gpu_chunks = num_gpu_chunks; } int get_num_gpu_chunks() { return _num_gpu_chunks; } + bool get_restrict_to_cylindrical_FOV() const; + void set_restrict_to_cylindrical_FOV(bool val); + protected: virtual void actual_forward_project(RelatedViewgrams& viewgrams, @@ -102,7 +101,7 @@ virtual void set_up( friend class ProjectorByBinPairUsingParallelproj; void set_helper(shared_ptr); bool _cuda_verbosity; - bool _use_truncation; + bool _restrict_to_cylindrical_FOV; int _num_gpu_chunks; }; diff --git a/src/include/stir/recon_buildblock/Parallelproj_projector/ProjectorByBinPairUsingParallelproj.h b/src/include/stir/recon_buildblock/Parallelproj_projector/ProjectorByBinPairUsingParallelproj.h index b08a2e058..a4d12e480 100644 --- a/src/include/stir/recon_buildblock/Parallelproj_projector/ProjectorByBinPairUsingParallelproj.h +++ b/src/include/stir/recon_buildblock/Parallelproj_projector/ProjectorByBinPairUsingParallelproj.h @@ -57,6 +57,9 @@ class ProjectorByBinPairUsingParallelproj : /// Set verbosity void set_verbosity(const bool verbosity); + bool get_restrict_to_cylindrical_FOV() const; + void set_restrict_to_cylindrical_FOV(bool val); + private: shared_ptr _helper; @@ -64,6 +67,8 @@ class ProjectorByBinPairUsingParallelproj : void initialise_keymap() override; bool post_processing() override; bool _verbosity; + bool _restrict_to_cylindrical_FOV; + bool _already_set_up; }; END_NAMESPACE_STIR diff --git a/src/recon_buildblock/Parallelproj_projector/BackProjectorByBinParallelproj.cxx b/src/recon_buildblock/Parallelproj_projector/BackProjectorByBinParallelproj.cxx index 77291a8d6..c9a2a7948 100644 --- a/src/recon_buildblock/Parallelproj_projector/BackProjectorByBinParallelproj.cxx +++ b/src/recon_buildblock/Parallelproj_projector/BackProjectorByBinParallelproj.cxx @@ -67,6 +67,7 @@ initialise_keymap() parser.add_start_key("Back Projector Using Parallelproj Parameters"); parser.add_stop_key("End Back Projector Using Parallelproj Parameters"); parser.add_key("verbosity", &_cuda_verbosity); + parser.add_key("restrict to cylindrical FOV", &_restrict_to_cylindrical_FOV); parser.add_key("num_gpu_chunks", &_num_gpu_chunks); } @@ -76,8 +77,22 @@ set_defaults() { _cuda_verbosity = true; _num_gpu_chunks = 1; + _restrict_to_cylindrical_FOV = true; } +bool +BackProjectorByBinParallelproj:: +get_restrict_to_cylindrical_FOV() const +{ + return this->_restrict_to_cylindrical_FOV; +} + +void +BackProjectorByBinParallelproj:: +set_restrict_to_cylindrical_FOV(bool val) +{ + this->_restrict_to_cylindrical_FOV = val; +} void BackProjectorByBinParallelproj::set_helper(shared_ptr helper) @@ -91,6 +106,10 @@ BackProjectorByBinParallelproj:: set_up(const shared_ptr& proj_data_info_sptr, const shared_ptr >& density_info_sptr) { +#ifdef STIR_TOF + if (proj_data_info_sptr->get_num_tof_poss() > 1) + error("STIR-ParallelProj interface does not support TOF data yet. Sorry!"); +#endif BackProjectorByBin::set_up(proj_data_info_sptr,density_info_sptr); check(*proj_data_info_sptr, *_density_sptr); _symmetries_sptr.reset(new TrivialDataSymmetriesForBins(proj_data_info_sptr)); @@ -201,10 +220,7 @@ get_output(DiscretisedDensity<3,float> &density) const // --------------------------------------------------------------- // std::copy(image_vec.begin(), image_vec.end(), density.begin_all()); - // After the back projection, we enforce a truncation outside of the FOV. - // This is because the parallelproj projector seems to have some trouble at the edges and this - // could cause some voxel values to spiral out of control. - //if (_use_truncation) + if (this->_restrict_to_cylindrical_FOV) { const float radius = p.get_proj_data_info_sptr()->get_scanner_sptr()->get_inner_ring_radius(); const float image_radius = _helper->voxsize[2]*_helper->imgdim[2]/2; diff --git a/src/recon_buildblock/Parallelproj_projector/ForwardProjectorByBinParallelproj.cxx b/src/recon_buildblock/Parallelproj_projector/ForwardProjectorByBinParallelproj.cxx index 020a84e25..73a0ae36f 100644 --- a/src/recon_buildblock/Parallelproj_projector/ForwardProjectorByBinParallelproj.cxx +++ b/src/recon_buildblock/Parallelproj_projector/ForwardProjectorByBinParallelproj.cxx @@ -42,7 +42,7 @@ ForwardProjectorByBinParallelproj::registered_name = "Parallelproj"; ForwardProjectorByBinParallelproj::ForwardProjectorByBinParallelproj() : - _cuda_verbosity(true), _use_truncation(true), _num_gpu_chunks(1) + _cuda_verbosity(true), _restrict_to_cylindrical_FOV(true), _num_gpu_chunks(1) { this->_already_set_up = false; this->_do_not_setup_helper = false; @@ -59,6 +59,7 @@ initialise_keymap() parser.add_start_key("Forward Projector Using Parallelproj Parameters"); parser.add_stop_key("End Forward Projector Using Parallelproj Parameters"); parser.add_key("verbosity", &_cuda_verbosity); + parser.add_key("restrict to cylindrical FOV", &_restrict_to_cylindrical_FOV); parser.add_key("num_gpu_chunks", &_num_gpu_chunks); } @@ -67,10 +68,24 @@ ForwardProjectorByBinParallelproj:: set_defaults() { _cuda_verbosity = true; - _use_truncation = true; + _restrict_to_cylindrical_FOV = true; _num_gpu_chunks = 1; } +bool +ForwardProjectorByBinParallelproj:: +get_restrict_to_cylindrical_FOV() const +{ + return this->_restrict_to_cylindrical_FOV; +} + +void +ForwardProjectorByBinParallelproj:: +set_restrict_to_cylindrical_FOV(bool val) +{ + this->_restrict_to_cylindrical_FOV = val; +} + void ForwardProjectorByBinParallelproj::set_helper(shared_ptr helper) @@ -84,6 +99,10 @@ ForwardProjectorByBinParallelproj:: set_up(const shared_ptr& proj_data_info_sptr, const shared_ptr >& density_info_sptr) { +#ifdef STIR_TOF + if (proj_data_info_sptr->get_num_tof_poss() > 1) + error("STIR-ParallelProj interface does not support TOF data yet. Sorry!"); +#endif ForwardProjectorByBin::set_up(proj_data_info_sptr,density_info_sptr); check(*proj_data_info_sptr, *_density_sptr); _symmetries_sptr.reset(new TrivialDataSymmetriesForBins(proj_data_info_sptr)); @@ -135,10 +154,7 @@ set_input(const DiscretisedDensity<3,float> & density) { ForwardProjectorByBin::set_input(density); - // Before forward projection, we enforce a truncation outside of the FOV. - // This is because the parallelproj projector seems to have some trouble at the edges and this - // could cause some voxel values to spiral out of control. - //if (_use_truncation) + if (this->_restrict_to_cylindrical_FOV) { const float radius = this->_proj_data_info_sptr->get_scanner_sptr()->get_inner_ring_radius(); const float image_radius = _helper->voxsize[2]*_helper->imgdim[2]/2; diff --git a/src/recon_buildblock/Parallelproj_projector/ProjectorByBinPairUsingParallelproj.cxx b/src/recon_buildblock/Parallelproj_projector/ProjectorByBinPairUsingParallelproj.cxx index 4627ccc84..22e9cf983 100644 --- a/src/recon_buildblock/Parallelproj_projector/ProjectorByBinPairUsingParallelproj.cxx +++ b/src/recon_buildblock/Parallelproj_projector/ProjectorByBinPairUsingParallelproj.cxx @@ -39,6 +39,7 @@ ProjectorByBinPairUsingParallelproj::initialise_keymap() parser.add_start_key("Projector Pair Using Parallelproj Parameters"); parser.add_stop_key("End Projector Pair Using Parallelproj Parameters"); parser.add_key("verbosity",&_verbosity); + parser.add_key("restrict to cylindrical FOV", &_restrict_to_cylindrical_FOV); } @@ -47,6 +48,8 @@ ProjectorByBinPairUsingParallelproj::set_defaults() { base_type::set_defaults(); this->set_verbosity(true); + this->set_restrict_to_cylindrical_FOV(true); + this->_already_set_up = false; } bool @@ -67,22 +70,48 @@ ProjectorByBinPairUsingParallelproj() set_defaults(); } +bool +ProjectorByBinPairUsingParallelproj:: +get_restrict_to_cylindrical_FOV() const +{ + return this->_restrict_to_cylindrical_FOV; +} + +void +ProjectorByBinPairUsingParallelproj:: +set_restrict_to_cylindrical_FOV(bool val) +{ + this->_already_set_up = this->_already_set_up && (this->_restrict_to_cylindrical_FOV == val); + this->_restrict_to_cylindrical_FOV = val; +} + Succeeded ProjectorByBinPairUsingParallelproj:: set_up(const shared_ptr& proj_data_info_sptr, const shared_ptr >& image_info_sptr) { - _helper = std::make_shared(*proj_data_info_sptr, *image_info_sptr); - dynamic_pointer_cast(this->forward_projector_sptr) - ->set_helper(_helper); - dynamic_pointer_cast(this->back_projector_sptr) - ->set_helper(_helper); + auto fwd_prj_downcast_sptr = + dynamic_pointer_cast(this->forward_projector_sptr); + if (!fwd_prj_downcast_sptr) + error("internal error: forward projector should be ParallelProj"); + + auto bck_prj_downcast_sptr = + dynamic_pointer_cast(this->back_projector_sptr); + if (!bck_prj_downcast_sptr) + error("internal error: back projector should be ParallelProj"); + + bck_prj_downcast_sptr->set_restrict_to_cylindrical_FOV(this->_restrict_to_cylindrical_FOV); + fwd_prj_downcast_sptr->set_restrict_to_cylindrical_FOV(this->_restrict_to_cylindrical_FOV); + this->_helper = std::make_shared(*proj_data_info_sptr, *image_info_sptr); + fwd_prj_downcast_sptr->set_helper(this->_helper); + bck_prj_downcast_sptr->set_helper(this->_helper); // the forward_projector->set_up etc will be called in the base class if (base_type::set_up(proj_data_info_sptr, image_info_sptr) != Succeeded::yes) return Succeeded::no; + this->_already_set_up = true; return Succeeded::yes; }