/[PAMELA software]/DarthVader/TrackerLevel2/src/TrkCore.cpp
ViewVC logotype

Contents of /DarthVader/TrackerLevel2/src/TrkCore.cpp

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1.35 - (show annotations) (download)
Fri Nov 28 09:00:37 2008 UTC (16 years, 3 months ago) by mocchiut
Branch: MAIN
CVS Tags: v6r01, v6r00
Changes since 1.34: +1 -1 lines
Event loss due to timesync bug fixed

1 /**
2 * \file TrkCore.cpp
3 * \author Elena Vannuccini
4 */
5 // .........................................................
6 // ROOT header files
7 // .........................................................
8 #include <TTree.h>
9 #include <TClassEdit.h>
10 #include <TObject.h>
11 #include <TList.h>
12 #include <TSystem.h>
13 #include <TSystemDirectory.h>
14 #include <TString.h>
15 #include <TFile.h>
16 #include <TClass.h>
17 #include <TCanvas.h>
18 #include <TH1.h>
19 #include <TH1F.h>
20 #include <TH2D.h>
21 #include <TLatex.h>
22 #include <TPad.h>
23 #include <TSQLServer.h>
24 #include <TSQLRow.h>
25 #include <TSQLResult.h>
26 #include <TBenchmark.h>
27
28 // .........................................................
29 // other general header files
30 // .........................................................
31 #include <fstream>
32 #include <iostream>
33 // .........................................................
34 // files in the inc directory
35 // .........................................................
36 #include <RunInfo.h>
37 #include <GLTables.h>
38 #include <TrkVerl2.h>
39 #include <TrkProcess.h>
40 //#include <TrkStruct.h>
41 #include <TrkLevel2.h>
42 #include <TrkLevel1.h>
43 #include <TrkLevel0.h>
44 #include <TrkHough.h>
45 // .........................................................
46 // YODA header files
47 // .........................................................
48 #include <PamelaRun.h>
49 //#include <RegistryEvent.h>
50 #include <physics/tracker/TrackerEvent.h>
51 #include <CalibTrk1Event.h>
52 #include <CalibTrk2Event.h>
53
54 // .........................................................
55 // namespaces
56 // .........................................................
57 using namespace std;
58 using namespace pamela;
59 // ================================================================================
60 //
61 //
62 // ================================================================================
63 /**
64 * \brief Tracker data reduction routine.
65 *
66 * It performs data reduction LEVEL0->LEVEL1->LEVEL2, producing LEVEL1 and or LEVEL2 output, in ROOT or HBOOK format.
67 * Input parameters:
68 * @param run id of the run to be processed (if run=0 a whole file is reprocessed)
69 * @param dbc pointer to BD server
70 * @param file1 output LEVEL1 file name (null if no LEVEL1 output is required)
71 * @param file2 output LEVEL2 file name (null if no LEVEL2 output is required)
72 * @param standalone (flag to run the program in standalone mode, that is without reading RunInfo)
73 */
74 //short int TrkCore(Int_t run, TSQLServer *dbc, TString file1, TString file2, Bool_t standalone)
75 //int TrkCore(UInt_t run, TFile *f2, TSQLServer *dbc, int ncustom, char*vcustom[]) // EMILIANO
76 int TrkCore(UInt_t run, TFile *f2, GL_TABLES *glt, int ncustom, char*vcustom[]) // EMILIANO
77 {
78 // ---------------------------
79 // Define some basic variables
80 // ---------------------------
81
82 TString filename = 0;
83 Long64_t nentries = 0LL;
84
85 // LEVEL0 input
86 TFile *f0 = 0;
87 TTree *physicsTree = 0;
88 TBranch *b_trk = 0;
89 TBranch *b_header = 0;
90 EventHeader *header = 0;
91 PscuHeader *pscu = 0;
92 TrkLevel0 *l0_event = new TrkLevel0();
93 // RunInfo
94 ItoRunInfo *runinfo = 0;
95 TArrayI *runlist = 0;
96 UInt_t from_run;
97 UInt_t to_run;
98 Bool_t reprocessing = false;
99 //LEVEL2 output - ROOT
100 TTree *t_level2 = 0;
101 TTree *t_clone = 0;
102 TrkLevel2 *l2_event = new TrkLevel2();
103 TrkLevel2 *l2_clone = new TrkLevel2();
104 //LEVEL1 output - ROOT
105 TrkLevel1 *l1_event = new TrkLevel1();
106 TrkHough *lh_event = new TrkHough();
107 // -----------------------
108 // -----------------------
109 // Handle input parameters
110 // (data)
111 //
112 // - open/create output files, determining the environment root/hbook from the estension
113 // - create the level1/level2 tree+branch/nt-uple
114 // -----------------------
115 // -----------------------
116 TrkProcess *p = new TrkProcess(run,f2);
117 if( p->HandleCustomPar(ncustom,vcustom) )return(p->ostatus);;
118 if( TrkParams::VerboseMode() )p->Dump();
119
120 // ===================================
121 // Open/Create level2 output file
122 // ===================================
123 if(p->get2){
124 //-------------------------------------------
125 // read RunInfo
126 //-------------------------------------------
127 if(!(p->standalone)){
128 // Open "RunInfo" tree and get list of run
129 runinfo = new ItoRunInfo(f2);
130 char *trkversion = TrkInfo(false);
131 int runinfo_error = runinfo->Update(run,"TRK",trkversion);
132 if( runinfo_error ) throw runinfo_error;
133 runlist = runinfo->GetRunList();
134 reprocessing = runinfo->IsReprocessing();//????
135 if(TrkParams::VerboseMode()){
136 cout << "#events "<< runinfo->GetFileEntries() << endl;// #events in the file
137 cout << "#runs "<< runinfo->GetRunEntries() << endl;// #runs in the file
138 cout << "1-st run "<< runlist->At(0) << endl;
139 cout << "last run "<< runlist->At(runinfo->GetRunEntries()-1) << endl;
140 cout << "1-st event "<< runinfo->GetFirstEntry() << endl;// first event of our run
141 cout << "last event+1 "<< runinfo->GetLastEntry() << endl;// first event AFTER the last event of our run
142 cout << "reprocessing "<< runinfo->IsReprocessing() << endl << endl;
143 };
144 };
145 //-------------------------------------------
146 //
147 //-------------------------------------------
148 // Take (if present) the old Tracker tree
149 t_clone = (TTree*)f2->Get("Tracker");
150 if( t_clone != NULL ) t_clone->SetBranchAddress("TrkLevel2",&l2_clone);
151 // Create NEW Tracker tree
152 t_level2 = new TTree("Tracker","PAMELA tracker level2 data ");
153 l2_event->Set(); // ****NBNBNBN*****
154 t_level2->Branch("TrkLevel2","TrkLevel2",&l2_event);
155 if(p->get1){
156 if(TrkParams::VerboseMode())cout << endl << "Requested LEVEL1 output" << endl;
157 l1_event->Set(); // ****NBNBNBN*****
158 t_level2->Branch("TrkLevel1","TrkLevel1",&l1_event);
159 }
160 if(p->geth){
161 if(TrkParams::VerboseMode())cout << endl << "Requested Hough-Transform output" << endl;
162 t_level2->Branch("TrkHough","TrkHough",&lh_event);
163 };
164 };
165
166
167 // -------------------------------------------
168 // define runs to be processed/reprocessed
169 // -------------------------------------------
170 if(run == 0){
171 // reprocessing ALL runs
172 if(p->standalone)throw -298; // reprocessing not implemented
173 from_run = runlist->At(0);
174 to_run = runlist->At(runinfo->GetRunEntries()-1);
175 }else{
176 // processing/reprocessing ONE single run
177 from_run = run;
178 to_run = run;
179 };
180
181 //
182 // init event counter
183 //
184 Int_t ev_count =0;
185 //
186 // create query-results objects
187 //
188 GL_RUN q1 = GL_RUN();
189 GL_TRK_CALIB q2 = GL_TRK_CALIB();
190 GL_ROOT q3 = GL_ROOT();
191 GL_PARAM q4 = GL_PARAM();
192
193 // ------------------------------------------------------------
194 // if reprocessing one run, copy all the events BEFORE the run
195 // ------------------------------------------------------------
196 if( !(p->standalone) ){
197 for(UInt_t i=0; i<runinfo->GetFirstEntry(); i++){
198 t_clone->GetEntry(i);
199 *l2_event = *l2_clone;
200 t_level2->Fill();
201 l2_event->Clear();
202 // COPY COPY COPY
203 };
204 };
205 // ------------------------------------------------------------
206 // ------------------------------------------------------------
207 // START LOOP OVER RUNS TO PROCESS/REPROCESS
208 // ------------------------------------------------------------
209 // ------------------------------------------------------------
210 for(UInt_t idRun = from_run; idRun <= to_run; idRun++){
211
212 if(TrkParams::VerboseMode()) cout << endl<<" ========================= Run: "<< idRun << endl;
213 UInt_t runheadtime = 0;
214 UInt_t runtrailtime = 0;
215 UInt_t runheadobt = 0;
216 UInt_t runtrailobt = 0;
217 UInt_t evfrom = 0;
218 UInt_t evto = 0;
219 UInt_t nevents = 0;
220 Int_t id_root_l0 =-1;
221 Int_t trk_calib_used = 0;
222
223 TString host = glt->CGetHost(); // EMILIANO
224 TString user = glt->CGetUser(); // EMILIANO
225 TString psw = glt->CGetPsw(); // EMILIANO
226 // TString host = "mysql://localhost/pamelaprod";
227 // TString user = "anonymous";
228 // TString psw = "";
229 // const char *pamdbhost=gSystem->Getenv("PAM_DBHOST");
230 // const char *pamdbuser=gSystem->Getenv("PAM_DBUSER");
231 // const char *pamdbpsw=gSystem->Getenv("PAM_DBPSW");
232 // if ( !pamdbhost ) pamdbhost = "";
233 // if ( !pamdbuser ) pamdbuser = "";
234 // if ( !pamdbpsw ) pamdbpsw = "";
235 // if ( strcmp(pamdbhost,"") ) host = pamdbhost;
236 // if ( strcmp(pamdbuser,"") ) user = pamdbuser;
237 // if ( strcmp(pamdbpsw,"") ) psw = pamdbpsw;
238 TSQLServer *dbc = TSQLServer::Connect(host.Data(),user.Data(),psw.Data());
239 stringstream myquery; // EMILIANO
240 myquery.str(""); // EMILIANO
241 myquery << "SET time_zone='+0:00'"; // EMILIANO
242 dbc->Query(myquery.str().c_str()); // EMILIANO
243 if(p->standalone){
244 // ==============================
245 // first query: retrieve run info
246 // ==============================
247 if (q1.Query_GL_RUN(idRun,dbc) )throw -50;
248 id_root_l0 = q1.ID_ROOT_L0;
249 runheadtime = q1.RUNHEADER_TIME;
250 runtrailtime = q1.RUNTRAILER_TIME;
251 runheadobt = q1.RUNHEADER_OBT;
252 runtrailobt = q1.RUNTRAILER_OBT;
253 evfrom = q1.EV_FROM;
254 evto = q1.EV_TO;
255 nevents = q1.NEVENTS;
256 trk_calib_used = q1.TRK_CALIB_USED;
257 }else{
258 // ==============================
259 // get run info from RunInfo tree
260 // ==============================
261 int runinfo_error = runinfo->GetRunInfo(idRun);
262 if( runinfo_error ) throw runinfo_error;
263 id_root_l0 = runinfo->ID_ROOT_L0;
264 runheadtime = runinfo->RUNHEADER_TIME;
265 runtrailtime = runinfo->RUNTRAILER_TIME;
266 runheadobt = runinfo->RUNHEADER_OBT;
267 runtrailobt = runinfo->RUNTRAILER_OBT;
268 evfrom = runinfo->EV_FROM;
269 evto = runinfo->EV_TO;
270 nevents = runinfo->NEVENTS;
271 trk_calib_used = runinfo->TRK_CALIB_USED;
272
273 };
274
275 //
276 if(TrkParams::VerboseMode()){
277 cout << "ROOT file ID "<< id_root_l0 << endl;
278 cout << "RunHeader time "<< runheadtime << endl;
279 cout << "RunTrailer time "<< runtrailtime << endl;
280 cout << " from event "<< evfrom << endl;
281 cout << " to event "<< evto << endl;
282 cout << " num. events "<< nevents << endl;
283 cout << "trk-calibration used "<< trk_calib_used << endl;
284 };
285 // ========================================================
286 // second query: search the LEVEL0 file that contains idRun
287 // ========================================================
288 TString lastfilename = filename;
289 if( q3.Query_GL_ROOT(id_root_l0,dbc) )throw -51;
290 filename = q3.PATH + q3.NAME;
291 // ========================================================
292 // Open the input LEVEL0 data file
293 // ========================================================
294 if(filename.CompareTo(lastfilename)){
295 if(!lastfilename.IsNull())f0->Close();
296 //if( debug ) cout << "Opening LEVEL0 file: "<< filename << endl;
297 if(TrkParams::VerboseMode()) cout << "Opening LEVEL0 file: "<< filename << endl;
298 FileStat_t t;
299 if( gSystem->GetPathInfo(filename.Data(),t) )throw -6;
300 f0 = new TFile(filename);
301 if ( !f0 ) throw -6;
302 physicsTree = (TTree*)f0->Get("Physics"); if(!physicsTree) throw -7;
303 b_header = physicsTree ->GetBranch("Header"); if(!b_header ) throw -8;
304 b_trk = physicsTree ->GetBranch("Tracker"); if(!b_trk ) throw -203;
305 l0_event->Set();
306 physicsTree->SetBranchAddress("Tracker" ,l0_event->GetPointerToTrackerEvent());
307 physicsTree->SetBranchAddress("Header",&header);
308
309 nentries = physicsTree->GetEntries();
310 if ( nentries < 1 && nevents ) throw -11; // EMILIANO: go on if the file is empty, why not? in principle we should not have any event to process but if the case the next line will take care of exiting; exit only if the file is empty but nevents is not zero
311 if ( nentries < (evto+1) && nevents > 0 ) throw -12; // EMILIANO: if NEVENTS = 0 and the file is empty everything is ok but due to a mistake in the
312 // variable EV_TO type (UInt_t instead of Int_t) that we don't want to change to avoid changing a lot of code, evto becomes inf hence
313 // the condition ( nentries < (evto+1) ) is satisfied and DV exit with error even if the error is only in the DB.
314
315 };
316
317 GL_TIMESYNC *dbtime = new GL_TIMESYNC(id_root_l0,"ID",dbc);
318
319
320 // =============================================================
321 // retrieve information about parameters to process LEVEL2
322 // =============================================================
323
324 TrkParams::Set(runinfo->GetGL_RUN(),dbc);
325 for(int i=0; i<p->npar; i++){
326 if(TrkParams::VerboseMode())cout<<" ((( force parameters from input path )))"<<endl;
327 TrkParams::Set(p->parpath[i],p->partype[i]);
328 }
329
330 TrkParams::Load();
331 if( !TrkParams::IsLoaded() )throw -52;
332
333 // =============================================================
334 // retrieve calibration file needed to reduce data
335 // =============================================================
336
337 TrkParams::SetCalib(runinfo->GetGL_RUN(),dbc);
338 TrkParams::LoadCalib( );
339 if( !TrkParams::CalibIsLoaded() )throw -53;
340
341 TBenchmark *reduction = new TBenchmark();
342 if(TrkParams::VerboseMode())reduction->Start("reduction");
343 Int_t ntrk=0;
344 // ====================================================
345 // start looping on events cointained in the data file
346 // ====================================================
347 if(dbc){
348 dbc->Close();
349 delete dbc;};
350 //
351 for (UInt_t re = evfrom+min(p->nskip,nevents); re < evfrom+nevents; re++){
352
353 ev_count++;
354
355
356 // if ( TrkParams::DebugMode() && re%100 == 0 && re > 0 ) cout << ".";
357
358 b_trk->GetEntry(re);
359 b_header->GetEntry(re);
360 pscu = header->GetPscuHeader();
361
362 if( TrkParams::DebugMode() )cout << ">>> "<<ev_count-1<<" @ OBT "<<pscu->GetOrbitalTime()<<endl;
363
364 if ( dbtime->DBabsTime(pscu->GetOrbitalTime()) > (runtrailtime+1) || dbtime->DBabsTime(pscu->GetOrbitalTime()) < (runheadtime-1)) {
365
366 if (TrkParams::VerboseMode()){
367 printf(" TrkCore - WARNING: event outside the run time window, skipping it\n");
368 cout << " OBT "<<pscu->GetOrbitalTime()<<" RUN "<<runheadobt<<"-"<<runtrailobt<<" ABS-time "<<dbtime->DBabsTime(pscu->GetOrbitalTime())<<" RUN "<<runheadtime<<"-"<<runtrailtime<<endl;
369 };
370 }else{
371 if ( TrkParams::DebugMode() )
372 printf("\n-------------------------------\nevent %d\n",re-evfrom);
373
374 p->ProcessEvent(l0_event);
375
376 // ----------------
377 // LEVEL1 output
378 // ----------------
379 if(p->get1){
380 if(p->ifroot1){ // root
381 l1_event->Clear();
382 // l1_event->SetFromLevel1Struct(&level1event_,p->full1);
383 l1_event->SetFromLevel1Struct(p->full1);
384 // t_level1->Fill();
385 }else{ // hbook
386 throw -299;
387 };
388 };
389 // ----------------
390 // HOUGH output
391 // ----------------
392 if(p->geth){
393 if(p->ifrooth){ // root
394 lh_event->Delete();
395 lh_event->SetFromHoughStruct(&houghevent_);
396 }else{ // hbook
397 throw -299;
398 };
399 };
400 // ----------------
401 // LEVEL2 output
402 // ----------------
403 if(p->get2){
404 l2_event->Clear();
405 if(p->get1){
406 l2_event->SetFromLevel2Struct(&level2event_,l1_event);//set references to level1
407 }else{
408 l2_event->SetFromLevel2Struct(&level2event_);
409 }
410 // l2_event->Dump();
411 t_level2->Fill();
412 if( l2_event->ntrk()>0 )ntrk++;
413 if(TrkParams::VerboseMode())l2_event->Dump();
414 };
415 };
416 }; // end loop on events
417 if(TrkParams::VerboseMode()){
418 cout << " >>> processed "<< ev_count <<" events"<< endl;
419 if(ev_count)cout << ntrk << " events with at least one track ("<<(Int_t)(100*ntrk/ev_count)<<"%)\n";
420 reduction->Show("reduction");
421 }
422 delete reduction;
423
424 delete dbtime;
425
426 }; // end loop on runs
427
428
429 // ------------------------------------------------------------
430 // if reprocessing one run, copy all the events AFTER the run
431 // ------------------------------------------------------------
432 if( !(p->standalone) ){
433 for(UInt_t i=runinfo->GetLastEntry()+1; i<runinfo->GetFileEntries(); i++){
434 t_clone->GetEntry(i);
435 *l2_event = *l2_clone;
436 t_level2->Fill();
437 l2_event->Clear();
438 // COPY COPY COPY
439 };
440 };
441 // ---------------
442 // close the files
443 // ---------------
444 if(p->get2){
445 if( t_clone )t_clone->Delete("all");//delete old tree from file
446 if( !(p->standalone) )runinfo->Close();
447 f2->Write("Tracker");
448 if( t_level2 )t_level2->Delete(); //delete new tree from memory
449
450 };
451
452 if(f0) if(f0->IsOpen()) f0->Close();
453 // if( f0_c->IsOpen() )f0_c->Close();
454
455 // lh_event->Delete();
456 l1_event->Delete();
457 l2_event->Delete();
458 l2_clone->Delete();
459
460 return(p->ostatus);
461 }
462

  ViewVC Help
Powered by ViewVC 1.1.23