source: trunk/Mars/mjobs/MJCalibrateSignal.cc@ 11691

Last change on this file since 11691 was 9519, checked in by tbretz, 16 years ago
*** empty log message ***
File size: 35.6 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): Thomas Bretz, 1/2004 <mailto:tbretz@astro.uni-wuerzburg.de>
19!
20! Copyright: MAGIC Software Development, 2000-2008
21!
22!
23\* ======================================================================== */
24
25/////////////////////////////////////////////////////////////////////////////
26//
27// MJCalibrateSignal
28//
29// This class is reading the output written by callisto. It calibrates
30// signal and time.
31//
32// The signal and time extractors are read from the callisto-output. In
33// pricipal you could overwrite these default using the resource file,
34// but this is NOT recommended!
35//
36/////////////////////////////////////////////////////////////////////////////
37#include "MJCalibrateSignal.h"
38
39#include <TEnv.h>
40#include <TFile.h>
41
42// Core
43#include "MLog.h"
44#include "MLogManip.h"
45
46#include "MDirIter.h"
47#include "MTaskList.h"
48#include "MParList.h"
49#include "MEvtLoop.h"
50
51#include "MStatusDisplay.h"
52
53// General containers
54#include "MGeomCam.h"
55#include "MBadPixelsCam.h"
56#include "MCalibConstCam.h"
57#include "MPedestalCam.h"
58#include "MArrivalTimeCam.h"
59
60// General histograms
61#include "MH3.h"
62#include "MHCamEvent.h"
63#include "MHVsTime.h"
64
65// Calibration containers
66#include "MCalibrationQECam.h"
67#include "MCalibrationBlindCam.h"
68#include "MCalibrationChargeCam.h"
69#include "MCalibrationRelTimeCam.h"
70#include "MCalibrationChargePINDiode.h"
71
72// Calibration histograms
73#include "MHCalibrationChargeCam.h"
74#include "MHCalibrationChargeBlindCam.h"
75#include "MHCalibrationChargePINDiode.h"
76#include "MHCalibrationRelTimeCam.h"
77
78// Tasks
79#include "MReadMarsFile.h"
80#include "MRawFileRead.h"
81#include "MTaskEnv.h"
82#include "MContinue.h"
83#include "MFillH.h"
84#include "MGeomApply.h"
85#include "MExtractTimeAndCharge.h"
86#include "MTriggerPatternDecode.h"
87#include "MCalibrationPatternDecode.h"
88#include "MCalibrationChargeCalc.h"
89#include "MCalibrationRelTimeCalc.h"
90#include "MCalibCalcFromPast.h"
91#include "MPedestalSubtract.h"
92#include "MPedCalcFromLoGain.h"
93#include "MCalibrateData.h"
94#include "MExtractPINDiode.h"
95#include "MExtractBlindPixel.h"
96#include "MCalibrateRelTimes.h"
97#include "MBadPixelsCalc.h"
98#include "MBadPixelsTreat.h"
99#include "MWriteRootFile.h"
100
101// Filter
102#include "MFTriggerPattern.h"
103#include "MFCosmics.h"
104#include "MFilterList.h"
105#include "MFDataPhrase.h"
106
107// Classes for writing movies
108#include "MFEvtNumber.h"
109#include "MMoviePrepare.h"
110#include "MMovieWrite.h"
111#include "MImgCleanStd.h"
112
113
114ClassImp(MJCalibrateSignal);
115
116using namespace std;
117
118// --------------------------------------------------------------------------
119//
120// Default constructor.
121//
122// Sets:
123// - fIsInterlaced to kTRUE
124// - fIsRelTimesUpdate to kFALSE
125// - fIsHiLoCalibration to kFALSE
126//
127MJCalibrateSignal::MJCalibrateSignal(const char *name, const char *title)
128 : fExtractor(0), fIsInterlaced(kTRUE), fIsRelTimesUpdate(kTRUE), fIsMovieMode(kFALSE), fIsTestMode(kFALSE)
129{
130 fName = name ? name : "MJCalibrateSignal";
131 fTitle = title ? title : "Tool to calibrate data";
132}
133
134MJCalibrateSignal::~MJCalibrateSignal()
135{
136 if (fExtractor)
137 delete fExtractor;
138}
139
140void MJCalibrateSignal::SetExtractor(const MExtractor *ext)
141{
142 if (fExtractor)
143 delete fExtractor;
144
145 fExtractor = ext ? (MExtractor*)ext->Clone() : NULL;
146}
147
148Bool_t MJCalibrateSignal::WriteResult() const
149{
150 if (IsNoStorage() || fIsMovieMode)
151 return kTRUE;
152
153 // FIXME: This is not nice because it will update the signal
154 // file always. Which might make the usage outside of
155 // callisto difficult.
156 TObjArray cont;
157 cont.Add(const_cast<TEnv*>(GetEnv()));
158 cont.Add(const_cast<MSequence*>(&fSequence));
159
160 TNamed cmdline("CommandLine", fCommandLine.Data());
161 cont.Add(&cmdline);
162
163 if (fDisplay)
164 {
165 TString title = "-- Calibrate Signal: ";
166 title += fSequence.GetSequence();
167 title += " --";
168 fDisplay->SetTitle(title, kFALSE);
169
170 cont.Add(fDisplay);
171 }
172
173 const TString name(Form("signal%08d.root", fSequence.GetSequence()));
174 return WriteContainer(cont, name, "UPDATE");
175}
176
177Bool_t MJCalibrateSignal::ReadCalibration(TObjArray &l, MBadPixelsCam &cam, MExtractor* &ext2, MExtractor* &ext3, TString &geom) const
178{
179 TString fname = Form("%s/calib%08d.root", fPathIn.Data(), fSequence.GetSequence());
180
181 *fLog << inf << "Reading from file: " << fname << endl;
182
183 TFile file(fname, "READ");
184 if (!file.IsOpen())
185 {
186 *fLog << err << dbginf << "ERROR - Could not open file " << fname << endl;
187 return kFALSE;
188 }
189
190 TObject *o = file.Get("ExtractSignal");
191 if (o && !o->InheritsFrom(MExtractor::Class()))
192 {
193 *fLog << err << dbginf << "ERROR - ExtractSignal read from " << fname << " doesn't inherit from MExtractor!" << endl;
194 return kFALSE;
195 }
196 ext3 = o ? (MExtractor*)o->Clone() : NULL;
197
198 o = file.Get("ExtractTime");
199 if (o && !o->InheritsFrom(MExtractor::Class()))
200 {
201 *fLog << err << dbginf << "ERROR - ExtractTime read from " << fname << " doesn't inherit from MExtractor!" << endl;
202 return kFALSE;
203 }
204 ext2 = o ? (MExtractor*)o->Clone() : NULL;
205 if (!ext3 && !ext2)
206 {
207 *fLog << err << dbginf << "ERROR - Neither ExtractSignal nor ExrtractTime found in " << fname << "!" << endl;
208 return kFALSE;
209 }
210
211 o = file.Get("MGeomCam");
212 if (o && !o->InheritsFrom(MGeomCam::Class()))
213 {
214 *fLog << err << dbginf << "ERROR - MGeomCam read from " << fname << " doesn't inherit from MGeomCam!" << endl;
215 return kFALSE;
216 }
217 geom = o ? o->ClassName() : "";
218
219 TObjArray cont(l);
220 cont.Add(&cam);
221 return ReadContainer(cont);
222}
223
224// --------------------------------------------------------------------------
225//
226// MJCalibration allows to setup several option by a resource file:
227// MJCalibrateSignal.RawData: yes,no
228//
229// For more details see the class description and the corresponding Getters
230//
231Bool_t MJCalibrateSignal::CheckEnvLocal()
232{
233 SetInterlaced(GetEnv("Interlaced", fIsInterlaced));
234 SetRelTimesUpdate(GetEnv("RelTimesUpdate", fIsRelTimesUpdate));
235 SetMovieMode(GetEnv("MovieMode", fIsMovieMode));
236
237 return MJCalib::CheckEnvLocal();
238}
239
240Bool_t MJCalibrateSignal::Process(MPedestalCam &pedcamab, MPedestalCam &pedcambias,
241 MPedestalCam &pedcamextr)
242{
243 if (!fSequence.IsValid())
244 {
245 *fLog << err << "ERROR - Sequence invalid..." << endl;
246 return kFALSE;
247 }
248
249 // --------------------------------------------------------------------------------
250
251 *fLog << inf;
252 fLog->Separator(GetDescriptor());
253 *fLog << "Calculate calibrated data from Sequence #";
254 *fLog << fSequence.GetSequence() << endl << endl;
255
256
257 if (!CheckEnv())
258 return kFALSE;
259
260 // --------------------------------------------------------------------------------
261
262 MDirIter iter;
263 if (fSequence.IsValid())
264 {
265 if (fSequence.GetRuns(iter, fIsTestMode?MSequence::kRawCal:MSequence::kRawDat)<=0)
266 return kFALSE;
267 }
268
269 // Read File
270 MHCalibrationChargeCam hchacam;
271 MHCalibrationChargeBlindCam hbndcam;
272 MHCalibrationChargePINDiode hpndiod;
273 MHCalibrationRelTimeCam hrelcam;
274 //MHCalibrationHiLoCam hilocam;
275 //MHCalibrationPulseTimeCam hpulcam;
276
277 hchacam.SetOscillations(kFALSE);
278 hbndcam.SetOscillations(kFALSE);
279 hrelcam.SetOscillations(kFALSE);
280
281 MCalibrationChargeCam calcam;
282 MCalibrationQECam qecam;
283 MCalibrationBlindCam bndcam;
284 MCalibrationChargePINDiode pind;
285 MCalibrationRelTimeCam tmcam;
286 //MCalibrationHiLoCam hilcam;
287 //MCalibrationPulseTimeCam pulcam;
288
289 MBadPixelsCam badpix;
290
291 TObjArray interlacedcont;
292 if (fIsInterlaced)
293 {
294 interlacedcont.Add(&hchacam);
295 if (IsUseBlindPixel())
296 interlacedcont.Add(&hbndcam);
297 if (IsUsePINDiode())
298 interlacedcont.Add(&hpndiod);
299 if (fIsRelTimesUpdate)
300 interlacedcont.Add(&hrelcam);
301 }
302
303 MExtractor *extractor1=fExtractor;
304 MExtractor *extractor2=0;
305 MExtractor *extractor3=0;
306 TString geom;
307
308 TObjArray calibcont;
309 calibcont.Add(&calcam);
310 calibcont.Add(&qecam);
311 calibcont.Add(&bndcam);
312 calibcont.Add(&tmcam);
313 if (IsUseBlindPixel())
314 calibcont.Add(&bndcam);
315
316 if (!ReadCalibration(calibcont, badpix, extractor2, extractor3, geom))
317 return kFALSE;
318
319 *fLog << all;
320 if (!geom.IsNull())
321 *fLog << inf << "Camera geometry found in file: " << geom << endl;
322 else
323 *fLog << inf << "No Camera geometry found using default <MGeomCamMagic>" << endl;
324
325 if (extractor3)
326 {
327 *fLog << underline << "Signal Extractor found in calibration file" << endl;
328 extractor3->Print();
329 *fLog << endl;
330 }
331 else
332 *fLog << inf << "No Signal Extractor: ExtractSignal in file." << endl;
333
334
335 if (extractor1)
336 {
337 *fLog << underline << "Modified Signal Extractor set by user." << endl;
338 extractor1->Print();
339 *fLog << endl;
340 }
341 else
342 {
343 *fLog << inf << "No modified Signal Extractor set by user... using ExtractSignal." << endl;
344 extractor1 = extractor3 ? (MExtractor*)extractor3->Clone() : 0;
345 }
346
347 if (extractor2)
348 {
349 *fLog << underline << "Time Extractor found in calibration file" << endl;
350 extractor2->Print();
351 *fLog << endl;
352 }
353 else
354 *fLog << inf << "No Time Extractor: ExtractTime in file." << endl;
355
356 // This is necessary for the case in which it is not in the files
357 MCalibConstCam constcam;
358
359 //MBadPixelsCam badcam;
360 constcam.SetBadPixels(&badpix);
361
362 // Setup Parlist
363 MParList plist;
364 plist.AddToList(this); // take care of fDisplay!
365 plist.AddToList(&badpix);
366 plist.AddToList(&constcam);
367 //plist.AddToList(&hilcam);
368 plist.AddToList(&calibcont); // Using AddToList(TObjArray *)
369 plist.AddToList(&interlacedcont); // Using AddToList(TObjArray *)
370 //plist.AddToList(&pulcam);
371
372 // Setup Tasklist
373 MTaskList tlist;
374 plist.AddToList(&tlist);
375
376 MReadMarsFile readmc("Events");
377 readmc.DisableAutoScheme();
378
379 MRawFileRead rawread(NULL);
380 rawread.SetForceMode(); // Ignore broken time-stamps
381
382 MRead *read = fSequence.IsMonteCarlo() ? (MRead*)&readmc : (MRead*)&rawread;
383 read->AddFiles(iter);
384
385 const TString fname(Form("s/(([0-9]+_)?(M[12]_)?[0-9.]+)_[CD]_(.*[.])(raw|raw[.]gz|root)$/%s\\/$1_Y_$4root/",
386 Esc(fPathOut).Data()));
387
388 // Skips MC which have no contents. This are precisely the
389 // events which fullfilled the MC Lvl1 trigger and an
390 // arbitrary cut (typically at 50phe) to speed up simulation
391 MContinue contmc("MRawEvtData.GetNumPixels<0.5", "ContEmptyMC");
392
393 //MPointingPosInterpolate pextr;
394 //pextr.AddFiles(&iter);
395
396 MGeomApply apply; // Only necessary to create geometry
397 if (!geom.IsNull())
398 apply.SetGeometry(geom);
399 //MBadPixelsMerge merge(&badpix);
400
401 // Make sure that pedcamab has the correct name
402 pedcamab.SetName("MPedestalFundamental");
403 pedcamextr.SetName("MPedestalFromExtractorRndm");
404 pedcambias.SetName("MPedestalFromExtractor");
405 plist.AddToList(&pedcamextr);
406 plist.AddToList(&pedcambias);
407 plist.AddToList(&pedcamab);
408
409 MArrivalTimeCam timecam;
410 plist.AddToList(&timecam);
411
412 // Check for interleaved events
413 MCalibrationPatternDecode caldec;
414 MTriggerPatternDecode decode;
415
416 MH3 hpat("MRawRunHeader.GetFileID", "MTriggerPattern.GetUnprescaled");
417 if (fSequence.IsMonteCarlo())
418 hpat.SetWeight("100./MMcRunHeader.fNumSimulatedShowers");
419 else
420 hpat.SetWeight("1./TMath::Max(MRawRunHeader.GetRunLength,1)");
421 hpat.SetName("TrigPat");
422 hpat.SetTitle("Rate of the trigger pattern [Hz];File Id;Trigger Pattern;Rate [Hz]");
423 hpat.InitLabels(MH3::kLabelsXY);
424 //hpat.DefineLabelsY("1=Lvl1;2=Cal;3=Cal;4=Lvl2;5=Cal;7=Cal;8=Ped;9=Ped+Trig;13=Ped+Trig;16=Pin;32=Sum");
425 hpat.DefaultLabelY("UNKNOWN");
426 hpat.DefineLabelY( 0, "0"); // 0: No pattern
427 hpat.DefineLabelY( 1, "Trig"); // 1: Lvl1
428 hpat.DefineLabelY( 2, "Cal"); // 2: Cal
429 hpat.DefineLabelY( 3, "Cal"); // Cal+Lvl1
430 hpat.DefineLabelY( 4, "Trig"); // 4: Lvl2
431 hpat.DefineLabelY( 5, "Trig"); // Lvl2+Lvl1
432 hpat.DefineLabelY( 7, "Cal"); // Lvl2+Cal+Lvl
433 hpat.DefineLabelY( 8, "Ped"); // 8: Ped
434 hpat.DefineLabelY( 9, "Ped+Trig"); // Ped+Lvl1
435 hpat.DefineLabelY(10, "Ped+Cal"); // Ped+Cal
436 hpat.DefineLabelY(12, "Ped+Trig"); // Ped+Lvl2
437 hpat.DefineLabelY(13, "Ped+Trig"); // Ped+Lvl2+Lvl1
438 hpat.DefineLabelY(16, "Pin"); // 16: Pin
439 hpat.DefineLabelY(32, "Sum"); // 32: Sum
440 hpat.DefineLabelY(33, "Trig"); // Sum+Lvl1
441 hpat.DefineLabelY(34, "Cal"); // Sum+Cal
442 hpat.DefineLabelY(35, "Cal"); // Sum+Cal+Lvl1
443 hpat.DefineLabelY(36, "Trig"); // Sum+Lvl2
444 hpat.DefineLabelY(37, "Trig"); // Sum+Lvl1+Lvl2
445 hpat.DefineLabelY(39, "Cal"); // Sum+Lvl2+Cal+Lvl1
446
447 MFillH fillpat(&hpat, "", "FillPattern");
448 fillpat.SetDrawOption("box");
449
450 // This will make that for data with version less than 5, where
451 // trigger patterns were not yet correct, all the events in the real
452 // data file will be processed. In any case there are no interleaved
453 // calibration events in such data, so this is fine.
454 // We allow only cosmics triggered events to pass (before prescaling)
455 MFTriggerPattern fcalped("SelectTrigEvts");
456 fcalped.SetInverted();
457 fcalped.SetDefault(kTRUE);
458 fcalped.DenyAll();
459 fcalped.AllowTriggerLvl1();
460 fcalped.AllowTriggerLvl2();
461 fcalped.AllowSumTrigger();
462 if (fIsTestMode)
463 fcalped.AllowCalibration();
464
465 // This will skip interleaved events with a cal- or ped-trigger
466 MContinue contcalped(&fcalped, "ContNonTrigger");
467
468 // Create the pedestal subtracted raw-data
469 MPedestalSubtract pedsub;
470 pedsub.SetPedestalCam(&pedcamab);
471
472 // Do signal and pedestal calculation
473 MPedCalcFromLoGain pedlo1("MPedCalcFundamental");
474 pedlo1.SetPedestalUpdate(kTRUE);
475 pedlo1.SetNamePedestalCamOut("MPedestalFundamental");
476
477 MPedCalcFromLoGain pedlo2("MPedCalcWithExtractorRndm");
478 pedlo2.SetPedestalUpdate(kTRUE);
479 pedlo2.SetRandomCalculation(kTRUE);
480 pedlo2.SetNamePedestalCamOut("MPedestalFromExtractorRndm");
481
482 MPedCalcFromLoGain pedlo3("MPedCalcWithExtractor");
483 pedlo3.SetPedestalUpdate(kTRUE);
484 pedlo3.SetRandomCalculation(kFALSE);
485 pedlo3.SetNamePedestalCamOut("MPedestalFromExtractor");
486
487 if (!extractor1)
488 {
489 *fLog << err << "ERROR - extractor1 == NULL" << endl;
490 return kFALSE;
491 }
492
493 // Setup to use the hi-gain extraction window in the lo-gain
494 // range (the start of the lo-gain range is added automatically
495 // by MPedCalcFromLoGain)
496 //
497 // The window size of the extractor is not yet initialized,
498 // so we have to stick to the extraction range
499 //
500 // Even if we would like to use a range comparable to the
501 // hi-gain extraction we use the lo-gain range to make
502 // sure that exclusions (eg. due to switching noise)
503 // are correctly handled.
504 //
505 pedlo1.SetRangeFromExtractor(*extractor1);
506
507 if (extractor1->InheritsFrom("MExtractTimeAndCharge"))
508 {
509 pedlo2.SetExtractor((MExtractTimeAndCharge*)extractor1);
510 pedlo3.SetExtractor((MExtractTimeAndCharge*)extractor1);
511 }
512 else
513 {
514 pedlo2.SetRangeFromExtractor(*extractor1);
515 pedlo3.SetRangeFromExtractor(*extractor1);
516 }
517
518 //------------------------------
519 //
520 // this is the filter to find pedestal events. For the Siegen FADC
521 // these are all events which are not calibration events because
522 // the pedestal is extracted from the lo-gain signal. For MUX
523 // data this are artifiially triggered events with the pedestal
524 // trigger flag, and for safty without Lvl1 or Lvl2 flag)
525 //
526 // For the time "before" the trigger pattern all events (Siegen FADC)
527 // can be considered to be pedestal, because this was also the time
528 // without artifially calibration events
529 //
530 // Deny or allow is done before prescaling.
531 //
532 MFTriggerPattern fped("SelectPedestals");
533 fped.SetDefault(kTRUE);
534 fped.DenyCalibration();
535 if (!extractor1->HasLoGain())
536 {
537 fped.DenyAll();
538 fped.RequirePedestal();
539 }
540
541 //------------------------------
542 //
543 // Apply a filter against cosmics (this is to make sure that the
544 // trigger system was working properly and no empty events survive)
545 // For every event 5% of the pixel must not be empty. In PostProcess
546 // an error is raised if more than 50% of the events were skipped.
547 //
548 // Since unsuitable pixels are not counted in the ideal case
549 // there is no empty pixel at all.
550 //
551 MFCosmics fcosmicscal;
552 fcosmicscal.SetNamePedestalCam("MPedestalFundamental"); //CORRECT?
553 fcosmicscal.SetMaxEmptyPixels(0.05);
554 fcosmicscal.SetMaxAcceptedFraction(0.5);
555
556 MContinue contcoscal(&fcosmicscal, "ContCosmicsCal");
557
558 //------------------------------
559 //
560 // Remove events which are too bright. These events could either
561 // be calibration events or events which produce many many
562 // saturating low-gains and thus spoil the dead pixel plot.
563 // These events don't seem to be physical events so we can fairly
564 // remove them (they are no analysable gammas anyway)
565 // So we remove all events which have less than 5% empty pixels.
566 //
567 MFCosmics fcosmicsbright;
568 fcosmicsbright.SetNamePedestalCam("MPedestalFundamental");
569 fcosmicsbright.SetMaxEmptyPixels(0.10);
570 fcosmicsbright.SetMinAcceptedFraction(0.9);
571
572 MContinue contbright(&fcosmicsbright, "ContBrightEvts");
573 contbright.SetInverted();
574
575 //------------------------------
576 //
577 // Thie signal extractors
578 //
579 MTaskEnv taskenv1("ExtractSignal");
580 MTaskEnv taskenv2("ExtractTime");
581 MTaskEnv taskenv3("ExtractInterlaced");
582 taskenv1.SetDefault(extractor1);
583 taskenv2.SetDefault(extractor2);
584 taskenv3.SetDefault(extractor3);
585
586 //
587 // This is new calibration to photo-electrons, hard-coded
588 // as decided at the Wuerzburg software meeting 26.01.05
589 //
590 MCalibrateData calib;
591 calib.SetSignalType(MCalibrateData::kPhe);
592 //calib.AddPedestal("Fundamental");
593 calib.AddPedestal("FromExtractor");
594 calib.AddPedestal("FromExtractorRndm");
595 calib.SetPedestalFlag(MCalibrateData::kEvent);
596
597 //----------------------------------------------------------
598
599 MExtractPINDiode pinext;
600 MExtractBlindPixel bldext;
601
602 // Execute for all events with the calibration trigger. If no
603 // trigger pattern is available do not execute it
604 // The selection is done before prescaling.
605 MFTriggerPattern fcalib("SelectCalEvts");
606 fcalib.SetDefault(kFALSE);
607 fcalib.DenyAll();
608 fcalib.RequireCalibration();
609 fcalib.AllowTriggerLvl1();
610 fcalib.AllowTriggerLvl2();
611 fcalib.AllowSumTrigger();
612
613 MCalibrationChargeCalc chcalc;
614 chcalc.SetContinousCalibration();
615 chcalc.SetExtractor(extractor3);
616
617 MCalibrationRelTimeCalc recalc;
618 MCalibCalcFromPast pacalc;
619
620 chcalc.SetPedestals(&pedcamextr);
621
622 pacalc.SetChargeCalc(&chcalc);
623 if (fIsRelTimesUpdate)
624 pacalc.SetRelTimeCalc(&recalc);
625 pacalc.SetCalibrate(&calib);
626
627 //
628 // Calibration histogramming
629 //
630 MFillH filpin(&hpndiod, "MExtractedSignalPINDiode", "FillPINDiode");
631 MFillH filbnd(&hbndcam, "MExtractedSignalBlindPixel", "FillBlindCam");
632 MFillH filcam(&hchacam, "MExtractedSignalCam", "FillChargeCam");
633 MFillH filtme(&hrelcam, "MArrivalTimeCam", "FillRelTime");
634 //MFillH filhil(&hilocam, "MExtractedSignalCam", "FillHiLoRatio");
635 //MFillH filpul(&hpulcam, "MRawEvtData", "FillPulseTime");
636 filpin.SetBit(MFillH::kDoNotDisplay);
637 filbnd.SetBit(MFillH::kDoNotDisplay);
638 filcam.SetBit(MFillH::kDoNotDisplay);
639 filtme.SetBit(MFillH::kDoNotDisplay);
640 //filhil.SetBit(MFillH::kDoNotDisplay);
641 //filpul.SetBit(MFillH::kDoNotDisplay);
642
643 MCalibrateRelTimes caltm;
644 MBadPixelsCalc bpcal;
645 MBadPixelsTreat treat;
646
647 //bpcal.SetNamePedPhotCam("MPedPhotFromExtractor");
648 bpcal.SetNamePedPhotCam("MPedPhotFromExtractorRndm");
649
650 //treat.AddNamePedPhotCam("MPedPhotFundamental");
651 treat.AddNamePedPhotCam("MPedPhotFromExtractor");
652 treat.AddNamePedPhotCam("MPedPhotFromExtractorRndm");
653 if (!extractor2 && !extractor1->InheritsFrom("MExtractTimeAndCharge"))
654 treat.SetProcessTimes(kFALSE);
655
656 MHCamEvent evt0( 0, "PedFLG", "Fundamental Pedestal from Lo Gain;;P [cnts/sl]");
657 MHCamEvent evt1( 2, "PedRmsFLG", "RMS from Extractor applied to ped.;;\\sigma_{p} [cnts/sl]");
658 MHCamEvent evt2( 0, "Extra'd", "Extracted Signal;;S [cnts/sl]");
659 // MHCamEvent evt3(4, "PedPhot", "Calibrated Pedestal;;P [phe]");
660 MHCamEvent evt4( 5, "PedRMS", "Calibrated RMS from Extractor applied to ped.;;\\sigma_{p} [phe]");
661 MHCamEvent evt5( 0, "Interp'd", "Interpolated Signal scaled with A/A_{0};;S [phe]");
662 MHCamEvent evt6(102, "Unsuitable", "Fraction of unsuitable events per Pixel;;[1]");
663 // MHCamEvent evt7( 6, "Times", "Calibrated Arrival Time;;T [fadc sl]");
664 MHCamEvent evt8( 0, "Conv", "Calibration Conv. Factors;;[phe/cnts]");
665 MHCamEvent evt9( 7, "PulsePos", "Pulse Position of cosmics (>50phe);;T [ns]");
666 MHCamEvent evtR( 4, "HiLoCal", "Hi-/Lo-Gain ratio;;Ratio");
667 MHCamEvent evtO( 7, "HiLoOff", "Lo-/Hi-Gain Offset;;Offset");
668 MHCamEvent evtC( 4, "CalPos", "Extracted pulse position of calibration pulses;;T [sl]");
669 //MHCamEvent evt2(0, "Extra'd", "Extracted Calibration Signal;;S [cnts/sl]");
670 evt2.SetErrorSpread(kFALSE);
671 evt5.SetErrorSpread(kFALSE);
672 evt6.SetErrorSpread(kFALSE);
673 evt6.SetThreshold();
674
675 MFillH fill0(&evt0, "MPedestalFundamental", "FillPedFLG");
676 MFillH fill1(&evt1, "MPedestalFromExtractorRndm", "FillPedRmsFLG");
677 MFillH fill2(&evt2, "MExtractedSignalCam", "FillExtracted");
678 // MFillH fill3(&evt3, "MPedPhotFundamental", "FillPedPhot");
679 MFillH fill4(&evt4, "MPedPhotFromExtractorRndm", "FillPedRMS");
680 MFillH fill5(&evt5, "MSignalCam", "FillInterpolated");
681 MFillH fill6(&evt6, "MBadPixelsCam", "FillUnsuitable");
682 // MFillH fill7(&evt7, "MSignalCam", "FillTimes");
683 MFillH fill8(&evt8, "MCalibConstCam", "FillConv");
684 MFillH fill9(&evt9, "MSignalCam", "FillPulsePos");
685 MFillH fillR(&evtR, "MExtractedSignalCam", "FillHiLoCal");
686 MFillH fillO(&evtO, "MArrivalTimeCam", "FillHiLoOff");
687 MFillH fillC(&evtC, "MArrivalTimeCam", "FillCalPos");
688 //MFillH fill2(&evt2, "MExtractedSignalCam", "FillExtractedSignal");
689
690 MHVsTime histbp("MBadPixelsCam.GetNumUnsuitable");
691 histbp.SetName("BadPixTm");
692 histbp.SetTitle("Number of unsuitable pixels;;N");
693 histbp.SetMinimum(0);
694
695 MHVsTime histdp("MSignalCam.GetNumPixelsUnmapped");
696 histdp.SetName("DeadPixTm");
697 histdp.SetTitle("Number of dead/unmapped pixels;;N");
698 histdp.SetMinimum(0);
699
700 // Task to fill the histogram
701 MFillH fillB(&histbp, "MTime", "FillBadPixTm");
702 MFillH fillD(&histdp, "MTime", "FillDeadPixTm");
703 fillB.SetNameTab("BadPixTm");
704 fillD.SetNameTab("DeadPixTm");
705
706 /*
707 MFillH fillP("MHPulseShape", "", "FillPulseShape");
708 fillP.SetNameTab("Pulse");
709 */
710
711 /*
712 MHVsTime hbadpix("MBadPixelsCam.GetNumUnsuitable");
713 hbadpix.SetNumEvents(50);
714 MFillH fillB(&hbadpix, "MTime");
715 */
716
717 MTaskEnv fillflorian("FinalFantasy");
718 fillflorian.SetDefault();
719
720 // The second rule is for the case reading raw-files!
721 MWriteRootFile write(2, fname, fOverwrite?"RECREATE":"NEW", "Calibrated Data");
722 // Run Header
723 write.AddContainer("MRawRunHeader", "RunHeaders");
724// write.AddContainer("MBadPixelsCam", "RunHeaders");
725 write.AddContainer("MGeomCam", "RunHeaders");
726 // Monte Carlo Headers
727 write.AddContainer("MMcRunHeader", "RunHeaders", kFALSE);
728 write.AddContainer("MMcFadcHeader", "RunHeaders", kFALSE);
729 write.AddContainer("MMcTrigHeader", "RunHeaders", kFALSE);
730 write.AddContainer("MMcConfigRunHeader", "RunHeaders", kFALSE);
731 write.AddContainer("MMcCorsikaRunHeader", "RunHeaders", kFALSE);
732 write.AddContainer("MCorsikaRunHeader", "RunHeaders", kFALSE);
733 // Monte Carlo
734 write.AddContainer("MMcEvt", "Events", kFALSE);
735 write.AddContainer("MMcTrig", "Events", kFALSE);
736 write.AddContainer("MCorsikaEvtHeader", "Events", kFALSE);
737 // Data tree
738 write.AddContainer("MSignalCam", "Events");
739 // write.AddContainer("MPedPhotFundamental", "Events");
740 write.AddContainer("MPedPhotFromExtractor", "Events");
741 write.AddContainer("MPedPhotFromExtractorRndm", "Events");
742 write.AddContainer("MTime", "Events", kFALSE);
743 write.AddContainer("MRawEvtHeader", "Events");
744 write.AddContainer("MTriggerPattern", "Events");
745
746 // Trees with slow-control information (obsolete
747 // if we don't read merpped root-files)
748 write.AddTree("Trigger", kFALSE);
749 write.AddTree("Drive", kFALSE);
750 write.AddTree("CC", kFALSE);
751 write.AddTree("Pyrometer", kFALSE);
752 write.AddTree("Currents", kFALSE);
753 write.AddTree("Camera", kFALSE);
754 // Slow-Control: Current-tree
755 write.AddContainer("MCameraDC", "Currents", kFALSE);
756 // Slow-Control: Camera-tree
757 write.AddContainer("MCameraAUX", "Camera", kFALSE);
758 write.AddContainer("MCameraCalibration", "Camera", kFALSE);
759 write.AddContainer("MCameraCooling", "Camera", kFALSE);
760 write.AddContainer("MCameraHV", "Camera", kFALSE);
761 write.AddContainer("MCameraLV", "Camera", kFALSE);
762 write.AddContainer("MCameraLids", "Camera", kFALSE);
763
764 // Write the special MC tree
765 MWriteRootFile writemc(2, fname, fOverwrite?"RECREATE":"NEW", "Calibrated Data");
766 writemc.SetName("WriteMC");
767 writemc.AddContainer("MMcEvtBasic", "OriginalMC", kFALSE);
768 if (fSequence.IsMonteCarlo())
769 writemc.AddCopySource("OriginalMC", kFALSE);
770
771 // Write the special calib tree
772 /*
773 MWriteRootFile writecal(2, fname, fOverwrite?"RECREATE":"NEW");
774 writecal.SetName("WriteCalib");
775 writecal.AddContainer("MBadPixelsCam", "Calib");
776 writecal.AddContainer("MCalibrationChargeCam", "Calib");
777 writecal.AddContainer("MCalibrationRelTimeCam", "Calib");
778 */
779
780 //-----------------------------------------------------------
781 // Build tasklist
782
783 MTaskList tlist2("AllEvents");
784
785 tlist2.AddToList(&caldec);
786 tlist2.AddToList(&decode);
787 tlist2.AddToList(&fillpat);
788 tlist2.AddToList(&apply);
789 //tlist2.AddToList(&merge);
790 tlist2.AddToList(&pedsub);
791
792 //-----------------------------------------------------------
793 // Pedestal extraction
794
795 MTaskList tlist3("PedEvents");
796 tlist3.SetFilter(&fped); // Deny events with cal-trigger
797
798 tlist2.AddToList(&fped); // If no lo-gain require ped-trigger
799 tlist2.AddToList(&tlist3); // and deny cosmics (lvl1/2) trigger
800
801 tlist3.AddToList(&pedlo1); // extract pedestal events
802 tlist3.AddToList(&pedlo2); // extract pedestal events
803 tlist3.AddToList(&pedlo3); // extract pedestal events
804 tlist3.AddToList(&fill0); // fill pedestal events
805 tlist3.AddToList(&fill1); // fill pedestal events
806
807 //-----------------------------------------------------------
808 // Calibration
809
810 MTaskList tlist4("CalEvents");
811 tlist4.SetFilter(&fcalib); // process only events with cal-trigger
812
813 //MFDataPhrase filcalco("MCalibrationConstCam.IsReadyToSave>0.5", "CalibConstFilter");
814 if (fIsInterlaced && !fIsTestMode)
815 {
816 // The task list is executed for all events with the calibration
817 // trigger
818 tlist2.AddToList(&fcalib); // MFTriggerPattern
819 tlist2.AddToList(&tlist4);
820
821 tlist4.AddToList(&taskenv3);
822 tlist4.AddToList(&contcoscal); // MContinue/ContCosmicsCal
823 if (IsUsePINDiode())
824 tlist4.AddToList(&pinext); // MExtractPINDiode
825 if (IsUseBlindPixel())
826 tlist4.AddToList(&bldext); // MExtractBlindPixel
827 tlist4.AddToList(&pacalc); // MCalibCalcFromPast
828 /*
829 tlist3.AddToList(&filcalco); // CalibConstFilter (IsReadyToSave)
830 fill8.SetFilter(&filcalco);
831 tlist3.AddToList(&fill8); // FillConvUpd
832 */
833
834 tlist4.AddToList(&filcam); // FillChargeCam
835 if (fIsRelTimesUpdate)
836 tlist4.AddToList(&filtme); // FillRelTime
837 if (IsUseBlindPixel())
838 tlist4.AddToList(&filbnd); // FillBlindCam
839 if (IsUsePINDiode())
840 tlist4.AddToList(&filpin); // FillPINDiode
841 tlist4.AddToList(&chcalc); // MCalibrationChargeCalc
842 if (fIsRelTimesUpdate)
843 tlist4.AddToList(&recalc); // MCalibrationRelTimeCam
844
845 tlist4.AddToList(&fillC); // FillCalPos
846
847 //tlist3.AddToList(&writecal); // MWriteRootFile
848 }
849
850 //-----------------------------------------------------------
851 // Cosmics extraction
852
853 // remove all events with a cal- or ped-trigger (no matter
854 // whether they have lvl1 or lvl2 or any other flag
855 tlist2.AddToList(&contcalped);
856
857 // Extract the signal
858 if (extractor1)
859 tlist2.AddToList(&taskenv1);
860
861 // remove all events which definitly don't have a signal
862 // using MFCosmics (ContCosmicsPed)
863 // tlist2.AddToList(&contcosped);
864
865 // Extract arrival time (if a dedicated extrator given)
866 if (extractor2)
867 tlist2.AddToList(&taskenv2);
868
869 // Check if we have an extremely bright event (remove them,
870 // they are presumably errornous events or calibration events)
871 if (!fIsTestMode)
872 tlist2.AddToList(&contbright);
873
874 /*
875 if (fIsHiLoCalibration)
876 {
877 plist.AddToList(&hilocam);
878 tlist2.AddToList(&filhil);
879 }
880
881 if (fIsPulsePosCheck)
882 {
883 plist.AddToList(&hpulcam);
884 tlist2.AddToList(&filpul);
885 }
886 */
887
888 tlist2.AddToList(&fill2);
889 tlist2.AddToList(&fill8); // FillConv
890 tlist2.AddToList(&calib); // MCalibrateData
891 if (extractor2 || extractor1->InheritsFrom("MExtractTimeAndCharge"))
892 tlist2.AddToList(&caltm);
893
894 tlist2.AddToList(&bpcal); // MBadPixelsCalc
895 tlist2.AddToList(&treat); // MBadPixelsTreat
896 tlist2.AddToList(&fill6);
897 // tlist2.AddToList(&fill3);
898 tlist2.AddToList(&fill4);
899 tlist2.AddToList(&fill5);
900 //if (extractor2 || extractor1->InheritsFrom("MExtractTimeAndCharge"))
901 // tlist2.AddToList(&fill7);
902 tlist2.AddToList(&fill9);
903 if (!fSequence.IsMonteCarlo())
904 {
905 tlist2.AddToList(&fillB);
906 tlist2.AddToList(&fillD);
907 }
908 if (extractor1->HasLoGain())
909 {
910 tlist2.AddToList(&fillR);
911 tlist2.AddToList(&fillO);
912 }
913
914 /*
915 MFillH fillC("MHCleaning", "", "FillClean");
916 tlist2.AddToList(&fillC);
917
918 //tlist2.AddToList(&fillP);
919 */
920
921 // ----- Start: Code for encoding movies -----
922
923 MMoviePrepare movprep;
924 MMovieWrite movwrite;
925 movprep.SetRangeFromExtractor(*extractor1);
926
927 //MFDataPhrase movfilt("MMovieData.fMax>150");
928 MFDataPhrase movfilt("MMovieData.fMax>5*MMovieData.fMedianPedestalRms", "MovieFilter");
929
930 MImgCleanStd movclean(8.5, 4.0);
931 movclean.SetMethod(MImgCleanStd::kAbsolute);
932
933 //movprep.SetFilter(&evtnum);
934 movclean.SetFilter(&movfilt);
935 movwrite.SetFilter(&movfilt);
936
937 MTaskList tlistmov("MovieEncoder");
938 tlistmov.AddToList(&movprep);
939 tlistmov.AddToList(&movfilt);
940 tlistmov.AddToList(&movclean);
941 tlistmov.AddToList(&movwrite);
942
943 MFEvtNumber evtnum;
944 //evtnum.SetSelector("ThetaSquared.fVal<0.04");
945 //evtnum.SetFileName("ganymed00000001.root");
946 tlistmov.SetFilter(&evtnum);
947
948 if (fIsMovieMode)
949 {
950 tlist2.AddToList(&evtnum);
951 tlist2.AddToList(&tlistmov);
952 }
953
954 // ----- End: Code for encoding movies -----
955
956 tlist2.AddToList(&fillflorian);
957
958 // Setup List for Drive-tree
959 //MPointingPosCalc pcalc;
960
961 // Now setup main tasklist
962 tlist.AddToList(read);
963
964 if (fSequence.IsMonteCarlo())
965 {
966 if (!fIsMovieMode && !HasNullOut())
967 tlist.AddToList(&writemc);
968 tlist.AddToList(&contmc);
969 }
970
971 //if (IsUseRootData())
972 // tlist2.AddToList(&pextr);
973 tlist.AddToList(&tlist2, fSequence.IsMonteCarlo() ? "Events" : "All");
974
975 //if (fSequence.IsMonteCarlo())
976 // tlist.AddToList(&pcalc, "Drive");
977
978 if (!fIsMovieMode && !HasNullOut())
979 tlist.AddToList(&write);
980
981 // Create and setup the eventloop
982 MEvtLoop evtloop(fName);
983 evtloop.SetParList(&plist);
984 evtloop.SetDisplay(fDisplay);
985 evtloop.SetLogStream(fLog);
986 if (!SetupEnv(evtloop))
987 return kFALSE;
988
989 // Execute first analysis
990 const Bool_t rc = evtloop.Eventloop(fMaxEvents);
991
992 // make sure owned object are deleted
993 if (extractor1 && extractor1!=fExtractor)
994 delete extractor1;
995 if (extractor2)
996 delete extractor2;
997 if (extractor3)
998 delete extractor3;
999
1000 // return if job failed
1001 if (!rc)
1002 {
1003 *fLog << err << GetDescriptor() << ": Failed." << endl;
1004 return kFALSE;
1005 }
1006
1007 // if everything went ok write and display result
1008 DisplayResult(plist);
1009
1010 /*
1011 if (fIsPixelCheck)
1012 {
1013 if (fIsPulsePosCheck)
1014 hpulcam[fCheckedPixId].DrawClone("");
1015
1016 //if (fIsHiLoCalibration)
1017 // hilocam[fCheckedPixId].DrawClone("");
1018 }
1019 interlacedcont.Add(&pulcam);
1020 */
1021
1022 //if (fIsHiLoCalibration)
1023 // interlacedcont.Add(&hilcam);
1024
1025 //if (fIsPulsePosCheck)
1026 // interlacedcont.Add(plist.FindObject("MHCalibrationPulseTimeCam"));
1027
1028 //if (fIsHiLoCalibration)
1029 // interlacedcont.Add(plist.FindObject("MHCalibrationHiLoCam"));
1030
1031 if (!WriteResult())
1032 return kFALSE;
1033
1034 // return if job went ok
1035 *fLog << all << GetDescriptor() << ": Done." << endl;
1036 *fLog << endl << endl;
1037
1038 return kTRUE;
1039}
1040
1041
1042void MJCalibrateSignal::DisplayResult(MParList &plist)
1043{
1044 /*
1045 if (!fDisplay || !fIsHiLoCalibration)
1046 return;
1047
1048 MCalibrationHiLoCam *hcam = (MCalibrationHiLoCam*)plist.FindObject("MCalibrationHiLoCam");
1049 MGeomCam *geom = (MGeomCam*)plist.FindObject("MGeomCam");
1050 if (!hcam || !geom)
1051 return;
1052
1053 // Create histograms to display
1054 MHCamera disp1(*geom, "HiLoConv", "Ratio Amplification HiGain vs. LoGain (Charges)");
1055 MHCamera disp2(*geom, "HiLoDiff", "Arrival Time Diff. HiGain vs. LoGain (Times)");
1056
1057 disp1.SetCamContent(*hcam, 0);
1058 disp1.SetCamError( *hcam, 1);
1059 disp2.SetCamContent(*hcam, 5);
1060 disp2.SetCamError( *hcam, 6);
1061
1062 disp1.SetYTitle("R [1]");
1063 disp2.SetYTitle("\\Delta T [FADC sl.]");
1064
1065 TCanvas &c1 = fDisplay->AddTab("HiLoConv");
1066 c1.Divide(2,3);
1067
1068 disp1.CamDraw(c1, 1, 2, 1);
1069 disp2.CamDraw(c1, 2, 2, 1);
1070*/
1071}
1072
Note: See TracBrowser for help on using the repository browser.