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

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