source: trunk/Mars/fact/processing/numevents.C

Last change on this file was 19268, checked in by tbretz, 3 years ago
Fixed another issue pointed out by the YIT compiler.
File size: 9.4 KB
Line 
1#ifndef __CLING__
2#include <iostream>
3
4#include <TCanvas.h>
5#include <TString.h>
6#include <TTree.h>
7#include <TFile.h>
8#include <TSystem.h>
9#include <TSQLResult.h>
10#include <TSQLRow.h>
11#include <TPRegexp.h>
12
13#include "MSQLMagic.h"
14#include "MAlphaFitter.h"
15#include "MHThetaSq.h"
16#include "MStatusArray.h"
17#endif
18
19using namespace std;
20
21int process(MSQLMagic &serv, TString gfile, TString sfile, TString night2, TString runid, TString table, Bool_t pernight, Int_t source)
22{
23 TString where;
24 TString vars;
25 //the variables to save into the db
26 Int_t fNumEvtsAfterCleaning = 0;
27 Int_t fNumEvtsAfterQualCuts = 0;
28 Int_t fNumEvtsAfterBgCuts = 0;
29 Float_t fNumBgEvts = 0;
30 Float_t fNumSigEvts = 0;
31 Float_t fNumExcEvts = 0;
32 Float_t fNumIslandsMean = 0;
33 Float_t fOnTimeAfterCuts = 0;
34
35 if (!pernight)
36 {
37 //check star file status
38 if (gSystem->AccessPathName(sfile))
39 {
40 cout << "ERROR - file " << sfile << " does not exist." << endl;
41 return 2;
42 }
43
44 //try to open the star file
45 TFile star_file(sfile, "READ");
46 if (!star_file.IsOpen())
47 {
48 cout << "ERROR - Could not open file " << sfile << endl;
49 return 2;
50 }
51
52 //Get the events tree
53 TTree *star_tree = (TTree *)star_file.Get("Events");
54 if (!star_tree)
55 {
56 cout << "ERROR - Could not read tree " << sfile << endl;
57 return 2;
58 }
59
60 //the number of events after cleaning
61 fNumEvtsAfterCleaning = star_tree->GetEntries();
62
63 //get the mean number of islands in the _I.root file on runbasis
64 star_tree->Draw("MImagePar.fNumIslands>>IslandHisto","","");
65 TH1F *HistJo = (TH1F*)gDirectory->Get("IslandHisto");
66 if (!HistJo)
67 {
68 cout << "ERROR - Reading of IslandHisto failed " << HistJo << endl;
69 return 2;
70 }
71 fNumIslandsMean = HistJo->GetMean();
72 }
73
74 //check ganymed file status
75 if (gSystem->AccessPathName(gfile))
76 {
77 cout << "ERROR - file " << gfile << " does not exist." << endl;
78 return 2;
79 }
80
81 //try to open the ganymed file
82 TFile ganymed_file(gfile, "READ");
83 if (!ganymed_file.IsOpen())
84 {
85 cout << "ERROR - Could not open file " << gfile << endl;
86 return 2;
87 }
88
89 //Get the status display contents
90 MStatusArray arr;
91 if (arr.Read()<=0)
92 {
93 cout << "ERROR - could not read MStatusArray, file " << gfile << endl;
94 return 2;
95 }
96
97 TH1D *vstime = (TH1D*)arr.FindObjectInCanvas("OnTime", "TH1D", "OnTime");
98 if (!vstime)
99 {
100 cout << "WARNING - Reading of Theta failed." << endl;
101 return 2;
102 }
103 fOnTimeAfterCuts = vstime->Integral();
104
105 //Get number of events after quality cuts
106 //from the number of entries in a histogram
107 TH1F *precut =
108 (TH1F*)arr.FindCanvas("PreCut")->FindObject("SizeSame");
109 if (!precut)
110 {
111 cout << "WARNING - Reading of PreCut canvas failed." << endl;
112 return 2;
113 }
114 fNumEvtsAfterQualCuts = precut->GetEntries();
115
116 //Get number of events after background rejection cuts
117 //from the number of entries in a histogram
118 TH1F *postcut =
119 (TH1F*)arr.FindCanvas("PostCut")->FindObject("SizeSame");
120 if (!postcut)
121 {
122 cout << "WARNING - Reading of PostCut canvas failed." << endl;
123 return 2;
124 }
125 fNumEvtsAfterBgCuts = postcut->GetEntries();
126
127 //Get signal, bg and excess events
128 MHThetaSq *halphaon = (MHThetaSq*)arr.FindObjectInCanvas("MHThetaSq", "Hist");
129 if (!halphaon)
130 {
131 cout << "WARNING - Reading of MHThetaSq failed, file " << gfile << endl;
132 return 2;
133 }
134 const MAlphaFitter *fit = &halphaon->GetAlphaFitter();
135 if (!fit)
136 {
137 cout << "WARNING - Reading of MAlphaFitter failed, file " << gfile << endl;
138 return 2;
139 }
140
141 //check if the excess events number is small and close to 0
142 //due to the computer precision error when subtracting sig and bg events
143 fNumBgEvts = fit->GetEventsBackground();
144 fNumSigEvts = fit->GetEventsSignal();
145 fNumExcEvts = fit->GetEventsExcess();
146 fNumExcEvts = fabs(fNumExcEvts)<1e-5 ? 0 : fNumExcEvts;
147
148 /*
149 cout<< fNumEvtsAfterCleaning<< " "
150 << fNumEvtsAfterQualCuts<< " "
151 << fNumEvtsAfterBgCuts << " "
152 << fNumBgEvts << " "
153 << fNumSigEvts << " "
154 << fNumExcEvts << endl;
155 */
156
157 //inserting or updating the information in the database
158 if (!pernight)
159 vars = Form("fRunID=%s, fNight=%s,", runid.Data(), night2.Data());
160 else
161 vars = Form("fNight=%s, fSourceKey=%d, ", night2.Data(), source);
162 vars += Form(" fNumEvtsAfterQualCuts=%d, "
163 " fNumEvtsAfterBgCuts=%d, "
164 " fNumBgEvts=%6.1f, "
165 " fNumSigEvts=%6.1f, "
166 " fNumExcEvts=%6.1f, "
167 " fOnTimeAfterCuts=%7.2f ",
168 fNumEvtsAfterQualCuts,
169 fNumEvtsAfterBgCuts,
170 fNumBgEvts,
171 fNumSigEvts,
172 fNumExcEvts,
173 fOnTimeAfterCuts
174 );
175 if (!pernight)
176 vars += Form(", fNumIslandsMean=%6.2f, fNumEvtsAfterCleaning=%d ",
177 fNumIslandsMean,
178 fNumEvtsAfterCleaning
179 );
180
181
182 if (pernight)
183 where = Form("fNight=%s AND fSourceKey=%d", night2.Data(), source);
184 else
185 where = Form("fRunID=%s AND fNight=%s", runid.Data(), night2.Data());
186
187 cout << "vars: " << vars << endl;
188 cout << "where: " << where << endl;
189
190 if (!serv.InsertUpdate(table, vars, where))
191 {
192 cout << "ERROR - insert/update to DB failed (vars:" << vars << " where: " << where << ")." << endl;
193 return 2;
194 }
195 return 1;
196}
197
198
199int numevents(Int_t night, TString inpath="/loc_data/analysis/", TString table= "AnalysisResultsRunLP", Bool_t dummy=kTRUE, Bool_t pernight=kFALSE, Int_t source=0)
200{
201
202 //connect to mysql server
203 //MSQLServer serv("sql.rc");
204 MSQLMagic serv("sql.rc");
205 if (!serv.IsConnected())
206 {
207 cout << "ERROR - Connection to database failed." << endl;
208 return 0;
209 }
210 serv.SetIsDummy(dummy);
211
212 //select runs, where star finished, format the list nicely
213 TString query;
214 query = " SELECT CONCAT(fNight, '_', LPAD(fRunID, 3, 0)), SUBSTRING(fNight, 1,4), SUBSTRING(fNight, 5,2), SUBSTRING(fNight, 7,2), ";
215 query += " fNight, fRunID FROM RunInfo";
216
217 //and only data runs
218 query+=" WHERE fRunTypeKEY=1 ";
219 query+=Form(" AND fNight=%d",night);
220 //query+=" AND fRunID=128 ";
221 if (pernight)
222 query+=" GROUP BY fNight";
223
224 cout << "Q: " << query << endl;
225
226 //post the query
227 TSQLResult *res = serv.Query(query);
228 if (!res)
229 {
230 cout << "ERROR - no result from query " << query << endl;
231 return 2;
232 }
233
234 //allocate variables to use in the loop
235 TString night2;
236 TString runid;
237 TString run;
238 TString year;
239 TString month;
240 TString day;
241 TString ganymed_fname;
242 TString star_fname;
243 Int_t rc=0;
244 Int_t numpb=0;
245
246 //loop over the data files
247 TSQLRow *row=0;
248
249 //cout the data names
250 cout<<"#evts_after_cleaning" << " "
251 <<"evts_after_qual-cuts" << " "
252 <<"evts_after_bg-cuts" << " "
253 <<"mean_evts_from_off_regions"<< " "
254 <<"evts_from_on_region" << " "
255 <<"calculated_excess_evts_in_file" << " " << endl;
256
257 while ((row=res->Next()))
258 {
259 //use the results from the query
260 night2=(*row)[4];
261 runid=(*row)[5];
262 run=(*row)[0];
263 year=(*row)[1];
264 month=(*row)[2];
265 day=(*row)[3];
266
267 // to be fixed: path for nightly ganymeds
268 if (!pernight)
269 ganymed_fname=inpath+"/ganymed_run/"+year+"/"+month+"/"+day+"/"+run+"-analysis.root";
270 else
271 if (!source)
272 cout << "ERROR - source " << source << " not valid. " << endl;
273 else
274 ganymed_fname=Form("%s/ganymed_night/%d/%s-analysis.root", inpath.Data(), source, night2.Data());
275 cout << "gf => " << ganymed_fname << endl;
276
277 if (!pernight)
278 {
279 star_fname=inpath+"/star/"+year+"/"+month+"/"+day+"/"+run+"_I.root";
280 cout << "sf => " << star_fname << endl;
281 }
282 rc=process(serv, ganymed_fname, star_fname, night2, runid, table, pernight, source);
283 if (rc>1)
284 numpb++;
285 //return rc;
286
287 }
288
289
290 if (numpb>0)
291 return numpb+1;//avoid that 1 (i.e. everything ok) is returned in case of problems
292 else
293 return 1;
294}
295
296int numevents(TString ganymedfile, TString starfile, TString table, Bool_t dummy=kTRUE)
297{
298
299 //connect to mysql server
300 MSQLMagic serv("sql.rc");
301 if (!serv.IsConnected())
302 {
303 cout << "ERROR - Connection to database failed." << endl;
304 return 0;
305 }
306 serv.SetIsDummy(dummy);
307 //get night and run from file name
308 TString basename=gSystem->BaseName(ganymedfile);
309 TPRegexp regexp("20[0-9][0-9][0-2][0-9][0-3][0-9]_[0-9][0-9][0-9]");
310 TString run = basename(regexp);
311 TPRegexp regexp2("20[0-9][0-9][0-2][0-9][0-3][0-9]");
312 TPRegexp regexp3("_[0-9][0-9][0-9]");
313 TPRegexp regexp4("[0-9][0-9][0-9]");
314 TString night = run(regexp2);
315 TString runnum = run(regexp3);
316 TString runid = runnum(regexp4);
317
318 return process(serv, ganymedfile, starfile, night, runid, table, false, 0);
319}
320
Note: See TracBrowser for help on using the repository browser.