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

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