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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1.6 - (hide annotations) (download)
Tue Sep 5 12:52:20 2006 UTC (18 years, 5 months ago) by pam-fi
Branch: MAIN
Changes since 1.5: +92 -99 lines
implemented class TrkLevel1

1 mocchiut 1.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     // .........................................................
27     // other general header files
28     // .........................................................
29     #include <fstream>
30     #include <iostream>
31     // .........................................................
32     // files in the inc directory
33     // .........................................................
34     #include <RunInfo.h>
35     #include <GLTables.h>
36     #include <TrkVerl2.h>
37     #include <TrkProcess.h>
38 pam-fi 1.2 //#include <TrkStruct.h>
39 mocchiut 1.1 #include <TrkLevel2.h>
40 pam-fi 1.6 #include <TrkLevel1.h>
41 mocchiut 1.1 #include <TrkLevel0.h>
42     // .........................................................
43     // YODA header files
44     // .........................................................
45     #include <PamelaRun.h>
46     #include <RegistryEvent.h>
47     #include <physics/tracker/TrackerEvent.h>
48     #include <CalibTrk1Event.h>
49     #include <CalibTrk2Event.h>
50 pam-fi 1.2
51 mocchiut 1.1 // .........................................................
52     // namespaces
53     // .........................................................
54     using namespace std;
55     using namespace pamela;
56     // ================================================================================
57     //
58     //
59     // ================================================================================
60     /**
61     * \brief Tracker data reduction routine.
62     *
63     * It performs data reduction LEVEL0->LEVEL1->LEVEL2, producing LEVEL1 and or LEVEL2 output, in ROOT or HBOOK format.
64     * Input parameters:
65     * @param run id of the run to be processed (if run=0 a whole file is reprocessed)
66     * @param dbc pointer to BD server
67     * @param file1 output LEVEL1 file name (null if no LEVEL1 output is required)
68     * @param file2 output LEVEL2 file name (null if no LEVEL2 output is required)
69     * @param standalone (flag to run the program in standalone mode, that is without reading RunInfo)
70     */
71     //short int TrkCore(Int_t run, TSQLServer *dbc, TString file1, TString file2, Bool_t standalone)
72     int TrkCore(ULong64_t run, TFile *f2, TSQLServer *dbc, int ncustom, char*vcustom[])
73     {
74     // ---------------------------
75     // Define some basic variables
76     // ---------------------------
77    
78     TString filename = 0;
79     TString trkcalibfile = 0;
80     Int_t last_trk_calib_used = 0;
81 pam-fi 1.2 Long64_t nevents = 0LL;
82    
83 mocchiut 1.1 // LEVEL0 input
84     TFile *f0 = 0;
85 pam-fi 1.2 TFile *f0_c = 0;
86 mocchiut 1.1 TTree *physicsTree = 0;
87     TBranch *b_registry = 0;
88     TBranch *b_trk = 0;
89     TBranch *b_header = 0;
90     RegistryEvent *reg = 0;
91     EventHeader *header = 0;
92     PscuHeader *pscu = 0;
93     TrkLevel0 *l0_event = 0;
94     // RunInfo
95     ItoRunInfo *runinfo = 0;
96     TArrayL *runlist = 0;
97     ULong64_t from_run;
98     ULong64_t to_run;
99     Bool_t reprocessing = false;
100     //LEVEL1 output - ROOT
101 pam-fi 1.6 TFile *f1;
102     TTree *t_level1 = 0;
103     TrkLevel1 *l1_event = new TrkLevel1();
104     // TrkLevel2Ref *l2_ref = new TrkLevel2Ref();
105 mocchiut 1.1 //LEVEL2 output - ROOT
106     TTree *t_level2 = 0;
107     TTree *t_clone = 0;
108     TrkLevel2 *l2_event = new TrkLevel2();
109     TrkLevel2 *l2_clone = new TrkLevel2();
110     // -----------------------
111     // -----------------------
112     // Handle input parameters
113     // (data)
114     //
115     // - open/create output files, determining the environment root/hbook from the estension
116     // - create the level1/level2 tree+branch/nt-uple
117     // -----------------------
118     // -----------------------
119     TrkProcess *p = new TrkProcess(run,f2);
120     p->HandleCustomPar(ncustom,vcustom);
121 pam-fi 1.2 if( p->DebugMode() )p->Dump();
122 mocchiut 1.1
123     // ===================================
124     // Open/Create level1 output file
125     // (NB output files other than LEVEL2 are
126     // created in a separate folder)
127     // ===================================
128     if(p->get1){
129 pam-fi 1.6 if(p->ifroot1){
130     // if(p->VerboseMode())
131     cout << endl << "Requested LEVEL1 output" << endl;
132     // Create NEW Tracker tree
133     t_level1 = new TTree("Clusters","PAMELA tracker level1 data ");
134     t_level1->Branch("TrkLevel1","TrkLevel1",&l1_event);
135     t_level1->SetDirectory(0);
136     }else{
137     throw -299;
138     };
139 mocchiut 1.1 };
140     // ===================================
141     // Open/Create level2 output file
142     // ===================================
143     if(p->get2){
144     // if(p->ifroot2){ // root
145 pam-fi 1.6 //-------------------------------------------
146     // read RunInfo
147     //-------------------------------------------
148     if(!(p->standalone)){
149     // Open "RunInfo" tree and get list of run
150     runinfo = new ItoRunInfo(f2);
151     char *trkversion = TrkInfo(false);
152     int runinfo_error = runinfo->Update(run,"TRK",trkversion);
153     if( runinfo_error ) throw runinfo_error;
154     runlist = runinfo->GetRunList();
155     reprocessing = runinfo->IsReprocessing();//????
156     if(p->VerboseMode()){
157     cout << "#events "<< runinfo->GetFileEntries() << endl;// #events in the file
158     cout << "#runs "<< runinfo->GetRunEntries() << endl;// #runs in the file
159     cout << "1-st run "<< runlist->At(0) << endl;
160     cout << "last run "<< runlist->At(runinfo->GetRunEntries()-1) << endl;
161     cout << "1-st event "<< runinfo->GetFirstEntry() << endl;// first event of our run
162     cout << "last event+1 "<< runinfo->GetLastEntry() << endl;// first event AFTER the last event of our run
163     cout << "reprocessing "<< runinfo->IsReprocessing() << endl << endl;
164     };
165     };
166     //-------------------------------------------
167     //
168     //-------------------------------------------
169     // Take (if present) the old Tracker tree
170     t_clone = (TTree*)f2->Get("Tracker");
171     if( t_clone != NULL ) t_clone->SetBranchAddress("TrkLevel2",&l2_clone);
172     // Create NEW Tracker tree
173     t_level2 = new TTree("Tracker","PAMELA tracker level2 data ");
174     t_level2->Branch("TrkLevel2","TrkLevel2",&l2_event);
175     // per ora lo commento... ma potrebbe essere una soluzione
176     /* if(p->get1){
177     t_level2->Branch("TrkLevel1","TrkLevel1",&l1_event);
178     t_level2->BranchRef();
179     };*/
180 mocchiut 1.1 // }else{ // hbook
181     // throw -299;
182     // }
183     };
184    
185    
186     // -------------------------------------------
187     // define runs to be processed/reprocessed
188     // -------------------------------------------
189     if(run == 0){
190     // reprocessing ALL runs
191     if(p->standalone)throw -298; // reprocessing not implemented
192     from_run = runlist->At(0);
193     to_run = runlist->At(runinfo->GetRunEntries()-1);
194     }else{
195     // processing/reprocessing ONE single run
196     from_run = run;
197     to_run = run;
198     };
199    
200     //
201     // init "expiration date" of calibrations and parameters
202     //
203     ULong64_t pedsig_time =0ULL;
204     ULong64_t B_time =0ULL;
205     ULong64_t mip_time =0ULL;
206     ULong64_t charge_time =0ULL;
207     ULong64_t eta_time =0ULL;
208     ULong64_t align_time =0ULL;
209     ULong64_t mask_time =0ULL;
210     //
211     // init event counter
212     //
213     Int_t ev_count =0;
214     //
215     // create query-results objects
216     //
217     GL_RUN q1 = GL_RUN();
218     GL_TRK_CALIB q2 = GL_TRK_CALIB();
219     GL_ROOT q3 = GL_ROOT();
220     GL_PARAM q4 = GL_PARAM();
221    
222     // ------------------------------------------------------------
223     // if reprocessing one run, copy all the events BEFORE the run
224     // ------------------------------------------------------------
225     if( !(p->standalone) ){
226 pam-fi 1.6 for(UInt_t i=0; i<runinfo->GetFirstEntry(); i++){
227     t_clone->GetEntry(i);
228     *l2_event = *l2_clone;
229     t_level2->Fill();
230     l2_event->Clear();
231     // COPY COPY COPY
232     };
233 mocchiut 1.1 };
234     // ------------------------------------------------------------
235     // ------------------------------------------------------------
236     // START LOOP OVER RUNS TO PROCESSED/REPROCESSED
237     // ------------------------------------------------------------
238     // ------------------------------------------------------------
239     for(ULong64_t idRun=from_run; idRun <= to_run; idRun++){
240    
241 pam-fi 1.2 if(p->VerboseMode()) cout << endl<<" ========================= Run: "<< idRun << endl;
242 mocchiut 1.1 ULong64_t runheadtime = 0ULL;
243     ULong64_t runtrailtime = 0ULL;
244     UInt_t evfrom = 0;
245     UInt_t evto = 0;
246     Int_t id_reg_run =-1;
247     Int_t trk_calib_used = 0;
248    
249     if(p->standalone){
250     // ==============================
251     // first query: retrieve run info
252     // ==============================
253     if (q1.Query_GL_RUN(idRun,dbc) )throw -50;
254     id_reg_run = q1.ID_REG_RUN;
255     runheadtime = q1.RUNHEADER_TIME;
256     runtrailtime = q1.RUNTRAILER_TIME;
257     evfrom = q1.EV_REG_PHYS_FROM;
258     evto = q1.EV_REG_PHYS_TO;
259     trk_calib_used = q1.TRK_CALIB_USED;
260     }else{
261     // ==============================
262     // get run info from RunInfo tree
263     // ==============================
264 pam-fi 1.3 int runinfo_error = runinfo->GetRunInfo(idRun);
265     if( runinfo_error ) throw runinfo_error;
266 mocchiut 1.1 id_reg_run = runinfo->ID_REG_RUN;
267     runheadtime = runinfo->RUNHEADER_TIME;
268     runtrailtime = runinfo->RUNTRAILER_TIME;
269     evfrom = runinfo->EV_REG_PHYS_FROM;
270     evto = runinfo->EV_REG_PHYS_TO;
271     trk_calib_used = runinfo->TRK_CALIB_USED;
272     };
273     //
274 pam-fi 1.2 if(p->VerboseMode()){
275 mocchiut 1.1 cout << "ROOT file ID "<< id_reg_run << endl;
276     cout << "RunHeader time "<< runheadtime << endl;
277     cout << "RunTrailer time "<< runtrailtime << endl;
278     cout << " from event "<< evfrom << endl;
279     cout << " to event "<< evto << endl;
280     cout << "trk-calibration used "<< trk_calib_used << endl;
281     };
282     // ========================================================
283     // second query: search the LEVEL0 file that contains idRun
284     // ========================================================
285     TString lastfilename = filename;
286     if( q3.Query_GL_ROOT(id_reg_run,dbc) )throw -51;
287     filename = q3.PATH + q3.NAME;
288     // ========================================================
289     // Open the input LEVEL0 data file
290     // ========================================================
291     if(filename.CompareTo(lastfilename)){
292     if(!lastfilename.IsNull())f0->Close();
293     //if( debug ) cout << "Opening LEVEL0 file: "<< filename << endl;
294 pam-fi 1.2 if(p->VerboseMode()) cout << "Opening LEVEL0 file: "<< filename << endl;
295     FileStat_t t;
296     if( gSystem->GetPathInfo(filename.Data(),t) )throw -6;
297 mocchiut 1.1 f0 = new TFile(filename);
298     if ( !f0 ) throw -6;
299     physicsTree = (TTree*)f0->Get("Physics"); if(!physicsTree) throw -7;
300     b_header = physicsTree ->GetBranch("Header"); if(!b_header ) throw -8;
301     b_registry = physicsTree ->GetBranch("Registry"); if(!b_registry ) throw -9;
302     b_trk = physicsTree ->GetBranch("Tracker"); if(!b_trk ) throw -203;
303     physicsTree->SetBranchAddress("Tracker" ,&l0_event);
304     physicsTree->SetBranchAddress("Registry",&reg);
305     physicsTree->SetBranchAddress("Header",&header);
306    
307     nevents = physicsTree->GetEntries();
308     if ( nevents < 1 ) throw -11;
309 pam-fi 1.3 if ( nevents < (evto+1)) throw -12;
310 mocchiut 1.1 };
311    
312     // =============================================================
313 pam-fi 1.2 // retrieve calibration file needed to reduce data
314 mocchiut 1.1 // =============================================================
315 pam-fi 1.2 // if run OBT is > last calibration "expiration date"
316     // - search for new calibration packet
317     // - load calibration parameters (full + truncated)
318 pam-fi 1.5 if(p->VerboseMode() )cout << "Full pedestals for cluster finding:";
319 pam-fi 1.2 if(runheadtime > pedsig_time){
320     if( q2.Query_GL_TRK_CALIB(runheadtime,dbc) )throw -53;
321     pedsig_time = q2.TO_TIME;
322     if(q2.EV_REG_CALIBTRK1 != q2.EV_REG_CALIBTRK2)
323     printf("WARNING!! ---> EV_REG_CALIBTRK1=%d it's different from EV_REG_CALIBTRK2=%d \n\n",q2.EV_REG_CALIBTRK1,q2.EV_REG_CALIBTRK2);
324     if( q3.Query_GL_ROOT(q2.ID_REG_CALIBTRK,dbc) )throw -51;
325 pam-fi 1.5 trkcalibfile = q3.PATH + q3.NAME;
326     if(p->VerboseMode()) cout << " >> Loading from LEVEL0 file: "<< trkcalibfile << endl;
327 pam-fi 1.2 if(trkcalibfile.CompareTo(filename)){
328 pam-fi 1.5 FileStat_t t;
329     if( gSystem->GetPathInfo(trkcalibfile.Data(),t) )throw -6;
330     f0_c=new TFile(trkcalibfile);
331     if ( !f0_c ) throw -6;
332 pam-fi 1.2 }else{
333 pam-fi 1.5 f0_c = f0;
334 pam-fi 1.2 };
335     if(p->VerboseMode()) {
336 pam-fi 1.5 cout << " calibration entries "<< q2.EV_REG_CALIBTRK1 << " " << q2.EV_REG_CALIBTRK2;
337 pam-fi 1.2 cout << " (from time "<< q2.FROM_TIME <<" to time "<< q2.TO_TIME <<")"<<endl;
338     };
339     pedsigbad_.FillACalibFrom(f0_c,q2.EV_REG_CALIBTRK1,q2.EV_REG_CALIBTRK2);
340 pam-fi 1.5 }else{
341     if(p->VerboseMode() )cout<< " already loaded. "<< endl;
342 pam-fi 1.2 };
343     if( p->DebugMode() ) for(int i=0; i<12; i++) cout << " DSP "<< i << " "<< pedsigbad_.pedestal[64][12][i] << endl;
344     // =============================================================
345     // retrieve calibration file needed to uncompress data
346     // =============================================================
347     // if the run was compressed using default calib
348     // load truncated pedestals from default
349     // otherwise reload them from on-line calibration
350     if(p->VerboseMode() )cout << "Truncated pedestals for uncompression:";
351 mocchiut 1.1 if(trk_calib_used==104 && last_trk_calib_used!=104){
352    
353 pam-fi 1.5 if(p->VerboseMode() )cout << " >> Loading default calibration " << endl;
354 pam-fi 1.2 if (q4.Query_GL_PARAM(runheadtime,"default calibration",dbc) )throw -52;
355     pedsigbad_.FillTCalibFrom(q4.PATH+q4.NAME);
356    
357     }else if( (trk_calib_used==1 | trk_calib_used==2) && last_trk_calib_used==104 ){
358    
359     if(p->VerboseMode()) cout << ">> Re-loading on-line calibration " << endl;
360     pedsigbad_.FillTCalibFrom(f0_c,q2.EV_REG_CALIBTRK1,q2.EV_REG_CALIBTRK2);
361     }else{
362     if(p->VerboseMode()) cout << ">> Using on-line calibration " << endl;
363     };
364     if( p->DebugMode() ) for(int i=0; i<12; i++) cout << " DSP "<< i << " "<< pedsigbad_.pedestal_t[64][12][i] << endl;
365     last_trk_calib_used = trk_calib_used;
366    
367 mocchiut 1.1 // =============================================================
368     // retrieve information about parameters to process LEVEL2
369     // =============================================================
370     //
371     // magnetic field
372     //
373 pam-fi 1.2 if(p->VerboseMode()) cout << "Loading parameter files : "<< endl;
374    
375 mocchiut 1.1 if(runheadtime > B_time){
376     if( q4.Query_GL_PARAM(runheadtime,"field",dbc) )throw -52;
377     B_time = q4.TO_TIME;
378     l2_event->LoadField(q4.PATH+q4.NAME);
379 pam-fi 1.2 // l2_event->LoadField(q4.PATH+q4.NAME);
380 mocchiut 1.1 if(path_.error) throw -215;
381     };
382     //
383     // mip conversion parameters
384     //
385     if(runheadtime > mip_time){
386     if( q4.Query_GL_PARAM(runheadtime,"trk mip",dbc) )throw -52;
387     mip_time = q4.TO_TIME;
388     path_.FillWith(q4.PATH+q4.NAME);
389     readmipparam_();
390     if(path_.error) throw -212;
391     };
392     //
393     // charge correlation parameters
394     //
395     if(runheadtime > charge_time){
396     if( q4.Query_GL_PARAM(runheadtime,"trk charge",dbc) )throw -52;
397     charge_time = q4.TO_TIME;
398     path_.FillWith(q4.PATH+q4.NAME);
399     readchargeparam_();
400     if(path_.error) throw -213;
401     };
402     //
403     // eta p.f.a. parameters
404     //
405     if(runheadtime > eta_time){
406     if( q4.Query_GL_PARAM(runheadtime,"trk pfa",dbc) )throw -52;
407     eta_time = q4.TO_TIME;
408     path_.FillWith(q4.PATH+q4.NAME);
409     readetaparam_();
410     if(path_.error) throw -214;
411     };
412     //
413     // alignment parameters
414     //
415     if(runheadtime > align_time){
416     if( q4.Query_GL_PARAM(runheadtime,"trk alignment",dbc) )throw -52;
417     align_time = q4.TO_TIME;
418     path_.FillWith(q4.PATH+q4.NAME);
419     readalignparam_();
420     if(path_.error) throw -211;
421     };
422     //
423     // viking mask
424     //
425     if(runheadtime > mask_time){
426     if( q4.Query_GL_PARAM(runheadtime,"trk mask",dbc) )throw -52;
427     mask_time = q4.TO_TIME;
428     path_.FillWith(q4.PATH+q4.NAME);
429     readvkmask_();
430     if(path_.error) throw -210;
431     };
432    
433     //
434     // ====================================================
435     // start looping on events cointained in the data file
436     // ====================================================
437     for (UInt_t re = evfrom; re <= evto; re++){
438    
439     ev_count++;
440     // if ( p->DEBUG && re%100 == 0 && re > 0 ) printf(" %iK \n",re/1000);
441 pam-fi 1.2 if ( p->DebugMode() && re%100 == 0 && re > 0 ) cout << ".";
442 mocchiut 1.1
443     b_registry->GetEntry(re);
444     b_trk->GetEntry(reg->event);
445     b_header->GetEntry(reg->event);
446     pscu = header->GetPscuHeader();
447     // pscu->Print();
448     // cout <<endl<< " *************** Pkt ID "<< (UINT32)pscu->Counter << endl;
449    
450 pam-fi 1.2 p->ProcessEvent(l0_event);
451 pam-fi 1.5
452 mocchiut 1.1 // ----------------
453     // LEVEL1 output
454     // ----------------
455 pam-fi 1.6 if(p->get1){
456     if(p->ifroot1){ // root
457     l1_event->Clear();
458     l1_event->SetFromLevel1Struct(&level1event_);
459     t_level1->Fill();
460     }else{ // hbook
461     };
462 mocchiut 1.1 };
463     // ----------------
464     // LEVEL2 output
465     // ----------------
466     if(p->get2){
467     // if(p->ifroot2){ // root
468 pam-fi 1.6 l2_event->Clear();
469     if(p->get1) l2_event->SetFromLevel2Struct(&level2event_,l1_event);//set references to level1
470     else l2_event->SetFromLevel2Struct(&level2event_);
471     t_level2->Fill();
472 mocchiut 1.1 // }else{ // hbook
473     // }
474     };
475    
476     }; // end loop on events
477 pam-fi 1.2 if(p->VerboseMode())cout << " >>> processed "<< ev_count <<" events"<< endl;
478 mocchiut 1.1 }; // end loop on runs
479     // ------------------------------------------------------------
480     // if reprocessing one run, copy all the events AFTER the run
481     // ------------------------------------------------------------
482     if( !(p->standalone) ){
483     for(UInt_t i=runinfo->GetLastEntry()+1; i<runinfo->GetFileEntries(); i++){
484     t_clone->GetEntry(i);
485     *l2_event = *l2_clone;
486     t_level2->Fill();
487     l2_event->Clear();
488     // COPY COPY COPY
489     };
490     };
491     // ---------------
492     // close the files
493     // ---------------
494 pam-fi 1.6 if(p->get1){
495     if(p->ifroot1){ // root
496     cout << "Creating file: "<<p->file1<< endl;
497     f1 = new TFile(p->file1,"RECREATE");
498     t_level1->Write();
499     // gDirectory->pwd();
500     // gDirectory->ls();
501     if( f1->IsOpen() ) f1->Close();
502     if( t_level1 )t_level1->Delete(); //delete new tree from memory
503     // f2->Write("Clusters");//non funziona
504    
505     }else{ // hbook
506     }
507     cout << "Done."<<endl;
508 mocchiut 1.1 };
509     if(p->get2){
510     // if(p->ifroot2){ // root
511 pam-fi 1.6 if( t_clone )t_clone->Delete("all");//delete old tree from file
512     if( !(p->standalone) )runinfo->Close();
513     // gDirectory->pwd();
514     // gDirectory->ls();
515     f2->Write("Tracker");
516     if( t_level2 )t_level2->Delete(); //delete new tree from memory
517 mocchiut 1.1
518     // }else{ // hbook
519     // }
520 pam-fi 1.6 };
521    
522     if( f0->IsOpen() ) f0->Close();
523 pam-fi 1.2 if( f0_c->IsOpen() )f0_c->Close();
524 pam-fi 1.6
525 pam-fi 1.2
526     return(p->ostatus);
527 mocchiut 1.1 }
528    

  ViewVC Help
Powered by ViewVC 1.1.23