![]() |
CASToR
1.0
Tomographic Reconstruction (PET/SPECT)
|
00001 00002 /* 00003 Implementation of class iOptimizerNEGML 00004 00005 - separators: X 00006 - doxygen: X 00007 - default initialization: X 00008 - CASTOR_DEBUG: 00009 - CASTOR_VERBOSE: 00010 */ 00011 00018 #include "iOptimizerNEGML.hh" 00019 #include "sOutputManager.hh" 00020 00021 // ===================================================================== 00022 // --------------------------------------------------------------------- 00023 // --------------------------------------------------------------------- 00024 // ===================================================================== 00025 00026 iOptimizerNEGML::iOptimizerNEGML() : vOptimizer() 00027 { 00028 // --------------------------- 00029 // Mandatory member parameters 00030 // --------------------------- 00031 00032 // Initial value at 1 00033 m_initialValue = 1.; 00034 // Only one backward image for NEGML 00035 m_nbBackwardImages = 1; 00036 // NEGML does not accept penalties 00037 //m_penaltyEnergyFunctionDerivativesOrder = 0; 00038 // NEGML is only compatible with histogram data 00039 m_listmodeCompatibility = false; 00040 m_histogramCompatibility = true; 00041 00042 // -------------------------- 00043 // Specific member parameters 00044 // -------------------------- 00045 00046 // This is the NEGML psi parameter 00047 m_psi = -1.; 00048 } 00049 00050 // ===================================================================== 00051 // --------------------------------------------------------------------- 00052 // --------------------------------------------------------------------- 00053 // ===================================================================== 00054 00055 iOptimizerNEGML::~iOptimizerNEGML() 00056 { 00057 } 00058 00059 // ===================================================================== 00060 // --------------------------------------------------------------------- 00061 // --------------------------------------------------------------------- 00062 // ===================================================================== 00063 00064 void iOptimizerNEGML::ShowHelpSpecific() 00065 { 00066 cout << "This optimizer is the NEGML algorithm from K. Van Slambrouck et al, IEEE TMI, Jan 2015, vol. 34, pp. 126-136." << endl; 00067 cout << "Subsets can be used. This implementation only consider the psi parameter, but not the alpha image design parameter" << endl; 00068 cout << "which is supposed to be 1 for all voxels. It implements equation 17 of the reference paper." << endl; 00069 cout << "The following options can be used (in this particular order when provided as a list):" << endl; 00070 cout << " initial image value: to set the uniform voxel value for the initial image" << endl; 00071 cout << " psi: to set the psi parameter that sets the transition from Poisson to Gaussian statistics (must be positive)." << endl; 00072 cout << " (if set to 0, then it is taken to infinity and implements equation 21 in the reference paper)." << endl; 00073 } 00074 00075 // ===================================================================== 00076 // --------------------------------------------------------------------- 00077 // --------------------------------------------------------------------- 00078 // ===================================================================== 00079 00080 int iOptimizerNEGML::ReadConfigurationFile(const string& a_configurationFile) 00081 { 00082 string key_word = ""; 00083 // Read the initial image value option 00084 key_word = "initial image value"; 00085 if (ReadDataASCIIFile(a_configurationFile, key_word, &m_initialValue, 1, KEYWORD_MANDATORY)) 00086 { 00087 Cerr("***** iOptimizerNEGML::ReadAndCheckConfigurationFile() -> Failed to get the '" << key_word << "' keyword !" << endl); 00088 return 1; 00089 } 00090 // Read the psi value option 00091 key_word = "psi"; 00092 if (ReadDataASCIIFile(a_configurationFile, key_word, &m_psi, 1, KEYWORD_MANDATORY)) 00093 { 00094 Cerr("***** iOptimizerNEGML::ReadAndCheckConfigurationFile() -> Failed to get the '" << key_word << "' keyword !" << endl); 00095 return 1; 00096 } 00097 // Normal end 00098 return 0; 00099 } 00100 00101 // ===================================================================== 00102 // --------------------------------------------------------------------- 00103 // --------------------------------------------------------------------- 00104 // ===================================================================== 00105 00106 int iOptimizerNEGML::ReadOptionsList(const string& a_optionsList) 00107 { 00108 // There are 2 floating point variables as options 00109 FLTNB options[2]; 00110 // Read them 00111 if (ReadStringOption(a_optionsList, options, 2, ",", "NEGML configuration")) 00112 { 00113 Cerr("***** iOptimizerNEGML::ReadAndCheckConfigurationFile() -> Failed to correctly read the list of options !" << endl); 00114 return 1; 00115 } 00116 // Affect options 00117 m_initialValue = options[0]; 00118 m_psi = options[1]; 00119 // Normal end 00120 return 0; 00121 } 00122 00123 // ===================================================================== 00124 // --------------------------------------------------------------------- 00125 // --------------------------------------------------------------------- 00126 // ===================================================================== 00127 00128 int iOptimizerNEGML::CheckSpecificParameters() 00129 { 00130 // Check that the psi value is positive (a zero value means we take psi to infinity) 00131 if (m_psi<0.) 00132 { 00133 Cerr("***** iOptimizerNEGML::Initialize() -> Provided psi value (" << m_psi << ") must be positive !" << endl); 00134 return 1; 00135 } 00136 // Normal end 00137 return 0; 00138 } 00139 00140 // ===================================================================== 00141 // --------------------------------------------------------------------- 00142 // --------------------------------------------------------------------- 00143 // ===================================================================== 00144 00145 int iOptimizerNEGML::InitializeSpecific() 00146 { 00147 // Verbose 00148 if (m_verbose>=1) 00149 { 00150 Cout("iOptimizerNEGML::Initialize() -> Use the NEGML optimizer" << endl); 00151 if (m_verbose>=2) 00152 { 00153 Cout(" --> Initial image value: " << m_initialValue << endl); 00154 if (m_psi>0.) Cout(" --> Psi: " << m_psi << endl); 00155 else if (m_psi==0.) Cout(" --> Psi taken to infinity" << endl); 00156 } 00157 } 00158 // Normal end 00159 return 0; 00160 } 00161 00162 // ===================================================================== 00163 // --------------------------------------------------------------------- 00164 // --------------------------------------------------------------------- 00165 // ===================================================================== 00166 00167 int iOptimizerNEGML::SensitivitySpecificOperations( FLTNB a_data, FLTNB a_forwardModel, FLTNB* ap_weight, 00168 FLTNB a_multiplicativeCorrections, FLTNB a_additiveCorrections, 00169 FLTNB a_quantificationFactor, oProjectionLine* ap_Line ) 00170 { 00171 // With psi taken to infinity, the weight is the forward projection of a uniform 1 image 00172 if (m_psi==0.) *ap_weight = ForwardProject(ap_Line); 00173 // Line weight here is a forward projection of alpha image divided by the max between psi and the forward model 00174 else *ap_weight = ForwardProject(ap_Line) / max( m_psi , a_forwardModel ); 00175 // That's all 00176 return 0; 00177 } 00178 00179 // ===================================================================== 00180 // --------------------------------------------------------------------- 00181 // --------------------------------------------------------------------- 00182 // ===================================================================== 00183 00184 int iOptimizerNEGML::DataSpaceSpecificOperations( FLTNB a_data, FLTNB a_forwardModel, FLTNB* ap_backwardValues, 00185 FLTNB a_multiplicativeCorrections, FLTNB a_additiveCorrections, 00186 FLTNB a_quantificationFactor, oProjectionLine* ap_Line ) 00187 { 00188 // Compute numerator 00189 FLTNB numerator = a_data - a_forwardModel; 00190 // Compute denominator that will be strictly positive 00191 FLTNB denominator = 1.; 00192 if (m_psi>0.) denominator = max( m_psi , a_forwardModel ); 00193 // Update backward values 00194 *ap_backwardValues = numerator / denominator; 00195 // That's all 00196 return 0; 00197 } 00198 00199 // ===================================================================== 00200 // --------------------------------------------------------------------- 00201 // --------------------------------------------------------------------- 00202 // ===================================================================== 00203 00204 int iOptimizerNEGML::ImageSpaceSpecificOperations( FLTNB a_currentImageValue, FLTNB* ap_newImageValue, 00205 FLTNB a_sensitivity, FLTNB* ap_correctionValues ) 00206 { 00207 // Compute image update factor 00208 FLTNB image_update_factor = *ap_correctionValues / a_sensitivity; 00209 // Update image 00210 *ap_newImageValue = a_currentImageValue + image_update_factor; 00211 // End 00212 return 0; 00213 } 00214 00215 // ===================================================================== 00216 // --------------------------------------------------------------------- 00217 // --------------------------------------------------------------------- 00218 // =====================================================================