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

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