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

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