source: trunk/MagicSoft/Mars/datacenter/macros/fillsignal.C@ 9052

Last change on this file since 9052 was 9005, checked in by tbretz, 16 years ago
*** empty log message ***
File size: 14.5 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!
18! Author(s): Thomas Bretz, 04/2005 <mailto:tbretz@astro.uni-wuerzburg.de>
19! Author(s): Daniela Dorner, 04/2005 <mailto:dorner@astro.uni-wuerzburg.de>
20!
21! Copyright: MAGIC Software Development, 2000-2008
22!
23!
24\* ======================================================================== */
25
26/////////////////////////////////////////////////////////////////////////////
27//
28// fillsignal.C
29// ============
30//
31// This macro is used to read the calibration-/callisto-output files
32// signal00000.root.
33//
34// From this file the mean pedrms, the mean signal and the pulse position
35// for the inner and outer camera is extracted and inserted into the database
36// in the table Calibration, where the results of callisto are stored.
37// The sequence number is extracted from the filename.
38//
39// Usage:
40// .x fillsignal.C("/magic/data/callisto/0004/00047421/signal00047421.root", kTRUE)
41//
42// The second argument is the 'dummy-mode'. If it is kTRUE dummy-mode is
43// switched on and nothing will be written into the database. This is usefull
44// for tests.
45//
46// The macro can also be run without ACLiC but this is a lot slower...
47//
48// Remark: Running it from the commandline looks like this:
49// root -q -l -b fillsignal.C+\(\"filename\"\,kFALSE\) 2>&1 | tee fillsignal.log
50//
51// Make sure, that database and password are corretly set in a resource
52// file called sql.rc and the resource file is found.
53//
54// Returns 2 in case of failure, 1 in case of success and 0 if the connection
55// to the database is not working.
56//
57/////////////////////////////////////////////////////////////////////////////
58#include <iostream>
59#include <iomanip>
60
61#include <TRegexp.h>
62
63#include <TH2.h>
64#include <TFile.h>
65#include <TGraph.h>
66#include <TSQLResult.h>
67
68#include "MSQLMagic.h"
69
70#include "MStatusArray.h"
71#include "MHCamera.h"
72#include "MHVsTime.h"
73
74#include "MCalibrationPulseTimeCam.h"
75#include "MCalibrationPix.h"
76
77using namespace std;
78
79int Process(MSQLMagic &serv, TString fname)
80{
81 TFile file(fname, "READ");
82 if (!file.IsOpen())
83 {
84 cout << "ERROR - Could not find file " << fname << endl;
85 return 2;
86 }
87
88 TString meanextpulpos("NULL");
89 TString rmsextpulpos("NULL");
90
91 MCalibrationPulseTimeCam *pt;
92 file.GetObject("MCalibrationPulseTimeCam", pt);
93 if (pt)
94 {
95 Double_t meanextpul = pt->GetAverageArea(0).GetHiGainMean();
96 Double_t rmsextpul = pt->GetAverageArea(0).GetHiGainRms();
97
98 if (meanextpul>=0 || rmsextpulpos>=0)
99 {
100 meanextpulpos.Form("%6.2f", meanextpul);
101 rmsextpulpos.Form("%6.2f", rmsextpul);
102 }
103 }
104
105
106 MStatusArray arr;
107 if (arr.Read()<=0)
108 {
109 cout << "ERROR - Reading of MStatusDisplay failed." << endl;
110 return 2;
111 }
112
113 MHCamera *cam = (MHCamera*)arr.FindObjectInCanvas("PedRMS;avg", "MHCamera", "PedRMS");
114 if (!cam)
115 {
116 cout << "ERROR - Reading of PedRMS;avg failed." << endl;
117 return 2;
118 }
119
120 MHCamera *cal = (MHCamera*)arr.FindObjectInCanvas("CalPos;avg", "MHCamera", "CalPos");
121 MHCamera *pul = (MHCamera*)arr.FindObjectInCanvas("PulsePos;avg", "MHCamera", "PulsePos");
122 if (!pul)
123 {
124 cout << "ERROR - Reading of PulsePos;avg failed." << endl;
125 return 2;
126 }
127/*
128 MHCamera *difflo = (MHCamera*)arr.FindObjectInCanvas("DiffLo;avg", "MHCamera", "DiffLo");
129 if (!difflo)
130 {
131 cout << "ERROR - Reading of DiffLo;avg failed." << endl;
132 return 2;
133 }
134 MHCamera *diffhi = (MHCamera*)arr.FindObjectInCanvas("DiffHi;avg", "MHCamera", "DiffHi");
135 if (!diffhi)
136 {
137 cout << "ERROR - Reading of DiffHi;avg failed." << endl;
138 return 2;
139 }
140*/
141
142 //get sequence number from the filename
143 TString sequence = fname(TRegexp("signal[0-9]+[.]root$"));
144 if (sequence.IsNull())
145 {
146 cout << "ERROR - Sequ# empty" << endl;
147 cout << "Sequ# " << sequence << endl;
148 return 2;
149 }
150
151 Int_t seq = atoi(sequence.Data()+6);
152
153 TString medpuloff("NULL");
154 TString devpuloff("NULL");
155 TString medhilocal("NULL");
156 TString devhilocal("NULL");
157
158 if (seq < 200000)
159 {
160 MHCamera *hilooff = (MHCamera*)arr.FindObjectInCanvas("HiLoOff;avg", "MHCamera", "HiLoOff");
161 if (!hilooff)
162 {
163 cout << "ERROR - Reading of HiLoOff failed." << endl;
164 return 2;
165 }
166
167 MHCamera *hilocal = (MHCamera*)arr.FindObjectInCanvas("HiLoCal;avg", "MHCamera", "HiLoCal");
168 if (!hilocal)
169 {
170 cout << "ERROR - Reading of HiLoCal failed." << endl;
171 return 2;
172 }
173
174 medpuloff.Form("%7.4f", hilooff->GetMedian());
175 devpuloff.Form("%7.4f", hilooff->GetDev());
176 medhilocal.Form("%6.2f", hilocal->GetMedian());
177 devhilocal.Form("%6.2f", hilocal->GetDev());
178 }
179
180 TArrayI inner(1);
181 inner[0] = 0;
182
183 TArrayI outer(1);
184 outer[0] = 1;
185
186 Int_t s0[] = { 1, 2, 3, 4, 5, 6 };
187
188 Stat_t meanrmsi = cam->GetMeanSectors(TArrayI(6, s0), inner);
189 Stat_t meanrmso = cam->GetMeanSectors(TArrayI(6, s0), outer);
190
191 if (meanrmsi<0 || meanrmso<0)
192 {
193 cout << "ERROR - MeanPedRMS inner or outer < 0 " << endl;
194 cout << "MeanPedRMS inner " << meanrmsi << endl;
195 cout << "MeanPedRMS outer " << meanrmso << endl;
196 return 2;
197 }
198
199 TString meanrmsinner=Form("%6.2f", meanrmsi);
200 TString meanrmsouter=Form("%6.2f", meanrmso);
201
202 cam = (MHCamera*)arr.FindObjectInCanvas("Interp'd;avg", "MHCamera", "Interp'd");
203 if (!cam)
204 {
205 cout << "ERROR - Reading of Interp'd;avg failed." << endl;
206 return 2;
207 }
208
209 Stat_t meansigi = cam->GetMeanSectors(TArrayI(6, s0), inner);
210 Stat_t meansigo = cam->GetMeanSectors(TArrayI(6, s0), outer);
211
212 if (meansigi<0 || meansigo<0)
213 {
214 cout << "ERRROR - MeanInterp'd inner or outer < 0 " << endl;
215 cout << "MeanInterp'd inner " << meansigi << endl;
216 cout << "MeanInterp'd outer " << meansigo << endl;
217 return 2;
218 }
219
220 TString meansiginner =Form("%6.2f", meansigi);
221 TString meansigouter =Form("%6.2f", meansigo);
222
223 TString calpos = cal ? Form("%5.1f", cal->GetMean()) : "NULL";
224
225 if (pul->GetMean()<0 || pul->GetRMS()<0)
226 {
227 cout << "ERROR - PulsePos'd mean or rms < 0 " << endl;
228 cout << "PulsePos'd mean " << pul->GetMean() << endl;
229 cout << "PulsePos'd rms " << pul->GetRMS() << endl;
230 return 2;
231 }
232
233 TString meanpulpos = Form("%6.2f", pul->GetMean());
234 TString rmspulpos = Form("%6.2f", pul->GetRMS());
235
236/*
237 Double_t meanhi = TMath::Nint(pulhi->GetMean()*100.)/100.;
238 Double_t rmshi = TMath::Nint(pulhi->GetRMS() *100.)/100.;
239
240 Double_t meanlo = TMath::Nint(pullo->GetMean()*100.)/100.;
241 Double_t rmslo = TMath::Nint(pullo->GetRMS() *100.)/100.;
242 pullo->Add(pullo, pulhi, 1, -1);
243 pullo->ResetBit(MHCamera::kProfile);
244
245 Double_t meanoff = TMath::Nint(pullo->GetMean()*100.)/100.;
246 Double_t rmsoff = TMath::Nint(pullo->GetRMS() *100.)/100.;
247
248 // USE MEDIAN INSTEAD? GetQuantiles(Int_t nprobSum, Double_t *q, const Double_t *probSum=0);
249
250 TString meanpulhi =Form("%6.2f", meanhi);
251 TString rmspulhi =Form("%6.2f", rmshi);
252
253 TString meanpullo =Form("%6.2f", meanlo);
254 TString rmspullo =Form("%6.2f", rmslo);
255 */
256
257 cam = (MHCamera*)arr.FindObjectInCanvas("Unsuitable;avg", "MHCamera", "Unsuitable");
258 if (!cam)
259 {
260 cout << "ERROR - Reading of Unsuitable;avg failed." << endl;
261 return 2;
262 }
263
264 Int_t unsuitable50 = cam->GetNumBinsAboveThreshold(0.50);
265 Int_t unsuitable01 = cam->GetNumBinsAboveThreshold(0.01);
266
267 TString unsuitablemax = "NULL";
268 TString deadmax = "NULL";
269
270 TGraph *gr = (TGraph*)arr.FindObjectInCanvas("BadPixTm", "TGraph", "BadPixTm");
271 if (gr)
272 {
273 const Int_t p = TMath::FloorNint(gr->GetN()*0.999);
274 unsuitablemax = Form("%d", TMath::Nint(TMath::KOrdStat(gr->GetN(), gr->GetY(), p)));
275 }
276
277 gr = (TGraph*)arr.FindObjectInCanvas("DeadPixTm", "TGraph", "DeadPixTm");
278 if (gr)
279 deadmax = Form("%d", TMath::Nint(TMath::MaxElement(gr->GetN(), gr->GetY())));
280
281 TString rateped = "NULL";
282 TString rateped2 = "NULL";
283 TString ratecal = "NULL";
284 TString ratetrig = "NULL";
285 TString ratesum = "NULL";
286 TString ratena = "NULL";
287 TString ratenull = "NULL";
288
289 TH2D *htp = (TH2D*)arr.FindObjectInCanvas("TrigPat", "TH2D", "TrigPat");
290 if (htp)
291 {
292 htp->ResetBit(TH1::kCanRebin);
293
294 Int_t iped = htp->GetYaxis()->FindBin("Ped");
295 Int_t iped2 = htp->GetYaxis()->FindBin("Ped+Trig");
296 Int_t ical = htp->GetYaxis()->FindBin("Cal");
297 Int_t itrig = htp->GetYaxis()->FindBin("Trig");
298 Int_t isum = htp->GetYaxis()->FindBin("Sum");
299 Int_t inull = htp->GetYaxis()->FindBin("0");
300 Int_t ina = htp->GetYaxis()->FindBin("UNKNOWN");
301
302 Int_t nx = htp->GetNbinsX();
303
304 rateped = iped <0 ? "NULL" : Form("%8.1f", htp->Integral(1, nx, iped, iped) / nx);
305 rateped2 = iped2<0 ? "NULL" : Form("%7.2f", htp->Integral(1, nx, iped2, iped2) / nx);
306 ratecal = ical <0 ? "NULL" : Form("%8.1f", htp->Integral(1, nx, ical, ical) / nx);
307 ratetrig = itrig<0 ? "NULL" : Form("%8.1f", htp->Integral(1, nx, itrig, itrig) / nx);
308 ratesum = isum <0 ? "NULL" : Form("%8.1f", htp->Integral(1, nx, isum, isum) / nx);
309 ratenull = inull<0 ? "NULL" : Form("%8.1f", htp->Integral(1, nx, inull, inull) / nx);
310 ratena = ina <0 ? "NULL" : Form("%7.2f", htp->Integral(1, nx, ina, ina) / nx);
311 }
312
313 // *****************************************************
314
315 cout << "Sequence #" << seq << endl;
316 cout << " Mean Ped RMS inner [phe] " << meanrmsinner << endl;
317 cout << " Mean Ped RMS outer [phe] " << meanrmsouter << endl;
318 cout << " Mean Signal inner [phe] " << meansiginner << endl;
319 cout << " Mean Signal outer [phe] " << meansigouter << endl;
320 cout << " Mean extracted PulsePos " << meanextpulpos << " +- " << rmsextpulpos << endl;
321 cout << " Mean calibrated PulsePos " << meanpulpos << " +- " << rmspulpos << endl;
322 cout << " Mean calib pulse pos " << calpos << endl;
323// cout << " Mean ext.HiGain PulsePos " << meanpulhi << " +- " << rmspulhi << endl;
324// cout << " Mean ext.LoGain PulsePos " << meanpullo << " +- " << rmspullo << endl;
325 cout << " Lo-Hi gain offset: " << medpuloff << " +- " << devpuloff << endl;
326 cout << " Hi/Lo gain ratio: " << medhilocal << " +- " << devhilocal << endl;
327 cout << " Unsuitable > 50%: " << setw(6) << unsuitable50 << endl;
328 cout << " Unsuitable > 1%: " << setw(6) << unsuitable01 << endl;
329 cout << " UnsuitableMax (99.9%) " << setw(6) << unsuitablemax << endl;
330 cout << " DeadMax " << setw(6) << deadmax << endl;
331 cout << " Rate Trigger [Hz] " << ratetrig << endl;
332 cout << " Rate SUM [Hz] " << ratesum << endl;
333 cout << " Rate Ped+Trigger [Hz] " << rateped2 << endl;
334 cout << " Rate Pedestal [Hz] " << rateped << endl;
335 cout << " Rate Calibration [Hz] " << ratecal << endl;
336 cout << " Rate 0 [Hz] " << ratenull << endl;
337 cout << " Rate UNKNOWN [Hz] " << ratena << endl;
338 cout << endl;
339
340 //build query and insert information into the database
341 // here only a update query is built, as this macro is exexuted after
342 // the macro fillcalib.C in the script fillcallisto
343 // and so the table Calibration is always updated
344 TString vars = Form(" fMeanPedRmsInner=%s, fMeanPedRmsOuter=%s, "
345 " fMeanSignalInner=%s, fMeanSignalOuter=%s, "
346 " fPulsePosMean=%s, fPulsePosRms=%s, "
347 " fPulsePosCheckMean=%s, fPulsePosCheckRms=%s, "
348 " fPulsePosCalib=%s, "
349 //" fPulsePosHiMean=%s, fPulsePosHiRms=%s, "
350 //" fPulsePosLoMean=%s, fPulsePosLoRms=%s, "
351 " fPulsePosOffMed=%s, fPulsePosOffDev=%s, "
352 " fHiLoGainRatioMed=%s, fHiLoGainRatioDev=%s, "
353 " fUnsuitable50=%d, fUnsuitable01=%d, "
354 " fUnsuitableMax=%s, fDeadMax=%s, "
355 " fRateTrigEvts=%s, fRateSumEvts=%s, "
356 " fRatePedTrigEvts=%s, fRatePedEvts=%s, "
357 " fRateCalEvts=%s, fRateNullEvts=%s, "
358 " fRateUnknownEvts=%s ",
359 meanrmsinner.Data(), meanrmsouter.Data(),
360 meansiginner.Data(), meansigouter.Data(),
361 meanpulpos.Data(), rmspulpos.Data(),
362 meanextpulpos.Data(), rmsextpulpos.Data(),
363 calpos.Data(),
364 //meanpulhi.Data(), rmspulhi.Data(),
365 //meanpullo.Data(), rmspullo.Data(),
366 medpuloff.Data(), devpuloff.Data(),
367 medhilocal.Data(), devhilocal.Data(),
368 unsuitable50, unsuitable01,
369 unsuitablemax.Data(), deadmax.Data(),
370 ratetrig.Data(), ratesum.Data(),
371 rateped2.Data(), rateped.Data(),
372 ratecal.Data(), ratenull.Data(),
373 ratena.Data()
374 );
375
376
377 return serv.Update("Calibration", vars, Form("fSequenceFirst=%d", seq)) ? 1 : 2;
378}
379
380int fillsignal(TString fname, Bool_t dummy=kTRUE)
381{
382 MSQLMagic serv("sql.rc");
383 if (!serv.IsConnected())
384 {
385 cout << "ERROR - Connection to database failed." << endl;
386 return 0;
387 }
388
389 cout << "fillsignal" << endl;
390 cout << "----------" << endl;
391 cout << endl;
392 cout << "Connected to " << serv.GetName() << endl;
393 cout << "File: " << fname << endl;
394 cout << endl;
395
396 serv.SetIsDummy(dummy);
397
398 //process file
399 return Process(serv, fname);
400}
Note: See TracBrowser for help on using the repository browser.