Skip to content

Commit

Permalink
simplified compute_ac_factors methods (C++/C/Python)
Browse files Browse the repository at this point in the history
  • Loading branch information
evgueni-ovtchinnikov committed May 31, 2024
1 parent 1b30b55 commit d2585d9
Show file tree
Hide file tree
Showing 5 changed files with 18 additions and 25 deletions.
12 changes: 4 additions & 8 deletions src/xSTIR/cSTIR/cstir.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -650,21 +650,17 @@ void* cSTIR_createPETAttenuationModel(const void* ptr_img, const void* ptr_am)
}

extern "C"
void* cSTIR_computeACF(const void* ptr_sino, const void* ptr_mumap, const void* ptr_am,
void* cSTIR_computeACF(const void* ptr_sino,
const void* ptr_att, void* ptr_af, void* ptr_acf)
{
try {
SPTR_FROM_HANDLE(STIRAcquisitionData, sptr_sino, ptr_sino);
SPTR_FROM_HANDLE(STIRImageData, sptr_mumap, ptr_mumap);
STIRAcquisitionData& sino = objectFromHandle<STIRAcquisitionData>(ptr_sino);
PETAttenuationModel& att = objectFromHandle<PETAttenuationModel>(ptr_att);
SPTR_FROM_HANDLE(STIRAcquisitionData, sptr_af, ptr_af);
SPTR_FROM_HANDLE(STIRAcquisitionData, sptr_acf, ptr_acf);
SPTR_FROM_HANDLE(AcqMod3DF, sptr_am, ptr_am);
SPTR_FROM_HANDLE(PETAttenuationModel, sptr_att, ptr_att);
PETAttenuationModel::compute_ac_factors(sptr_sino, sptr_mumap,
sptr_am, sptr_att, sptr_af, sptr_acf);
PETAttenuationModel::compute_ac_factors(sino, att, sptr_af, sptr_acf);
HANDLE_FROM_SPTR(STIRAcquisitionData, sptr_af, ptr_af);
HANDLE_FROM_SPTR(STIRAcquisitionData, sptr_acf, ptr_acf);
//HANDLE_FROM_SPTR(PETAttenuationModel, sptr_att, ptr_att);
return (void*) new DataHandle;
}
CATCH;
Expand Down
4 changes: 2 additions & 2 deletions src/xSTIR/cSTIR/include/sirf/STIR/cstir.h
Original file line number Diff line number Diff line change
Expand Up @@ -79,7 +79,7 @@ extern "C" {
void* cSTIR_createPETAttenuationModel
(const void* ptr_img, const void* ptr_am);
void* cSTIR_computeACF
(const void* ptr_sino, const void* ptr_mumap, const void* ptr_am, const void* ptr_att, void* ptr_acf, void* ptr_iacf);
(const void* ptr_sino, const void* ptr_att, void* ptr_acf, void* ptr_iacf);
void* cSTIR_chainPETAcquisitionSensitivityModels
(const void* ptr_first, const void* ptr_second);
void* cSTIR_setupAcquisitionSensitivityModel(void* ptr_sm, void* ptr_ad);
Expand All @@ -96,7 +96,7 @@ extern "C" {
int subset_num, int num_subsets);
void* cSTIR_acquisitionModelBwdReplace(void* ptr_am, void* ptr_ad,
int subset_num, int num_subsets, void* ptr_im);
void* cSTIR_get_MatrixInfo(void* ptr);
void* cSTIR_get_MatrixInfo(void* ptr);

// Acquisition Model Matrix
void* cSTIR_setupSPECTUBMatrix(const void* h_smx, const void* h_acq, const void* h_img);
Expand Down
15 changes: 3 additions & 12 deletions src/xSTIR/cSTIR/include/sirf/STIR/stir_x.h
Original file line number Diff line number Diff line change
Expand Up @@ -557,25 +557,16 @@ The actual algorithm is described in
*/
static void compute_ac_factors(
// input arguments
std::shared_ptr<STIRAcquisitionData> acq_templ_sptr,
std::shared_ptr<STIRImageData> mu_map_sptr,
std::shared_ptr<PETAcquisitionModel>& am_sptr,
std::shared_ptr<PETAttenuationModel>& asm_sptr,
STIRAcquisitionData& acq_templ,
PETAttenuationModel& acq_sens_mod,
// output arguments
std::shared_ptr<STIRAcquisitionData>& af_sptr,
std::shared_ptr<STIRAcquisitionData>& acf_sptr)
{
//PETAcquisitionModelUsingRayTracingMatrix acq_mod;
PETAcquisitionModel& acq_mod = *am_sptr;
acq_mod.set_up(acq_templ_sptr, mu_map_sptr);
PETAttenuationModel& acq_sens_mod = *asm_sptr;
acq_sens_mod.set_up(acq_templ_sptr->get_exam_info_sptr(),
acq_templ_sptr->get_proj_data_info_sptr()->create_shared_clone());
af_sptr = acq_templ_sptr->clone();
af_sptr = acq_templ.clone();
af_sptr->fill(1.0);
acf_sptr = af_sptr->clone();
acq_sens_mod.unnormalise(*af_sptr);
//acq_sens_mod.normalise(*acf_sptr);
acf_sptr->inv(0, *af_sptr);
}

Expand Down
8 changes: 6 additions & 2 deletions src/xSTIR/cSTIR/tests/test7.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -88,10 +88,14 @@ int main()
std::cout << "===== sinograms norm: " << sinograms_sptr->norm() << '\n';
std::cout << "===== randoms norm: " << randoms_sptr->norm() << '\n';

std::shared_ptr<PETAttenuationModel> att_sptr = std::shared_ptr<PETAttenuationModel>(new PETAttenuationModel(mu_map, am));
am.set_up(sinograms_sptr, mu_map_sptr);
STIRAcquisitionData& sinograms = *sinograms_sptr;
PETAttenuationModel att(mu_map, am);
att.set_up(sinograms.get_exam_info_sptr(),
sinograms.get_proj_data_info_sptr()->create_shared_clone());
std::shared_ptr<STIRAcquisitionData> af_sptr; // attenuation factor
std::shared_ptr<STIRAcquisitionData> acf_sptr; // the inverse of the above
PETAttenuationModel::compute_ac_factors(sinograms_sptr, mu_map_sptr, am_sptr, att_sptr, af_sptr, acf_sptr);
PETAttenuationModel::compute_ac_factors(sinograms, att, af_sptr, acf_sptr);
std::cout << "===== norm of the attenuation factor: " << af_sptr->norm() << '\n';
std::cout << "===== norm of the attenuation correction factor: " << acf_sptr->norm() << '\n';

Expand Down
4 changes: 3 additions & 1 deletion src/xSTIR/pSTIR/STIR.py
Original file line number Diff line number Diff line change
Expand Up @@ -1744,7 +1744,9 @@ def compute_attenuation_factors(sinograms, mu_map):
attn = AcquisitionSensitivityModel(mu_map, am)
af = AcquisitionData(sinograms)
acf = AcquisitionData(sinograms)
try_calling(pystir.cSTIR_computeACF(sinograms.handle, mu_map.handle, am.handle, attn.handle, af.handle, acf.handle))
am.set_up(sinograms, mu_map)
attn.set_up(sinograms)
try_calling(pystir.cSTIR_computeACF(sinograms.handle, attn.handle, af.handle, acf.handle))
return af, acf

def __del__(self):
Expand Down

0 comments on commit d2585d9

Please sign in to comment.