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