source: trunk/Mars/fact/analysis/callisto.C@ 17876

Last change on this file since 17876 was 17876, checked in by tbretz, 11 years ago
Removed some old commented code
File size: 28.7 KB
Line 
1#include "MLogManip.h"
2
3int callisto(const char *seqfile="seq/2012/01/23/20120123_023.seq", const char *outpath = "output")
4{
5 // ======================================================
6
7 // true: Display correctly mapped pixels in the camera displays
8 // but the value-vs-index plot is in software/spiral indices
9 // false: Display pixels in hardware/linear indices,
10 // but the order is the camera display is distorted
11 bool usemap = true;
12
13 // map file to use (get that from La Palma!)
14 const char *map = usemap ? "/home/fact/FACT++/FACTmap111030.txt" : NULL;
15
16 Bool_t maximum = kTRUE;
17
18 const char *lp_template = maximum ?
19 "template-lp-extractor-maximum.root" :
20 "template-lp-extractor-leading-edge.root";
21
22 const char *pulse_template = "template-pulse.root";
23
24 // ------------------------------------------------------
25
26 bool use_delays=false;
27
28 int spike_removal=3;
29
30 // ------------------------------------------------------
31
32 // Extraction range in slices. It will always(!) contain the full range
33 // of integration
34 const int first_slice = 20; // 10ns
35 const int last_slice = 250; // 125ns
36
37 // Note that rise and fall time mean different things whether you use IntegralFixed or IntegraRel:
38 //
39 // IntegralFixed:
40 // * fRiseTime: Number of slices left from arrival time
41 // * fFallTime: Number of slices right from arrival time
42 // IntegralRel:
43 // * fRiseTime: Number of slices left from maximum time
44 // * fFallTime: Number of slices right from maximum time
45 //
46 const int rise_time_cal = maximum ? 40 : 10; // was 13; 5ns
47 const int fall_time_cal = maximum ? 120 : 160; // was 23; 80ns
48
49 const int rise_time_dat = maximum ? 10 : 2; // was 13; was 10; 1ns
50 const int fall_time_dat = maximum ? 40 : 48; // was 23; was 40; 24ns
51
52 // Extraction type: Extract integral and half leading edge
53
54 const int type = maximum ? (MExtralgoSpline::kAmplitudeRel) : (MExtralgoSpline::kIntegralFixed);
55 //const int type = MExtralgoSpline::kIntegralFixed;
56
57
58 const double heighttm = 0.5; // IntegralAbs { 1.5pe * 9.6mV/pe } / IntegralRel { 0.5 }
59
60 Long_t max = 0; // All
61 Long_t max0 = max; // Time marker
62 Long_t max1 = max; // Light pulser
63 //Long_t max2 = 3000; // Calibration ratio
64 Long_t max3 = max; // Pedestal Rndm
65 Long_t max4 = max; // Pedestal Ext
66 Long_t max5 = max; // Data
67
68 // ======================================================
69
70 if (map && gSystem->AccessPathName(map, kFileExists))
71 {
72 gLog << err << "ERROR - Cannot access mapping file '" << map << "'" << endl;
73 return 1;
74 }
75
76 // The sequence file which defines the files for the analysis
77 MSequence seq(seqfile);
78 if (!seq.IsValid())
79 {
80 gLog << err << "ERROR - Sequence '" << seqfile << "' invalid!" << endl;
81 return 2;
82 }
83
84 // --------------------------------------------------------------------------------
85
86 gLog.Separator("Callisto");
87 gLog << all << "Calibrate data of sequence '" << seq.GetFileName() << "'" << endl;
88 gLog << endl;
89
90 // ------------------------------------------------------
91
92 ostringstream drsname;
93 drsname << gSystem->DirName(seqfile) << "/";
94 drsname << seq.GetNight().GetNightAsInt() << "_";
95 drsname << Form("%03d", seq.GetDrsSequence()) << ".drs.seq";
96
97 MSequence drs(drsname.str().c_str());
98 if (!drs.IsValid())
99 {
100 gLog << err << "ERROR - DRS sequence invalid!" << endl;
101 return 3;
102 }
103
104 gLog << all << "DRS sequence file: " << drsname.str() << '\n' << endl;
105
106 TString drsfile = seq.GetFileName(0, MSequence::kRawDrs);
107 if (drsfile.IsNull())
108 {
109 gLog << err << "No DRS file available in sequence." << endl;
110 return 4;
111 }
112
113 TString timfile = drs.GetFileName(0, MSequence::kFitsDat);
114 TString drs1024 = drs.GetFileName(0, MSequence::kFitsDrs);
115 TString pedfile = seq.GetFileName(0, MSequence::kFitsPed);
116 TString calfile = seq.GetFileName(0, MSequence::kFitsCal);
117
118 gLog << all;
119 gLog << "DRS calib 300: " << drsfile << '\n';
120 gLog << "DRS calib 1024: " << drs1024 << "\n\n";
121
122 MDrsCalibration drscalib300;
123 if (!drscalib300.ReadFits(drsfile.Data()))
124 return 5;
125
126 MDrsCalibration drscalib1024;
127 if (!drscalib1024.ReadFits(drs1024.Data()))
128 return 6;
129
130 gLog << all;
131 gLog << "Time calibration : " << timfile << '\n';
132 gLog << "Pedestal file: " << pedfile << '\n';
133 gLog << "Light Pulser file: " << calfile << '\n' << endl;
134
135 // ------------------------------------------------------
136
137 MDirIter iter;
138 if (seq.GetRuns(iter, MSequence::kFitsDat)<=0)
139 {
140 gLog << err << "ERROR - Sequence valid but without files." << endl;
141 return 7;
142 }
143 iter.Print();
144
145 // ======================================================
146
147 MStatusArray arrt, arrp;
148
149 TFile ft(lp_template);
150 if (arrt.Read()<=0)
151 {
152 gLog << err << "ERROR - Reading LP template from " << lp_template << endl;
153 return 100;
154 }
155
156 MHCamera *lpref = (MHCamera*)arrt.FindObjectInCanvas("ExtCalSig;avg", "MHCamera", "Cam");
157 if (!lpref)
158 {
159 gLog << err << "ERROR - LP Template not found in " << lp_template << endl;
160 return 101;
161 }
162 lpref->SetDirectory(0);
163
164 MHCamera *gain = (MHCamera*)arrt.FindObjectInCanvas("gain", "MHCamera", "Gain");
165 if (!gain)
166 {
167 gLog << err << "ERROR - Gain not found in " << lp_template << endl;
168 return 101;
169 }
170 gain->SetDirectory(0);
171
172 TFile fp(pulse_template);
173 if (arrp.Read()<=0)
174 {
175 gLog << err << "ERROR - Reading Pulse template from " << pulse_template << endl;
176 return 102;
177 }
178
179 TH1F *hpulse = (TH1F*)arrp.FindObjectInCanvas("hPixelEdgeMean0_0", "TH1F", "cgpPixelPulses0");
180 if (!hpulse)
181 {
182 gLog << err << "ERROR - Pulse Template not found in " << pulse_template << endl;
183 return 103;
184 }
185 hpulse->SetDirectory(0);
186
187 // ======================================================
188
189 MStatusDisplay *d = new MStatusDisplay;
190
191 MBadPixelsCam badpixels;
192 badpixels.InitSize(1440);
193 badpixels[ 424].SetUnsuitable(MBadPixelsPix::kUnsuitable);
194 badpixels[ 583].SetUnsuitable(MBadPixelsPix::kUnsuitable);
195 badpixels[ 830].SetUnsuitable(MBadPixelsPix::kUnsuitable);
196 badpixels[ 923].SetUnsuitable(MBadPixelsPix::kUnsuitable);
197 badpixels[1208].SetUnsuitable(MBadPixelsPix::kUnsuitable);
198 badpixels[1399].SetUnsuitable(MBadPixelsPix::kUnsuitable);
199
200 // Twin pixel
201 // 113
202 // 115
203 // 354
204 // 423
205 // 1195
206 // 1393
207
208 MDrsCalibrationTime timecam;
209
210 // Plot the trigger pattern rates vs. run-number
211 MH3 hrate("MRawRunHeader.GetFileID", "MRawEvtHeader.GetTriggerID&0xff00");
212 hrate.SetWeight("1./TMath::Max(MRawRunHeader.GetRunLength,1)");
213 hrate.SetName("Rate");
214 hrate.SetTitle("Event rate [Hz];File Id;Trigger Type;");
215 hrate.InitLabels(MH3::kLabelsXY);
216 hrate.DefineLabelY( 0, "Data"); // What if TriggerID==0 already???
217 hrate.DefineLabelY(0x100, "Cal");
218 hrate.DefineLabelY(0x400, "Ped");
219 // hrate.DefaultLabelY("ERROR");
220
221 Bool_t isinteg =
222 (type&MExtralgoSpline::kIntegral) ||
223 (type&MExtralgoSpline::kFixedWidth) ||
224 (type&MExtralgoSpline::kDynWidth)
225 ? kTRUE : kFALSE;
226
227 gStyle->SetOptFit(kTRUE);
228
229 // ======================================================
230
231 gLog << endl;
232 gLog.Separator("Processing DRS timing calibration run");
233
234 MTaskList tlist0;
235
236 MParList plist0;
237 plist0.AddToList(&tlist0);
238 plist0.AddToList(&drscalib1024);
239 plist0.AddToList(&badpixels);
240 plist0.AddToList(&timecam);
241
242 MEvtLoop loop0("DetermineTimeCal");
243 loop0.SetDisplay(d);
244 loop0.SetParList(&plist0);
245
246 // ------------------ Setup the tasks ---------------
247
248 MRawFitsRead read0(timfile);
249
250 MContinue cont0("MRawEvtHeader.GetTriggerID!=33792", "SelectTim");
251
252 MGeomApply apply0;
253
254 MDrsCalibApply drsapply0;
255 drsapply0.SetRemoveSpikes(spike_removal);
256
257 MFillH fill0("MHDrsCalibrationTime");
258 fill0.SetNameTab("DeltaT");
259
260 tlist0.AddToList(&read0);
261 tlist0.AddToList(&apply0);
262 tlist0.AddToList(&drsapply0);
263 tlist0.AddToList(&cont0);
264 tlist0.AddToList(&fill0);
265
266 if (!loop0.Eventloop(max0))
267 return 8;
268
269 if (!loop0.GetDisplay())
270 return 9;
271
272 /*
273 MHDrsCalibrationT *t = (MHDrsCalibrationT*)plist4.FindObject("MHDrsCalibrationT");
274 t->SetDisplay(d);
275 t->PlotAll();
276 */
277
278 // ======================================================
279
280 gLog << endl;
281 gLog.Separator("Processing external light pulser run");
282
283 MTaskList tlist1;
284
285 MParList plist1;
286 plist1.AddToList(&tlist1);
287 plist1.AddToList(&drscalib300);
288 plist1.AddToList(&badpixels);
289 plist1.AddToList(&timecam);
290
291 MEvtLoop loop1("DetermineCalConst");
292 loop1.SetDisplay(d);
293 loop1.SetParList(&plist1);
294
295 // ------------------ Setup the tasks ---------------
296
297 MRawFitsRead read1;
298 read1.LoadMap(map);
299 read1.AddFile(calfile);
300
301 MContinue cont1("(MRawEvtHeader.GetTriggerID&0xff00)!=0x100", "SelectCal");
302
303 MGeomApply apply1;
304
305 MDrsCalibApply drsapply1;
306 drsapply1.SetRemoveSpikes(spike_removal);
307
308 // ---
309
310 MExtractTimeAndChargeSpline extractor1b("ExtractPulse");
311 extractor1b.SetRange(first_slice, last_slice);
312 extractor1b.SetRiseTimeHiGain(rise_time_cal);
313 extractor1b.SetFallTimeHiGain(fall_time_cal);
314 extractor1b.SetHeightTm(heighttm);
315 extractor1b.SetChargeType(type);
316 extractor1b.SetSaturationLimit(600000);
317 extractor1b.SetNoiseCalculation(kFALSE);
318
319 MExtractTimeAndChargeSpline extractor1c("ExtractAmplitude");
320 extractor1c.SetRange(first_slice, last_slice);
321 extractor1c.SetChargeType(MExtralgoSpline::kAmplitude);
322 extractor1c.SetSaturationLimit(600000);
323 extractor1c.SetNoiseCalculation(kFALSE);
324 extractor1c.SetNameSignalCam("Amplitude");
325 extractor1c.SetNameTimeCam("AmplitudePos");
326
327 // ---
328
329 MHCamEvent evt1a(5, "CalRatio", "Ratio per slice between integrated signal and amplitude;; r [1/n]");
330 evt1a.SetNameSub("Amplitude", kTRUE);
331 MFillH fill1a(&evt1a, "MExtractedSignalCam", "FillRatio");
332 fill1a.SetDrawOption("gaus");
333
334 MParameterD ratio1a;
335 ratio1a.SetVal(1./(fall_time_cal+rise_time_cal));
336 fill1a.SetWeight(&ratio1a);
337
338 // ---
339
340 MHCamEvent evt1f(0, "ExtCalSig", "Extracted calibration signal;;S [mV·sl]");
341 MHCamEvent evt1g(4, "ExtCalTm", "Extracted arrival times;;T [sl]");
342 MHCamEvent evt1h(6, "CalCalTm", "Calibrated arrival times;;T [sl]");
343
344 MHSectorVsTime hist1rmsb("ExtSigVsTm");
345 MHSectorVsTime hist1tmb("CalTmVsTm");
346 hist1rmsb.SetTitle("Extracted calibration vs event number;;S [mV·sl]");
347 hist1rmsb.SetType(0);
348 hist1tmb.SetTitle("Extracted arrival time vs event number;;T [sl]");
349 //hist1tmb.SetType(4);
350 hist1tmb.SetType(6);
351
352 MFillH fill1f(&evt1f, "MExtractedSignalCam", "FillExtSig");
353 MFillH fill1g(&evt1g, "MArrivalTimeCam", "FillExtTm");
354 MFillH fill1h(&evt1h, "MSignalCam", "FillCalTm");
355 MFillH fill1r(&hist1rmsb, "MExtractedSignalCam", "FillExtSigVsTm");
356 //MFillH fill1j(&hist1tmb, "MArrivalTimeCam", "FillExtTmVsTm");
357 MFillH fill1j(&hist1tmb, "MSignalCam", "FillCalTmVsTm");
358
359 fill1f.SetDrawOption("gaus");
360 fill1h.SetDrawOption("gaus");
361
362 // ---
363
364 MCalibrateDrsTimes calctm1a("CalibrateCalEvents");
365 calctm1a.SetNameUncalibrated("UncalibratedTimes");
366
367 MBadPixelsTreat treat1;
368 treat1.SetProcessPedestalRun(kFALSE);
369 treat1.SetProcessPedestalEvt(kFALSE);
370
371 // ---
372
373 MHCamEvent evt1c(6, "ExtCalTmShift", "Relative extracted arrival time of calibration pulse (w.r.t. event-median);;\\Delta T [ns]");
374 MHCamEvent evt1d(6, "CalCalTmShift", "Relative calibrated arrival time of calibration pulse (w.r.t. event-median);;\\Delta T [ns]");
375
376 evt1c.SetMedianShift();
377 evt1d.SetMedianShift();
378
379 MFillH fill1c(&evt1c, "UncalibratedTimes", "FillExtCalTm");
380 MFillH fill1d(&evt1d, "MSignalCam", "FillCalCalTm");
381 fill1d.SetDrawOption("gaus");
382
383 // ------------------ Setup eventloop and run analysis ---------------
384
385 tlist1.AddToList(&read1);
386 tlist1.AddToList(&apply1);
387 tlist1.AddToList(&drsapply1);
388 tlist1.AddToList(&cont1);
389 tlist1.AddToList(&extractor1b);
390 if (isinteg)
391 {
392 tlist1.AddToList(&extractor1c);
393 tlist1.AddToList(&fill1a);
394 }
395 tlist1.AddToList(&calctm1a);
396 tlist1.AddToList(&treat1);
397 tlist1.AddToList(&fill1f);
398 tlist1.AddToList(&fill1g);
399 tlist1.AddToList(&fill1h);
400 tlist1.AddToList(&fill1r);
401 tlist1.AddToList(&fill1j);
402 tlist1.AddToList(&fill1c);
403 tlist1.AddToList(&fill1d);
404
405 if (!loop1.Eventloop(max1))
406 return 10;
407
408 if (!loop1.GetDisplay())
409 return 11;
410
411 if (use_delays)
412 timecam.SetDelays(*evt1h.GetHist());
413
414 // ========================= Result ==================================
415
416 Double_t avgS = evt1f.GetHist()->GetMean();
417 Double_t medS = evt1f.GetHist()->GetMedian();
418 Double_t rmsS = evt1f.GetHist()->GetRMS();
419 Double_t maxS = evt1f.GetHist()->GetMaximum();
420
421 MArrayF der1(hpulse->GetNbinsX());
422 MArrayF der2(hpulse->GetNbinsX());
423
424 MExtralgoSpline spline(hpulse->GetArray()+1, hpulse->GetNbinsX(),
425 der1.GetArray(), der2.GetArray());
426 spline.SetRiseFallTime(rise_time_dat, fall_time_dat);
427 spline.SetExtractionType(type);
428 spline.SetHeightTm(heighttm);
429
430 spline.Extract(hpulse->GetMaximumBin()-1);
431
432 // The pulser signal is most probably around 400mV/9.5mV
433 // IntegraFixed 2/48 corresponds to roughly 215mV*50slices
434 Double_t scale = 1./spline.GetSignal();
435
436 MArrayD calib(1440);
437 for (int i=0; i<1440; i++)
438 {
439 Double_t g = gain->GetBinContent(i+1)>0.5 ? gain->GetBinContent(i+1) : 1;
440 if (evt1f.GetHist()->GetBinContent(i+1)>0 && !badpixels[i].IsUnsuitable())
441 calib[i] = lpref->GetBinContent(i+1) / evt1f.GetHist()->GetBinContent(i+1) / g;
442 }
443
444 gROOT->SetSelectedPad(0);
445 d->AddTab("PulseTemp");
446 gPad->SetGrid();
447 hpulse->SetNameTitle("Pulse", "Single p.e. pulse template");
448 hpulse->SetDirectory(0);
449 hpulse->SetLineColor(kBlack);
450 hpulse->DrawCopy();
451
452 TAxis *ax = hpulse->GetXaxis();
453
454 Double_t w = hpulse->GetBinWidth(1);
455 Double_t T = w*(spline.GetTime()+0.5) +ax->GetXmin();
456 Double_t H = w*(hpulse->GetMaximumBin()+0.5)+ax->GetXmin();
457
458 TLine line;
459 line.SetLineColor(kRed);
460 line.DrawLine(T-rise_time_dat*w, spline.GetHeight(),
461 T+fall_time_dat*w, spline.GetHeight());
462 line.DrawLine(T, spline.GetHeight()/4, T, 3*spline.GetHeight()/4);
463 line.DrawLine(T-rise_time_dat*w, 0,
464 T-rise_time_dat*w, spline.GetHeight());
465 line.DrawLine(T+fall_time_dat*w, 0,
466 T+fall_time_dat*w, spline.GetHeight());
467
468 TGraph gg;
469 for (int ix=1; ix<=hpulse->GetNbinsX(); ix++)
470 for (int i=0; i<10; i++)
471 {
472 Double_t x = hpulse->GetBinLowEdge(ix)+i*hpulse->GetBinWidth(ix)/10.;
473 gg.SetPoint(gg.GetN(), x+w/2, spline.EvalAt(ix-1+i/10.));
474 }
475
476 gg.SetLineColor(kBlue);
477 gg.SetMarkerColor(kBlue);
478 gg.SetMarkerStyle(kFullDotMedium);
479 gg.DrawClone("L");
480
481 gROOT->SetSelectedPad(0);
482 d->AddTab("CalConst");
483 MGeomCamFACT fact;
484 MHCamera hcalco(fact);
485 hcalco.SetName("CalConst");
486 hcalco.SetTitle(Form("Relative calibration constant [%.0f/pe]", 1./scale));
487 hcalco.SetCamContent(calib);
488 hcalco.SetAllUsed();
489 //hcalco.Scale(scale);
490 hcalco.DrawCopy();
491
492 // ======================================================
493
494 gLog << endl;
495 gLog.Separator("Extracting random pedestal");
496
497 MTaskList tlist3;
498
499 MParList plist3;
500 plist3.AddToList(&tlist3);
501 plist3.AddToList(&drscalib300);
502 plist3.AddToList(&badpixels);
503 plist3.AddToList(&timecam);
504
505 MEvtLoop loop3("DetermineRndmPed");
506 loop3.SetDisplay(d);
507 loop3.SetParList(&plist3);
508
509 // ------------------ Setup the tasks ---------------
510
511 MRawFitsRead read3;
512 read3.LoadMap(map);
513 read3.AddFile(pedfile);
514
515 MFillH fill3a(&hrate);
516
517 MContinue cont3("(MRawEvtHeader.GetTriggerID&0xff00)!=0x400", "SelectPed");
518
519 MGeomApply apply3;
520
521 MDrsCalibApply drsapply3;
522 drsapply3.SetRemoveSpikes(spike_removal);
523
524 //---
525
526 MExtractTimeAndChargeSpline extractor3;
527 extractor3.SetRange(first_slice, last_slice);
528 extractor3.SetRiseTimeHiGain(rise_time_dat);
529 extractor3.SetFallTimeHiGain(fall_time_dat);
530 extractor3.SetHeightTm(heighttm);
531 extractor3.SetChargeType(type);
532 extractor3.SetSaturationLimit(600000);
533 extractor3.SetNoiseCalculation(kTRUE);
534
535 MCalibrateFact conv3;
536 conv3.SetScale(scale);
537 conv3.SetCalibConst(calib);
538
539 MBadPixelsTreat treat3;
540 treat3.SetProcessPedestalRun(kFALSE);
541 treat3.SetProcessPedestalEvt(kFALSE);
542 treat3.SetProcessTimes(kFALSE);
543
544 MHCamEvent evt3b(0, "PedRdm","Interpolated random pedestal;;Signal [~phe]");
545 //evt2b.SetErrorSpread(kFALSE);
546
547 MFillH fill3b(&evt3b, "MSignalCam", "FillPedRdm");
548 fill3b.SetDrawOption("gaus");
549
550 // ------------------ Setup eventloop and run analysis ---------------
551
552 tlist3.AddToList(&read3);
553 tlist3.AddToList(&apply3);
554 tlist3.AddToList(&drsapply3);
555 tlist3.AddToList(&cont3);
556 tlist3.AddToList(&extractor3);
557// tlist3.AddToList(&fill3a);
558 tlist3.AddToList(&conv3);
559 tlist3.AddToList(&treat3);
560 tlist3.AddToList(&fill3b);
561
562 if (!loop3.Eventloop(max3))
563 return 14;
564
565 if (!loop3.GetDisplay())
566 return 15;
567
568 // ======================================================
569
570 gLog << endl;
571 gLog.Separator("Extracting pedestal");
572
573 MTaskList tlist4;
574
575 MParList plist4;
576 plist4.AddToList(&tlist4);
577 plist4.AddToList(&drscalib300);
578 plist4.AddToList(&badpixels);
579 plist4.AddToList(&timecam);
580
581 MEvtLoop loop4("DetermineExtractedPed");
582 loop4.SetDisplay(d);
583 loop4.SetParList(&plist4);
584
585 // ------------------ Setup the tasks ---------------
586
587 MRawFitsRead read4;
588 read4.LoadMap(map);
589 read4.AddFile(pedfile);
590
591 MContinue cont4("(MRawEvtHeader.GetTriggerID&0xff00)!=0x400", "SelectPed");
592
593 MGeomApply apply4;
594
595 MDrsCalibApply drsapply4;
596 drsapply4.SetRemoveSpikes(spike_removal);
597
598 MExtractTimeAndChargeSpline extractor4;
599 extractor4.SetRange(first_slice, last_slice);
600 extractor4.SetRiseTimeHiGain(rise_time_dat);
601 extractor4.SetFallTimeHiGain(fall_time_dat);
602 extractor4.SetHeightTm(heighttm);
603 extractor4.SetChargeType(type);
604 extractor4.SetSaturationLimit(600000);
605 extractor4.SetNoiseCalculation(kFALSE);
606
607 MCalibrateFact conv4;
608 conv4.SetScale(scale);
609 conv4.SetCalibConst(calib);
610
611 MBadPixelsTreat treat4;
612 treat4.SetProcessPedestalRun(kFALSE);
613 treat4.SetProcessPedestalEvt(kFALSE);
614
615 MHCamEvent evt4b(0, "PedExt","Interpolated extracted pedestal;;Signal [~phe]");
616 //evt4b.SetErrorSpread(kFALSE);
617
618 MFillH fill4b(&evt4b, "MSignalCam", "FillPedExt");
619 fill4b.SetDrawOption("gaus");
620
621 // ------------------ Setup eventloop and run analysis ---------------
622
623 tlist4.AddToList(&read4);
624 tlist4.AddToList(&apply4);
625 tlist4.AddToList(&drsapply4);
626 tlist4.AddToList(&cont4);
627 tlist4.AddToList(&extractor4);
628// tlist4.AddToList(&fill4a);
629 tlist4.AddToList(&conv4);
630 tlist4.AddToList(&treat4);
631 tlist4.AddToList(&fill4b);
632
633 if (!loop4.Eventloop(max4))
634 return 15;
635
636 if (!loop4.GetDisplay())
637 return 16;
638
639 // ===================================================================
640
641 gLog << endl;
642 gLog.Separator("Extracting and calibration data");
643
644 MTaskList tlist5;
645
646 MParList plist5;
647 plist5.AddToList(&tlist5);
648 plist5.AddToList(&drscalib300);
649 plist5.AddToList(&badpixels);
650 plist5.AddToList(&timecam);
651
652 MEvtLoop loop5("CalibratingData");
653 loop5.SetDisplay(d);
654 loop5.SetParList(&plist5);
655
656 // ------------------ Setup the tasks ---------------
657
658 MRawFitsRead read5;
659 read5.LoadMap(map);
660 read5.AddFiles(iter);
661
662 MFillH fill5a(&hrate);
663
664 MGeomApply apply5;
665
666 MDrsCalibApply drsapply5;
667 drsapply5.SetRemoveSpikes(spike_removal);
668
669 MFDataPhrase filterdat("(MRawEvtHeader.GetTriggerID&0xff00)==0", "SelectDat");
670 MFDataPhrase filtercal("(MRawEvtHeader.GetTriggerID&0xff00)==0x100", "SelectCal");
671 MFDataPhrase filterped("(MRawEvtHeader.GetTriggerID&0xff00)==0x400", "SelectPed");
672 MFDataPhrase filterncl("(MRawEvtHeader.GetTriggerID&0xff00)!=0x100", "SelectNonCal");
673
674 //MContinue cont4("MRawEvtHeader.GetTriggerID!=4", "SelectData");
675
676 // ---
677
678 MExtractTimeAndChargeSpline extractor5dat;
679 extractor5dat.SetRange(first_slice, last_slice);
680 extractor5dat.SetRiseTimeHiGain(rise_time_dat);
681 extractor5dat.SetFallTimeHiGain(fall_time_dat);
682 extractor5dat.SetHeightTm(heighttm);
683 extractor5dat.SetChargeType(type);
684 extractor5dat.SetSaturationLimit(600000);
685 extractor5dat.SetNoiseCalculation(kFALSE);
686
687 MExtractTimeAndChargeSpline extractor5cal;
688 extractor5cal.SetRange(first_slice, last_slice);
689 extractor5cal.SetRiseTimeHiGain(rise_time_cal);
690 extractor5cal.SetFallTimeHiGain(fall_time_cal);
691 extractor5cal.SetHeightTm(heighttm);
692 extractor5cal.SetChargeType(type);
693 extractor5cal.SetSaturationLimit(600000);
694 extractor5cal.SetNoiseCalculation(kFALSE);
695
696 MExtractTimeAndChargeSpline extractor5tm("ExtractTM");
697 extractor5tm.SetRange(last_slice, 294);
698 extractor5tm.SetRiseTimeHiGain(1);
699 extractor5tm.SetFallTimeHiGain(1);
700 extractor5tm.SetHeightTm(0.5);
701 extractor5tm.SetChargeType(MExtralgoSpline::kAmplitudeRel);
702 extractor5tm.SetSaturationLimit(600000);
703 extractor5tm.SetNoiseCalculation(kFALSE);
704 extractor5tm.SetNameSignalCam("TimeMarkerAmplitude");
705 extractor5tm.SetNameTimeCam("TimeMarkerTime");
706
707 extractor5dat.SetFilter(&filterncl);
708 extractor5cal.SetFilter(&filtercal);
709 //extractor4tm.SetFilter(&filtercal);
710
711 // ---
712 MCalibrateFact conv5;
713 conv5.SetScale(scale);
714 conv5.SetCalibConst(calib);
715
716 MCalibrateDrsTimes calctm5;
717 calctm5.SetNameUncalibrated("UncalibratedTimes");
718
719 MCalibrateDrsTimes calctm5tm("CalibrateTimeMarker");
720 calctm5tm.SetNameArrivalTime("TimeMarkerTime");
721 calctm5tm.SetNameUncalibrated("UncalTimeMarker");
722 calctm5tm.SetNameCalibrated("TimeMarker");
723 calctm5tm.SetTimeMarker();
724 //calctm4tm.SetFilter(&filtercal);
725
726 MHCamEvent evt5m(6, "ExtTm", "Extracted arrival times of calibration pulse;;\\Delta T [ns]");
727 MHCamEvent evt5n(6, "CalTm", "Calibrated arrival times of calibration pulse;;\\Delta T [ns]");
728 MHCamEvent evt5q(6, "ExtTmShift", "Relative extracted arrival times of calibration pulse (w.r.t. event-median);;\\Delta T [ns]");
729 MHCamEvent evt5r(6, "CalTmShift", "Relative calibrated arrival times of calibration pulse (w.r.t. event-median);;\\Delta T [ns]");
730 MHCamEvent evt5s(6, "ExtTM", "Extracted absolute time marker position;;T [sl]");
731 MHCamEvent evt5t(6, "CalTM", "Calibrated absolute time marker position;;T [ns]");
732 MHCamEvent evt5u(6, "ExtTMshift", "Relative extracted time marker position (w.r.t. event-median);;\\Delta T [ns]");
733 MHCamEvent evt5v(6, "CalTMshift", "Relative calibrated time marker position (w.r.t. event-median);;\\Delta T [ns]");
734 MHCamEvent evt5w(6, "ExtDiff", "Difference between extracted arrival time of time marker and calibration pulse;;\\Delta T [ns]");
735 MHCamEvent evt5x(6, "CalDiff", "Difference between calibrated arrival time of time marker and calibration pulse;;\\Delta T [ns]");
736
737 evt5w.SetNameSub("UncalibratedTimes");
738 evt5x.SetNameSub("MSignalCam");
739
740 evt5q.SetMedianShift();
741 evt5r.SetMedianShift();
742 evt5u.SetMedianShift();
743 evt5v.SetMedianShift();
744 //evt4w.SetMedianShift();
745 //evt4x.SetMedianShift();
746
747 MFillH fill5m(&evt5m, "UncalibratedTimes", "FillExtTm");
748 MFillH fill5n(&evt5n, "MSignalCam", "FillCalTm");
749 MFillH fill5q(&evt5q, "UncalibratedTimes", "FillExtTmShift");
750 MFillH fill5r(&evt5r, "MSignalCam" , "FillCalTmShift");
751 MFillH fill5s(&evt5s, "UncalTimeMarker", "FillExtTM");
752 MFillH fill5t(&evt5t, "TimeMarker", "FillCalTM");
753 MFillH fill5u(&evt5u, "UncalTimeMarker", "FillExtTMshift");
754 MFillH fill5v(&evt5v, "TimeMarker", "FillCalTMshift");
755 MFillH fill5w(&evt5w, "UncalTimeMarker", "FillExtDiff");
756 MFillH fill5x(&evt5x, "TimeMarker", "FillCalDiff");
757
758 fill5m.SetDrawOption("gaus");
759 fill5n.SetDrawOption("gaus");
760 fill5q.SetDrawOption("gaus");
761 fill5r.SetDrawOption("gaus");
762 //fill5s.SetDrawOption("gaus");
763 //fill5t.SetDrawOption("gaus");
764 //fill5u.SetDrawOption("gaus");
765 //fill5v.SetDrawOption("gaus");
766 //fill5w.SetDrawOption("gaus");
767 //fill5x.SetDrawOption("gaus");
768
769
770 MBadPixelsTreat treat5;
771 treat5.SetProcessPedestalRun(kFALSE);
772 treat5.SetProcessPedestalEvt(kFALSE);
773
774 MHSectorVsTime hist5cal("CalVsTm");
775 MHSectorVsTime hist5ped("PedVsTm");
776 hist5cal.SetTitle("Median calibrated calibration signal vs event number;;Signal [~phe]");
777 hist5ped.SetTitle("Median calibrated pedestal signal vs event number;;Signal [~phe]");
778 hist5cal.SetType(0);
779 hist5ped.SetType(0);
780 hist5cal.SetMinimum(0);
781 hist5ped.SetMinimum(0);
782 hist5cal.SetUseMedian();
783 hist5ped.SetUseMedian();
784 hist5cal.SetNameTime("MTime");
785 hist5ped.SetNameTime("MTime");
786
787 MFillH fill5cal(&hist5cal, "MSignalCam", "FillCalVsTm");
788 MFillH fill5ped(&hist5ped, "MSignalCam", "FillPedVsTm");
789 fill5cal.SetFilter(&filtercal);
790 fill5ped.SetFilter(&filterped);
791
792 MHCamEvent evt5b(0, "ExtSig", "Extracted signal;;S [mV·sl]");
793 MHCamEvent evt5c(0, "CalSig", "Calibrated and interpolated signal;;S [~phe]");
794 MHCamEvent evt5d(4, "ExtSigTm", "Extracted time;;T [sl]");
795 MHCamEvent evt5e(6, "CalSigTm", "Calibrated and interpolated time;;T [ns]");
796
797 MFillH fill5b(&evt5b, "MExtractedSignalCam", "FillExtSig");
798 MFillH fill5c(&evt5c, "MSignalCam", "FillCalSig");
799 MFillH fill5d(&evt5d, "MArrivalTimeCam", "FillExtTm");
800 MFillH fill5e(&evt5e, "MSignalCam", "FillCalTm");
801
802 fill5c.SetDrawOption("gaus");
803 fill5d.SetDrawOption("gaus");
804 fill5e.SetDrawOption("gaus");
805
806 /*
807 fill4b.SetFilter(&filterdat);
808 fill4c.SetFilter(&filterdat);
809 fill4d.SetFilter(&filterdat);
810 fill4e.SetFilter(&filterdat);
811 */
812
813 //MFSoftwareTrigger swtrig;
814 //MContinue contsw(&swtrig, "FilterSwTrigger", "Software trigger");
815 //contsw.SetInverted();
816
817 const TString fname(Form("s/([0-9]+_[0-9]+)[.]fits([.][fg]z)?$/%s\\/$1_C.root/",
818 MJob::Esc(outpath).Data()));
819
820 // The second rule is for the case reading raw-files!
821 MWriteRootFile write5(2, fname, "RECREATE", "Calibrated Data");
822 write5.AddContainer("MRawRunHeader", "RunHeaders");
823 write5.AddContainer("MGeomCam", "RunHeaders");
824 write5.AddContainer("MSignalCam", "Events");
825 write5.AddContainer("MTime", "Events");
826 write5.AddContainer("MRawEvtHeader", "Events");
827 //write.AddContainer("MTriggerPattern", "Events");
828
829 // ------------------ Setup histograms and fill tasks ----------------
830
831 MContinue test;
832 test.SetFilter(&filterncl);
833
834 MTaskList tlist5tm;
835 tlist5tm.AddToList(&extractor5tm);
836 tlist5tm.AddToList(&calctm5tm);
837 tlist5tm.AddToList(&fill5m);
838 tlist5tm.AddToList(&fill5n);
839 tlist5tm.AddToList(&fill5q);
840 tlist5tm.AddToList(&fill5r);
841 //tlist5tm.AddToList(&fill5s);
842 //tlist5tm.AddToList(&fill5t);
843 tlist5tm.AddToList(&fill5u);
844 tlist5tm.AddToList(&fill5v);
845 tlist5tm.AddToList(&fill5w);
846 tlist5tm.AddToList(&fill5x);
847 tlist5tm.SetFilter(&filtercal);
848
849 MTaskList tlist5dat;
850 tlist5dat.AddToList(&fill5b);
851 tlist5dat.AddToList(&fill5c);
852 tlist5dat.AddToList(&fill5d);
853 tlist5dat.AddToList(&fill5e);
854 tlist5dat.SetFilter(&filterdat);
855
856 tlist5.AddToList(&read5);
857 tlist5.AddToList(&apply5);
858 tlist5.AddToList(&drsapply5);
859 tlist5.AddToList(&filterncl);
860 //tlist5.AddToList(&test);
861 tlist5.AddToList(&filterdat);
862 tlist5.AddToList(&filtercal);
863 tlist5.AddToList(&filterped);
864 tlist5.AddToList(&fill5a);
865 tlist5.AddToList(&extractor5dat);
866 tlist5.AddToList(&extractor5cal);
867 tlist5.AddToList(&calctm5);
868 tlist5.AddToList(&tlist5tm);
869 tlist5.AddToList(&conv5);
870 tlist5.AddToList(&treat5);
871 tlist5.AddToList(&fill5ped);
872 tlist5.AddToList(&fill5cal);
873 tlist5.AddToList(&tlist5dat);
874 tlist5.AddToList(&write5);
875
876 if (!loop5.Eventloop(max4))
877 return 18;
878
879 if (!loop5.GetDisplay())
880 return 19;
881
882 TString title = "-- Calibrated signal #";
883 title += seq.GetSequence();
884 title += " (";
885 title += drsfile;
886 title += ") --";
887 d->SetTitle(title, kFALSE);
888
889 TString path;
890 path += Form("%s/20%6d_%03d-calibration.root", outpath,
891 seq.GetSequence()/1000, seq.GetSequence()%1000);
892
893 d->SaveAs(path);
894
895 return 0;
896}
Note: See TracBrowser for help on using the repository browser.