1 | /* ======================================================================== *\
|
---|
2 | !
|
---|
3 | ! *
|
---|
4 | ! * This file is part of MARS, the MAGIC Analysis and Reconstruction
|
---|
5 | ! * Software. It is distributed to you in the hope that it can be a useful
|
---|
6 | ! * and timesaving tool in analysing Data of imaging Cerenkov telescopes.
|
---|
7 | ! * It is distributed WITHOUT ANY WARRANTY.
|
---|
8 | ! *
|
---|
9 | ! * Permission to use, copy, modify and distribute this software and its
|
---|
10 | ! * documentation for any purpose is hereby granted without fee,
|
---|
11 | ! * provided that the above copyright notice appear in all copies and
|
---|
12 | ! * that both that copyright notice and this permission notice appear
|
---|
13 | ! * in supporting documentation. It is provided "as is" without express
|
---|
14 | ! * or implied warranty.
|
---|
15 | ! *
|
---|
16 | !
|
---|
17 | !
|
---|
18 | ! Author(s): David Paneque 12/2003 <mailto:dpaneque@mppmu.mpg.de>
|
---|
19 | !
|
---|
20 | !
|
---|
21 | ! Copyright: MAGIC Software Development, 2000-2003
|
---|
22 | !
|
---|
23 | !
|
---|
24 | \* ======================================================================== */
|
---|
25 |
|
---|
26 | /////////////////////////////////////////////////////////////////////////////
|
---|
27 | // //
|
---|
28 | // MFindSupercutsONOFFThetaLoop //
|
---|
29 | // //
|
---|
30 | // Class for optimizing the parameters of the supercuts //
|
---|
31 | // Using ON and OFF data runing over data which are subdivided in
|
---|
32 | // several Theta ranges. //
|
---|
33 | // //
|
---|
34 | // //
|
---|
35 | /////////////////////////////////////////////////////////////////////////////
|
---|
36 |
|
---|
37 |
|
---|
38 | #include "MFindSupercutsONOFFThetaLoop.h"
|
---|
39 |
|
---|
40 |
|
---|
41 | #include <math.h> // fabs
|
---|
42 |
|
---|
43 | #include <TFile.h>
|
---|
44 | #include <TArrayD.h>
|
---|
45 | #include <TMinuit.h>
|
---|
46 | #include <TCanvas.h>
|
---|
47 | #include <TStopwatch.h>
|
---|
48 | #include <TVirtualFitter.h>
|
---|
49 |
|
---|
50 | #include <iostream.h>
|
---|
51 | #include <fstream.h>
|
---|
52 | #include <stdlib.h>
|
---|
53 | #include <stdio.h>
|
---|
54 | #include <string.h>
|
---|
55 |
|
---|
56 | #include "MFindSupercutsONOFF.h"
|
---|
57 | #include "MBinning.h"
|
---|
58 | #include "MContinue.h"
|
---|
59 | #include "MSupercuts.h"
|
---|
60 | #include "MSupercutsCalcONOFF.h"
|
---|
61 | #include "MDataElement.h"
|
---|
62 | #include "MDataMember.h"
|
---|
63 |
|
---|
64 |
|
---|
65 | #include <TPostScript.h>
|
---|
66 |
|
---|
67 | #include "MEvtLoop.h"
|
---|
68 | #include "MFCT1SelFinal.h"
|
---|
69 | #include "MF.h"
|
---|
70 | #include "MFEventSelector.h"
|
---|
71 | #include "MFEventSelector2.h"
|
---|
72 | #include "MFillH.h"
|
---|
73 | //#include "MGeomCamCT1Daniel.h"
|
---|
74 | //#include "MGeomCamCT1.h"
|
---|
75 | #include "MGeomCamMagic.h" // only for magic
|
---|
76 | #include "MFRandomSplit.h"
|
---|
77 | #include "MH3.h"
|
---|
78 | #include "MHCT1Supercuts.h"
|
---|
79 | #include "MHFindSignificance.h" // To be removed at some point...
|
---|
80 | #include "MHFindSignificanceONOFF.h"
|
---|
81 | #include "MHMatrix.h"
|
---|
82 | #include "MHOnSubtraction.h"
|
---|
83 | #include "MDataValue.h"
|
---|
84 | // #include "MDataString.h"
|
---|
85 |
|
---|
86 | #include "MLog.h"
|
---|
87 | #include "MLogManip.h"
|
---|
88 | #include "MMatrixLoop.h"
|
---|
89 | #include "MMinuitInterface.h"
|
---|
90 | #include "MParList.h"
|
---|
91 | #include "MProgressBar.h"
|
---|
92 | #include "MReadMarsFile.h"
|
---|
93 | #include "MReadTree.h"
|
---|
94 | #include "MTaskList.h"
|
---|
95 |
|
---|
96 |
|
---|
97 | ClassImp(MFindSupercutsONOFF);
|
---|
98 |
|
---|
99 | using namespace std;
|
---|
100 |
|
---|
101 | // --------------------------------------------------------------------------
|
---|
102 | //
|
---|
103 | // Default constructor.
|
---|
104 | //
|
---|
105 | MFindSupercutsONOFFThetaLoop::MFindSupercutsONOFFThetaLoop(const char *name, const char *title)
|
---|
106 | {
|
---|
107 | fName = name ? name : "MFindSupercutsONOFF";
|
---|
108 | fTitle = title ? title : "Optimizer of the supercuts";
|
---|
109 |
|
---|
110 | //---------------------------
|
---|
111 |
|
---|
112 |
|
---|
113 | // Cuts (low and up) in variable ThetaOrig.fVal
|
---|
114 | // The default is not cut, i.e. all values (0-1) are taken
|
---|
115 | // fThetaOrig.fVal is measured in Radians; thus 1 = 57 degrees.
|
---|
116 |
|
---|
117 | fThetaMin = 0.0;
|
---|
118 | fThetaMax = 1.0;
|
---|
119 |
|
---|
120 | fAlphaSig = 20; // By default, signal is expected in alpha<20 for CT1
|
---|
121 |
|
---|
122 | // By default, bkg region is set to alpha range 30-90 for CT1
|
---|
123 | fAlphaBkgMin = 30;
|
---|
124 | fAlphaBkgMax = 90;
|
---|
125 |
|
---|
126 | fActualCosThetaBinCenter = 0;
|
---|
127 |
|
---|
128 | fNormFactorTrainHist = NULL;
|
---|
129 | //fNormFactorTrainHist -> SetDirectory(0);
|
---|
130 | fNormFactorTestHist = NULL;
|
---|
131 | //fNormFactorTestHist -> SetDirectory(0);
|
---|
132 | fSigmaLiMaTrainHist = NULL;
|
---|
133 | //fSigmaLiMaTrainHist -> SetDirectory(0);
|
---|
134 | fSigmaLiMaTestHist = NULL;
|
---|
135 | //fSigmaLiMaTestHist -> SetDirectory(0);
|
---|
136 | fNexTrainHist = NULL;
|
---|
137 | //fSigmaLiMaTestHist -> SetDirectory(0);
|
---|
138 | fNexTestHist = NULL;
|
---|
139 | //fNexTestHist -> SetDirectory(0);
|
---|
140 |
|
---|
141 |
|
---|
142 | fSuccessfulThetaBinsHist = NULL;
|
---|
143 |
|
---|
144 | /*
|
---|
145 | fNexVSAlphaSigTrain = NULL;
|
---|
146 | fNexVSAlphaSigTest = NULL;
|
---|
147 |
|
---|
148 | fSignificanceVSAlphaSigTrain = NULL;
|
---|
149 | fSignificanceVSAlphaSigTest = NULL;
|
---|
150 | */
|
---|
151 |
|
---|
152 |
|
---|
153 | fNEvtsInTrainMatrixONHist = NULL;
|
---|
154 | fNEvtsInTestMatrixONHist = NULL;
|
---|
155 | fNEvtsInTrainMatrixOFFHist = NULL;
|
---|
156 | fNEvtsInTestMatrixOFFHist = NULL;
|
---|
157 |
|
---|
158 |
|
---|
159 | fTrainMatrixONFilenameVector = NULL;
|
---|
160 | fTestMatrixONFilenameVector= NULL;
|
---|
161 | fTrainMatrixOFFFilenameVector= NULL;
|
---|
162 | fTestMatrixOFFFilenameVector= NULL;
|
---|
163 |
|
---|
164 |
|
---|
165 | fOptSCParamFilenameVector = NULL;
|
---|
166 |
|
---|
167 | fThetaRangeStringVector = NULL;
|
---|
168 |
|
---|
169 | fPsFilename = NULL;
|
---|
170 |
|
---|
171 |
|
---|
172 |
|
---|
173 | fTuneNormFactor = kTRUE; // Norm factors will be corrected using the total amount of OFF events before cuts and the estimated excess events
|
---|
174 |
|
---|
175 | // fNormFactorTrain = fNormFactorTrain - Ngammas/EventsInTrainMatrixOFF
|
---|
176 |
|
---|
177 | fGammaEfficiency = 0.5; // Fraction of gammas that remain after cuts
|
---|
178 | // Quantity that will have to be determined with MC, yet for the
|
---|
179 | // time being I set it to 0.5
|
---|
180 |
|
---|
181 | // boolean variables are set to kFALSE as default
|
---|
182 |
|
---|
183 | fReadMatricesFromFile = kFALSE;
|
---|
184 | fOptimizeParameters = kFALSE;
|
---|
185 | fTestParameters = kFALSE;
|
---|
186 |
|
---|
187 | // Fraction of Train/Test ON/OFF samples is set to 0.5 as default
|
---|
188 |
|
---|
189 | fWhichFractionTrain = 0.5; // number <= 1; specifying fraction of ON Train events
|
---|
190 | fWhichFractionTest = 0.5; // number <= 1; specifying fraction of ON Test events
|
---|
191 |
|
---|
192 | fWhichFractionTrainOFF = 0.5; // number <= 1; specifying fraction of OFF Train events
|
---|
193 | fWhichFractionTestOFF = 0.5; // number <= 1; specifying fraction of OFF Test events
|
---|
194 |
|
---|
195 |
|
---|
196 |
|
---|
197 | // use quantities computed from the fits
|
---|
198 | // The variable allows the user to NOT use these quantities when there is not
|
---|
199 | // enough statistics and fit not always is possible
|
---|
200 | fUseFittedQuantities = kTRUE;
|
---|
201 |
|
---|
202 |
|
---|
203 |
|
---|
204 | // Boolean variable that controls wether the optimization of the
|
---|
205 | // parameters (MMinuitInterface::CallMinuit(..) in function FindParams(..))
|
---|
206 | // takes place or not. kTRUE will skip such optimization.
|
---|
207 | // This variable is useful to test the optmized parameters (previously found
|
---|
208 | // and stored in root file) on the TRAIN sample.
|
---|
209 |
|
---|
210 | fSkipOptimization = kFALSE;
|
---|
211 |
|
---|
212 | // Boolean variable that allows the user to write the initial parameters
|
---|
213 | // into the root file that will be used to store the optimum cuts.
|
---|
214 | // If fUseInitialSCParams = kTRUE , parameters are written.
|
---|
215 | // In this way, the initial SC parameters can be applied on the data (train/test)
|
---|
216 |
|
---|
217 | // The initial parameters are ONLY written to the root file if
|
---|
218 | // there is NO SC params optimization, i.e., if variable
|
---|
219 | // fSkipOptimization = kTRUE;
|
---|
220 |
|
---|
221 | // The default value is obviously kFALSE.
|
---|
222 |
|
---|
223 | fUseInitialSCParams = kTRUE;
|
---|
224 |
|
---|
225 |
|
---|
226 |
|
---|
227 | // fInitSCPar and fInitSCParSteps initialized to zero
|
---|
228 |
|
---|
229 | fInitSCPar.Set(0);
|
---|
230 |
|
---|
231 | fInitSCParSteps.Set(0);
|
---|
232 |
|
---|
233 | }
|
---|
234 |
|
---|
235 | // --------------------------------------------------------------------------
|
---|
236 | //
|
---|
237 | // Default destructor.
|
---|
238 | //
|
---|
239 | MFindSupercutsONOFFThetaLoop::~MFindSupercutsONOFFThetaLoop()
|
---|
240 | {
|
---|
241 |
|
---|
242 | *fLog << "Destructor of MFindSupercutsONOFFThetaLoop is called"
|
---|
243 | << endl;
|
---|
244 |
|
---|
245 | fPsFilename = NULL;
|
---|
246 |
|
---|
247 | delete fNormFactorTrainHist;
|
---|
248 | delete fNormFactorTestHist;
|
---|
249 | delete fSigmaLiMaTrainHist;
|
---|
250 | delete fSigmaLiMaTestHist;
|
---|
251 | delete fNexTrainHist;
|
---|
252 | delete fNexTestHist;
|
---|
253 |
|
---|
254 |
|
---|
255 |
|
---|
256 | delete fSuccessfulThetaBinsHist;
|
---|
257 |
|
---|
258 |
|
---|
259 | delete fNEvtsInTrainMatrixONHist;
|
---|
260 | delete fNEvtsInTestMatrixONHist;
|
---|
261 | delete fNEvtsInTrainMatrixOFFHist;
|
---|
262 | delete fNEvtsInTestMatrixOFFHist;
|
---|
263 |
|
---|
264 | delete [] fTrainMatrixONFilenameVector;
|
---|
265 | delete [] fTestMatrixONFilenameVector;
|
---|
266 | delete [] fTrainMatrixOFFFilenameVector;
|
---|
267 | delete [] fTestMatrixOFFFilenameVector;
|
---|
268 |
|
---|
269 | delete [] fOptSCParamFilenameVector;
|
---|
270 |
|
---|
271 | delete [] fThetaRangeStringVector;
|
---|
272 |
|
---|
273 |
|
---|
274 | *fLog << "Destructor of MFindSupercutsONOFFThetaLoop finished "
|
---|
275 | << "successfully."
|
---|
276 | << endl;
|
---|
277 |
|
---|
278 | }
|
---|
279 |
|
---|
280 |
|
---|
281 | // --------------------------------------------------------------------------
|
---|
282 | //
|
---|
283 | // Function that sets the name of the PostScript file where alpha distributions
|
---|
284 | // for the different Theta bins will be stored.
|
---|
285 | // It also initializes
|
---|
286 |
|
---|
287 |
|
---|
288 |
|
---|
289 | void MFindSupercutsONOFFThetaLoop::SetPostScriptFile (TPostScript* PsFile)
|
---|
290 | {
|
---|
291 | fPsFilename = PsFile;
|
---|
292 |
|
---|
293 | *fLog << "MFindSupercutsONOFF : Results (alpha distributions with excess and significances) will be stored in PostScript file "
|
---|
294 | << fPsFilename -> GetName() << endl;
|
---|
295 |
|
---|
296 | }
|
---|
297 |
|
---|
298 |
|
---|
299 | Bool_t MFindSupercutsONOFFThetaLoop::SetGammaEfficiency (Double_t gammaeff)
|
---|
300 | {
|
---|
301 | if (gammaeff <= 0.0 || gammaeff > 1)
|
---|
302 | {
|
---|
303 | *fLog << " ****** ERROR *********" << endl
|
---|
304 | << "MFindSupercutsONOFFThetaLoop :: SetGammaEfficiency; " << endl
|
---|
305 | << "The requested Gamma efficiency " << gammaeff << " can not occur" << endl
|
---|
306 | << "Member data fGammaEfficiency can not be set properly." << endl
|
---|
307 | << "Default value for fGammaEfficiency, which is "
|
---|
308 | << fGammaEfficiency << " will remain." << endl;
|
---|
309 |
|
---|
310 |
|
---|
311 | return kFALSE;
|
---|
312 | }
|
---|
313 |
|
---|
314 | fGammaEfficiency = gammaeff;
|
---|
315 |
|
---|
316 | return kTRUE;
|
---|
317 |
|
---|
318 | }
|
---|
319 |
|
---|
320 |
|
---|
321 | Bool_t MFindSupercutsONOFFThetaLoop::ReadSCParamsFromAsciiFile(const char* filename,
|
---|
322 | Int_t Nparams)
|
---|
323 | {
|
---|
324 | const Int_t NColumns = 2;
|
---|
325 |
|
---|
326 | Double_t tmp_vector[NColumns];
|
---|
327 | Int_t counter = 0;
|
---|
328 | Double_t tmpdouble = 0.0;
|
---|
329 |
|
---|
330 |
|
---|
331 | // File is read and data stored in data_vector
|
---|
332 |
|
---|
333 | ifstream infile (filename, ios::in);
|
---|
334 |
|
---|
335 | if ( !infile )
|
---|
336 | {
|
---|
337 | cout << "Input file could not be opened" << endl;
|
---|
338 | return kFALSE;
|
---|
339 | // exit (1);
|
---|
340 | }
|
---|
341 |
|
---|
342 | counter = 0;
|
---|
343 | while (infile >> tmp_vector [0])
|
---|
344 | {
|
---|
345 | for (Int_t i = 1; i < NColumns; i++)
|
---|
346 | {
|
---|
347 | infile >> tmp_vector[i];
|
---|
348 | }
|
---|
349 | counter++;
|
---|
350 | }
|
---|
351 |
|
---|
352 | infile.close();
|
---|
353 |
|
---|
354 |
|
---|
355 | // Length of TArray vectors is set to the length of the
|
---|
356 | // vectors defined in ascii file.
|
---|
357 | // If such length does not match with the specified
|
---|
358 | // by variable Nparams, the function complains and
|
---|
359 | // aborts.
|
---|
360 |
|
---|
361 | if (counter != Nparams)
|
---|
362 | {
|
---|
363 | *fLog << "MFindSupercutsONOFFThetaLoop::ReadSCParamsFromAsciiFile" << endl
|
---|
364 | << "Length of vectors in ascii file " << filename << " is " << counter
|
---|
365 | << " , which does not match with the expected length " << Nparams << endl
|
---|
366 | << " fInitSCPar and fInitSCParSteps are NOT retrieved from ascii file" << endl;
|
---|
367 |
|
---|
368 | return kFALSE;
|
---|
369 | }
|
---|
370 |
|
---|
371 | fInitSCPar.Set(counter);
|
---|
372 | fInitSCParSteps.Set(counter);
|
---|
373 |
|
---|
374 |
|
---|
375 | // Vectors are filled
|
---|
376 |
|
---|
377 | ifstream infile2 (filename, ios::in);
|
---|
378 |
|
---|
379 | if ( !infile2 )
|
---|
380 | {
|
---|
381 | cout << "Input file could not be opened" << endl;
|
---|
382 | return kFALSE;
|
---|
383 | // exit (1);
|
---|
384 | }
|
---|
385 |
|
---|
386 | for (Int_t i = 0; i < fInitSCPar.GetSize(); i++)
|
---|
387 | {
|
---|
388 |
|
---|
389 | infile2 >> tmpdouble;
|
---|
390 | fInitSCPar[i] = tmpdouble;
|
---|
391 |
|
---|
392 | infile2 >> tmpdouble;
|
---|
393 | fInitSCParSteps[i] = tmpdouble;
|
---|
394 |
|
---|
395 | }
|
---|
396 |
|
---|
397 |
|
---|
398 | // Print the vectors read
|
---|
399 |
|
---|
400 | *fLog << "***** INITIAL SC PARAMETERS AND STEPS TAKEN FROM ASCII FILE ********* " << endl;
|
---|
401 | for (Int_t i = 0; i < fInitSCPar.GetSize(); i++)
|
---|
402 | {
|
---|
403 | cout << fInitSCPar[i] << setw(20) << fInitSCParSteps[i] << endl;
|
---|
404 | }
|
---|
405 |
|
---|
406 | //
|
---|
407 |
|
---|
408 | // File close
|
---|
409 | infile2.close();
|
---|
410 |
|
---|
411 | return kTRUE;
|
---|
412 | }
|
---|
413 |
|
---|
414 |
|
---|
415 | Bool_t MFindSupercutsONOFFThetaLoop::SetInitSCPar (TArrayD &d)
|
---|
416 | {
|
---|
417 | *fLog << "MFindSupercutsONOFFThetaLoop; Initial SC parameters are "
|
---|
418 | << "set using function SetInitSCPar (TArrayD &d)" << endl;
|
---|
419 |
|
---|
420 | fInitSCPar = d;
|
---|
421 |
|
---|
422 | return kTRUE;
|
---|
423 | }
|
---|
424 |
|
---|
425 |
|
---|
426 | Bool_t MFindSupercutsONOFFThetaLoop::SetInitSCParSteps (TArrayD &d)
|
---|
427 | {
|
---|
428 |
|
---|
429 | *fLog << "MFindSupercutsONOFFThetaLoop; Initial SC parameters STEPS are "
|
---|
430 | << "set using function SetInitSCParSteps (TArrayD &d)" << endl;
|
---|
431 |
|
---|
432 | fInitSCParSteps = d;
|
---|
433 |
|
---|
434 | return kTRUE;
|
---|
435 | }
|
---|
436 |
|
---|
437 |
|
---|
438 |
|
---|
439 | Bool_t MFindSupercutsONOFFThetaLoop::SetAlphaSig(Double_t alphasig)
|
---|
440 | {
|
---|
441 | // check that alpha is within the limits 0-90
|
---|
442 | if (alphasig <= 0 || alphasig > 90)
|
---|
443 | {
|
---|
444 | *fLog << "MFindSupercutsONOFFThetaLoop ::SetAlphaSig; "
|
---|
445 | << "value " << alphasig << " is not within the the "
|
---|
446 | << "logical limits of alpha; 0-90" << endl;
|
---|
447 | return kFALSE;
|
---|
448 | }
|
---|
449 |
|
---|
450 |
|
---|
451 | fAlphaSig = alphasig;
|
---|
452 |
|
---|
453 | return kTRUE;
|
---|
454 | }
|
---|
455 |
|
---|
456 | Bool_t MFindSupercutsONOFFThetaLoop::SetAlphaBkgMin(Double_t alpha)
|
---|
457 | {
|
---|
458 | // check that alpha is within the limits 0-90
|
---|
459 | if (alpha <= 0 || alpha >= 90)
|
---|
460 | {
|
---|
461 | *fLog << "MFindSupercutsONOFFThetaLoop::SetAlphaBkgMin; "
|
---|
462 | << "value " << alpha << " is not within the the "
|
---|
463 | << "logical limits of alpha; 0-90" << endl;
|
---|
464 | return kFALSE;
|
---|
465 | }
|
---|
466 |
|
---|
467 |
|
---|
468 | fAlphaBkgMin = alpha;
|
---|
469 |
|
---|
470 | return kTRUE;
|
---|
471 | }
|
---|
472 |
|
---|
473 |
|
---|
474 | Bool_t MFindSupercutsONOFFThetaLoop::SetAlphaBkgMax(Double_t alpha)
|
---|
475 | {
|
---|
476 | // check that alpha is within the limits 0-90
|
---|
477 | if (alpha <= 0 || alpha > 90.001)
|
---|
478 | {
|
---|
479 | *fLog << "MFindSupercutsONOFFThetaLoop ::SetAlphaBkgMax; "
|
---|
480 | << "value " << alpha << " is not within the the "
|
---|
481 | << "logical limits of alpha; 0-90" << endl;
|
---|
482 | return kFALSE;
|
---|
483 | }
|
---|
484 |
|
---|
485 |
|
---|
486 | fAlphaBkgMax = alpha;
|
---|
487 |
|
---|
488 | return kTRUE;
|
---|
489 | }
|
---|
490 |
|
---|
491 |
|
---|
492 | // Function that checks that the values of the member data
|
---|
493 | // fAlphaSig, fAlphaBkgMin and fAlphaBkgMax make sense
|
---|
494 | // (ie, fAlphaSig < fAlphaBkgMin < fAlphaBkgMax)
|
---|
495 |
|
---|
496 | Bool_t MFindSupercutsONOFFThetaLoop::CheckAlphaSigBkg()
|
---|
497 | {
|
---|
498 |
|
---|
499 | if (fAlphaSig > fAlphaBkgMin)
|
---|
500 | {
|
---|
501 | *fLog << "MFindSupercutsONOFFThetaLoop ::CheckAlphaSigBkg(); " << endl
|
---|
502 | << "fAlphaSig > fAlphaBkgMin, which should not occur..." << endl
|
---|
503 | << "fAlphaSig = " << fAlphaSig << ", fAlphaBkgMin = " << fAlphaBkgMin
|
---|
504 | << endl;
|
---|
505 |
|
---|
506 | return kFALSE;
|
---|
507 | }
|
---|
508 |
|
---|
509 | if (fAlphaBkgMax < fAlphaBkgMin)
|
---|
510 | {
|
---|
511 | *fLog << "MFindSupercutsONOFFThetaLoop::CheckAlphaSigBkg(); " << endl
|
---|
512 | << "fAlphaBkgMin > fAlphaBkgMax, which should not occur..." << endl
|
---|
513 | << "fAlphaBkgMin = " << fAlphaBkgMin << ", fAlphaBkgMax = " << fAlphaBkgMax
|
---|
514 | << endl;
|
---|
515 |
|
---|
516 | return kFALSE;
|
---|
517 | }
|
---|
518 |
|
---|
519 | return kTRUE;
|
---|
520 |
|
---|
521 | }
|
---|
522 |
|
---|
523 |
|
---|
524 |
|
---|
525 |
|
---|
526 | // Function that sets the value of fCosThetaRangeVector and
|
---|
527 | // fCosThetaBinCenterVector
|
---|
528 |
|
---|
529 | Bool_t MFindSupercutsONOFFThetaLoop::SetCosThetaRangeVector(
|
---|
530 | const TArrayD &d)
|
---|
531 | {
|
---|
532 | if (d.GetSize() < 2)
|
---|
533 | {
|
---|
534 | *fLog << err
|
---|
535 | << "MFindSupercutsONOFFThetaLoop::SetCosThetaRangeVector: "
|
---|
536 | << endl
|
---|
537 | << "Size of d is smaller than 2... "
|
---|
538 | << "fCosThetaRangeVector can not be defined properly"
|
---|
539 | << endl;
|
---|
540 | return kFALSE;
|
---|
541 | }
|
---|
542 |
|
---|
543 | fCosThetaRangeVector.Set(d.GetSize());
|
---|
544 | fCosThetaRangeVector = d;
|
---|
545 |
|
---|
546 |
|
---|
547 | // Give a "soft warning" if Vector has increasing values
|
---|
548 | // of CosThetas; which means decreasing values of Thetas
|
---|
549 | /*
|
---|
550 | if (fCosThetaRangeVector[0] <= fCosThetaRangeVector[1])
|
---|
551 | {
|
---|
552 | *fLog << "MFindSupercutsONOFFThetaLoop::SetCosThetaRangeVector; "
|
---|
553 | << endl
|
---|
554 | << "Values in vector fCosThetaRangeVector are in increasing "
|
---|
555 | << "order, i.e., values in theta will be in decreasing order"
|
---|
556 | << endl
|
---|
557 | << "DO NOT forget it !!" << endl;
|
---|
558 | }
|
---|
559 |
|
---|
560 | */
|
---|
561 |
|
---|
562 | // fCosThetaBinCenterVector is defined
|
---|
563 |
|
---|
564 | Int_t dim = fCosThetaRangeVector.GetSize() - 1;
|
---|
565 | fCosThetaBinCenterVector.Set(dim);
|
---|
566 |
|
---|
567 | for (Int_t i = 0; i < dim; i++)
|
---|
568 | {
|
---|
569 | fCosThetaBinCenterVector[i] =
|
---|
570 | (fCosThetaRangeVector[i+1]+fCosThetaRangeVector[i])/2;
|
---|
571 | }
|
---|
572 |
|
---|
573 | return kTRUE;
|
---|
574 | }
|
---|
575 |
|
---|
576 | // Function that sets the range of Theta (fThetaMin, fThetaMax)
|
---|
577 | // where data will be used in the iteration
|
---|
578 | // thetabin. This theta range is used when geting data from
|
---|
579 | // file and filling matrices. If matrices are already built
|
---|
580 | // and just read from root files, this theta range is useless.
|
---|
581 |
|
---|
582 |
|
---|
583 | // **** IMPORTANT *********
|
---|
584 | // thetabin is an integer value that goes from 1 to
|
---|
585 | // fCosThetaBinCenterVector.GetSize()
|
---|
586 | // ************************
|
---|
587 |
|
---|
588 | Bool_t MFindSupercutsONOFFThetaLoop::SetThetaRange (Int_t thetabin)
|
---|
589 | {
|
---|
590 | // Check that iteration thetabin is smaller than
|
---|
591 | // components of fCosThetaBinCenterVector
|
---|
592 |
|
---|
593 | if (fCosThetaRangeVector.GetSize() < 2 ||
|
---|
594 | fCosThetaBinCenterVector.GetSize() < 1)
|
---|
595 |
|
---|
596 | {
|
---|
597 | *fLog << "MFindSupercutsONOFFThetaLoop::SetThetaRange; " << endl
|
---|
598 | << "Dimension of vectors fCosThetaRangeVector or/and "
|
---|
599 | << "fCosThetaBinCenterVector are smaller than 2 and 1 "
|
---|
600 | << "respectively. Range for theta can not be defined." << endl;
|
---|
601 |
|
---|
602 | return kFALSE;
|
---|
603 | }
|
---|
604 |
|
---|
605 | if (thetabin < 1 ||
|
---|
606 | thetabin > fCosThetaBinCenterVector.GetSize())
|
---|
607 | {
|
---|
608 | *fLog << "MFindSupercutsONOFFThetaLoop::SetThetaRange; " << endl
|
---|
609 | << "thetabin value is " << thetabin
|
---|
610 | << " Is is outside the possible values given by "
|
---|
611 | << "vector fCosThetaRangeVector. "
|
---|
612 | << "Range for theta can not be defined." << endl;
|
---|
613 |
|
---|
614 | return kFALSE;
|
---|
615 | }
|
---|
616 |
|
---|
617 |
|
---|
618 | fThetaMin = fCosThetaRangeVector[thetabin-1];
|
---|
619 | fThetaMax = fCosThetaRangeVector[thetabin];
|
---|
620 |
|
---|
621 |
|
---|
622 | // Now costheta values have to be converted to Radians
|
---|
623 | // Function TMath::ACos(x) gives angle in Radians
|
---|
624 |
|
---|
625 |
|
---|
626 | fThetaMin = TMath::ACos(fThetaMin);
|
---|
627 | fThetaMax = TMath::ACos(fThetaMax);
|
---|
628 |
|
---|
629 | // If fThetaMin > fThetaMax means that fCosThetaRangeVector
|
---|
630 | // is defined with increasing values of CosTheta, which
|
---|
631 | // means decreasing values of Theta
|
---|
632 |
|
---|
633 | // In such case, values of fThetaMin and fThetaMax are swapped
|
---|
634 |
|
---|
635 | if (fThetaMin > fThetaMax)
|
---|
636 | {
|
---|
637 | Double_t tmp = fThetaMin;
|
---|
638 | fThetaMin = fThetaMax;
|
---|
639 | fThetaMax = tmp;
|
---|
640 | }
|
---|
641 |
|
---|
642 | // Info for user
|
---|
643 | *fLog << "-------------------------------------------------"<< endl
|
---|
644 | << "MFindSupercutsONOFFThetaLoop::SetThetaRange; " << endl
|
---|
645 | << "Theta range for this iteration is set to "
|
---|
646 | << fThetaMin << "-" << fThetaMax << " Radians." << endl
|
---|
647 | << "-------------------------------------------------"<< endl;
|
---|
648 |
|
---|
649 | return kTRUE;
|
---|
650 |
|
---|
651 |
|
---|
652 | }
|
---|
653 |
|
---|
654 |
|
---|
655 | Bool_t MFindSupercutsONOFFThetaLoop::SetNormFactorTrainHist()
|
---|
656 | {
|
---|
657 | if (fCosThetaRangeVector.GetSize() <= 1)
|
---|
658 | {
|
---|
659 | *fLog << err << "MFindSupercutsONOFFThetaLoop::SetNormFactorTrainHist: "
|
---|
660 | << endl
|
---|
661 | << "Size of fCosThetaRangeVector is <=1 ... fCosThetaRangeVector must exist and have, at least, 2 components defining a single bin"
|
---|
662 | << endl;
|
---|
663 | return kFALSE;
|
---|
664 | }
|
---|
665 |
|
---|
666 |
|
---|
667 | Int_t HistoNbins = fCosThetaRangeVector.GetSize() - 1;
|
---|
668 | Double_t* BinsVector = fCosThetaRangeVector.GetArray();
|
---|
669 |
|
---|
670 | fNormFactorTrainHist = new TH1F ("NormFactorTrainHist",
|
---|
671 | "NormFactorTrainHist",
|
---|
672 | HistoNbins,BinsVector );
|
---|
673 |
|
---|
674 | return kTRUE;
|
---|
675 |
|
---|
676 | }
|
---|
677 |
|
---|
678 |
|
---|
679 | Bool_t MFindSupercutsONOFFThetaLoop::SetNormFactorTestHist()
|
---|
680 | {
|
---|
681 | if (fCosThetaRangeVector.GetSize() <= 1)
|
---|
682 | {
|
---|
683 | *fLog << err << "MFindSupercutsONOFFThetaLoop::SetNormFactorTestHist: "
|
---|
684 | << endl
|
---|
685 | << "Size of fCosThetaRangeVector is <=1 ... fCosThetaRangeVector must exist and have, at least, 2 components defining a single bin"
|
---|
686 | << endl;
|
---|
687 | return kFALSE;
|
---|
688 | }
|
---|
689 |
|
---|
690 |
|
---|
691 | Int_t HistoNbins = fCosThetaRangeVector.GetSize() - 1;
|
---|
692 | Double_t* BinsVector = fCosThetaRangeVector.GetArray();
|
---|
693 |
|
---|
694 | fNormFactorTestHist = new TH1F ("NormFactorTestHist",
|
---|
695 | "NormFactorTestHist",
|
---|
696 | HistoNbins, BinsVector);
|
---|
697 |
|
---|
698 | return kTRUE;
|
---|
699 | }
|
---|
700 |
|
---|
701 |
|
---|
702 | Bool_t MFindSupercutsONOFFThetaLoop::SetSigmaLiMaTrainHist()
|
---|
703 | {
|
---|
704 | if (fCosThetaRangeVector.GetSize() <= 1)
|
---|
705 | {
|
---|
706 | *fLog << err << "MFindSupercutsONOFFThetaLoop::SetSigmaLiMaTrainHist: "
|
---|
707 | << endl
|
---|
708 | << "Size of fCosThetaRangeVector is <=1 ... fCosThetaRangeVector must exist and have, at least, 2 components defining a single bin"
|
---|
709 | << endl;
|
---|
710 | return kFALSE;
|
---|
711 | }
|
---|
712 |
|
---|
713 |
|
---|
714 | // At some point, definition of bins should be able to
|
---|
715 | // contain variable bins...
|
---|
716 |
|
---|
717 | Int_t HistoNbins = fCosThetaRangeVector.GetSize() - 1;
|
---|
718 |
|
---|
719 | Double_t* BinsVector = fCosThetaRangeVector.GetArray();
|
---|
720 |
|
---|
721 | fSigmaLiMaTrainHist = new TH1F ("SigmaLiMaTrainHist",
|
---|
722 | "SigmaLiMaTrainHist",
|
---|
723 | HistoNbins, BinsVector);
|
---|
724 |
|
---|
725 | return kTRUE;
|
---|
726 | }
|
---|
727 |
|
---|
728 |
|
---|
729 |
|
---|
730 | Bool_t MFindSupercutsONOFFThetaLoop::SetSigmaLiMaTestHist()
|
---|
731 | {
|
---|
732 | if (fCosThetaRangeVector.GetSize() <= 1)
|
---|
733 | {
|
---|
734 | *fLog << err << "MFindSupercutsONOFFThetaLoop::SetSigmaLiMaTestHist: "
|
---|
735 | << endl
|
---|
736 | << "Size of fCosThetaRangeVector is <=1 ... fCosThetaRangeVector must exist and have, at least, 2 components defining a single bin"
|
---|
737 | << endl;
|
---|
738 | return kFALSE;
|
---|
739 | }
|
---|
740 |
|
---|
741 |
|
---|
742 |
|
---|
743 | Int_t HistoNbins = fCosThetaRangeVector.GetSize() - 1;
|
---|
744 |
|
---|
745 | Double_t* BinsVector = fCosThetaRangeVector.GetArray();
|
---|
746 | fSigmaLiMaTestHist = new TH1F ("SigmaLiMaTestHist",
|
---|
747 | "SigmaLiMaTestHist",
|
---|
748 | HistoNbins, BinsVector);
|
---|
749 |
|
---|
750 | return kTRUE;
|
---|
751 |
|
---|
752 | }
|
---|
753 |
|
---|
754 |
|
---|
755 | Bool_t MFindSupercutsONOFFThetaLoop::SetNexTrainHist()
|
---|
756 | {
|
---|
757 | if (fCosThetaRangeVector.GetSize() <= 1)
|
---|
758 | {
|
---|
759 | *fLog << err << "MFindSupercutsONOFFThetaLoop::SetNexTrainHist: "
|
---|
760 | << endl
|
---|
761 | << "Size of fCosThetaRangeVector is <=1 ... fCosThetaRangeVector must exist and have, at least, 2 components defining a single bin"
|
---|
762 | << endl;
|
---|
763 | return kFALSE;
|
---|
764 | }
|
---|
765 |
|
---|
766 |
|
---|
767 | Int_t HistoNbins = fCosThetaRangeVector.GetSize() - 1;
|
---|
768 |
|
---|
769 | Double_t* BinsVector = fCosThetaRangeVector.GetArray();
|
---|
770 |
|
---|
771 | fNexTrainHist = new TH1F ("NexTrainHist",
|
---|
772 | "NexTrainHist",
|
---|
773 | HistoNbins, BinsVector);
|
---|
774 |
|
---|
775 | return kTRUE;
|
---|
776 |
|
---|
777 | }
|
---|
778 |
|
---|
779 |
|
---|
780 | Bool_t MFindSupercutsONOFFThetaLoop::SetNexTestHist()
|
---|
781 | {
|
---|
782 | if (fCosThetaRangeVector.GetSize() <= 1)
|
---|
783 | {
|
---|
784 | *fLog << err << "MFindSupercutsONOFFThetaLoop::SetNexTestHist: "
|
---|
785 | << endl
|
---|
786 | << "Size of fCosThetaRangeVector is <=1 ... fCosThetaRangeVector must exist and have, at least, 2 components defining a single bin"
|
---|
787 | << endl;
|
---|
788 | return kFALSE;
|
---|
789 | }
|
---|
790 |
|
---|
791 |
|
---|
792 | Int_t HistoNbins = fCosThetaRangeVector.GetSize() - 1;
|
---|
793 | Double_t* BinsVector = fCosThetaRangeVector.GetArray();
|
---|
794 |
|
---|
795 | fNexTestHist = new TH1F ("NexTestHist",
|
---|
796 | "NexTestHist",
|
---|
797 | HistoNbins, BinsVector);
|
---|
798 |
|
---|
799 | return kTRUE;
|
---|
800 |
|
---|
801 | }
|
---|
802 |
|
---|
803 | Bool_t MFindSupercutsONOFFThetaLoop::SetSuccessfulThetaBinsHist()
|
---|
804 | {
|
---|
805 | if (fCosThetaRangeVector.GetSize() <= 1)
|
---|
806 | {
|
---|
807 | *fLog << err << "MFindSupercutsONOFFThetaLoop::SetSuccessfulThetaBinsHist(): "
|
---|
808 | << endl
|
---|
809 | << "Size of fCosThetaRangeVector is <=1 ... "
|
---|
810 | << "fCosThetaRangeVector must exist and have >= 2 components (defining a single bin)"
|
---|
811 | << endl;
|
---|
812 | return kFALSE;
|
---|
813 | }
|
---|
814 |
|
---|
815 |
|
---|
816 | Int_t HistoNbins = fCosThetaRangeVector.GetSize() - 1;
|
---|
817 | Double_t* BinsVector = fCosThetaRangeVector.GetArray();
|
---|
818 |
|
---|
819 | fSuccessfulThetaBinsHist = new TH1F ("SuccessfulThetaBinsHist",
|
---|
820 | "SuccessfulThetaBinsHist",
|
---|
821 | HistoNbins, BinsVector);
|
---|
822 |
|
---|
823 | return kTRUE;
|
---|
824 |
|
---|
825 | }
|
---|
826 |
|
---|
827 |
|
---|
828 |
|
---|
829 |
|
---|
830 | Bool_t MFindSupercutsONOFFThetaLoop::SetNexSigmaLiMaNormFactorNEvtsTrainTestHist()
|
---|
831 | {
|
---|
832 | if (fCosThetaRangeVector.GetSize() <= 1)
|
---|
833 | {
|
---|
834 | *fLog << err
|
---|
835 | << "MFindSupercutsONOFFThetaLoop::SetNexSigmaLiMaNormFactorTrainTestHist: "
|
---|
836 | << endl
|
---|
837 | << "Size of fCosThetaRangeVector is <=1 ... "
|
---|
838 | << "fCosThetaRangeVector must exist and have, at least, "
|
---|
839 | << "2 components defining a single bin"
|
---|
840 | << endl;
|
---|
841 | return kFALSE;
|
---|
842 | }
|
---|
843 |
|
---|
844 |
|
---|
845 | Int_t HistoNbins = fCosThetaRangeVector.GetSize() - 1;
|
---|
846 |
|
---|
847 | Double_t* BinsVector = fCosThetaRangeVector.GetArray();
|
---|
848 |
|
---|
849 | char* x_axis_title = {"Cos(Theta)"};
|
---|
850 |
|
---|
851 | fNormFactorTrainHist = new TH1F ("NormFactorTrainHist",
|
---|
852 | "NormFactorTrainHist",
|
---|
853 | HistoNbins, BinsVector);
|
---|
854 |
|
---|
855 |
|
---|
856 | fNormFactorTrainHist -> SetXTitle(x_axis_title);
|
---|
857 | fNormFactorTrainHist -> SetYTitle("Normalization Factor (Non/Noff)");
|
---|
858 |
|
---|
859 |
|
---|
860 | fNormFactorTestHist = new TH1F ("NormFactorTestHist",
|
---|
861 | "NormFactorTestHist",
|
---|
862 | HistoNbins, BinsVector);
|
---|
863 |
|
---|
864 | fNormFactorTestHist -> SetXTitle(x_axis_title);
|
---|
865 | fNormFactorTestHist -> SetYTitle("Normalization Factor (Non/Noff)");
|
---|
866 |
|
---|
867 | fSigmaLiMaTrainHist = new TH1F ("SigmaLiMaTrainHist",
|
---|
868 | "SigmaLiMaTrainHist",
|
---|
869 | HistoNbins, BinsVector);
|
---|
870 |
|
---|
871 | fSigmaLiMaTrainHist -> SetXTitle(x_axis_title);
|
---|
872 | fSigmaLiMaTrainHist -> SetYTitle("Significance");
|
---|
873 |
|
---|
874 |
|
---|
875 | fSigmaLiMaTestHist = new TH1F ("SigmaLiMaTestHist",
|
---|
876 | "SigmaLiMaTestHist",
|
---|
877 | HistoNbins, BinsVector);
|
---|
878 |
|
---|
879 | fSigmaLiMaTestHist -> SetXTitle(x_axis_title);
|
---|
880 | fSigmaLiMaTestHist -> SetYTitle("Significance");
|
---|
881 |
|
---|
882 |
|
---|
883 |
|
---|
884 | fNexTrainHist = new TH1F ("NexTrainHist",
|
---|
885 | "NexTrainHist",
|
---|
886 | HistoNbins, BinsVector);
|
---|
887 |
|
---|
888 | fNexTrainHist -> SetXTitle(x_axis_title);
|
---|
889 | fNexTrainHist -> SetYTitle("Excess Events");
|
---|
890 |
|
---|
891 |
|
---|
892 | fNexTestHist = new TH1F ("NexTestHist",
|
---|
893 | "NexTestHist",
|
---|
894 | HistoNbins, BinsVector);
|
---|
895 |
|
---|
896 | fNexTestHist -> SetXTitle(x_axis_title);
|
---|
897 | fNexTestHist -> SetYTitle("Excess Events");
|
---|
898 |
|
---|
899 |
|
---|
900 |
|
---|
901 | fNEvtsInTrainMatrixONHist = new TH1F("NEvtsInTrainMatrixONHist",
|
---|
902 | "NEvtsInTrainMatrixONHist",
|
---|
903 | HistoNbins, BinsVector);
|
---|
904 |
|
---|
905 |
|
---|
906 | fNEvtsInTrainMatrixONHist -> SetXTitle(x_axis_title);
|
---|
907 | fNEvtsInTrainMatrixONHist -> SetYTitle("Number of ON Events");
|
---|
908 |
|
---|
909 |
|
---|
910 | fNEvtsInTestMatrixONHist = new TH1F("NEvtsInTestMatrixONHist",
|
---|
911 | "NEvtsInTestMatrixONHist",
|
---|
912 | HistoNbins, BinsVector);
|
---|
913 |
|
---|
914 | fNEvtsInTestMatrixONHist -> SetXTitle(x_axis_title);
|
---|
915 | fNEvtsInTestMatrixONHist -> SetYTitle("Number of ON Events");
|
---|
916 |
|
---|
917 |
|
---|
918 | fNEvtsInTrainMatrixOFFHist = new TH1F("NEvtsInTrainMatrixOFFHist",
|
---|
919 | "NEvtsInTrainMatrixOFFHist",
|
---|
920 | HistoNbins, BinsVector);
|
---|
921 |
|
---|
922 | fNEvtsInTrainMatrixOFFHist -> SetXTitle(x_axis_title);
|
---|
923 | fNEvtsInTrainMatrixOFFHist -> SetYTitle("Number of OFF Events");
|
---|
924 |
|
---|
925 |
|
---|
926 | fNEvtsInTestMatrixOFFHist = new TH1F("NEvtsInTestMatrixOFFHist",
|
---|
927 | "NEvtsInTestMatrixOFFHist",
|
---|
928 | HistoNbins, BinsVector);
|
---|
929 |
|
---|
930 | fNEvtsInTestMatrixOFFHist -> SetXTitle(x_axis_title);
|
---|
931 | fNEvtsInTestMatrixOFFHist -> SetYTitle("Number of OFF Events");
|
---|
932 |
|
---|
933 |
|
---|
934 | // Histograms references are removed from the current directory
|
---|
935 |
|
---|
936 | /* IT DOES NOT WORK !!! THESE STATEMENTS PRODUCE SEGMENTATION FAULT.
|
---|
937 | fNormFactorTrainHist -> SetDirectory(NULL);
|
---|
938 |
|
---|
939 | fNormFactorTestHist -> SetDirectory(NULL);
|
---|
940 |
|
---|
941 | fSigmaLiMaTrainHist -> SetDirectory(NULL);
|
---|
942 |
|
---|
943 | fSigmaLiMaTestHist -> SetDirectory(NULL);
|
---|
944 |
|
---|
945 | fSigmaLiMaTestHist -> SetDirectory(NULL);
|
---|
946 |
|
---|
947 | fNexTestHist -> SetDirectory(NULL);
|
---|
948 | */
|
---|
949 |
|
---|
950 |
|
---|
951 | return kTRUE;
|
---|
952 |
|
---|
953 | }
|
---|
954 |
|
---|
955 | void MFindSupercutsONOFFThetaLoop::WriteSuccessfulThetaBinsHistToFile()
|
---|
956 | {
|
---|
957 | *fLog << "MFindSupercutsONOFFThetaLoop : Writing histogram "
|
---|
958 | << " 'SuccessfulThetaBinsHist'"
|
---|
959 | << " (if contents != 0) into root file "
|
---|
960 | << fAlphaDistributionsRootFilename << endl;
|
---|
961 |
|
---|
962 | TFile rootfile(fAlphaDistributionsRootFilename, "UPDATE",
|
---|
963 | "Alpha Distributions for several Theta bins");
|
---|
964 |
|
---|
965 |
|
---|
966 | if (fSuccessfulThetaBinsHist)
|
---|
967 | {
|
---|
968 | if (fSuccessfulThetaBinsHist->GetEntries() > 0.5)
|
---|
969 | fSuccessfulThetaBinsHist -> Write();
|
---|
970 | }
|
---|
971 |
|
---|
972 |
|
---|
973 | rootfile.Close();
|
---|
974 |
|
---|
975 | }
|
---|
976 |
|
---|
977 | void MFindSupercutsONOFFThetaLoop::WriteNexSigmaLiMaNormFactorNEvtsTrainTestHistToFile()
|
---|
978 | {
|
---|
979 |
|
---|
980 | *fLog << "MFindSupercutsONOFFThetaLoop : Writing histograms 'NexSigmaLiMaNormFactorNEvtsTrainTest'"
|
---|
981 | << " (if they have contents != 0) into root file "
|
---|
982 | << fAlphaDistributionsRootFilename << endl;
|
---|
983 |
|
---|
984 | TFile rootfile(fAlphaDistributionsRootFilename, "UPDATE",
|
---|
985 | "Alpha Distributions for several Theta bins");
|
---|
986 |
|
---|
987 |
|
---|
988 | if (fNormFactorTrainHist)
|
---|
989 | {
|
---|
990 | if (fNormFactorTrainHist->GetEntries() > 0.5)
|
---|
991 | fNormFactorTrainHist -> Write();
|
---|
992 | }
|
---|
993 |
|
---|
994 | if (fNormFactorTestHist)
|
---|
995 | {
|
---|
996 | if (fNormFactorTestHist->GetEntries() > 0.5)
|
---|
997 | fNormFactorTestHist -> Write();
|
---|
998 | }
|
---|
999 |
|
---|
1000 | if (fSigmaLiMaTrainHist)
|
---|
1001 | {
|
---|
1002 | if (fSigmaLiMaTrainHist->GetEntries() > 0.5)
|
---|
1003 | fSigmaLiMaTrainHist -> Write();
|
---|
1004 | }
|
---|
1005 |
|
---|
1006 | if (fSigmaLiMaTestHist)
|
---|
1007 | {
|
---|
1008 | if (fSigmaLiMaTestHist->GetEntries() > 0.5)
|
---|
1009 | fSigmaLiMaTestHist -> Write();
|
---|
1010 | }
|
---|
1011 |
|
---|
1012 | if (fNexTrainHist)
|
---|
1013 | {
|
---|
1014 | if (fNexTrainHist->GetEntries() > 0.5)
|
---|
1015 | fNexTrainHist -> Write();
|
---|
1016 | }
|
---|
1017 |
|
---|
1018 | if (fNexTestHist)
|
---|
1019 | {
|
---|
1020 | if (fNexTestHist->GetEntries() > 0.5)
|
---|
1021 | fNexTestHist -> Write();
|
---|
1022 | }
|
---|
1023 |
|
---|
1024 |
|
---|
1025 | if (fNEvtsInTrainMatrixONHist)
|
---|
1026 | {
|
---|
1027 | if (fNEvtsInTrainMatrixONHist->GetEntries() > 0.5)
|
---|
1028 | fNEvtsInTrainMatrixONHist -> Write();
|
---|
1029 | }
|
---|
1030 |
|
---|
1031 | if (fNEvtsInTestMatrixONHist)
|
---|
1032 | {
|
---|
1033 | if (fNEvtsInTestMatrixONHist->GetEntries() > 0.5)
|
---|
1034 | fNEvtsInTestMatrixONHist -> Write();
|
---|
1035 | }
|
---|
1036 |
|
---|
1037 | if (fNEvtsInTrainMatrixOFFHist)
|
---|
1038 | {
|
---|
1039 | if (fNEvtsInTrainMatrixOFFHist->GetEntries() > 0.5)
|
---|
1040 | fNEvtsInTrainMatrixOFFHist -> Write();
|
---|
1041 | }
|
---|
1042 |
|
---|
1043 |
|
---|
1044 | if (fNEvtsInTestMatrixOFFHist)
|
---|
1045 | {
|
---|
1046 | if (fNEvtsInTestMatrixOFFHist->GetEntries() > 0.5)
|
---|
1047 | fNEvtsInTestMatrixOFFHist -> Write();
|
---|
1048 | }
|
---|
1049 |
|
---|
1050 |
|
---|
1051 | rootfile.Close();
|
---|
1052 |
|
---|
1053 | }
|
---|
1054 |
|
---|
1055 | void MFindSupercutsONOFFThetaLoop::SetFractionTrainTestOnOffEvents(
|
---|
1056 | Double_t fontrain,
|
---|
1057 | Double_t fontest,
|
---|
1058 | Double_t fofftrain,
|
---|
1059 | Double_t fofftest)
|
---|
1060 | {
|
---|
1061 |
|
---|
1062 | fWhichFractionTrain = fontrain;
|
---|
1063 | fWhichFractionTest = fontest;
|
---|
1064 | fWhichFractionTrainOFF = fofftrain;
|
---|
1065 | fWhichFractionTestOFF = fofftest;
|
---|
1066 |
|
---|
1067 |
|
---|
1068 | }
|
---|
1069 |
|
---|
1070 | void MFindSupercutsONOFFThetaLoop::SetReadMatricesFromFile (Bool_t b)
|
---|
1071 | {
|
---|
1072 |
|
---|
1073 | *fLog << "---------------------------------------------------" << endl
|
---|
1074 | << "MFindSupercutsONOFFThetaLoop: " << endl;
|
---|
1075 |
|
---|
1076 | if (b)
|
---|
1077 | { *fLog << "Matrices are read from root files." << endl;}
|
---|
1078 | else
|
---|
1079 | { *fLog << "Matrices are produced and stored into root files."
|
---|
1080 | << endl;
|
---|
1081 | }
|
---|
1082 |
|
---|
1083 |
|
---|
1084 | fReadMatricesFromFile = b;
|
---|
1085 | }
|
---|
1086 |
|
---|
1087 | Bool_t MFindSupercutsONOFFThetaLoop::SetAlphaDistributionsRootFilename(
|
---|
1088 | const TString &name)
|
---|
1089 | {
|
---|
1090 |
|
---|
1091 | if (fPathForFiles.IsNull())
|
---|
1092 | {
|
---|
1093 | *fLog << err
|
---|
1094 | << "MFindSupercutsONOFFThetaLoop::SetAlphaDistributionsRootFilename; "
|
---|
1095 | << endl
|
---|
1096 | << "fPathForFiles is not defined yet, hence "
|
---|
1097 | << "name for rootfile containing alpha distributions "
|
---|
1098 | << "can not be defined properly" << endl;
|
---|
1099 |
|
---|
1100 | return kFALSE;
|
---|
1101 | }
|
---|
1102 |
|
---|
1103 |
|
---|
1104 | fAlphaDistributionsRootFilename = (fPathForFiles);
|
---|
1105 | fAlphaDistributionsRootFilename += (name);
|
---|
1106 |
|
---|
1107 |
|
---|
1108 | *fLog << "MFindSupercutsONOFFThetaLoop: fAlphaDistributionsRootFilename = "
|
---|
1109 | << fAlphaDistributionsRootFilename << endl;
|
---|
1110 |
|
---|
1111 |
|
---|
1112 | return kTRUE;
|
---|
1113 |
|
---|
1114 |
|
---|
1115 | }
|
---|
1116 |
|
---|
1117 | // Function to set Names of root files containing matrices and
|
---|
1118 | // optimizedparameters and also fThetaRangeStringVector vector
|
---|
1119 |
|
---|
1120 | Bool_t MFindSupercutsONOFFThetaLoop::SetALLNames()
|
---|
1121 | {
|
---|
1122 | // Names for files only possible if CosThetaRangeVector
|
---|
1123 | // is defined properly
|
---|
1124 | if (fCosThetaRangeVector.GetSize() < 2)
|
---|
1125 | {
|
---|
1126 | *fLog << err
|
---|
1127 | << "MFindSupercutsONOFFThetaLoop::SetALLNames; "
|
---|
1128 | << endl
|
---|
1129 | << "Size of fCosThetaRangeVector is smaller than 2... "
|
---|
1130 | << "fCosThetaRangeVector is not defined properly, "
|
---|
1131 | << "and thus, names for rootfiles containing data matrices "
|
---|
1132 | << "and optimized parameters for the different theta "
|
---|
1133 | << "ranges can not defined properly." << endl;
|
---|
1134 |
|
---|
1135 | return kFALSE;
|
---|
1136 | }
|
---|
1137 |
|
---|
1138 | if (fPathForFiles.IsNull())
|
---|
1139 | {
|
---|
1140 | *fLog << err
|
---|
1141 | << "MFindSupercutsONOFFThetaLoop::Setnames; "
|
---|
1142 | << endl
|
---|
1143 | << "fPathForFiles is not defined yet, hence "
|
---|
1144 | << "names for rootfiles containing data matrices "
|
---|
1145 | << "and optimized parameters for the different theta "
|
---|
1146 | << "ranges can not defined properly." << endl;
|
---|
1147 |
|
---|
1148 | return kFALSE;
|
---|
1149 | }
|
---|
1150 |
|
---|
1151 |
|
---|
1152 | // Name of rootfiles for Supercuts parameters follow the
|
---|
1153 | // follwing structure
|
---|
1154 |
|
---|
1155 | // Path + "OptSCParametersONOFF" + "ThetaRange"
|
---|
1156 | // + int {(fThetaMin_fThetaMax)*1000} + ".root"
|
---|
1157 |
|
---|
1158 |
|
---|
1159 | // Name of rootfiles for Matrices follow the follwing structure
|
---|
1160 |
|
---|
1161 | // Path + "ON(OFF)Data" + "Train(Test)" + "Matrix" + Fraction
|
---|
1162 | // + "ThetaRange" + int {(fThetaMin_fThetaMax)*1000} + ".root"
|
---|
1163 |
|
---|
1164 |
|
---|
1165 | Int_t tmpdim = fCosThetaRangeVector.GetSize() - 1;
|
---|
1166 | const Int_t VectorDimension = tmpdim;
|
---|
1167 |
|
---|
1168 | fThetaRangeStringVector = new TString[VectorDimension];
|
---|
1169 |
|
---|
1170 | fOptSCParamFilenameVector = new TString[VectorDimension];
|
---|
1171 |
|
---|
1172 | fTrainMatrixONFilenameVector = new TString[VectorDimension];
|
---|
1173 | fTestMatrixONFilenameVector = new TString[VectorDimension];
|
---|
1174 | fTrainMatrixOFFFilenameVector = new TString[VectorDimension];
|
---|
1175 | fTestMatrixOFFFilenameVector = new TString[VectorDimension];
|
---|
1176 |
|
---|
1177 |
|
---|
1178 | Int_t ThetaMinTimes1000 = 0;
|
---|
1179 | Int_t ThetaMaxTimes1000 = 0;
|
---|
1180 | Int_t FractionONTrainTimes100 = 0;
|
---|
1181 | Int_t FractionONTestTimes100 = 0;
|
---|
1182 | Int_t FractionOFFTrainTimes100 = 0;
|
---|
1183 | Int_t FractionOFFTestTimes100 = 0;
|
---|
1184 |
|
---|
1185 | char ThetaRangeString[100];
|
---|
1186 |
|
---|
1187 | for (Int_t i = 0; i < VectorDimension; i++)
|
---|
1188 | {
|
---|
1189 | if (!SetThetaRange(i+1))
|
---|
1190 | {
|
---|
1191 | *fLog << "MFindSupercutsONOFFThetaLoop::SetALLNames; "
|
---|
1192 | << endl
|
---|
1193 | << "Values for fThetaMin and fThetaMax could NOT "
|
---|
1194 | << "be computed for theta bin "
|
---|
1195 | << (i+1) << "; SetNames can not go on successfully."
|
---|
1196 | << endl;
|
---|
1197 |
|
---|
1198 | return kFALSE;
|
---|
1199 | }
|
---|
1200 |
|
---|
1201 |
|
---|
1202 | ThetaMinTimes1000 = int(fThetaMin*1000);
|
---|
1203 | ThetaMaxTimes1000 = int(fThetaMax*1000);
|
---|
1204 | sprintf(ThetaRangeString, "ThetaRange%d_%dmRad",
|
---|
1205 | ThetaMinTimes1000, ThetaMaxTimes1000);
|
---|
1206 |
|
---|
1207 | fThetaRangeStringVector[i] = (ThetaRangeString);
|
---|
1208 |
|
---|
1209 | // tmp
|
---|
1210 | *fLog << "ThetaRangeString = " << fThetaRangeStringVector[i] << endl;
|
---|
1211 | // endtmp
|
---|
1212 |
|
---|
1213 | Double_t tmpdouble = 0.0;
|
---|
1214 | tmpdouble = fWhichFractionTrain*100 - int(fWhichFractionTrain*100);
|
---|
1215 |
|
---|
1216 | FractionONTrainTimes100 = tmpdouble < 0.5 ? int(fWhichFractionTrain*100):
|
---|
1217 | int(fWhichFractionTrain*100) + 1;
|
---|
1218 |
|
---|
1219 |
|
---|
1220 | tmpdouble = fWhichFractionTest*100 - int(fWhichFractionTest*100);
|
---|
1221 |
|
---|
1222 | FractionONTestTimes100 = tmpdouble < 0.5 ? int(fWhichFractionTest*100):
|
---|
1223 | int(fWhichFractionTest*100) + 1;
|
---|
1224 |
|
---|
1225 | tmpdouble = fWhichFractionTrainOFF*100 - int(fWhichFractionTrainOFF*100);
|
---|
1226 |
|
---|
1227 | FractionOFFTrainTimes100 = tmpdouble < 0.5 ? int(fWhichFractionTrainOFF*100):
|
---|
1228 | int(fWhichFractionTrainOFF*100) + 1;
|
---|
1229 |
|
---|
1230 | tmpdouble = fWhichFractionTestOFF*100 - int(fWhichFractionTestOFF*100);
|
---|
1231 | FractionOFFTestTimes100 = tmpdouble < 0.5 ? int(fWhichFractionTestOFF*100):
|
---|
1232 | int(fWhichFractionTestOFF*100) + 1;
|
---|
1233 |
|
---|
1234 |
|
---|
1235 |
|
---|
1236 | // File for SC parameters
|
---|
1237 | fOptSCParamFilenameVector[i] = (fPathForFiles);
|
---|
1238 | fOptSCParamFilenameVector[i] += ("OptSCParametersONOFF");
|
---|
1239 | fOptSCParamFilenameVector[i] += (fThetaRangeStringVector[i]);
|
---|
1240 | fOptSCParamFilenameVector[i] += (".root");
|
---|
1241 |
|
---|
1242 |
|
---|
1243 |
|
---|
1244 | // File for Data Matrices
|
---|
1245 |
|
---|
1246 | fTrainMatrixONFilenameVector[i] = (fPathForFiles);
|
---|
1247 | fTrainMatrixONFilenameVector[i] += ("ONDataTrainMatrixFraction");
|
---|
1248 | fTrainMatrixONFilenameVector[i] += (FractionONTrainTimes100);
|
---|
1249 | fTrainMatrixONFilenameVector[i] += (fThetaRangeStringVector[i]);
|
---|
1250 | fTrainMatrixONFilenameVector[i] += (".root");
|
---|
1251 |
|
---|
1252 | fTestMatrixONFilenameVector[i] = (fPathForFiles);
|
---|
1253 | fTestMatrixONFilenameVector[i] += ("ONDataTestMatrixFraction");
|
---|
1254 | fTestMatrixONFilenameVector[i] += (FractionONTestTimes100);
|
---|
1255 | fTestMatrixONFilenameVector[i] += (fThetaRangeStringVector[i]);
|
---|
1256 | fTestMatrixONFilenameVector[i] += (".root");
|
---|
1257 |
|
---|
1258 |
|
---|
1259 | fTrainMatrixOFFFilenameVector[i] = (fPathForFiles);
|
---|
1260 | fTrainMatrixOFFFilenameVector[i] += ("OFFDataTrainMatrixFraction");
|
---|
1261 | fTrainMatrixOFFFilenameVector[i] += (FractionOFFTrainTimes100);
|
---|
1262 | fTrainMatrixOFFFilenameVector[i] += (fThetaRangeStringVector[i]);
|
---|
1263 | fTrainMatrixOFFFilenameVector[i] += (".root");
|
---|
1264 |
|
---|
1265 | fTestMatrixOFFFilenameVector[i] = (fPathForFiles);
|
---|
1266 | fTestMatrixOFFFilenameVector[i] += ("OFFDataTestMatrixFraction");
|
---|
1267 | fTestMatrixOFFFilenameVector[i] += (FractionOFFTestTimes100);
|
---|
1268 | fTestMatrixOFFFilenameVector[i] += (fThetaRangeStringVector[i]);
|
---|
1269 | fTestMatrixOFFFilenameVector[i] += (".root");
|
---|
1270 |
|
---|
1271 | }
|
---|
1272 | // Info concerning names is printed
|
---|
1273 |
|
---|
1274 | *fLog << "--------------------------------------------------" << endl
|
---|
1275 | << "MFindSupercutsONOFFThetaLoop::Setnames; "
|
---|
1276 | << endl
|
---|
1277 | << "Names for root files containing Opt SC Param. and Data "
|
---|
1278 | << "matrices for the different theta bins are the "
|
---|
1279 | << "following ones: " << endl
|
---|
1280 | << "MeanCosTheta, OptSCParm, ONDataTrainMatrix, ..." << endl;
|
---|
1281 |
|
---|
1282 | for (Int_t i = 0; i < VectorDimension; i++)
|
---|
1283 | {
|
---|
1284 | *fLog << fCosThetaBinCenterVector[i] << ", " << endl
|
---|
1285 | << fOptSCParamFilenameVector[i] << ", "<< endl
|
---|
1286 | << fTrainMatrixONFilenameVector[i] << ", " << endl
|
---|
1287 | << fTestMatrixONFilenameVector[i] << ", " << endl
|
---|
1288 | << fTrainMatrixOFFFilenameVector[i] << ", " << endl
|
---|
1289 | << fTestMatrixOFFFilenameVector[i] << ", " << endl
|
---|
1290 | << endl << endl;
|
---|
1291 | }
|
---|
1292 |
|
---|
1293 | *fLog << "--------------------------------------------------" << endl;
|
---|
1294 |
|
---|
1295 | return kTRUE;
|
---|
1296 |
|
---|
1297 | }
|
---|
1298 |
|
---|
1299 |
|
---|
1300 |
|
---|
1301 | Bool_t MFindSupercutsONOFFThetaLoop::LoopOverThetaRanges()
|
---|
1302 | {
|
---|
1303 |
|
---|
1304 |
|
---|
1305 | *fLog << "--------------------------------------------------------------" << endl
|
---|
1306 | << "MFindSupercutsONOFFThetaLoop::LoopOverThetaRanges function starts "
|
---|
1307 | << endl;
|
---|
1308 |
|
---|
1309 |
|
---|
1310 |
|
---|
1311 | if (fDataONRootFilename.IsNull())
|
---|
1312 | {
|
---|
1313 | *fLog << "MFindSupercutsONOFF::LoopOverThetaRanges; root file for ON data is not defined... "
|
---|
1314 | << "aborting..."
|
---|
1315 | << endl;
|
---|
1316 | return kFALSE;
|
---|
1317 | }
|
---|
1318 |
|
---|
1319 |
|
---|
1320 | if (fDataOFFRootFilename.IsNull())
|
---|
1321 | {
|
---|
1322 | *fLog << "MFindSupercutsONOFF:: LoopOverThetaRanges; root file for OFF data is not defined... aborting"
|
---|
1323 | << endl;
|
---|
1324 | return kFALSE;
|
---|
1325 | }
|
---|
1326 |
|
---|
1327 |
|
---|
1328 | if (fHadronnessName.IsNull())
|
---|
1329 | {
|
---|
1330 | *fLog << "MFindSupercutsONOFF::LoopOverThetaRanges; hadronness name for ON data is not defined... aborting"
|
---|
1331 | << endl;
|
---|
1332 | return kFALSE;
|
---|
1333 | }
|
---|
1334 |
|
---|
1335 |
|
---|
1336 | if (fHadronnessNameOFF.IsNull())
|
---|
1337 | {
|
---|
1338 | *fLog << "MFindSupercutsONOFF::LoopOverThetaRanges; hadronness name for OFF data is not defined... aborting"
|
---|
1339 | << endl;
|
---|
1340 | return kFALSE;
|
---|
1341 | }
|
---|
1342 |
|
---|
1343 |
|
---|
1344 | if (fCosThetaRangeVector.GetSize() < 2)
|
---|
1345 |
|
---|
1346 | {
|
---|
1347 | *fLog << "MFindSupercutsONOFFThetaLoop::LoopOverThetaRanges; " << endl
|
---|
1348 | << "Dimension of vector fCosThetaRangeVector "
|
---|
1349 | << "is smaller than 2"
|
---|
1350 | << "Theta ranges are not well defined... aborting ..." << endl;
|
---|
1351 |
|
---|
1352 | return kFALSE;
|
---|
1353 | }
|
---|
1354 |
|
---|
1355 |
|
---|
1356 | // Check if pointers to root file names (opt SC parameters and matrices)
|
---|
1357 | // are empty
|
---|
1358 |
|
---|
1359 |
|
---|
1360 | if (!fTrainMatrixONFilenameVector ||
|
---|
1361 | !fTestMatrixONFilenameVector ||
|
---|
1362 | !fTrainMatrixOFFFilenameVector ||
|
---|
1363 | !fTestMatrixOFFFilenameVector ||
|
---|
1364 | !fOptSCParamFilenameVector)
|
---|
1365 | {
|
---|
1366 | *fLog << "MFindSupercutsONOFFThetaLoop::LoopOverThetaRanges; " << endl
|
---|
1367 | << "Pointers to root files names for Opt SC parameters and data "
|
---|
1368 | << "matrices are not well defined. Aborting... " << endl;
|
---|
1369 |
|
---|
1370 | return kFALSE;
|
---|
1371 | }
|
---|
1372 |
|
---|
1373 |
|
---|
1374 | if (fAlphaDistributionsRootFilename.IsNull())
|
---|
1375 | {
|
---|
1376 | *fLog << "MFindSupercutsONOFFThetaLoop::LoopOverThetaRanges; " << endl
|
---|
1377 | << "Root filename where to store alpha distributions "
|
---|
1378 | << "(fAlphaDistributionsRootFilename) is not defined." << endl
|
---|
1379 | << "Loop over theta angles can not go on..." << endl;
|
---|
1380 |
|
---|
1381 | return kTRUE;
|
---|
1382 | }
|
---|
1383 |
|
---|
1384 |
|
---|
1385 |
|
---|
1386 |
|
---|
1387 | if (!fPsFilename)
|
---|
1388 | {
|
---|
1389 | *fLog << "MFindSupercutsONOFFThetaLoop::LoopOverThetaRanges; " << endl
|
---|
1390 | << "Object from Class TPostScript is not defined. "
|
---|
1391 | << "Plots will be stored in a a file named 'PsFilesTmp.ps' located "
|
---|
1392 | << "in current directory" << endl;
|
---|
1393 |
|
---|
1394 | // fPsFilename = new TPostScript("PsFilesTmp.ps" ,111);
|
---|
1395 |
|
---|
1396 |
|
---|
1397 | }
|
---|
1398 |
|
---|
1399 |
|
---|
1400 |
|
---|
1401 |
|
---|
1402 |
|
---|
1403 |
|
---|
1404 | // All drinks are in house... the party starts !!
|
---|
1405 |
|
---|
1406 | Int_t ThetaLoopIterations = fCosThetaRangeVector.GetSize() - 1;
|
---|
1407 |
|
---|
1408 | *fLog << "looping over the following Cos Theta Angle Ranges : "<< endl;
|
---|
1409 | for (Int_t i = 0; i < ThetaLoopIterations; i++)
|
---|
1410 | {
|
---|
1411 | *fLog << i+1 << " -> "
|
---|
1412 | << fCosThetaRangeVector[i] << "-" << fCosThetaRangeVector[i+1]
|
---|
1413 | << endl;
|
---|
1414 | }
|
---|
1415 |
|
---|
1416 | *fLog << "--------------------------------------------------------------------"
|
---|
1417 | << endl;
|
---|
1418 |
|
---|
1419 |
|
---|
1420 |
|
---|
1421 | // Let's create histograms where SigmaLiMa, Nex and Normfactors will be stored
|
---|
1422 |
|
---|
1423 | if (!SetNexSigmaLiMaNormFactorNEvtsTrainTestHist())
|
---|
1424 | {
|
---|
1425 | *fLog << "MFindSupercutsONOFFThetaLoop::LoopOverThetaRanges; " << endl
|
---|
1426 | << "Histograms where to store Significances, Nex and Normalization "
|
---|
1427 | << "factors, NEvts for all theta bins could not be created... Aborting."
|
---|
1428 | << endl;
|
---|
1429 |
|
---|
1430 | return kFALSE;
|
---|
1431 | }
|
---|
1432 |
|
---|
1433 |
|
---|
1434 | // Let's create histogram where successfully optimized theta bins will be stored
|
---|
1435 |
|
---|
1436 | if (!SetSuccessfulThetaBinsHist())
|
---|
1437 | {
|
---|
1438 | *fLog << "MFindSupercutsONOFFThetaLoop::LoopOverThetaRanges; " << endl
|
---|
1439 | << "Histogram where to store successfully optimized theta bins "
|
---|
1440 | << " could not be created... Aborting."
|
---|
1441 | << endl;
|
---|
1442 |
|
---|
1443 | return kFALSE;
|
---|
1444 | }
|
---|
1445 |
|
---|
1446 |
|
---|
1447 |
|
---|
1448 |
|
---|
1449 |
|
---|
1450 | Double_t NormFactor = 0.0;
|
---|
1451 | Double_t Nex = 0.0;
|
---|
1452 | Double_t SigmaLiMa = 0.0;
|
---|
1453 | Int_t NEvtsON = 0;
|
---|
1454 | Int_t NEvtsOFF = 0;
|
---|
1455 |
|
---|
1456 | Bool_t MinimizationSuccessful = kTRUE;
|
---|
1457 | Double_t MinimizationValue = 1.0; // Value quantifying the minimization, and stored in
|
---|
1458 | // histogram fSuccessfulThetaBinsHist. For the time being it is 0.0 (minimization failed).
|
---|
1459 | // and 1.0 (minimization succeeded). Yet in future it may take other values quantifying
|
---|
1460 | // other aspects...
|
---|
1461 |
|
---|
1462 |
|
---|
1463 | for (Int_t i = 0; i < ThetaLoopIterations; i++)
|
---|
1464 | {
|
---|
1465 |
|
---|
1466 | // Get the Theta range that will be used in this iteration
|
---|
1467 | // Remember that argument of the function is iteration number,
|
---|
1468 | // and it starts at 1...
|
---|
1469 |
|
---|
1470 | if (!SetThetaRange(i+1))
|
---|
1471 | {
|
---|
1472 | *fLog << "MFindSupercutsONOFFThetaLoop::LoopOverThetaRanges; "
|
---|
1473 | << endl
|
---|
1474 | << "Values for fThetaMin and fThetaMax could NOT "
|
---|
1475 | << "be computed for theta bin "
|
---|
1476 | << (i+1) << "; LoopOverThetaRanges is aborting."
|
---|
1477 | << endl;
|
---|
1478 |
|
---|
1479 | return kFALSE;
|
---|
1480 | }
|
---|
1481 |
|
---|
1482 |
|
---|
1483 |
|
---|
1484 | // Object of class MFindSupercutsONOFF will be created and used
|
---|
1485 |
|
---|
1486 |
|
---|
1487 | MFindSupercutsONOFF FindSupercuts("MFindSupercutsONOFF",
|
---|
1488 | "Optimizer for the supercuts");
|
---|
1489 |
|
---|
1490 |
|
---|
1491 | // Set Theta range that will be used
|
---|
1492 | FindSupercuts.SetThetaRange(fThetaMin, fThetaMax);
|
---|
1493 |
|
---|
1494 |
|
---|
1495 | // Alphamax where signal is expected is set
|
---|
1496 | FindSupercuts.SetAlphaSig(fAlphaSig);
|
---|
1497 |
|
---|
1498 | // Bkg alpha region is set
|
---|
1499 | FindSupercuts.SetAlphaBkgMin(fAlphaBkgMin);
|
---|
1500 | FindSupercuts.SetAlphaBkgMax(fAlphaBkgMax);
|
---|
1501 |
|
---|
1502 | // alpha bkg and signal region set in object FindSupercuts
|
---|
1503 | // are re-checked in order to be sure that make sense
|
---|
1504 |
|
---|
1505 | FindSupercuts.CheckAlphaSigBkg();
|
---|
1506 |
|
---|
1507 |
|
---|
1508 | // bining for alpha plots is set
|
---|
1509 |
|
---|
1510 | FindSupercuts.SetAlphaPlotBinining(fNAlphaBins, fAlphaBinLow, fAlphaBinUp);
|
---|
1511 |
|
---|
1512 |
|
---|
1513 |
|
---|
1514 |
|
---|
1515 | // Bool variable NormFactorFromAlphaBkg is set in FindSupercuts
|
---|
1516 | // in order to decide the way in which mormalization factors
|
---|
1517 | // for TRAIN and TEST samples are computed
|
---|
1518 |
|
---|
1519 |
|
---|
1520 | FindSupercuts.SetNormFactorFromAlphaBkg(fNormFactorFromAlphaBkg);
|
---|
1521 |
|
---|
1522 |
|
---|
1523 | // Define size range of events used to fill the data matrices
|
---|
1524 |
|
---|
1525 | FindSupercuts.SetSizeRange(fSizeCutLow, fSizeCutUp);
|
---|
1526 |
|
---|
1527 |
|
---|
1528 |
|
---|
1529 |
|
---|
1530 | // Boolean variable that controls wether the optimization of the
|
---|
1531 | // parameters (MMinuitInterface::CallMinuit(..) in function FindParams(..))
|
---|
1532 | // takes place or not. kTRUE will skip such optimization.
|
---|
1533 | // This variable is useful to test the optmized parameters (previously found
|
---|
1534 | // and stored in root file) on the TRAIN sample.
|
---|
1535 |
|
---|
1536 | FindSupercuts.SetSkipOptimization(fSkipOptimization);
|
---|
1537 |
|
---|
1538 | // Boolean variable that allows the user to write the initial parameters
|
---|
1539 | // into the root file that will be used to store the optimum cuts.
|
---|
1540 | // If fUseInitialSCParams = kTRUE , parameters are written.
|
---|
1541 | // In this way, the initial SC parameters can be applied on the data (train/test)
|
---|
1542 |
|
---|
1543 | FindSupercuts.SetUseInitialSCParams(fUseInitialSCParams);
|
---|
1544 |
|
---|
1545 |
|
---|
1546 | // Set boolean variable that controls wether the theta variable
|
---|
1547 | // is used or not in the computation of the dynamica cuts
|
---|
1548 |
|
---|
1549 | FindSupercuts.SetVariableNotUseTheta(fNotUseTheta);
|
---|
1550 |
|
---|
1551 | // Set boolean variable that controls wether cuts are dynamical
|
---|
1552 | // or static.
|
---|
1553 |
|
---|
1554 | FindSupercuts.SetVariableUseStaticCuts(fUseStaticCuts);
|
---|
1555 |
|
---|
1556 |
|
---|
1557 | FindSupercuts.SetUseFittedQuantities(fUseFittedQuantities);
|
---|
1558 |
|
---|
1559 |
|
---|
1560 | // fGammaEfficiency is set in object FindSupercuts
|
---|
1561 |
|
---|
1562 | FindSupercuts.SetGammaEfficiency (fGammaEfficiency);
|
---|
1563 |
|
---|
1564 |
|
---|
1565 | // Filename with SuperCuts parameters is set
|
---|
1566 |
|
---|
1567 | FindSupercuts.SetFilenameParam(fOptSCParamFilenameVector[i]);
|
---|
1568 |
|
---|
1569 | // Hadroness containers are set
|
---|
1570 |
|
---|
1571 | FindSupercuts.SetHadronnessName(fHadronnessName);
|
---|
1572 | FindSupercuts.SetHadronnessNameOFF(fHadronnessNameOFF);
|
---|
1573 |
|
---|
1574 | // Rootfilename where to store alpha distribution histograms is set
|
---|
1575 |
|
---|
1576 | FindSupercuts.SetAlphaDistributionsRootFilename(fAlphaDistributionsRootFilename);
|
---|
1577 |
|
---|
1578 |
|
---|
1579 | // Set the names of all trees
|
---|
1580 | // used to store the supercuts applied to ON/OFF Train/Test samples
|
---|
1581 |
|
---|
1582 | // In future, names might be controlled by user, but for the time
|
---|
1583 | // being I'll make things simple, and names are set to some defaults
|
---|
1584 | // taking into account the theta range used
|
---|
1585 |
|
---|
1586 | // Containers will be stored in rootfile specified by variable
|
---|
1587 | // fAlphaDistributionsRootFilename
|
---|
1588 |
|
---|
1589 | FindSupercuts.SetSupercutsAppliedTreeNames();
|
---|
1590 |
|
---|
1591 |
|
---|
1592 |
|
---|
1593 |
|
---|
1594 |
|
---|
1595 | // Object of the class TPostScritp is set...
|
---|
1596 | // BUT THAT'S NOT WORKING PROPERLY !!!!!
|
---|
1597 |
|
---|
1598 | // FindSupercuts.SetPostScriptFile (fPsFilename);
|
---|
1599 |
|
---|
1600 | // ********************************************************
|
---|
1601 | // Due to the failure of the use of object TPostScript
|
---|
1602 | // to make a Ps document with all plots, I decided to use the
|
---|
1603 | // standard way (SaveAs(filename.ps)) to store plots related
|
---|
1604 | // to alpha distributions
|
---|
1605 | // for ON and OFF and BEFORE and AFTER cuts (VERY IMPORTANT
|
---|
1606 | // TO COMPUTE EFFICIENCIES IN CUTS)
|
---|
1607 |
|
---|
1608 | // Psfilename is set inside function LoopOverThetaRanges()
|
---|
1609 | // and given to the object MFindSupercutsONOFF created
|
---|
1610 | // within this loop.
|
---|
1611 |
|
---|
1612 | // This will have to be removed as soon as the TPostScript
|
---|
1613 | // solutions works...
|
---|
1614 | // ********************************************************
|
---|
1615 |
|
---|
1616 | TString PsFilenameString = (fPathForFiles);
|
---|
1617 | PsFilenameString += ("SC");
|
---|
1618 | PsFilenameString += (fThetaRangeStringVector[i]);
|
---|
1619 | // Name to be completed inside object MFindSupercutsONOFF
|
---|
1620 |
|
---|
1621 | FindSupercuts.SetPsFilenameString(PsFilenameString);
|
---|
1622 |
|
---|
1623 |
|
---|
1624 | if (fReadMatricesFromFile)
|
---|
1625 | {
|
---|
1626 | *fLog << "MFindSupercutsONOFFThetaLoop::LoopOverThetaRanges; " << endl
|
---|
1627 | << "Reading Data matrices from root files... " << endl;
|
---|
1628 |
|
---|
1629 | FindSupercuts.ReadMatrix(fTrainMatrixONFilenameVector[i],
|
---|
1630 | fTestMatrixONFilenameVector[i]);
|
---|
1631 | FindSupercuts.ReadMatrixOFF(fTrainMatrixOFFFilenameVector[i],
|
---|
1632 | fTestMatrixOFFFilenameVector[i]);
|
---|
1633 |
|
---|
1634 | }
|
---|
1635 |
|
---|
1636 | else
|
---|
1637 | {
|
---|
1638 | *fLog << "MFindSupercutsONOFFThetaLoop::LoopOverThetaRanges; " << endl
|
---|
1639 | << "Reading Data from files " << fDataONRootFilename
|
---|
1640 | << " and " << fDataOFFRootFilename << "(for ON and OFF data "
|
---|
1641 | << "respectively "
|
---|
1642 | << " and producing data matrices that will be stored "
|
---|
1643 | << "in root files. " << endl;
|
---|
1644 |
|
---|
1645 | // ON data
|
---|
1646 |
|
---|
1647 | FindSupercuts.DefineTrainTestMatrixThetaRange(
|
---|
1648 | fDataONRootFilename,
|
---|
1649 | fWhichFractionTrain,
|
---|
1650 | fWhichFractionTest,
|
---|
1651 | fThetaMin, fThetaMax,
|
---|
1652 | fTrainMatrixONFilenameVector[i],
|
---|
1653 | fTestMatrixONFilenameVector[i]);
|
---|
1654 |
|
---|
1655 | // OFF data
|
---|
1656 |
|
---|
1657 | FindSupercuts.DefineTrainTestMatrixOFFThetaRange(
|
---|
1658 | fDataOFFRootFilename,
|
---|
1659 | fWhichFractionTrainOFF,
|
---|
1660 | fWhichFractionTestOFF,
|
---|
1661 | fThetaMin, fThetaMax,
|
---|
1662 | fTrainMatrixOFFFilenameVector[i],
|
---|
1663 | fTestMatrixOFFFilenameVector[i]);
|
---|
1664 |
|
---|
1665 | }
|
---|
1666 |
|
---|
1667 |
|
---|
1668 | MinimizationSuccessful = kTRUE;
|
---|
1669 |
|
---|
1670 | if (fOptimizeParameters)
|
---|
1671 | {
|
---|
1672 | MinimizationSuccessful =
|
---|
1673 | FindSupercuts.FindParams("",fInitSCPar, fInitSCParSteps);
|
---|
1674 |
|
---|
1675 | if (MinimizationSuccessful)
|
---|
1676 | {
|
---|
1677 | // Minimization was successful
|
---|
1678 | // NormFactor, Nex, SigmaLiMa and NEvts Train histograms are updated
|
---|
1679 |
|
---|
1680 | *fLog << "MFindSupercutsONOFFThetaLoop::LoopOverThetaRanges; "
|
---|
1681 | << endl
|
---|
1682 | << "Minimization was successful" << endl
|
---|
1683 | << "Updating values for NormFactor, Nex, SigmaLiMa and NEvts "
|
---|
1684 | << "Train histograms for " << fThetaRangeStringVector[i] << endl;
|
---|
1685 |
|
---|
1686 |
|
---|
1687 | NormFactor = FindSupercuts.GetNormFactorTrain();
|
---|
1688 | SigmaLiMa = FindSupercuts.GetSigmaLiMaTrain();
|
---|
1689 | Nex = FindSupercuts.GetNexTrain();
|
---|
1690 |
|
---|
1691 | NEvtsON = FindSupercuts.GetMatrixTrain()->GetM().GetNrows();
|
---|
1692 | NEvtsOFF = FindSupercuts.GetMatrixTrainOFF()->GetM().GetNrows();
|
---|
1693 |
|
---|
1694 | // MinimizationValue is updated...
|
---|
1695 |
|
---|
1696 | MinimizationValue = 1.0;
|
---|
1697 |
|
---|
1698 | }
|
---|
1699 | else
|
---|
1700 | {
|
---|
1701 | // Minimization was NOT successful
|
---|
1702 | // NormFactor, Nex, SigmaLiMa and NEvts Train histograms are updated
|
---|
1703 |
|
---|
1704 | *fLog << "MFindSupercutsONOFFThetaLoop::LoopOverThetaRanges; "
|
---|
1705 | << endl
|
---|
1706 | << "Minimization was NOT successful" << endl
|
---|
1707 | << "Updating values for NormFactor, and NEvts "
|
---|
1708 | << "Train histograms for " << fThetaRangeStringVector[i] << endl
|
---|
1709 | << "Nex and SigmaLiMa are set to ZERO" << endl;
|
---|
1710 |
|
---|
1711 |
|
---|
1712 | NormFactor = FindSupercuts.GetNormFactorTrain();
|
---|
1713 |
|
---|
1714 | NEvtsON = FindSupercuts.GetMatrixTrain()->GetM().GetNrows();
|
---|
1715 | NEvtsOFF = FindSupercuts.GetMatrixTrainOFF()->GetM().GetNrows();
|
---|
1716 |
|
---|
1717 | SigmaLiMa = 0.0;
|
---|
1718 | Nex = 0.0;
|
---|
1719 |
|
---|
1720 | // MinimizationValue is updated...
|
---|
1721 |
|
---|
1722 | MinimizationValue = 0.0;
|
---|
1723 |
|
---|
1724 | }
|
---|
1725 |
|
---|
1726 |
|
---|
1727 |
|
---|
1728 | cout << "---- Train Histograms will be updated with the following numbers: ------ "
|
---|
1729 | << endl
|
---|
1730 | << "Mean CosTheta for this Thetabin = "
|
---|
1731 | << fCosThetaBinCenterVector[i] << endl
|
---|
1732 | << "Minimization status (1.0 is successful, 0.0 is UNsuccessful) = "
|
---|
1733 | << MinimizationValue << endl
|
---|
1734 | << "Number of ON Events in Train sample used = "
|
---|
1735 | << NEvtsON << endl
|
---|
1736 | << "Number of OFF Events in Train sample used = "
|
---|
1737 | << NEvtsOFF << endl
|
---|
1738 | << "Normalization factor (Non/Noff) for Train sample = "
|
---|
1739 | << NormFactor << endl
|
---|
1740 | << "Excess events and SigmaLiMa: "
|
---|
1741 | << Nex << "; " << SigmaLiMa << endl << endl;
|
---|
1742 |
|
---|
1743 |
|
---|
1744 |
|
---|
1745 | fNormFactorTrainHist -> Fill(fCosThetaBinCenterVector[i],
|
---|
1746 | NormFactor);
|
---|
1747 |
|
---|
1748 | fNexTrainHist -> Fill(fCosThetaBinCenterVector[i],
|
---|
1749 | Nex);
|
---|
1750 |
|
---|
1751 | fSigmaLiMaTrainHist -> Fill(fCosThetaBinCenterVector[i],
|
---|
1752 | SigmaLiMa);
|
---|
1753 |
|
---|
1754 |
|
---|
1755 | fNEvtsInTrainMatrixONHist -> Fill(fCosThetaBinCenterVector[i],
|
---|
1756 | NEvtsON);
|
---|
1757 |
|
---|
1758 | fNEvtsInTrainMatrixOFFHist -> Fill(fCosThetaBinCenterVector[i],
|
---|
1759 | NEvtsOFF);
|
---|
1760 |
|
---|
1761 |
|
---|
1762 | fSuccessfulThetaBinsHist -> Fill(fCosThetaBinCenterVector[i],
|
---|
1763 | MinimizationValue);
|
---|
1764 |
|
---|
1765 | cout << "Histograms (TRAIN) updated successfully " << endl;
|
---|
1766 |
|
---|
1767 |
|
---|
1768 |
|
---|
1769 |
|
---|
1770 | }
|
---|
1771 |
|
---|
1772 |
|
---|
1773 | if (fTestParameters && MinimizationSuccessful)
|
---|
1774 | {
|
---|
1775 |
|
---|
1776 | // Cuts are applied to test sample
|
---|
1777 |
|
---|
1778 | FindSupercuts.TestParamsOnTestSample();
|
---|
1779 |
|
---|
1780 | // NormFactor, Nex and SigmaLiMa Test histograms are updated
|
---|
1781 |
|
---|
1782 | *fLog << "MFindSupercutsONOFFThetaLoop::LoopOverThetaRanges; " << endl
|
---|
1783 | << "Updating values for NormFactor, Nex, SigmaLiMa and NEvts "
|
---|
1784 | << "Test histograms for " << fThetaRangeStringVector[i] << endl;
|
---|
1785 |
|
---|
1786 |
|
---|
1787 | NormFactor = FindSupercuts.GetNormFactorTest();
|
---|
1788 | SigmaLiMa = FindSupercuts.GetSigmaLiMaTest();
|
---|
1789 | Nex = FindSupercuts.GetNexTest();
|
---|
1790 |
|
---|
1791 | NEvtsON = FindSupercuts.GetMatrixTest()->GetM().GetNrows();
|
---|
1792 | NEvtsOFF = FindSupercuts.GetMatrixTestOFF()->GetM().GetNrows();
|
---|
1793 |
|
---|
1794 |
|
---|
1795 | cout << "---- Test Histograms will be updated with the following numbers: ------ "
|
---|
1796 | << endl
|
---|
1797 | << "Mean CosTheta for this Thetabin = "
|
---|
1798 | << fCosThetaBinCenterVector[i] << endl
|
---|
1799 | << "Number of ON Events in Test sample used = "
|
---|
1800 | << NEvtsON << endl
|
---|
1801 | << "Number of OFF Events in Test sample used = "
|
---|
1802 | << NEvtsOFF << endl
|
---|
1803 | << "Normalization factor (Non/Noff) for Test sample = "
|
---|
1804 | << NormFactor << endl
|
---|
1805 | << "Excess events and SigmaLiMa: "
|
---|
1806 | << Nex << "; " << SigmaLiMa << endl << endl;
|
---|
1807 |
|
---|
1808 |
|
---|
1809 |
|
---|
1810 |
|
---|
1811 | fNormFactorTestHist -> Fill(fCosThetaBinCenterVector[i],
|
---|
1812 | NormFactor);
|
---|
1813 |
|
---|
1814 | fNexTestHist -> Fill(fCosThetaBinCenterVector[i],
|
---|
1815 | Nex);
|
---|
1816 |
|
---|
1817 | fSigmaLiMaTestHist -> Fill(fCosThetaBinCenterVector[i],
|
---|
1818 | SigmaLiMa);
|
---|
1819 |
|
---|
1820 |
|
---|
1821 | fNEvtsInTestMatrixONHist -> Fill(fCosThetaBinCenterVector[i],
|
---|
1822 | NEvtsON);
|
---|
1823 |
|
---|
1824 | fNEvtsInTestMatrixOFFHist -> Fill(fCosThetaBinCenterVector[i],
|
---|
1825 | NEvtsOFF);
|
---|
1826 |
|
---|
1827 |
|
---|
1828 | cout << "Histograms (TEST) updated successfully " << endl;
|
---|
1829 |
|
---|
1830 |
|
---|
1831 |
|
---|
1832 |
|
---|
1833 | }
|
---|
1834 |
|
---|
1835 | }
|
---|
1836 |
|
---|
1837 |
|
---|
1838 | // PsFile is closed
|
---|
1839 |
|
---|
1840 |
|
---|
1841 | //*fLog << "Closing PsFile ";
|
---|
1842 | //*fLog << fPsFilename -> GetName() << endl;
|
---|
1843 |
|
---|
1844 | //fPsFilename -> Close();
|
---|
1845 |
|
---|
1846 |
|
---|
1847 |
|
---|
1848 | if (fOptimizeParameters || fTestParameters)
|
---|
1849 | {
|
---|
1850 | // Histograms cotanining SigmaLiMa, Nex and NormFactors and NEvts are written
|
---|
1851 | // in the same root file where alpha distributions where stored
|
---|
1852 | // i.e. fAlphaDistributionsRootFilename
|
---|
1853 |
|
---|
1854 | *fLog << "MFindSupercutsONOFFThetaLoop::LoopOverThetaRanges; " << endl
|
---|
1855 | << "Executing WriteNexSigmaLiMaNormFactorNEvtsTrainTestHistToFile()"
|
---|
1856 | << endl;
|
---|
1857 |
|
---|
1858 | WriteNexSigmaLiMaNormFactorNEvtsTrainTestHistToFile();
|
---|
1859 |
|
---|
1860 |
|
---|
1861 | if (fOptimizeParameters)
|
---|
1862 | { // Histogram containing minimization status for all theta bins defined
|
---|
1863 | // in vector fCosThetaVector is stored in file fAlphaDistributionsRootFilename
|
---|
1864 |
|
---|
1865 | *fLog << "MFindSupercutsONOFFThetaLoop::LoopOverThetaRanges; " << endl
|
---|
1866 | << "Executing WriteSuccessfulThetaBinsHistToFile()"
|
---|
1867 | << endl;
|
---|
1868 |
|
---|
1869 | WriteSuccessfulThetaBinsHistToFile();
|
---|
1870 |
|
---|
1871 | }
|
---|
1872 |
|
---|
1873 | }
|
---|
1874 |
|
---|
1875 |
|
---|
1876 | return kTRUE;
|
---|
1877 | }
|
---|
1878 |
|
---|
1879 |
|
---|
1880 |
|
---|
1881 |
|
---|
1882 |
|
---|
1883 | // Function that loops over the alpha distributions (ON-OFF)
|
---|
1884 | // stored in root file defined by fAlphaDistributionsRootFilename
|
---|
1885 | // and computes the significance and Nex (using MHFindSignificanceONOFF::FindSigma)
|
---|
1886 | // for several cuts in alpha (0-fAlphaSig; in bins defined for alpha distributions
|
---|
1887 | // by user).
|
---|
1888 |
|
---|
1889 | // It initializes the histograms (memeber data), fills them and store them
|
---|
1890 | // in root file defined by fAlphaDistributionsRootFilename. A single histogram
|
---|
1891 | // for each theta bin.
|
---|
1892 |
|
---|
1893 |
|
---|
1894 |
|
---|
1895 | Bool_t MFindSupercutsONOFFThetaLoop::ComputeNexSignificanceVSAlphaSig()
|
---|
1896 | {
|
---|
1897 |
|
---|
1898 |
|
---|
1899 | *fLog << endl << endl
|
---|
1900 | << "****************************************************************************" << endl
|
---|
1901 | << "Computing number of excess events and significance vs " << endl
|
---|
1902 | << " the cut in alpha parameter (fAlphaSig) using function " << endl
|
---|
1903 | << " MFindSupercutsONOFFThetaLoop::ComputeNexSignificanceVSAlphaSig() " << endl
|
---|
1904 | << "****************************************************************************" << endl
|
---|
1905 | << endl;
|
---|
1906 |
|
---|
1907 |
|
---|
1908 | // Check that all elements needed by the function are available
|
---|
1909 |
|
---|
1910 |
|
---|
1911 | if (fAlphaDistributionsRootFilename.IsNull())
|
---|
1912 | {
|
---|
1913 | *fLog << err
|
---|
1914 | << "MFindSupercutsONOFFThetaLoop::ComputeNexSignificanceVSAlphaSig; "
|
---|
1915 | << endl
|
---|
1916 | << "fAlphaDistributionsRootFilename is not defined yet, hence "
|
---|
1917 | << "histograms containing alpha distributions for the different theta "
|
---|
1918 | << "ranges can not be retrieved. " << endl;
|
---|
1919 |
|
---|
1920 | return kFALSE;
|
---|
1921 | }
|
---|
1922 |
|
---|
1923 |
|
---|
1924 |
|
---|
1925 | if (fCosThetaRangeVector.GetSize() < 2)
|
---|
1926 | {
|
---|
1927 | *fLog << err
|
---|
1928 | << "MFindSupercutsONOFFThetaLoop::ComputeNexSignificanceVSAlphaSig; "
|
---|
1929 | << endl
|
---|
1930 | << "Size of fCosThetaRangeVector is smaller than 2... "
|
---|
1931 | << "fCosThetaRangeVector is not defined properly, "
|
---|
1932 | << "and thus, number of theta loops and names "
|
---|
1933 | << " of the histograms containing alpha distributions "
|
---|
1934 | << "for the different theta "
|
---|
1935 | << "ranges can not defined properly." << endl;
|
---|
1936 |
|
---|
1937 | return kFALSE;
|
---|
1938 | }
|
---|
1939 |
|
---|
1940 |
|
---|
1941 |
|
---|
1942 | // It seems that i have all I need for the time being. Let's start the party...
|
---|
1943 |
|
---|
1944 |
|
---|
1945 | Int_t tmpint = 0; // variable that will be used to store temporally integers
|
---|
1946 | Double_t tmpdouble = 0.0; // variable that will be used to store temporally doubles
|
---|
1947 |
|
---|
1948 | tmpint = fCosThetaRangeVector.GetSize() - 1;
|
---|
1949 | const Int_t NThetaBins = tmpint;
|
---|
1950 |
|
---|
1951 |
|
---|
1952 | // Definition of the bins for Nex and Significance vs Alpha histograms
|
---|
1953 | // The bins in Alpha will start with ZERO and will end
|
---|
1954 | // with 5 bins further after fAlphaSig.
|
---|
1955 |
|
---|
1956 | Double_t AlphaBinWidth = double((fAlphaBinUp-fAlphaBinLow)/fNAlphaBins);
|
---|
1957 | Double_t FirstBinLow = 0.0;
|
---|
1958 | tmpint = int((fAlphaSig/AlphaBinWidth) + 5.0);
|
---|
1959 |
|
---|
1960 | while ((FirstBinLow + AlphaBinWidth*tmpint) > fAlphaBkgMin)
|
---|
1961 | {
|
---|
1962 | tmpint--;
|
---|
1963 | }
|
---|
1964 |
|
---|
1965 | *fLog << "Information about the binning of the histograms Nex and Significance vs alphasig "
|
---|
1966 | << endl
|
---|
1967 | << "FirstBinLow, AlphaBinWidth, NBins : "
|
---|
1968 | << FirstBinLow << ", " << AlphaBinWidth << ", " << tmpint-1 << endl;
|
---|
1969 |
|
---|
1970 |
|
---|
1971 |
|
---|
1972 | const Int_t NexSignVSAlphaXBinsVectorDimension = tmpint;
|
---|
1973 |
|
---|
1974 |
|
---|
1975 | Double_t* XBinsVector = new Double_t [NexSignVSAlphaXBinsVectorDimension];
|
---|
1976 |
|
---|
1977 | for (Int_t i = 0; i < NexSignVSAlphaXBinsVectorDimension; i++)
|
---|
1978 | {
|
---|
1979 | XBinsVector[i] = FirstBinLow + i*AlphaBinWidth;
|
---|
1980 | // cout << XBinsVector[i] << endl;
|
---|
1981 | }
|
---|
1982 |
|
---|
1983 |
|
---|
1984 |
|
---|
1985 | // Axis labels for Nex and Significance vs Alpha histograms
|
---|
1986 |
|
---|
1987 | TString x_axis_label = ("Cut in alpha [\\circ]");
|
---|
1988 | TString y_axis_label_nex = ("Excess Events");
|
---|
1989 | TString y_axis_label_significance = ("Significance (LiMa formula 17)");
|
---|
1990 |
|
---|
1991 |
|
---|
1992 | // Arguments for function MHFindSignificanceONOFF::FindSigmaONOFF()
|
---|
1993 | // are defined in here. The same arguments (obviously) will be used
|
---|
1994 | // in all theta bins and TRAIN and TEST samples.
|
---|
1995 |
|
---|
1996 | // variable alphasig will be defined and varied withing the theta loop
|
---|
1997 |
|
---|
1998 | Double_t alphasig = 0.0;
|
---|
1999 |
|
---|
2000 | // Minimum value of alpha for bkg region in ON data
|
---|
2001 | const Double_t alphabkgmin = fAlphaBkgMin;
|
---|
2002 |
|
---|
2003 | // Maximum value of alpha for bkg region in ON data
|
---|
2004 | const Double_t alphabkgmax = fAlphaBkgMax;
|
---|
2005 |
|
---|
2006 | // degree of polynomial function used to fit ON data in background region
|
---|
2007 | const Int_t degree = 2;
|
---|
2008 |
|
---|
2009 | // degree of polynomial function used to fit OFF data in all alpha region (0-90)
|
---|
2010 | const Int_t degreeOFF = 2;
|
---|
2011 |
|
---|
2012 | const Bool_t drawpoly = kFALSE;
|
---|
2013 | const Bool_t fitgauss = kFALSE;
|
---|
2014 | const Bool_t print = kFALSE;
|
---|
2015 |
|
---|
2016 | const Bool_t saveplots = kFALSE; // No plots are saved
|
---|
2017 |
|
---|
2018 | const Bool_t RebinHistograms = kTRUE;
|
---|
2019 | const Bool_t ReduceDegree = kFALSE;
|
---|
2020 |
|
---|
2021 |
|
---|
2022 |
|
---|
2023 | Double_t Nex = 0.0;
|
---|
2024 | Double_t Significance = 0.0;
|
---|
2025 |
|
---|
2026 |
|
---|
2027 | // Names of histograms containing Nex and Significance vs alphasig
|
---|
2028 |
|
---|
2029 | TString NexVSAlphaName = ("NexVSAlpha"); // Names to be completed in theta loop
|
---|
2030 | TString SignificanceVSAlphaName = ("SignificanceVSAlpha");// Names to be completed in theta loop
|
---|
2031 |
|
---|
2032 |
|
---|
2033 |
|
---|
2034 | // Names of histograms to be retrieved from root file
|
---|
2035 | // fAlphaDistributionsRootFilename
|
---|
2036 |
|
---|
2037 | TString NormFactorTrainHistoName = ("NormFactorTrainHist");
|
---|
2038 | TString NormFactorTestHistoName = ("NormFactorTestHist");
|
---|
2039 |
|
---|
2040 | TString SuccessfulThetaBinsHistoName = ("SuccessfulThetaBinsHist");
|
---|
2041 |
|
---|
2042 |
|
---|
2043 |
|
---|
2044 | // Histogram containing information about the success of the optimization
|
---|
2045 | // for the several theta bins. The content of those bins where optimization was not
|
---|
2046 | // successful, are set to ZERO.
|
---|
2047 |
|
---|
2048 |
|
---|
2049 | TH1F* SuccessfulThetaBinsHisto;
|
---|
2050 |
|
---|
2051 |
|
---|
2052 |
|
---|
2053 | // Histograms to store the normalization factors
|
---|
2054 |
|
---|
2055 | TH1F* NormFactorTrainHisto;
|
---|
2056 | TH1F* NormFactorTestHisto;
|
---|
2057 |
|
---|
2058 |
|
---|
2059 |
|
---|
2060 | // Histograms that will be used to store the alpha distributions
|
---|
2061 | // retrieved from fAlphaDistributionsRootFilename, and that will be use as arguments in
|
---|
2062 | // function MHFindSignificanceONOFF::FindSigma
|
---|
2063 |
|
---|
2064 | // The same histograms (pointer to histograms actually) will be used for TRAIN and TEST
|
---|
2065 | // and for the different ThetaNBinsUsed.
|
---|
2066 |
|
---|
2067 | TH1F* AlphaONAfterCuts;
|
---|
2068 | TH1F* AlphaOFFAfterCuts;
|
---|
2069 |
|
---|
2070 |
|
---|
2071 |
|
---|
2072 |
|
---|
2073 | // Vector containing the histogram index of the successfully
|
---|
2074 | // optimized theta bins.
|
---|
2075 | // It is filled with the contents of the histogram
|
---|
2076 | // fSuccessfulThetaBinsHist retrieved from the root
|
---|
2077 | // file fAlphaDistributionsRootFilename
|
---|
2078 | // Remember that histogram index START AT 1 !!!
|
---|
2079 | TArrayI SuccessfulThetaBinsVector;
|
---|
2080 |
|
---|
2081 |
|
---|
2082 |
|
---|
2083 |
|
---|
2084 |
|
---|
2085 |
|
---|
2086 |
|
---|
2087 | // Open root file for retrieving information
|
---|
2088 |
|
---|
2089 | *fLog << "MFindSupercutsONOFFThetaLoop::ComputeNexSignificanceVSAlphaSig()" << endl
|
---|
2090 | << "Opening root file " << fAlphaDistributionsRootFilename << endl;
|
---|
2091 |
|
---|
2092 | TFile rootfile (fAlphaDistributionsRootFilename, "READ");
|
---|
2093 |
|
---|
2094 | // Retrieving SuccessfulThetaBinsHisto from root file.
|
---|
2095 |
|
---|
2096 | *fLog << "MFindSupercutsONOFFThetaLoop::ComputeNexSignificanceVSAlphaSig()" << endl
|
---|
2097 | << "Retrieving SuccessfulThetaBinsHisto from root file and filling vector "
|
---|
2098 | << " SuccessfulThetaBinsVector" << endl;
|
---|
2099 |
|
---|
2100 |
|
---|
2101 | cout << "Getting histogram with name " << SuccessfulThetaBinsHistoName << endl;
|
---|
2102 |
|
---|
2103 | SuccessfulThetaBinsHisto = (TH1F*) rootfile.Get(SuccessfulThetaBinsHistoName);
|
---|
2104 |
|
---|
2105 |
|
---|
2106 |
|
---|
2107 | // Check that bins in this histo correspond with the
|
---|
2108 | // CosTheta intervals defined by vector fCosThetaRangeVector
|
---|
2109 |
|
---|
2110 |
|
---|
2111 | if (SuccessfulThetaBinsHisto -> GetNbinsX() != NThetaBins)
|
---|
2112 | {
|
---|
2113 | *fLog << "MFindSupercutsONOFFThetaLoop::ComputeNexSignificanceVSAlphaSig()"
|
---|
2114 | << endl
|
---|
2115 | << "Number of theta bins defined by vector fCosThetaRangeVector ("
|
---|
2116 | << NThetaBins << ") does not match with the bins in histogram "
|
---|
2117 | << SuccessfulThetaBinsHistoName << "("
|
---|
2118 | << SuccessfulThetaBinsHisto -> GetNbinsX() << ")" << endl
|
---|
2119 | << "Function execution will be ABORTED." << endl;
|
---|
2120 | return kFALSE;
|
---|
2121 | }
|
---|
2122 |
|
---|
2123 |
|
---|
2124 | // Filling vector SuccessfulThetaBinsVector
|
---|
2125 |
|
---|
2126 | tmpint = 0;
|
---|
2127 | for (Int_t i = 0; i < NThetaBins; i++)
|
---|
2128 | {
|
---|
2129 | if ((SuccessfulThetaBinsHisto -> GetBinContent(i+1)) > 0.5)
|
---|
2130 | {
|
---|
2131 | tmpint++;
|
---|
2132 | }
|
---|
2133 | }
|
---|
2134 |
|
---|
2135 | SuccessfulThetaBinsVector.Set(tmpint);
|
---|
2136 |
|
---|
2137 | tmpint = 0;
|
---|
2138 | for (Int_t i = 0; i < NThetaBins; i++)
|
---|
2139 | {
|
---|
2140 | if ((SuccessfulThetaBinsHisto -> GetBinContent(i+1)) > 0.5)
|
---|
2141 | {
|
---|
2142 | if(tmpint < SuccessfulThetaBinsVector.GetSize())
|
---|
2143 | {
|
---|
2144 | SuccessfulThetaBinsVector[tmpint] = i+1;
|
---|
2145 | tmpint++;
|
---|
2146 | }
|
---|
2147 | else
|
---|
2148 | {
|
---|
2149 | *fLog << "MFindSupercutsONOFFThetaLoop::ComputeNexSignificanceVSAlphaSig()"
|
---|
2150 | << endl
|
---|
2151 | << "Problem when filling vector 'SuccessfulThetaBinsVector'"
|
---|
2152 | << endl
|
---|
2153 | << "Function execution will be ABORTED." << endl;
|
---|
2154 | return kFALSE;
|
---|
2155 | }
|
---|
2156 |
|
---|
2157 |
|
---|
2158 | }
|
---|
2159 | }
|
---|
2160 |
|
---|
2161 | // variable defining the theta bins where optimization was successful, and
|
---|
2162 | // hence, alpha distributions stored in root file are meaningful (up to some extent)
|
---|
2163 |
|
---|
2164 | tmpint = SuccessfulThetaBinsVector.GetSize();
|
---|
2165 | const Int_t SuccessfulNThetaBins = tmpint;
|
---|
2166 |
|
---|
2167 | // TEMP
|
---|
2168 | // cout << "const Int_t SuccessfulNThetaBins = " << SuccessfulNThetaBins << endl;
|
---|
2169 | // ENDTEMP
|
---|
2170 |
|
---|
2171 | // Vectors of pointers to histograms are defined to store the histograms
|
---|
2172 | // containing Nex and significance for TRAIN and TEST sample for all the
|
---|
2173 | // successful theta bins
|
---|
2174 |
|
---|
2175 | // Histograms will be initialized in the loop over theta bins
|
---|
2176 |
|
---|
2177 | TH1F* NexVSAlphaSigTrainVector[SuccessfulNThetaBins];
|
---|
2178 | TH1F* NexVSAlphaSigTestVector[SuccessfulNThetaBins];
|
---|
2179 |
|
---|
2180 |
|
---|
2181 | TH1F* SignificanceVSAlphaSigTrainVector[SuccessfulNThetaBins];
|
---|
2182 | TH1F* SignificanceVSAlphaSigTestVector[SuccessfulNThetaBins];
|
---|
2183 |
|
---|
2184 |
|
---|
2185 | // ***************************************************
|
---|
2186 | // HISTOGRAMS FOR TRAIN SAMPLE ARE COMPUTED
|
---|
2187 | // ***************************************************
|
---|
2188 |
|
---|
2189 |
|
---|
2190 |
|
---|
2191 | if (fOptimizeParameters)
|
---|
2192 | { // Computation for TRAIN sample
|
---|
2193 |
|
---|
2194 |
|
---|
2195 | *fLog << endl << endl
|
---|
2196 | << "****************************************************************************" << endl
|
---|
2197 | << "Computing Nex and significance vs cut in alpha (alphasig) for TRAIN sample" << endl
|
---|
2198 | << "****************************************************************************" << endl
|
---|
2199 | << endl << endl;
|
---|
2200 |
|
---|
2201 |
|
---|
2202 | // Retrieving NormFactorTrainHisto from root file.
|
---|
2203 |
|
---|
2204 | cout << "Getting histogram with name " << NormFactorTrainHistoName << endl;
|
---|
2205 |
|
---|
2206 | NormFactorTrainHisto = (TH1F*) rootfile.Get(NormFactorTrainHistoName);
|
---|
2207 |
|
---|
2208 |
|
---|
2209 | // Check that bins in this histo correspond with the
|
---|
2210 | // CosTheta intervals defined by vector fCosThetaRangeVector
|
---|
2211 |
|
---|
2212 | ////////////////// START CHECK //////////////////////////////////////
|
---|
2213 |
|
---|
2214 | if (NormFactorTrainHisto -> GetNbinsX() != NThetaBins)
|
---|
2215 | {
|
---|
2216 | *fLog << "MFindSupercutsONOFFThetaLoop::ComputeNexSignificanceVSAlphaSig()"
|
---|
2217 | << endl
|
---|
2218 | << "Number of theta bins defined by vector fCosThetaRangeVector ("
|
---|
2219 | << NThetaBins << ") does not match with the bins in histogram "
|
---|
2220 | << NormFactorTrainHistoName << "("
|
---|
2221 | << NormFactorTrainHisto -> GetNbinsX() << ")" << endl
|
---|
2222 | << "Function execution will be ABORTED." << endl;
|
---|
2223 | return kFALSE;
|
---|
2224 | }
|
---|
2225 |
|
---|
2226 | // histo test
|
---|
2227 | /*
|
---|
2228 | cout << "NORMFACTORTRAINHIST Test: BinCenter and Value" << endl;
|
---|
2229 | for (Int_t k = 0; k < NThetaBins; k++)
|
---|
2230 | {
|
---|
2231 |
|
---|
2232 | cout << NormFactorTrainHisto -> GetBinCenter(k+1)
|
---|
2233 | << "; " << NormFactorTrainHisto -> GetBinContent(k+1)<< endl;
|
---|
2234 | }
|
---|
2235 | */
|
---|
2236 |
|
---|
2237 | for (Int_t i = 0; i < SuccessfulThetaBinsVector.GetSize(); i++)
|
---|
2238 | {
|
---|
2239 | tmpdouble = NormFactorTrainHisto -> GetBinCenter(SuccessfulThetaBinsVector[i]);
|
---|
2240 | if (fCosThetaRangeVector[SuccessfulThetaBinsVector[i]-1] > tmpdouble ||
|
---|
2241 | fCosThetaRangeVector[SuccessfulThetaBinsVector[i]] < tmpdouble)
|
---|
2242 | {
|
---|
2243 | *fLog << "MFindSupercutsONOFFThetaLoop::ComputeNexSignificanceVSAlphaSig()"
|
---|
2244 | << endl
|
---|
2245 | << "Bins defined by vector fCosThetaRangeVector "
|
---|
2246 | << "do not match with the bins of histogram "
|
---|
2247 | << NormFactorTrainHistoName << endl
|
---|
2248 | << "CosTheta: " << fCosThetaRangeVector[SuccessfulThetaBinsVector[i]-1]
|
---|
2249 | << "-" << fCosThetaRangeVector[SuccessfulThetaBinsVector[i]] << endl
|
---|
2250 | << "NormFactorHist Bin: " << tmpdouble << endl
|
---|
2251 | << "Function execution will be ABORTED." << endl;
|
---|
2252 |
|
---|
2253 | return kFALSE;
|
---|
2254 | }
|
---|
2255 |
|
---|
2256 | }
|
---|
2257 |
|
---|
2258 | ////////////////// END CHECK //////////////////////////////////////
|
---|
2259 |
|
---|
2260 |
|
---|
2261 |
|
---|
2262 |
|
---|
2263 | // Check that the theta range string (needed to compute names for alpha histograms)
|
---|
2264 | // is well defined
|
---|
2265 |
|
---|
2266 | for (Int_t i = 0; i < SuccessfulThetaBinsVector.GetSize(); i++)
|
---|
2267 | {
|
---|
2268 | if (fThetaRangeStringVector[SuccessfulThetaBinsVector[i]-1].IsNull())
|
---|
2269 | {
|
---|
2270 | *fLog << "MFindSupercutsONOFFThetaLoop::ComputeNexSignificanceVSAlphaSig()"
|
---|
2271 | << endl
|
---|
2272 | << "Component " << SuccessfulThetaBinsVector[i]-1
|
---|
2273 | << " of fThetaRangeStringVector is EMPTY"
|
---|
2274 | << endl
|
---|
2275 | << "Function execution will be ABORTED." << endl;
|
---|
2276 |
|
---|
2277 | return kFALSE;
|
---|
2278 |
|
---|
2279 |
|
---|
2280 | }
|
---|
2281 | }
|
---|
2282 |
|
---|
2283 |
|
---|
2284 | // *************************************************************
|
---|
2285 | // *************************************************************
|
---|
2286 |
|
---|
2287 | // Normfactors and successfultheta bins are ok.
|
---|
2288 | // I can now compute the Nex and Sigma for train sample (ON-OFF)
|
---|
2289 |
|
---|
2290 | // *************************************************************
|
---|
2291 | // *************************************************************
|
---|
2292 |
|
---|
2293 |
|
---|
2294 |
|
---|
2295 |
|
---|
2296 | // Loop over the several successful theta bins doing the following
|
---|
2297 |
|
---|
2298 | // 0) Initializing histograms Nex and Significance (with names!!)
|
---|
2299 | // 1) Retrieving alpha histograms (ON-OFF), and normalization factor
|
---|
2300 | // 2) Computing Nex and Signficance vs alpha and filling histograms
|
---|
2301 |
|
---|
2302 |
|
---|
2303 | *fLog << endl
|
---|
2304 | << "Starting loop over the successfully optimized theta ranges " << endl
|
---|
2305 | << endl;
|
---|
2306 |
|
---|
2307 | Double_t NormalizationFactor = 0.0;
|
---|
2308 |
|
---|
2309 | for (Int_t i = 0; i < SuccessfulThetaBinsVector.GetSize(); i++)
|
---|
2310 | {
|
---|
2311 |
|
---|
2312 | *fLog << fThetaRangeStringVector[SuccessfulThetaBinsVector[i]-1] << endl << endl;
|
---|
2313 |
|
---|
2314 | // 0) Initializing histograms Nex and Significance (with names!!)
|
---|
2315 |
|
---|
2316 | TString NexVSAlphaNameTmp = (NexVSAlphaName);
|
---|
2317 | TString SignificanceVSAlphaNameTmp = (SignificanceVSAlphaName);
|
---|
2318 |
|
---|
2319 | NexVSAlphaNameTmp += "Train";
|
---|
2320 | NexVSAlphaNameTmp += fThetaRangeStringVector[SuccessfulThetaBinsVector[i]-1];
|
---|
2321 |
|
---|
2322 | SignificanceVSAlphaNameTmp += "Train";
|
---|
2323 | SignificanceVSAlphaNameTmp += fThetaRangeStringVector[SuccessfulThetaBinsVector[i]-1];
|
---|
2324 |
|
---|
2325 |
|
---|
2326 | // TEMP
|
---|
2327 | // cout << NexVSAlphaNameTmp << endl;
|
---|
2328 | // cout << SignificanceVSAlphaNameTmp << endl;
|
---|
2329 | // ENDTEMP
|
---|
2330 |
|
---|
2331 |
|
---|
2332 |
|
---|
2333 |
|
---|
2334 | NexVSAlphaSigTrainVector[i] = new TH1F (NexVSAlphaNameTmp.Data(),
|
---|
2335 | NexVSAlphaNameTmp.Data(),
|
---|
2336 | NexSignVSAlphaXBinsVectorDimension - 1,
|
---|
2337 | XBinsVector);
|
---|
2338 |
|
---|
2339 | // It is important to unlink the the created histograms from
|
---|
2340 | // the current directory (defined now by "rootfile"),
|
---|
2341 | // so that we can do whatever we want to
|
---|
2342 | // with them regardless of such directory is opened or closed
|
---|
2343 | // The reference of the histogram is removed from the current
|
---|
2344 | // directory by means of the member function "SetDirectory(0)"
|
---|
2345 |
|
---|
2346 |
|
---|
2347 | NexVSAlphaSigTrainVector[i] -> SetDirectory(0);
|
---|
2348 |
|
---|
2349 |
|
---|
2350 | SignificanceVSAlphaSigTrainVector[i] = new TH1F (SignificanceVSAlphaNameTmp.Data(),
|
---|
2351 | SignificanceVSAlphaNameTmp.Data(),
|
---|
2352 | NexSignVSAlphaXBinsVectorDimension - 1,
|
---|
2353 | XBinsVector);
|
---|
2354 |
|
---|
2355 | SignificanceVSAlphaSigTrainVector[i] -> SetDirectory(0);
|
---|
2356 |
|
---|
2357 |
|
---|
2358 |
|
---|
2359 | // axis are defined
|
---|
2360 |
|
---|
2361 | NexVSAlphaSigTrainVector[i] -> SetXTitle(x_axis_label.Data());
|
---|
2362 | NexVSAlphaSigTrainVector[i] -> SetYTitle(y_axis_label_nex.Data());
|
---|
2363 |
|
---|
2364 | SignificanceVSAlphaSigTrainVector[i] ->SetXTitle(x_axis_label.Data());
|
---|
2365 | SignificanceVSAlphaSigTrainVector[i] ->SetYTitle(y_axis_label_significance.Data());
|
---|
2366 |
|
---|
2367 |
|
---|
2368 | // 1) Retrieving alpha histograms (ON-OFF) and normalization factor
|
---|
2369 |
|
---|
2370 | NormalizationFactor =
|
---|
2371 | NormFactorTrainHisto -> GetBinContent(SuccessfulThetaBinsVector[i]);
|
---|
2372 |
|
---|
2373 |
|
---|
2374 | TString TmpAlphaHistoName = ("TrainAlphaAfterCuts");
|
---|
2375 | TmpAlphaHistoName += fThetaRangeStringVector[SuccessfulThetaBinsVector[i]-1];
|
---|
2376 |
|
---|
2377 | *fLog << "Getting Histogram with name " << endl
|
---|
2378 | << TmpAlphaHistoName << endl;
|
---|
2379 |
|
---|
2380 | AlphaONAfterCuts = (TH1F*) rootfile.Get(TmpAlphaHistoName);
|
---|
2381 |
|
---|
2382 |
|
---|
2383 | TString TmpAlphaOFFHistoName = ("TrainAlphaOFFAfterCuts");
|
---|
2384 | TmpAlphaOFFHistoName += fThetaRangeStringVector[SuccessfulThetaBinsVector[i]-1];
|
---|
2385 |
|
---|
2386 | *fLog << "Getting Histogram with name " << endl
|
---|
2387 | << TmpAlphaOFFHistoName << endl;
|
---|
2388 |
|
---|
2389 | AlphaOFFAfterCuts = (TH1F*) rootfile.Get(TmpAlphaOFFHistoName);
|
---|
2390 |
|
---|
2391 |
|
---|
2392 | // 2) Computing Nex and Signficance vs alpha and filling histograms
|
---|
2393 | // Significance is found using MHFindSignificanceONOFF::FindSigmaONOFF
|
---|
2394 |
|
---|
2395 |
|
---|
2396 | *fLog << endl << endl
|
---|
2397 | << "Starting loop to compute Nex and Significance vs alphasig " << endl
|
---|
2398 | << "for TRAIN sample and "
|
---|
2399 | << fThetaRangeStringVector[SuccessfulThetaBinsVector[i]-1]
|
---|
2400 | << " :" << endl;
|
---|
2401 |
|
---|
2402 |
|
---|
2403 | for (Int_t j = 1; j < NexSignVSAlphaXBinsVectorDimension; j++)
|
---|
2404 | {
|
---|
2405 |
|
---|
2406 |
|
---|
2407 | alphasig = XBinsVector[j];
|
---|
2408 |
|
---|
2409 | MHFindSignificanceONOFF findsig;
|
---|
2410 | findsig.SetRebin(RebinHistograms);
|
---|
2411 | findsig.SetReduceDegree(ReduceDegree);
|
---|
2412 | findsig.SetUseFittedQuantities(fUseFittedQuantities);
|
---|
2413 |
|
---|
2414 | // Dummy name for psfile
|
---|
2415 | TString psfilename = ("DummyPs");
|
---|
2416 |
|
---|
2417 | findsig.FindSigmaONOFF(AlphaONAfterCuts, AlphaOFFAfterCuts,
|
---|
2418 | NormalizationFactor,
|
---|
2419 | alphabkgmin, alphabkgmax,
|
---|
2420 | degree, degreeOFF,
|
---|
2421 | alphasig,
|
---|
2422 | drawpoly, fitgauss, print,
|
---|
2423 | saveplots, psfilename);
|
---|
2424 |
|
---|
2425 |
|
---|
2426 | Significance = double(findsig.GetSignificance());
|
---|
2427 |
|
---|
2428 | // so far there is not a variable that contains
|
---|
2429 | // NexONOFF or NexONOFFFitted (in class MHFindSignificanceONOFF)
|
---|
2430 | // according to the value of the variable fUseFittedQuantities.
|
---|
2431 | // As soon as I have some time I will implement it (with the
|
---|
2432 | // corresponging GET member function) and will change
|
---|
2433 | // MFindSignificanceONOFF accordingly to make use of it.
|
---|
2434 |
|
---|
2435 | // For the time being, I will survive with an "if statement".
|
---|
2436 |
|
---|
2437 |
|
---|
2438 | if(fUseFittedQuantities)
|
---|
2439 | {
|
---|
2440 | Nex = double(findsig.GetNexONOFFFitted());
|
---|
2441 | }
|
---|
2442 | else
|
---|
2443 | {
|
---|
2444 | Nex = double(findsig.GetNexONOFF());
|
---|
2445 | }
|
---|
2446 |
|
---|
2447 |
|
---|
2448 |
|
---|
2449 | *fLog << endl << "Cut in alpha, Nex and significance : "
|
---|
2450 | << alphasig << ", " << Nex << ", " << Significance << endl << endl;
|
---|
2451 |
|
---|
2452 |
|
---|
2453 | // updating histograms
|
---|
2454 |
|
---|
2455 | tmpdouble = alphasig - AlphaBinWidth/2.0;
|
---|
2456 |
|
---|
2457 | NexVSAlphaSigTrainVector[i] -> Fill(tmpdouble, Nex);
|
---|
2458 | SignificanceVSAlphaSigTrainVector[i] -> Fill(tmpdouble, Significance);
|
---|
2459 |
|
---|
2460 |
|
---|
2461 |
|
---|
2462 |
|
---|
2463 |
|
---|
2464 |
|
---|
2465 |
|
---|
2466 | }
|
---|
2467 |
|
---|
2468 |
|
---|
2469 | // Dynamic memory allocated in this loop is released
|
---|
2470 |
|
---|
2471 |
|
---|
2472 | *fLog << "Memory allocated by temporal alpha histograms "
|
---|
2473 | << "will be released" << endl;
|
---|
2474 |
|
---|
2475 | delete AlphaONAfterCuts;
|
---|
2476 | delete AlphaOFFAfterCuts;
|
---|
2477 |
|
---|
2478 | }
|
---|
2479 |
|
---|
2480 |
|
---|
2481 | }
|
---|
2482 |
|
---|
2483 |
|
---|
2484 |
|
---|
2485 | // ***************************************************
|
---|
2486 | // HISTOGRAMS FOR TEST SAMPLE ARE COMPUTED
|
---|
2487 | // ***************************************************
|
---|
2488 |
|
---|
2489 |
|
---|
2490 |
|
---|
2491 |
|
---|
2492 | if(fTestParameters)
|
---|
2493 | {// computation for TEST sample
|
---|
2494 |
|
---|
2495 | *fLog << endl << endl
|
---|
2496 | << "**************************************************************************" << endl
|
---|
2497 | << "Computing Nex and significance vs cut in alpha (alphasig) for TEST sample" << endl
|
---|
2498 | << "**************************************************************************" << endl
|
---|
2499 | << endl << endl;
|
---|
2500 |
|
---|
2501 | // Retrieving NormFactorTestHisto from root file.
|
---|
2502 |
|
---|
2503 | *fLog << "Getting histogram with name " << NormFactorTestHistoName << endl;
|
---|
2504 |
|
---|
2505 | NormFactorTestHisto = (TH1F*) rootfile.Get(NormFactorTestHistoName);
|
---|
2506 |
|
---|
2507 |
|
---|
2508 | // Check that bins in this histo correspond with the
|
---|
2509 | // CosTheta intervals defined by vector fCosThetaRangeVector
|
---|
2510 |
|
---|
2511 | ////////////////// START CHECK //////////////////////////////////////
|
---|
2512 |
|
---|
2513 | if (NormFactorTestHisto -> GetNbinsX() != NThetaBins)
|
---|
2514 | {
|
---|
2515 | *fLog << "MFindSupercutsONOFFThetaLoop::ComputeNexSignificanceVSAlphaSig()"
|
---|
2516 | << endl
|
---|
2517 | << "Number of theta bins defined by vector fCosThetaRangeVector ("
|
---|
2518 | << NThetaBins << ") does not match with the bins in histogram "
|
---|
2519 | << NormFactorTestHistoName << "("
|
---|
2520 | << NormFactorTestHisto -> GetNbinsX() << ")" << endl
|
---|
2521 | << "Function execution will be ABORTED." << endl;
|
---|
2522 | return kFALSE;
|
---|
2523 | }
|
---|
2524 |
|
---|
2525 | // histo test
|
---|
2526 | /*
|
---|
2527 | cout << "NORMFACTORTRAINHIST Test: BinCenter and Value" << endl;
|
---|
2528 | for (Int_t k = 0; k < NThetaBins; k++)
|
---|
2529 | {
|
---|
2530 |
|
---|
2531 | cout << NormFactorTestHisto -> GetBinCenter(k+1)
|
---|
2532 | << "; " << NormFactorTestHisto -> GetBinContent(k+1)<< endl;
|
---|
2533 | }
|
---|
2534 | */
|
---|
2535 |
|
---|
2536 | for (Int_t i = 0; i < SuccessfulThetaBinsVector.GetSize(); i++)
|
---|
2537 | {
|
---|
2538 | tmpdouble = NormFactorTestHisto -> GetBinCenter(SuccessfulThetaBinsVector[i]);
|
---|
2539 | if (fCosThetaRangeVector[SuccessfulThetaBinsVector[i]-1] > tmpdouble ||
|
---|
2540 | fCosThetaRangeVector[SuccessfulThetaBinsVector[i]] < tmpdouble)
|
---|
2541 | {
|
---|
2542 | *fLog << "MFindSupercutsONOFFThetaLoop::ComputeNexSignificanceVSAlphaSig()"
|
---|
2543 | << endl
|
---|
2544 | << "Bins defined by vector fCosThetaRangeVector "
|
---|
2545 | << "do not match with the bins of histogram "
|
---|
2546 | << NormFactorTestHistoName << endl
|
---|
2547 | << "CosTheta: " << fCosThetaRangeVector[SuccessfulThetaBinsVector[i]-1]
|
---|
2548 | << "-" << fCosThetaRangeVector[SuccessfulThetaBinsVector[i]] << endl
|
---|
2549 | << "NormFactorHist Bin: " << tmpdouble << endl
|
---|
2550 | << "Function execution will be ABORTED." << endl;
|
---|
2551 |
|
---|
2552 | return kFALSE;
|
---|
2553 | }
|
---|
2554 |
|
---|
2555 | }
|
---|
2556 |
|
---|
2557 | ////////////////// END CHECK //////////////////////////////////////
|
---|
2558 |
|
---|
2559 |
|
---|
2560 |
|
---|
2561 |
|
---|
2562 | // Check that the theta range string (needed to compute names for alpha histograms)
|
---|
2563 | // is well defined
|
---|
2564 |
|
---|
2565 | for (Int_t i = 0; i < SuccessfulThetaBinsVector.GetSize(); i++)
|
---|
2566 | {
|
---|
2567 | if (fThetaRangeStringVector[SuccessfulThetaBinsVector[i]-1].IsNull())
|
---|
2568 | {
|
---|
2569 | *fLog << "MFindSupercutsONOFFThetaLoop::ComputeNexSignificanceVSAlphaSig()"
|
---|
2570 | << endl
|
---|
2571 | << "Component " << SuccessfulThetaBinsVector[i]-1
|
---|
2572 | << " of fThetaRangeStringVector is EMPTY"
|
---|
2573 | << endl
|
---|
2574 | << "Function execution will be ABORTED." << endl;
|
---|
2575 |
|
---|
2576 | return kFALSE;
|
---|
2577 |
|
---|
2578 |
|
---|
2579 | }
|
---|
2580 | }
|
---|
2581 |
|
---|
2582 |
|
---|
2583 | // *************************************************************
|
---|
2584 | // *************************************************************
|
---|
2585 |
|
---|
2586 | // Normfactors and successfultheta bins are ok.
|
---|
2587 | // I can now compute the Nex and Sigma for test sample (ON-OFF)
|
---|
2588 |
|
---|
2589 | // *************************************************************
|
---|
2590 | // *************************************************************
|
---|
2591 |
|
---|
2592 | // Loop over the several successful theta bins doing the following
|
---|
2593 |
|
---|
2594 | // 0) Initializing histograms Nex and Significance (with names!!)
|
---|
2595 | // 1) Retrieving alpha histograms (ON-OFF), and normalization factor
|
---|
2596 | // 2) Computing Nex and Signficance vs alpha and filling histograms
|
---|
2597 |
|
---|
2598 |
|
---|
2599 | *fLog << endl
|
---|
2600 | << "Starting loop over the successfully optimized theta ranges " << endl
|
---|
2601 | << endl;
|
---|
2602 |
|
---|
2603 |
|
---|
2604 | Double_t NormalizationFactor = 0.0;
|
---|
2605 |
|
---|
2606 | for (Int_t i = 0; i < SuccessfulThetaBinsVector.GetSize(); i++)
|
---|
2607 | {
|
---|
2608 |
|
---|
2609 | *fLog << fThetaRangeStringVector[SuccessfulThetaBinsVector[i]-1] << endl << endl;
|
---|
2610 |
|
---|
2611 | // 0) Initializing histograms Nex and Significance (with names!!)
|
---|
2612 |
|
---|
2613 | TString NexVSAlphaNameTmp = (NexVSAlphaName);
|
---|
2614 | TString SignificanceVSAlphaNameTmp = (SignificanceVSAlphaName);
|
---|
2615 |
|
---|
2616 | NexVSAlphaNameTmp += "Test";
|
---|
2617 | NexVSAlphaNameTmp += fThetaRangeStringVector[SuccessfulThetaBinsVector[i]-1];
|
---|
2618 |
|
---|
2619 | SignificanceVSAlphaNameTmp += "Test";
|
---|
2620 | SignificanceVSAlphaNameTmp += fThetaRangeStringVector[SuccessfulThetaBinsVector[i]-1];
|
---|
2621 |
|
---|
2622 |
|
---|
2623 | // TEMP
|
---|
2624 | // cout << NexVSAlphaNameTmp << endl;
|
---|
2625 | // cout << SignificanceVSAlphaNameTmp << endl;
|
---|
2626 | // ENDTEMP
|
---|
2627 |
|
---|
2628 |
|
---|
2629 |
|
---|
2630 |
|
---|
2631 | NexVSAlphaSigTestVector[i] = new TH1F (NexVSAlphaNameTmp.Data(),
|
---|
2632 | NexVSAlphaNameTmp.Data(),
|
---|
2633 | NexSignVSAlphaXBinsVectorDimension - 1,
|
---|
2634 | XBinsVector);
|
---|
2635 |
|
---|
2636 | // It is important to unlink the the created histograms from
|
---|
2637 | // the current directory (defined now by "rootfile"),
|
---|
2638 | // so that we can do whatever we want to
|
---|
2639 | // with them regardless of such directory is opened or closed
|
---|
2640 | // The reference of the histogram is removed from the current
|
---|
2641 | // directory by means of the member function "SetDirectory(0)"
|
---|
2642 |
|
---|
2643 |
|
---|
2644 | NexVSAlphaSigTestVector[i] -> SetDirectory(0);
|
---|
2645 |
|
---|
2646 |
|
---|
2647 | SignificanceVSAlphaSigTestVector[i] = new TH1F (SignificanceVSAlphaNameTmp.Data(),
|
---|
2648 | SignificanceVSAlphaNameTmp.Data(),
|
---|
2649 | NexSignVSAlphaXBinsVectorDimension - 1,
|
---|
2650 | XBinsVector);
|
---|
2651 |
|
---|
2652 | SignificanceVSAlphaSigTestVector[i] -> SetDirectory(0);
|
---|
2653 |
|
---|
2654 |
|
---|
2655 |
|
---|
2656 | // axis are defined
|
---|
2657 |
|
---|
2658 | NexVSAlphaSigTestVector[i] -> SetXTitle(x_axis_label.Data());
|
---|
2659 | NexVSAlphaSigTestVector[i] -> SetYTitle(y_axis_label_nex.Data());
|
---|
2660 |
|
---|
2661 | SignificanceVSAlphaSigTestVector[i] ->SetXTitle(x_axis_label.Data());
|
---|
2662 | SignificanceVSAlphaSigTestVector[i] ->SetYTitle(y_axis_label_significance.Data());
|
---|
2663 |
|
---|
2664 |
|
---|
2665 | // 1) Retrieving alpha histograms (ON-OFF) and normalization factor
|
---|
2666 |
|
---|
2667 | NormalizationFactor =
|
---|
2668 | NormFactorTestHisto -> GetBinContent(SuccessfulThetaBinsVector[i]);
|
---|
2669 |
|
---|
2670 |
|
---|
2671 | TString TmpAlphaHistoName = ("TestAlphaAfterCuts");
|
---|
2672 | TmpAlphaHistoName += fThetaRangeStringVector[SuccessfulThetaBinsVector[i]-1];
|
---|
2673 |
|
---|
2674 | *fLog << "Getting Histogram with name " << endl
|
---|
2675 | << TmpAlphaHistoName << endl;
|
---|
2676 |
|
---|
2677 | AlphaONAfterCuts = (TH1F*) rootfile.Get(TmpAlphaHistoName);
|
---|
2678 |
|
---|
2679 |
|
---|
2680 | TString TmpAlphaOFFHistoName = ("TestAlphaOFFAfterCuts");
|
---|
2681 | TmpAlphaOFFHistoName += fThetaRangeStringVector[SuccessfulThetaBinsVector[i]-1];
|
---|
2682 |
|
---|
2683 | *fLog << "Getting Histogram with name " << endl
|
---|
2684 | << TmpAlphaOFFHistoName << endl;
|
---|
2685 |
|
---|
2686 | AlphaOFFAfterCuts = (TH1F*) rootfile.Get(TmpAlphaOFFHistoName);
|
---|
2687 |
|
---|
2688 |
|
---|
2689 | // 2) Computing Nex and Signficance vs alpha and filling histograms
|
---|
2690 | // Significance is found using MHFindSignificanceONOFF::FindSigmaONOFF
|
---|
2691 |
|
---|
2692 |
|
---|
2693 | *fLog << endl << endl
|
---|
2694 | << "Starting loop to compute Nex and Significance vs alphasig " << endl
|
---|
2695 | << "for TRAIN sample and "
|
---|
2696 | << fThetaRangeStringVector[SuccessfulThetaBinsVector[i]-1]
|
---|
2697 | << " :" << endl;
|
---|
2698 |
|
---|
2699 |
|
---|
2700 | for (Int_t j = 1; j < NexSignVSAlphaXBinsVectorDimension; j++)
|
---|
2701 | {
|
---|
2702 |
|
---|
2703 |
|
---|
2704 | alphasig = XBinsVector[j];
|
---|
2705 |
|
---|
2706 | MHFindSignificanceONOFF findsig;
|
---|
2707 | findsig.SetRebin(RebinHistograms);
|
---|
2708 | findsig.SetReduceDegree(ReduceDegree);
|
---|
2709 | findsig.SetUseFittedQuantities(fUseFittedQuantities);
|
---|
2710 |
|
---|
2711 | // Dummy name for psfile
|
---|
2712 | TString psfilename = ("DummyPs");
|
---|
2713 |
|
---|
2714 | findsig.FindSigmaONOFF(AlphaONAfterCuts, AlphaOFFAfterCuts,
|
---|
2715 | NormalizationFactor,
|
---|
2716 | alphabkgmin, alphabkgmax,
|
---|
2717 | degree, degreeOFF,
|
---|
2718 | alphasig,
|
---|
2719 | drawpoly, fitgauss, print,
|
---|
2720 | saveplots, psfilename);
|
---|
2721 |
|
---|
2722 |
|
---|
2723 |
|
---|
2724 | Significance = double(findsig.GetSignificance());
|
---|
2725 |
|
---|
2726 | // so far there is not a variable that contains
|
---|
2727 | // NexONOFF or NexONOFFFitted (in class MHFindSignificanceONOFF)
|
---|
2728 | // according to the value of the variable fUseFittedQuantities.
|
---|
2729 | // As soon as I have some time I will implement it (with the
|
---|
2730 | // corresponging GET member function) and will change
|
---|
2731 | // MFindSignificanceONOFF accordingly to make use of it.
|
---|
2732 |
|
---|
2733 | // For the time being, I will survive with an "if statement".
|
---|
2734 |
|
---|
2735 |
|
---|
2736 | if(fUseFittedQuantities)
|
---|
2737 | {
|
---|
2738 | Nex = double(findsig.GetNexONOFFFitted());
|
---|
2739 | }
|
---|
2740 | else
|
---|
2741 | {
|
---|
2742 | Nex = double(findsig.GetNexONOFF());
|
---|
2743 | }
|
---|
2744 |
|
---|
2745 |
|
---|
2746 | *fLog << endl << "Cut in alpha, Nex and significance : "
|
---|
2747 | << alphasig << ", " << Nex << ", " << Significance << endl << endl;
|
---|
2748 |
|
---|
2749 |
|
---|
2750 | // updating histograms
|
---|
2751 |
|
---|
2752 | tmpdouble = alphasig - AlphaBinWidth/2.0;
|
---|
2753 |
|
---|
2754 | NexVSAlphaSigTestVector[i] -> Fill(tmpdouble, Nex);
|
---|
2755 | SignificanceVSAlphaSigTestVector[i] -> Fill(tmpdouble, Significance);
|
---|
2756 |
|
---|
2757 |
|
---|
2758 |
|
---|
2759 |
|
---|
2760 |
|
---|
2761 |
|
---|
2762 |
|
---|
2763 | }
|
---|
2764 |
|
---|
2765 |
|
---|
2766 | // Dynamic memory allocated in this loop is released
|
---|
2767 |
|
---|
2768 |
|
---|
2769 | *fLog << "Memory allocated by temporal alpha histograms "
|
---|
2770 | << "will be released" << endl;
|
---|
2771 |
|
---|
2772 | delete AlphaONAfterCuts;
|
---|
2773 | delete AlphaOFFAfterCuts;
|
---|
2774 |
|
---|
2775 | }
|
---|
2776 |
|
---|
2777 |
|
---|
2778 |
|
---|
2779 |
|
---|
2780 |
|
---|
2781 |
|
---|
2782 |
|
---|
2783 |
|
---|
2784 | }
|
---|
2785 |
|
---|
2786 | rootfile.Close();
|
---|
2787 |
|
---|
2788 |
|
---|
2789 |
|
---|
2790 | // *************************************************************
|
---|
2791 | // Histograms Nex and Significance vs alpha are written into
|
---|
2792 | // root file fAlphaDistributionsRootFilename
|
---|
2793 |
|
---|
2794 |
|
---|
2795 |
|
---|
2796 |
|
---|
2797 | *fLog <<"Histograms containing Nex and significance vs alphasig "
|
---|
2798 | << " are saved into root file " << endl
|
---|
2799 | << fAlphaDistributionsRootFilename << endl;
|
---|
2800 |
|
---|
2801 |
|
---|
2802 | TFile rootfiletowrite (fAlphaDistributionsRootFilename, "UPDATE",
|
---|
2803 | "Alpha Distributions for several Theta bins");
|
---|
2804 |
|
---|
2805 |
|
---|
2806 | for (Int_t i = 0; i < SuccessfulThetaBinsVector.GetSize(); i++)
|
---|
2807 | {
|
---|
2808 |
|
---|
2809 | if (fOptimizeParameters)
|
---|
2810 | {// Histograms for train sample are written to root file
|
---|
2811 |
|
---|
2812 | if(NexVSAlphaSigTrainVector[i])
|
---|
2813 | {
|
---|
2814 | if(NexVSAlphaSigTrainVector[i]-> GetEntries() > 0.5)
|
---|
2815 | NexVSAlphaSigTrainVector[i]-> Write();
|
---|
2816 |
|
---|
2817 | }
|
---|
2818 |
|
---|
2819 | if(SignificanceVSAlphaSigTrainVector[i])
|
---|
2820 | {
|
---|
2821 | if(SignificanceVSAlphaSigTrainVector[i]-> GetEntries() > 0.5)
|
---|
2822 | SignificanceVSAlphaSigTrainVector[i]-> Write();
|
---|
2823 |
|
---|
2824 | }
|
---|
2825 | }
|
---|
2826 |
|
---|
2827 |
|
---|
2828 |
|
---|
2829 | if (fTestParameters)
|
---|
2830 | {// Histograms for test sample are written to root file
|
---|
2831 |
|
---|
2832 | if(NexVSAlphaSigTestVector[i])
|
---|
2833 | {
|
---|
2834 | if(NexVSAlphaSigTestVector[i]-> GetEntries() > 0.5)
|
---|
2835 | NexVSAlphaSigTestVector[i]-> Write();
|
---|
2836 |
|
---|
2837 | }
|
---|
2838 |
|
---|
2839 | if(SignificanceVSAlphaSigTestVector[i])
|
---|
2840 | {
|
---|
2841 | if(SignificanceVSAlphaSigTestVector[i]-> GetEntries() > 0.5)
|
---|
2842 | SignificanceVSAlphaSigTestVector[i]-> Write();
|
---|
2843 |
|
---|
2844 | }
|
---|
2845 | }
|
---|
2846 |
|
---|
2847 |
|
---|
2848 |
|
---|
2849 | }
|
---|
2850 |
|
---|
2851 |
|
---|
2852 | rootfiletowrite.Close();
|
---|
2853 |
|
---|
2854 |
|
---|
2855 |
|
---|
2856 |
|
---|
2857 | // Dynamic memory allocated is released
|
---|
2858 |
|
---|
2859 | delete XBinsVector;
|
---|
2860 |
|
---|
2861 | // function finished successfully... hopefully...
|
---|
2862 | return kTRUE;
|
---|
2863 |
|
---|
2864 | }
|
---|
2865 |
|
---|
2866 |
|
---|
2867 |
|
---|
2868 | // Function that gets the histograms with the alpha distributions
|
---|
2869 | // (for all the theta bins specified by fCosThetaRangeVector)
|
---|
2870 | // stored in fAlphaDistributionsRootFilename, and combine them
|
---|
2871 | // (correcting OFF histograms with the normalization factors stored
|
---|
2872 | // in NormFactorTrain or NormFactorTest) to get one single
|
---|
2873 | // Alpha distribution for ON and another one for OFF.
|
---|
2874 | // Then these histograms are given as arguments to
|
---|
2875 | // the function MHFindSignificanceONOFF::FindSigmaONOFF,
|
---|
2876 | // (Object of this class is created) to compute the
|
---|
2877 | // Overall Excess events and significance, that will be
|
---|
2878 | // stored in variables fOverallNexTrain/Test and fOverallSigmaLiMaTrain/Test
|
---|
2879 |
|
---|
2880 |
|
---|
2881 |
|
---|
2882 |
|
---|
2883 | Bool_t MFindSupercutsONOFFThetaLoop::ComputeOverallSignificance(Bool_t CombineTrainData,
|
---|
2884 | Bool_t CombineTestData)
|
---|
2885 | {
|
---|
2886 |
|
---|
2887 | // First of all let's check that we have the proper the "drinks" for the party
|
---|
2888 |
|
---|
2889 | if (fAlphaDistributionsRootFilename.IsNull())
|
---|
2890 | {
|
---|
2891 | *fLog << "MFindSupercutsONOFFThetaLoop::ComputeOverallSignificance()" << endl
|
---|
2892 | << "Root file containing all alpha distributions "
|
---|
2893 | <<" (fAlphaDistributionsRootFilename) is not well defined." << endl
|
---|
2894 | << "Execution of the function is ABORTED" << endl;
|
---|
2895 | return kFALSE;
|
---|
2896 | }
|
---|
2897 |
|
---|
2898 |
|
---|
2899 | if (fCosThetaRangeVector.GetSize() < 2)
|
---|
2900 |
|
---|
2901 | {
|
---|
2902 | *fLog << "MFindSupercutsONOFFThetaLoop::ComputeOverallSignificance()" << endl
|
---|
2903 | << "Dimension of vector fCosThetaRangeVector is smaller than 2."
|
---|
2904 | << "So, vector fCosThetaRangeVector is NOT properly defined" << endl
|
---|
2905 | << "Function execution will be ABORTED" << endl;
|
---|
2906 |
|
---|
2907 | return kFALSE;
|
---|
2908 | }
|
---|
2909 |
|
---|
2910 |
|
---|
2911 | TArrayI SuccessfulThetaBinsVector;
|
---|
2912 | // Vector containing the histogram index of the successfully
|
---|
2913 | // optimized theta bins.
|
---|
2914 | // It is filled with the contents of the histogram
|
---|
2915 | // fSuccessfulThetaBinsHist retrieved from the root
|
---|
2916 | // file fAlphaDistributionsRootFilename
|
---|
2917 | // Remember that histogram index START AT 1 !!!
|
---|
2918 |
|
---|
2919 | TH1F* CombinedAlphaTrainON = NULL;
|
---|
2920 | TString CombinedAlphaTrainONName = ("CombinedAlphaTrainON");
|
---|
2921 | TH1F* CombinedAlphaTrainOFF = NULL;
|
---|
2922 | TString CombinedAlphaTrainOFFName = ("CombinedAlphaTrainOFF");
|
---|
2923 |
|
---|
2924 | TH1F* CombinedAlphaTestON = NULL;
|
---|
2925 | TString CombinedAlphaTestONName = ("CombinedAlphaTestON");
|
---|
2926 |
|
---|
2927 | TH1F* CombinedAlphaTestOFF = NULL;
|
---|
2928 | TString CombinedAlphaTestOFFName = ("CombinedAlphaTestOFF");
|
---|
2929 |
|
---|
2930 |
|
---|
2931 |
|
---|
2932 |
|
---|
2933 |
|
---|
2934 | Int_t NAlphaBins = 0; // Amount of bins of alpha histograms is expected to
|
---|
2935 | // be equal for all alpha histograms
|
---|
2936 |
|
---|
2937 | // Vectors storing the computed error for each of the NAlphaBins
|
---|
2938 | // Only 2 vectors are needed for OFF Train and OFF Test sample, whose
|
---|
2939 | // alpha distributions are corrected with normalization factors.
|
---|
2940 |
|
---|
2941 | // The ON sample is just added (no normalization factors are applied)
|
---|
2942 | // and therefore, the error of the bin is the assumed by default
|
---|
2943 | // (i.e Sqrt(BinContent)
|
---|
2944 |
|
---|
2945 | TArrayD ErrorInCombinedAlphaTrainOFF;
|
---|
2946 | TArrayD ErrorInCombinedAlphaTestOFF;
|
---|
2947 |
|
---|
2948 | TArrayD CombinedAlphaTrainOFFWithoutNormalization;
|
---|
2949 | TArrayD CombinedAlphaTestOFFWithoutNormalization;
|
---|
2950 |
|
---|
2951 |
|
---|
2952 |
|
---|
2953 | const Double_t SmallQuantity = 0.01; // This quantity will be used as the
|
---|
2954 | // maximum difference between quantities (double) that are expected to be equal,
|
---|
2955 | // as the center of the bins in the several alpha distributions
|
---|
2956 |
|
---|
2957 | TH1F* NormFactorTrainHisto;
|
---|
2958 | TString NormFactorTrainHistoName = ("NormFactorTrainHist");
|
---|
2959 | TH1F* NormFactorTestHisto;
|
---|
2960 | TString NormFactorTestHistoName = ("NormFactorTestHist");
|
---|
2961 | TString SuccessfulThetaBinsHistoName = ("SuccessfulThetaBinsHist");
|
---|
2962 |
|
---|
2963 |
|
---|
2964 | TH1F* SuccessfulThetaBinsHisto;
|
---|
2965 |
|
---|
2966 | TH1F* TmpTH1FHisto;
|
---|
2967 |
|
---|
2968 |
|
---|
2969 | // Retrieve number of theta bins that have to be combined
|
---|
2970 | // from TArrayD fCosThetaRangeVector
|
---|
2971 |
|
---|
2972 | Int_t tmpdim = fCosThetaRangeVector.GetSize() - 1;
|
---|
2973 | const Int_t NThetaBins = tmpdim;
|
---|
2974 |
|
---|
2975 | // Variables to store the mean normalization factor of the
|
---|
2976 | // combined OFF sample.
|
---|
2977 | // MeanNormFactor = Sum (NormFactor_i*Noff_i)/Sum (Noff_i)
|
---|
2978 |
|
---|
2979 |
|
---|
2980 | Double_t MeanNormFactorTrain = 0.0;
|
---|
2981 | Double_t MeanNormFactorTest = 0.0;
|
---|
2982 |
|
---|
2983 | Double_t TotalNumberOFFTrain = 0.0;
|
---|
2984 | Double_t TotalNumberOFFTest = 0.0;
|
---|
2985 |
|
---|
2986 |
|
---|
2987 |
|
---|
2988 | Double_t EventCounter = 0.0;// variable used to count events used to fill histograms
|
---|
2989 | // it is double (and not Int) because the contents of an histogram bin might not
|
---|
2990 | // be an integer value... due to normalization factors !!
|
---|
2991 | Int_t tmpint = 0;
|
---|
2992 | Double_t tmpdouble = 0.0; // tmp variable double
|
---|
2993 | Double_t tmpdouble2 = 0.0; // tmp variable double
|
---|
2994 | Double_t tmperror = 0.0; // quantity used to store
|
---|
2995 | // temporaly the bin error in OFF histogram
|
---|
2996 | Double_t tmpnormfactor = 0.0; // variable used to store normalization factor
|
---|
2997 | // of the theta bin i (used for TRAIN and TEST sample)
|
---|
2998 |
|
---|
2999 |
|
---|
3000 |
|
---|
3001 |
|
---|
3002 | // Open root file
|
---|
3003 |
|
---|
3004 |
|
---|
3005 |
|
---|
3006 |
|
---|
3007 | *fLog << "MFindSupercutsONOFFThetaLoop::ComputeOverallSignificance()" << endl
|
---|
3008 | << "Opening root file " << fAlphaDistributionsRootFilename << endl;
|
---|
3009 |
|
---|
3010 |
|
---|
3011 |
|
---|
3012 | TFile rootfile (fAlphaDistributionsRootFilename, "READ");
|
---|
3013 |
|
---|
3014 | // Retrieving SuccessfulThetaBinsHisto from root file.
|
---|
3015 |
|
---|
3016 |
|
---|
3017 | *fLog << "MFindSupercutsONOFFThetaLoop::ComputeOverallSignificance()" << endl
|
---|
3018 | << "Retrieving SuccessfulThetaBinsHisto from root file and filling vector "
|
---|
3019 | << " SuccessfulThetaBinsVector" << endl;
|
---|
3020 |
|
---|
3021 |
|
---|
3022 | cout << "Getting histogram " << SuccessfulThetaBinsHistoName << endl;
|
---|
3023 |
|
---|
3024 | SuccessfulThetaBinsHisto = (TH1F*) rootfile.Get(SuccessfulThetaBinsHistoName);
|
---|
3025 |
|
---|
3026 |
|
---|
3027 |
|
---|
3028 | // Check that bins in this histo correspond with the
|
---|
3029 | // CosTheta intervals defined by vector fCosThetaRangeVector
|
---|
3030 |
|
---|
3031 |
|
---|
3032 | if (SuccessfulThetaBinsHisto -> GetNbinsX() != NThetaBins)
|
---|
3033 | {
|
---|
3034 | *fLog << "MFindSupercutsONOFFThetaLoop::ComputeOverallSignificance()"
|
---|
3035 | << endl
|
---|
3036 | << "Number of theta bins defined by vector fCosThetaRangeVector ("
|
---|
3037 | << NThetaBins << ") does not match with the bins in histogram "
|
---|
3038 | << SuccessfulThetaBinsHistoName << "("
|
---|
3039 | << SuccessfulThetaBinsHisto -> GetNbinsX() << ")" << endl
|
---|
3040 | << "Function execution will be ABORTED." << endl;
|
---|
3041 | return kFALSE;
|
---|
3042 | }
|
---|
3043 |
|
---|
3044 |
|
---|
3045 | // Filling vector SuccessfulThetaBinsVector
|
---|
3046 |
|
---|
3047 | tmpint = 0;
|
---|
3048 | for (Int_t i = 0; i < NThetaBins; i++)
|
---|
3049 | {
|
---|
3050 | if ((SuccessfulThetaBinsHisto -> GetBinContent(i+1)) > 0.5)
|
---|
3051 | {
|
---|
3052 | tmpint++;
|
---|
3053 | }
|
---|
3054 | }
|
---|
3055 |
|
---|
3056 | SuccessfulThetaBinsVector.Set(tmpint);
|
---|
3057 |
|
---|
3058 | tmpint = 0;
|
---|
3059 | for (Int_t i = 0; i < NThetaBins; i++)
|
---|
3060 | {
|
---|
3061 | if ((SuccessfulThetaBinsHisto -> GetBinContent(i+1)) > 0.5)
|
---|
3062 | {
|
---|
3063 | if(tmpint < SuccessfulThetaBinsVector.GetSize())
|
---|
3064 | {
|
---|
3065 | SuccessfulThetaBinsVector[tmpint] = i+1;
|
---|
3066 | tmpint++;
|
---|
3067 | }
|
---|
3068 | else
|
---|
3069 | {
|
---|
3070 | *fLog << "MFindSupercutsONOFFThetaLoop::ComputeOverallSignificance()"
|
---|
3071 | << endl
|
---|
3072 | << "Problem when filling vector 'SuccessfulThetaBinsVector'"
|
---|
3073 | << endl
|
---|
3074 | << "Function execution will be ABORTED." << endl;
|
---|
3075 | return kFALSE;
|
---|
3076 | }
|
---|
3077 |
|
---|
3078 |
|
---|
3079 | }
|
---|
3080 | }
|
---|
3081 |
|
---|
3082 |
|
---|
3083 |
|
---|
3084 |
|
---|
3085 |
|
---|
3086 | // ********* HISTOGRAMS FOR TRAIN SAMPLE WILL BE FILLED *********///
|
---|
3087 |
|
---|
3088 | if (CombineTrainData)
|
---|
3089 | {
|
---|
3090 |
|
---|
3091 | // Retrieving NormFactorTrainHisto from root file.
|
---|
3092 |
|
---|
3093 | NormFactorTrainHisto = (TH1F*) rootfile.Get(NormFactorTrainHistoName);
|
---|
3094 | // NormFactorTrainHisto-> DrawCopy();
|
---|
3095 |
|
---|
3096 | // Check that bins in this histo correspond with the
|
---|
3097 | // CosTheta intervals defined by vector fCosThetaRangeVector
|
---|
3098 |
|
---|
3099 |
|
---|
3100 | if (NormFactorTrainHisto -> GetNbinsX() != NThetaBins)
|
---|
3101 | {
|
---|
3102 | *fLog << "MFindSupercutsONOFFThetaLoop::ComputeOverallSignificance()"
|
---|
3103 | << endl
|
---|
3104 | << "Number of theta bins defined by vector fCosThetaRangeVector ("
|
---|
3105 | << NThetaBins << ") does not match with the bins in histogram "
|
---|
3106 | << NormFactorTrainHistoName << "("
|
---|
3107 | << NormFactorTrainHisto -> GetNbinsX() << ")" << endl
|
---|
3108 | << "Function execution will be ABORTED." << endl;
|
---|
3109 | return kFALSE;
|
---|
3110 | }
|
---|
3111 |
|
---|
3112 | // histo test
|
---|
3113 | /*
|
---|
3114 | cout << "NORMFACTORTRAINHIST TEST: BinCenter and Value" << endl;
|
---|
3115 | for (Int_t k = 0; k < NThetaBins; k++)
|
---|
3116 | {
|
---|
3117 |
|
---|
3118 | cout << NormFactorTrainHisto -> GetBinCenter(k+1)
|
---|
3119 | << "; " << NormFactorTrainHisto -> GetBinContent(k+1)<< endl;
|
---|
3120 | }
|
---|
3121 | */
|
---|
3122 |
|
---|
3123 | for (Int_t i = 0; i < SuccessfulThetaBinsVector.GetSize(); i++)
|
---|
3124 | {
|
---|
3125 | tmpdouble = NormFactorTrainHisto -> GetBinCenter(SuccessfulThetaBinsVector[i]);
|
---|
3126 | if (fCosThetaRangeVector[SuccessfulThetaBinsVector[i]-1] > tmpdouble ||
|
---|
3127 | fCosThetaRangeVector[SuccessfulThetaBinsVector[i]] < tmpdouble)
|
---|
3128 | {
|
---|
3129 | *fLog << "MFindSupercutsONOFFThetaLoop::ComputeOverallSignificance()"
|
---|
3130 | << endl
|
---|
3131 | << "Bins defined by vector fCosThetaRangeVector "
|
---|
3132 | << "do not match with the bins of histogram "
|
---|
3133 | << NormFactorTrainHistoName << endl
|
---|
3134 | << "CosTheta: " << fCosThetaRangeVector[SuccessfulThetaBinsVector[i]-1]
|
---|
3135 | << "-" << fCosThetaRangeVector[SuccessfulThetaBinsVector[i]] << endl
|
---|
3136 | << "NormFactorHist Bin: " << tmpdouble << endl
|
---|
3137 | << "Function execution will be ABORTED." << endl;
|
---|
3138 |
|
---|
3139 | return kFALSE;
|
---|
3140 | }
|
---|
3141 |
|
---|
3142 | }
|
---|
3143 |
|
---|
3144 |
|
---|
3145 | // Let's summ up all ON alpha distributions
|
---|
3146 |
|
---|
3147 | EventCounter = 0;
|
---|
3148 |
|
---|
3149 | for (Int_t i = 0; i < SuccessfulThetaBinsVector.GetSize(); i++)
|
---|
3150 | {
|
---|
3151 | if (fThetaRangeStringVector[SuccessfulThetaBinsVector[i]-1].IsNull())
|
---|
3152 | {
|
---|
3153 | *fLog << "MFindSupercutsONOFFThetaLoop::ComputeOverallSignificance()"
|
---|
3154 | << endl
|
---|
3155 | << "Component " << SuccessfulThetaBinsVector[i]-1
|
---|
3156 | << " of fThetaRangeStringVector is EMPTY"
|
---|
3157 | << endl
|
---|
3158 | << "Function execution will be ABORTED." << endl;
|
---|
3159 | return kFALSE;
|
---|
3160 |
|
---|
3161 |
|
---|
3162 | }
|
---|
3163 |
|
---|
3164 |
|
---|
3165 |
|
---|
3166 | TString TmpHistoName = ("TrainAlphaAfterCuts");
|
---|
3167 | TmpHistoName += fThetaRangeStringVector[SuccessfulThetaBinsVector[i]-1];
|
---|
3168 |
|
---|
3169 | *fLog << "Getting Histogram with name " << endl
|
---|
3170 | << TmpHistoName << endl;
|
---|
3171 |
|
---|
3172 | if (i == 0)
|
---|
3173 | {
|
---|
3174 | CombinedAlphaTrainON = (TH1F*) rootfile.Get(TmpHistoName);
|
---|
3175 | NAlphaBins = CombinedAlphaTrainON -> GetNbinsX();
|
---|
3176 |
|
---|
3177 | // Name of histogram is set
|
---|
3178 | CombinedAlphaTrainON -> SetName(CombinedAlphaTrainONName);
|
---|
3179 |
|
---|
3180 | // update EventCounter
|
---|
3181 | EventCounter += CombinedAlphaTrainON -> GetEntries();
|
---|
3182 |
|
---|
3183 | // tmp
|
---|
3184 | cout << "NAlphaBins in histo CombinedAlphaTrainON: "
|
---|
3185 | << NAlphaBins << endl;
|
---|
3186 | // endtmp
|
---|
3187 | }
|
---|
3188 | else
|
---|
3189 | {
|
---|
3190 | TmpTH1FHisto = (TH1F*) rootfile.Get(TmpHistoName);
|
---|
3191 |
|
---|
3192 | // update EventCounter
|
---|
3193 | EventCounter += TmpTH1FHisto -> GetEntries();
|
---|
3194 |
|
---|
3195 | for (Int_t j = 1; j <= NAlphaBins; j++)
|
---|
3196 | {
|
---|
3197 | // At some time alpha bins might not be the same for
|
---|
3198 | // the different alpha distributions. That's why
|
---|
3199 | // I will check it before summing the alpha histo contents
|
---|
3200 | tmpdouble = CombinedAlphaTrainON->GetBinCenter(j) -
|
---|
3201 | TmpTH1FHisto->GetBinCenter(j);
|
---|
3202 | tmpdouble = TMath::Abs(tmpdouble);
|
---|
3203 |
|
---|
3204 |
|
---|
3205 | if (tmpdouble > SmallQuantity)
|
---|
3206 | {
|
---|
3207 | *fLog << "MFindSupercutsONOFFThetaLoop::ComputeOverallSignificance()"
|
---|
3208 | << endl
|
---|
3209 | << "Bins among the several alpha ON histograms "
|
---|
3210 | << "for TRAIN sample do not match"
|
---|
3211 | << endl
|
---|
3212 | << "Function execution will be ABORTED." << endl;
|
---|
3213 | return kFALSE;
|
---|
3214 | }
|
---|
3215 |
|
---|
3216 | tmpdouble = CombinedAlphaTrainON->GetBinContent(j) +
|
---|
3217 | TmpTH1FHisto->GetBinContent(j);
|
---|
3218 |
|
---|
3219 |
|
---|
3220 | CombinedAlphaTrainON-> SetBinContent(j,tmpdouble);
|
---|
3221 |
|
---|
3222 | }
|
---|
3223 |
|
---|
3224 | // Dynamic memory allocated for TmpTH1FHisto is released
|
---|
3225 |
|
---|
3226 |
|
---|
3227 | //tmp
|
---|
3228 | cout << "Memory allocated by temporal histogram TmpTH1FHisto "
|
---|
3229 | << "will be released" << endl;
|
---|
3230 | //endtmp
|
---|
3231 | delete TmpTH1FHisto;
|
---|
3232 |
|
---|
3233 | }
|
---|
3234 |
|
---|
3235 |
|
---|
3236 | }
|
---|
3237 |
|
---|
3238 | // Entries of histogram CombinedAlphaTrainON is set to
|
---|
3239 | // the events counted in EventCounter
|
---|
3240 |
|
---|
3241 | tmpint = int (EventCounter);
|
---|
3242 | CombinedAlphaTrainON->SetEntries(tmpint);
|
---|
3243 |
|
---|
3244 | // Let's summ up all OFF alpha distributions
|
---|
3245 |
|
---|
3246 | EventCounter = 0;
|
---|
3247 |
|
---|
3248 | for (Int_t i = 0; i < SuccessfulThetaBinsVector.GetSize(); i++)
|
---|
3249 | {
|
---|
3250 | if (fThetaRangeStringVector[SuccessfulThetaBinsVector[i]-1].IsNull())
|
---|
3251 | {
|
---|
3252 | *fLog << "MFindSupercutsONOFFThetaLoop::ComputeOverallSignificance()"
|
---|
3253 | << endl
|
---|
3254 | << "Component " << SuccessfulThetaBinsVector[i]-1
|
---|
3255 | << " of fThetaRangeStringVector is EMPTY"
|
---|
3256 | << endl
|
---|
3257 | << "Function execution will be ABORTED." << endl;
|
---|
3258 | return kFALSE;
|
---|
3259 |
|
---|
3260 |
|
---|
3261 | }
|
---|
3262 |
|
---|
3263 | TString TmpHistoName = ("TrainAlphaOFFAfterCuts");
|
---|
3264 | TmpHistoName += fThetaRangeStringVector[SuccessfulThetaBinsVector[i]-1];
|
---|
3265 |
|
---|
3266 | // Normalization constant for this theta bin is stored now in tmpnormfactor
|
---|
3267 | tmpnormfactor = NormFactorTrainHisto -> GetBinContent(SuccessfulThetaBinsVector[i]);
|
---|
3268 |
|
---|
3269 | // temp
|
---|
3270 | cout << "Normalization factor for bin "
|
---|
3271 | << SuccessfulThetaBinsVector[i] << " is "
|
---|
3272 | << tmpnormfactor << endl;
|
---|
3273 |
|
---|
3274 | // endtemp
|
---|
3275 |
|
---|
3276 | if (i == 0)
|
---|
3277 | {
|
---|
3278 | CombinedAlphaTrainOFF = (TH1F*) rootfile.Get(TmpHistoName);
|
---|
3279 | NAlphaBins = CombinedAlphaTrainOFF -> GetNbinsX();
|
---|
3280 |
|
---|
3281 |
|
---|
3282 | // Name of histogram is set
|
---|
3283 | CombinedAlphaTrainOFF -> SetName(CombinedAlphaTrainOFFName);
|
---|
3284 |
|
---|
3285 | // Event counter is updated taking into account the
|
---|
3286 | // normalization factor
|
---|
3287 |
|
---|
3288 | EventCounter += tmpnormfactor * CombinedAlphaTrainOFF -> GetEntries();
|
---|
3289 |
|
---|
3290 | // Contribution to total number of OFF events and mean normfactor
|
---|
3291 | // of this theta bin is computed
|
---|
3292 |
|
---|
3293 | TotalNumberOFFTrain = CombinedAlphaTrainOFF -> GetEntries();
|
---|
3294 | MeanNormFactorTrain = tmpnormfactor *
|
---|
3295 | CombinedAlphaTrainOFF -> GetEntries();
|
---|
3296 |
|
---|
3297 |
|
---|
3298 |
|
---|
3299 | // ****
|
---|
3300 | // Dimension of ErrorInCombinedAlphaTrainOFF is set
|
---|
3301 | ErrorInCombinedAlphaTrainOFF.Set(NAlphaBins);
|
---|
3302 |
|
---|
3303 | CombinedAlphaTrainOFFWithoutNormalization.Set(NAlphaBins);
|
---|
3304 |
|
---|
3305 | // Histogram now is normalized (to the respective ON data)
|
---|
3306 | // by using the normalization constants stored in
|
---|
3307 | // histogram NormFactorTrainHisto
|
---|
3308 |
|
---|
3309 | // Error of the normalized bins is computed and stored in
|
---|
3310 | // vector ErrorInCombinedAlphaTrainOFF
|
---|
3311 |
|
---|
3312 | for (Int_t j = 1; j <= NAlphaBins; j++)
|
---|
3313 | {
|
---|
3314 |
|
---|
3315 | // Number of events (without normalization) are
|
---|
3316 | // stored in vector CombinedAlphaTrainOFFWithoutNormalization
|
---|
3317 |
|
---|
3318 | CombinedAlphaTrainOFFWithoutNormalization[j-1] =
|
---|
3319 | CombinedAlphaTrainOFF -> GetBinContent(j);
|
---|
3320 |
|
---|
3321 | // Bin content is set
|
---|
3322 | tmpdouble2 = tmpnormfactor *
|
---|
3323 | CombinedAlphaTrainOFF -> GetBinContent(j);
|
---|
3324 |
|
---|
3325 | CombinedAlphaTrainOFF -> SetBinContent(j,tmpdouble2);
|
---|
3326 |
|
---|
3327 | // (Bin Error)^2 is computed and stored in
|
---|
3328 | // ErrorInCombinedAlphaTrainOFF
|
---|
3329 | // Bin error = Sqrt(BinContent) X Normalization Factor
|
---|
3330 |
|
---|
3331 | tmperror = tmpdouble2*tmpnormfactor;
|
---|
3332 | ErrorInCombinedAlphaTrainOFF[j-1] = tmperror;
|
---|
3333 |
|
---|
3334 | }
|
---|
3335 |
|
---|
3336 | }
|
---|
3337 | else
|
---|
3338 | {
|
---|
3339 | TmpTH1FHisto = (TH1F*) rootfile.Get(TmpHistoName);
|
---|
3340 |
|
---|
3341 |
|
---|
3342 | // Event counter is updated taking into account the
|
---|
3343 | // normalization factor
|
---|
3344 |
|
---|
3345 | EventCounter += tmpnormfactor * TmpTH1FHisto -> GetEntries();
|
---|
3346 |
|
---|
3347 | // Contribution to total number of OFF events and mean normfactor
|
---|
3348 | // of this theta bin is computed
|
---|
3349 |
|
---|
3350 | TotalNumberOFFTrain += TmpTH1FHisto -> GetEntries();
|
---|
3351 | MeanNormFactorTrain += tmpnormfactor *
|
---|
3352 | TmpTH1FHisto -> GetEntries();
|
---|
3353 |
|
---|
3354 |
|
---|
3355 |
|
---|
3356 | for (Int_t j = 1; j <= NAlphaBins; j++)
|
---|
3357 | {
|
---|
3358 | // At some time alpha bins might not be the same for
|
---|
3359 | // the different alpha distributions. That's why
|
---|
3360 | // I will check it before summing the alpha histo contents
|
---|
3361 | tmpdouble2 = CombinedAlphaTrainOFF->GetBinCenter(j) -
|
---|
3362 | TmpTH1FHisto->GetBinCenter(j);
|
---|
3363 | tmpdouble2 = TMath::Abs(tmpdouble2);
|
---|
3364 |
|
---|
3365 | if (tmpdouble2 > SmallQuantity)
|
---|
3366 | {
|
---|
3367 | *fLog << "MFindSupercutsONOFFThetaLoop::ComputeOverallSignificance()"
|
---|
3368 | << endl
|
---|
3369 | << "Bins among the several alpha OFF histograms "
|
---|
3370 | << "for TRAIN sample do not match"
|
---|
3371 | << endl
|
---|
3372 | << "Function execution will be ABORTED." << endl;
|
---|
3373 | return kFALSE;
|
---|
3374 | }
|
---|
3375 |
|
---|
3376 |
|
---|
3377 | // Number of events (without normalization) are
|
---|
3378 | // added in vector CombinedAlphaTrainOFFWithoutNormalization
|
---|
3379 |
|
---|
3380 | CombinedAlphaTrainOFFWithoutNormalization[j-1] +=
|
---|
3381 | TmpTH1FHisto -> GetBinContent(j);
|
---|
3382 |
|
---|
3383 | // Histogram bin contents must be normalized (to the respective ON data)
|
---|
3384 | // by using the normalization constants stored in
|
---|
3385 | // histogram NormFactorTrainHisto before being added
|
---|
3386 | // to CombinedAlphaTrainOFF
|
---|
3387 |
|
---|
3388 |
|
---|
3389 | tmpdouble2 = CombinedAlphaTrainOFF->GetBinContent(j) +
|
---|
3390 | (TmpTH1FHisto->GetBinContent(j) * tmpnormfactor);
|
---|
3391 |
|
---|
3392 |
|
---|
3393 | CombinedAlphaTrainOFF-> SetBinContent(j,tmpdouble2);
|
---|
3394 |
|
---|
3395 |
|
---|
3396 | // (Bin Error)^2 is computed and added to
|
---|
3397 | // ErrorInCombinedAlphaTrainOFF
|
---|
3398 | // Bin error = Sqrt(BinContent) X Normalization Factor
|
---|
3399 |
|
---|
3400 | tmperror = TmpTH1FHisto->GetBinContent(j) * tmpnormfactor;
|
---|
3401 | ErrorInCombinedAlphaTrainOFF[j-1] += tmperror;
|
---|
3402 |
|
---|
3403 | }
|
---|
3404 |
|
---|
3405 | // Dynamic memory allocated for TmpTH1FHisto is released
|
---|
3406 |
|
---|
3407 | delete TmpTH1FHisto;
|
---|
3408 |
|
---|
3409 | }
|
---|
3410 |
|
---|
3411 |
|
---|
3412 |
|
---|
3413 |
|
---|
3414 |
|
---|
3415 | }
|
---|
3416 |
|
---|
3417 |
|
---|
3418 |
|
---|
3419 |
|
---|
3420 | // Mean Normalization factor is computed for the combined OFF histogram.
|
---|
3421 | // This factor will be used when computing the significance via
|
---|
3422 | // MHFindSignificance::FindSigmaONOFF().
|
---|
3423 |
|
---|
3424 | // The bin contents (and errors) of teh combined OFF TRAIN histo will be
|
---|
3425 | // normalized (divided) by teh mean norm constant, so that
|
---|
3426 | // the histogram for OFF data given as argument in FindSigmaONOFF
|
---|
3427 | // is equivalent to a non-normalized histogram.
|
---|
3428 |
|
---|
3429 | MeanNormFactorTrain = MeanNormFactorTrain/TotalNumberOFFTrain;
|
---|
3430 |
|
---|
3431 |
|
---|
3432 |
|
---|
3433 | // Sqrt(Contents) is performed in ErrorInCombinedAlphaTrainOFF,
|
---|
3434 | // in order to compute the real error in the bin content of histogram
|
---|
3435 | // CombinedAlphaTrainOFF.
|
---|
3436 |
|
---|
3437 | // Then error and contents (normalized with mean normfctor) are
|
---|
3438 | // set to histogram bin
|
---|
3439 |
|
---|
3440 | for (Int_t j = 1; j <= NAlphaBins; j++)
|
---|
3441 | {
|
---|
3442 | tmpdouble2 =
|
---|
3443 | (CombinedAlphaTrainOFF -> GetBinContent(j))/MeanNormFactorTrain;
|
---|
3444 |
|
---|
3445 | CombinedAlphaTrainOFF -> SetBinContent(j,tmpdouble2);
|
---|
3446 |
|
---|
3447 | ErrorInCombinedAlphaTrainOFF[j-1] =
|
---|
3448 | TMath::Sqrt(ErrorInCombinedAlphaTrainOFF[j-1])/MeanNormFactorTrain;
|
---|
3449 | // Proper error treatment when adding quantities to histogram bin
|
---|
3450 |
|
---|
3451 | CombinedAlphaTrainOFF -> SetBinError(j,ErrorInCombinedAlphaTrainOFF[j-1]);
|
---|
3452 |
|
---|
3453 | // ****** TMP **************
|
---|
3454 | /*
|
---|
3455 | // Error computation assuming a given (no normalized) bin content
|
---|
3456 |
|
---|
3457 | tmpdouble = CombinedAlphaTrainOFF -> GetBinContent(j);
|
---|
3458 |
|
---|
3459 | if (CombinedAlphaTrainOFFWithoutNormalization[j-1] < SmallQuantity)
|
---|
3460 | {tmpnormfactor = 0;}
|
---|
3461 | else
|
---|
3462 | {tmpnormfactor = tmpdouble/CombinedAlphaTrainOFFWithoutNormalization[j-1];}
|
---|
3463 |
|
---|
3464 | tmperror = TMath::Sqrt(CombinedAlphaTrainOFFWithoutNormalization[j-1]) *
|
---|
3465 | tmpnormfactor;
|
---|
3466 |
|
---|
3467 | CombinedAlphaTrainOFF -> SetBinError(j,tmperror);
|
---|
3468 |
|
---|
3469 | */
|
---|
3470 | // ****** END TMP **************
|
---|
3471 |
|
---|
3472 | // tmp
|
---|
3473 | /*
|
---|
3474 | cout << CombinedAlphaTrainOFF -> GetBinContent(j) << " +/- "
|
---|
3475 | << CombinedAlphaTrainOFF -> GetBinError(j)
|
---|
3476 | << " Number Events without normalization: "
|
---|
3477 | << CombinedAlphaTrainOFFWithoutNormalization[j-1] << endl;
|
---|
3478 | */
|
---|
3479 | // endtmp
|
---|
3480 |
|
---|
3481 | }
|
---|
3482 |
|
---|
3483 |
|
---|
3484 | // Entries of histogram CombinedAlphaTrainON is set to
|
---|
3485 | // the events counted during histogram filling
|
---|
3486 |
|
---|
3487 | EventCounter = EventCounter/MeanNormFactorTrain;
|
---|
3488 | tmpint = int (EventCounter);
|
---|
3489 | CombinedAlphaTrainOFF->SetEntries(tmpint);
|
---|
3490 |
|
---|
3491 |
|
---|
3492 | cout << "SILLY INFO FOR TRAIN SAMPLE: Total number of events Before nomralization and "
|
---|
3493 | << " re-re normalized" << endl
|
---|
3494 | << TotalNumberOFFTrain << " - " << EventCounter << endl;
|
---|
3495 |
|
---|
3496 | }
|
---|
3497 |
|
---|
3498 | // ********* COMBINED HISTOGRAMS FOR TRAIN SAMPLE ALREADY FILLED ********* ///
|
---|
3499 |
|
---|
3500 |
|
---|
3501 |
|
---|
3502 |
|
---|
3503 | // ********* HISTOGRAMS FOR TEST SAMPLE WILL BE FILLED ********* ///
|
---|
3504 |
|
---|
3505 |
|
---|
3506 | if (CombineTestData)
|
---|
3507 | {
|
---|
3508 |
|
---|
3509 | // Retrieving NormFactorTestHisto from root file.
|
---|
3510 |
|
---|
3511 | NormFactorTestHisto = (TH1F*) rootfile.Get(NormFactorTestHistoName);
|
---|
3512 | // NormFactorTestHisto-> DrawCopy();
|
---|
3513 |
|
---|
3514 | // Check that bins in this histo correspond with the
|
---|
3515 | // CosTheta intervals defined by vector fCosThetaRangeVector
|
---|
3516 |
|
---|
3517 |
|
---|
3518 | if (NormFactorTestHisto -> GetNbinsX() != NThetaBins)
|
---|
3519 | {
|
---|
3520 | *fLog << "MFindSupercutsONOFFThetaLoop::ComputeOverallSignificance()"
|
---|
3521 | << endl
|
---|
3522 | << "Number of theta bins defined by vector fCosThetaRangeVector ("
|
---|
3523 | << NThetaBins << ") does not match with the bins in histogram "
|
---|
3524 | << NormFactorTestHistoName << "("
|
---|
3525 | << NormFactorTestHisto -> GetNbinsX() << ")" << endl
|
---|
3526 | << "Function execution will be ABORTED." << endl;
|
---|
3527 | return kFALSE;
|
---|
3528 | }
|
---|
3529 |
|
---|
3530 | // histo test
|
---|
3531 | /*
|
---|
3532 | cout << "NORMFACTORTRAINHIST TEST: BinCenter and Value" << endl;
|
---|
3533 | for (Int_t k = 0; k < NThetaBins; k++)
|
---|
3534 | {
|
---|
3535 |
|
---|
3536 | cout << NormFactorTestHisto -> GetBinCenter(k+1)
|
---|
3537 | << "; " << NormFactorTestHisto -> GetBinContent(k+1)<< endl;
|
---|
3538 | }
|
---|
3539 | */
|
---|
3540 |
|
---|
3541 |
|
---|
3542 |
|
---|
3543 | for (Int_t i = 0; i < SuccessfulThetaBinsVector.GetSize(); i++)
|
---|
3544 | {
|
---|
3545 | tmpdouble = NormFactorTestHisto -> GetBinCenter(SuccessfulThetaBinsVector[i]);
|
---|
3546 | if (fCosThetaRangeVector[SuccessfulThetaBinsVector[i]-1] > tmpdouble ||
|
---|
3547 | fCosThetaRangeVector[SuccessfulThetaBinsVector[i]] < tmpdouble)
|
---|
3548 | {
|
---|
3549 | *fLog << "MFindSupercutsONOFFThetaLoop::ComputeOverallSignificance()"
|
---|
3550 | << endl
|
---|
3551 | << "Bins defined by vector fCosThetaRangeVector "
|
---|
3552 | << "do not match with the bins of histogram "
|
---|
3553 | << NormFactorTestHistoName << endl
|
---|
3554 | << "CosTheta: " << fCosThetaRangeVector[SuccessfulThetaBinsVector[i]-1]
|
---|
3555 | << "-" << fCosThetaRangeVector[SuccessfulThetaBinsVector[i]] << endl
|
---|
3556 | << "NormFactorHist Bin: " << tmpdouble << endl
|
---|
3557 | << "Function execution will be ABORTED." << endl;
|
---|
3558 |
|
---|
3559 | return kFALSE;
|
---|
3560 | }
|
---|
3561 |
|
---|
3562 | }
|
---|
3563 |
|
---|
3564 |
|
---|
3565 | // Let's summ up all ON alpha distributions
|
---|
3566 |
|
---|
3567 | // Event counter (counts events used in the histogram filling) is
|
---|
3568 | // initialized to zero
|
---|
3569 |
|
---|
3570 | EventCounter = 0;
|
---|
3571 |
|
---|
3572 | for (Int_t i = 0; i < SuccessfulThetaBinsVector.GetSize(); i++)
|
---|
3573 | {
|
---|
3574 | if (fThetaRangeStringVector[SuccessfulThetaBinsVector[i]-1].IsNull())
|
---|
3575 | {
|
---|
3576 | *fLog << "MFindSupercutsONOFFThetaLoop::ComputeOverallSignificance()"
|
---|
3577 | << endl
|
---|
3578 | << "Component " << SuccessfulThetaBinsVector[i]-1
|
---|
3579 | << " of fThetaRangeStringVector is EMPTY"
|
---|
3580 | << endl
|
---|
3581 | << "Function execution will be ABORTED." << endl;
|
---|
3582 | return kFALSE;
|
---|
3583 |
|
---|
3584 |
|
---|
3585 | }
|
---|
3586 |
|
---|
3587 |
|
---|
3588 |
|
---|
3589 | TString TmpHistoName = ("TestAlphaAfterCuts");
|
---|
3590 | TmpHistoName += fThetaRangeStringVector[SuccessfulThetaBinsVector[i]-1];
|
---|
3591 |
|
---|
3592 | *fLog << "Getting Histogram with name " << endl
|
---|
3593 | << TmpHistoName << endl;
|
---|
3594 |
|
---|
3595 | if (i == 0)
|
---|
3596 | {
|
---|
3597 | CombinedAlphaTestON = (TH1F*) rootfile.Get(TmpHistoName);
|
---|
3598 | NAlphaBins = CombinedAlphaTestON -> GetNbinsX();
|
---|
3599 |
|
---|
3600 | // Name of histogram is set
|
---|
3601 | CombinedAlphaTestON -> SetName(CombinedAlphaTestONName);
|
---|
3602 |
|
---|
3603 | // EventCounter is updated
|
---|
3604 |
|
---|
3605 | EventCounter += CombinedAlphaTestON -> GetEntries();
|
---|
3606 |
|
---|
3607 | // tmp
|
---|
3608 | cout << "NAlphaBins in histo CombinedAlphaTestON: "
|
---|
3609 | << NAlphaBins << endl;
|
---|
3610 | // endtmp
|
---|
3611 | }
|
---|
3612 | else
|
---|
3613 | {
|
---|
3614 | TmpTH1FHisto = (TH1F*) rootfile.Get(TmpHistoName);
|
---|
3615 |
|
---|
3616 | // EventCounter is updated
|
---|
3617 |
|
---|
3618 | EventCounter += TmpTH1FHisto -> GetEntries();
|
---|
3619 |
|
---|
3620 | for (Int_t j = 1; j <= NAlphaBins; j++)
|
---|
3621 | {
|
---|
3622 | // At some time alpha bins might not be the same for
|
---|
3623 | // the different alpha distributions. That's why
|
---|
3624 | // I will check it before summing the alpha histo contents
|
---|
3625 | tmpdouble = CombinedAlphaTestON->GetBinCenter(j) -
|
---|
3626 | TmpTH1FHisto->GetBinCenter(j);
|
---|
3627 | tmpdouble = TMath::Abs(tmpdouble);
|
---|
3628 |
|
---|
3629 |
|
---|
3630 | if (tmpdouble > SmallQuantity)
|
---|
3631 | {
|
---|
3632 | *fLog << "MFindSupercutsONOFFThetaLoop::ComputeOverallSignificance()"
|
---|
3633 | << endl
|
---|
3634 | << "Bins among the several alpha ON histograms "
|
---|
3635 | << "for TEST sample do not match"
|
---|
3636 | << endl
|
---|
3637 | << "Function execution will be ABORTED." << endl;
|
---|
3638 | return kFALSE;
|
---|
3639 | }
|
---|
3640 |
|
---|
3641 | tmpdouble = CombinedAlphaTestON->GetBinContent(j) +
|
---|
3642 | TmpTH1FHisto->GetBinContent(j);
|
---|
3643 |
|
---|
3644 |
|
---|
3645 | CombinedAlphaTestON-> SetBinContent(j,tmpdouble);
|
---|
3646 |
|
---|
3647 | }
|
---|
3648 |
|
---|
3649 | // Dynamic memory allocated for TmpTH1FHisto is released
|
---|
3650 |
|
---|
3651 |
|
---|
3652 | //tmp
|
---|
3653 | cout << "Memory allocated by temporal histogram TmpTH1FHisto "
|
---|
3654 | << "will be released" << endl;
|
---|
3655 | //endtmp
|
---|
3656 | delete TmpTH1FHisto;
|
---|
3657 |
|
---|
3658 | }
|
---|
3659 |
|
---|
3660 |
|
---|
3661 |
|
---|
3662 |
|
---|
3663 | }
|
---|
3664 |
|
---|
3665 | // Entries of histogram CombinedAlphaTestON is set to
|
---|
3666 | // the events counted during histogram filling
|
---|
3667 | tmpint = int (EventCounter);
|
---|
3668 | CombinedAlphaTestON->SetEntries(tmpint);
|
---|
3669 |
|
---|
3670 |
|
---|
3671 | // Let's summ up all OFF alpha distributions
|
---|
3672 |
|
---|
3673 | EventCounter = 0;
|
---|
3674 |
|
---|
3675 | for (Int_t i = 0; i < SuccessfulThetaBinsVector.GetSize(); i++)
|
---|
3676 | {
|
---|
3677 | if (fThetaRangeStringVector[SuccessfulThetaBinsVector[i]-1].IsNull())
|
---|
3678 | {
|
---|
3679 | *fLog << "MFindSupercutsONOFFThetaLoop::ComputeOverallSignificance()"
|
---|
3680 | << endl
|
---|
3681 | << "Component " << SuccessfulThetaBinsVector[i]-1
|
---|
3682 | << " of fThetaRangeStringVector is EMPTY"
|
---|
3683 | << endl
|
---|
3684 | << "Function execution will be ABORTED." << endl;
|
---|
3685 | return kFALSE;
|
---|
3686 |
|
---|
3687 |
|
---|
3688 | }
|
---|
3689 |
|
---|
3690 | TString TmpHistoName = ("TestAlphaOFFAfterCuts");
|
---|
3691 | TmpHistoName += fThetaRangeStringVector[SuccessfulThetaBinsVector[i]-1];
|
---|
3692 |
|
---|
3693 | // Normalization constant for this theta bin is stored now in tmpnormfactor
|
---|
3694 | tmpnormfactor = NormFactorTestHisto -> GetBinContent(SuccessfulThetaBinsVector[i]);
|
---|
3695 |
|
---|
3696 | // temp
|
---|
3697 | cout << "Normalization factor for bin "
|
---|
3698 | << SuccessfulThetaBinsVector[i] << " is "
|
---|
3699 | << tmpnormfactor << endl;
|
---|
3700 |
|
---|
3701 | // endtemp
|
---|
3702 |
|
---|
3703 | if (i == 0)
|
---|
3704 | {
|
---|
3705 | CombinedAlphaTestOFF = (TH1F*) rootfile.Get(TmpHistoName);
|
---|
3706 | NAlphaBins = CombinedAlphaTestOFF -> GetNbinsX();
|
---|
3707 |
|
---|
3708 |
|
---|
3709 | // Name of histogram is set
|
---|
3710 | CombinedAlphaTestOFF -> SetName(CombinedAlphaTestOFFName);
|
---|
3711 |
|
---|
3712 | // EventCounter is updated taking into account
|
---|
3713 | // normalization factor
|
---|
3714 |
|
---|
3715 | EventCounter += tmpnormfactor*CombinedAlphaTestOFF -> GetEntries();
|
---|
3716 |
|
---|
3717 |
|
---|
3718 | // Contribution to total number of OFF events and mean normfactor
|
---|
3719 | // of this theta bin is computed
|
---|
3720 |
|
---|
3721 | TotalNumberOFFTest = CombinedAlphaTestOFF -> GetEntries();
|
---|
3722 | MeanNormFactorTest = tmpnormfactor *
|
---|
3723 | CombinedAlphaTestOFF -> GetEntries();
|
---|
3724 |
|
---|
3725 |
|
---|
3726 | // Dimension of ErrorInCombinedAlphaTrainOFF is set
|
---|
3727 | ErrorInCombinedAlphaTestOFF.Set(NAlphaBins);
|
---|
3728 |
|
---|
3729 | CombinedAlphaTestOFFWithoutNormalization.Set(NAlphaBins);
|
---|
3730 |
|
---|
3731 |
|
---|
3732 | // Histogram now is normalized (to the respective ON data)
|
---|
3733 | // by using the normalization constants stored in
|
---|
3734 | // histogram NormFactorTestHisto
|
---|
3735 |
|
---|
3736 |
|
---|
3737 | for (Int_t j = 1; j <= NAlphaBins; j++)
|
---|
3738 | {
|
---|
3739 |
|
---|
3740 | // Vector containing number of events without
|
---|
3741 | // normalization is filled
|
---|
3742 |
|
---|
3743 | CombinedAlphaTestOFFWithoutNormalization[j-1] =
|
---|
3744 | CombinedAlphaTestOFF -> GetBinContent(j);
|
---|
3745 |
|
---|
3746 |
|
---|
3747 |
|
---|
3748 | tmpdouble2 = tmpnormfactor *
|
---|
3749 | CombinedAlphaTestOFF -> GetBinContent(j);
|
---|
3750 |
|
---|
3751 | CombinedAlphaTestOFF -> SetBinContent(j,tmpdouble2);
|
---|
3752 |
|
---|
3753 | // (Bin Error)^2 is computed and stored in
|
---|
3754 | // ErrorInCombinedAlphaTestOFF
|
---|
3755 | // Bin error = Sqrt(BinContent) X Normalization Factor
|
---|
3756 |
|
---|
3757 | tmperror = tmpdouble2*tmpnormfactor;
|
---|
3758 | ErrorInCombinedAlphaTestOFF[j-1] = tmperror;
|
---|
3759 |
|
---|
3760 |
|
---|
3761 | }
|
---|
3762 |
|
---|
3763 | }
|
---|
3764 | else
|
---|
3765 | {
|
---|
3766 | TmpTH1FHisto = (TH1F*) rootfile.Get(TmpHistoName);
|
---|
3767 |
|
---|
3768 |
|
---|
3769 | // EventCounter is updated taking into account
|
---|
3770 | // normalization factor
|
---|
3771 |
|
---|
3772 | EventCounter += tmpnormfactor* TmpTH1FHisto-> GetEntries();
|
---|
3773 |
|
---|
3774 | // Contribution to total number of OFF events and mean normfactor
|
---|
3775 | // of this theta bin is computed
|
---|
3776 |
|
---|
3777 |
|
---|
3778 | TotalNumberOFFTest += TmpTH1FHisto -> GetEntries();
|
---|
3779 | MeanNormFactorTest += tmpnormfactor *
|
---|
3780 | TmpTH1FHisto -> GetEntries();
|
---|
3781 |
|
---|
3782 |
|
---|
3783 |
|
---|
3784 | for (Int_t j = 1; j <= NAlphaBins; j++)
|
---|
3785 | {
|
---|
3786 | // At some time alpha bins might not be the same for
|
---|
3787 | // the different alpha distributions. That's why
|
---|
3788 | // I will check it before summing the alpha histo contents
|
---|
3789 | tmpdouble2 = CombinedAlphaTestOFF->GetBinCenter(j) -
|
---|
3790 | TmpTH1FHisto->GetBinCenter(j);
|
---|
3791 | tmpdouble2 = TMath::Abs(tmpdouble2);
|
---|
3792 |
|
---|
3793 | if (tmpdouble2 > SmallQuantity)
|
---|
3794 | {
|
---|
3795 | *fLog << "MFindSupercutsONOFFThetaLoop::ComputeOverallSignificance()"
|
---|
3796 | << endl
|
---|
3797 | << "Bins among the several alpha OFF histograms "
|
---|
3798 | << "for TRAIN sample do not match"
|
---|
3799 | << endl
|
---|
3800 | << "Function execution will be ABORTED." << endl;
|
---|
3801 | return kFALSE;
|
---|
3802 | }
|
---|
3803 |
|
---|
3804 |
|
---|
3805 | // Vector containing number of events without
|
---|
3806 | // normalization is updated
|
---|
3807 |
|
---|
3808 | // Vector containing number of events without
|
---|
3809 | // normalization is filled
|
---|
3810 |
|
---|
3811 | CombinedAlphaTestOFFWithoutNormalization[j-1] +=
|
---|
3812 | TmpTH1FHisto -> GetBinContent(j);
|
---|
3813 |
|
---|
3814 |
|
---|
3815 |
|
---|
3816 |
|
---|
3817 |
|
---|
3818 |
|
---|
3819 | // Histogram bin contents must be normalized (to the respective ON data)
|
---|
3820 | // by using the normalization constants stored in
|
---|
3821 | // histogram NormFactorTestHisto before being added
|
---|
3822 | // to CombinedAlphaTestOFF
|
---|
3823 |
|
---|
3824 |
|
---|
3825 | tmpdouble2 = CombinedAlphaTestOFF->GetBinContent(j) +
|
---|
3826 | (TmpTH1FHisto->GetBinContent(j) * tmpnormfactor);
|
---|
3827 |
|
---|
3828 |
|
---|
3829 | CombinedAlphaTestOFF-> SetBinContent(j,tmpdouble2);
|
---|
3830 |
|
---|
3831 |
|
---|
3832 | // (Bin Error)^2 is computed and added to
|
---|
3833 | // ErrorInCombinedAlphaTestOFF
|
---|
3834 | // Bin error = Sqrt(BinContent) X Normalization Factor
|
---|
3835 |
|
---|
3836 | tmperror = TmpTH1FHisto->GetBinContent(j) * tmpnormfactor;
|
---|
3837 | ErrorInCombinedAlphaTestOFF[j-1] += tmperror;
|
---|
3838 |
|
---|
3839 |
|
---|
3840 | }
|
---|
3841 |
|
---|
3842 | // Dynamic memory allocated for TmpTH1FHisto is released
|
---|
3843 |
|
---|
3844 | delete TmpTH1FHisto;
|
---|
3845 |
|
---|
3846 | }
|
---|
3847 |
|
---|
3848 |
|
---|
3849 |
|
---|
3850 |
|
---|
3851 |
|
---|
3852 | }
|
---|
3853 |
|
---|
3854 |
|
---|
3855 |
|
---|
3856 | // Mean Normalization factor is computed for the combined OFF histogram.
|
---|
3857 | // This factor will be used when computing the significance via
|
---|
3858 | // MHFindSignificance::FindSigmaONOFF().
|
---|
3859 |
|
---|
3860 | // The bin contents (and errors) of teh combined OFF TEST histo will be
|
---|
3861 | // normalized (divided) by teh mean norm constant, so that
|
---|
3862 | // the histogram for OFF data given as argument in FindSigmaONOFF
|
---|
3863 | // is equivalent to a non-normalized histogram.
|
---|
3864 |
|
---|
3865 | MeanNormFactorTest = MeanNormFactorTest/TotalNumberOFFTest;
|
---|
3866 |
|
---|
3867 |
|
---|
3868 |
|
---|
3869 | // Sqrt(Contents) is performed in ErrorInCombinedAlphaTestOFF,
|
---|
3870 | // in order to compute the real error in the bin content of histogram
|
---|
3871 | // CombinedAlphaTestOFF.
|
---|
3872 |
|
---|
3873 | // Then error and contents (normalized with mean normfctor) are
|
---|
3874 | // set to histogram bin
|
---|
3875 |
|
---|
3876 |
|
---|
3877 |
|
---|
3878 |
|
---|
3879 |
|
---|
3880 | for (Int_t j = 1; j <= NAlphaBins; j++)
|
---|
3881 | {
|
---|
3882 |
|
---|
3883 | tmpdouble2 =
|
---|
3884 | (CombinedAlphaTestOFF -> GetBinContent(j))/MeanNormFactorTest;
|
---|
3885 |
|
---|
3886 | CombinedAlphaTestOFF -> SetBinContent(j,tmpdouble2);
|
---|
3887 |
|
---|
3888 | ErrorInCombinedAlphaTestOFF[j-1] =
|
---|
3889 | TMath::Sqrt(ErrorInCombinedAlphaTestOFF[j-1])/MeanNormFactorTest;
|
---|
3890 |
|
---|
3891 | // Proper error treatment for the bin contents
|
---|
3892 | CombinedAlphaTestOFF -> SetBinError(j,ErrorInCombinedAlphaTestOFF[j-1]);
|
---|
3893 |
|
---|
3894 |
|
---|
3895 | // ****** TMP **************
|
---|
3896 | /*
|
---|
3897 | // Error computation assuming a given (no normalized) bin content
|
---|
3898 |
|
---|
3899 | tmpdouble = CombinedAlphaTestOFF -> GetBinContent(j);
|
---|
3900 |
|
---|
3901 | if (CombinedAlphaTestOFFWithoutNormalization[j-1] < SmallQuantity)
|
---|
3902 | {tmpnormfactor = 0;}
|
---|
3903 | else
|
---|
3904 | {tmpnormfactor = tmpdouble/CombinedAlphaTestOFFWithoutNormalization[j-1];}
|
---|
3905 |
|
---|
3906 | tmperror = (TMath::Sqrt(CombinedAlphaTestOFFWithoutNormalization[j-1])) * tmpnormfactor;
|
---|
3907 |
|
---|
3908 | CombinedAlphaTestOFF -> SetBinError(j,tmperror);
|
---|
3909 |
|
---|
3910 | */
|
---|
3911 |
|
---|
3912 | // ****** ENDTMP **************
|
---|
3913 |
|
---|
3914 | // tmp
|
---|
3915 | // cout << CombinedAlphaTestOFF -> GetBinContent(j) << " +/- "
|
---|
3916 | // << CombinedAlphaTestOFF -> GetBinError(j) << endl;
|
---|
3917 |
|
---|
3918 | // endtmp
|
---|
3919 |
|
---|
3920 | }
|
---|
3921 |
|
---|
3922 |
|
---|
3923 |
|
---|
3924 | // Entries of histogram CombinedAlphaTestOFF is set to
|
---|
3925 | // the events counted during histogram filling
|
---|
3926 |
|
---|
3927 | EventCounter = EventCounter/MeanNormFactorTest;
|
---|
3928 | tmpint = int(EventCounter);
|
---|
3929 | CombinedAlphaTestOFF->SetEntries(tmpint);
|
---|
3930 |
|
---|
3931 |
|
---|
3932 | cout << "SILLY INFO FOR TEST SAMPLE: Total number of events Before nomralization and "
|
---|
3933 | << " re-re normalized" << endl
|
---|
3934 | << TotalNumberOFFTest << " - " << EventCounter << endl;
|
---|
3935 |
|
---|
3936 | /*
|
---|
3937 | // Sumw2() test
|
---|
3938 |
|
---|
3939 | cout << "******** TEST ********* " << endl;
|
---|
3940 | cout << "Sumw2() is executed: " << endl;
|
---|
3941 | CombinedAlphaTestOFF -> Sumw2();
|
---|
3942 | for (Int_t j = 1; j <= NAlphaBins; j++)
|
---|
3943 | {
|
---|
3944 | // tmp
|
---|
3945 | cout << CombinedAlphaTestOFF -> GetBinContent(j) << " +/- "
|
---|
3946 | << CombinedAlphaTestOFF -> GetBinError(j) << endl;
|
---|
3947 |
|
---|
3948 | // endtmp
|
---|
3949 |
|
---|
3950 | }
|
---|
3951 |
|
---|
3952 | */
|
---|
3953 |
|
---|
3954 |
|
---|
3955 |
|
---|
3956 | }
|
---|
3957 |
|
---|
3958 |
|
---|
3959 |
|
---|
3960 | // ********* COMBINED HISTOGRAMS FOR TEST SAMPLE ALREADY FILLED ********* ///
|
---|
3961 |
|
---|
3962 |
|
---|
3963 | /*
|
---|
3964 | CombinedAlphaTrainON-> DrawCopy();
|
---|
3965 | gPad -> SaveAs("tmpONTrain.ps");
|
---|
3966 |
|
---|
3967 | CombinedAlphaTrainOFF-> DrawCopy();
|
---|
3968 | gPad -> SaveAs("tmpOFFTrain.ps");
|
---|
3969 |
|
---|
3970 |
|
---|
3971 | CombinedAlphaTestON-> DrawCopy();
|
---|
3972 | gPad -> SaveAs("tmpONTest.ps");
|
---|
3973 |
|
---|
3974 | CombinedAlphaTestOFF-> DrawCopy();
|
---|
3975 | gPad -> SaveAs("tmpOFFTest.ps");
|
---|
3976 |
|
---|
3977 | */
|
---|
3978 |
|
---|
3979 |
|
---|
3980 | // **** COMPUTATION OF SIGNIFICANCES FOR COMBINED ALPHA HISTOGRAMS ****** ///
|
---|
3981 |
|
---|
3982 | if (SuccessfulThetaBinsVector.GetSize() >= 1)
|
---|
3983 | {
|
---|
3984 | // There is, at least, ONE theta bin for which optimization of
|
---|
3985 | // supercuts was SUCCESSFUL, and thus, it is worth to
|
---|
3986 | // compute significance of signal
|
---|
3987 |
|
---|
3988 |
|
---|
3989 | // Significance is found using MHFindSignificanceONOFF::FindSigmaONOFF
|
---|
3990 |
|
---|
3991 |
|
---|
3992 |
|
---|
3993 | // Maximum value of alpha below which signal is expected
|
---|
3994 | const Double_t alphasig = fAlphaSig;
|
---|
3995 |
|
---|
3996 | // Minimum value of alpha for bkg region in ON data
|
---|
3997 | const Double_t alphabkgmin = fAlphaBkgMin;
|
---|
3998 |
|
---|
3999 | // Maximum value of alpha for bkg region in ON data
|
---|
4000 | const Double_t alphabkgmax = fAlphaBkgMax;
|
---|
4001 |
|
---|
4002 | // degree of polynomial function used to fit ON data in background region
|
---|
4003 | const Int_t degree = 2;
|
---|
4004 |
|
---|
4005 | // degree of polynomial function used to fit OFF data in all alpha region (0-90)
|
---|
4006 | const Int_t degreeOFF = 2;
|
---|
4007 |
|
---|
4008 | const Bool_t drawpoly = kTRUE;
|
---|
4009 | const Bool_t fitgauss = kTRUE;
|
---|
4010 | const Bool_t print = kTRUE;
|
---|
4011 |
|
---|
4012 | const Bool_t saveplots = kTRUE; // Save final alpha plot for ON and OFF
|
---|
4013 | // (with Nex and significance) in Psfile
|
---|
4014 |
|
---|
4015 |
|
---|
4016 |
|
---|
4017 |
|
---|
4018 | MHFindSignificanceONOFF findsig;
|
---|
4019 | findsig.SetRebin(kTRUE);
|
---|
4020 | findsig.SetReduceDegree(kFALSE);
|
---|
4021 |
|
---|
4022 |
|
---|
4023 |
|
---|
4024 | if (CombineTrainData)
|
---|
4025 | {
|
---|
4026 |
|
---|
4027 | // Name of psfile where Alpha plot will be stored
|
---|
4028 | TString psfilename = ("TRAINSampleCombined");
|
---|
4029 |
|
---|
4030 | findsig.FindSigmaONOFF(CombinedAlphaTrainON, CombinedAlphaTrainOFF,
|
---|
4031 | MeanNormFactorTrain,
|
---|
4032 | alphabkgmin, alphabkgmax,
|
---|
4033 | degree, degreeOFF,
|
---|
4034 | alphasig,
|
---|
4035 | drawpoly, fitgauss, print,
|
---|
4036 | saveplots, psfilename);
|
---|
4037 |
|
---|
4038 |
|
---|
4039 | fOverallSigmaLiMaTrain = double(findsig.GetSignificance());
|
---|
4040 | fOverallNexTrain = double(findsig.GetNexONOFF());
|
---|
4041 |
|
---|
4042 | *fLog << "MFindSupercutsONOFFThetaLoop::ComputeOverallSignificance" << endl
|
---|
4043 | << "After performing the combined analysis for TRAIN sample it was found "
|
---|
4044 | << "the following Nex and Significance: " << fOverallNexTrain
|
---|
4045 | << ", " << fOverallSigmaLiMaTrain << endl;
|
---|
4046 |
|
---|
4047 |
|
---|
4048 | }
|
---|
4049 |
|
---|
4050 |
|
---|
4051 | if (CombineTestData)
|
---|
4052 | {
|
---|
4053 |
|
---|
4054 | // Name of psfile where Alpha plot will be stored
|
---|
4055 | TString psfilename = ("TESTSampleCombined");
|
---|
4056 |
|
---|
4057 | findsig.FindSigmaONOFF(CombinedAlphaTestON, CombinedAlphaTestOFF,
|
---|
4058 | MeanNormFactorTest,
|
---|
4059 | alphabkgmin, alphabkgmax,
|
---|
4060 | degree, degreeOFF,
|
---|
4061 | alphasig,
|
---|
4062 | drawpoly, fitgauss, print,
|
---|
4063 | saveplots, psfilename);
|
---|
4064 |
|
---|
4065 | fOverallSigmaLiMaTest = double(findsig.GetSignificance());
|
---|
4066 | fOverallNexTest = double(findsig.GetNexONOFF());
|
---|
4067 |
|
---|
4068 | *fLog << "MFindSupercutsONOFFThetaLoop::ComputeOverallSignificance" << endl
|
---|
4069 | << "After performing the combined analysis for TEST sample it was found "
|
---|
4070 | << "the following Nex and Significance: " << fOverallNexTest
|
---|
4071 | << ", " << fOverallSigmaLiMaTest << endl;
|
---|
4072 |
|
---|
4073 | }
|
---|
4074 |
|
---|
4075 |
|
---|
4076 |
|
---|
4077 |
|
---|
4078 | }
|
---|
4079 |
|
---|
4080 |
|
---|
4081 | return kTRUE;
|
---|
4082 | }
|
---|
4083 |
|
---|
4084 |
|
---|
4085 |
|
---|
4086 |
|
---|
4087 |
|
---|
4088 |
|
---|
4089 |
|
---|
4090 |
|
---|