source: trunk/MagicSoft/Mars/macros/pedphotcalc.C@ 3421

Last change on this file since 3421 was 3298, checked in by gaug, 21 years ago
*** empty log message ***
File size: 11.1 KB
Line 
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 "MCalibrationChargeCam.h"
37#include "MHCalibrationChargeCam.h"
38#include "MCalibrationChargePINDiode.h"
39#include "MHCalibrationChargePINDiode.h"
40#include "MCalibrationChargeBlindPix.h"
41#include "MHCalibrationChargeBlindPix.h"
42#include "MCalibrationChargeCalc.h"
43#include "MExtractedSignalCam.h"
44#include "MExtractSignal.h"
45#include "MExtractedSignalPINDiode.h"
46#include "MExtractPINDiode.h"
47#include "MExtractedSignalBlindPixel.h"
48#include "MExtractBlindPixel.h"
49#include "MCerPhotEvt.h"
50#include "MCalibrate.h"
51#include "MPedPhotCalc.h"
52#include "MPedPhotCam.h"
53#include "MHCamera.h"
54
55#include "TTimer.h"
56#include "TString.h"
57#include "TCanvas.h"
58
59#include <iostream>
60
61
62const TString pedfile="/disc02/Data/rootdata/Miscellaneous/2003_12_19/20031218_03522_P_Park_E.root";
63const TString calfile="/disc02/Data/rootdata/Miscellaneous/2003_12_19/20031218_03527_C_Park_E.root";
64
65void pedphotcalc(TString pedname=pedfile,TString calname=calfile)
66{
67
68 /*
69 This macro is an example of the use of MPedPhotCalc. It computes
70 the pedestal mean and rms from pedestal files undergoing the
71 signal extraction + calibration chain, in units of photons. It
72 produces plots of the relevant computed quantities.
73 */
74
75 /************************************/
76 /* FIRST LOOP: PEDESTAL COMPUTATION */
77 /************************************/
78
79 MParList plist;
80 MTaskList tlist;
81 plist.AddToList(&tlist);
82
83 // containers
84 MPedestalCam pedcam;
85 plist.AddToList(&pedcam);
86
87 //tasks
88 MReadMarsFile read("Events", pedname);
89 MGeomApply geomapl;
90 MPedCalcPedRun pedcalc;
91 MGeomCamMagic geomcam;
92
93 read.DisableAutoScheme();
94
95 tlist.AddToList(&read);
96 tlist.AddToList(&geomapl);
97 tlist.AddToList(&pedcalc);
98
99 // Create and execute the event looper
100 MEvtLoop evtloop;
101 evtloop.SetParList(&plist);
102
103 cout << "*************************" << endl;
104 cout << "** COMPUTING PEDESTALS **" << endl;
105 cout << "*************************" << endl;
106 if (!evtloop.Eventloop())
107 return;
108
109 tlist.PrintStatistics();
110
111 /****************************/
112 /* SECOND LOOP: CALIBRATION */
113 /****************************/
114
115 MParList plist2;
116 MTaskList tlist2;
117 plist2.AddToList(&tlist2);
118
119 // containers
120 MCalibrationChargeCam calcam;
121 MExtractedSignalCam sigcam;
122
123 MCalibrationChargeCam chargecam;
124 MCalibrationChargePINDiode chargepin;
125 MCalibrationChargeBlindPixel chargeblind;
126 chargeblind.SetColor(MCalibrationChargeBlindPix::kECT1);
127 chargepin.SetColor( MCalibrationChargePINDiode::kECT1);
128 //
129 plist2.AddToList(&geomcam);
130 plist2.AddToList(&pedcam);
131 plist2.AddToList(&calcam);
132 plist2.AddToList(&sigcam);
133 plist2.AddToList(&chargecam);
134 plist2.AddToList(&chargepin);
135 plist2.AddToList(&chargeblind);
136
137 //tasks
138 MReadMarsFile read2("Events", calname);
139 read2.DisableAutoScheme();
140
141 MExtractPINDiode pincalc;
142 MExtractBlindPixel blindcalc;
143 MExtractSignal sigcalc;
144 MCalibrationChargeCalc calcalc;
145
146 MFillH fillpin("MHCalibrationChargePINDiode", "MExtractedSignalPINDiode");
147 MFillH fillblind("MHCalibrationChargeBlindPixel", "MExtractedSignalBlindPixel");
148 MFillH fillcam("MHCalibrationChargeCam" , "MExtractedSignalCam");
149 //
150 // Apply a filter against cosmics
151 // (was directly in MCalibrationChargeCalc in earlier versions)
152 //
153 MFCosmics cosmics;
154 MContinue cont(&cosmics);
155
156 tlist2.AddToList(&read2);
157 tlist2.AddToList(&geomapl);
158 tlist2.AddToList(&blindcalc);
159 tlist2.AddToList(&pincalc);
160 tlist2.AddToList(&sigcalc);
161 //
162 // In case, you want to skip the cosmics rejection,
163 // uncomment the next line
164 //
165 tlist2.AddToList(&cont);
166 //
167 tlist2.AddToList(&fillpin);
168 tlist2.AddToList(&fillblind);
169 tlist2.AddToList(&fillcam);
170 tlist2.AddToList(&calcalc);
171
172 // Execute second loop (calibration)
173 MEvtLoop evtloop2;
174 evtloop2.SetParList(&plist2);
175
176 cout << "***********************************" << endl;
177 cout << "** COMPUTING CALIBRATION FACTORS **" << endl;
178 cout << "***********************************" << endl;
179 if (!evtloop2.Eventloop())
180 return;
181
182 tlist2.PrintStatistics();
183
184
185 /************************************************************************/
186 /* THIRD LOOP: PEDESTAL COMPUTATION USING EXTRACTED SIGNAL (IN PHOTONS) */
187 /************************************************************************/
188
189 // Create an empty Parameter List and an empty Task List
190 MParList plist3;
191 MTaskList tlist3;
192 plist3.AddToList(&tlist3);
193
194 // containers
195 MCerPhotEvt photcam;
196 MPedPhotCam pedphotcam;
197 MExtractedSignalCam extedsig;
198
199 plist3.AddToList(&geomcam);
200 plist3.AddToList(&pedcam);
201 plist3.AddToList(&calcam);
202 plist3.AddToList(&photcam);
203 plist3.AddToList(&extedsig);
204 plist3.AddToList(&pedphotcam);
205
206 //tasks
207 MReadMarsFile read3("Events", pedname);
208 MCalibrate photcalc;
209 MPedPhotCalc pedphotcalc;
210
211 read3.DisableAutoScheme();
212
213 tlist3.AddToList(&read3);
214 tlist3.AddToList(&geomapl);
215 tlist3.AddToList(&sigcalc);
216 tlist3.AddToList(&photcalc);
217 tlist3.AddToList(&pedphotcalc);
218
219 // Create and execute eventloop
220 MEvtLoop evtloop3;
221 evtloop3.SetParList(&plist3);
222
223 cout << "*************************************************************" << endl;
224 cout << "** COMPUTING PEDESTALS USING EXTRACTED SIGNAL (IN PHOTONS) **" << endl;
225 cout << "*************************************************************" << endl;
226
227 if (!evtloop3.Eventloop())
228 return;
229 tlist3.PrintStatistics();
230
231
232 /**********************/
233 /* PRODUCE NICE PLOTS */
234 /**********************/
235
236 const UShort_t npix = 577;
237
238 // declare histograms
239 // pedestals
240 TH1F* pedhist = new TH1F("pedhist","Pedestal",npix,0,npix);
241 TH1F* pedrmshist = new TH1F("pedrmshist","Pedestal RMS",npix,0,npix);
242 TH1F* peddist = new TH1F("peddist","Pedestal",100,0,20);
243 TH1F* pedrmsdist = new TH1F("pedrmsdist","Pedestal RMS",100,0,15);
244
245 // calibration factors
246 TH1F* calhist = new TH1F("calhist","Calibration factors",npix,0,npix);
247 TH1F* caldist = new TH1F("caldist","Calibration factors",100,0,1);
248
249 // pedestal signals
250 TH1F* pedphothist = new TH1F("pedphothist","Pedestal Signal",npix,0,npix);
251 TH1F* pedphotrmshist = new TH1F("pedphotrmshist","Pedestal Signal RMS",npix,0,npix);
252 TH1F* pedphotdist = new TH1F("pedphotdist","Pedestal Signal",100,-0.4,0.4);
253 TH1F* pedphotrmsdist = new TH1F("pedphotrmsdist","Pedestal Signal RMS",100,0,25);
254
255 // fill histograms
256 for(int i=0;i<npix;i++)
257 {
258 const Float_t ped = pedcam[i].GetPedestal();
259 const Float_t pedrms = pedcam[i].GetPedestalRms();
260 const Float_t cal = calcam[i].GetMeanConversionBlindPixelMethod();
261 const Float_t calerr = calcam[i].GetErrorConversionBlindPixelMethod();
262 const Float_t pedphot = pedphotcam[i].GetMean();
263 const Float_t pedphotrms = pedphotcam[i].GetRms();
264
265 pedhist->Fill(i,ped);
266 peddist->Fill(ped);
267 pedrmshist->Fill(i,pedrms);
268 pedrmsdist->Fill(pedrms);
269
270 calhist->Fill(i,cal);
271 caldist->Fill(cal);
272
273 pedphothist->Fill(i,pedphot);
274 pedphotdist->Fill(pedphot);
275 pedphotrmshist->Fill(i,pedphotrms);
276 pedphotrmsdist->Fill(pedphotrms);
277 }
278
279 // Draw
280 gROOT->Reset();
281 gStyle->SetCanvasColor(0);
282 TCanvas* canvas = new TCanvas("canvas","pedphotcalc.C", 0, 100, 650, 800);
283 canvas->SetBorderMode(0); // Delete the Canvas' border line
284 canvas->cd();
285 canvas->SetBorderMode(0);
286
287 canvas->Divide(2,4);
288
289 // draw pedestal histo
290 canvas->cd(1);
291 gPad->cd();
292 gPad->SetBorderMode(0);
293
294 pedhist->SetStats(kFALSE);
295 pedhist->GetXaxis()->SetTitle("Pixel SW Id");
296 pedhist->GetYaxis()->SetTitle("Pedestal (ADC counts)");
297 pedrmshist->SetStats(kFALSE);
298 pedrmshist->SetLineColor(2);
299 pedhist->Draw();
300 pedrmshist->Draw("same");
301
302 TLegend* leg1 = new TLegend(.14,.68,.39,.88);
303 leg1->SetHeader("");
304 leg1->AddEntry(pedhist,"Pedestal","L");
305 leg1->AddEntry(pedrmshist,"Pedestal RMS","L");
306 leg1->SetFillColor(0);
307 leg1->SetLineColor(0);
308 leg1->SetBorderSize(0);
309 leg1->Draw();
310
311 // draw pedestal distribution
312 canvas->cd(2);
313 gPad->cd();
314 gPad->SetBorderMode(0);
315 peddist->GetXaxis()->SetTitle("Pedestal (ADC counts)");
316 pedrmsdist->SetLineColor(2);
317 peddist->Draw();
318 pedrmsdist->Draw("same");
319
320
321 TLegend* leg2 = new TLegend(.14,.68,.39,.88);
322 leg2->SetHeader("");
323 leg2->AddEntry(pedhist,"Pedestal","L");
324 leg2->AddEntry(pedrmshist,"Pedestal RMS","L");
325 leg2->SetFillColor(0);
326 leg2->SetLineColor(0);
327 leg2->SetBorderSize(0);
328 leg2->Draw();
329
330 // draw calibration histo
331 canvas->cd(3);
332 gPad->cd();
333 gPad->SetBorderMode(0);
334 calhist->GetXaxis()->SetTitle("Pixel SW Id");
335 calhist->SetMaximum(1);
336 calhist->SetMinimum(0);
337 calhist->GetYaxis()->SetTitle("Calibration factor (photon/ADC count)");
338 calhist->SetStats(kFALSE);
339 calhist->Draw();
340
341 // draw calibration distribution
342 canvas->cd(4);
343 gPad->cd();
344 gPad->SetBorderMode(0);
345 caldist->GetXaxis()->SetTitle("Calibration factor (photon/ADC count)");
346 caldist->Draw();
347
348 // draw pedestal signal histo
349 canvas->cd(5);
350 gPad->cd();
351 gPad->SetBorderMode(0);
352 pedphothist->GetXaxis()->SetTitle("Pixel SW Id");
353 pedphothist->GetYaxis()->SetTitle("Pedestal signal (photons)");
354 pedphothist->SetStats(kFALSE);
355 pedphothist->Draw();
356
357 // draw pedestal signal distribution
358 canvas->cd(6);
359 gPad->cd();
360 gPad->SetBorderMode(0);
361 pedphotdist->GetXaxis()->SetTitle("Pedestal signal (photons)");
362 pedphotdist->Draw();
363
364 // draw pedestal signal rms histo
365 canvas->cd(7);
366 gPad->cd();
367 gPad->SetBorderMode(0);
368 pedphotrmshist->GetXaxis()->SetTitle("Pixel SW Id");
369 pedphotrmshist->GetYaxis()->SetTitle("Pedestal signal rms (photons)");
370 pedphotrmshist->SetStats(kFALSE);
371 pedphotrmshist->Draw();
372
373 // draw pedestal signal rms distribution
374 canvas->cd(8);
375 gPad->cd();
376 gPad->SetBorderMode(0);
377 pedphotrmsdist->GetXaxis()->SetTitle("Pedestal signal rms (photons)");
378 pedphotrmsdist->Draw();
379
380 return;
381}
382
Note: See TracBrowser for help on using the repository browser.