![]() |
CASToR
1.0
Tomographic Reconstruction (PET/SPECT)
|
00001 00002 /* 00003 Implementation of class iOptimizerLandweber 00004 00005 - separators: X 00006 - doxygen: X 00007 - default initialization: X 00008 - CASTOR_DEBUG: 00009 - CASTOR_VERBOSE: 00010 */ 00011 00018 #include "iOptimizerLandweber.hh" 00019 #include "sOutputManager.hh" 00020 00021 // ===================================================================== 00022 // --------------------------------------------------------------------- 00023 // --------------------------------------------------------------------- 00024 // ===================================================================== 00025 00026 iOptimizerLandweber::iOptimizerLandweber() : vOptimizer() 00027 { 00028 // --------------------------- 00029 // Mandatory member parameters 00030 // --------------------------- 00031 00032 // Initial value at 0 00033 m_initialValue = 0.; 00034 // Only one backward image for Landweber 00035 m_nbBackwardImages = 1; 00036 // Landweber does not accept penalties 00037 //m_penaltyEnergyFunctionDerivativesOrder = 0; 00038 // Landweber is not compatible with listmode data 00039 m_listmodeCompatibility = false; 00040 // Landweber is compatible with histogram data only 00041 m_histogramCompatibility = true; 00042 00043 // -------------------------- 00044 // Specific member parameters 00045 // -------------------------- 00046 00047 m_relaxationFactor = -1.; 00048 } 00049 00050 // ===================================================================== 00051 // --------------------------------------------------------------------- 00052 // --------------------------------------------------------------------- 00053 // ===================================================================== 00054 00055 iOptimizerLandweber::~iOptimizerLandweber() 00056 { 00057 } 00058 00059 // ===================================================================== 00060 // --------------------------------------------------------------------- 00061 // --------------------------------------------------------------------- 00062 // ===================================================================== 00063 00064 void iOptimizerLandweber::ShowHelpSpecific() 00065 { 00066 cout << "This optimizer implements the standard Landweber algorithm for least-squares optimization." << endl; 00067 cout << "No penalty can be used with this algorithm." << endl; 00068 cout << "Be aware that the relaxation parameter is not automatically set, so it often requires some" << endl; 00069 cout << "trials and errors to find an optimal setting. Also, remember that this algorithm is particularly" << endl; 00070 cout << "slow to converge." << endl; 00071 cout << "The following options can be used (in this particular order when provided as a list):" << endl; 00072 cout << " initial image value: to set the uniform voxel value for the initial image" << endl; 00073 cout << " relaxation factor: to set the relaxation factor applied to the update" << endl; 00074 } 00075 00076 // ===================================================================== 00077 // --------------------------------------------------------------------- 00078 // --------------------------------------------------------------------- 00079 // ===================================================================== 00080 00081 int iOptimizerLandweber::ReadConfigurationFile(const string& a_configurationFile) 00082 { 00083 string key_word = ""; 00084 // Read the initial image value option 00085 key_word = "initial image value"; 00086 if (ReadDataASCIIFile(a_configurationFile, key_word, &m_initialValue, 1, KEYWORD_MANDATORY)) 00087 { 00088 Cerr("***** iOptimizerLandweber::ReadAndCheckConfigurationFile() -> Failed to get the '" << key_word << "' keyword !" << endl); 00089 return 1; 00090 } 00091 // Read the relaxation factor option 00092 key_word = "relaxation factor"; 00093 if (ReadDataASCIIFile(a_configurationFile, key_word, &m_relaxationFactor, 1, KEYWORD_MANDATORY)) 00094 { 00095 Cerr("***** iOptimizerLandweber::ReadAndCheckConfigurationFile() -> Failed to get the '" << key_word << "' keyword !" << endl); 00096 return 1; 00097 } 00098 // Normal end 00099 return 0; 00100 } 00101 00102 // ===================================================================== 00103 // --------------------------------------------------------------------- 00104 // --------------------------------------------------------------------- 00105 // ===================================================================== 00106 00107 int iOptimizerLandweber::ReadOptionsList(const string& a_optionsList) 00108 { 00109 // There are 2 floating point variables as options 00110 FLTNB options[4]; 00111 00112 // Read them 00113 if (ReadStringOption(a_optionsList, options, 2, ",", "Landweber configuration")) 00114 { 00115 Cerr("***** iOptimizerLandweber::ReadAndCheckConfigurationFile() -> Failed to correctly read the list of options !" << endl); 00116 return 1; 00117 } 00118 // Affect options 00119 m_initialValue = options[0]; 00120 m_relaxationFactor = options[1]; 00121 00122 // Normal end 00123 return 0; 00124 } 00125 00126 // ===================================================================== 00127 // --------------------------------------------------------------------- 00128 // --------------------------------------------------------------------- 00129 // ===================================================================== 00130 00131 int iOptimizerLandweber::CheckSpecificParameters() 00132 { 00133 // Check that relaxation factor value is strictly positive 00134 if (m_relaxationFactor<=0.) 00135 { 00136 Cerr("***** iOptimizerLandweber->Initialize() -> Provided relaxation factor (" << m_relaxationFactor << ") must be strictly positive !" << endl); 00137 return 1; 00138 } 00139 // Normal end 00140 return 0; 00141 } 00142 00143 // ===================================================================== 00144 // --------------------------------------------------------------------- 00145 // --------------------------------------------------------------------- 00146 // ===================================================================== 00147 00148 int iOptimizerLandweber::InitializeSpecific() 00149 { 00150 // Verbose 00151 if (m_verbose>=1) 00152 { 00153 Cout("iOptimizerLandweber::Initialize() -> Use the Landweber algorithm" << endl); 00154 if (m_verbose>=2) 00155 { 00156 Cout(" --> Initial image value: " << m_initialValue << endl); 00157 Cout(" --> Relaxation factor: " << m_relaxationFactor << endl); 00158 } 00159 } 00160 // Normal end 00161 return 0; 00162 } 00163 00164 // ===================================================================== 00165 // --------------------------------------------------------------------- 00166 // --------------------------------------------------------------------- 00167 // ===================================================================== 00168 00169 int iOptimizerLandweber::SensitivitySpecificOperations( FLTNB a_data, FLTNB a_forwardModel, FLTNB* ap_weight, 00170 FLTNB a_multiplicativeCorrections, FLTNB a_additiveCorrections, 00171 FLTNB a_quantificationFactor, oProjectionLine* ap_Line ) 00172 { 00173 // Line weight here is simply 1 00174 *ap_weight = 1.; 00175 // That's all 00176 return 0; 00177 } 00178 00179 // ===================================================================== 00180 // --------------------------------------------------------------------- 00181 // --------------------------------------------------------------------- 00182 // ===================================================================== 00183 00184 int iOptimizerLandweber::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 // Simply subtract the model from the data 00189 *ap_backwardValues = (a_data - a_forwardModel); 00190 // That's all 00191 return 0; 00192 } 00193 00194 // ===================================================================== 00195 // --------------------------------------------------------------------- 00196 // --------------------------------------------------------------------- 00197 // ===================================================================== 00198 00199 int iOptimizerLandweber::ImageSpaceSpecificOperations( FLTNB a_currentImageValue, FLTNB* ap_newImageValue, 00200 FLTNB a_sensitivity, FLTNB* ap_correctionValues ) 00201 { 00202 // Compute image update factor 00203 FLTNB image_update_factor = *ap_correctionValues * m_relaxationFactor / a_sensitivity; 00204 // Update image 00205 *ap_newImageValue = a_currentImageValue + image_update_factor; 00206 // End 00207 return 0; 00208 } 00209 00210 // ===================================================================== 00211 // --------------------------------------------------------------------- 00212 // --------------------------------------------------------------------- 00213 // =====================================================================