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

Last change on this file since 3198 was 3188, checked in by gaug, 21 years ago
*** empty log message ***
File size: 10.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 "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
53const TString pedfile="/disc02/Data/rootdata/Miscellaneous/2003_12_19/20031218_03522_P_Park_E.root";
54const TString calfile="/disc02/Data/rootdata/Miscellaneous/2003_12_19/20031218_03527_C_Park_E.root";
55
56void 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
Note: See TracBrowser for help on using the repository browser.