source: trunk/MagicSoft/Mars/datacenter/macros/fillstar.C@ 8203

Last change on this file since 8203 was 8203, checked in by tbretz, 18 years ago
*** empty log message ***
File size: 11.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, 05/2005 <mailto:tbretz@astro.uni-wuerzburg.de>
19! Author(s): Daniela Dorner, 05/2005 <mailto:dorner@astro.uni-wuerzburg.de>
20!
21! Copyright: MAGIC Software Development, 2000-2006
22!
23!
24\* ======================================================================== */
25
26/////////////////////////////////////////////////////////////////////////////
27//
28// fillstar.C
29// ==========
30//
31// This macro is used to read the star-output files.
32// These files are automatically called star00000000.root.
33// From these files the muon parameters (psf, muon number, ratio, muon rate),
34// the rate, the number of islands, the effective ontime, the maximum
35// humidity and a parameter describing the inhomogeneity are extracted from
36// the status display.
37// The sequence number is extracted from the filename.
38//
39// Usage:
40// .x fillstar.C("/magic/data/star/0004/00047421/star00047421.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 fillstar.C+\(\"filename\"\,kFALSE\) 2>&1 | tee fillstar.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 <TH1.h>
65#include <TH2.h>
66#include <TGraph.h>
67#include <TProfile.h>
68#include <TFile.h>
69#include <TSQLResult.h>
70#include <TSQLRow.h>
71
72#include "MSQLServer.h"
73
74#include "MHMuonPar.h"
75#include "MStatusArray.h"
76#include "MGeomCamMagic.h"
77#include "MBadPixelsCam.h"
78
79using namespace std;
80
81// --------------------------------------------------------------------------
82//
83// Checks whether an entry is already existing
84//
85Bool_t ExistStr(MSQLServer &serv, const char *column, const char *table, Int_t test)
86{
87 TString query(Form("SELECT %s FROM %s WHERE %s='%d'", column, table, column, test));
88 TSQLResult *res = serv.Query(query);
89 if (!res)
90 return kFALSE;
91
92 TSQLRow *row;
93
94 Bool_t rc = kFALSE;
95 while ((row=res->Next()))
96 {
97 if ((*row)[0])
98 {
99 rc = kTRUE;
100 break;
101 }
102 }
103
104 delete res;
105
106 return rc;
107}
108
109
110int Process(MSQLServer &serv, TString fname, Bool_t dummy)
111{
112 TFile file(fname, "READ");
113 if (!file.IsOpen())
114 {
115 cout << "ERROR - Could not find file " << fname << endl;
116 return 2;
117 }
118
119
120 MStatusArray arr;
121 if (arr.Read()<=0)
122 {
123 cout << "ERROR - Reading of MStatusDisplay failed." << endl;
124 return 2;
125 }
126
127 TH2F *hcog = (TH2F*)arr.FindObjectInCanvas("Center", "TH2F", "MHHillas");
128 if (!hcog)
129 {
130 cout << "WARNING - Reading of MHHillas failed." << endl;
131 return 2;
132 }
133
134 MHMuonPar *hmuon = (MHMuonPar*)arr.FindObjectInCanvas("MHMuonPar", "MHMuonPar", "MHMuonPar");
135 if (!hmuon)
136 {
137 cout << "WARNING - Reading of MHMuon failed." << endl;
138 return 2;
139 }
140
141 Double_t val[6];
142 for (int x=1; x<hcog->GetNbinsX(); x++)
143 for (int y=1; y<hcog->GetNbinsY(); y++)
144 {
145 Stat_t px = hcog->GetXaxis()->GetBinCenter(x);
146 Stat_t py = hcog->GetYaxis()->GetBinCenter(y);
147 Int_t i = (TMath::Nint(3*TMath::ATan2(px,py)/TMath::Pi())+6)%6;
148 val[i] += hcog->GetBinContent(x, y);
149 }
150
151 Double_t inhom = TMath::RMS(6, val)*6/hcog->GetEntries()*100;
152 inhom = TMath::Nint(inhom*10)/10.;
153 TString inhomogen = Form("%5.1f", inhom);
154
155 Float_t mw = hmuon->GetMeanWidth();
156 Float_t psf = 58.5480*mw - 37.3829;
157 psf = TMath::Nint(psf*10)/10.;
158 if (psf<0)
159 psf=0;
160 TString PSF = Form("%5.1f", psf);
161 Int_t num = (int)hmuon->GetEntries();
162
163 Float_t ratiodatamc = (hmuon->GetMeanSize()/9327.)*100;
164 TString ratio = Form("%5.1f", ratiodatamc);
165
166 TH1 *h = (TH1*)arr.FindObjectInCanvas("Islands", "TH1F", "MHImagePar");
167 if (!h)
168 {
169 cout << "WARNING - Reading of Islands failed." << endl;
170 return 2;
171 }
172
173 Float_t quality = h->GetMean();
174 quality = TMath::Nint(quality*100)/100.;
175 TString islands = Form("%6.2f", quality);
176
177 h = (TH1*)arr.FindObjectInCanvas("EffOnTheta", "TH1D", "EffOnTime");
178 if (!h)
179 {
180 cout << "WARNING - Reading of EffOnTheta failed." << endl;
181 return 2;
182 }
183
184 Float_t effon = h->Integral();
185 Float_t mrate = num/effon;
186 mrate = TMath::Nint(mrate*100)/100.;
187 if (mrate<0)
188 mrate=0;
189 TString muonrate = Form("%6.2f", mrate);
190 Int_t effontime = TMath::Nint(effon);
191
192 h = (TH1*)arr.FindObjectInCanvas("ProjDeltaT", "TH1D", "EffOnTime");
193 if (!h)
194 {
195 cout << "WARNING - Reading of ProjDeltaT failed." << endl;
196 return 2;
197 }
198
199 Int_t numevents = (int)h->GetEntries();
200 Int_t datarate = (int)(numevents/effon);
201
202 TGraph *g = (TGraph*)arr.FindObjectInCanvas("Humidity", "TGraph", "MHWeather");
203 if (!g)
204 {
205 cout << "WARNING - Reading of Humidity failed." << endl;
206 return 2;
207 }
208
209 Double_t max = TMath::MaxElement(g->GetN(), g->GetY());
210 TString maxhum = Form("%6.1f", max);
211
212
213 g = (TGraph*)arr.FindObjectInCanvas("NumStars", "TGraph", "MHPointing");
214 if (!g)
215 {
216 cout << "WARNING - Reading of NumStars failed." << endl;
217// return 2;
218 }
219
220 Double_t numstarmed = g ? TMath::Median(g->GetN(), g->GetY()) : -1;
221 TString numstarsmed = Form("%5.1f", numstarmed);
222 Double_t numstarrms = g ? g->GetRMS(2) : -1;
223 TString numstarsrms = Form("%5.1f", numstarrms);
224
225 g = (TGraph*)arr.FindObjectInCanvas("NumStarsCor", "TGraph", "MHPointing");
226 if (!g)
227 {
228 cout << "WARNING - Reading of NumStarsCor failed." << endl;
229// return 2;
230 }
231
232 Double_t numcormed = g ? TMath::Median(g->GetN(), g->GetY()) : -1;
233 TString numcorsmed = Form("%5.1f", numcormed);
234 Double_t numcorrms = g ? g->GetRMS(2) : -1;
235 TString numcorsrms = Form("%5.1f", numcorrms);
236
237 g = (TGraph*)arr.FindObjectInCanvas("Brightness", "TGraph", "MHPointing");
238 if (!g)
239 {
240 cout << "WARNING - Reading of SkyBrightness failed." << endl;
241// return 2;
242 }
243
244 Double_t brightnessmed = g ? TMath::Median(g->GetN(), g->GetY()) : -1;
245 TString skybrightnessmed = Form("%5.1f", brightnessmed);
246 Double_t brightnessrms = g ? g->GetRMS(2) : -1;
247 TString skybrightnessrms = Form("%5.1f", brightnessrms);
248
249
250 TString sequence = fname(TRegexp("star[0-9]+[.]root$"));
251 if (sequence.IsNull())
252 {
253 cout << "WARNING - Could not determine sequence# from filename: " << fname << endl;
254 return 2;
255 }
256
257 Int_t seq = atoi(sequence.Data()+5);
258
259 cout << "Sequence #" << seq << endl;
260 cout << " Inhomogeneity " << inhomogen << endl;
261 cout << " PSF [mm] " << PSF << endl;
262 cout << " Island Quality " << islands << endl;
263 cout << " Ratio [%] " << ratio << endl;
264 cout << " Muon Number " << num << endl;
265 cout << " EffOnTime [s] " << effontime << endl;
266 cout << " Muon Rate [Hz] " << muonrate << endl;
267 cout << " # of Events " << numevents << endl;
268 cout << " Rate after ImgCl [Hz] " << datarate << endl;
269 cout << " Maximum Humidity [%] " << maxhum << endl;
270 cout << " Number of Stars " << numstarsmed << " +/- " << numstarsrms << endl;
271 cout << " Number of cor. Stars " << numcormed << " +/- " << numcorrms << endl;
272 cout << " Skybrightness " << skybrightnessmed << " +/- " << skybrightnessrms << endl;
273
274 TString query;
275 if (!ExistStr(serv, "fSequenceFirst", "Star", seq))
276 {
277 query = Form("INSERT Star SET"
278 " fSequenceFirst=%d,"
279 " fMeanNumberIslands=%s, "
280 " fRatio=%s, "
281 " fMuonNumber=%d, "
282 " fEffOnTime=%d, "
283 " fMuonRate=%s, "
284 " fPSF=%s, "
285 " fDataRate=%d, "
286 " fMaxHumidity=%s ,"
287 " fNumStarsMed=%s ,"
288 " fNumStarsRMS=%s ,"
289 " fNumStarsCorMed=%s ,"
290 " fNumStarsCorRMS=%s ,"
291 " fBrightnessMed=%s ,"
292 " fBrightnessRMS=%s ,"
293 " fInhomogeneity=%s ",
294 seq, islands.Data(), ratio.Data(),
295 num, effontime,
296 muonrate.Data(), PSF.Data(),
297 datarate, maxhum.Data(),
298 numstarsmed.Data(), numstarsrms.Data(),
299 numcorsmed.Data(), numcorsrms.Data(),
300 skybrightnessmed.Data(), skybrightnessrms.Data(),
301 inhomogen.Data());
302 }
303 else
304 {
305 query = Form("UPDATE Star SET"
306 " fMeanNumberIslands=%s, "
307 " fRatio=%s, "
308 " fMuonNumber=%d, "
309 " fEffOnTime=%d, "
310 " fMuonRate=%s, "
311 " fPSF=%s, "
312 " fDataRate=%d, "
313 " fMaxHumidity=%s, "
314 " fNumStarsMed=%s ,"
315 " fNumStarsRMS=%s ,"
316 " fNumStarsCorMed=%s ,"
317 " fNumStarsCorRMS=%s ,"
318 " fBrightnessMed=%s ,"
319 " fBrightnessRMS=%s ,"
320 " fInhomogeneity=%s "
321 " WHERE fSequenceFirst=%d ",
322 islands.Data(), ratio.Data(),
323 num, effontime,
324 muonrate.Data(), PSF.Data(),
325 datarate, maxhum.Data(),
326 numstarsmed.Data(), numstarsrms.Data(),
327 numcorsmed.Data(), numcorsrms.Data(),
328 skybrightnessmed.Data(), skybrightnessrms.Data(),
329 inhomogen.Data(), seq);
330 }
331
332 cout << "Q: " << query << endl;
333
334 if (dummy)
335 return 1;
336
337 TSQLResult *res = serv.Query(query);
338 if (!res)
339 {
340 cout << "ERROR - Query failed: " << query << endl;
341 return 2;
342 }
343 delete res;
344
345 return 1;
346}
347
348int fillstar(TString fname, Bool_t dummy=kTRUE)
349{
350 TEnv env("sql.rc");
351
352 MSQLServer serv(env);
353 if (!serv.IsConnected())
354 {
355 cout << "ERROR - Connection to database failed." << endl;
356 return 0;
357 }
358
359 cout << "fillstar" << endl;
360 cout << "---------" << endl;
361 cout << endl;
362 cout << "Connected to " << serv.GetName() << endl;
363 cout << "File: " << fname << endl;
364 cout << endl;
365
366 return Process(serv, fname, dummy);
367}
Note: See TracBrowser for help on using the repository browser.