source: trunk/Mars/fact/analysis/mc/callisto.C@ 17699

Last change on this file since 17699 was 17699, checked in by tbretz, 11 years ago
Improved the calling convention. There is really no need to give so many arguments. Pipe some Monte Carlo values through. I wonder what special MC macros are for which omit all the Monte Carlo data in their output.
File size: 21.7 KB
Line 
1#include <sstream>
2#include <iostream>
3
4#include "MLog.h"
5#include "MLogManip.h"
6
7#include "TH1F.h"
8#include "TFile.h"
9#include "TStyle.h"
10#include "TGraph.h"
11#include "TLine.h"
12
13#include "../mcore/DrsCalib.h"
14#include "MDrsCalibration.h"
15#include "MExtralgoSpline.h"
16#include "MSequence.h"
17#include "MStatusArray.h"
18#include "MHCamera.h"
19#include "MJob.h"
20#include "MWriteRootFile.h"
21#include "MHCamera.h"
22#include "MBadPixelsCam.h"
23#include "MBadPixelsPix.h"
24#include "MDirIter.h"
25#include "MTaskList.h"
26#include "MFDataPhrase.h"
27#include "MArrayF.h"
28#include "MBadPixelsTreat.h"
29#include "MCalibrateDrsTimes.h"
30#include "MHSectorVsTime.h"
31#include "MHCamEvent.h"
32#include "MExtractTimeAndChargeSpline.h"
33#include "MFillH.h"
34#include "MDrsCalibApply.h"
35#include "MGeomApply.h"
36#include "MContinue.h"
37#include "MRawFitsRead.h"
38#include "MEvtLoop.h"
39#include "MParList.h"
40#include "MStatusDisplay.h"
41#include "MDrsCalibrationTime.h"
42#include "MH3.h"
43#include "MGeomCamFACT.h"
44#include "MCalibrateFact.h"
45#include "MParameters.h"
46#include "MWriteAsciiFile.h"
47
48using namespace std;
49
50/* Maybe you wanna use this macro like this:
51 *
52 * 0.) ---- call root ----
53 * root -b
54 *
55 * 1.) ---- compile the stuff ----
56 * .L fact/analysis/callisto_buildable_no_sequence_file.C++
57 * <read a lot of warnings>
58 *
59 * 2.) ---- you can call it then ----
60 * Therefore you need to specify all the paths ... see below.
61 *
62 * When you wanna call the stuff directly from the bash make sure to
63 * escape the bracets and quotes correctly.
64 *
65 * your can do:
66 * root -b -q callisto_buildable_no_sequence_file.C++'("path1","path2",...)'
67 * or:
68 * root -b -q callisto_buildable_no_sequence_file.C++(\"path1\",\"$HOME\",...)
69 * using bash enviroment variables like $HOME is not possible in the upper variant.
70 */
71
72int callisto(const TString drsfile="test300samples.drs.fits",
73 const TString pedfile="00000001.001_P_MonteCarlo000_Events.fits",
74 const TString datfile="00000003.387_D_MonteCarlo010_Events.fits",
75 TString outfile = "",
76 TString displayfile = "", TString displaytitle = "")
77{
78
79 // ======================================================
80
81 if (displaytitle.IsNull())
82 displaytitle = gSystem->BaseName(datfile);
83
84 FileStat_t fstat;
85 int rc = gSystem->GetPathInfo(outfile, fstat);
86 bool isdir = !rc || R_ISDIR(fstat.fMode);
87
88 const char *buf = gSystem->ConcatFileName(outfile, "callisto.root");
89 outfile = buf;
90 delete [] buf;
91
92 if (displayfile.IsNull())
93 {
94 displayfile = outfile;
95 displayfile.Insert(displayfile.Last('.'), "-display");
96 }
97 else
98 {
99 if (isdir && gSystem->DirName(displayfile)==TString("."))
100 {
101 buf = gSystem->ConcatFileName(outfile, displayfile);
102 displayfile = buf;
103 delete [] buf;
104 }
105 }
106
107 // ======================================================
108
109 // true: Display correctly mapped pixels in the camera displays
110 // but the value-vs-index plot is in software/spiral indices
111 // false: Display pixels in hardware/linear indices,
112 // but the order is the camera display is distorted
113 bool usemap = true;
114
115 // map file to use (get that from La Palma!)
116 const char *map = usemap ? "TestForThomas/FACT/FACTmap111030.txt" : NULL;
117
118 Bool_t maximum = kTRUE;
119
120 //const char *lp_template = maximum ?
121 // "/cm/shared/apps/fact/Mars_svn_LP/template-lp-extractor-maximum.root" :
122 // "/cm/shared/apps/fact/Mars_svn_LP/template-lp-extractor-leading-edge.root";
123
124 const char *pulse_template = "TestForThomas/FACT/template-pulse.root";
125
126 // ------------------------------------------------------
127
128 // Calib: 51 / 90 / 197 (20% TH)
129 // Data: 52 / 64 / 104 (20% TH)
130
131 // Extraction range in slices. It will always(!) contain the full range
132 // of integration
133 const int first_slice = 20; // 10ns
134 const int last_slice = 250; // 125ns
135
136 // Note that rise and fall time mean different things whether you use IntegralFixed or IntegraRel:
137 //
138 // IntegralFixed:
139 // * fRiseTime: Number of slices left from arrival time
140 // * fFallTime: Number of slices right from arrival time
141 // IntegralRel:
142 // * fRiseTime: Number of slices left from maximum time
143 // * fFallTime: Number of slices right from maximum time
144 //
145 const int rise_time_cal = maximum ? 40 : 10; // was 13; 5ns
146 const int fall_time_cal = maximum ? 120 : 160; // was 23; 80ns
147
148 const int rise_time_dat = maximum ? 10 : 2; // was 13; was 10; 1ns
149 const int fall_time_dat = maximum ? 40 : 48; // was 23; was 40; 24ns
150
151 // Extraction type: Extract integral and half leading edge
152
153 const MExtralgoSpline::ExtractionType_t type = maximum ? (MExtralgoSpline::kIntegralRel) : (MExtralgoSpline::kIntegralFixed);
154 //const int type = MExtralgoSpline::kIntegralFixed;
155
156
157 const double heighttm = 0.5; // IntegralAbs { 1.5pe * 9.6mV/pe } / IntegralRel { 0.5 }
158
159 Long_t max = 0; // All
160 Long_t max3 = max; // Pedestal Rndm
161 Long_t max4 = max; // Pedestal Ext
162
163 // ======================================================
164
165 if (map && gSystem->AccessPathName(map, kFileExists))
166 {
167 gLog << err << "ERROR - Cannot access mapping file '" << map << "'" << endl;
168 return 1;
169 }
170
171 gLog.Separator("Callisto");
172 gLog << all;
173 gLog << "Data File: " << datfile << '\n';
174 gLog << "DRS calib 300: " << drsfile << endl;;
175
176 MDrsCalibration drscalib300;
177 if (!drscalib300.ReadFits(drsfile.Data())) {
178 gLog << err << "ERROR - Cannot access drscallib300 file '" << drsfile << "'" << endl;
179 return 5;
180 }
181 gLog << all;
182 gLog << "Pedestal file: " << pedfile << '\n';
183 gLog << "Output file: " << outfile << '\n';
184 gLog << "Display file: " << displayfile << '\n';
185 gLog << "Display title: " << displaytitle << endl;
186
187 // ------------------------------------------------------
188 MStatusArray arrt, arrp;
189
190 // TFile ft(lp_template);
191 // if (arrt.Read()<=0)
192 // {
193 // gLog << err << "ERROR - Reading LP template from " << lp_template << endl;
194 // return 100;
195 // }
196
197 // MHCamera *lpref = (MHCamera*)arrt.FindObjectInCanvas("ExtCalSig;avg", "MHCamera", "Cam");
198 // if (!lpref)
199 // {
200 // gLog << err << "ERROR - LP Template not found in " << lp_template << endl;
201 // return 101;
202 // }
203 // lpref->SetDirectory(0);
204
205 // MHCamera *gain = (MHCamera*)arrt.FindObjectInCanvas("gain", "MHCamera", "Gain");
206 // if (!gain)
207 // {
208 // gLog << err << "ERROR - Gain not found in " << lp_template << endl;
209 // return 101;
210 // }
211 // gain->SetDirectory(0);
212
213 TFile fp(pulse_template);
214 if (arrp.Read()<=0)
215 {
216 gLog << err << "ERROR - Reading Pulse template from " << pulse_template << endl;
217 return 102;
218 }
219
220 TH1F *hpulse = (TH1F*)arrp.FindObjectInCanvas("hPixelEdgeMean0_0", "TH1F", "cgpPixelPulses0");
221 if (!hpulse)
222 {
223 gLog << err << "ERROR - Pulse Template not found in " << pulse_template << endl;
224 return 103;
225 }
226 hpulse->SetDirectory(0);
227 // ======================================================
228
229 MStatusDisplay *d = new MStatusDisplay;
230
231 MBadPixelsCam badpixels;
232 badpixels.InitSize(1440);
233 /*
234 badpixels[ 424].SetUnsuitable(MBadPixelsPix::kUnsuitable);
235 badpixels[ 583].SetUnsuitable(MBadPixelsPix::kUnsuitable);
236 badpixels[ 830].SetUnsuitable(MBadPixelsPix::kUnsuitable);
237 badpixels[ 923].SetUnsuitable(MBadPixelsPix::kUnsuitable);
238 badpixels[1208].SetUnsuitable(MBadPixelsPix::kUnsuitable);
239 badpixels[1399].SetUnsuitable(MBadPixelsPix::kUnsuitable);
240 */
241 // Twin pixel
242 // 113
243 // 115
244 // 354
245 // 423
246 // 1195
247 // 1393
248
249 MDrsCalibrationTime timecam;
250
251 // Plot the trigger pattern rates vs. run-number
252 MH3 hrate("MRawRunHeader.GetFileID", "MRawEvtHeader.GetTriggerID&0xff00");
253 hrate.SetWeight("1./TMath::Max(MRawRunHeader.GetRunLength,1)");
254 hrate.SetName("Rate");
255 hrate.SetTitle("Event rate [Hz];File Id;Trigger Type;");
256 hrate.InitLabels(MH3::kLabelsXY);
257 hrate.DefineLabelY( 0, "Data"); // What if TriggerID==0 already???
258 hrate.DefineLabelY(0x100, "Cal");
259 hrate.DefineLabelY(0x400, "Ped");
260 // hrate.DefaultLabelY("ERROR");
261 gStyle->SetOptFit(kTRUE);
262
263
264 // ========================= Result ==================================
265
266 //~ Double_t avgS = evt1f.GetHist()->GetMean();
267 //~ Double_t medS = evt1f.GetHist()->GetMedian();
268 //~ Double_t rmsS = evt1f.GetHist()->GetRMS();
269 //~ Double_t maxS = evt1f.GetHist()->GetMaximum();
270
271 MArrayF der1(hpulse->GetNbinsX());
272 MArrayF der2(hpulse->GetNbinsX());
273
274 MExtralgoSpline spline(hpulse->GetArray()+1, hpulse->GetNbinsX(),
275 der1.GetArray(), der2.GetArray());
276 spline.SetRiseFallTime(rise_time_dat, fall_time_dat);
277 spline.SetExtractionType(type);
278 spline.SetHeightTm(heighttm);
279
280 spline.Extract(hpulse->GetMaximumBin()-1);
281
282 // The pulser signal is most probably around 400mV/9.5mV
283 // IntegraFixed 2/48 corresponds to roughly 215mV*50slices
284 Double_t scale = 1./spline.GetSignal();
285
286 MArrayD calib(1440);
287 for (int i=0; i<1440; i++)
288 calib[i] =1.;
289
290 gROOT->SetSelectedPad(0);
291 d->AddTab("PulseTemp");
292 gPad->SetGrid();
293 hpulse->SetNameTitle("Pulse", "Single p.e. pulse template");
294 hpulse->SetDirectory(0);
295 hpulse->SetLineColor(kBlack);
296 hpulse->DrawCopy();
297
298 TAxis *ax = hpulse->GetXaxis();
299
300 Double_t w = hpulse->GetBinWidth(1);
301 Double_t T = w*(spline.GetTime()+0.5) +ax->GetXmin();
302 //~ Double_t H = w*(hpulse->GetMaximumBin()+0.5)+ax->GetXmin();
303
304 TLine line;
305 line.SetLineColor(kRed);
306 line.DrawLine(T-rise_time_dat*w, spline.GetHeight(),
307 T+fall_time_dat*w, spline.GetHeight());
308 line.DrawLine(T, spline.GetHeight()/4, T, 3*spline.GetHeight()/4);
309 line.DrawLine(T-rise_time_dat*w, 0,
310 T-rise_time_dat*w, spline.GetHeight());
311 line.DrawLine(T+fall_time_dat*w, 0,
312 T+fall_time_dat*w, spline.GetHeight());
313
314 TGraph gg;
315 for (int ix=1; ix<=hpulse->GetNbinsX(); ix++)
316 for (int i=0; i<10; i++)
317 {
318 Double_t x = hpulse->GetBinLowEdge(ix)+i*hpulse->GetBinWidth(ix)/10.;
319 gg.SetPoint(gg.GetN(), x+w/2, spline.EvalAt(ix-1+i/10.));
320 }
321
322 gg.SetLineColor(kBlue);
323 gg.SetMarkerColor(kBlue);
324 gg.SetMarkerStyle(kFullDotMedium);
325 gg.DrawClone("L");
326
327 gROOT->SetSelectedPad(0);
328 d->AddTab("CalConst");
329 MGeomCamFACT fact;
330 MHCamera hcalco(fact);
331 hcalco.SetName("CalConst");
332 hcalco.SetTitle(Form("Relative calibration constant [%.0f/pe]", 1./scale));
333 hcalco.SetCamContent(calib);
334 hcalco.SetAllUsed();
335 //hcalco.Scale(scale);
336 hcalco.DrawCopy();
337
338 // ======================================================
339
340 gLog << endl;
341 gLog.Separator("Extracting random pedestal");
342
343 MTaskList tlist3;
344
345 MParList plist3;
346 plist3.AddToList(&tlist3);
347 plist3.AddToList(&drscalib300);
348 plist3.AddToList(&badpixels);
349 plist3.AddToList(&timecam);
350
351 MEvtLoop loop3("DetermineRndmPed");
352 loop3.SetDisplay(d);
353 loop3.SetParList(&plist3);
354
355 // ------------------ Setup the tasks ---------------
356
357 MRawFitsRead read3;
358 read3.LoadMap(map);
359 read3.AddFile(pedfile);
360
361 MFillH fill3a(&hrate);
362
363 MContinue cont3("(MRawEvtHeader.GetTriggerID&0xff00)!=0x400", "SelectPed");
364
365 MGeomApply apply3;
366
367 MDrsCalibApply drsapply3;
368
369 //---
370
371 MExtractTimeAndChargeSpline extractor3;
372 extractor3.SetRange(first_slice, last_slice);
373 extractor3.SetRiseTimeHiGain(rise_time_dat);
374 extractor3.SetFallTimeHiGain(fall_time_dat);
375 extractor3.SetHeightTm(heighttm);
376 extractor3.SetChargeType(type);
377 extractor3.SetSaturationLimit(600000);
378 extractor3.SetNoiseCalculation(kTRUE);
379
380// MHCamEvent evt2a(0, "PedRdm", "Extracted Pedestal Signal;;S");
381
382// MFillH fill2a(&evt2a, "MExtractedSignalCam", "FillPedRndm");
383
384 // Use this for data, but not for calibration events
385// evt2a.SetErrorSpread(kFALSE);
386
387 /*
388 MCalibrateData conv3;
389 conv3.SetCalibrationMode(MCalibrateData::kNone);
390 conv3.SetPedestalFlag(MCalibrateData::kNo);
391 conv3.SetCalibConvMinLimit(0);
392 conv3.SetCalibConvMaxLimit(10000);
393 conv3.SetScaleFactor(scale);
394 */
395
396 MCalibrateFact conv3;
397 conv3.SetScale(scale);
398 conv3.SetCalibConst(calib);
399
400 MBadPixelsTreat treat3;
401 treat3.SetProcessPedestalRun(kFALSE);
402 treat3.SetProcessPedestalEvt(kFALSE);
403 treat3.SetProcessTimes(kFALSE);
404
405 MHCamEvent evt3b(0, "PedRdm","Interpolated random pedestal;;Signal [~phe]");
406 //evt2b.SetErrorSpread(kFALSE);
407
408 MFillH fill3b(&evt3b, "MSignalCam", "FillPedRdm");
409 fill3b.SetDrawOption("gaus");
410
411 // ------------------ Setup eventloop and run analysis ---------------
412
413 tlist3.AddToList(&read3);
414 tlist3.AddToList(&apply3);
415 tlist3.AddToList(&drsapply3);
416 tlist3.AddToList(&cont3);
417 tlist3.AddToList(&extractor3);
418// tlist3.AddToList(&fill3a);
419 tlist3.AddToList(&conv3);
420 tlist3.AddToList(&treat3);
421 tlist3.AddToList(&fill3b);
422
423 if (!loop3.Eventloop(max3))
424 return 14;
425
426 if (!loop3.GetDisplay())
427 return 15;
428
429 // ======================================================
430
431 gLog << endl;
432 gLog.Separator("Extracting pedestal");
433
434 MTaskList tlist4;
435
436 MParList plist4;
437 plist4.AddToList(&tlist4);
438 plist4.AddToList(&drscalib300);
439 plist4.AddToList(&badpixels);
440 plist4.AddToList(&timecam);
441
442 MEvtLoop loop4("DetermineExtractedPed");
443 loop4.SetDisplay(d);
444 loop4.SetParList(&plist4);
445
446 // ------------------ Setup the tasks ---------------
447
448 MRawFitsRead read4;
449 read4.LoadMap(map);
450 read4.AddFile(pedfile);
451
452 MContinue cont4("(MRawEvtHeader.GetTriggerID&0xff00)!=0x400", "SelectPed");
453
454 MGeomApply apply4;
455
456 MDrsCalibApply drsapply4;
457
458 MExtractTimeAndChargeSpline extractor4;
459 extractor4.SetRange(first_slice, last_slice);
460 extractor4.SetRiseTimeHiGain(rise_time_dat);
461 extractor4.SetFallTimeHiGain(fall_time_dat);
462 extractor4.SetHeightTm(heighttm);
463 extractor4.SetChargeType(type);
464 extractor4.SetSaturationLimit(600000);
465 extractor4.SetNoiseCalculation(kFALSE);
466
467 // MHCamEvent evt3a(0, "PedExt", "Extracted Pedestal Signal;;S");
468
469 // MFillH fill3a(&evt3a, "MExtractedSignalCam", "FillPedExt");
470
471 // Use this for data, but not for calibration events
472// evt3a.SetErrorSpread(kFALSE);
473/*
474 MCalibrateData conv4;
475 conv4.SetCalibrationMode(MCalibrateData::kNone);
476 conv4.SetPedestalFlag(MCalibrateData::kNo);
477 conv4.SetCalibConvMinLimit(0);
478 conv4.SetCalibConvMaxLimit(10000);
479 conv4.SetScaleFactor(scale);
480*/
481 MCalibrateFact conv4;
482 conv4.SetScale(scale);
483 conv4.SetCalibConst(calib);
484
485 MBadPixelsTreat treat4;
486 treat4.SetProcessPedestalRun(kFALSE);
487 treat4.SetProcessPedestalEvt(kFALSE);
488
489 MHCamEvent evt4b(0, "PedExt","Interpolated extracted pedestal;;Signal [~phe]");
490 //evt4b.SetErrorSpread(kFALSE);
491
492 MFillH fill4b(&evt4b, "MSignalCam", "FillPedExt");
493 fill4b.SetDrawOption("gaus");
494
495 // ------------------ Setup eventloop and run analysis ---------------
496
497 tlist4.AddToList(&read4);
498 tlist4.AddToList(&apply4);
499 tlist4.AddToList(&drsapply4);
500 tlist4.AddToList(&cont4);
501 tlist4.AddToList(&extractor4);
502// tlist4.AddToList(&fill4a);
503 tlist4.AddToList(&conv4);
504 tlist4.AddToList(&treat4);
505 tlist4.AddToList(&fill4b);
506
507 if (!loop4.Eventloop(max4))
508 return 15;
509
510 if (!loop4.GetDisplay())
511 return 16;
512
513 // ===================================================================
514
515 gLog << endl;
516 gLog.Separator("Extracting and calibration data");
517
518 MTaskList tlist5;
519
520 MParList plist5;
521 plist5.AddToList(&tlist5);
522 plist5.AddToList(&drscalib300);
523 plist5.AddToList(&badpixels);
524 plist5.AddToList(&timecam);
525
526 MEvtLoop loop5("CalibratingData");
527 loop5.SetDisplay(d);
528 loop5.SetParList(&plist5);
529
530 // ------------------ Setup the tasks ---------------
531
532 MRawFitsRead read5;
533 read5.LoadMap(map);
534 read5.AddFile(datfile);
535
536 MFillH fill5a(&hrate);
537
538 MGeomApply apply5;
539
540 MDrsCalibApply drsapply5;
541
542 MFDataPhrase filterdat("(MRawEvtHeader.GetTriggerID&0xff00)==0", "SelectDat");
543 MFDataPhrase filtercal("(MRawEvtHeader.GetTriggerID&0xff00)==0x100", "SelectCal");
544 MFDataPhrase filterped("(MRawEvtHeader.GetTriggerID&0xff00)==0x400", "SelectPed");
545 MFDataPhrase filterncl("(MRawEvtHeader.GetTriggerID&0xff00)!=0x100", "SelectNonCal");
546
547 //MContinue cont4("MRawEvtHeader.GetTriggerID!=4", "SelectData");
548
549 // ---
550
551 MExtractTimeAndChargeSpline extractor5dat;
552 extractor5dat.SetRange(first_slice, last_slice);
553 extractor5dat.SetRiseTimeHiGain(rise_time_dat);
554 extractor5dat.SetFallTimeHiGain(fall_time_dat);
555 extractor5dat.SetHeightTm(heighttm);
556 extractor5dat.SetChargeType(type);
557 extractor5dat.SetSaturationLimit(600000);
558 extractor5dat.SetNoiseCalculation(kFALSE);
559
560 MExtractTimeAndChargeSpline extractor5cal;
561 extractor5cal.SetRange(first_slice, last_slice);
562 extractor5cal.SetRiseTimeHiGain(rise_time_cal);
563 extractor5cal.SetFallTimeHiGain(fall_time_cal);
564 extractor5cal.SetHeightTm(heighttm);
565 extractor5cal.SetChargeType(type);
566 extractor5cal.SetSaturationLimit(600000);
567 extractor5cal.SetNoiseCalculation(kFALSE);
568
569 MExtractTimeAndChargeSpline extractor5tm("ExtractTM");
570 extractor5tm.SetRange(last_slice, 294);
571 extractor5tm.SetRiseTimeHiGain(1);
572 extractor5tm.SetFallTimeHiGain(1);
573 extractor5tm.SetHeightTm(0.5);
574 extractor5tm.SetChargeType(MExtralgoSpline::kAmplitudeRel);
575 extractor5tm.SetSaturationLimit(600000);
576 extractor5tm.SetNoiseCalculation(kFALSE);
577 extractor5tm.SetNameSignalCam("TimeMarkerAmplitude");
578 extractor5tm.SetNameTimeCam("TimeMarkerTime");
579
580 extractor5dat.SetFilter(&filterncl);
581 extractor5cal.SetFilter(&filtercal);
582 //extractor4tm.SetFilter(&filtercal);
583
584 // ---
585/*
586 MCalibrateData conv5;
587 conv5.SetCalibrationMode(MCalibrateData::kNone);
588 conv5.SetPedestalFlag(MCalibrateData::kNo);
589 conv5.SetCalibConvMinLimit(0);
590 conv5.SetCalibConvMaxLimit(10000);
591 conv5.SetScaleFactor(scale);
592*/
593 MCalibrateFact conv5;
594 conv5.SetScale(scale);
595 conv5.SetCalibConst(calib);
596
597 MCalibrateDrsTimes calctm5;
598 calctm5.SetNameUncalibrated("UncalibratedTimes");
599
600 MCalibrateDrsTimes calctm5tm("CalibrateTimeMarker");
601 calctm5tm.SetNameArrivalTime("TimeMarkerTime");
602 calctm5tm.SetNameUncalibrated("UncalTimeMarker");
603 calctm5tm.SetNameCalibrated("TimeMarker");
604 calctm5tm.SetTimeMarker();
605 //calctm4tm.SetFilter(&filtercal);
606
607 MBadPixelsTreat treat5;
608 treat5.SetProcessPedestalRun(kFALSE);
609 treat5.SetProcessPedestalEvt(kFALSE);
610
611 MHCamEvent evt5b(0, "ExtSig", "Extracted signal;;S [mV·sl]");
612 MHCamEvent evt5c(0, "CalSig", "Calibrated and interpolated signal;;S [~phe]");
613 MHCamEvent evt5d(4, "ExtSigTm", "Extracted time;;T [sl]");
614 MHCamEvent evt5e(6, "CalSigTm", "Calibrated and interpolated time;;T [ns]");
615
616 MFillH fill5b(&evt5b, "MExtractedSignalCam", "FillExtSig");
617 MFillH fill5c(&evt5c, "MSignalCam", "FillCalSig");
618 MFillH fill5d(&evt5d, "MArrivalTimeCam", "FillExtTm");
619 MFillH fill5e(&evt5e, "MSignalCam", "FillCalTm");
620
621 fill5c.SetDrawOption("gaus");
622 fill5d.SetDrawOption("gaus");
623 fill5e.SetDrawOption("gaus");
624
625 /*
626 fill4b.SetFilter(&filterdat);
627 fill4c.SetFilter(&filterdat);
628 fill4d.SetFilter(&filterdat);
629 fill4e.SetFilter(&filterdat);
630 */
631
632 //MFSoftwareTrigger swtrig;
633 //MContinue contsw(&swtrig, "FilterSwTrigger", "Software trigger");
634 //contsw.SetInverted();
635
636 // The second rule is for the case reading raw-files!
637
638 MWriteRootFile write5(outfile, "RECREATE", "Calibrated Data", 2);
639 write5.AddContainer("MRawRunHeader", "RunHeaders");
640 write5.AddContainer("MGeomCam", "RunHeaders");
641 write5.AddContainer("MMcCorsikaRunHeader", "RunHeaders", kFALSE);
642 write5.AddContainer("MCorsikaRunHeader", "RunHeaders", kFALSE);
643 write5.AddContainer("MMcRunHeader", "RunHeaders", kFALSE);
644
645 // Common events
646 write5.AddContainer("MCorsikaEvtHeader", "Events", kFALSE);
647 write5.AddContainer("MMcEvt", "Events", kFALSE);
648 write5.AddContainer("IncidentAngle", "Events", kFALSE);
649 write5.AddContainer("MPointingPos", "Events", kFALSE);
650 write5.AddContainer("MSignalCam", "Events");
651 write5.AddContainer("MTime", "Events");
652 write5.AddContainer("MRawEvtHeader", "Events");
653 //write.AddContainer("MTriggerPattern", "Events");
654
655 // ------------------ Setup histograms and fill tasks ----------------
656
657 MContinue test;
658 test.SetFilter(&filterncl);
659
660 MTaskList tlist5tm;
661 tlist5tm.AddToList(&extractor5tm);
662 tlist5tm.AddToList(&calctm5tm);
663 tlist5tm.SetFilter(&filtercal);
664
665 MTaskList tlist5dat;
666 tlist5dat.AddToList(&fill5b);
667 tlist5dat.AddToList(&fill5c);
668 tlist5dat.AddToList(&fill5d);
669 tlist5dat.AddToList(&fill5e);
670 tlist5dat.SetFilter(&filterdat);
671
672 tlist5.AddToList(&read5);
673 tlist5.AddToList(&apply5);
674 tlist5.AddToList(&drsapply5);
675 tlist5.AddToList(&filterncl);
676 //tlist5.AddToList(&test);
677 tlist5.AddToList(&filterdat);
678 tlist5.AddToList(&filtercal);
679 tlist5.AddToList(&filterped);
680 tlist5.AddToList(&fill5a);
681 tlist5.AddToList(&extractor5dat);
682 tlist5.AddToList(&extractor5cal);
683 tlist5.AddToList(&calctm5);
684 tlist5.AddToList(&tlist5tm);
685 tlist5.AddToList(&conv5);
686 tlist5.AddToList(&treat5);
687 tlist5.AddToList(&tlist5dat);
688 tlist5.AddToList(&write5);
689
690 if (!loop5.Eventloop(max4))
691 return 18;
692
693 if (!loop5.GetDisplay())
694 return 19;
695
696 d->SetTitle(displaytitle, kFALSE);
697 d->SaveAs(displayfile);
698
699 return 0;
700}
Note: See TracBrowser for help on using the repository browser.