source: tags/Mars-V0.8.5/macros/calibration.C

Last change on this file was 4348, checked in by gaug, 21 years ago
*** empty log message ***
File size: 14.8 KB
Line 
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): Markus Gaug, 11/2003 <mailto:markus@ifae.es>
19!
20! Copyright: MAGIC Software Development, 2000-2004
21!
22!
23\* ======================================================================== */
24/////////////////////////////////////////////////////////////////////////////
25//
26// calibration.C
27//
28// Needs as arguments the run number of a calibration file ("*_C_*.root") and
29// the run number of the corresponding pedestal file ("*_P_*.root").
30//
31// The TString inpath has to be set correctly.
32//
33// The macro searches for the pulser colour which corresponds to the calibration
34// run number. If the run number is smaller than 20000, pulser colour "CT1"
35// is assumed, otherwise, it searches for the strings "green", "blue", "uv" or
36// "ct1" in the filenames. If no colour or multiple colours are found, the
37// execution is aborted.
38//
39// The container MBadPixelsCam is created and followed during the execution of the
40// rest of the macro.
41//
42// A first loop over the pedestal file is performed using the class MJPedestal
43//
44// The container MCalibrationQECam is created and followed during the execution of the
45// rest of the macro.
46//
47// A loop over the calibration files is performed using the class MJCalibration.
48// The results are displayed using the MJCalibration::SetNormalDisplay() facility,
49// but other displays can easily be uncommented.
50// The call to MJCalibration skips the relative time calibration, which can be
51// uncommented as well.
52//
53// Last, a third loop is performed over the calibration file again in order to
54// "calibrate" it and test the resulting outcome.
55//
56/////////////////////////////////////////////////////////////////////////////
57#include "MJPedestal.h"
58#include "MJCalibration.h"
59#include "MJExtractSignal.h"
60#include "MJExtractCalibTest.h"
61#include "MExtractFixedWindowPeakSearch.h"
62#include "MExtractSlidingWindow.h"
63#include "MExtractFixedWindow.h"
64#include "MExtractFixedWindowSpline.h"
65#include "MExtractAmplitudeSpline.h"
66#include "MExtractTimeHighestIntegral.h"
67#include "MExtractTimeFastSpline.h"
68#include "MRunIter.h"
69#include "MStatusDisplay.h"
70#include "MCalibrationQECam.h"
71#include "MHCalibrationTestCam.h"
72#include "MHCalibrationTestPix.h"
73#include "MBadPixelsCam.h"
74#include "MArgs.h"
75#include "MArray.h"
76#include "MLog.h"
77#include "MParContainer.h"
78
79#include "TStyle.h"
80#include "TObject.h"
81#include "TObjectTable.h"
82#include "TSystem.h"
83
84#include <fstream>
85
86using namespace std;
87
88static TString outpath = "./";
89static TString inpath = "/home/rootdata/Calib/2004_05_24/";
90static TString badfile = "";
91//
92// the default pedestal run for the calibration
93//
94static const Int_t pedrun = 26569;
95//
96// the default start calibration run
97//
98static const Int_t calrun1 = 26568;
99//
100// the default last calibration run (if 0, only one run is taken, otherwise consecutive runs
101// between calrun1 and calrun2)
102//
103static const Int_t calrun2 = 0;
104//
105// A switch to output debugging information about Objects use
106//
107static Bool_t debug = kFALSE;
108//
109// A switch to use the Blind Pixel
110//
111static Bool_t blindpix = kTRUE;
112//
113// A switch to use the PIN Diode
114//
115static Bool_t pindiode = kFALSE;
116//
117// Tell if you want to calibrate times:
118//
119static Bool_t useTimes = kTRUE;
120//
121// Tell if you want to use the display:
122//
123static Bool_t useDisplay = kTRUE;
124//
125// Tell if you want to test the result afterwards
126//
127static Bool_t useTest = kTRUE;
128//
129Int_t calibration(const Int_t prun=pedrun,
130 const Int_t crun1=calrun1, const Int_t crun2=calrun2)
131{
132
133 if (debug)
134 TObject::SetObjectStat(kTRUE);
135
136 //
137 // Choose the signal Extractor:
138 //
139 // MExtractFixedWindowPeakSearch extractor;
140 // MExtractSlidingWindow extractor;
141 MExtractFixedWindowSpline extractor;
142 // MExtractAmplitudeSpline extractor;
143 // MExtractTimeAndChargeSpline extractor;
144 //
145 // Set Ranges or Windows
146 //
147 extractor.SetRange(1,14,2,13);
148 // extractor.SetWindows(8,8);
149
150 //
151 // Choose the arrival time Extractor:
152 //
153 // MExtractTimeHighestIntegral timeext;
154 MExtractTimeFastSpline timeext;
155 //
156 // Set Ranges or Windows
157 //
158 timeext.SetRange(0,7,3,8);
159
160 MRunIter pruns;
161 MRunIter cruns;
162
163 pruns.AddRun(prun,inpath);
164
165 if (crun2==0)
166 cruns.AddRun(crun1,inpath);
167 else
168 cruns.AddRuns(crun1,crun2,inpath);
169
170 gStyle->SetOptStat(1111);
171 gStyle->SetOptFit();
172 gStyle->SetTitleSize(0.1,"u");
173 gStyle->SetLineWidth(1);
174
175 MStatusDisplay *display = NULL;
176
177 if (useDisplay)
178 {
179 display = new MStatusDisplay;
180 display->SetUpdateTime(3000);
181 display->Resize(850,700);
182 }
183
184 /************************************/
185 /* FIRST LOOP: PEDESTAL COMPUTATION */
186 /************************************/
187
188 MCalibrationQECam qecam;
189 MBadPixelsCam badcam;
190 //
191 // If you want to exclude pixels from the beginning, read
192 // an ascii-file with the corr. pixel numbers (see MBadPixelsCam)
193 //
194 if (!badfile.IsNull())
195 {
196 ifstream f(badfile.Data());
197 badcam.AsciiRead(f);
198 f.close();
199 }
200
201 MJPedestal pedloop;
202 pedloop.SetExtractor(&extractor);
203 pedloop.SetInput(&pruns);
204 pedloop.SetOutputPath(outpath.Data());
205 if (useDisplay)
206 pedloop.SetDisplay(display);
207 pedloop.SetBadPixels(badcam);
208
209 if (!pedloop.Process())
210 return 1;
211
212 /****************************************/
213 /* SECOND LOOP: CALIBRATION COMPUTATION */
214 /****************************************/
215
216 MJCalibration calloop;
217
218 if (debug)
219 calloop.SetDebug();
220 //
221 // If you want to run the data-check on RAW DATA!!!, choose:
222 // calloop.SetDataCheck();
223 //
224 // If you want to see the data-check plots only, choose:
225 calloop.SetDataCheckDisplay();
226 //
227 // For everything, you have ever dreamed of, choose:
228 // calloop.SetFullDisplay();
229
230 //
231 // If you want to calibrate the times as well, choose:
232 //
233 calloop.SetRelTimeCalibration(useTimes);
234 calloop.SetExtractor(&extractor);
235 calloop.SetTimeExtractor(&timeext);
236 calloop.SetInput(&cruns);
237 calloop.SetOutputPath(outpath.Data());
238 if (useDisplay)
239 calloop.SetDisplay(display);
240 calloop.SetUseBlindPixel(blindpix);
241 calloop.SetUsePINDiode(pindiode);
242 calloop.SetQECam(qecam);
243 calloop.SetBadPixels(pedloop.GetBadPixels());
244
245 if (!calloop.Process(pedloop.GetPedestalCam()))
246 return 2;
247
248 //
249 // The next lines are the use the Print() function and have
250 // all the results as ascii-tables:
251 //
252 if (debug)
253 {
254 MCalibrationChargeCam &chargecam = calloop.GetCalibrationCam();
255 MCalibrationQECam &nqecam = calloop.GetQECam();
256 MBadPixelsCam &badbad = calloop.GetBadPixels();
257 chargecam.Print();
258 nqecam.Print();
259 badbad.Print();
260 }
261
262 gLog << endl;
263 gLog << "Mean number of photons from pulser Inner pixels (F-Factor Method): "
264 << calloop.GetCalibrationCam().GetNumPhotonsFFactorMethod()
265 << " +- " << calloop.GetCalibrationCam().GetNumPhotonsFFactorMethodErr() << endl;
266 gLog << endl;
267
268 /********************************************************************/
269 /* THIRD LOOP: APPLY CALIBRATION TO THE CALIBRATION FILES FOR TESTS */
270 /********************************************************************/
271
272 if (useTest)
273 {
274
275 MJExtractCalibTest testloop;
276
277 testloop.SetExtractor(&extractor);
278 testloop.SetTimeExtractor(&timeext);
279 testloop.SetInput(&cruns);
280 testloop.SetOutputPath(outpath);
281 if (useDisplay)
282 testloop.SetDisplay(display);
283 testloop.SetBadPixels(calloop.GetBadPixels());
284
285 if (!testloop.ProcessD(pedloop.GetPedestalCam(),calloop.GetCalibrationCam(),calloop.GetQECam()))
286 return 3;
287
288 if (useTimes)
289 if (!testloop.ProcessT(pedloop.GetPedestalCam(),calloop.GetRelTimeCam()))
290 return 4;
291
292 gLog << endl;
293 gLog << "Mean equiv. number of photons from cascades per mm^2 Inner pixels: "
294 << testloop.GetTestCam().GetMeanMeanPhotPerArea(0)
295 << " +- " << testloop.GetTestCam().GetRmsMeanPhotPerArea(0) << endl;
296 gLog << "Sigma equiv. number of photons from cascades per mm^2 Inner pixels: "
297 << testloop.GetTestCam().GetMeanSigmaPhotPerArea(0)
298 << " +- " << testloop.GetTestCam().GetRmsSigmaPhotPerArea(0) << endl;
299
300 gLog << endl;
301 gLog << "Mean equiv. number of photons from cascades per mm^2 Outer pixels: "
302 << testloop.GetTestCam().GetMeanMeanPhotPerArea(1)
303 << " +- " << testloop.GetTestCam().GetRmsMeanPhotPerArea(1) << endl;
304 gLog << "Sigma equiv. number of photons from cascades per mm^2 Outer pixels: "
305 << testloop.GetTestCam().GetMeanSigmaPhotPerArea(1)
306 << " +- " << testloop.GetTestCam().GetRmsSigmaPhotPerArea(1) << endl;
307 gLog << endl;
308
309 }
310
311 /********************************************************************/
312 /* FOURTH LOOP: APPLY CALIBRATION TO THE PEDESTAL FILES */
313 /********************************************************************/
314
315 MJExtractSignal pedphotloop;
316
317 pedphotloop.SetExtractor(&extractor);
318 pedphotloop.SetTimeExtractor(&timeext);
319 pedphotloop.SetInput(&pruns);
320 pedphotloop.SetOutputPath(outpath);
321 if (useDisplay)
322 pedphotloop.SetDisplay(display);
323 pedphotloop.SetBadPixels(calloop.GetBadPixels());
324
325 if (!pedphotloop.ProcessP(pedloop.GetPedestalCam(),calloop.GetCalibrationCam(),calloop.GetQECam()))
326 return 5;
327
328
329 if (debug)
330 TObject::SetObjectStat(kFALSE);
331
332 //
333 // Debugging at the end:
334 //
335 if (debug)
336 gObjectTable->Print();
337
338 //
339 // List of useful containers:
340 //
341/*
342 MPedestalCam &pedcam = pedloop.GetPedestalCam();
343 MCalibrationChargeCam &chargecam = calloop.GetCalibrationCam();
344 MCalibrationQECam &qecam = calloop.GetCalibrationCam();
345 MBadPixelsCam &badcam = calloop.GetBadPixels();
346 MHCalibrationTestCam &testcam = testloop.GetTestCam();
347 MHCalibrationTestTimeCam &testtime = testloop.GetTestTimeCam();
348*/
349
350 //
351 // List of interesting plots:
352 //
353/*
354 testcam[200].DrawClone("fourierevents");
355 testcam.GetAverageHiGainArea(0).DrawClone();
356 testcam.GetAverageLoGainArea(0).DrawClone();
357 testcam.GetAverageHiGainArea(1).DrawClone();
358 testcam.GetAverageLoGainArea(1).DrawClone();
359
360 testcam.GetAverageHiGainSector(1).DrawClone();
361 testcam.GetAverageLoGainSector(1).DrawClone();
362 testcam.GetAverageHiGainSector(2).DrawClone();
363 testcam.GetAverageLoGainSector(2).DrawClone();
364 testcam.GetAverageHiGainSector(3).DrawClone();
365 testcam.GetAverageLoGainSector(3).DrawClone();
366 testcam.GetAverageHiGainSector(4).DrawClone();
367 testcam.GetAverageLoGainSector(4).DrawClone();
368 testcam.GetAverageHiGainSector(5).DrawClone();
369 testcam.GetAverageLoGainSector(5).DrawClone();
370 testcam.GetAverageHiGainSector(6).DrawClone();
371 testcam.GetAverageLoGainSector(6).DrawClone();
372*/
373
374 return 0;
375
376}
377
378static void Usage()
379{
380 gLog << endl;
381 gLog << "Usage:" << endl;
382 gLog << endl;
383 gLog << " calibration [ped.run nr.] [first cal.run nr.] [last cal.run nr.]" << endl ;
384 gLog << endl;
385 gLog << " ped.run.nr: Run number of the pedestal file." << endl;
386 gLog << " first cal.run nr.: Run number of the first calibration file." << endl;
387 gLog << " last cal.run nr.: Run number of the last calibration file." << endl;
388 gLog << endl;
389 gLog << "All calibration runs between (first cal.run nr.) and (last cal.run nr.) will be used" << endl;
390 gLog << "If last.cal.run.nr is 0 (default), only one calibration run is taken" << endl;
391 gLog << endl;
392 gLog << "Additional Options: " << endl;
393 gLog << " --inpath=# Find the data in inpath" << endl;
394 gLog << " --outpath=# Write the output containers to outpath" << endl;
395 gLog << " --badfile=# Use the file # to exclude pixels from the beginning" << endl;
396 gLog << " --debug Use the TObjectTable for debugging " << endl;
397 gLog << " and write out the pixels as ascii tables" << endl;
398 gLog << " --useTimes Calibrate the relative arrival times" << endl;
399 gLog << " --useTest Use the class MJExtractCalibTest to test the calibration on itself" << endl;
400 gLog << " --skipBlindPix Skip the blind pixel calibration" << endl;
401 gLog << " --skipPINDiode Skip the PIN Diode calibration" << endl;
402}
403
404
405int main(int argc, char **argv)
406{
407 //
408 // Evaluate arguments
409 //
410 MArgs arg(argc, argv);
411
412 if (arg.HasOnly("-?") || arg.HasOnly("-h") || arg.HasOnly("--help"))
413 {
414 Usage();
415 return -1;
416 }
417
418 debug = arg.HasOnlyAndRemove("--debug") || arg.HasOnlyAndRemove("-d");
419 useTimes = arg.HasOnlyAndRemove("--useTimes") || arg.HasOnlyAndRemove("-t");
420 useTest = arg.HasOnlyAndRemove("--useTest") || arg.HasOnlyAndRemove("-e");
421 blindpix = !(arg.HasOnlyAndRemove("--skipBlindPix"));
422 pindiode = !(arg.HasOnlyAndRemove("--skipPINDiode"));
423
424 if (arg.HasOption("--inpath="))
425 inpath = arg.GetStringAndRemove("--inpath=");
426
427 if (arg.HasOption("--outpath="))
428 outpath = arg.GetStringAndRemove("--outpath=");
429
430 if (arg.HasOption("--badfile="))
431 badfile = arg.GetStringAndRemove("--badfile=");
432
433 if (gSystem->AccessPathName(badfile,kFileExists))
434 {
435 gLog << "WARNING: the bad pixels file '" << badfile.Data() << "' doesn't exist." << endl;
436 badfile = "";
437 }
438
439 // check for the right usage of the program
440 //
441 if (arg.GetNumArguments()>4)
442 {
443 Usage();
444 return -1;
445 }
446
447 //
448 // Initialize Non-GUI (batch) mode
449 //
450 gROOT->SetBatch();
451
452 //
453 // Switch off the display
454 //
455 useDisplay = kFALSE;
456
457
458 //
459 // check for the arguments
460 //
461 Int_t pedr = 0;
462 Int_t calr1 = 0;
463 Int_t calr2 = 0;
464
465 const Int_t nargs = arg.GetNumArguments();
466
467 if (nargs>=3)
468 {
469 pedr = arg.GetArgumentInt(0);
470 calr1 = arg.GetArgumentInt(1);
471 calr2 = arg.GetArgumentInt(2);
472 return calibration(pedr,calr1,calr2);
473 }
474
475 if (nargs>=2)
476 {
477 pedr = arg.GetArgumentInt(0);
478 calr1 = arg.GetArgumentInt(1);
479 return calibration(pedr,calr1);
480 }
481
482 if (nargs>=1)
483 {
484 pedr = arg.GetArgumentInt(0);
485 gLog << "PEDR: " << pedr << endl;
486 return calibration(pedr);
487 }
488
489 return calibration();
490}
Note: See TracBrowser for help on using the repository browser.