CASToR  1.0
Tomographic Reconstruction (PET/SPECT)
oDynamicDataManager.cc
Go to the documentation of this file.
00001 
00002 /*
00003   Implementation of class oDynamicDataManager
00004 
00005   - separators: X
00006   - doxygen: X
00007   - default initialization: X
00008   - CASTOR_DEBUG: none
00009   - CASTOR_VERBOSE: X
00010 */
00011 
00012 
00021 #include "oDynamicDataManager.hh"
00022 #include "oImageDimensionsAndQuantification.hh"
00023 
00024 // =====================================================================
00025 // ---------------------------------------------------------------------
00026 // ---------------------------------------------------------------------
00027 // =====================================================================
00028 /*
00029   Constructor oDynamicDataManager
00030 */
00031 oDynamicDataManager::oDynamicDataManager()
00032 {
00033   mp_ID = NULL;
00034   m_verbose = -1;
00035   m_nbTimeFrames = 1;
00036   mp_currentFrameIndex = NULL;
00037   m_respGatingFlag = false;
00038   m_rMotionCorrFlag = false;
00039   m_nbRespGates = 1;
00040   m2p_nbEventsPerRespGate = NULL;
00041   m2p_indexLastEventRespGate = NULL;
00042   mp_currentRespGateIndex = NULL;
00043   m_cardGatingFlag = false;
00044   m_cMotionCorrFlag = false;
00045   m_nbCardGates = 1;
00046   m2p_nbEventsPerCardGate = NULL;
00047   m2p_indexLastEventCardGate = NULL;
00048   mp_currentCardGateIndex = NULL;
00049   m_pMotionCorrFlag = false;
00050   m_nbPMotionTriggers = 0;
00051   mp_listPMotionTriggers = NULL;
00052   mp_currentPMotionIndex = NULL;
00053 }
00054 
00055 // =====================================================================
00056 // ---------------------------------------------------------------------
00057 // ---------------------------------------------------------------------
00058 // =====================================================================
00059 /*
00060   Destructor oDynamicDataManager
00061 */
00062 oDynamicDataManager::~oDynamicDataManager()
00063 {
00064   for (int fr=0; fr<m_nbTimeFrames; fr++)
00065   {
00066     if (m2p_nbEventsPerRespGate && m2p_nbEventsPerRespGate[fr]) 
00067       delete m2p_nbEventsPerRespGate[fr];
00068     if (m2p_nbEventsPerCardGate && m2p_nbEventsPerCardGate[fr]) 
00069       delete m2p_nbEventsPerCardGate[fr];
00070     if (m2p_indexLastEventRespGate && m2p_indexLastEventRespGate[fr]) 
00071       delete m2p_indexLastEventRespGate[fr];
00072     if (m2p_indexLastEventCardGate && m2p_indexLastEventCardGate[fr]) 
00073       delete m2p_indexLastEventCardGate[fr];
00074   }
00075 
00076   if (m2p_nbEventsPerRespGate!=NULL) delete m2p_nbEventsPerRespGate;
00077   if (m2p_nbEventsPerCardGate!=NULL) delete m2p_nbEventsPerCardGate;
00078   if (m2p_indexLastEventRespGate!= NULL) delete m2p_indexLastEventRespGate;
00079   if (m2p_indexLastEventCardGate!= NULL) delete m2p_indexLastEventCardGate;
00080   if (mp_currentFrameIndex!=NULL) delete mp_currentFrameIndex;
00081   if (mp_currentRespGateIndex!=NULL) delete mp_currentRespGateIndex;
00082   if (mp_currentCardGateIndex!=NULL) delete mp_currentCardGateIndex;
00083   if (mp_currentPMotionIndex!=NULL) delete mp_currentPMotionIndex;
00084   if (mp_listPMotionTriggers!=NULL) delete mp_listPMotionTriggers;
00085 }
00086 
00087 // =====================================================================
00088 // ---------------------------------------------------------------------
00089 // ---------------------------------------------------------------------
00090 // =====================================================================
00091 /*
00092   \fn InitDynamicData
00093   \param a_nbRespGates
00094   \param a_nbCardGates
00095   \param a_pathTo4DDataFile : path to an ASCII file containing dynamic metadata regarding the acquisition
00096   \param a_rmCorrFlag : indicate whether respiratory motion correction is enabled (1) or disabled (0) 
00097   \param a_cmCorrFlag : indicate whether cardiac motion correction is enabled (1) or disabled (0) 
00098   \param a_pmCorrFlag : indicate whether involuntary patient motion correction is enabled (1) or disabled (0) 
00099   \brief Main function for instanciation and initialization of the member variables and arrays. Call the specific initialization function depending of the type of dataset.
00100   \todo  warning in the documentation that Respiratory motion correction is not supported for cardiac-gated data in the current implementation
00101   \return 0 if success, positive value otherwise.
00102 */
00103 int oDynamicDataManager::InitDynamicData( int a_nbRespGates, int a_nbCardGates, const string& a_pathTo4DDataFile,
00104                                           int a_rmCorrFlag, int a_cmCorrFlag, int a_pmCorrFlag)
00105 {
00106   if (m_verbose>=3) Cout("oDynamicDataManager::InitDynamicData() -> Initialize dynamic data management" << endl);
00107 
00108   // Initialization
00109   m_nbTimeFrames = mp_ID->GetNbTimeFrames();
00110   m_nbRespGates = a_nbRespGates;
00111   m_nbCardGates = a_nbCardGates;
00112   if (m_nbRespGates > 1) m_respGatingFlag = true;
00113   if (m_nbCardGates > 1) m_cardGatingFlag = true;
00114   m_rMotionCorrFlag = a_rmCorrFlag;
00115   m_cMotionCorrFlag = a_cmCorrFlag;
00116   m_pMotionCorrFlag = a_pmCorrFlag;
00117   
00118   // Instanciation & initialization for current indices (used during the loop on events to know which frame/gate the event belongs to)
00119   mp_currentFrameIndex     = new int[mp_ID->GetNbThreadsForProjection()];
00120   mp_currentRespGateIndex  = new int[mp_ID->GetNbThreadsForProjection()];
00121   mp_currentCardGateIndex  = new int[mp_ID->GetNbThreadsForProjection()];
00122   mp_currentPMotionIndex = new int[mp_ID->GetNbThreadsForProjection()];
00123   for (int th=0 ; th<mp_ID->GetNbThreadsForProjection() ; th++)
00124   {
00125     mp_currentFrameIndex    [0] = 0;
00126     mp_currentRespGateIndex [0] = 0;
00127     mp_currentCardGateIndex [0] = 0;
00128     mp_currentPMotionIndex[0] = 0;
00129   }
00130   
00131   // Instanciation & initialization of number of events per frame for gated reconstructions
00132   m2p_nbEventsPerRespGate = new int64_t*[m_nbTimeFrames];
00133   m2p_nbEventsPerCardGate = new int64_t*[m_nbTimeFrames];
00134   m2p_indexLastEventRespGate = new int64_t*[m_nbTimeFrames];
00135   m2p_indexLastEventCardGate = new int64_t*[m_nbTimeFrames];
00136   
00137   for (int fr=0; fr<m_nbTimeFrames; fr++)
00138   {
00139     m2p_nbEventsPerRespGate[fr] = new int64_t[m_nbRespGates];
00140     m2p_nbEventsPerCardGate[fr] = new int64_t[m_nbRespGates*m_nbCardGates];
00141     m2p_indexLastEventRespGate[fr] = new int64_t[m_nbRespGates];
00142     m2p_indexLastEventCardGate[fr] = new int64_t[m_nbRespGates*m_nbCardGates];
00143   
00144     for (int rg=0; rg<m_nbRespGates; rg++)
00145     {
00146       m2p_nbEventsPerRespGate[fr][rg] = -1;
00147       m2p_indexLastEventRespGate[fr][rg] = -1;
00148       
00149       for (int cg=0; cg<m_nbCardGates; cg++)
00150       {
00151         m2p_nbEventsPerCardGate[fr][rg*m_nbCardGates+cg] = -1;
00152         m2p_indexLastEventCardGate[fr][rg*m_nbCardGates+cg] = -1;
00153       }
00154     }
00155   }
00156 
00157   // Throw error if cardiac gating is enabled alongside respiratory motion correction (not available)
00158   // TODO : warning in the documentation that Respiratory motion correction is not supported for cardiac-gated data in the current implementation
00159   if (m_rMotionCorrFlag && m_cardGatingFlag)
00160   {
00161     Cerr("***** oDynamicDataManager::InitDynamicData()-> Respiratory motion correction is enabled for cardiac-gated data. This is not supported in the current implementation  !" << endl);
00162     return 1;
00163   }
00164     
00165   // Initialization for respiratory/cardiac gated reconstruction
00166   // Note : for Analytic Projection, gating flag could be enabled in order to project an image containing several gates, but no description file is required
00167   //        InitDynamicDataGating() is restricted to reconstruction
00168   //        Check on the existence of a_pathTo4DDataFile during reconstruction is performed onnly in Grecon
00169   if ((m_respGatingFlag || m_cardGatingFlag) && !a_pathTo4DDataFile.empty() )
00170   {
00171     if(m_verbose >=3) Cout("oDynamicDataManager::InitDynamicData()-> Initializing data for gated reconstruction... " << endl);
00172       
00173     if (InitDynamicDataGating(a_pathTo4DDataFile))
00174     {
00175       Cerr("***** oDynamicDataManager::InitDynamicData()-> A problem occured during the dynamic gating initialization !" << endl);
00176       return 1;
00177     }
00178   }
00179 
00180   // Initialization with involuntary patient motion corrected reconstruction
00181   if (m_pMotionCorrFlag && InitDynamicDataPatientMotion(a_pathTo4DDataFile) )
00182   {
00183     Cerr("***** oDynamicDataManager::InitDynamicData()-> A problem occured during the patient involuntary motion correction initialization !" << endl);
00184     return 1;
00185   }
00186 
00187   // Some feedback
00188   if (m_verbose>=3) 
00189   {
00190     if(m_respGatingFlag)  Cout("oDynamicDataManager::InitDynamicData()-> Respiratory gating enabled" << endl);
00191     if(m_rMotionCorrFlag) Cout("oDynamicDataManager::InitDynamicData()-> Respiratory image-based motion correction enabled" << endl);
00192     if(m_cardGatingFlag)  Cout("oDynamicDataManager::InitDynamicData()-> Cardiac gating enabled" << endl);
00193     if(m_cMotionCorrFlag) Cout("oDynamicDataManager::InitDynamicData()-> Cardiac image-based motion correction enabled" << endl);
00194     if(m_pMotionCorrFlag) Cout("oDynamicDataManager::InitDynamicData()-> Involuntary image-based patient motion correction enabled" << endl);
00195   }
00196     
00197   // End
00198   return 0;
00199 }
00200 
00201 // =====================================================================
00202 // ---------------------------------------------------------------------
00203 // ---------------------------------------------------------------------
00204 // =====================================================================
00205 /*
00206   \fn InitDynamicDataGating
00207   \param a_pathToFile : path to an ASCII file containing dynamic metadata regarding the acquisition
00208   \brief Initialisation of arrays containing informations about the data splitting and respiratory/cardiac gated reconstruction
00209   \todo check the reconstruction with cardiac gated reconstruction. Some adjustements will be required for double gating
00210   \return 0 if success, positive value otherwise.
00211 */
00212 int oDynamicDataManager::InitDynamicDataGating(const string& a_pathToFile)
00213 {
00214   if(m_verbose >=3) Cout("oDynamicDataManager::InitDynamicDataGating() ... " << endl);
00215     
00216   // Respiratory gating enabled
00217   if (m_respGatingFlag)
00218   {
00219     // Read information about respiratory gating (number of events per respiratory gate per frame)
00220     if ( ReadDataASCIIFile(a_pathToFile, 
00221                           "respiratory_gates_data_splitting", 
00222                            m2p_nbEventsPerRespGate, 
00223                            m_nbRespGates, 
00224                            m_nbTimeFrames, 
00225                            KEYWORD_MANDATORY) )
00226     {
00227       Cerr("***** oDynamicDataManager::InitDynamicDataGating() -> Error while trying to read 'respiratory_gates_data_splitting' in file '" << a_pathToFile << "' !" << endl);
00228       return 1;
00229     }
00230     
00231     // Get the index of the last event of each respiratory gate using the number of events inside each gate
00232     uint64_t event_number_sum = 0;
00233     if (m_verbose>=3) Cout("oDynamicDataManager::InitDynamicDataGating() :" << endl);
00234             
00235     for (int fr=0; fr<m_nbTimeFrames; fr++)
00236       for (int rg=0; rg<m_nbRespGates; rg++)
00237       { 
00238         m2p_indexLastEventRespGate[fr][rg] += m2p_nbEventsPerRespGate[fr][rg] + event_number_sum;
00239         event_number_sum += m2p_nbEventsPerRespGate[fr][rg]; 
00240         if (m_verbose>=3) Cout("Number of events for frame #" << fr << ", respiratory gate #" << rg << " = " << m2p_nbEventsPerRespGate[fr][rg] << endl);
00241       } 
00242   }
00243   
00244   
00245   // Cardiac gating enabled
00246   if (m_cardGatingFlag)
00247   {
00248     // Read information about cardiac gating (number of events per cardiac gate per respiratory gates times frames)
00249     if ( ReadDataASCIIFile(a_pathToFile, 
00250                           "cardiac_gates_data_splitting", 
00251                            m2p_nbEventsPerCardGate,
00252                            m_nbCardGates, 
00253                            m_nbTimeFrames*m_nbRespGates, 
00254                            KEYWORD_MANDATORY) )
00255     {
00256       Cerr("***** oDynamicDataManager::InitDynamicDataGating() -> Error while trying to read 'cardiac_gates_data_splitting' in file '" << a_pathToFile << "' !" << endl);
00257       return 1;
00258     }
00259       
00260     // Get the index of the last event of each cardiac gate using the number of events inside each gate
00261     uint64_t event_number_sum = 0;
00262     if (m_verbose >=3) Cout("oDynamicDataManager::InitDynamicDataGating() :" << endl);
00263         
00264     for (int fr=0; fr<m_nbTimeFrames; fr++)
00265       for (int g=0; g<m_nbRespGates*m_nbCardGates; g++)
00266       { 
00267         m2p_indexLastEventCardGate[fr][g] += m2p_nbEventsPerCardGate[fr][g] + event_number_sum;
00268         event_number_sum += m2p_nbEventsPerCardGate[fr][g]; 
00269         if (m_verbose>=3) Cout("Number of events for frame #" << fr << ", cardiac gate #" << g << " = " << m2p_nbEventsPerCardGate[fr][g] << endl);
00270       } 
00271   }
00272   
00273   // End
00274   return 0;
00275 }
00276 
00277 // =====================================================================
00278 // ---------------------------------------------------------------------
00279 // ---------------------------------------------------------------------
00280 // =====================================================================
00281 /*
00282   \fn InitDynamicDataPatientMotion
00283   \param a_pathToFile : path to an ASCII file containing dynamic metadata regarding the acquisition
00284   \brief Initialisation of involuntary patient motion correction information, if any.
00285   \todo Implementation currently IN PROGRESS. Will have to define how motion subset information is provided as it depends of each datafile, if the acquisition contain several bed steps.
00286   \return 0 if success, positive value otherwise.
00287 */
00288 int oDynamicDataManager::InitDynamicDataPatientMotion(const string& a_pathToFile)
00289 {
00290   if(m_verbose >=3) Cout("oDynamicDataManager::InitDynamicDataPatientMotion() ... " << endl);
00291   Cerr("!!!!! WARNING oDynamicDataManager::InitDynamicDataGating() -> Still WiP !" << endl);
00292   // First get the number of time triggers into the dynamic file
00293   if (ReadDataASCIIFile(a_pathToFile, 
00294                        "nb_involuntary_motion_triggers", 
00295                         &m_nbPMotionTriggers, 
00296                         1, 
00297                         KEYWORD_MANDATORY))
00298     {
00299       Cerr("***** oDynamicDataManager::InitDynamicDataPatientMotion() -> Involuntary motion correction enabled but could" << endl
00300         << " not find the number of triggers in the dynamic file '" << a_pathToFile << "' !" << endl);
00301       return 1;
00302     }
00303   
00304   // Allocate time triggers and read them from the dynamic file
00305   mp_listPMotionTriggers = (uint32_t*)malloc(m_nbPMotionTriggers*sizeof(uint32_t));
00306 
00307   if (ReadDataASCIIFile(a_pathToFile, 
00308                       "list_involuntary_motion_triggers", 
00309                         mp_listPMotionTriggers, 
00310                         m_nbPMotionTriggers, 
00311                         KEYWORD_MANDATORY))
00312   {
00313     Cerr("***** oDynamicDataManager::InitDynamicDataPatientMotion() -> Involuntary motion correction enabled but list of triggers" << endl
00314       << "                                                             not found in the dynamic file '" << a_pathToFile << "' !" << endl);
00315     return 1;
00316   }
00317   // End
00318   return 0;
00319 }
00320 
00321 // =====================================================================
00322 // ---------------------------------------------------------------------
00323 // ---------------------------------------------------------------------
00324 // =====================================================================
00325 /*
00326   \fn SetDynamicSpecificQuantificationFactors
00327   \param FLTNB** a2p_quantificationFactors : 2 dimensional [timeframe][gate] set of quantitative factors to update
00328   \brief Compute gate-specific quantificative factors using the number of events within each gate, and update the quantitative factors passed in argument
00329   \return 0 if success, positive value otherwise.
00330 */
00331 int oDynamicDataManager::SetDynamicSpecificQuantificationFactors(FLTNB** a2p_quantificationFactors)
00332 {
00333   if(m_verbose >=3) Cout("oDynamicDataManager::SetDynamicSpecificQuantificationFactors() ... " << endl);
00334     
00335   // COMPUTE GATE-SPECIFIC QUANTITATIVE FACTORS
00336   if (m_nbRespGates>1 || m_nbCardGates>1)
00337   {
00338     for(int fr=0 ; fr<m_nbTimeFrames ; fr++)
00339     {
00340       // We have cardiac gates and don't correct for motion -> quantification factors required for cardiac gates
00341       if (m_cardGatingFlag && !m_cMotionCorrFlag)
00342       {
00343         uint64_t total_events_in_frame = 0;
00344         for(int g=0 ; g<m_nbRespGates*m_nbCardGates ; g++) total_events_in_frame += m2p_nbEventsPerCardGate[fr][g];
00345 
00346         if (m_verbose>=3) Cout("oDynamicDataManager::SetDynamicSpecificQuantificationFactors() -> Cardiac gating correction factors :" << endl);
00347         for(int g=0 ; g<m_nbRespGates*m_nbCardGates ; g++)
00348         {
00349           a2p_quantificationFactors[fr][g] *= (FLTNB)total_events_in_frame/(FLTNB)m2p_nbEventsPerCardGate[fr][g];
00350           if (m_verbose>=3) Cout("Frame #" << fr << ", cardiac gate #" << g << " = " << a2p_quantificationFactors[fr][g] << endl);
00351         }
00352       }
00353       
00354       // We have resp gates and don't correct for resp motion (and no independent cardiac gate reconstruction) -> quantification factors required for cardiac gates
00355       else if(m_respGatingFlag && !m_rMotionCorrFlag)
00356       {
00357         uint64_t total_events_in_frame = 0;
00358         for(int rg=0 ; rg<m_nbRespGates ; rg++) total_events_in_frame += m2p_nbEventsPerRespGate[fr][rg];
00359 
00360         if (m_verbose>=3) Cout("oDynamicDataManager::SetDynamicSpecificQuantificationFactors() -> Respiratory gating correction factors :" << endl);
00361         for(int g=0 ; g<m_nbRespGates*m_nbCardGates ; g++)
00362         {
00363           int rg = int(g/m_nbCardGates);
00364           a2p_quantificationFactors[fr][g] *= (FLTNB)total_events_in_frame/(FLTNB)m2p_nbEventsPerRespGate[fr][rg];
00365           if (m_verbose>=3) Cout("Frame #" << fr << ", gate #" << g << " = " << a2p_quantificationFactors[fr][g] << endl);
00366         }
00367       }
00368     }
00369   }
00370   
00371   // End
00372   return 0;
00373 }
00374 
00375 // =====================================================================
00376 // ---------------------------------------------------------------------
00377 // ---------------------------------------------------------------------
00378 // =====================================================================
00385 int oDynamicDataManager::CheckParameters(int64_t a_nbEvents)
00386 {
00387   if (m_verbose>=3) Cout("oDynamicDataManager::CheckParameters() -> Check parameters for dynamic data settings" << endl);
00388     
00389   // Check image dimensions
00390   if (mp_ID==NULL)
00391   {
00392     Cerr("***** oDynamicDataManager::CheckParameters() -> No image dimensions provided !" << endl);
00393     return 1;
00394   }
00395   // Check resp gates
00396   if (m_nbRespGates<0)
00397   {
00398     Cerr("***** oDynamicDataManager::CheckParameters() -> Wrong number of respiratory gates !" << endl);
00399     return 1;
00400   }
00401   // Check card gates
00402   if (m_nbCardGates<0)
00403   {
00404     Cerr("***** oDynamicDataManager::CheckParameters() -> Wrong number of respiratory gates !" << endl);
00405     return 1;
00406   }
00407   // Check involuntary motion triggers
00408   if (m_nbPMotionTriggers<0)
00409   {
00410     Cerr("***** oDynamicDataManager::CheckParameters() -> Wrong number of involuntary motion subsets provided !" << endl);
00411     return 1;
00412   }
00413   // Check verbosity
00414   if (m_verbose<0)
00415   {
00416     Cerr("***** oDynamicDataManager::CheckParameters() -> Wrong verbosity level provided !" << endl);
00417     return 1;
00418   }
00419   // Feedback regarding the enabled/disabled options for the reconstruction
00420   if (m_verbose>=2)
00421   {
00422     if (m_respGatingFlag) Cout("oDynamicDataManager::CheckParameters() -> Respiratory gating is enabled" << endl);
00423     if (m_rMotionCorrFlag) Cout("oDynamicDataManager::CheckParameters() -> Respiratory motion correction enabled" << endl);
00424     if (m_cardGatingFlag) Cout("oDynamicDataManager::CheckParameters() -> Cardiac gating is enabled" << endl);
00425     if (m_cMotionCorrFlag) Cout("oDynamicDataManager::CheckParameters() -> Cardiac motion correction is enabled" << endl);
00426     if (m_pMotionCorrFlag) Cout("oDynamicDataManager::CheckParameters() -> Involuntary motion correction is enabled" << endl);
00427   }
00428 
00429   // Check consistency between number of events in the datafile and the potential splitting of the dynamic data into respiratory or cardiac gates (if any gating is enabled)
00430   if (m_respGatingFlag)
00431   {
00432     int last_fr = m_nbTimeFrames-1;
00433     int last_rg = m_nbRespGates-1;
00434     if (m2p_indexLastEventRespGate[last_fr][last_rg]+1 != a_nbEvents)
00435     {
00436       Cerr("***** oDynamicDataManager::CheckParameters() -> Problem while checking consistency of dynamic data !" << endl
00437         << "                                                The number of events in the datafile (" << a_nbEvents
00438         << ") is different from the total number of events in respiratory gates (" << m2p_indexLastEventRespGate[last_fr][last_rg] << ") as initialized in the gating file !" << endl);
00439       return 1;
00440     }
00441   }
00442   if (m_cardGatingFlag)
00443   {
00444     int last_fr = m_nbTimeFrames-1;
00445     int last_rg = m_nbCardGates-1;
00446     if (m2p_indexLastEventCardGate[last_fr][last_rg]+1 != a_nbEvents)
00447     {
00448       Cerr("***** oDynamicDataManager::CheckParameters() -> Problem while checking consistency of dynamic data !" << endl
00449         << "                                                The number of events in the datafile (" << a_nbEvents
00450         << ") is different to the total number of events in cardiac gates (" << m2p_indexLastEventRespGate[last_fr][last_rg] << ") as initialized in the gating file !" << endl);
00451       return 1;
00452     }
00453   }
00454   
00455   // Normal end
00456   return 0;
00457 }
00458 
00459 // =====================================================================
00460 // ---------------------------------------------------------------------
00461 // ---------------------------------------------------------------------
00462 // =====================================================================
00467 void oDynamicDataManager::ResetCurrentDynamicIndices()
00468 {
00469   #ifdef CASTOR_VERBOSE
00470   if (m_verbose >=4) Cout("oDynamicDataManager::ResetCurrentDynamicIndices() -> Reset indices" << endl);
00471   #endif
00472   
00473   for (int th=0; th<mp_ID->GetNbThreadsForProjection(); th++)
00474   {
00475     // We initialize the frame index to -1 because at the beginning of the events loop, we don't know yet
00476     // if we are in the first frame (the user can reconstruct only a part of the whole datafile).
00477     mp_currentFrameIndex[th] = -1;
00478     mp_currentPMotionIndex[th] = 0;
00479     mp_currentRespGateIndex[th] = 0;
00480     mp_currentCardGateIndex[th] = 0;
00481   }
00482 }
00483 
00484 // =====================================================================
00485 // ---------------------------------------------------------------------
00486 // ---------------------------------------------------------------------
00487 // =====================================================================
00488 /*
00489   \fn DynamicSwitch
00490   \param a_currentEventIndex : Index of the current event in the iterative reconstruction loop
00491   \param a_currentTime : Timestamp of the current event in the iterative reconstruction loop
00492   \param a_bed
00493   \param a_th
00494   \brief This function is called in the reconstruction event loop. It is used to check if the current event belongs to a new respiratory/cardiac/involuntary motion gate
00495          Increment the index reporting the current relative gate in this case
00496   \todo  Check implementation for cardiac gating and involuntary patient motion correction
00497   \return 0 if nothing to do,
00498           1 if the current timestamp is inferior to the start time of the first frame,
00499           2 if a deformation is required (gate has changed and motion correction is enabled)
00500 */
00501 int oDynamicDataManager::DynamicSwitch(int64_t a_currentEventIndex, uint32_t a_currentTime, int a_bed, int a_th)
00502 {
00503   #ifdef CASTOR_VERBOSE
00504   if (m_verbose >=4) Cout("oDynamicDataManager::DynamicSwitch() -> Enter function" << endl);
00505   #endif
00506   
00507   // The value that will be returned
00508   int return_value = DYNAMIC_SWITCH_NOTHING;
00509 
00510   // -------------------------------------------------------------------------------------------------
00511   // Step 1: involuntary motion management
00512   // -------------------------------------------------------------------------------------------------
00513 
00514   // Only do this if the involuntary motion correction is enabled (meaning we must proceed to a deformation)
00515   if (m_pMotionCorrFlag)
00516   {
00517     // Search if we passed one of the next motion triggers (starting at current index)
00518     for (int mt=mp_currentPMotionIndex[a_th]; mt<m_nbPMotionTriggers; mt++)
00519     {
00520       // If we passed this trigger, set the return value to DEFORMATION. However, we continue the loop
00521       // in the case where we passed multiple triggers.
00522       if (a_currentTime>=mp_listPMotionTriggers[mt])
00523       {
00524         mp_currentPMotionIndex[a_th] = mt+1;
00525         #ifdef CASTOR_VERBOSE
00526         if (m_verbose >=4) Cout("oDynamicDataManager::DynamicSwitch() -> Thread " << a_th << ", gate for patient motion correction switched to " << mp_currentPMotionIndex[a_th] << endl
00527                              << "                                        at event #" << a_currentEventIndex << ", timestamp = " << a_currentTime << endl);
00528         #endif
00529         return_value = DYNAMIC_SWITCH_DEFORMATION;
00530       }
00531     }
00532   }
00533 
00534   // -------------------------------------------------------------------------------------------------
00535   // Step 2: frame management
00536   // -------------------------------------------------------------------------------------------------
00537 
00538   // Special case if the frame number is still -1
00539   if (mp_currentFrameIndex[a_th]==-1)
00540   {
00541     // We check if we are before the first frame, in this case we directly return a CONTINUE signal to go to the next event
00542     if (a_currentTime<mp_ID->GetFrameTimeStartInMs(a_bed,0)) return DYNAMIC_SWITCH_CONTINUE;
00543     // Otherwise, we now are at least in the first frame (which will be updated right after this if section)
00544     else mp_currentFrameIndex[a_th] = 0;
00545   }
00546 
00547   // A boolean to know later if the frame index has changed
00548   bool frame_has_changed = false;
00549 
00550   // Now we search for the first frame index in which the event belongs, starting from the current frame index. Note that we do that
00551   // this way (not simply incrementing the frame index) because we want to deal with the case where the next event managed by this
00552   // thread could possibly skip multiple frames at once.
00553   for (int fr=mp_currentFrameIndex[a_th]; fr<m_nbTimeFrames; fr++)
00554   {
00555     // If the current time is less than the time stop of the frame 'fr', then the event belongs to this frame
00556     if (a_currentTime<mp_ID->GetFrameTimeStopInMs(a_bed,fr))
00557     {
00558       // Test if the frame has actually changed
00559       if (mp_currentFrameIndex[a_th]!=fr)
00560       {
00561         // Set the boolean to true
00562         frame_has_changed = true;
00563         // Update the current frame index
00564         mp_currentFrameIndex[a_th] = fr;
00565       }
00566       break;
00567     }
00568   }
00569 
00570   #ifdef CASTOR_VERBOSE
00571   if (frame_has_changed && m_verbose >=4)
00572     Cout("oDynamicDataManager::DynamicSwitch() -> Thread " << a_th << ", frame switched to " << mp_currentFrameIndex[a_th] << endl);
00573   #endif
00574 
00575   // -------------------------------------------------------------------------------------------------
00576   // Step 3: respiratory gating management
00577   // -------------------------------------------------------------------------------------------------
00578 
00579   // A boolean to know later if the respiratory index has changed
00580   bool resp_gate_has_changed = false;
00581 
00582   // Do this only if respiratory gating is enabled
00583   if (m_respGatingFlag)
00584   {
00585     // Test if the frame index has changed
00586     if (frame_has_changed)
00587     {
00588       // Reset the respiratory gate index
00589       mp_currentRespGateIndex[a_th] = 0;
00590       // No deformation signal need to be sent, as the forward/backward image matrices for the next frame will be in the reference position as required
00591     }
00592     // For this frame, search the first gate (from the current gate index) for which the provided event index is below the event-index-stop
00593     for (int rg=mp_currentRespGateIndex[a_th]; rg<m_nbRespGates; rg++)
00594     {
00595       // If the current event index is below the last event of this gate, then the event belongs to this gate
00596       // (We won't enter again in the if statement due to the flag setting to true)
00597       if (a_currentEventIndex<=m2p_indexLastEventRespGate[mp_currentFrameIndex[a_th]][rg] && resp_gate_has_changed == false)
00598       {
00599         // Test if the gate has actually changed
00600         if (mp_currentRespGateIndex[a_th]!=rg)
00601         {
00602           // Update the current gate index
00603           mp_currentRespGateIndex[a_th] = rg;
00604           // Verbose
00605           #ifdef CASTOR_VERBOSE
00606           if (m_verbose>=4) Cout("oDynamicDataManager::DynamicSwitch() -> Thread " << a_th << ", respiratory gate switch to " << mp_currentRespGateIndex[a_th]
00607                                                                                    << " on event " << a_currentEventIndex << endl
00608                                                                                    << "current frame : " << mp_currentFrameIndex[a_th] << endl
00609                                                                                    << "current respiratory gate index " << mp_currentRespGateIndex[a_th] << endl);
00610           #endif
00611           // If motion correction is enabled, then we should return a DEFORMATION signal
00612           if (m_rMotionCorrFlag) return_value = DYNAMIC_SWITCH_DEFORMATION;
00613         }
00614         // Set the boolean to true
00615         resp_gate_has_changed = true;
00616       }
00617     }
00618   }
00619   
00620   // -------------------------------------------------------------------------------------------------
00621   // Step 4: cardiac gating management
00622   // -------------------------------------------------------------------------------------------------
00623 
00624   // A boolean to know later if the respiratory index has changed
00625   bool card_gate_has_changed = false;
00626   
00627   // Do this only if cardiac gating is enabled
00628   if (m_cardGatingFlag)
00629   {
00630     // Test if the frame or respiratory index have changed
00631     if (frame_has_changed || resp_gate_has_changed)
00632     {
00633       // Reset the cardiac gate index
00634       mp_currentCardGateIndex[a_th] = 0;
00635       // Thus if one apply cardiac motion correction, we should return a DEFORMATION signal
00636       if (m_cMotionCorrFlag) return_value = DYNAMIC_SWITCH_DEFORMATION;
00637     }
00638     // For this frame and respiratory gate, search the first gate (from the current gate index) for which the provided event index is below the event-index-stop
00639     for (int cg=mp_currentCardGateIndex[a_th]; cg<m_nbCardGates; cg++)
00640     {
00641       // If the current event index is below the event-index-stop of this gate, then the event belongs to this gate
00642       if (a_currentEventIndex<m2p_indexLastEventCardGate[mp_currentFrameIndex[a_th]][mp_currentRespGateIndex[a_th]*m_nbCardGates+cg]  && card_gate_has_changed == false)
00643       {
00644         // Test if the gate has actually changed
00645         if (mp_currentCardGateIndex[a_th]!=cg)
00646         {
00647           // Update the current gate index
00648           mp_currentCardGateIndex[a_th] = cg;
00649           // Verbose
00650           #ifdef CASTOR_VERBOSE
00651           if (m_verbose>=4) Cout("oDynamicDataManager::DynamicSwitch() -> Thread " << a_th << ", cardiac gate switch to " <<  mp_currentCardGateIndex[a_th]
00652                                                                                    << " on event " << a_currentEventIndex << endl
00653                                                                                    << "current frame : " << mp_currentFrameIndex[a_th] << endl);
00654           #endif
00655           // If motion correction is enabled, then we should return a DEFORMATION signal
00656           if (m_cMotionCorrFlag) return_value = DYNAMIC_SWITCH_DEFORMATION;
00657         }
00658         // Set the boolean to true
00659         card_gate_has_changed = true;
00660       }
00661     }
00662   }
00663   
00664   // -------------------------------------------------------------------------------------------------
00665   // Step 5: just return the value !
00666   // -------------------------------------------------------------------------------------------------
00667 
00668   // Return the status of the dynamic switch
00669   return return_value;
00670 }
 All Classes Files Functions Variables Typedefs Defines