1 | /* ======================================================================== *\
|
---|
2 | !
|
---|
3 | ! *
|
---|
4 | ! * This file is part of MARS, the MAGIC Analysis and Reconstruction
|
---|
5 | ! * Software. It is distributed to you in the hope that it can be a useful
|
---|
6 | ! * and timesaving tool in analysing Data of imaging Cerenkov telescopes.
|
---|
7 | ! * It is distributed WITHOUT ANY WARRANTY.
|
---|
8 | ! *
|
---|
9 | ! * Permission to use, copy, modify and distribute this software and its
|
---|
10 | ! * documentation for any purpose is hereby granted without fee,
|
---|
11 | ! * provided that the above copyright notice appear in all copies and
|
---|
12 | ! * that both that copyright notice and this permission notice appear
|
---|
13 | ! * in supporting documentation. It is provided "as is" without express
|
---|
14 | ! * or implied warranty.
|
---|
15 | ! *
|
---|
16 | !
|
---|
17 | ! Author(s): Josep Flix, 01/2004 <mailto:jflix@ifae.es>
|
---|
18 | ! Javier Rico, 01/2004 <mailto:jrico@ifae.es>
|
---|
19 | !
|
---|
20 | ! (based on bootcampstandardanalysis.C by Javier López)
|
---|
21 | !
|
---|
22 | !
|
---|
23 | ! Copyright: MAGIC Software Development, 2000-2004
|
---|
24 | !
|
---|
25 | !
|
---|
26 | \* ======================================================================== */
|
---|
27 |
|
---|
28 | #include "MParList.h"
|
---|
29 | #include "MTaskList.h"
|
---|
30 | #include "MPedestalCam.h"
|
---|
31 | #include "MReadMarsFile.h"
|
---|
32 | #include "MGeomApply.h"
|
---|
33 | #include "MPedCalcPedRun.h"
|
---|
34 | #include "MGeomCamMagic.h"
|
---|
35 | #include "MEvtLoop.h"
|
---|
36 | #include "MCalibrationCam.h"
|
---|
37 | #include "MExtractedSignalCam.h"
|
---|
38 | #include "MExtractSignal.h"
|
---|
39 | #include "MCalibrationCalc.h"
|
---|
40 | #include "MCerPhotEvt.h"
|
---|
41 | #include "MCalibrate.h"
|
---|
42 | #include "MPedPhotCalc.h"
|
---|
43 | #include "MPedPhotCam.h"
|
---|
44 | #include "MHCamera.h"
|
---|
45 |
|
---|
46 | #include "TTimer.h"
|
---|
47 | #include "TString.h"
|
---|
48 | #include "TCanvas.h"
|
---|
49 |
|
---|
50 | #include <iostream>
|
---|
51 |
|
---|
52 |
|
---|
53 | const TString pedfile="/disc02/Data/rootdata/Miscellaneous/2003_12_19/20031218_03522_P_Park_E.root";
|
---|
54 | const TString calfile="/disc02/Data/rootdata/Miscellaneous/2003_12_19/20031218_03527_C_Park_E.root";
|
---|
55 |
|
---|
56 | void pedphotcalc(TString pedname=pedfile,TString calname=calfile)
|
---|
57 | {
|
---|
58 |
|
---|
59 | /*
|
---|
60 | This macro is an example of the use of MPedPhotCalc. It computes
|
---|
61 | the pedestal mean and rms from pedestal files undergoing the
|
---|
62 | signal extraction + calibration chain, in units of photons. It
|
---|
63 | produces plots of the relevant computed quantities.
|
---|
64 | */
|
---|
65 |
|
---|
66 | /************************************/
|
---|
67 | /* FIRST LOOP: PEDESTAL COMPUTATION */
|
---|
68 | /************************************/
|
---|
69 |
|
---|
70 | MParList plist;
|
---|
71 | MTaskList tlist;
|
---|
72 | plist.AddToList(&tlist);
|
---|
73 |
|
---|
74 | // containers
|
---|
75 | MPedestalCam pedcam;
|
---|
76 | plist.AddToList(&pedcam);
|
---|
77 |
|
---|
78 | //tasks
|
---|
79 | MReadMarsFile read("Events", pedname);
|
---|
80 | MGeomApply geomapl;
|
---|
81 | MPedCalcPedRun pedcalc;
|
---|
82 | MGeomCamMagic geomcam;
|
---|
83 |
|
---|
84 | read.DisableAutoScheme();
|
---|
85 |
|
---|
86 | tlist.AddToList(&read);
|
---|
87 | tlist.AddToList(&geomapl);
|
---|
88 | tlist.AddToList(&pedcalc);
|
---|
89 |
|
---|
90 | // Create and execute the event looper
|
---|
91 | MEvtLoop evtloop;
|
---|
92 | evtloop.SetParList(&plist);
|
---|
93 |
|
---|
94 | cout << "*************************" << endl;
|
---|
95 | cout << "** COMPUTING PEDESTALS **" << endl;
|
---|
96 | cout << "*************************" << endl;
|
---|
97 | if (!evtloop.Eventloop())
|
---|
98 | return;
|
---|
99 |
|
---|
100 | tlist.PrintStatistics();
|
---|
101 |
|
---|
102 | /****************************/
|
---|
103 | /* SECOND LOOP: CALIBRATION */
|
---|
104 | /****************************/
|
---|
105 |
|
---|
106 | MParList plist2;
|
---|
107 | MTaskList tlist2;
|
---|
108 | plist2.AddToList(&tlist2);
|
---|
109 |
|
---|
110 | // containers
|
---|
111 | MCalibrationCam calcam;
|
---|
112 | MExtractedSignalCam sigcam;
|
---|
113 |
|
---|
114 | plist2.AddToList(&geomcam);
|
---|
115 | plist2.AddToList(&pedcam);
|
---|
116 | plist2.AddToList(&calcam);
|
---|
117 | plist2.AddToList(&sigcam);
|
---|
118 |
|
---|
119 | //tasks
|
---|
120 | MReadMarsFile read2("Events", calname);
|
---|
121 | read2.DisableAutoScheme();
|
---|
122 |
|
---|
123 | MExtractSignal sigcalc;
|
---|
124 | MCalibrationCalc calcalc;
|
---|
125 |
|
---|
126 | //
|
---|
127 | // Apply a filter against cosmics
|
---|
128 | // (was directly in MCalibrationCalc in earlier versions)
|
---|
129 | //
|
---|
130 | MFCosmics cosmics;
|
---|
131 | MContinue cont(&cosmics);
|
---|
132 |
|
---|
133 | //
|
---|
134 | // As long as we don't have digital modules,
|
---|
135 | // we have to set the color of the pulser LED by hand
|
---|
136 | //
|
---|
137 | calcalc.SetPulserColor(MCalibrationCalc::kECT1);
|
---|
138 |
|
---|
139 | tlist2.AddToList(&read2);
|
---|
140 | tlist2.AddToList(&geomapl);
|
---|
141 | tlist2.AddToList(&sigcalc);
|
---|
142 | //
|
---|
143 | // In case, you want to skip the cosmics rejection,
|
---|
144 | // uncomment the next line
|
---|
145 | //
|
---|
146 | tlist2.AddToList(&cont);
|
---|
147 | //
|
---|
148 | tlist2.AddToList(&calcalc);
|
---|
149 |
|
---|
150 | // Execute second loop (calibration)
|
---|
151 | MEvtLoop evtloop2;
|
---|
152 | evtloop2.SetParList(&plist2);
|
---|
153 |
|
---|
154 | cout << "***********************************" << endl;
|
---|
155 | cout << "** COMPUTING CALIBRATION FACTORS **" << endl;
|
---|
156 | cout << "***********************************" << endl;
|
---|
157 | if (!evtloop2.Eventloop())
|
---|
158 | return;
|
---|
159 |
|
---|
160 | tlist2.PrintStatistics();
|
---|
161 |
|
---|
162 |
|
---|
163 | /************************************************************************/
|
---|
164 | /* THIRD LOOP: PEDESTAL COMPUTATION USING EXTRACTED SIGNAL (IN PHOTONS) */
|
---|
165 | /************************************************************************/
|
---|
166 |
|
---|
167 | // Create an empty Parameter List and an empty Task List
|
---|
168 | MParList plist3;
|
---|
169 | MTaskList tlist3;
|
---|
170 | plist3.AddToList(&tlist3);
|
---|
171 |
|
---|
172 | // containers
|
---|
173 | MCerPhotEvt photcam;
|
---|
174 | MPedPhotCam pedphotcam;
|
---|
175 | MExtractedSignalCam extedsig;
|
---|
176 |
|
---|
177 | plist3.AddToList(&geomcam);
|
---|
178 | plist3.AddToList(&pedcam);
|
---|
179 | plist3.AddToList(&calcam);
|
---|
180 | plist3.AddToList(&photcam);
|
---|
181 | plist3.AddToList(&extedsig);
|
---|
182 | plist3.AddToList(&pedphotcam);
|
---|
183 |
|
---|
184 | //tasks
|
---|
185 | MReadMarsFile read3("Events", pedname);
|
---|
186 | MCalibrate photcalc;
|
---|
187 | MPedPhotCalc pedphotcalc;
|
---|
188 |
|
---|
189 | read3.DisableAutoScheme();
|
---|
190 |
|
---|
191 | tlist3.AddToList(&read3);
|
---|
192 | tlist3.AddToList(&geomapl);
|
---|
193 | tlist3.AddToList(&sigcalc);
|
---|
194 | tlist3.AddToList(&photcalc);
|
---|
195 | tlist3.AddToList(&pedphotcalc);
|
---|
196 |
|
---|
197 | // Create and execute eventloop
|
---|
198 | MEvtLoop evtloop3;
|
---|
199 | evtloop3.SetParList(&plist3);
|
---|
200 |
|
---|
201 | cout << "*************************************************************" << endl;
|
---|
202 | cout << "** COMPUTING PEDESTALS USING EXTRACTED SIGNAL (IN PHOTONS) **" << endl;
|
---|
203 | cout << "*************************************************************" << endl;
|
---|
204 |
|
---|
205 | if (!evtloop3.Eventloop())
|
---|
206 | return;
|
---|
207 | tlist3.PrintStatistics();
|
---|
208 |
|
---|
209 |
|
---|
210 | /**********************/
|
---|
211 | /* PRODUCE NICE PLOTS */
|
---|
212 | /**********************/
|
---|
213 |
|
---|
214 | const UShort_t npix = 577;
|
---|
215 |
|
---|
216 | // declare histograms
|
---|
217 | // pedestals
|
---|
218 | TH1F* pedhist = new TH1F("pedhist","Pedestal",npix,0,npix);
|
---|
219 | TH1F* pedrmshist = new TH1F("pedrmshist","Pedestal RMS",npix,0,npix);
|
---|
220 | TH1F* peddist = new TH1F("peddist","Pedestal",100,0,20);
|
---|
221 | TH1F* pedrmsdist = new TH1F("pedrmsdist","Pedestal RMS",100,0,15);
|
---|
222 |
|
---|
223 | // calibration factors
|
---|
224 | TH1F* calhist = new TH1F("calhist","Calibration factors",npix,0,npix);
|
---|
225 | TH1F* caldist = new TH1F("caldist","Calibration factors",100,0,1);
|
---|
226 |
|
---|
227 | // pedestal signals
|
---|
228 | TH1F* pedphothist = new TH1F("pedphothist","Pedestal Signal",npix,0,npix);
|
---|
229 | TH1F* pedphotrmshist = new TH1F("pedphotrmshist","Pedestal Signal RMS",npix,0,npix);
|
---|
230 | TH1F* pedphotdist = new TH1F("pedphotdist","Pedestal Signal",100,-0.4,0.4);
|
---|
231 | TH1F* pedphotrmsdist = new TH1F("pedphotrmsdist","Pedestal Signal RMS",100,0,25);
|
---|
232 |
|
---|
233 | // fill histograms
|
---|
234 | for(int i=0;i<npix;i++)
|
---|
235 | {
|
---|
236 | const Float_t ped = pedcam[i].GetPedestal();
|
---|
237 | const Float_t pedrms = pedcam[i].GetPedestalRms();
|
---|
238 | const Float_t cal = calcam[i].GetMeanConversionBlindPixelMethod();
|
---|
239 | const Float_t calerr = calcam[i].GetErrorConversionBlindPixelMethod();
|
---|
240 | const Float_t pedphot = pedphotcam[i].GetMean();
|
---|
241 | const Float_t pedphotrms = pedphotcam[i].GetRms();
|
---|
242 |
|
---|
243 | pedhist->Fill(i,ped);
|
---|
244 | peddist->Fill(ped);
|
---|
245 | pedrmshist->Fill(i,pedrms);
|
---|
246 | pedrmsdist->Fill(pedrms);
|
---|
247 |
|
---|
248 | calhist->Fill(i,cal);
|
---|
249 | caldist->Fill(cal);
|
---|
250 |
|
---|
251 | pedphothist->Fill(i,pedphot);
|
---|
252 | pedphotdist->Fill(pedphot);
|
---|
253 | pedphotrmshist->Fill(i,pedphotrms);
|
---|
254 | pedphotrmsdist->Fill(pedphotrms);
|
---|
255 | }
|
---|
256 |
|
---|
257 | // Draw
|
---|
258 | gROOT->Reset();
|
---|
259 | gStyle->SetCanvasColor(0);
|
---|
260 | TCanvas* canvas = new TCanvas("canvas","pedphotcalc.C", 0, 100, 650, 800);
|
---|
261 | canvas->SetBorderMode(0); // Delete the Canvas' border line
|
---|
262 | canvas->cd();
|
---|
263 | canvas->SetBorderMode(0);
|
---|
264 |
|
---|
265 | canvas->Divide(2,4);
|
---|
266 |
|
---|
267 | // draw pedestal histo
|
---|
268 | canvas->cd(1);
|
---|
269 | gPad->cd();
|
---|
270 | gPad->SetBorderMode(0);
|
---|
271 |
|
---|
272 | pedhist->SetStats(kFALSE);
|
---|
273 | pedhist->GetXaxis()->SetTitle("Pixel SW Id");
|
---|
274 | pedhist->GetYaxis()->SetTitle("Pedestal (ADC counts)");
|
---|
275 | pedrmshist->SetStats(kFALSE);
|
---|
276 | pedrmshist->SetLineColor(2);
|
---|
277 | pedhist->Draw();
|
---|
278 | pedrmshist->Draw("same");
|
---|
279 |
|
---|
280 | TLegend* leg1 = new TLegend(.14,.68,.39,.88);
|
---|
281 | leg1->SetHeader("");
|
---|
282 | leg1->AddEntry(pedhist,"Pedestal","L");
|
---|
283 | leg1->AddEntry(pedrmshist,"Pedestal RMS","L");
|
---|
284 | leg1->SetFillColor(0);
|
---|
285 | leg1->SetLineColor(0);
|
---|
286 | leg1->SetBorderSize(0);
|
---|
287 | leg1->Draw();
|
---|
288 |
|
---|
289 | // draw pedestal distribution
|
---|
290 | canvas->cd(2);
|
---|
291 | gPad->cd();
|
---|
292 | gPad->SetBorderMode(0);
|
---|
293 | peddist->GetXaxis()->SetTitle("Pedestal (ADC counts)");
|
---|
294 | pedrmsdist->SetLineColor(2);
|
---|
295 | peddist->Draw();
|
---|
296 | pedrmsdist->Draw("same");
|
---|
297 |
|
---|
298 |
|
---|
299 | TLegend* leg2 = new TLegend(.14,.68,.39,.88);
|
---|
300 | leg2->SetHeader("");
|
---|
301 | leg2->AddEntry(pedhist,"Pedestal","L");
|
---|
302 | leg2->AddEntry(pedrmshist,"Pedestal RMS","L");
|
---|
303 | leg2->SetFillColor(0);
|
---|
304 | leg2->SetLineColor(0);
|
---|
305 | leg2->SetBorderSize(0);
|
---|
306 | leg2->Draw();
|
---|
307 |
|
---|
308 | // draw calibration histo
|
---|
309 | canvas->cd(3);
|
---|
310 | gPad->cd();
|
---|
311 | gPad->SetBorderMode(0);
|
---|
312 | calhist->GetXaxis()->SetTitle("Pixel SW Id");
|
---|
313 | calhist->SetMaximum(1);
|
---|
314 | calhist->SetMinimum(0);
|
---|
315 | calhist->GetYaxis()->SetTitle("Calibration factor (photon/ADC count)");
|
---|
316 | calhist->SetStats(kFALSE);
|
---|
317 | calhist->Draw();
|
---|
318 |
|
---|
319 | // draw calibration distribution
|
---|
320 | canvas->cd(4);
|
---|
321 | gPad->cd();
|
---|
322 | gPad->SetBorderMode(0);
|
---|
323 | caldist->GetXaxis()->SetTitle("Calibration factor (photon/ADC count)");
|
---|
324 | caldist->Draw();
|
---|
325 |
|
---|
326 | // draw pedestal signal histo
|
---|
327 | canvas->cd(5);
|
---|
328 | gPad->cd();
|
---|
329 | gPad->SetBorderMode(0);
|
---|
330 | pedphothist->GetXaxis()->SetTitle("Pixel SW Id");
|
---|
331 | pedphothist->GetYaxis()->SetTitle("Pedestal signal (photons)");
|
---|
332 | pedphothist->SetStats(kFALSE);
|
---|
333 | pedphothist->Draw();
|
---|
334 |
|
---|
335 | // draw pedestal signal distribution
|
---|
336 | canvas->cd(6);
|
---|
337 | gPad->cd();
|
---|
338 | gPad->SetBorderMode(0);
|
---|
339 | pedphotdist->GetXaxis()->SetTitle("Pedestal signal (photons)");
|
---|
340 | pedphotdist->Draw();
|
---|
341 |
|
---|
342 | // draw pedestal signal rms histo
|
---|
343 | canvas->cd(7);
|
---|
344 | gPad->cd();
|
---|
345 | gPad->SetBorderMode(0);
|
---|
346 | pedphotrmshist->GetXaxis()->SetTitle("Pixel SW Id");
|
---|
347 | pedphotrmshist->GetYaxis()->SetTitle("Pedestal signal rms (photons)");
|
---|
348 | pedphotrmshist->SetStats(kFALSE);
|
---|
349 | pedphotrmshist->Draw();
|
---|
350 |
|
---|
351 | // draw pedestal signal rms distribution
|
---|
352 | canvas->cd(8);
|
---|
353 | gPad->cd();
|
---|
354 | gPad->SetBorderMode(0);
|
---|
355 | pedphotrmsdist->GetXaxis()->SetTitle("Pedestal signal rms (photons)");
|
---|
356 | pedphotrmsdist->Draw();
|
---|
357 |
|
---|
358 | return;
|
---|
359 | }
|
---|
360 |
|
---|