source: trunk/MagicSoft/Mars/mjobs/MJExtractCalibTest.cc@ 6073

Last change on this file since 6073 was 5942, checked in by mazin, 20 years ago
*** empty log message ***
File size: 14.4 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, 04/2004 <mailto:markus@ifae.es>
19!
20! Copyright: MAGIC Software Development, 2000-2004
21!
22!
23\* ======================================================================== */
24
25/////////////////////////////////////////////////////////////////////////////
26//
27// MJExtractCalibTest
28//
29// If the flag SetDataCheckDisplay() is set, only the most important distributions
30// are displayed.
31// Otherwise, (default: SetNormalDisplay()), a good selection of plots is given
32//
33/////////////////////////////////////////////////////////////////////////////
34#include "MJExtractCalibTest.h"
35
36#include <TFile.h>
37#include <TStyle.h>
38#include <TCanvas.h>
39#include <TSystem.h>
40
41#include "MLog.h"
42#include "MLogManip.h"
43
44#include "MRunIter.h"
45#include "MParList.h"
46#include "MTaskList.h"
47#include "MTaskEnv.h"
48#include "MEvtLoop.h"
49
50#include "MHCamera.h"
51
52#include "MPedestalCam.h"
53#include "MBadPixelsCam.h"
54#include "MBadPixelsTreat.h"
55#include "MBadPixelsCalc.h"
56#include "MBadPixelsMerge.h"
57#include "MCerPhotEvt.h"
58#include "MArrivalTime.h"
59#include "MCalibrationChargeCam.h"
60#include "MCalibrationRelTimeCam.h"
61#include "MCalibrationQECam.h"
62#include "MCalibrationTestCam.h"
63#include "MCalibrationTestCalc.h"
64#include "MHCamEvent.h"
65#include "MHCalibrationTestCam.h"
66
67#include "MReadMarsFile.h"
68#include "MGeomApply.h"
69#include "MExtractTimeAndChargeSlidingWindow.h"
70#include "MExtractor.h"
71#include "MExtractTime.h"
72#include "MExtractTimeFastSpline.h"
73#include "MFCosmics.h"
74#include "MContinue.h"
75#include "MFillH.h"
76#include "MCalibrateData.h"
77#include "MCalibrateRelTimes.h"
78
79#include "MTriggerPattern.h"
80#include "MTriggerPatternDecode.h"
81#include "MFTriggerPattern.h"
82
83#include "MStatusDisplay.h"
84
85ClassImp(MJExtractCalibTest);
86
87using namespace std;
88// --------------------------------------------------------------------------
89//
90// Default constructor.
91//
92// Sets fUseCosmicsFilter to kTRUE, fRuns to 0, fExtractor to NULL, fTimeExtractor to NULL
93// fDisplay to kNormalDisplay
94//
95MJExtractCalibTest::MJExtractCalibTest(const char *name, const char *title)
96 : fUseCosmicsFilter(kTRUE), fRuns(NULL), fExtractor(NULL), fTimeExtractor(NULL),
97 fDisplayType(kNormalDisplay), fGeometry("MGeomCamMagic")
98{
99 fName = name ? name : "MJExtractCalibTest";
100 fTitle = title ? title : "Tool to extract, calibrate and test signals from a file";
101}
102
103
104void MJExtractCalibTest::DisplayResult(MParList &plist)
105{
106 if (!fDisplay)
107 return;
108
109 //
110 // Update display
111 //
112 TString title = fDisplay->GetTitle();
113 title += "-- Extraction-Calibration-Test ";
114 title += fRuns->GetRunsAsString();
115 title += " --";
116 fDisplay->SetTitle(title);
117
118 //
119 // Get container from list
120 //
121 MGeomCam &geomcam = *(MGeomCam*) plist.FindObject("MGeomCam");
122 MHCalibrationTestCam &testcam = *(MHCalibrationTestCam*)plist.FindObject("MHCalibrationTestCam");
123
124 // Create histograms to display
125 MHCamera disp1 (geomcam, "Test;Photons", "Mean of calibrated Photons");
126 MHCamera disp2 (geomcam, "Test;SigmaPhotons", "Sigma of calibrated photons");
127 MHCamera disp3 (geomcam, "Test;PhotonsPerArea", "Equiv. Cherenkov Photons per Area");
128 MHCamera disp4 (geomcam, "Test;SigmaPhotPerArea", "Sigma equiv. Cher. Photons per Area");
129 MHCamera disp5 (geomcam, "Test;Phot", "Calibrated Photons");
130 MHCamera disp6 (geomcam, "Test;PhotPerArea", "Calibrated Photons per Area");
131 MHCamera disp7 (geomcam, "Test;NotInterpolate", "Not interpolated pixels");
132 MHCamera disp8 (geomcam, "Test;DeviatingPhots", "Deviating Number Photons");
133
134 // Fitted charge means and sigmas
135 disp1.SetCamContent(testcam, 0);
136 disp1.SetCamError( testcam, 1);
137 disp2.SetCamContent(testcam, 2);
138 disp2.SetCamError( testcam, 3);
139 disp3.SetCamContent(testcam, 7);
140 disp3.SetCamError( testcam, 8);
141 disp4.SetCamContent(testcam, 9);
142 disp4.SetCamError( testcam, 10);
143
144 disp5.SetCamContent(fTestCam, 0);
145 disp5.SetCamError( fTestCam, 1);
146 disp6.SetCamContent(fTestCam, 2);
147 disp6.SetCamError( fTestCam, 3);
148 disp7.SetCamError( fTestCam, 4);
149
150 disp8.SetCamError( fBadPixels, 22);
151
152
153 disp1.SetYTitle("Photons");
154 disp2.SetYTitle("\\sigma_{phot}");
155 disp3.SetYTitle("Photons per Area [mm^{-2}]");
156 disp4.SetYTitle("\\sigma_{phot} per Area [mm^{-2}]");
157
158 disp5.SetYTitle("Photons");
159 disp6.SetYTitle("Photons per Area [mm^{-2}]");
160 disp7.SetYTitle("[1]");
161 disp8.SetYTitle("[1]");
162
163 gStyle->SetOptStat(1111);
164 gStyle->SetOptFit();
165
166 if (fDisplayType == kNormalDisplay)
167 {
168
169 TCanvas &c = fDisplay->AddTab("TestCharges");
170 c.Divide(4,4);
171
172 disp1.CamDraw(c, 1, 4, 2, 1);
173 disp2.CamDraw(c, 2, 4, 2, 1);
174 disp3.CamDraw(c, 3, 4, 1, 1);
175 disp4.CamDraw(c, 4, 4, 2, 1);
176 }
177
178 TCanvas &c2 = fDisplay->AddTab("TestResult");
179 c2.Divide(2,4);
180
181 disp5.CamDraw(c2, 1, 2, 2, 1);
182 disp6.CamDraw(c2, 2, 2, 2, 1);
183
184 TCanvas &c3 = fDisplay->AddTab("TestDefects");
185 c3.Divide(2,2);
186
187 disp7.CamDraw(c3, 1, 2, 0);
188 disp8.CamDraw(c3, 2, 2, 0);
189
190 return;
191
192}
193
194
195void MJExtractCalibTest::DisplayResultT(MParList &plist)
196{
197 if (!fDisplay)
198 return;
199
200 //
201 // Update display
202 //
203 TString title = fDisplay->GetTitle();
204 title += "-- Extraction-Calibration-Test-Time";
205 title += fRuns->GetRunsAsString();
206 title += " --";
207 fDisplay->SetTitle(title);
208
209 //
210 // Get container from list
211 //
212 MGeomCam &geomcam = *(MGeomCam*)plist.FindObject("MGeomCam");
213
214 // Create histograms to display
215 MHCamera disp1 (geomcam, "Test;Arr.Times", "Mean of calibrated Arr.Times");
216 MHCamera disp2 (geomcam, "Test;SigmaArr.Times", "Sigma of calibrated Arr.Times");
217
218 // Fitted charge means and sigmas
219 disp1.SetCamContent(fTestTimeCam, 0);
220 disp1.SetCamError( fTestTimeCam, 1);
221 disp2.SetCamContent(fTestTimeCam, 2);
222 disp2.SetCamError( fTestTimeCam, 3);
223
224 disp1.SetYTitle("Mean Arr.Times [FADC units]");
225 disp2.SetYTitle("\\sigma_{t} [FADC units]");
226
227 gStyle->SetOptStat(1111);
228 gStyle->SetOptFit();
229
230 TCanvas &c = fDisplay->AddTab("TestTimes");
231 c.Divide(2,4);
232
233 disp1.CamDraw(c, 1, 2, 5, 1);
234 disp2.CamDraw(c, 2, 2, 5, 1);
235
236 return;
237
238}
239
240
241const char* MJExtractCalibTest::GetOutputFile() const
242{
243
244 if (fSequence.IsValid())
245 return Form("%s/test%08d.root", (const char*)fPathOut, fSequence.GetSequence());
246
247 if (!fRuns)
248 return "";
249
250 return Form("%s/%s-Test.root", (const char*)fPathOut, (const char*)fRuns->GetRunsAsFileName());
251
252}
253
254
255Bool_t MJExtractCalibTest::Process(MPedestalCam &pedcam,
256 MCalibrationChargeCam &calcam, MCalibrationQECam &qecam,
257 MCalibrationRelTimeCam &reltime,
258 Byte_t filetype)
259{
260 // const TString fname = GetOutputFile();
261
262// if (gSystem->AccessPathName(fname, kFileExists))
263 return ProcessFile(pedcam,calcam,qecam,reltime,filetype);
264
265 // return kTRUE;
266}
267
268Bool_t MJExtractCalibTest::ProcessFile(MPedestalCam &pedcam,
269 MCalibrationChargeCam &calcam, MCalibrationQECam &qecam,
270 MCalibrationRelTimeCam &relcam,
271 Byte_t filetype)
272{
273 if (!fRuns)
274 {
275 *fLog << err << "No Runs choosen... abort." << endl;
276 return kFALSE;
277 }
278 if (fRuns->GetNumRuns() != fRuns->GetNumEntries())
279 {
280 *fLog << err << "Number of files found doesn't match number of runs... abort." << endl;
281 return kFALSE;
282 }
283
284 *fLog << inf;
285 fLog->Separator(GetDescriptor());
286 *fLog << "Calculate MExtractedSignalCam from Runs " << fRuns->GetRunsAsString() << endl;
287 *fLog << endl;
288
289 MBadPixelsCam badcam;
290
291 // Setup Lists
292 MParList plist;
293 plist.AddToList(this); // take care of fDisplay!
294 plist.AddToList(&fTestCam);
295 plist.AddToList(&fTestTimeCam);
296 plist.AddToList(&badcam);
297 plist.AddToList(&pedcam);
298 plist.AddToList(&calcam);
299 plist.AddToList(&qecam);
300 plist.AddToList(&relcam);
301
302 MCerPhotEvt cerphot;
303 MPedPhotCam pedphot;
304 MHCalibrationTestCam testcam;
305
306 plist.AddToList(&cerphot);
307 plist.AddToList(&pedphot);
308 plist.AddToList(&testcam);
309
310 pedcam.SetName("MPedestalFundamental");
311
312 MTaskList tlist;
313 plist.AddToList(&tlist);
314
315 // Setup Task-lists
316 MReadMarsFile read("Events");
317 read.DisableAutoScheme();
318 static_cast<MRead&>(read).AddFiles(*fRuns);
319
320 // Check for interleaved events
321 MTriggerPatternDecode decode;
322 MFTriggerPattern fcalib;
323 fcalib.DenyCalibration();
324 MContinue conttp(&fcalib, "ContTrigPattern");
325
326 MGeomApply apply; // Only necessary to craete geometry
327 apply.SetGeometry(fGeometry);
328 MBadPixelsMerge merge(&fBadPixels);
329
330 MExtractTimeAndChargeSlidingWindow extrsw;
331 MExtractTimeFastSpline extime;
332 extime.SetPedestals(&pedcam);
333
334 MTaskEnv taskenv1("ExtractSignal");
335 MTaskEnv taskenv2("ExtractTime");
336
337 if (fExtractor)
338 {
339 fExtractor->SetPedestals(&pedcam);
340 taskenv1.SetDefault(fExtractor);
341 }
342 if (fTimeExtractor)
343 {
344 fTimeExtractor->SetPedestals(&pedcam);
345 taskenv2.SetDefault(fTimeExtractor);
346 }
347 else
348 {
349 extrsw.SetPedestals(&pedcam);
350 extrsw.SetWindowSize(8,8);
351 taskenv2.SetDefault(&extrsw);
352 *fLog << warn << GetDescriptor()
353 << ": No extractor has been chosen, take default MExtractTimeAndChargeSlidingWindow " << endl;
354 }
355
356 MCalibrateData photcalc;
357 MCalibrateRelTimes caltimes;
358 if (filetype==3) // MC file
359 {
360 photcalc.SetCalibrationMode(MCalibrateData::kFfactor);
361 photcalc.SetPedestalFlag(MCalibrateData::kRun);
362 photcalc.AddPedestal("MPedestalCam", "MPedPhotFundamental");
363 }
364 else
365 {
366 photcalc.SetCalibrationMode(MCalibrateData::kFfactor);
367 photcalc.AddPedestal("Fundamental");
368 photcalc.SetPedestalFlag(MCalibrateData::kEvent);
369 photcalc.SetSignalType(MCalibrateData::kPhot);
370 }
371
372 MBadPixelsCalc badcalc;
373 MBadPixelsTreat badtreat;
374 badtreat.SetProcessTimes(kFALSE);
375
376 badcalc.SetNamePedPhotCam("MPedPhotFundamental");
377 //badtreat.SetUseInterpolation();
378 badtreat.AddNamePedPhotCam("MPedPhotFundamental");
379
380 MCalibrationTestCalc testcalc;
381
382 if (!fSequence.IsValid())
383 {
384 testcalc.SetOutputPath(fPathOut);
385 testcalc.SetOutputFile(Form("%s-TestCalibStat.txt",(const char*)fRuns->GetRunsAsFileName()));
386 }
387
388 MHCamEvent evt0(0,"Signal", "Un-Calibrated Signal;;S [FADC cnts]" );
389 MHCamEvent evt1(3,"CalSig", "Calibrated and Interpolated Signal;;S [\\gamma]");
390 MHCamEvent evt2(0,"Times" , "Arrival Time;;T [slice]");
391
392 MFillH fill0(&evt0, "MExtractedSignalCam", "FillUncalibrated");
393 MFillH fill1(&evt1, "MCerPhotEvt", "FillCalibrated");
394 MFillH fill2(&evt2, "MArrivalTime","FillTimes");
395
396 MFillH fillcam("MHCalibrationTestCam", "MCerPhotEvt");
397 fillcam.SetNameTab("Test");
398 MFillH filltme("MHCalibrationTestTimeCam", "MArrivalTime");
399 filltme.SetNameTab("TestTime");
400
401 MFCosmics cosmics;
402 cosmics.SetNamePedestalCam("MPedestalFundamental");
403 MContinue contcos(&cosmics,"ContCosmics");
404
405 tlist.AddToList(&read);
406 tlist.AddToList(&decode);
407 tlist.AddToList(&apply);
408 tlist.AddToList(&merge);
409 tlist.AddToList(&conttp);
410 tlist.AddToList(&taskenv1);
411 if (!fExtractor->InheritsFrom("MExtractTimeAndCharge"))
412 tlist.AddToList(&taskenv2);
413 tlist.AddToList(&contcos);
414 tlist.AddToList(&fill0);
415 tlist.AddToList(&photcalc);
416 tlist.AddToList(&caltimes);
417 tlist.AddToList(&badcalc);
418 tlist.AddToList(&badtreat);
419 tlist.AddToList(&fill1);
420 tlist.AddToList(&fill2);
421 tlist.AddToList(&fillcam);
422 tlist.AddToList(&filltme);
423 tlist.AddToList(&testcalc);
424
425 // Create and setup the eventloop
426 MEvtLoop evtloop(fName);
427 evtloop.SetParList(&plist);
428 evtloop.SetDisplay(fDisplay);
429 evtloop.SetLogStream(fLog);
430
431 // Execute first analysis
432 if (!evtloop.Eventloop())
433 {
434 *fLog << err << GetDescriptor() << ": Failed." << endl;
435 return kFALSE;
436 }
437
438 tlist.PrintStatistics();
439
440 DisplayResult(plist);
441
442 if (!WriteResult())
443 return kFALSE;
444
445 *fLog << inf << GetDescriptor() << ": Done." << endl;
446
447 return kTRUE;
448}
449
450Bool_t MJExtractCalibTest::ReadPedPhotCam()
451{
452
453 const TString fname = GetOutputFile();
454
455 if (gSystem->AccessPathName(fname, kFileExists))
456 {
457 *fLog << err << "Input file " << fname << " doesn't exist." << endl;
458 return kFALSE;
459 }
460
461 *fLog << inf << "Reading from file: " << fname << endl;
462
463 TFile file(fname, "READ");
464 if (fPedPhotCam.Read()<=0)
465 {
466 *fLog << "Unable to read MPedPhotCam from " << fname << endl;
467 return kFALSE;
468 }
469
470 if (file.FindKey("MBadPixelsCam"))
471 {
472 MBadPixelsCam bad;
473 if (bad.Read()<=0)
474 {
475 *fLog << "Unable to read MBadPixelsCam from " << fname << endl;
476 return kFALSE;
477 }
478 fBadPixels.Merge(bad);
479 }
480
481 if (fDisplay /*&& !fDisplay->GetCanvas("Pedestals")*/) // FIXME!
482 fDisplay->Read();
483
484 return kTRUE;
485}
486
487Bool_t MJExtractCalibTest::WriteResult()
488{
489
490 if (fPathOut.IsNull())
491 return kTRUE;
492
493 const TString oname(GetOutputFile());
494
495 *fLog << inf << "Writing to file: " << oname << endl;
496
497 TFile file(oname, "UPDATE");
498
499 if (fDisplay && fDisplay->Write()<=0)
500 {
501 *fLog << err << "Unable to write MStatusDisplay to " << oname << endl;
502 return kFALSE;
503 }
504
505 if (fPedPhotCam.Write()<=0)
506 {
507 *fLog << err << "Unable to write MPedPhotCam to " << oname << endl;
508 return kFALSE;
509 }
510
511 if (fTestCam.Write()<=0)
512 {
513 *fLog << err << "Unable to write MCalibrationTestCam to " << oname << endl;
514 return kFALSE;
515 }
516
517 if (fTestTimeCam.Write()<=0)
518 {
519 *fLog << err << "Unable to write MCalibrationTestCam to " << oname << endl;
520 return kFALSE;
521 }
522
523 return kTRUE;
524
525}
526
527
528
529Bool_t MJExtractCalibTest::CheckEnv()
530{
531 if (HasEnv("DataCheckDisplay"))
532 fDisplayType = GetEnv("DataCheckDisplay", kFALSE) ? kDataCheckDisplay : kNormalDisplay;
533
534 SetOverwrite(GetEnv("Overwrite", fOverwrite));
535
536 return MJob::CheckEnv();
537}
Note: See TracBrowser for help on using the repository browser.