diff --git a/documentation/release_5.2.htm b/documentation/release_5.2.htm index b16aa7624..3a6abc3e2 100644 --- a/documentation/release_5.2.htm +++ b/documentation/release_5.2.htm @@ -39,6 +39,9 @@

New functionality

  • Global variables in SPECTUB have been substituted by class members, such that multiple SPECTUB projectors can be used.
    PR #1169.
  • +
  • Global variables in PinholeSPECTUB have been substituted by class members, such that multiple PinholeSPECTUB projectors can be used. +
    PR #1212. +
  • Scatter estimation is now smoothed in axial direction for BlocksOnCylindrical scanners.
    PR #1172.
  • diff --git a/src/include/stir/recon_buildblock/PinholeSPECTUB_Tools.h b/src/include/stir/recon_buildblock/PinholeSPECTUB_Tools.h index c90b48fda..c135a3c4f 100644 --- a/src/include/stir/recon_buildblock/PinholeSPECTUB_Tools.h +++ b/src/include/stir/recon_buildblock/PinholeSPECTUB_Tools.h @@ -384,43 +384,43 @@ namespace SPECTUB_mph //... functions from wmtools_SPECT.cpp ......................................... - void wm_alloc ( int * Nitems); // allocate wm + void wm_alloc ( int * Nitems, wm_da_type &wm, wmh_mph_type &wmh ); // allocate wm - void free_wm () ; // delete wm + // void free_wm () ; // delete wm - void write_wm_FC_mph (); // write double array weight matrix + // void write_wm_FC_mph (); // write double array weight matrix - void write_wm_hdr_mph (); // write header of a matrix + // void write_wm_hdr_mph (); // write header of a matrix - void write_wm_STIR_mph (); // write matrix in STIR format + // void write_wm_STIR_mph (); // write matrix in STIR format - void read_prj_params_mph (); // read ring parameters from a file + void read_prj_params_mph ( wmh_mph_type &wmh ); // read ring parameters from a file - void read_coll_params_mph (); // read collimator parameters from a file + void read_coll_params_mph ( wmh_mph_type &wmh ); // read collimator parameters from a file - void which_hole(); + // void which_hole(); - void fill_pcf (); // fill precalculated functions + void fill_pcf ( wmh_mph_type &wmh, pcf_type &pcf ); // fill precalculated functions - void free_pcf (); // fill precalculated functions + // void free_pcf (); // fill precalculated functions void calc_cumsum ( discrf2d_type *f ); - void generate_msk_mph ( bool *msk_3d, float *att ); // create a boolean mask for wm (no weights outside the msk) + void generate_msk_mph ( bool *msk_3d, float *att, wmh_mph_type &wmh ); // create a boolean mask for wm (no weights outside the msk) - void read_msk_file_mph ( bool * msk ); // read mask from a file + // void read_msk_file_mph ( bool * msk ); // read mask from a file std::string wm_SPECT_read_value_1d ( std::ifstream * stream1, char DELIMITER ); - void wm_SPECT_read_hvalues_mph ( std::ifstream * stream1, char DELIMITER, int * nh, bool do_cyl ); + void wm_SPECT_read_hvalues_mph ( std::ifstream * stream1, char DELIMITER, int * nh, bool do_cyl, wmh_mph_type &wmh ); - void read_att_map_mph ( float *attmap ); // read attenuation map from a file + // void read_att_map_mph ( float *attmap ); // read attenuation map from a file diff --git a/src/include/stir/recon_buildblock/PinholeSPECTUB_Weight3d.h b/src/include/stir/recon_buildblock/PinholeSPECTUB_Weight3d.h index d8ec97b2e..2e1b08e5a 100644 --- a/src/include/stir/recon_buildblock/PinholeSPECTUB_Weight3d.h +++ b/src/include/stir/recon_buildblock/PinholeSPECTUB_Weight3d.h @@ -18,26 +18,19 @@ namespace SPECTUB_mph { -/* Global variables: the matrix etc TODO - Need to be defined elsewhere - - A lot of the functions below modify these variables. -*/ - //... global variables .............................................. - - extern wmh_mph_type wmh; - extern wm_da_type wm; - extern pcf_type pcf; void wm_calculation_mph ( bool do_estim, const int kOS, psf2d_type *psf2d_bin, - psf2d_type *psf_subs , - psf2d_type *psf2d_aux , + psf2d_type *psf_subs, + psf2d_type *psf2d_aux, psf2d_type *kern, float *attmap, bool *msk_3d, - int *Nitems ); + int *Nitems, + wmh_mph_type &wmh, + wm_da_type &wm, + pcf_type &pcf ); //... geometric component ............................................ @@ -49,21 +42,21 @@ namespace SPECTUB_mph //bool check_zang_obl( lor_type * l, voxel_type * vox, hole_type * h); - void voxel_projection_mph ( lor_type * l, voxel_type * v, hole_type * h ); + void voxel_projection_mph ( lor_type * l, voxel_type * v, hole_type * h, wmh_mph_type &wmh ); - void fill_psfi ( psf2d_type * kern ); + void fill_psfi ( psf2d_type * kern, wmh_mph_type &wmh ); void downsample_psf ( psf2d_type * psf_in, psf2d_type * psf_out, int factor, bool do_calc ); void psf_convol( psf2d_type * psf1, psf2d_type * psf_aux, psf2d_type * psf2, bool do_calc ); - float bresenh_f( int i1, int j1, int i2, int j2, float ** f, int imax, int jmax, float dcr ); + float bresenh_f( int i1, int j1, int i2, int j2, float ** f, int imax, int jmax, float dcr, wmh_mph_type &wmh, pcf_type &pcf ); - void fill_psf_geo ( psf2d_type *psf2d, lor_type *l, discrf2d_type *f, int factor, bool do_calc ); + void fill_psf_geo ( psf2d_type *psf2d, lor_type *l, discrf2d_type *f, int factor, bool do_calc, wmh_mph_type &wmh ); - void fill_psf_depth ( psf2d_type *psf2d, lor_type *l, discrf2d_type *f, int factor, bool do_calc ); + void fill_psf_depth ( psf2d_type *psf2d, lor_type *l, discrf2d_type *f, int factor, bool do_calc, wmh_mph_type &wmh, pcf_type &pcf ); void psf_convol ( psf2d_type * psf2d, psf2d_type * psf_aux, psf2d_type * kern ); @@ -71,7 +64,7 @@ namespace SPECTUB_mph //... attenuation................................................... - float calc_att_mph ( bin_type bin, voxel_type vox, float *attmap ); + float calc_att_mph ( bin_type bin, voxel_type vox, float *attmap, wmh_mph_type &wmh ); int comp_dist ( float dx, float dy, float dz, float dlast ); diff --git a/src/include/stir/recon_buildblock/ProjMatrixByBinPinholeSPECTUB.h b/src/include/stir/recon_buildblock/ProjMatrixByBinPinholeSPECTUB.h index cd2ccff3f..44d3ad2c2 100644 --- a/src/include/stir/recon_buildblock/ProjMatrixByBinPinholeSPECTUB.h +++ b/src/include/stir/recon_buildblock/ProjMatrixByBinPinholeSPECTUB.h @@ -137,6 +137,9 @@ class ProjMatrixByBinPinholeSPECTUB : //! Default constructor (calls set_defaults()) ProjMatrixByBinPinholeSPECTUB(); + // disable copy-constructor as currently unsafe to copy due to bare pointers + ProjMatrixByBinPinholeSPECTUB(const ProjMatrixByBinPinholeSPECTUB&) = delete; + //! Destructor (deallocates UB SPECT memory) ~ProjMatrixByBinPinholeSPECTUB(); @@ -244,6 +247,10 @@ class ProjMatrixByBinPinholeSPECTUB : bool already_setup; + mutable SPECTUB_mph::wmh_mph_type wmh; // weight matrix header. + mutable SPECTUB_mph::wm_da_type wm; // double array weight matrix structure. + mutable SPECTUB_mph::pcf_type pcf; // pre-calculated functions + virtual void calculate_proj_matrix_elems_for_one_bin(ProjMatrixElemsForOneBin&) const; diff --git a/src/recon_buildblock/PinholeSPECTUB_Tools.cxx b/src/recon_buildblock/PinholeSPECTUB_Tools.cxx index 18b695081..866fac2ad 100644 --- a/src/recon_buildblock/PinholeSPECTUB_Tools.cxx +++ b/src/recon_buildblock/PinholeSPECTUB_Tools.cxx @@ -27,8 +27,6 @@ #include #include "stir/info.h" #include "stir/error.h" -#include "stir/warning.h" -#include "stir/error.h" #include #include #include @@ -63,17 +61,12 @@ float dg2rd = boost::math::constants::pi() / (float)180. ; #define DELIMITER1 '#' //delimiter character in input parameter text file #define DELIMITER2 '%' //delimiter character in input parameter text file -//... global variables .............................................. - -extern wmh_mph_type wmh; -extern wm_da_type wm; -extern pcf_type pcf; //============================================================================= //=== wm_alloc ============================================================= //============================================================================= -void wm_alloc( int * Nitems) +void wm_alloc( int * Nitems, wm_da_type &wm, wmh_mph_type &wmh ) { //... double array wm.val and wm.col ..................................................... @@ -279,7 +272,7 @@ void write_wm_STIR_mph() //=== precalculated functions =============================================== //============================================================================== -void fill_pcf() +void fill_pcf( wmh_mph_type &wmh, pcf_type &pcf ) { //... distribution function for a round shape hole ................. @@ -576,7 +569,7 @@ void read_prj_params_mph() }*/ -void read_prj_params_mph() +void read_prj_params_mph( wmh_mph_type &wmh ) { string token; detel_type d; @@ -729,7 +722,7 @@ void read_prj_params_mph() //=== read collimator params mph =============================================== //============================================================================== -void read_coll_params_mph( ) +void read_coll_params_mph( wmh_mph_type &wmh ) { string token; vector param; @@ -757,9 +750,9 @@ void read_coll_params_mph( ) // wmh.collim.holes = new hole_type [ wmh.collim.Nht ]; int nh = 0; - if ( wmh.collim.model == "cyl" ) wm_SPECT_read_hvalues_mph( &stream1, DELIMITER, &nh, true ); + if ( wmh.collim.model == "cyl" ) wm_SPECT_read_hvalues_mph( &stream1, DELIMITER, &nh, true, wmh ); else{ - if ( wmh.collim.model == "pol" ) wm_SPECT_read_hvalues_mph( &stream1, DELIMITER, &nh, false ); + if ( wmh.collim.model == "pol" ) wm_SPECT_read_hvalues_mph( &stream1, DELIMITER, &nh, false, wmh ); else error_wmtools_SPECT_mph ( 334, 0, wmh.collim.model ); } @@ -790,7 +783,7 @@ void read_coll_params_mph( ) //======== wm_SPECT_read_hvalues_mph ============================== //===================================================================== -void wm_SPECT_read_hvalues_mph( ifstream * stream1 , char DELIMITER, int * nh, bool do_cyl ) +void wm_SPECT_read_hvalues_mph( ifstream * stream1 , char DELIMITER, int * nh, bool do_cyl, wmh_mph_type &wmh ) { size_t pos1, pos2, pos3; @@ -994,7 +987,7 @@ void wm_SPECT_read_hvalues_mph( ifstream * stream1 , char DELIMITER, int * nh, b //=== generate_msk_mph ======================================================== //============================================================================= -void generate_msk_mph ( bool *msk_3d, float *attmap ) +void generate_msk_mph ( bool *msk_3d, float *attmap, wmh_mph_type &wmh ) { // bool do_save_resulting_msk = true; diff --git a/src/recon_buildblock/PinholeSPECTUB_Weight3d.cxx b/src/recon_buildblock/PinholeSPECTUB_Weight3d.cxx index 7004c6f87..f649fa6d5 100644 --- a/src/recon_buildblock/PinholeSPECTUB_Weight3d.cxx +++ b/src/recon_buildblock/PinholeSPECTUB_Weight3d.cxx @@ -66,7 +66,10 @@ void wm_calculation_mph ( bool do_calc, psf2d_type *kern, float *attmap, bool *msk_3d, - int *Nitems ) + int *Nitems, + wmh_mph_type &wmh, + wm_da_type &wm, + pcf_type &pcf ) { voxel_type vox; // structure with voxel information bin_type bin; // structure with bin information @@ -162,7 +165,7 @@ void wm_calculation_mph ( bool do_calc, //...vector voxel-hole, angles and distances............................... - voxel_projection_mph ( &l , &vox , h ); + voxel_projection_mph ( &l , &vox , h, wmh ); //... hole shape ....................................... @@ -173,16 +176,16 @@ void wm_calculation_mph ( bool do_calc, if ( wmh.do_subsamp ){ - if ( wmh.do_depth ) fill_psf_depth ( psf_subs, &l, f, wmh.subsamp, do_calc ); + if ( wmh.do_depth ) fill_psf_depth ( psf_subs, &l, f, wmh.subsamp, do_calc, wmh, pcf ); - else fill_psf_geo ( psf_subs, &l, f, wmh.subsamp, do_calc ); + else fill_psf_geo ( psf_subs, &l, f, wmh.subsamp, do_calc, wmh ); if ( wmh.do_psfi ) psf_convol ( psf_subs, psf_aux, kern, do_calc ); downsample_psf( psf_subs, psf_bin, wmh.subsamp, do_calc ); } - else fill_psf_geo ( psf_bin, &l, f, 1, do_calc ); + else fill_psf_geo ( psf_bin, &l, f, 1, do_calc, wmh ); //... calculus of simple attenuation ............................. @@ -194,7 +197,7 @@ void wm_calculation_mph ( bool do_calc, bin.y = d->y0 + l.x1d_l * d->sinth; bin.z = d->z0 + l.z1d_l ; - coeff_att = calc_att_mph( bin, vox, attmap ); + coeff_att = calc_att_mph( bin, vox, attmap, wmh ); } } @@ -236,7 +239,7 @@ void wm_calculation_mph ( bool do_calc, bin.y = d->ybin0 + (float) ib * d->incy ; bin.z = d->zbin0 + (float) jb * wmh.prj.thcm ; - coeff_att = calc_att_mph( bin, vox, attmap ); + coeff_att = calc_att_mph( bin, vox, attmap, wmh ); } //... calculus and storage of the weight............ @@ -269,7 +272,7 @@ void wm_calculation_mph ( bool do_calc, //=== fill_psfi ============================================================ //========================================================================== -void fill_psfi( psf2d_type * kern ) +void fill_psfi( psf2d_type * kern, wmh_mph_type &wmh ) { //float K0 = (float)0.39894228040143 / wmh.prj.sgm_i ; //Normalization factor: 1/sqrt(2*M_PI)/sigma float K0 = (1.0f/boost::math::constants::root_two_pi()) / wmh.prj.sgm_i ; //Normalization factor: 1/sqrt(2*M_PI)/sigma @@ -385,7 +388,7 @@ bool check_zang_par( voxel_type * v, hole_type * h ){ //=== voxel_projection ===================================================== //========================================================================== -void voxel_projection_mph ( lor_type * l, voxel_type * v, hole_type * h ) +void voxel_projection_mph ( lor_type * l, voxel_type * v, hole_type * h, wmh_mph_type &wmh ) { //...vector voxel-hole, angles and distances............................... @@ -444,7 +447,7 @@ void voxel_projection_mph ( lor_type * l, voxel_type * v, hole_type * h ) //=== fill_psf_geo ========================================================= //========================================================================== -void fill_psf_geo ( psf2d_type * psf, lor_type *l, discrf2d_type *f, int factor, bool do_calc ) +void fill_psf_geo ( psf2d_type * psf, lor_type *l, discrf2d_type *f, int factor, bool do_calc, wmh_mph_type &wmh ) { psf->xc = l->x1d_l + wmh.prj.FOVxcmd2; // x distance of center of PSF to the begin of the FOVcm psf->zc = l->z1d_l + wmh.prj.FOVzcmd2; // z distance of center of PSF to the begin of the FOVcm @@ -516,7 +519,7 @@ void fill_psf_geo ( psf2d_type * psf, lor_type *l, discrf2d_type *f, int factor, //=== fill_psf_depth ========================================================== //============================================================================= -void fill_psf_depth( psf2d_type *psf, lor_type *l, discrf2d_type *f, int factor, bool do_calc ) +void fill_psf_depth( psf2d_type *psf, lor_type *l, discrf2d_type *f, int factor, bool do_calc, wmh_mph_type &wmh, pcf_type &pcf ) { float xc_d = l->x1d_l + wmh.prj.FOVxcmd2; // x distance of center of PSF from the begin of the FOVcm @@ -609,7 +612,7 @@ void fill_psf_depth( psf2d_type *psf, lor_type *l, discrf2d_type *f, int factor, if_d = if0_d + i * incxf_d ; if_dc = if0_dc + i * incxf_dc ; - v = bresenh_f( if_d, jf_d, if_dc, jf_dc, f->val, f->i_max, f->j_max, dcr ); + v = bresenh_f( if_d, jf_d, if_dc, jf_dc, f->val, f->i_max, f->j_max, dcr, wmh, pcf ); psf->val [ j ][ i ] += v ; psf->val [ j - 1 ][ i - 1 ] += v ; @@ -628,7 +631,7 @@ void fill_psf_depth( psf2d_type *psf, lor_type *l, discrf2d_type *f, int factor, if_d = if0_d ; if_dc = if0_dc ; - v = bresenh_f( if_d, jf_d, if_dc, jf_dc, f->val, f->i_max, f->j_max, dcr ); + v = bresenh_f( if_d, jf_d, if_dc, jf_dc, f->val, f->i_max, f->j_max, dcr, wmh, pcf ); psf->val [ j ][ 0 ] += v ; psf->val [ j - 1 ][ 0 ] -= v ; @@ -636,7 +639,7 @@ void fill_psf_depth( psf2d_type *psf, lor_type *l, discrf2d_type *f, int factor, if_d = if0_d + psf->dimx * incxf_d ; if_dc = if0_dc + psf->dimx * incxf_d ; - v = bresenh_f( if_d, jf_d, if_dc, jf_dc, f->val, f->i_max, f->j_max, dcr ); + v = bresenh_f( if_d, jf_d, if_dc, jf_dc, f->val, f->i_max, f->j_max, dcr, wmh, pcf ); psf->val [ j ][ psf->dimx - 1 ] -= v ; psf->val [ j - 1 ][ psf->dimx - 1 ] += v ; @@ -652,7 +655,7 @@ void fill_psf_depth( psf2d_type *psf, lor_type *l, discrf2d_type *f, int factor, if_d = if0_d + i * incxf_d ; if_dc = if0_dc + i * incxf_dc ; - v = bresenh_f( if_d, jf_d, if_dc, jf_dc, f->val, f->i_max, f->j_max, dcr ); + v = bresenh_f( if_d, jf_d, if_dc, jf_dc, f->val, f->i_max, f->j_max, dcr, wmh, pcf ); psf->val [ 0 ][ i ] += v ; psf->val [ 0 ][ i - 1 ] -= v ; @@ -660,7 +663,7 @@ void fill_psf_depth( psf2d_type *psf, lor_type *l, discrf2d_type *f, int factor, jf_d = jf0_d + psf->dimz * inczf_d ; jf_dc = jf0_dc + psf->dimz * inczf_d ; - v = bresenh_f( if_d, jf_d, if_dc, jf_dc, f->val, f->i_max, f->j_max, dcr ); + v = bresenh_f( if_d, jf_d, if_dc, jf_dc, f->val, f->i_max, f->j_max, dcr, wmh, pcf ); psf->val [ psf->dimz - 1 ][ i ] -= v ; psf->val [ psf->dimz - 1 ][ i - 1 ] += v ; @@ -674,7 +677,7 @@ void fill_psf_depth( psf2d_type *psf, lor_type *l, discrf2d_type *f, int factor, jf_d = jf0_d ; jf_dc = jf0_dc ; - v = bresenh_f( if_d, jf_d, if_dc, jf_dc, f->val, f->i_max, f->j_max, dcr ); + v = bresenh_f( if_d, jf_d, if_dc, jf_dc, f->val, f->i_max, f->j_max, dcr, wmh, pcf ); psf->val [ 0 ][ 0 ] += v ; psf->sum += v ; @@ -684,7 +687,7 @@ void fill_psf_depth( psf2d_type *psf, lor_type *l, discrf2d_type *f, int factor, if_d = if0_d + psf->dimx * incxf_d ; if_dc = if0_dc + psf->dimx * incxf_d ; - v = bresenh_f( if_d, jf_d, if_dc, jf_dc, f->val, f->i_max, f->j_max, dcr ); + v = bresenh_f( if_d, jf_d, if_dc, jf_dc, f->val, f->i_max, f->j_max, dcr, wmh, pcf ); psf->val [ 0 ][ psf->dimx - 1 ] -= v ; psf->sum -= v ; @@ -694,7 +697,7 @@ void fill_psf_depth( psf2d_type *psf, lor_type *l, discrf2d_type *f, int factor, jf_d = jf0_d + psf->dimz * inczf_d ; jf_dc = jf0_dc + psf->dimz * inczf_d ; - v = bresenh_f( if_d, jf_d, if_dc, jf_dc, f->val, f->i_max, f->j_max, dcr ); + v = bresenh_f( if_d, jf_d, if_dc, jf_dc, f->val, f->i_max, f->j_max, dcr, wmh, pcf ); psf->val [ psf->dimz - 1 ][ psf->dimx - 1 ] += v ; psf->sum += v ; @@ -704,7 +707,7 @@ void fill_psf_depth( psf2d_type *psf, lor_type *l, discrf2d_type *f, int factor, if_d = if0_d ; if_dc = if0_dc ; - v = bresenh_f( if_d, jf_d, if_dc, jf_dc, f->val, f->i_max, f->j_max, dcr ); + v = bresenh_f( if_d, jf_d, if_dc, jf_dc, f->val, f->i_max, f->j_max, dcr, wmh, pcf ); psf->val [ psf->dimz - 1 ][ 0 ] -= v ; psf->sum -= v ; @@ -837,7 +840,7 @@ void psf_convol( psf2d_type * psf, psf2d_type * psf_aux, psf2d_type * kern, bool //=== bresenh_f ============================================================ //========================================================================== -float bresenh_f( int i1, int j1, int i2, int j2, float ** f , int imax, int jmax, float dcr ) +float bresenh_f( int i1, int j1, int i2, int j2, float ** f , int imax, int jmax, float dcr, wmh_mph_type &wmh, pcf_type &pcf ) { int er; //the error term @@ -932,7 +935,7 @@ float bresenh_f( int i1, int j1, int i2, int j2, float ** f , int imax, int jmax //=== calc_att_mph ============================================================= //============================================================================= -float calc_att_mph( bin_type bin, voxel_type vox, float * attmap ) +float calc_att_mph( bin_type bin, voxel_type vox, float * attmap, wmh_mph_type &wmh ) { float dx, dy, dz; float dlast_x, dlast_y, dlast_z, dlast; diff --git a/src/recon_buildblock/ProjMatrixByBinPinholeSPECTUB.cxx b/src/recon_buildblock/ProjMatrixByBinPinholeSPECTUB.cxx index 2ff20f51f..3c258b16b 100644 --- a/src/recon_buildblock/ProjMatrixByBinPinholeSPECTUB.cxx +++ b/src/recon_buildblock/ProjMatrixByBinPinholeSPECTUB.cxx @@ -56,8 +56,6 @@ #include "stir/is_null_ptr.h" #include "stir/Coordinate3D.h" #include "stir/info.h" -#include "stir/warning.h" -#include "stir/error.h" #include "stir/CPUTimer.h" #ifdef STIR_OPENMP #include "stir/num_threads.h" @@ -83,22 +81,12 @@ //void wm_inputs_mph( char ** argv, int argc ); //void read_inputs_mph( vector param ); -//.. global variables .................................................................... - -namespace SPECTUB_mph -{ - wmh_mph_type wmh; // weight matrix header. Global variable - wm_da_type wm; // double array weight matrix structure. Global variable - pcf_type pcf; // pre-calculated functions -} - using namespace std; using namespace SPECTUB_mph; using std::string; START_NAMESPACE_STIR - const char * const ProjMatrixByBinPinholeSPECTUB::registered_name = "Pinhole SPECT UB"; @@ -503,8 +491,6 @@ set_up( } #endif - using namespace SPECTUB_mph; - std::stringstream info_stream; const VoxelsOnCartesianGrid * image_info_ptr = @@ -675,12 +661,12 @@ set_up( //... files with complementary information ................. - read_prj_params_mph(); - read_coll_params_mph(); + read_prj_params_mph( wmh ); + read_coll_params_mph( wmh ); //... precalculated functions ................ - fill_pcf(); + fill_pcf( wmh, pcf ); //... other variables ......................... @@ -736,7 +722,7 @@ set_up( //... to generate mask.......................................................... msk_3d = new bool [ wmh.vol.Nvox ]; - //generate_msk_mph( msk_3d, attmap ); + //generate_msk_mph( msk_3d, attmap, wmh ); // generate mask from mask file if (!is_null_ptr(mask_image_sptr)) @@ -754,7 +740,7 @@ set_up( // we do this to avoid using its own read_msk_file_mph wmh.do_msk_file = false; wmh.do_msk_att = true; - generate_msk_mph( msk_3d, mask_from_file ); + generate_msk_mph( msk_3d, mask_from_file, wmh ); delete[] mask_from_file; } @@ -763,13 +749,13 @@ set_up( { if (is_null_ptr(attmap)) error("No attenuation image set, so cannot compute the mask image from it."); - generate_msk_mph( msk_3d, attmap ); + generate_msk_mph( msk_3d, attmap, wmh ); } // generate mask from object radius else { wmh.do_msk_file = false; - generate_msk_mph( msk_3d, (float *)0 ); + generate_msk_mph( msk_3d, (float *)0, wmh ); } //... initialize psf2d in bins .................................................. @@ -807,7 +793,7 @@ set_up( for ( int i = 0 ; i < kern.max_dimz ; i++ ) kern.val[ i ] = new float [ kern.max_dimx ]; - fill_psfi ( &kern ); + fill_psfi ( &kern, wmh ); psf_subs.max_dimx += kern.max_dimx - 1 ; psf_subs.max_dimz += kern.max_dimz - 1 ; @@ -867,7 +853,7 @@ set_up( // size estimation for (int kOS = 0 ; kOS < wmh.prj.NOS ; kOS++) { - wm_calculation_mph ( false , kOS, &psf_bin, &psf_subs, &psf_aux, &kern, attmap, msk_3d, Nitems[ kOS ] ); + wm_calculation_mph ( false , kOS, &psf_bin, &psf_subs, &psf_aux, &kern, attmap, msk_3d, Nitems[ kOS ], wmh, wm, pcf ); } info(boost::format("Done estimating size of matrix. Execution time, CPU %1% s") % timer.value(), 2); @@ -889,8 +875,6 @@ delete_PinholeSPECTUB_arrays() if (!this->already_setup) return; - using namespace SPECTUB_mph; - //... freeing matrix memory.................................... delete [] wm.val; @@ -953,8 +937,6 @@ void ProjMatrixByBinPinholeSPECTUB:: compute_one_subset(const int kOS) const { - using namespace SPECTUB_mph; - CPUTimer timer; timer.start(); @@ -971,11 +953,11 @@ compute_one_subset(const int kOS) const //... memory allocation for wm float arrays .................................................... - wm_alloc( Nitems[kOS] ); + wm_alloc( Nitems[kOS], wm, wmh ); //... wm calculation ............................................................................... - wm_calculation_mph ( true, kOS, &psf_bin, &psf_subs, &psf_aux, &kern, attmap, msk_3d, Nitems[kOS] ); + wm_calculation_mph ( true, kOS, &psf_bin, &psf_subs, &psf_aux, &kern, attmap, msk_3d, Nitems[kOS], wmh, wm, pcf ); info(boost::format("Weight matrix calculation done, CPU %1% s") % timer.value(), 2);