Bug Summary

File:programs/AmplitudeAnalysis/fit/fit.cc
Location:line 182, column 18
Description:The right operand of '<' is a garbage value

Annotated Source Code

1#include <iostream>
2#include <fstream>
3#include <sstream>
4#include <complex>
5#include <string>
6#include <vector>
7#include <utility>
8#include <map>
9
10#include "TSystem.h"
11
12#include "AMPTOOLS_DATAIO/ROOTDataReader.h"
13#include "AMPTOOLS_DATAIO/ROOTDataReaderBootstrap.h"
14#include "AMPTOOLS_DATAIO/ROOTDataReaderWithTCut.h"
15#include "AMPTOOLS_DATAIO/ROOTDataReaderTEM.h"
16#include "AMPTOOLS_AMPS/TwoPSAngles.h"
17#include "AMPTOOLS_AMPS/TwoPSHelicity.h"
18#include "AMPTOOLS_AMPS/TwoPiAngles.h"
19#include "AMPTOOLS_AMPS/TwoPiAngles_amp.h"
20#include "AMPTOOLS_AMPS/TwoPiWt_primakoff.h"
21#include "AMPTOOLS_AMPS/TwoPiWt_sigma.h"
22#include "AMPTOOLS_AMPS/TwoPiW_brokenetas.h"
23#include "AMPTOOLS_AMPS/TwoPitdist.h"
24#include "AMPTOOLS_AMPS/TwoPiNC_tdist.h"
25#include "AMPTOOLS_AMPS/TwoPiEtas_tdist.h"
26#include "AMPTOOLS_AMPS/TwoPiAngles_primakoff.h"
27#include "AMPTOOLS_AMPS/ThreePiAngles.h"
28#include "AMPTOOLS_AMPS/ThreePiAnglesSchilling.h"
29#include "AMPTOOLS_AMPS/TwoPiAnglesRadiative.h"
30#include "AMPTOOLS_AMPS/Zlm.h"
31#include "AMPTOOLS_AMPS/BreitWigner.h"
32#include "AMPTOOLS_AMPS/BreitWigner3body.h"
33#include "AMPTOOLS_AMPS/b1piAngAmp.h"
34#include "AMPTOOLS_AMPS/omegapiAngAmp.h"
35#include "AMPTOOLS_AMPS/Uniform.h"
36#include "AMPTOOLS_AMPS/polCoef.h"
37#include "AMPTOOLS_AMPS/dblRegge.h"
38#include "AMPTOOLS_AMPS/dblReggeMod.h"
39#include "AMPTOOLS_AMPS/omegapi_amplitude.h"
40#include "AMPTOOLS_AMPS/Vec_ps_refl.h"
41#include "AMPTOOLS_AMPS/Piecewise.h"
42
43#include "MinuitInterface/MinuitMinimizationManager.h"
44#include "IUAmpTools/AmpToolsInterface.h"
45#include "IUAmpTools/FitResults.h"
46#include "IUAmpTools/ConfigFileParser.h"
47#include "IUAmpTools/ConfigurationInfo.h"
48
49using std::complex;
50using namespace std;
51
52double runSingleFit(ConfigurationInfo* cfgInfo, bool useMinos, int maxIter, string seedfile) {
53 AmpToolsInterface ati( cfgInfo );
54
55 cout << "LIKELIHOOD BEFORE MINIMIZATION: " << ati.likelihood() << endl;
56
57 MinuitMinimizationManager* fitManager = ati.minuitMinimizationManager();
58 fitManager->setMaxIterations(maxIter);
59
60 if( useMinos ){
61
62 fitManager->minosMinimization();
63 }
64 else{
65
66 fitManager->migradMinimization();
67 }
68
69 bool fitFailed =
70 ( fitManager->status() != 0 && fitManager->eMatrixStatus() != 3 );
71
72 if( fitFailed ){
73 cout << "ERROR: fit failed use results with caution..." << endl;
74 return 1e6;
75 }
76
77 cout << "LIKELIHOOD AFTER MINIMIZATION: " << ati.likelihood() << endl;
78
79 ati.finalizeFit();
80
81 if( seedfile.size() != 0 && !fitFailed ){
82 ati.fitResults()->writeSeed( seedfile );
83 }
84
85 return ati.likelihood();
86}
87
88void runRndFits(ConfigurationInfo* cfgInfo, bool useMinos, int maxIter, string seedfile, int numRnd, double maxFraction) {
89 AmpToolsInterface ati( cfgInfo );
90 string fitName = cfgInfo->fitName();
91
92 cout << "LIKELIHOOD BEFORE MINIMIZATION: " << ati.likelihood() << endl;
93
94 MinuitMinimizationManager* fitManager = ati.minuitMinimizationManager();
95 fitManager->setMaxIterations(maxIter);
96
97 vector< vector<string> > parRangeKeywords = cfgInfo->userKeywordArguments("parRange");
98
99 // keep track of best fit (mininum log-likelihood)
100 double minLL = 0;
101 int minFitTag = -1;
102
103 for(int i=0; i<numRnd; i++) {
104 cout << endl << "###############################" << endl;
105 cout << "FIT " << i << " OF " << numRnd << endl;
106 cout << endl << "###############################" << endl;
107
108 // randomize parameters
109 ati.randomizeProductionPars(maxFraction);
110 for(size_t ipar=0; ipar<parRangeKeywords.size(); ipar++) {
111 ati.randomizeParameter(parRangeKeywords[ipar][0], atof(parRangeKeywords[ipar][1].c_str()), atof(parRangeKeywords[ipar][2].c_str()));
112 }
113
114 if(useMinos)
115 fitManager->minosMinimization();
116 else
117 fitManager->migradMinimization();
118
119 bool fitFailed = (fitManager->status() != 0 && fitManager->eMatrixStatus() != 3);
120
121 if( fitFailed )
122 cout << "ERROR: fit failed use results with caution..." << endl;
123
124 cout << "LIKELIHOOD AFTER MINIMIZATION: " << ati.likelihood() << endl;
125
126 ati.finalizeFit(to_string(i));
127
128 if( seedfile.size() != 0 && !fitFailed ){
129 string seedfile_rand = seedfile + Form("_%d.txt", i);
130 ati.fitResults()->writeSeed( seedfile_rand );
131 }
132
133 // update best fit
134 if( !fitFailed && ati.likelihood() < minLL ) {
135 minLL = ati.likelihood();
136 minFitTag = i;
137 }
138 }
139
140 // print best fit results
141 if(minFitTag < 0) cout << "ALL FITS FAILED!" << endl;
142 else {
143 cout << "MINIMUM LIKELIHOOD FROM " << minFitTag << " of " << numRnd << " RANDOM PRODUCTION PARS = " << minLL << endl;
144 gSystem->Exec(Form("cp %s_%d.fit %s.fit", fitName.data(), minFitTag, fitName.data()));
145 if( seedfile.size() != 0 )
146 gSystem->Exec(Form("cp %s_%d.txt %s.txt", seedfile.data(), minFitTag, seedfile.data()));
147 }
148}
149
150void runParScan(ConfigurationInfo* cfgInfo, bool useMinos, int maxIter, string seedfile, string parScan) {
151 double minVal, maxVal, stepSize;
152 int steps;
7
'steps' declared without an initial value
153
154 vector< vector<string> > parScanKeywords = cfgInfo->userKeywordArguments("parScan");
155
156 if(parScanKeywords.size()==0) {
8
Taking false branch
157 cout << "No parScan keyword found in configuration file. Set up at least one parameter for scanning! Aborting." << endl;
158 return;
159 } else {
160 for(size_t ipar=0; ipar<parScanKeywords.size(); ipar++) {
9
Loop condition is true. Entering loop body
11
Loop condition is false. Execution continues on line 172
161 if(parScanKeywords[ipar][0]==parScan) {
10
Taking false branch
162 minVal = atof(parScanKeywords[ipar][1].c_str());
163 maxVal = atof(parScanKeywords[ipar][2].c_str());
164 stepSize = atof(parScanKeywords[ipar][3].c_str());
165 steps = trunc((maxVal-minVal)/stepSize)+1;
166 break;
167 } else
168 cout << "Skipping configuration to scan " << parScanKeywords[ipar][0] << "since scanning of " << parScan << " was requested..." << endl;
169 }
170 }
171
172 AmpToolsInterface ati( cfgInfo );
173
174 string fitName = cfgInfo->fitName();
175 cout << "LIKELIHOOD BEFORE MINIMIZATION: " << ati.likelihood() << endl;
176
177 ParameterManager* parMgr = ati.parameterManager();
178 MinuitMinimizationManager* fitManager = ati.minuitMinimizationManager();
179 fitManager->setMaxIterations(maxIter);
180
181
182 for(int i=0; i<steps; i++) {
12
The right operand of '<' is a garbage value
183 cout << endl << "###############################" << endl;
184 cout << "FIT " << i << " OF " << steps << endl;
185 cout << endl << "###############################" << endl;
186
187 // set parameter to be scanned
188 vector<ParameterInfo*> parInfoVec = cfgInfo->parameterList();
189
190 auto parItr = parInfoVec.begin();
191 for( ; parItr != parInfoVec.end(); ++parItr ) {
192 if( (**parItr).parName() == parScan ) break;
193 }
194
195 if( parItr == parInfoVec.end() ){
196 cout << "ERROR: request to scan nonexistent parameter: " << parScan << endl;
197 return;
198 }
199
200 // set and fix parameter for scan
201 double value = minVal + i*stepSize;
202 parMgr->setAmpParameter( parScan, value );
203
204 cfgInfo->setFitName(fitName + "_scan");
205
206 if(useMinos)
207 fitManager->minosMinimization();
208 else
209 fitManager->migradMinimization();
210
211 bool fitFailed = (fitManager->status() != 0 && fitManager->eMatrixStatus() != 3);
212
213 if( fitFailed )
214 cout << "ERROR: fit failed use results with caution..." << endl;
215
216 cout << "LIKELIHOOD AFTER MINIMIZATION: " << ati.likelihood() << endl;
217
218 ati.finalizeFit(to_string(i));
219
220 if( seedfile.size() != 0 && !fitFailed ){
221 string seedfile_scan = seedfile + Form("_scan_%d.txt", i);
222 ati.fitResults()->writeSeed( seedfile_scan );
223 }
224 }
225}
226
227int main( int argc, char* argv[] ){
228
229 // set default parameters
230
231 bool useMinos = false;
232
233 string configfile;
234 string seedfile;
235 string scanPar;
236 int numRnd = 0;
237 int maxIter = 10000;
238
239 // parse command line
240
241 for (int i = 1; i < argc; i++){
1
Assuming 'i' is >= 'argc'
2
Loop condition is false. Execution continues on line 272
242
243 string arg(argv[i]);
244
245 if (arg == "-c"){
246 if ((i+1 == argc) || (argv[i+1][0] == '-')) arg = "-h";
247 else configfile = argv[++i]; }
248 if (arg == "-s"){
249 if ((i+1 == argc) || (argv[i+1][0] == '-')) arg = "-h";
250 else seedfile = argv[++i]; }
251 if (arg == "-r"){
252 if ((i+1 == argc) || (argv[i+1][0] == '-')) arg = "-h";
253 else numRnd = atoi(argv[++i]); }
254 if (arg == "-m"){
255 if ((i+1 == argc) || (argv[i+1][0] == '-')) arg = "-h";
256 else maxIter = atoi(argv[++i]); }
257 if (arg == "-n") useMinos = true;
258 if (arg == "-p"){
259 if ((i+1 == argc) || (argv[i+1][0] == '-')) arg = "-h";
260 else scanPar = argv[++i]; }
261 if (arg == "-h"){
262 cout << endl << " Usage for: " << argv[0] << endl << endl;
263 cout << " -n \t\t\t\t\t use MINOS instead of MIGRAD" << endl;
264 cout << " -c <file>\t\t\t\t config file" << endl;
265 cout << " -s <output file>\t\t\t for seeding next fit based on this fit (optional)" << endl;
266 cout << " -r <int>\t\t\t Perform <int> fits each seeded with random parameters" << endl;
267 cout << " -p <parameter> \t\t\t\t Perform a scan of given parameter. Stepsize, min, max are to be set in cfg file" << endl;
268 cout << " -m <int>\t\t\t Maximum number of fit iterations" << endl;
269 exit(1);}
270 }
271
272 if (configfile.size() == 0){
3
Taking false branch
273 cout << "No config file specified" << endl;
274 exit(1);
275 }
276
277 ConfigFileParser parser(configfile);
278 ConfigurationInfo* cfgInfo = parser.getConfigurationInfo();
279 cfgInfo->display();
280
281 AmpToolsInterface::registerAmplitude( BreitWigner() );
282 AmpToolsInterface::registerAmplitude( BreitWigner3body() );
283 AmpToolsInterface::registerAmplitude( TwoPSAngles() );
284 AmpToolsInterface::registerAmplitude( TwoPSHelicity() );
285 AmpToolsInterface::registerAmplitude( TwoPiAngles() );
286 AmpToolsInterface::registerAmplitude( TwoPiAngles_amp() );
287 AmpToolsInterface::registerAmplitude( TwoPiAngles_primakoff() );
288 AmpToolsInterface::registerAmplitude( TwoPiWt_primakoff() );
289 AmpToolsInterface::registerAmplitude( TwoPiWt_sigma() );
290 AmpToolsInterface::registerAmplitude( TwoPitdist() );
291 AmpToolsInterface::registerAmplitude( ThreePiAngles() );
292 AmpToolsInterface::registerAmplitude( ThreePiAnglesSchilling() );
293 AmpToolsInterface::registerAmplitude( TwoPiAnglesRadiative() );
294 AmpToolsInterface::registerAmplitude( Zlm() );
295 AmpToolsInterface::registerAmplitude( b1piAngAmp() );
296 AmpToolsInterface::registerAmplitude( omegapiAngAmp() );
297 AmpToolsInterface::registerAmplitude( polCoef() );
298 AmpToolsInterface::registerAmplitude( Uniform() );
299 AmpToolsInterface::registerAmplitude( dblRegge() );
300 AmpToolsInterface::registerAmplitude( omegapi_amplitude() );
301 AmpToolsInterface::registerAmplitude( Vec_ps_refl() );
302 AmpToolsInterface::registerAmplitude( Piecewise() );
303
304 AmpToolsInterface::registerDataReader( ROOTDataReader() );
305 AmpToolsInterface::registerDataReader( ROOTDataReaderBootstrap() );
306 AmpToolsInterface::registerDataReader( ROOTDataReaderWithTCut() );
307 AmpToolsInterface::registerDataReader( ROOTDataReaderTEM() );
308
309 if(numRnd==0){
4
Taking true branch
310 if(scanPar=="")
5
Taking false branch
311 runSingleFit(cfgInfo, useMinos, maxIter, seedfile);
312 else
313 runParScan(cfgInfo, useMinos, maxIter, seedfile, scanPar);
6
Calling 'runParScan'
314 } else {
315 runRndFits(cfgInfo, useMinos, maxIter, seedfile, numRnd, 0.5);
316 }
317
318 return 0;
319}
320
321