/[PAMELA software]/PamelaLevel2/src/PamLevel2.cpp
ViewVC logotype

Contents of /PamelaLevel2/src/PamLevel2.cpp

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1.39 - (show annotations) (download)
Mon May 14 11:57:36 2007 UTC (17 years, 9 months ago) by pam-fi
Branch: MAIN
Changes since 1.38: +60 -15 lines
fixed bug during recovery of level0 event

1 #include <PamLevel2.h>
2 //--------------------------------------
3 //
4 //
5 //--------------------------------------
6 /**
7 * Default constructor
8 */
9 PamTrack::PamTrack(){
10 trk_track = 0;
11 calo_track = 0;
12 tof_track = 0;
13 candeleteobj = 0;
14 };
15 //--------------------------------------
16 //
17 //
18 //--------------------------------------
19 /**
20 * Constructor
21 */
22 PamTrack::PamTrack(TrkTrack* t, CaloTrkVar* c, ToFTrkVar* o){
23
24 trk_track = 0;
25 calo_track = 0;
26 tof_track = 0;
27 // if(t)trk_track = new TrkTrack(*t);
28 // if(c)calo_track = new CaloTrkVar(*c);
29 // if(o)tof_track = new ToFTrkVar(*o);
30 if(t)trk_track = t;
31 if(c)calo_track = c;
32 if(o)tof_track = o;
33 candeleteobj = 0;
34
35 };
36 PamTrack::PamTrack(const PamTrack& track){
37
38 TrkTrack *t = track.trk_track;
39 CaloTrkVar *c = track.calo_track;
40 ToFTrkVar *o = track.tof_track;
41
42 trk_track = 0;
43 calo_track = 0;
44 tof_track = 0;
45 if(t)trk_track = new TrkTrack(*t);
46 if(c)calo_track = new CaloTrkVar(*c);
47 if(o)tof_track = new ToFTrkVar(*o);
48 candeleteobj = 1;
49
50 }
51 void PamTrack::Clear(){
52
53 // cout << "PamTrack::Clear() "<<candeleteobj<<endl;
54 if(candeleteobj){
55 if(trk_track) trk_track->TrkTrack::Clear();
56 if(calo_track) calo_track->CaloTrkVar::Clear();//???
57 if(tof_track) tof_track->ToFTrkVar::Clear();//???
58 }else{
59 trk_track = 0;
60 calo_track = 0;
61 tof_track = 0;
62 }
63 }
64 void PamTrack::Delete(){
65 // cout << "PamTrack::Delete() "<<candeleteobj<<endl;
66 if(candeleteobj){
67 if(trk_track) {
68 trk_track->TrkTrack::Clear();
69 delete trk_track;
70 }
71 if(calo_track){
72 calo_track->CaloTrkVar::Clear();//???
73 delete calo_track;
74 }
75 if(tof_track){
76 tof_track->ToFTrkVar::Clear();//???
77 delete tof_track;
78 }
79 }else{
80 Clear();
81 }
82 }
83 //--------------------------------------
84 //
85 //
86 //--------------------------------------
87 /**
88 * Default Constructor
89 */
90 PamLevel2::PamLevel2(){
91 Initialize();
92 };
93
94 /**
95 * Constructor with given dir, list and detectors
96 */
97 PamLevel2::PamLevel2(TString ddir,TString list,TString detlist){
98 Initialize();
99 GetPamTree(GetListOfLevel2Files(ddir,list),detlist);
100 GetRunTree(GetListOfLevel2Files(ddir,list));
101 };
102
103 /**
104 * Constructor with given dir and list
105 */
106 PamLevel2::PamLevel2(TString ddir,TString list){
107 Initialize();
108 GetPamTree(GetListOfLevel2Files(ddir,list),"");
109 GetRunTree(GetListOfLevel2Files(ddir,list));
110 };
111
112
113 void PamLevel2::Initialize(){
114
115 // trk2_obj = TrkLevel2::GetTrkLevel2();
116 // trk1_obj = TrkLevel1::GetTrkLevel1();
117 // trkh_obj = TrkHough::GetTrkHough();
118 // calo1_obj = CaloLevel1::GetCaloLevel1();
119 // calo2_obj = CaloLevel2::GetCaloLevel2();
120 // tof_obj = ToFLevel2::GetToFLevel2();
121 // trig_obj = TrigLevel2::GetTrigLevel2();
122 // s4_obj = S4Level2::GetS4Level2();
123 // nd_obj = NDLevel2::GetNDLevel2();
124 // ac_obj = AcLevel2::GetAcLevel2();
125 // orb_obj = OrbitalInfo::GetOrbitalInfo();
126
127 h0_obj = 0;
128 trk0_obj = 0;
129 calo0_obj = 0;
130
131 trk2_obj = 0;
132 trk1_obj = 0;
133 trkh_obj = 0;
134 calo1_obj = 0;
135 calo2_obj = 0;
136 tof_obj = 0;
137 trig_obj = 0;
138 s4_obj = 0;
139 nd_obj = 0;
140 ac_obj = 0;
141 orb_obj = 0;
142
143 run_obj = 0;//new GL_RUN();
144 soft_obj = 0;// Emiliano
145 irun = -1LL;
146 runfirstentry = 0LL;
147 runlastentry = 0LL;
148
149 l0_file = NULL;
150 l0_tree = NULL;
151 iroot = -1;
152 dbc = 0;
153
154 prevshift = 0;
155
156 run_tree = NULL;
157 run_tree_clone = NULL;
158 sel_tree = NULL;
159 sel_tree_clone = NULL;
160
161 irunentry = -1LL;
162 pam_tree = NULL;
163 for(Int_t i=0; i<NCLONES; i++ )pam_tree_clone[i]=NULL;
164
165 host = "mysql://localhost/pamelaprod";
166 user = "anonymous";
167 psw = "";
168 const char *pamdbhost=gSystem->Getenv("PAM_DBHOST");
169 const char *pamdbuser=gSystem->Getenv("PAM_DBUSER");
170 const char *pamdbpsw=gSystem->Getenv("PAM_DBPSW");
171 if ( !pamdbhost ) pamdbhost = "";
172 if ( !pamdbuser ) pamdbuser = "";
173 if ( !pamdbpsw ) pamdbpsw = "";
174 if ( strcmp(pamdbhost,"") ) host = pamdbhost;
175 if ( strcmp(pamdbuser,"") ) user = pamdbuser;
176 if ( strcmp(pamdbpsw,"") ) psw = pamdbpsw;
177
178
179
180 // sorted_tracks = 0;//new TRefArray();
181
182 CAL0 = false;
183 CAL1 = true;
184 CAL2 = true;
185 TRK2 = true;
186 TRK1 = false;
187 TRKh = false;
188 TRKh = false;
189 TRG = true;
190 TOF = true;
191 TOF0 = false;
192 S4 = true;
193 ND = true;
194 AC = true;
195 ORB = true;
196
197 RUN = true;
198
199 SELLI = -1;
200
201 tsorted=0;
202 timage=0;
203
204 howtosort = "+CAL+TOF";
205 //howtosort = "+TOF";
206 sortthr = 100.;
207
208 };
209 /**
210 * Delete the event (act as Dtor)
211 */
212 void PamLevel2::Delete(){
213
214 if(run_obj)delete run_obj;
215 if(soft_obj)delete soft_obj; //Emiliano
216
217 // cout << "void PamLevel2::Clear()"<<endl;
218 if(h0_obj) delete h0_obj;
219 if(trk0_obj) delete trk0_obj;
220 if(calo0_obj) delete calo0_obj;
221 if(trk1_obj) delete trk1_obj;
222 if(trk2_obj) delete trk2_obj;
223 if(trkh_obj) delete trkh_obj;
224 if(calo1_obj)delete calo1_obj;
225 if(calo2_obj)delete calo2_obj;
226 if(tof_obj) delete tof_obj;
227 if(trig_obj) delete trig_obj;
228 if(s4_obj) delete s4_obj;
229 if(nd_obj) delete nd_obj;
230 if(ac_obj) delete ac_obj;
231 if(orb_obj) delete orb_obj;
232
233 if(tsorted){
234 tsorted->Delete();
235 delete tsorted;
236 }
237 if(timage){
238 timage->Delete();
239 delete timage;
240 }
241
242 if(dbc){
243 dbc->Close();
244 delete dbc;
245 }
246
247 if(l0_file)l0_file->Close();
248 // if(pam_tree)pam_tree->Delete();;
249
250 if ( pam_tree ){
251 //
252 // we have first to find which chains we have to delete, then delete the main chain and only after delete the friend chains. Any other order causes a segfault...
253 //
254 TList *temp = pam_tree->GetListOfFriends();
255 TList *contents = new TList; // create chain friend list
256 contents->SetOwner();
257 TIter next(temp);
258 TChain *questo = 0;
259 while ( (questo = (TChain*) next()) ){
260 TString name = questo->GetName();
261 contents->Add((TChain*)gROOT->FindObject(name.Data()));// add object to the list
262 };
263 //
264 // deleting the main chain
265 //
266 pam_tree->Delete();
267 //
268 // deleting the friends...
269 //
270 TIter next2(contents);
271 TChain *questa = 0;
272 while ( questa = (TChain*)next2() ){
273 TString name = questa->GetName();
274 questa->Delete();
275 questa=NULL;
276 };
277 //
278 };
279 pam_tree = NULL;
280
281 if(run_tree)run_tree->Delete();;
282 if(sel_tree)sel_tree->Delete();;
283 for(Int_t i=0; i<NCLONES; i++ )if(pam_tree_clone[i])pam_tree_clone[i]->Delete();;
284 if(run_tree_clone)run_tree_clone->Delete();;
285 if(sel_tree_clone)sel_tree_clone->Delete();;
286
287
288 };
289
290 /**
291 * Clear the event (NB! does not deallocate objects)
292 */
293 void PamLevel2::Clear(){
294
295 // cout << "void PamLevel2::Clear()"<<endl;
296
297 //
298 // This method is called once for every entry but RunInfo and SoftInfo do not change until the next run so we cannot clear them here unless we don't
299 // want to load them for each event even if they are the same...
300 //
301 // if(run_obj)delete run_obj;
302 // if(run_obj) run_obj->Clear(); // Emiliano: Do not deallocate run_obj here, it will give segmentation fault! call clear instead
303 // if(soft_obj) soft_obj->Clear();
304
305 if(h0_obj) h0_obj->Clear();
306 // if(trk0_obj) trk0_obj->Clear();
307 if(calo0_obj) calo0_obj->Clear();
308 if(trk1_obj) trk1_obj->Clear();
309 if(trk2_obj) trk2_obj->Clear();
310 if(trkh_obj) trkh_obj->Clear();
311 if(calo1_obj)calo1_obj->Clear();
312 if(calo2_obj)calo2_obj->Clear();
313 if(tof_obj) tof_obj->Clear();
314 if(trig_obj) trig_obj->Clear();
315 if(s4_obj) s4_obj->Clear();
316 if(nd_obj) nd_obj->Clear();
317 if(ac_obj) ac_obj->Clear();
318 if(orb_obj) orb_obj->Clear();
319
320 // if(sorted_tracks)sorted_tracks->Clear();
321 // sorted_tracks.Clear();
322
323 if(tsorted){
324 tsorted->Delete();
325 }
326 if(timage){
327 timage->Delete();
328 }
329 };
330
331 void PamLevel2::Reset(){
332 //
333 // First of all clear everything
334 //
335 Clear();
336 //
337 // close and reset chains and pointers
338 //
339 if ( pam_tree ){
340 //
341 // we have first to find which chains we have to delete, then delete the main chain and only after delete the friend chains. Any other order causes a segfault...
342 //
343 TList *temp = pam_tree->GetListOfFriends();
344 TList *contents = new TList; // create chain friend list
345 contents->SetOwner();
346 TIter next(temp);
347 TChain *questo = 0;
348 while ( (questo = (TChain*) next()) ){
349 TString name = questo->GetName();
350 contents->Add((TChain*)gROOT->FindObject(name.Data()));// add object to the list
351 };
352 //
353 // deleting the main chain
354 //
355 pam_tree->Delete();
356 //
357 // deleting the friends...
358 //
359 TIter next2(contents);
360 TChain *questa = 0;
361 while ( questa = (TChain*)next2() ){
362 TString name = questa->GetName();
363 questa->Delete();
364 questa=NULL;
365 };
366 //
367 };
368 pam_tree = NULL;
369 //
370 if(run_tree)run_tree->Delete();;
371 run_tree = NULL;
372 if(sel_tree)sel_tree->Delete();;
373 sel_tree = NULL;
374 //
375 // Close file
376 //
377 if(l0_file)l0_file->Close("R");
378 l0_file = NULL;
379 //
380 h0_obj = 0;
381 trk0_obj = 0;
382 calo0_obj = 0;
383 //
384 trk2_obj = 0;
385 trk1_obj = 0;
386 trkh_obj = 0;
387 calo1_obj = 0;
388 calo2_obj = 0;
389 tof_obj = 0;
390 trig_obj = 0;
391 s4_obj = 0;
392 nd_obj = 0;
393 ac_obj = 0;
394 orb_obj = 0;
395 //
396 // Reset run pointers
397 //
398 run_obj = 0;//new GL_RUN();
399 soft_obj = 0;// Emiliano
400 irun = -1;
401 runfirstentry = 0ULL;
402 runlastentry = 0ULL;
403 //
404 };
405
406 Bool_t PamLevel2::IsGood(){
407 Bool_t goodev=true;
408 // if(trk2_obj && trk2_obj->UnpackError() != 0 ) goodev = false;
409 if(calo2_obj && calo2_obj->good != 1) goodev = false;
410 if(tof_obj && tof_obj->unpackError != 0) goodev = false;
411 if(trig_obj && trig_obj->unpackError != 0) goodev = false;
412 if(s4_obj && s4_obj->unpackError != 0) goodev = false;
413 if(nd_obj && nd_obj->unpackError != 0) goodev = false;
414 if(ac_obj && ac_obj->unpackError != 255) goodev = false;
415 // if(orb_obj)
416 return goodev;
417 };
418
419 //--------------------------------------
420 //
421 //
422 //--------------------------------------
423 void *PamLevel2::GetPointerTo(const char* c ){
424
425 TString objname = c;
426
427 if(!objname.CompareTo("TrkLevel1")) {
428 if(!trk1_obj){
429 trk1_obj = new TrkLevel1();
430 trk1_obj->Set();
431 }
432 return &trk1_obj;
433 };
434 if(!objname.CompareTo("TrkLevel2")) {
435 if(!trk2_obj){
436 trk2_obj = new TrkLevel2();
437 trk2_obj->Set();
438 }
439 return &trk2_obj;
440 };
441 if(!objname.CompareTo("TrkHough")) {
442 if(!trkh_obj) trkh_obj = new TrkHough();
443 return &trkh_obj;
444 };
445 if(!objname.CompareTo("CaloLevel1")) {
446 if(!calo1_obj) calo1_obj = new CaloLevel1();
447 return &calo1_obj;
448 };
449 if(!objname.CompareTo("CaloLevel2")) {
450 if(!calo2_obj){
451 calo2_obj = new CaloLevel2();
452 calo2_obj->Set();
453 };
454 return &calo2_obj;
455 };
456 if(!objname.CompareTo("ToFLevel2")) {
457 if(!tof_obj){
458 tof_obj = new ToFLevel2();
459 tof_obj->Set();
460 }
461 return &tof_obj;
462 };
463 if(!objname.CompareTo("TrigLevel2")) {
464 if(!trig_obj) trig_obj = new TrigLevel2();
465 return &trig_obj;
466 };
467 if(!objname.CompareTo("S4Level2")) {
468 if(!s4_obj) s4_obj = new S4Level2();
469 return &s4_obj;
470 };
471 if(!objname.CompareTo("NDLevel2")) {
472 if(!nd_obj) nd_obj = new NDLevel2();
473 return &nd_obj;
474 };
475 if(!objname.CompareTo("AcLevel2")) {
476 if(!ac_obj) ac_obj = new AcLevel2();
477 return &ac_obj;
478 };
479 if(!objname.CompareTo("OrbitalInfo")){
480 if(!orb_obj) orb_obj = new OrbitalInfo();
481 return &orb_obj;
482 };
483
484 if(!objname.CompareTo("RunInfo"))return &run_obj;
485
486 if(!objname.CompareTo("SoftInfo"))return &soft_obj; // Emiliano
487
488 return NULL;
489 };
490 //--------------------------------------
491 //
492 //
493 //--------------------------------------
494 /**
495 * Retrieves the calorimeter track matching the seqno-th tracker stored track.
496 * (If seqno = -1 retrieves the self-trigger calorimeter track)
497 */
498 CaloTrkVar *PamLevel2::GetCaloStoredTrack(int seqno){
499
500 if( !calo2_obj )return 0;
501
502 if( calo2_obj->CaloLevel2::ntrk()==0 ){
503 cout << "PamLevel2::GetCaloStoredTrack(int) : requested tracker SeqNo "<< seqno <<" but no Calorimeter tracks are stored"<<endl;
504 return NULL;
505 };
506
507 CaloTrkVar *c = 0;
508 Int_t it_calo=0;
509
510 do{
511 c = calo2_obj->CaloLevel2::GetCaloTrkVar(it_calo);
512 it_calo++;
513 } while( c && seqno != c->trkseqno && it_calo < calo2_obj->CaloLevel2::ntrk());
514
515 if(!c || seqno != c->trkseqno){
516 c = 0;
517 if(seqno!=-1)cout << "PamLevel2::GetCaloStoredTrack(int) : requested tracker SeqNo "<< seqno <<" does not match Calorimeter stored tracks"<<endl;
518 };
519 return c;
520
521 };
522 //--------------------------------------
523 //
524 //
525 //--------------------------------------
526 /**
527 * Retrieves the ToF track matching the seqno-th tracker stored track.
528 * (If seqno = -1 retrieves the tracker-independent tof track)
529 */
530 ToFTrkVar *PamLevel2::GetToFStoredTrack(int seqno){
531
532 if( !tof_obj )return 0;
533
534 if( tof_obj->ToFLevel2::ntrk()==0 ){
535 cout << "PamLevel2::GetToFStoredTrack(int) : requested tracker SeqNo "<< seqno <<" but no ToF tracks are stored"<<endl;
536 return NULL;
537 };
538
539 ToFTrkVar *c = 0;
540 Int_t it_tof=0;
541
542 do{
543 c = tof_obj->ToFLevel2::GetToFTrkVar(it_tof);
544 it_tof++;
545 } while( c && seqno != c->trkseqno && it_tof < tof_obj->ToFLevel2::ntrk());
546
547 if(!c || seqno != c->trkseqno){
548 c = 0;
549 if(seqno!=-1)cout << "PamLevel2::GetToFStoredTrack(int) : requested tracker SeqNo "<< seqno <<" does not match ToF stored tracks"<<endl;
550 };
551 return c;
552
553 };
554
555 //--------------------------------------
556 //
557 //
558 //--------------------------------------
559 /**
560 * Give the pamela track associated to a tracker track, retrieving related calorimeter and tof track information.
561 */
562 PamTrack* PamLevel2::GetPamTrackAlong(TrkTrack* t){
563
564 cout <<"PamLevel2::GetPamTrackAlong(TrkTrack* t) **obsolete** "<<endl;
565 cout <<"(if you use it, remember to delete the PamTrack object)"<<endl;
566
567 CaloTrkVar *c = 0;
568 ToFTrkVar *o = 0;
569
570 if(CAL2) c = GetCaloStoredTrack(t->GetSeqNo());
571 if(TOF) o = GetToFStoredTrack(t->GetSeqNo());
572
573 // if(t && c && o)track = new PamTrack(t,c,o);
574 PamTrack *track = new PamTrack(t,c,o);
575
576 return track;
577
578 };
579 //--------------------------------------
580 //
581 //
582 //--------------------------------------
583 /**
584 * Retrieves the it-th stored track.
585 * It override TrkLevel2::GetTrack(int it).
586 * @param itrk Track number, ranging from 0 to GetNTracks().
587 */
588
589 PamTrack* PamLevel2::GetStoredTrack(Int_t itrk){
590
591 cout <<"PamLevel2::GetStoredTrack(Int_t itrk) **to-be-updated** "<<endl;
592 cout <<"for the moment, better use separately the methods: TrkLevel2::GetStoredTrack(seqno) CaloLevel2::GetCaloTrkVar(Int_t notrack) ToFLevel2::GetToFTrkVar(Int_t notrack)"<<endl;
593 cout <<"(if you use it, remember to delete the PamTrack object)"<<endl;
594 PamTrack *track = 0;
595
596 if( itrk >=0 && itrk < trk2_obj->TrkLevel2::ntrk() ){
597
598 TrkTrack *t = trk2_obj->TrkLevel2::GetStoredTrack(itrk);
599 track = GetPamTrackAlong(t);
600
601
602
603 }else{
604 cout << "PamLevel2::GetStoredTrack(int) : tracker track SeqNo "<< itrk <<" does not exist (GetNTracks() = "<<trk2_obj->TrkLevel2::GetNTracks()<<")"<<endl;
605 };
606
607 return track;
608
609 }
610 //--------------------------------------
611 //
612
613 /**
614 * Sort physical (tracker) tracks and stores them in the TRefArray (of TrkTrack objects) which pointer is PamLevel2::sorted_tracks. Left here as backward compatibility method.
615 **/
616 void PamLevel2::SortTracks(TString how){
617 printf(" WARNING! obsolete, use SortTracks() and SetSortingMethod(TString) instead! \n Setting sorting method to %s \n",how.Data());
618 howtosort = how;
619 SortTracks();
620 };
621
622 //
623 //--------------------------------------
624 /**
625 * Sort physical (tracker) tracks and stores them in the TRefArray (of TrkTrack objects) which pointer is PamLevel2::sorted_tracks.
626 * @param how String to set the sorting cryterium (es: "CAL" or "TRK+CAL+TOF" ecc...).
627 * Sorting cryteria:
628 * TRK: lower chi**2
629 * CAL: lower Y spatial residual on the first calorimeter plane
630 * TOF: bigger numebr of hit PMTs along the track, on S12 S21 S32 (where paddles are along the Y axis).
631 * The default sorting cryterium is "TOF+CAL".
632 *
633 * The total number of physical tracks is always given by GetNTracks() and the it-th physical track can be retrieved by means of the methods GetTrack(int it) and GetTrack(int it, TString how).
634 */
635 void PamLevel2::SortTracks(){
636 TString how = howtosort;
637
638 // cout <<" PamLevel2::SortTracks(TString how) "<<endl;
639 if( !trk2_obj ){
640 cout << "void PamLevel2::SortTracks(): TrkLevel2 not loaded !!!";
641 return;
642 };
643 // cout << "call SortTracs() "<<endl;
644 //Save current Object count
645 Int_t ObjectNumber = TProcessID::GetObjectCount();
646
647 // cout << "ObjectNumber "<<ObjectNumber <<endl;
648
649 // if(!sorted_tracks)sorted_tracks = new TRefArray();
650 // sorted_tracks->Clear();
651 // sorted_tracks.Clear();
652
653 if(!tsorted)tsorted = new TClonesArray("PamTrack",trk2_obj->GetNTracks());
654 tsorted->Delete();
655 TClonesArray &ttsorted = *tsorted;
656 if(!timage)timage = new TClonesArray("PamTrack",trk2_obj->GetNTracks());
657 timage->Delete();
658 TClonesArray &ttimage = *timage;
659
660
661 // loop over the tracks sorted by the tracker
662 Bool_t use_TRK = how.Contains("TRK", TString::kIgnoreCase);
663 Bool_t use_CAL = how.Contains("CAL", TString::kIgnoreCase);
664 Bool_t use_TOF = how.Contains("TOF", TString::kIgnoreCase);
665
666 if( !CAL2 && use_CAL) use_CAL = false;
667 if( !TOF && use_TOF) use_TOF = false;
668
669 if( !TRK2 ){
670 // cout << "SortTracks() : without tracker does not work!!! (not yet)" << endl;
671 return;
672 };
673
674 // cout << "use_CAL "<<use_CAL<<" use_TOF "<<use_TOF<<" use_TRK "<<use_TRK <<endl;
675
676 for(Int_t i=0; i < trk2_obj->TrkLevel2::GetNTracks(); i++){
677
678 TrkTrack *ts = 0;
679 CaloTrkVar *cs = 0;
680 ToFTrkVar *os = 0;
681
682 // get tracker tracks
683 TrkTrack *tp = trk2_obj->TrkLevel2::GetTrack(i); //tracker
684 CaloTrkVar *cp = GetCaloStoredTrack(tp->GetSeqNo());
685 ToFTrkVar *op = GetToFStoredTrack(tp->GetSeqNo());
686
687 TrkTrack *ti = 0; //tracker (image)
688 CaloTrkVar *ci = 0;
689 ToFTrkVar *oi = 0;
690 // cout << "trk track n. "<<i << " "<<hex<< tp <<dec<< endl;
691 // if track has an image, check image selection
692 if(tp->HasImage()){
693
694 ti = trk2_obj->TrkLevel2::GetTrackImage(i); //tracker (image)
695 ci = GetCaloStoredTrack(ti->GetSeqNo());
696 oi = GetToFStoredTrack(ti->GetSeqNo());
697
698 // cout << "its image "<<i << " "<<hex<< ti <<dec<< endl;
699
700 //assign starting scores
701 Int_t tp_score = 0; //main track sorted by the tracker
702 Int_t ti_score = 0; //image track
703
704 // ------------------------
705 // calorimeter check
706 // ------------------------
707 // check the Y spatial residual on the first calorimeter plane
708 // (cut on calorimeter variables from Emiliano)
709 if( use_CAL && !calo2_obj ){
710 cout << "void PamLevel2::SortTracks(): howtosort= "<<how<<" but CaloLevel2 not loaded !!!";
711 return;
712 };
713 if(
714 use_CAL &&
715 calo2_obj->npcfit[1] > 5 && //no. of fit planes on Y view
716 calo2_obj->varcfit[1] < 1000. && //fit variance on Y view
717 cp && ci &&
718 true){
719
720
721 Float_t resy_p = cp->tbar[0][1] - calo2_obj->cbar[0][1]; if(resy_p < 0)resy_p= - resy_p;
722 Float_t resy_i = ci->tbar[0][1] - calo2_obj->cbar[0][1]; if(resy_i < 0)resy_i= - resy_i;
723
724 if(resy_p <= resy_i) tp_score++;
725 else ti_score++;
726
727 // cout << "CALO "<<tp_score<<ti_score<<endl;
728
729 };
730 // ------------------------
731 // TOF check
732 // ------------------------
733 // check the number of hit pmts along the track
734 // on S12 S21 and S32, where paddles are parallel to Y axis
735 if( use_TOF && !tof_obj ){
736 cout << "void PamLevel2::SortTracks(): howtosort= "<<how<<" but ToFLevel2 not loaded !!!";
737 return;
738 };
739 if( use_TOF && op && oi ){
740
741 //
742 Float_t sen = 0.;
743 for (Int_t ih=0; ih < op->npmtadc; ih++){
744 Int_t pl = tof_obj->GetPlaneIndex( (op->pmtadc).At(ih) );
745 if ( pl == 2 || pl == 3 || pl == 4 || pl == 5 ) sen += (op->dedx).At(ih);
746 };
747 for (Int_t ih=0; ih < oi->npmtadc; ih++){
748 Int_t pl = tof_obj->GetPlaneIndex( (oi->pmtadc).At(ih) );
749 if ( pl == 2 || pl == 3 || pl == 4 || pl == 5 ) sen += (oi->dedx).At(ih);
750 };
751 //
752 if ( sen >= sortthr ){
753 //printf(" IS A NUCLEUS! en = %f \n",sen);
754 //
755 // is a nucleus use a different algorithm
756 //
757 Int_t nz = 6; Float_t zin[6]; // << define TOF z-coordinates
758 for(Int_t ip=0; ip<nz; ip++)
759 zin[ip] = tof_obj->ToFLevel2::GetZTOF(tof_obj->ToFLevel2::GetToFPlaneID(ip)); // << read ToF plane z-coordinates
760 Trajectory *tr = new Trajectory(nz,zin);
761 //
762 Int_t nphit_p =0;
763 Int_t nphit_i =0;
764 Float_t enhit_p = 0.;
765 Float_t enhit_i = 0.;
766 //
767 for (Int_t ih=0; ih < op->npmtadc; ih++){
768 Int_t pl = tof_obj->GetPlaneIndex( (op->pmtadc).At(ih) );
769 if(pl == 1 || pl == 2 || pl == 5){
770 nphit_p++;
771 enhit_p += (op->dedx).At(ih);
772 };
773 };
774 //
775 tp->DoTrack2(tr);
776 //
777 if ( fabs(tr->y[0]-oi->ytofpos[0]) < 2. ){
778 for (Int_t ih=0; ih < op->npmtadc; ih++){
779 Int_t pl = tof_obj->GetPlaneIndex( (op->pmtadc).At(ih) );
780 if(pl == 0){
781 nphit_p++;
782 enhit_p += (op->dedx).At(ih);
783 };
784 };
785 };
786 if ( fabs(tr->y[3]-oi->ytofpos[1]) < 2. ){
787 for (Int_t ih=0; ih < op->npmtadc; ih++){
788 Int_t pl = tof_obj->GetPlaneIndex( (op->pmtadc).At(ih) );
789 if(pl == 3){
790 nphit_p++;
791 enhit_p += (op->dedx).At(ih);
792 };
793 };
794 };
795 if ( fabs(tr->y[4]-oi->ytofpos[2]) < 2. ){
796 for (Int_t ih=0; ih < op->npmtadc; ih++){
797 Int_t pl = tof_obj->GetPlaneIndex( (op->pmtadc).At(ih) );
798 if(pl == 4){
799 nphit_p++;
800 enhit_p += (op->dedx).At(ih);
801 };
802 };
803 };
804
805 for (Int_t ih=0; ih < oi->npmtadc; ih++){
806 Int_t pl = tof_obj->GetPlaneIndex( (oi->pmtadc).At(ih) );
807 if(pl == 1 || pl == 2 || pl == 5){
808 nphit_i++;
809 enhit_i += (op->dedx).At(ih);
810 };
811 };
812 //
813 ti->DoTrack2(tr);
814 //
815 if ( fabs(tr->y[0]-oi->ytofpos[0]) < 2. ){
816 for (Int_t ih=0; ih < oi->npmtadc; ih++){
817 Int_t pl = tof_obj->GetPlaneIndex( (oi->pmtadc).At(ih) );
818 if(pl == 0){
819 nphit_i++;
820 enhit_i += (op->dedx).At(ih);
821 };
822 };
823 };
824 if ( fabs(tr->y[3]-oi->ytofpos[1]) < 2. ){
825 for (Int_t ih=0; ih < oi->npmtadc; ih++){
826 Int_t pl = tof_obj->GetPlaneIndex( (oi->pmtadc).At(ih) );
827 if(pl == 3){
828 nphit_i++;
829 enhit_i += (op->dedx).At(ih);
830 };
831 };
832 };
833 if ( fabs(tr->y[4]-oi->ytofpos[2]) < 2. ){
834 for (Int_t ih=0; ih < oi->npmtadc; ih++){
835 Int_t pl = tof_obj->GetPlaneIndex( (oi->pmtadc).At(ih) );
836 if(pl == 4){
837 nphit_i++;
838 enhit_i += (op->dedx).At(ih);
839 };
840 };
841 };
842
843
844 if(
845 use_TOF &&
846 (nphit_p+nphit_i) !=0 &&
847 true){
848
849 // printf(" seqno %i nphit_p %i nphit_i %i enhit_p %f enhit_i %f \n",trk2_obj->TrkLevel2::GetSeqNo(i),nphit_p,nphit_i,enhit_p,enhit_i);
850 // printf(" score p %i score i %i \n",tp_score,ti_score);
851 // if( enhit_p > enhit_i ) tp_score++;
852 // if( nphit_p >= nphit_i && enhit_p > enhit_i ) tp_score++;
853 if ( nphit_p > nphit_i ) tp_score++;
854 if ( nphit_p < nphit_i ) ti_score++;
855 if ( nphit_p == nphit_i ){
856 if ( enhit_p > enhit_i ) tp_score++;
857 else ti_score++;
858 };
859 // printf(" dopo score p %i score i %i \n",tp_score,ti_score);
860 };
861 delete tr;
862 //
863 } else {
864 //
865 // NOT a NUCLEUS
866 //
867 //printf(" NOT a NUCLEUS! en = %f \n",sen);
868
869 Int_t nphit_p =0;
870 Int_t nphit_i =0;
871
872
873 /* cout << "track: npmtadc "<< op->npmtadc << endl;
874 cout << "track: npmttdc "<< op->npmttdc << endl;
875 cout << "image: npmtadc "<< oi->npmtadc << endl;
876 cout << "image: npmttdc "<< oi->npmttdc << endl;*/
877
878 for (Int_t ih=0; ih < op->npmtadc; ih++){
879 Int_t pl = tof_obj->GetPlaneIndex( (op->pmtadc).At(ih) );
880 if(pl == 1 || pl == 2 || pl == 5)nphit_p++;
881 };
882
883 for (Int_t ih=0; ih < oi->npmtadc; ih++){
884 Int_t pl = tof_obj->GetPlaneIndex( (oi->pmtadc).At(ih) );
885 if(pl == 1 || pl == 2 || pl == 5)nphit_i++;
886 };
887
888 if(
889 use_TOF &&
890 (nphit_p+nphit_i) !=0 &&
891 true){
892
893 if( nphit_p >= nphit_i) tp_score++;
894 else ti_score++;
895 };
896 };
897 // cout << "TOF "<<tp_score<<ti_score<<endl;
898 };
899 if(tp_score == ti_score) use_TRK = true;
900 // ------------------------
901 // tracker check
902 // ------------------------
903 // chi**2 difference is not always large enough to distinguish among
904 // the real track and its image.
905 // Tracker check will be applied always when calorimeter and tof information is ambiguous.
906 if(use_TRK){
907 if( tp->chi2 > 0 && tp->chi2 < ti->chi2 ) tp_score++ ;
908 else if( ti->chi2 > 0 && ti->chi2 < tp->chi2 ) ti_score++ ;
909 // cout << "TRK "<<tp_score<<ti_score<<endl;
910 };
911
912 // ------------------------
913 // the winner is....
914 // ------------------------
915 if (tp_score > ti_score) {
916 // ts = tp;//the track sorted by the tracker!!
917 // cs = cp;
918 // os = op;
919 // cout << "It's primary!" << endl;
920
921 }else if (tp_score < ti_score) {
922
923 // cout << "It's image!" << endl;
924
925 ts = ti;//its image!!
926 cs = ci;
927 os = oi;
928
929 ti = tp;//its image!!
930 ci = cp;
931 oi = op;
932
933 tp = ts;//its image!!
934 cp = cs;
935 op = os;
936
937
938 }else {
939 // ts = tp;
940 // cs = cp;
941 // os = op;
942 // cout << "Warning - track image ambiguity not solved" << endl;
943 // cout << ts->GetNtot() << " "<< ts->chi2 << " " << npcfit[1] << " "<< nphit_p << endl;
944 };
945
946 }else{
947 // ts = tp;
948 // cs = cp;
949 // os = op;
950 };
951
952 // cout <<" SortTracks() "<<i<<" -- "<<ts<<endl;
953 // sorted_tracks->Add(ts);//save the track in the sorted array
954 // sorted_tracks.Add(ts);//save the track in the sorted array
955 // sorted_tracks.Add(tp);//save the track in the sorted array
956 // cout << "SortTracks:: sorted_tracks->Add(it) "<<i<<" "<<ts<<endl;
957 // cout<<"o "<<tp<<endl;
958 // cout<<"o "<<cp<<endl;
959 // cout<<"o "<<op<<endl;
960 new(ttsorted[i]) PamTrack(tp,cp,op);
961 new(ttimage[i]) PamTrack(ti,ci,oi);
962 };
963
964 if( tsorted->GetEntries() != trk2_obj->GetNTracks() ){
965 cout << "void PamLevel2::SortTracks(): tsorted->GetEntries() "<<tsorted->GetEntries()<<" != trk2_obj->GetNTracks() = "<<trk2_obj->GetNTracks() <<endl;
966 tsorted->Delete(); tsorted=0;
967 timage->Delete(); timage=0;
968 }
969
970 //Restore Object count
971 //To save space in the table keeping track of all referenced objects
972 //We reset the object count to what it was at the beginning of the event.
973 TProcessID::SetObjectCount(ObjectNumber);
974
975 };
976 //--------------------------------------
977 //
978 //
979 //--------------------------------------
980 /**
981 * This method overrides TrkLevel2::GetTracks(), where sorting is done by decreasing number of fit points and increasing chi^2.
982 * PamLevel2::GetTracks() keeps the same track order given by TrkLevel2::GetTracks(), but checks image selection by using calorimeter and ToF tracking information.
983 */
984 // TRefArray *PamLevel2::GetTracks(){
985
986 // // *-*-*-*-*-*-*-*-*-*-*-*-*
987 // SortTracks("+CAL+TOF");
988 // // *-*-*-*-*-*-*-*-*-*-*-*-*
989
990 // // return sorted_tracks;
991 // return &sorted_tracks;
992
993 // };
994 TClonesArray *PamLevel2::GetTracks(){
995
996 // *-*-*-*-*-*-*-*-*-*-*-*-*
997 SortTracks();
998 // *-*-*-*-*-*-*-*-*-*-*-*-*
999
1000 return tsorted;
1001
1002 };
1003 //--------------------------------------
1004 //
1005 //
1006 //--------------------------------------
1007 /**
1008 * Retrieves the it-th Pamela "physical" track.
1009 * It override TrkLevel2::GetTrack(int it).
1010 * @param it Track number, ranging from 0 to GetNTracks().
1011 */
1012 PamTrack *PamLevel2::GetTrack(int it){
1013
1014 // if(!trk2_obj) return 0;
1015
1016 // // *-*-*-*+-*-*-*-*-*-*-*-*-*
1017 // SortTracks("+CAL+TOF");
1018 // // *-*-*-*-*-*-*-*-*-*-*-*-*
1019 // // if(!sorted_tracks)return 0;
1020 // if(sorted_tracks.IsEmpty())return 0;
1021
1022 // PamTrack *track = 0;
1023
1024 // if( it >=0 && it < trk2_obj->TrkLevel2::GetNTracks() && it<sorted_tracks.GetEntriesFast() ){
1025 // TrkTrack *t = (TrkTrack*)sorted_tracks.At(it);
1026 // track = GetPamTrackAlong(t);
1027 // }else{
1028 // cout << "PamLevel2::GetTrack(int) : tracker track SeqNo "<< it <<" does not exist (GetNTracks() = "<<trk2_obj->TrkLevel2::GetNTracks()<<")"<<endl;
1029 // };
1030
1031 // return track;
1032
1033 // cout << "PamLevel2::GetTrack(int it) "<<endl;
1034 // *-*-*-*-*-*-*-*-*-*-*-*-*
1035 SortTracks();
1036 // *-*-*-*-*-*-*-*-*-*-*-*-*
1037 if(!tsorted)return 0;
1038 if(!tsorted->GetEntries())return 0;
1039
1040 PamTrack *track = 0;
1041
1042 if( it >=0 && it < trk2_obj->TrkLevel2::GetNTracks()){
1043 // TrkTrack *t = (TrkTrack*)sorted_tracks.At(it);
1044 // track = GetPamTrackAlong(t);
1045 TClonesArray &t = *(tsorted);
1046 track = (PamTrack*)t[it];
1047
1048 }else{
1049 cout << "PamLevel2::GetTrack(int) : tracker track SeqNo "<< it <<" does not exist (GetNTracks() = "<<trk2_obj->TrkLevel2::GetNTracks()<<")"<<endl;
1050 };
1051
1052 return track;
1053
1054 };
1055
1056 //--------------------------------------
1057 //
1058 //
1059 //--------------------------------------
1060 /**
1061 * Retrieves (if present) the image of the it-th Pamela "physical" track, sorted by the method PamLevel2::SortTracks().
1062 * @param it Track number, ranging from 0 to GetNTracks().
1063 */
1064 PamTrack *PamLevel2::GetTrackImage(int it){
1065
1066 // if(!trk2_obj) return 0;
1067
1068 // // *-*-*-*-*-*-*-*-*-*-*-*-*
1069 // SortTracks("+CAL+TOF");
1070 // // *-*-*-*-*-*-*-*-*-*-*-*-*
1071 // // if(!sorted_tracks)return 0;
1072 // if(sorted_tracks.IsEmpty())return 0;
1073
1074 // PamTrack *image = 0;
1075
1076 // if( it >=0 && it < trk2_obj->TrkLevel2::GetNTracks() ){
1077 // TrkTrack *temp = (TrkTrack*)sorted_tracks.At(it);
1078 // if( temp->HasImage() ){
1079 // TrkTrack *t = trk2_obj->TrkLevel2::GetStoredTrack(temp->GetImageSeqNo());
1080 // image = GetPamTrackAlong(t);
1081 // }else{
1082 // cout <<"PamLevel2::GetTrackImage(int) : Track SeqNo "<<it<<" does not have image"<<endl;
1083 // };
1084 // }else{
1085 // cout << "PamLevel2::GetTrackImage(int) : Tracker track SeqNo "<< it <<" does not exist (GetNTracks() = "<<trk2_obj->TrkLevel2::GetNTracks()<<")"<<endl;
1086 // };
1087
1088 // return image;
1089
1090
1091 // *-*-*-*-*-*-*-*-*-*-*-*-*
1092 SortTracks();
1093 // *-*-*-*-*-*-*-*-*-*-*-*-*
1094 if(!timage)return 0;
1095 if(!timage->GetEntries())return 0;
1096
1097 PamTrack *image = 0;
1098
1099 if( it >=0 && it < trk2_obj->TrkLevel2::GetNTracks() ){
1100 TClonesArray &t = *(tsorted);
1101 PamTrack *temp = (PamTrack*)t[it];
1102 if( temp->GetTrkTrack()->HasImage() ){
1103 TClonesArray &t = *(timage);
1104 image = (PamTrack*)t[it];
1105 }else{
1106 // cout <<"PamLevel2::GetTrackImage(int) : Track SeqNo "<<it<<" does not have image"<<endl;
1107 };
1108 }else{
1109 cout << "PamLevel2::GetTrackImage(int) : Tracker track SeqNo "<< it <<" does not exist (GetNTracks() = "<<trk2_obj->TrkLevel2::GetNTracks()<<")"<<endl;
1110 };
1111
1112 return image;
1113 }
1114
1115 //--------------------------------------
1116 //
1117 //
1118 //--------------------------------------
1119 /**
1120 * Get the Pamela detector trees in a single file and make them friends.
1121 * @param f TFile pointer
1122 * @param detlist String to select trees to be included
1123 * @return Pointer to a TTree
1124 */
1125 TTree *PamLevel2::GetPamTree(TFile *f, TString detlist ){
1126
1127 // if( !detlist.IsNull() ) SetWhichTrees(detlist);
1128 // else GetWhichTrees(f);
1129 if ( pam_tree ){
1130 printf("WARNING: TTree *PamLevel2::GetPamTree(TFile *fl, TString detlist) -- pam_tree already exists!\n ");
1131 return pam_tree;
1132 };
1133 //
1134
1135 cout << "TTree *PamLevel2::GetPamTree(TFile *f, TString detlist ) -- obsolte "<<endl;
1136
1137 SetWhichTrees(detlist);
1138
1139 TTree *Trout =0;
1140
1141 TString fname = f->GetName();
1142 if(!CheckLevel2File(fname))return NULL;
1143
1144 // UInt_t *found=0;
1145
1146 cout<< "GetPamTree(TFile*,TString): detector list --> ";
1147 if(TRK1)cout<<"TRK1 ";
1148 if(TRK2)cout<<"TRK2 ";
1149 if(TRKh)cout<<"TRKH ";
1150 if(CAL1)cout<<"CAL1 ";
1151 if(CAL2)cout<<"CAL2 ";
1152 if(TOF)cout<<"TOF ";
1153 if(TRG)cout<<"TRG ";
1154 if(AC)cout<<"AC ";
1155 if(ND)cout<<"ND ";
1156 if(S4)cout<<"S4 ";
1157 if(ORB)cout<<"ORB ";
1158 cout << endl;
1159 if(SELLI)cout<<">>> Found selection-list <<<"<<endl;
1160
1161 f->cd();
1162
1163 // Tracker
1164 TTree *T = (TTree*)f->Get("Tracker");
1165 if(T && (TRK2||TRK1||TRKh)) {
1166 if(TRK2)T->SetBranchAddress("TrkLevel2", GetPointerTo("TrkLevel2"));
1167 // else T->SetBranchStatus("TrkLevel2",0,found);
1168 if(TRK2)cout << "Tracker : set branch address TrkLevel2"<<endl;
1169 if(TRK1)T->SetBranchAddress("TrkLevel1", GetPointerTo("TrkLevel1"));
1170 // else T->SetBranchStatus("TrkLevel1",0,found);
1171 if(TRK1)cout << "Tracker : set branch address TrkLevel1"<<endl;
1172 if(TRKh)T->SetBranchAddress("TrkHough", GetPointerTo("TrkHough"));
1173 // else T->SetBranchStatus("TrkHough",0,found);
1174 if(TRKh)cout << "Tracker : set branch address TrkHough"<<endl;
1175 if(!Trout)Trout=T;
1176 else Trout->AddFriend(T);
1177 }else{
1178 cout << "Tracker : missing tree"<<endl;
1179 };
1180 // Calorimeter
1181 TTree *C = (TTree*)f->Get("Calorimeter");
1182 if(C && (CAL2||CAL1)) {
1183 if(CAL2)C->SetBranchAddress("CaloLevel2", GetPointerTo("CaloLevel2"));
1184 // else C->SetBranchStatus("CaloLevel2",0,found);
1185 if(CAL2)cout << "Calorimeter : set branch address CaloLevel2"<<endl;
1186 if(CAL1)C->SetBranchAddress("CaloLevel1", GetPointerTo("CaloLevel1"));
1187 // else C->SetBranchStatus("CaloLevel1",0,found);
1188 if(CAL1)cout << "Calorimeter : set branch address CaloLevel1"<<endl;
1189 if(!Trout)Trout=C;
1190 else Trout->AddFriend(C);
1191 }else{
1192 cout << "Calorimeter : missing tree"<<endl;
1193 };
1194
1195 // ToF
1196 TTree *O = (TTree*)f->Get("ToF");
1197 if(O && TOF) {
1198 O->SetBranchAddress("ToFLevel2", GetPointerTo("ToFLevel2"));
1199 cout << "ToF : set branch address ToFLevel2"<<endl;
1200 if(!Trout)Trout=O;
1201 else Trout->AddFriend(O);
1202 }else{
1203 cout << "ToF : missing tree"<<endl;
1204 };
1205 // Trigger
1206 TTree *R = (TTree*)f->Get("Trigger");
1207 if(R && TRG) {
1208 R->SetBranchAddress("TrigLevel2", GetPointerTo("TrigLevel2"));
1209 cout << "Trigger : set branch address TrigLevel2"<<endl;
1210 if(!Trout)Trout=R;
1211 else Trout->AddFriend(R);
1212 }else{
1213 cout << "Trigger : missing tree"<<endl;
1214 };
1215 // S4
1216 TTree *S = (TTree*)f->Get("S4");
1217 if(S && S4) {
1218 S->SetBranchAddress("S4Level2", GetPointerTo("S4Level2"));
1219 cout << "S4 : set branch address S4Level2"<<endl;
1220 if(!Trout)Trout=S;
1221 else Trout->AddFriend(S);
1222 }else{
1223 cout << "S4 : missing tree"<<endl;
1224 };
1225 // Neutron Detector
1226 TTree *N = (TTree*)f->Get("NeutronD");
1227 if(N && ND) {
1228 N->SetBranchAddress("NDLevel2", GetPointerTo("NDLevel2"));
1229 cout << "NeutronD : set branch address NDLevel2"<<endl;
1230 if(!Trout)Trout=N;
1231 else Trout->AddFriend(N);
1232 }else{
1233 cout << "NeutronD : missing tree"<<endl;
1234 };
1235 // Anticounters
1236 TTree *A = (TTree*)f->Get("Anticounter");
1237 if(A && AC) {
1238 A->SetBranchAddress("AcLevel2", GetPointerTo("AcLevel2"));
1239 cout << "Anticounter : set branch address AcLevel2"<<endl;
1240 if(!Trout)Trout=A;
1241 else Trout->AddFriend(A);
1242 }else{
1243 cout << "Anticounter : missing tree"<<endl;
1244 };
1245 // Orbital Info
1246 TTree *B = (TTree*)f->Get("OrbitalInfo");
1247 if(B && ORB) {
1248 B->SetBranchAddress("OrbitalInfo", GetPointerTo("OrbitalInfo"));
1249 cout << "OrbitalInfo : set branch address OrbitalInfo"<<endl;
1250 if(!Trout)Trout=B;
1251 else Trout->AddFriend(B);
1252 }else{
1253 cout << "OrbitalInfo : missing tree"<<endl;
1254 };
1255
1256 TTree *L = (TTree*)f->Get("SelectionList");
1257 if(L && SELLI==1) {
1258 // L->SetBranchAddress("RunEntry",&irun);
1259 // cout << "SelectionList: set branch address RunEntry"<<endl;
1260 // L->SetBranchAddress("EventEntry",&irunentry);
1261 // cout << "SelectionList: set branch address EventEntry"<<endl;
1262 // if(!Trout)Trout=O;
1263 // else Trout->AddFriend("SelectionList");
1264 cout << " TTree *PamLevel2::GetPamTree(TFile, TString) >>> SelectionList not implemente!!!!"<<endl;
1265 sel_tree = 0;
1266 }else{
1267 cout << "SelectionList : missing tree"<<endl;
1268 };
1269
1270 cout<<endl<<" Number of entries: "<<Trout->GetEntries()<<endl<<endl;
1271
1272 pam_tree = (TChain*)Trout;
1273
1274 return Trout;
1275
1276 }
1277 //--------------------------------------
1278 //
1279 //
1280 //--------------------------------------
1281 /**
1282 * Get list of Level2 files.
1283 * @param ddir Level2 data directory.
1284 * @param flisttxt Name of txt file containing file list.
1285 * @return Pointer to a TList of TSystemFiles
1286 * If no input file list is given , all the Level2 files inside the directory are processed.
1287 */
1288 TList* PamLevel2::GetListOfLevel2Files(TString ddir, TString flisttxt = ""){
1289
1290 TString wdir = gSystem->WorkingDirectory();
1291
1292 // if(ddir=="")ddir = wdir;
1293 // TSystemDirectory *datadir = new TSystemDirectory(gSystem->BaseName(ddir),ddir);
1294 if ( ddir != ""){
1295 cout << "Level2 data directory : "<< ddir << endl;
1296 } else {
1297 cout << "Level2 data directory not given as input: trying to evaluate from list or taking working directory " << endl;
1298 };
1299 TList *contents = new TList; // create output list
1300 contents->SetOwner();
1301
1302 // char *fullpath;
1303 // const char *fullpath;
1304
1305 // if no input file list is given:
1306 if ( flisttxt != "" ){
1307
1308 // if( !gSystem->IsFileInIncludePath(flisttxt,&fullpath) ){
1309 // if( !(fullpath = gSystem->FindFile("./",flisttxt)) ){
1310 // cout <<"File "<<flisttxt<<" not found"<<endl;
1311 // return 0;
1312 // }
1313 // flisttxt = fullpath;
1314 if ( !flisttxt.EndsWith(".root") ){
1315
1316 flisttxt = gSystem->ConcatFileName(gSystem->DirName(flisttxt),gSystem->BaseName(flisttxt));
1317
1318 if( !gSystem->ChangeDirectory(ddir) ){
1319 cout << "Cannot change directory : "<<ddir<<endl;
1320 return 0;
1321 }
1322
1323 cout <<"Input file list : " << flisttxt <<endl;
1324 ifstream in;
1325 in.open(flisttxt, ios::in); //open input file list
1326 if(!in.good()){
1327 cout <<" ERROR opening the file "<<endl;
1328 gSystem->ChangeDirectory(wdir); // back to the working directory
1329 return 0;
1330 }
1331 int line=0;
1332 while (1) {
1333 TString file;
1334 in >> file;
1335 if (!in.good()) break;
1336 line++;
1337 // cout <<"(1) " << file << endl;
1338 if(file.IsNull()){
1339 cout << "-- list interrupted at line "<<line <<endl;
1340 break;
1341 }
1342 if(file.Contains("#"))file = file(0,file.First("#"));
1343 // cout <<"(2) " << file << endl;
1344 // if( gSystem->IsFileInIncludePath(file,&fullpath) ){
1345 // if( (fullpath = gSystem->FindFile(ddir,file)) ){
1346 if( file.EndsWith(".root") ){
1347 TString filedir;
1348 if (ddir != ""){
1349 filedir = ddir; // take the input dir
1350 } else {
1351 gSystem->ChangeDirectory(wdir); // back to the working directory
1352 filedir = gSystem->DirName(file); // this will take the path if exist in the list otherwise it will return automatically the working dir
1353 };
1354 char *fullpath = gSystem->ConcatFileName(gSystem->DirName(filedir),gSystem->BaseName(file));
1355 contents->Add(new TSystemFile(fullpath,gSystem->DirName(fullpath)));// add file to the list
1356 delete fullpath;
1357 }
1358 // }else{
1359 // if(file.Data()!="")cout << "File: "<<file<<" ---> missing "<< endl;
1360 // };
1361 };
1362 in.close();
1363 } else {
1364 if(flisttxt.Contains("#"))flisttxt = flisttxt(0,flisttxt.First("#"));
1365 char *fullpath = gSystem->ConcatFileName(gSystem->DirName(flisttxt),gSystem->BaseName(flisttxt));
1366 contents->Add(new TSystemFile(fullpath,gSystem->DirName(fullpath)));// add file to the list
1367 delete fullpath;
1368 };
1369 }else{
1370
1371 cout << "No input file list given."<<endl;
1372 cout << "Check for existing root files."<<endl;
1373 // cout << "Warking directory: "<< gSystem->WorkingDirectory()<< endl;
1374 if ( ddir == "" ){
1375 ddir = wdir;
1376 cout << "Level2 data directory : "<< ddir << endl;
1377 };
1378
1379 TSystemDirectory *datadir = new TSystemDirectory(gSystem->BaseName(ddir),ddir);
1380 TList *temp = datadir->GetListOfFiles();
1381 // temp->Print();
1382 // cout << "*************************************" << endl;
1383
1384 TIter next(temp);
1385 TSystemFile *questo = 0;
1386
1387
1388 while ( (questo = (TSystemFile*) next()) ) {
1389 TString name = questo-> GetName();
1390 if( name.EndsWith(".root") ){
1391 // const char *fullpath = gSystem->FindFile(ddir,name);
1392 // char *fullpath;
1393 // gSystem->IsFileInIncludePath(name,&fullpath);
1394 char *fullpath = gSystem->ConcatFileName(gSystem->DirName(ddir),gSystem->BaseName(name));
1395 contents->Add(new TSystemFile(fullpath,gSystem->DirName(fullpath)));
1396 delete fullpath;
1397 };
1398 }
1399 delete temp;
1400 delete datadir;
1401
1402 };
1403 gSystem->ChangeDirectory(wdir); // back to the working directory
1404 // cout << endl << "Selected files:" << endl;
1405 // contents->Print();
1406 cout << contents->GetEntries()<<" files \n";
1407 // cout << endl;
1408 // cout << "Working directory: "<< gSystem->WorkingDirectory()<< endl;
1409 return contents;
1410 };
1411 //--------------------------------------
1412 //
1413 //
1414 //--------------------------------------
1415 /**
1416 * Get the Pamela detector chains from a list of files and make them friends.
1417 * @param fl Pointer to a TList of TSystemFiles
1418 * @param detlist String to select trees to be included
1419 * @return Pointer to a TChain
1420 */
1421 TChain *PamLevel2::GetPamTree(TList *fl, TString detlist ){
1422 //
1423 //
1424 //
1425 if ( pam_tree ){
1426 printf("WARNING: TChain *PamLevel2::GetPamTree(TList *fl, TString detlist) -- pam_tree already exists!\n ");
1427 return pam_tree;
1428 };
1429 //
1430
1431 TChain *Trout =0;
1432
1433 // if( !detlist.IsNull() )SetWhichTrees(detlist);
1434 SetWhichTrees(detlist);
1435
1436 cout<< "GetPamTree(TList*,TString): input detector list --> ";
1437 if(TRK1)cout<<"TRK1 ";
1438 if(TRK2)cout<<"TRK2 ";
1439 if(TRKh)cout<<"TRKH ";
1440 if(CAL1)cout<<"CAL1 ";
1441 if(CAL2)cout<<"CAL2 ";
1442 if(TOF)cout<<"TOF ";
1443 if(TRG)cout<<"TRG ";
1444 if(AC)cout<<"AC ";
1445 if(ND)cout<<"ND ";
1446 if(S4)cout<<"S4 ";
1447 if(ORB)cout<<"ORB ";
1448 cout << endl;
1449
1450 TChain *T = 0;
1451 TChain *C = 0;
1452 TChain *O = 0;
1453 TChain *R = 0;
1454 TChain *S = 0;
1455 TChain *N = 0;
1456 TChain *A = 0;
1457 TChain *B = 0;
1458
1459 TChain *L = 0;
1460
1461 if(TRK2||TRK1||TRKh) T = new TChain("Tracker");
1462 if(CAL2||CAL1) C = new TChain("Calorimeter");
1463 if(TOF) O = new TChain("ToF");
1464 if(TRG) R = new TChain("Trigger");
1465 if(S4) S = new TChain("S4");
1466 if(ND) N = new TChain("NeutronD");
1467 if(AC) A = new TChain("Anticounter");
1468 if(ORB) B = new TChain("OrbitalInfo");
1469 L = new TChain("SelectionList");
1470
1471 // loop over files and create chains
1472 TIter next(fl);
1473 TSystemFile *questo = 0;
1474 while ( (questo = (TSystemFile*) next()) ) {
1475 TString name = questo->GetName();
1476 cout << "File: "<< name << endl;
1477 if( CheckLevel2File(name) ){
1478 if(TRK2||TRK1||TRKh) T->Add(name);
1479 if(CAL1||CAL2) C->Add(name);
1480 if(TOF) O->Add(name);
1481 if(TRG) R->Add(name);
1482 if(S4) S->Add(name);
1483 if(ND) N->Add(name);
1484 if(AC) A->Add(name);
1485 if(ORB) B->Add(name);
1486 if(SELLI==1) L->Add(name);
1487 };
1488 };
1489
1490 cout << "done chain \n";
1491 cout << "+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+" <<endl;
1492
1493 // UInt_t *found=0;
1494 // Tracker
1495 if(T && (TRK2||TRK1||TRKh)) {
1496 if(TRK2)T->SetBranchAddress("TrkLevel2", GetPointerTo("TrkLevel2"));
1497 // else T->SetBranchStatus("TrkLevel2",0,found);
1498 if(TRK2)cout << "Tracker : set branch address TrkLevel2"<<endl;
1499 if(TRK1)T->SetBranchAddress("TrkLevel1", GetPointerTo("TrkLevel1"));
1500 // else T->SetBranchStatus("TrkLevel1",0,found);
1501 if(TRK1)cout << "Tracker : set branch address TrkLevel1"<<endl;
1502 if(TRKh)T->SetBranchAddress("TrkHough", GetPointerTo("TrkHough"));
1503 // else T->SetBranchStatus("TrkHough",0,found);
1504 if(TRKh)cout << "Tracker : set branch address TrkHough"<<endl;
1505 if(!Trout)Trout=T;
1506 else Trout->AddFriend("Tracker");
1507 }else{
1508 cout << "Tracker : missing tree"<<endl;
1509 };
1510 // Calorimeter
1511 if(C && (CAL2||CAL1)) {
1512 if(CAL2)C->SetBranchAddress("CaloLevel2", GetPointerTo("CaloLevel2"));
1513 // else C->SetBranchStatus("CaloLevel2",0,found);
1514 if(CAL2)cout << "Calorimeter : set branch address CaloLevel2"<<endl;
1515 if(CAL1)C->SetBranchAddress("CaloLevel1", GetPointerTo("CaloLevel1"));
1516 // else C->SetBranchStatus("CaloLevel1",0,found);
1517 if(CAL1)cout << "Calorimeter : set branch address CaloLevel1"<<endl;
1518 if(!Trout)Trout=C;
1519 else Trout->AddFriend("Calorimeter");
1520 }else{
1521 cout << "Calorimeter : missing tree"<<endl;
1522 };
1523 // ToF
1524 if(O && TOF) {
1525 O->SetBranchAddress("ToFLevel2", GetPointerTo("ToFLevel2"));
1526 cout << "ToF : set branch address ToFLevel2"<<endl;
1527 if(!Trout)Trout=O;
1528 else Trout->AddFriend("ToF");
1529 }else{
1530 cout << "ToF : missing tree"<<endl;
1531 };
1532 // Trigger
1533 if(R && TRG) {
1534 R->SetBranchAddress("TrigLevel2", GetPointerTo("TrigLevel2"));
1535 cout << "Trigger : set branch address TrigLevel2"<<endl;
1536 if(!Trout)Trout=O;
1537 else Trout->AddFriend("Trigger");
1538 }else{
1539 cout << "Trigger : missing tree"<<endl;
1540 };
1541 // S4
1542 if(S && S4) {
1543 S->SetBranchAddress("S4Level2", GetPointerTo("S4Level2"));
1544 cout << "S4 : set branch address S4Level2"<<endl;
1545 if(!Trout)Trout=O;
1546 else Trout->AddFriend("S4");
1547 }else{
1548 cout << "S4 : missing tree"<<endl;
1549 };
1550 // Neutron Detector
1551 if(N && ND) {
1552 N->SetBranchAddress("NDLevel2", GetPointerTo("NDLevel2"));
1553 cout << "NeutronD : set branch address NDLevel2"<<endl;
1554 if(!Trout)Trout=O;
1555 else Trout->AddFriend("NeutronD");
1556 }else{
1557 cout << "NeutronD : missing tree"<<endl;
1558 };
1559 // Anticounters
1560 if(A && AC) {
1561 A->SetBranchAddress("AcLevel2", GetPointerTo("AcLevel2"));
1562 cout << "Anticounter : set branch address AcLevel2"<<endl;
1563 if(!Trout)Trout=O;
1564 else Trout->AddFriend("Anticounter");
1565 }else{
1566 cout << "Anticounter : missing tree"<<endl;
1567 };
1568 // Orbital Info
1569 if(B && ORB) {
1570 B->SetBranchAddress("OrbitalInfo", GetPointerTo("OrbitalInfo"));
1571 cout << "OrbitalInfo : set branch address OrbitalInfo"<<endl;
1572 if(!Trout)Trout=O;
1573 else Trout->AddFriend("OrbitalInfo");
1574 }else{
1575 cout << "OrbitalInfo : missing tree"<<endl;
1576 };
1577 cout << "+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+" <<endl;
1578
1579 // Selection List
1580 if(L && SELLI==1) {
1581 cout<<">>> Found selection-list <<<"<<endl;
1582 L->SetBranchAddress("RunEntry",&irun);
1583 cout << "SelectionList: set branch address RunEntry"<<endl;
1584 L->SetBranchAddress("EventEntry",&irunentry);
1585 cout << "SelectionList: set branch address EventEntry"<<endl;
1586 sel_tree = L;
1587 // if(!Trout)Trout=O;
1588 // else Trout->AddFriend("SelectionList");
1589 cout << "----------------------------------------------------" <<endl;
1590 }else{
1591 // cout << "SelectionList : missing tree"<<endl;
1592 if(L)L->Delete();
1593 };
1594
1595 // cout<<endl<<" Number of entries: "<<Trout->GetEntries()<<endl<<endl;
1596
1597 pam_tree = Trout;
1598
1599 return Trout;
1600 }
1601
1602
1603
1604 //--------------------------------------
1605 //
1606 //
1607 //--------------------------------------
1608 /**
1609 * Set branch addresses for Pamela friend trees
1610 */
1611 void PamLevel2::SetBranchAddress(TTree *t){
1612
1613 TRK2 = TRK2 & t->GetBranchStatus("TrkLevel2");
1614 TRK1 = TRK1 & t->GetBranchStatus("TrkLevel1");
1615 TRKh = TRKh & t->GetBranchStatus("TrkHough");
1616 CAL2 = CAL2 & t->GetBranchStatus("CaloLevel2");
1617 CAL1 = CAL1 & t->GetBranchStatus("CaloLevel1");
1618 TOF = TOF & t->GetBranchStatus("ToFLevel2");
1619 TRG = TRG & t->GetBranchStatus("TrigLevel2");
1620 S4 = S4 & t->GetBranchStatus("S4Level2");
1621 ND = ND & t->GetBranchStatus("NDLevel2");
1622 AC = AC & t->GetBranchStatus("AcLevel2");
1623 ORB = ORB & t->GetBranchStatus("OrbitalInfo");
1624
1625
1626 // Tracker
1627 if(TRK1) {
1628 t->SetBranchAddress("TrkLevel1", GetPointerTo("TrkLevel2"));
1629 cout << "Tracker : set branch address TrkLevel1"<<endl;
1630 };
1631 if(TRK2) {
1632 t->SetBranchAddress("TrkLevel2", GetPointerTo("TrkLevel1"));
1633 cout << "Tracker : set branch address TrkLevel2"<<endl;
1634 };
1635 if(TRKh) {
1636 t->SetBranchAddress("TrkHough", GetPointerTo("TrkHough"));
1637 cout << "Tracker : set branch address TrkHough"<<endl;
1638 };
1639
1640 // Calorimeter
1641 if(CAL1) {
1642 t->SetBranchAddress("CaloLevel1", GetPointerTo("CaloLevel1"));
1643 cout << "Calorimeter : set branch address CaloLevel1"<<endl;
1644 };
1645 if(CAL2) {
1646 t->SetBranchAddress("CaloLevel2", GetPointerTo("CaloLevel2"));
1647 cout << "Calorimeter : set branch address CaloLevel2"<<endl;
1648 };
1649
1650 // ToF
1651 if(TOF) {
1652 t->SetBranchAddress("ToFLevel2", GetPointerTo("ToFLevel2"));
1653 cout << "ToF : set branch address ToFLevel2"<<endl;
1654 };
1655 // Trigger
1656 if(TRG) {
1657 t->SetBranchAddress("TrigLevel2", GetPointerTo("TrigLevel2"));
1658 cout << "Trigger : set branch address TrigLevel2"<<endl;
1659 };
1660 // S4
1661 if(S4) {
1662 t->SetBranchAddress("S4Level2", GetPointerTo("S4Level2"));
1663 cout << "S4 : set branch address S4Level2"<<endl;
1664 };
1665 // Neutron Detector
1666 if(ND) {
1667 t->SetBranchAddress("NDLevel2", GetPointerTo("NDLevel2"));
1668 cout << "NeutronD : set branch address NDLevel2"<<endl;
1669 };
1670 // Anticounters
1671 if(AC) {
1672 t->SetBranchAddress("AcLevel2", GetPointerTo("AcLevel2"));
1673 cout << "Anticounter : set branch address AcLevel2"<<endl;
1674 };
1675 // OrbitalInfo
1676 if(ORB) {
1677 t->SetBranchAddress("OrbitalInfo", GetPointerTo("OrbitalInfo"));
1678 cout << "OrbitalInfo : set branch address OrbitalInfo"<<endl;
1679 };
1680 // SelectionList
1681 if(SELLI==1) {
1682 t->SetBranchAddress("RunEntry", &irun);
1683 cout << "SelectionList: set branch address RunEntry"<<endl;
1684 t->SetBranchAddress("EventEntry", &irunentry);
1685 cout << "SelectionList: set branch address EventEntry"<<endl;
1686 };
1687
1688 }
1689 /**
1690 * Set branch addresses for Pamela friend trees
1691 */
1692 void PamLevel2::SetBranchAddress(TChain *t){
1693
1694 TRK2 = TRK2 & t->GetBranchStatus("TrkLevel2");
1695 TRK1 = TRK1 & t->GetBranchStatus("TrkLevel1");
1696 TRKh = TRKh & t->GetBranchStatus("TrkHough");
1697 CAL1 = CAL1 & t->GetBranchStatus("CaloLevel1");
1698 CAL2 = CAL2 & t->GetBranchStatus("CaloLevel2");
1699 TOF = TOF & t->GetBranchStatus("ToFLevel2");
1700 TRG = TRG & t->GetBranchStatus("TrigLevel2");
1701 S4 = S4 & t->GetBranchStatus("S4Level2");
1702 ND = ND & t->GetBranchStatus("NDLevel2");
1703 AC = AC & t->GetBranchStatus("AcLevel2");
1704 ORB = ORB & t->GetBranchStatus("OrbitalInfo");
1705
1706 // Tracker
1707 if(TRK2) {
1708 t->SetBranchAddress("TrkLevel2", GetPointerTo("TrkLevel2"));
1709 cout << "Tracker : set branch address TrkLevel2"<<endl;
1710 };
1711 if(TRK1) {
1712 t->SetBranchAddress("TrkLevel1", GetPointerTo("TrkLevel1"));
1713 cout << "Tracker : set branch address TrkLevel1"<<endl;
1714 };
1715 if(TRKh) {
1716 t->SetBranchAddress("TrkHough", GetPointerTo("TrkHough"));
1717 cout << "Tracker : set branch address TrkHough"<<endl;
1718 };
1719
1720 // Calorimeter
1721 if(CAL2) {
1722 t->SetBranchAddress("CaloLevel2", GetPointerTo("CaloLevel2"));
1723 cout << "Calorimeter : set branch address CaloLevel2"<<endl;
1724 };
1725 if(CAL1) {
1726 t->SetBranchAddress("CaloLevel1", GetPointerTo("CaloLevel1"));
1727 cout << "Calorimeter : set branch address CaloLevel1"<<endl;
1728 };
1729
1730 // ToF
1731 if(TOF) {
1732 t->SetBranchAddress("ToFLevel2", GetPointerTo("ToFLevel2"));
1733 cout << "ToF : set branch address ToFLevel2"<<endl;
1734 };
1735 // Trigger
1736 if(TRG) {
1737 t->SetBranchAddress("TrigLevel2", GetPointerTo("TrigLevel2"));
1738 cout << "Trigger : set branch address TrigLevel2"<<endl;
1739 };
1740 // S4
1741 if(S4) {
1742 t->SetBranchAddress("S4Level2", GetPointerTo("S4Level2"));
1743 cout << "S4 : set branch address S4Level2"<<endl;
1744 };
1745 // Neutron Detector
1746 if(ND) {
1747 t->SetBranchAddress("NDLevel2", GetPointerTo("NDLevel2"));
1748 cout << "NeutronD : set branch address NDLevel2"<<endl;
1749 };
1750 // Anticounters
1751 if(AC) {
1752 t->SetBranchAddress("AcLevel2", GetPointerTo("AcLevel2"));
1753 cout << "Anticounter : set branch address AcLevel2"<<endl;
1754 };
1755 // OrbitalInfo
1756 if(ORB) {
1757 t->SetBranchAddress("OrbitalInfo", GetPointerTo("OrbitalInfo"));
1758 cout << "OrbitalInfo : set branch address OrbitalInfo"<<endl;
1759 };
1760 // SelectionList
1761 if(SELLI==1) {
1762 t->SetBranchAddress("RunEntry", &irun);
1763 cout << "SelectionList: set branch address RunEntry"<<endl;
1764 t->SetBranchAddress("EventEntry", &irunentry);
1765 cout << "SelectionList: set branch address EventEntry"<<endl;
1766 };
1767
1768 }
1769
1770
1771 //--------------------------------------
1772 //
1773 //
1774 //--------------------------------------
1775 /**
1776 * Get the Run tree chain from a list of files.
1777 * @param fl Pointer to a TList of TSystemFiles
1778 * @return Pointer to a TChain
1779 */
1780 TChain *PamLevel2::GetRunTree(TList *fl){
1781 //
1782 //
1783 //
1784 if ( run_tree ){
1785 printf("WARNING: TChain *PamLevel2::GetRunTree(TList *fl) -- run_tree already exists!\n ");
1786 return run_tree;
1787 };
1788 //
1789 TChain *R = new TChain("Run");
1790
1791 // loop over files and create chains
1792 TIter next(fl);
1793 TSystemFile *questo = 0;
1794 while ( (questo = (TSystemFile*) next()) ) {
1795 TString name = questo->GetName();
1796 // cout << "File: "<< name << endl;
1797 if( CheckLevel2File(name) ){
1798 R->Add(name);
1799 };
1800 }
1801
1802 if(R->GetNtrees()){
1803
1804 R->SetBranchAddress("RunInfo", GetPointerTo("RunInfo"));
1805 cout << "Run : set branch address RunInfo"<<endl;
1806 R->SetBranchAddress("SoftInfo", GetPointerTo("SoftInfo")); // Emiliano
1807 cout << "Software : set branch address SoftInfo"<<endl; // Emiliano
1808 }else{
1809 delete R;
1810 R=0;
1811 }
1812
1813 run_tree = R;
1814
1815 return R;
1816
1817 }
1818 //--------------------------------------
1819 //
1820 //
1821 //--------------------------------------
1822 /**
1823 * Get the Run tree from a file.
1824 * @param f Pointer to a TFile
1825 * @return Pointer to a TTree
1826 */
1827 TTree *PamLevel2::GetRunTree(TFile *f){
1828 if ( run_tree ){
1829 printf("WARNING: TTree *PamLevel2::GetRunTree(TFile *f) -- run_tree already exists!\n ");
1830 return run_tree;
1831 };
1832
1833
1834 cout << "TTree *PamLevel2::GetRunTree(TFile *f) -- obsolte "<<endl;
1835
1836 TTree *T = (TTree*)f->Get("Run");
1837
1838 if(T){
1839 T->SetBranchAddress("RunInfo", GetPointerTo("RunInfo"));
1840 cout << "Run : set branch address RunInfo"<<endl;
1841 T->SetBranchAddress("SoftInfo", GetPointerTo("SoftInfo")); // Emiliano
1842 cout << "Software : set branch address SoftInfo"<<endl; // Emiliano
1843 }
1844
1845 run_tree = (TChain*)T;
1846
1847 return T;
1848
1849 }
1850 /**
1851 * Update the runinfo informations (to be used to have Run infos event by event basis)
1852 * @param run Pointer to the chain/tree which contains run infos
1853 * @return true if a new run has been read, false if it is still the same run
1854 */
1855 Bool_t PamLevel2::UpdateRunInfo(TChain *run, Long64_t iev){
1856 //
1857 // check if we have already called once GetEntry, if not call it
1858 //
1859 cout << "Bool_t PamLevel2::UpdateRunInfo(TChain *run, Long64_t iev) --- NON FUNZIONA BENISSIMO.... "<<endl;
1860 if(!run){
1861 cout << "Bool_t PamLevel2::UpdateRunInfo(TChain *run, ULong64_t iev) -- ERROR -- missing RunInfo tree "<<endl;
1862 return(false);
1863 }
1864 if ( run->GetEntries() <= 0 ) return(false);
1865 //
1866
1867 // Int_t oldrun = irun;
1868 Long64_t oldrun = irun;
1869
1870 // --------------------------------------
1871 // if it is a full file (not preselected)
1872 // --------------------------------------
1873 if(SELLI==0){
1874
1875 //
1876 // the absolute time is necessary to relate the event with the run
1877 //
1878 if( !GetOrbitalInfo() ){
1879 cout << "Bool_t PamLevel2::UpdateRunInfo(TChain *run, ULong64_t iev) -- ERROR -- missing OrbitalInfo "<<endl;
1880 return(false);
1881 }
1882
1883 //
1884 // the first time the routine is called, set run search from the beginning
1885 //
1886 if ( irun < 0LL ){
1887 irun = 0LL;
1888 run->GetEntry(irun);
1889 runfirstentry = 0LL;
1890 runlastentry += (Long64_t)(this->GetRunInfo()->NEVENTS);
1891 if ( (Long64_t)(this->GetRunInfo()->NEVENTS) > 0LL ) runlastentry -= 1LL;
1892 };
1893 //
1894 if ( irun == run->GetEntries()-1LL &&
1895 !(GetOrbitalInfo()->absTime >= GetRunInfo()->RUNHEADER_TIME &&
1896 GetOrbitalInfo()->absTime <= GetRunInfo()->RUNTRAILER_TIME)
1897 ){
1898 irun = -1LL;
1899 runfirstentry = 0LL;
1900 runlastentry = -1LL;
1901 };
1902 // modificato il controllo sull'aggiornamento del run, per evitare problemi
1903 // dovuti agli eventi annidati (NB! NEVENTS conta anche questi!!)
1904 //
1905 bool fromfirst = true;
1906 //
1907 while ( !(GetOrbitalInfo()->absTime >= GetRunInfo()->RUNHEADER_TIME && GetOrbitalInfo()->absTime <= GetRunInfo()->RUNTRAILER_TIME) && irun < run->GetEntries()-1LL ){
1908 // while ( iev > (runfirstentry+(ULong64_t)(this->GetRunInfo()->NEVENTS-1)) && irun < run->GetEntries() ){
1909 irun++;
1910 run->GetEntry(irun);
1911 runfirstentry = runlastentry;
1912 if ( (Long64_t)(this->GetRunInfo()->NEVENTS) > 0LL ) runfirstentry += 1LL;
1913 runlastentry += (Long64_t)(this->GetRunInfo()->NEVENTS);
1914 // cout << " ))))) UPDATE RUN INFO ((((( @iev "<<iev<<" run "<<GetRunInfo()->ID<<" irun "<<irun<<endl;
1915 // cout << "runfirstentry "<<runfirstentry<<endl;
1916 // printf(" iev %llu %u %llu \n",iev,this->GetRunInfo()->NEVENTS,(ULong64_t)(runfirstentry+(ULong64_t)(this->GetRunInfo()->NEVENTS)));
1917 // printf(" abstime %u trailertime %u \n",GetOrbitalInfo()->absTime,GetRunInfo()->RUNTRAILER_TIME);
1918 // printf(" IDRUN %u \n",GetRunInfo()->ID);
1919 //
1920 // prevshift = 0;
1921 //
1922 if ( irun == (Long64_t)(run->GetEntries()-1LL) && fromfirst && !(GetOrbitalInfo()->absTime >= GetRunInfo()->RUNHEADER_TIME && GetOrbitalInfo()->absTime <= GetRunInfo()->RUNTRAILER_TIME)){
1923 printf(" resetting irun (it should NOT happen!!!)\n");
1924 fromfirst = false;
1925 irun = 0;
1926 run->GetEntry(irun);
1927 runfirstentry = 0ULL;
1928 runlastentry += (Long64_t)(this->GetRunInfo()->NEVENTS);
1929 if ( (Long64_t)(this->GetRunInfo()->NEVENTS) > 0LL ) runlastentry -= 1LL;
1930 };
1931 //
1932 };
1933 //
1934 if (
1935 !(GetOrbitalInfo()->absTime >= GetRunInfo()->RUNHEADER_TIME &&
1936 GetOrbitalInfo()->absTime <= GetRunInfo()->RUNTRAILER_TIME)
1937 ) {
1938 printf(" Something very wrong here: cannot find RUN containing absolute time %u \n",GetOrbitalInfo()->absTime);
1939 return false;
1940 }
1941 //
1942 if ( irun == oldrun || irun >= run->GetEntries() ) return(false);
1943 //
1944 // printf(" iev %llu irun %i nevents %u 1st %llu last %llu \n",iev,irun,this->GetRunInfo()->NEVENTS,(ULong64_t)runfirstentry,(ULong64_t)runlastentry);
1945 //
1946 prevshift = 0;
1947 cout << " ))))) UPDATE RUN INFO ((((( @iev "<<iev<<" run "<<GetRunInfo()->ID<<" irun "<<irun<<endl;
1948 // cout << "runfirstentry "<<runfirstentry<<endl;
1949 return(true);
1950 };
1951 // ----------------------------------------------------
1952 // if it is a preselected file (there is SelectionList)
1953 // NBNB - the event tree MUST be read first
1954 // ----------------------------------------------------
1955 if(SELLI==1){
1956 sel_tree->GetEntry(iev);
1957 // cout << irun << " "<< irunentry << endl;
1958 if(irun != oldrun){
1959 run->GetEntry(irun);
1960 cout << " ))))) UPDATE RUN INFO ((((( @iev "<<iev<<" run "<<GetRunInfo()->ID<<" irun "<<irun<<endl;
1961 prevshift = 0;
1962 return true;
1963 }
1964 return false;
1965 }
1966
1967 return false;
1968 //
1969 };
1970
1971 Bool_t PamLevel2::UpdateRunInfo(Long64_t iev){
1972
1973 if(!run_tree){
1974 cout << " Bool_t PamLevel2::UpdateRunInfo(ULong64_t iev) -- ERROR -- run tree not loaded"<<endl;
1975 return false;
1976 }
1977 if ( run_tree->GetEntries() <= 0 ) {
1978 cout << " Bool_t PamLevel2::UpdateRunInfo(ULong64_t iev) -- ERROR -- run tree is empty"<<endl;
1979 return(false);
1980 }
1981
1982 Int_t oldrun = irun; // store current run index
1983
1984 // +-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-
1985 // if it is a full file (not preselected)
1986 // +-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-
1987 if(SELLI==0){
1988
1989 // ---------------------------------------------------------------
1990 // the absolute time is necessary to relate the event with the run
1991 // ---------------------------------------------------------------
1992 if( !GetOrbitalInfo() ){
1993 cout << "Bool_t PamLevel2::UpdateRunInfo(Long64_t iev) -- ERROR -- missing OrbitalInfo "<<endl;
1994 return(false);
1995 }
1996
1997 // -----------------------------------------------------------------------
1998 // the first time the routine is called, set run search from the beginning
1999 // -----------------------------------------------------------------------
2000 if ( irun < 0 ){
2001 irun = 0LL;
2002 irunentry = 0;
2003 prevshift = 0;
2004 run_tree->GetEntry(irun);
2005 };
2006 //
2007 bool fromfirst = true; // first loop over runs
2008
2009 // Bool_t hasfrag = false;
2010 // if( GetRunInfo()->ID_RUN_FRAG!=0 && GetRunInfo()->ID_RUN_FRAG != GetRunInfo()->ID)hasfrag=true;
2011 // ULong64_t fragid = GetRunInfo()->ID_RUN_FRAG;
2012
2013 // ------------------------------------------------------
2014 // loop over runs to find the one that contains the event
2015 // ------------------------------------------------------
2016 while (
2017 (
2018 (
2019 !(GetOrbitalInfo()->absTime >= GetRunInfo()->RUNHEADER_TIME && // check on absolute time (s)
2020 GetOrbitalInfo()->absTime <= GetRunInfo()->RUNTRAILER_TIME) &&
2021 !(GetOrbitalInfo()->OBT >= GetRunInfo()->RUNHEADER_OBT && // additional check on OBT (ms)
2022 GetOrbitalInfo()->OBT <= GetRunInfo()->RUNTRAILER_OBT)
2023 )
2024 || GetRunInfo()->NEVENTS==0
2025 || !(irunentry < GetRunInfo()->NEVENTS-1-prevshift)
2026 )
2027 && irun < run_tree->GetEntries() ){
2028
2029 irun++;
2030 // ------------------------------------
2031 // if the end of run tree is reached...
2032 // ------------------------------------
2033 if( irun == run_tree->GetEntries() ){
2034 if(!fromfirst){
2035 // -----------------------------------------------------
2036 // if it happened already once and the run was not found
2037 // ---> exit with error
2038 // -----------------------------------------------------
2039 cout << "Bool_t PamLevel2::UpdateRunInfo(Long64_t iev) -- ERROR -- event entry #"<<iev<<" does not belong to any run (should not happen)" <<endl;
2040 return false;
2041 }
2042 // -----------------------------------------
2043 // ...otherwise repeat search from beginning
2044 // -----------------------------------------
2045 cout << "Bool_t PamLevel2::UpdateRunInfo(Long64_t iev) -- WARNING -- reached end of run tree. searchin again from beginning " <<endl;
2046 fromfirst = false;
2047 irun = 0LL;
2048 runfirstentry = 0LL;
2049 }
2050 // -------------------------------------------------------------------
2051 // save the index of the first entry of the run, relative to pam_tree,
2052 // and read a new run
2053 // -------------------------------------------------------------------
2054 if(irun>0)runfirstentry += (GetRunInfo()->NEVENTS)-prevshift;
2055 irunentry = 0;
2056 prevshift = 0;
2057 run_tree->GetEntry(irun);
2058 if(GetRunInfo()->RUNHEADER_OBT>GetRunInfo()->RUNTRAILER_OBT ){
2059 cout << "Bool_t PamLevel2::UpdateRunInfo(Long64_t iev) -- WARNING -- irun "<<irun<<" has RUNHEADER_OBT>=RUNTRAILER_OBT " <<endl;
2060 cout << " (NB!! in this case some events are assigned to a wrong run)"<<endl;
2061 }
2062 // if(hasfrag && fragid != GetRunInfo()->ID){
2063 // cout << "... where is the next fragment ??"<<endl;
2064 // }
2065 };
2066
2067
2068 // --------------------------------------
2069 // if there was no need to update the run
2070 // ---> exit with FALSE
2071 // --------------------------------------
2072 if ( irun == oldrun ) return(false);
2073
2074 // --------------------------------------
2075 // ... otherwise
2076 // ---> exit with TRUE
2077 // --------------------------------------
2078 cout << endl << " ))))) UPDATE RUN INFO ((((( @iev "<<iev<<" run "<<GetRunInfo()->ID<<" irun "<<irun<<endl;
2079 // ----------------------------------------------------
2080 // then check if the run has a fragment
2081 // in this case we have to switch to the next fragment
2082 // when the end of the first fragment is reached
2083 // ----------------------------------------------------
2084 if(
2085 GetRunInfo()->ID_RUN_FRAG != 0 &&
2086 // GetRunInfo()->ID_RUN_FRAG != GetRunInfo()->ID &&
2087 true ){
2088 cout << "* fragment *"<<endl;
2089 }
2090
2091 return(true);
2092 };
2093 // +-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-
2094 // if it is a preselected file (there is SelectionList)
2095 // +-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-
2096 if(SELLI==1){
2097 sel_tree->GetEntry(iev);
2098 if(irun != oldrun){
2099 run_tree->GetEntry(irun);
2100 cout << endl << " ))))) UPDATE RUN INFO ((((( @iev "<<iev<<" run "<<GetRunInfo()->ID<<" irun "<<irun<<endl;
2101 prevshift = 0;
2102 return true;
2103 }
2104 return false;
2105 }
2106
2107 return false;
2108 //
2109 };
2110 /**
2111 * Update the runinfo informations (to be used to have Run infos event by event basis)
2112 * @param run Pointer to the chain/tree which contains run infos
2113 * @return true if a new run has been read, false if it is still the same run
2114 */
2115 Bool_t PamLevel2::UpdateRunInfo(TTree *run, Long64_t iev){
2116 return(UpdateRunInfo((TChain*)run,iev));
2117 };
2118
2119 //--------------------------------------
2120 //
2121 //
2122 //--------------------------------------
2123 /**
2124 * Set which trees shoul be analysed
2125 * @param detlist TString containing the sequence of trees required
2126 */
2127 void PamLevel2::SetWhichTrees(TString detlist){
2128
2129 if(detlist.IsNull() || detlist.Contains("+ALL", TString::kIgnoreCase)){
2130 CAL0 = false;
2131 CAL1 = true;
2132 CAL2 = true;
2133 TRK2 = true;
2134 TRK1 = false;
2135 TRKh = false;
2136 TRK0 = false;
2137 TRG = true;
2138 TOF = true;
2139 TOF0 = false;
2140 S4 = true;
2141 ND = true;
2142 AC = true;
2143 ORB = true;
2144 }else if( detlist.Contains("-ALL", TString::kIgnoreCase) ){
2145 CAL0 = false;
2146 CAL1 = false;
2147 CAL2 = false;
2148 TRK2 = false;
2149 TRK1 = false;
2150 TRKh = false;
2151 TRK0 = false;
2152 TRG = false;
2153 TOF = false;
2154 TOF0 = false;
2155 S4 = false;
2156 ND = false;
2157 AC = false;
2158 ORB = false;
2159 };
2160
2161 // -------------------------------------------------------------------------
2162 if( detlist.Contains("CAL1", TString::kIgnoreCase) ){
2163 if ( detlist.Contains("-CAL1", TString::kIgnoreCase) )CAL1=false;
2164 if ( detlist.Contains("+CAL1", TString::kIgnoreCase) )CAL1=true;
2165 };
2166
2167 if( detlist.Contains("CAL0", TString::kIgnoreCase) ){
2168 if ( detlist.Contains("-CAL0", TString::kIgnoreCase) )CAL0=false;
2169 if ( detlist.Contains("+CAL0", TString::kIgnoreCase) )CAL0=true;
2170 };
2171
2172 if( detlist.Contains("CAL2", TString::kIgnoreCase)){
2173 if ( detlist.Contains("-CAL2", TString::kIgnoreCase) )CAL2=false;
2174 if ( detlist.Contains("+CAL2", TString::kIgnoreCase) )CAL2=true;
2175 };
2176
2177 if( detlist.Contains("+CAL", TString::kIgnoreCase) && !CAL1 && !CAL2 )CAL2=true;
2178 if( detlist.Contains("-CAL", TString::kIgnoreCase) && CAL1 && CAL2 ){
2179 CAL2=false;
2180 CAL1=false;
2181 }
2182 // -------------------------------------------------------------------------
2183 if( detlist.Contains("TRK0", TString::kIgnoreCase) ){
2184 if ( detlist.Contains("-TRK0", TString::kIgnoreCase) )TRK0=false;
2185 if ( detlist.Contains("+TRK0", TString::kIgnoreCase) )TRK0=true;
2186 };
2187
2188 if( detlist.Contains("TRK1", TString::kIgnoreCase) ){
2189 if ( detlist.Contains("-TRK1", TString::kIgnoreCase) )TRK1=false;
2190 if ( detlist.Contains("+TRK1", TString::kIgnoreCase) )TRK1=true;
2191 };
2192
2193 if( detlist.Contains("TRK2", TString::kIgnoreCase)){
2194 if ( detlist.Contains("-TRK2", TString::kIgnoreCase) )TRK2=false;
2195 if ( detlist.Contains("+TRK2", TString::kIgnoreCase) )TRK2=true;
2196 };
2197
2198 if( detlist.Contains("TRKh", TString::kIgnoreCase)){
2199 if ( detlist.Contains("-TRKh", TString::kIgnoreCase) )TRKh=false;
2200 if ( detlist.Contains("+TRKh", TString::kIgnoreCase) )TRKh=true;
2201 };
2202
2203 if( detlist.Contains("+TRK", TString::kIgnoreCase) && !TRK1 && !TRK2 && !TRKh )TRK2=true;
2204 if( detlist.Contains("-TRK", TString::kIgnoreCase) && TRK1 && TRK2 && TRKh){
2205 TRK2=false;
2206 TRK1=false;
2207 TRKh=false;
2208 }
2209 // -------------------------------------------------------------------------
2210
2211 if( detlist.Contains("-TRG", TString::kIgnoreCase) )TRG = false;
2212 else if( detlist.Contains("+TRG", TString::kIgnoreCase) )TRG = true;
2213
2214 if( detlist.Contains("-TOF", TString::kIgnoreCase) )TOF = false;
2215 else if( detlist.Contains("+TOF", TString::kIgnoreCase) )TOF = true;
2216
2217 if( detlist.Contains("-TOF0", TString::kIgnoreCase) )TOF0 = false;
2218 else if( detlist.Contains("+TOF0", TString::kIgnoreCase) )TOF0 = true;
2219
2220 if( detlist.Contains("-S4", TString::kIgnoreCase) )S4 = false;
2221 else if( detlist.Contains("+S4", TString::kIgnoreCase) )S4 = true;
2222
2223 if( detlist.Contains("-ND", TString::kIgnoreCase) )ND = false;
2224 else if( detlist.Contains("+ND", TString::kIgnoreCase) )ND = true;
2225
2226 if( detlist.Contains("-AC", TString::kIgnoreCase) )AC = false;
2227 else if( detlist.Contains("+AC", TString::kIgnoreCase) )AC = true;
2228
2229 if( detlist.Contains("-ORB", TString::kIgnoreCase) )ORB = false;
2230 else if( detlist.Contains("+ORB", TString::kIgnoreCase) )ORB = true;
2231
2232 // cout<< "Set detector list --> ";
2233 // if(TRK1)cout<<"TRK1 ";
2234 // if(TRK2)cout<<"TRK2 ";
2235 // if(TRKh)cout<<"TRKH ";
2236 // if(CAL1)cout<<"CAL1 ";
2237 // if(CAL2)cout<<"CAL2 ";
2238 // if(TOF)cout<<"TOF ";
2239 // if(TRG)cout<<"TRG ";
2240 // if(AC)cout<<"AC ";
2241 // if(ND)cout<<"ND ";
2242 // if(S4)cout<<"S4 ";
2243 // if(ORB)cout<<"ORB ";
2244 // cout << endl;
2245
2246 };
2247
2248
2249 /**
2250 * Set tree/branch detector flags from the content of a tree
2251 */
2252 void PamLevel2::GetWhichTrees(TFile* f){
2253
2254
2255
2256 cout << "void PamLevel2::GetWhichTrees(TFile* f) --- WARNING!! --- ...potrebbe non funzionare "<<endl;
2257 // -----------
2258 // reset flags
2259 // -----------
2260 CAL1 = false;
2261 CAL2 = false;
2262 TRK2 = false;
2263 TRK1 = false;
2264 TRKh = false;
2265 TRG = false;
2266 TOF = false;
2267 S4 = false;
2268 ND = false;
2269 AC = false;
2270 ORB = false;
2271
2272 RUN = false;
2273
2274 cout << "Checking file: "<<f->GetName()<<endl;
2275 if( !f || f->IsZombie() ){
2276 cout << "File: "<< f->GetName() <<" Non valid root file"<< endl;
2277 return;
2278 }
2279
2280 TList *lk = f->GetListOfKeys();
2281 TIter next(lk);
2282 TKey *key =0;
2283
2284 Int_t nev = 0;
2285
2286 while( (key = (TKey*)next()) ){
2287
2288 if( !strcmp(key->GetName(),"Run" ) )RUN = true;
2289
2290 //=========================================================
2291 if( !strcmp(key->GetName(),"Trigger" ) ){
2292 TRG = true;
2293 Int_t nevt = ((TTree*)f->Get("Trigger"))->GetEntries();
2294 if( nev && nevt!=nev){
2295 cout << "File: "<< f->GetName() <<" Trigger tree has "<<nevt<<" events instead of "<<nev<< endl;
2296 TRG = false;
2297 }else nev=nevt;
2298 }
2299 //=========================================================
2300 if( !strcmp(key->GetName(),"ToF" ) ){
2301 TOF = true;
2302 Int_t nevt = ((TTree*)f->Get("ToF"))->GetEntries();
2303 if( nev && nevt!=nev){
2304 cout << "File: "<< f->GetName() <<" ToF tree has "<<nevt<<" events instead of "<<nev<< endl;
2305 TOF = false;
2306 }else nev=nevt;
2307 }
2308 //=========================================================
2309 if( !strcmp(key->GetName(),"S4" ) ){
2310 S4 = true;
2311 Int_t nevt = ((TTree*)f->Get("S4"))->GetEntries();
2312 if( nev && nevt!=nev){
2313 cout << "File: "<< f->GetName() <<" S4 tree has "<<nevt<<" events instead of "<<nev<< endl;
2314 S4 = false;
2315 }else nev=nevt;
2316 }
2317 //=========================================================
2318
2319 if( !strcmp(key->GetName(),"NeutronD" ) ){
2320 ND = true;
2321 Int_t nevt = ((TTree*)f->Get("NeutronD"))->GetEntries();
2322 if( nev && nevt!=nev){
2323 cout << "File: "<< f->GetName() <<"NeutronD tree has "<<nevt<<" events instead of "<<nev<< endl;
2324 ND =false;
2325 }else nev=nevt;
2326 }
2327 //=========================================================
2328 if( !strcmp(key->GetName(),"Anticounter") ){
2329 AC = true;
2330 Int_t nevt = ((TTree*)f->Get("Anticounter"))->GetEntries();
2331 if( nev && nevt!=nev){
2332 cout << "File: "<< f->GetName() <<" Anticounter tree has "<<nevt<<" events instead of "<<nev<< endl;
2333 AC =false;
2334 }else nev=nevt;
2335 }
2336 //=========================================================
2337 if( !strcmp(key->GetName(),"OrbitalInfo") ){
2338 ORB = true;
2339 Int_t nevt = ((TTree*)f->Get("OrbitalInfo"))->GetEntries();
2340 if( nev && nevt!=nev){
2341 cout << "File: "<< f->GetName() <<" OrbitalInfo tree has "<<nevt<<" events instead of "<<nev<< endl;
2342 ORB = false;
2343 }else nev=nevt;
2344 }
2345 //=========================================================
2346 if( !strcmp(key->GetName(),"Tracker" ) ){
2347 TTree *T = (TTree*)f->Get("Tracker");
2348 for(Int_t i=0; i<T->GetListOfBranches()->GetEntries(); i++){
2349 TString name = T->GetListOfBranches()->At(i)->GetName();
2350 if( !name.CompareTo("TrkLevel1") )TRK1=true;
2351 if( !name.CompareTo("TrkLevel2") )TRK2=true;
2352 if( !name.CompareTo("TrkHough") )TRKh=true;
2353 };
2354 Int_t nevt = T->GetEntries();
2355 if( nev && nevt!=nev){
2356 cout << "File: "<< f->GetName() <<" Tracker tree has "<<nevt<<" events instead of "<<nev<< endl;
2357 TRK1 = false;
2358 TRK2 = false;
2359 TRKh = false;
2360 }else nev=nevt;
2361 T->Delete();
2362 };
2363 //=========================================================
2364 if( !strcmp(key->GetName(),"Calorimeter" ) ){
2365 TTree *T = (TTree*)f->Get("Calorimeter");
2366 for(Int_t i=0; i<T->GetListOfBranches()->GetEntries(); i++){
2367 TString name = T->GetListOfBranches()->At(i)->GetName();
2368 if( !name.CompareTo("CaloLevel1") )CAL1=true;
2369 if( !name.CompareTo("CaloLevel2") )CAL2=true;
2370 };
2371 Int_t nevt = T->GetEntries();
2372 if( nev && nevt!=nev){
2373 cout << "File: "<< f->GetName() <<" Calorimeter tree has "<<nevt<<" events instead of "<<nev<< endl;
2374 CAL1 = false;
2375 CAL2 = false;
2376 }else nev=nevt;
2377 T->Delete();
2378 };
2379
2380 };
2381
2382 delete lk;
2383
2384 // cout<< "Get detector list from input file --> ";
2385 // if(TRK1)cout<<"TRK1 ";
2386 // if(TRK2)cout<<"TRK2 ";
2387 // if(TRKh)cout<<"TRKH ";
2388 // if(CAL1)cout<<"CAL1 ";
2389 // if(CAL2)cout<<"CAL2 ";
2390 // if(TOF)cout<<"TOF ";
2391 // if(TRG)cout<<"TRG ";
2392 // if(AC)cout<<"AC ";
2393 // if(ND)cout<<"ND ";
2394 // if(S4)cout<<"S4 ";
2395 // if(ORB)cout<<"ORB ";
2396 // cout << endl;
2397
2398 return ;
2399
2400 };
2401
2402
2403 //--------------------------------------
2404 //
2405 //
2406 //--------------------------------------
2407 /**
2408 * Check if a file contains selected Pamela Level2 trees.
2409 * @param name File name
2410 * @return true if the file is ok.
2411 */
2412 Bool_t PamLevel2::CheckLevel2File(TString name){
2413
2414 Bool_t CAL1__ok = false;
2415 Bool_t CAL2__ok = false;
2416 Bool_t TRK2__ok = false;
2417 Bool_t TRK1__ok = false;
2418 Bool_t TRKh__ok = false;
2419 Bool_t TRG__ok = false;
2420 Bool_t TOF__ok = false;
2421 Bool_t S4__ok = false;
2422 Bool_t ND__ok = false;
2423 Bool_t AC__ok = false;
2424 Bool_t ORB__ok = false;
2425
2426 Bool_t RUN__ok = false;
2427
2428 Bool_t SELLI__ok = false;
2429
2430 cout << "Checking file: "<<name<<endl;
2431 TFile *f = new TFile(name.Data());
2432 if( !f || f->IsZombie() ){
2433 cout << "File: "<< f->GetName() <<" discarded ---- Non valid root file"<< endl; return false;
2434 }
2435 // cout << "Get list of keys: "<<f<<endl;
2436 TList *lk = f->GetListOfKeys();
2437 // lk->Print();
2438 TIter next(lk);
2439 TKey *key =0;
2440
2441 Int_t nev = 0;
2442
2443 while( (key = (TKey*)next()) ){
2444
2445 // cout << key->GetName() << endl;
2446 // cout << key->GetName() << ""<<key->GetClassName()<<endl;
2447 // cout << " Get tree: " << f->Get(key->GetName())<<endl;
2448 // nev_previous = nev;
2449 // cout << " n.entries "<< nev <<endl;
2450 // if( key->GetClassName()=="TTree" && nev_previous && nev != nev_previous ){
2451 // nev = ((TTree*)f->Get(key->GetName()))->GetEntries();
2452 // cout << "File: "<< f->GetName() <<" discarded ---- "<< key->GetName() << " tree: n.entries does not match "<<nev<<" "<<nev_previous<< endl;
2453 // return false;
2454 // };
2455
2456
2457 if( !strcmp(key->GetName(),"Run" ) )RUN__ok = true;
2458
2459 //=========================================================
2460 if( !strcmp(key->GetName(),"SelectionList" ) ){
2461 SELLI__ok = true;
2462 if(SELLI==1){
2463 Int_t nevt = ((TTree*)f->Get("SelectionList"))->GetEntries();
2464 if( nev && nevt!=nev){
2465 cout << "File: "<< f->GetName() <<" discarded ---- SelectionList tree has "<<nevt<<" events instead of "<<nev<< endl;
2466 return false;
2467 }
2468 nev=nevt;
2469 }
2470 }
2471
2472 //=========================================================
2473 if( !strcmp(key->GetName(),"Trigger" ) ){
2474 TRG__ok = true;
2475 if(TRG){
2476 Int_t nevt = ((TTree*)f->Get("Trigger"))->GetEntries();
2477 if( nev && nevt!=nev){
2478 cout << "File: "<< f->GetName() <<" discarded ---- Trigger tree has "<<nevt<<" events instead of "<<nev<< endl;
2479 return false;
2480 }
2481 nev=nevt;
2482 }
2483 }
2484 //=========================================================
2485 if( !strcmp(key->GetName(),"ToF" ) ){
2486 TOF__ok = true;
2487 if(TOF){
2488 Int_t nevt = ((TTree*)f->Get("ToF"))->GetEntries();
2489 if( nev && nevt!=nev){
2490 cout << "File: "<< f->GetName() <<" discarded ---- ToF tree has "<<nevt<<" events instead of "<<nev<< endl;
2491 return false;
2492 }
2493 nev=nevt;
2494 }
2495 }
2496 //=========================================================
2497 if( !strcmp(key->GetName(),"S4" ) ){
2498 S4__ok = true;
2499 if(S4){
2500 Int_t nevt = ((TTree*)f->Get("S4"))->GetEntries();
2501 if( nev && nevt!=nev){
2502 cout << "File: "<< f->GetName() <<" discarded ---- S4 tree has "<<nevt<<" events instead of "<<nev<< endl;
2503 return false;
2504 }
2505 nev=nevt;
2506 }
2507 }
2508 //=========================================================
2509
2510 if( !strcmp(key->GetName(),"NeutronD" ) ){
2511 ND__ok = true;
2512 if(ND){
2513 Int_t nevt = ((TTree*)f->Get("NeutronD"))->GetEntries();
2514 if( nev && nevt!=nev){
2515 cout << "File: "<< f->GetName() <<" discarded ---- NeutronD tree has "<<nevt<<" events instead of "<<nev<< endl;
2516 return false;
2517 }
2518 nev=nevt;
2519 }
2520 }
2521 //=========================================================
2522 if( !strcmp(key->GetName(),"Anticounter") ){
2523 AC__ok = true;
2524 if(AC){
2525 Int_t nevt = ((TTree*)f->Get("Anticounter"))->GetEntries();
2526 if( nev && nevt!=nev){
2527 cout << "File: "<< f->GetName() <<" discarded ---- Anticounter tree has "<<nevt<<" events instead of "<<nev<< endl;
2528 return false;
2529 }
2530 nev=nevt;
2531 }
2532 }
2533 //=========================================================
2534 if( !strcmp(key->GetName(),"OrbitalInfo") ){
2535 ORB__ok = true;
2536 if(ORB){
2537 Int_t nevt = ((TTree*)f->Get("OrbitalInfo"))->GetEntries();
2538 if( nev && nevt!=nev){
2539 cout << "File: "<< f->GetName() <<" discarded ---- OrbitalInfo tree has "<<nevt<<" events instead of "<<nev<< endl;
2540 return false;
2541 }
2542 nev=nevt;
2543 }
2544 }
2545 //=========================================================
2546 if( !strcmp(key->GetName(),"Tracker" ) ){
2547 TTree *T = (TTree*)f->Get("Tracker");
2548 if(TRK1||TRK2||TRKh){
2549 Int_t nevt = T->GetEntries();
2550 if( nev && nevt!=nev){
2551 cout << "File: "<< f->GetName() <<" discarded ---- Tracker tree has "<<nevt<<" events instead of "<<nev<< endl;
2552 return false;
2553 }
2554 nev=nevt;
2555 }
2556 for(Int_t i=0; i<T->GetListOfBranches()->GetEntries(); i++){
2557 TString name = T->GetListOfBranches()->At(i)->GetName();
2558 if( !name.CompareTo("TrkLevel1") )TRK1__ok=true;
2559 if( !name.CompareTo("TrkLevel2") )TRK2__ok=true;
2560 if( !name.CompareTo("TrkHough") )TRKh__ok=true;
2561 };
2562 T->Delete();
2563 };
2564 //=========================================================
2565 if( !strcmp(key->GetName(),"Calorimeter" ) ){
2566 TTree *T = (TTree*)f->Get("Calorimeter");
2567 if(CAL1||CAL2){
2568 Int_t nevt = T->GetEntries();
2569 if( nev && nevt!=nev){
2570 cout << "File: "<< f->GetName() <<" discarded ---- Calorimeter tree has "<<nevt<<" events instead of "<<nev<< endl;
2571 return false;
2572 }
2573 nev=nevt;
2574 }
2575 for(Int_t i=0; i<T->GetListOfBranches()->GetEntries(); i++){
2576 TString name = T->GetListOfBranches()->At(i)->GetName();
2577 if( !name.CompareTo("CaloLevel1") )CAL1__ok=true;
2578 if( !name.CompareTo("CaloLevel2") )CAL2__ok=true;
2579 };
2580 T->Delete();
2581 };
2582
2583 };
2584
2585 if( SELLI==-1 )SELLI = (Int_t)SELLI__ok;
2586 if( SELLI==0 && SELLI__ok ){
2587 cout << "File: "<< f->GetName() <<" discarded ---- found SelectionList (it is not a full-event file)" << endl;
2588 return false;
2589 }
2590 if( SELLI==1 && !SELLI__ok ){
2591 cout << "File: "<< f->GetName() <<" discarded ---- SelectionList missing" << endl;
2592 return false;
2593 }
2594
2595 // cout << "SELLI "<<SELLI<<endl;
2596
2597 // cout<< "CheckLevel2File(TString): detector list --> ";
2598 // if(TRK1__ok)cout<<"TRK1 ";
2599 // if(TRK2__ok)cout<<"TRK2 ";
2600 // if(TRKh__ok)cout<<"TRKH ";
2601 // if(CAL1__ok)cout<<"CAL1 ";
2602 // if(CAL2__ok)cout<<"CAL2 ";
2603 // if(TOF__ok)cout<<"TOF ";
2604 // if(TRG__ok)cout<<"TRG ";
2605 // if(AC__ok)cout<<"AC ";
2606 // if(ND__ok)cout<<"ND ";
2607 // if(S4__ok)cout<<"S4 ";
2608 // if(ORB__ok)cout<<"ORB ";
2609 // cout << endl;
2610
2611
2612 if(TRK2 && TRK1__ok)TRK1=1;
2613 // ----------------------------------------------------------------------------
2614 // NOTA
2615 // se c'e` il level1, lo devo necessarimente leggere.
2616 // infatti (non ho capito perche`) i cluster vengono letti e allocati in memoria
2617 // comunque, ma non vengono disallocati da PamLevel2::Clear()
2618 // ----------------------------------------------------------------------------
2619
2620
2621 if(!RUN__ok) {
2622 cout << "File: "<< f->GetName() <<" *WARNING* ---- Missing RunInfo tree"<< endl;
2623 // return false;
2624 };
2625
2626 if(CAL1 && !CAL1__ok){
2627 cout << "File: "<< f->GetName() <<" discarded ---- Missing CaloLevel1 branch"<< endl;
2628 return false;
2629 };
2630 if(CAL2 && !CAL2__ok){
2631 cout << "File: "<< f->GetName() <<" discarded ---- Missing CaloLevel2 branch"<< endl;
2632 return false;
2633 };
2634 if(TRK2 && !TRK2__ok){
2635 cout << "File: "<< f->GetName() <<" discarded ---- Missing TrkLevel2 branch"<< endl;
2636 return false;
2637 };
2638 if(TRK1 && !TRK1__ok){
2639 cout << "File: "<< f->GetName() <<" discarded ---- Missing TrkLevel1 branch"<< endl;
2640 return false;
2641 };
2642 if(TRKh && !TRKh__ok){
2643 cout << "File: "<< f->GetName() <<" discarded ---- Missing TrkHough branch"<< endl;
2644 return false;
2645 };
2646 if(ORB && !ORB__ok){
2647 cout << "File: "<< f->GetName() <<" discarded ---- Missing ORB tree"<< endl;
2648 return false;
2649 };
2650 if(AC && !AC__ok){
2651 cout << "File: "<< f->GetName() <<" discarded ---- Missing AC tree"<< endl;
2652 return false;
2653 };
2654 if(S4 && !S4__ok){
2655 cout << "File: "<< f->GetName() <<" discarded ---- Missing S4 tree"<< endl;
2656 return false;
2657 };
2658 if(TOF && !TOF__ok){
2659 cout << "File: "<< f->GetName() <<" discarded ---- Missing ToF tree"<< endl;
2660 return false;
2661 };
2662
2663 if(ND && !ND__ok){
2664 cout << "File: "<< f->GetName() <<" discarded ---- Missing ND tree"<< endl;
2665 return false;
2666 };
2667 if(TRG && !TRG__ok){
2668 cout << "File: "<< f->GetName() <<" discarded ---- Missing Trigger tree"<< endl;
2669 return false;
2670 };
2671
2672
2673 // lk->Delete();
2674 // delete lk;
2675 f->Close();
2676
2677 // cout<< "CheckLevel2File(TString): detector list --> ";
2678 // if(TRK1)cout<<"TRK1 ";
2679 // if(TRK2)cout<<"TRK2 ";
2680 // if(TRKh)cout<<"TRKH ";
2681 // if(CAL1)cout<<"CAL1 ";
2682 // if(CAL2)cout<<"CAL2 ";
2683 // if(TOF)cout<<"TOF ";
2684 // if(TRG)cout<<"TRG ";
2685 // if(AC)cout<<"AC ";
2686 // if(ND)cout<<"ND ";
2687 // if(S4)cout<<"S4 ";
2688 // if(ORB)cout<<"ORB ";
2689 // cout << endl;
2690
2691 return true;
2692
2693 };
2694
2695
2696 /**
2697 * Create clone-trees
2698 */
2699 void PamLevel2::CreateCloneTrees0( TChain *fChain, TFile *ofile ){
2700
2701 cout << "+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+" <<endl;
2702 cout << "Create clones of PAMELA trees "<<endl;
2703
2704 Int_t i=0;
2705 pam_tree_clone[i] = fChain->GetTree()->CloneTree(0);
2706 TString name = pam_tree_clone[i]->GetName();
2707 name.Append("_clone");
2708 // pam_tree_clone[i]->SetName(name.Data());
2709 cout << pam_tree_clone[i]->GetName() <<endl;
2710 i++;
2711
2712 TList *li = fChain->GetListOfFriends();
2713 TIter next(li);
2714 TFriendElement* T_friend=0;
2715 ofile->cd();
2716 while( (T_friend = (TFriendElement*)next()) ){
2717 // cout<<T_friend->IsA()->GetName()<<" "<<T_friend->GetName()<<hex << T_friend->GetTree() << dec<<endl;
2718 // cout<<T_friend->GetTree()->GetName()<< endl;
2719 pam_tree_clone[i] = T_friend->GetTree()->CloneTree(0);
2720 pam_tree_clone[i]->SetAutoSave(1000000);
2721 name = pam_tree_clone[i]->GetName();
2722 name.Append("_clone");
2723 // pam_tree_clone[i]->SetName(name.Data());
2724 cout << pam_tree_clone[i]->GetName() << endl;
2725 i++;
2726 }
2727
2728 delete li;
2729
2730 cout << "+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+" <<endl;
2731
2732 }
2733
2734 /**
2735 * Create clone-trees
2736 */
2737 void PamLevel2::CreateCloneTrees(TFile *ofile){
2738
2739
2740 // if the pointer is null, create a default file
2741
2742 if(!ofile){
2743 cout << "void PamLevel2::CreateCloneTrees(TFile*) -- WARNING -- Creating file: clone-tree.root "<<endl;
2744 ofile = new TFile("clone-tree.root","recreate");
2745 }
2746
2747 ofile->cd();
2748
2749 cout << "+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+" <<endl;
2750 cout << "Create new PAMELA trees "<<endl;
2751
2752
2753 run_tree_clone = new TTree("Run","PAMELA Level2 data from the GL_RUN table ");
2754 run_tree_clone->Branch("RunInfo","GL_RUN",GetPointerTo("RunInfo"));
2755 cout << "Run : branch RunInfo"<<endl;
2756 run_tree_clone->Branch("SoftInfo","SoftInfo",GetPointerTo("SoftInfo"));
2757 cout << "Run : branch SoftInfo"<<endl;
2758 // ------------------
2759 // replicate run tree
2760 // ------------------
2761 cout << "----------------------------------------------------"<<endl;
2762 cout << "irun\t | RUN\t NEVENTS\t absolute time"<<endl;
2763 for (Int_t i=0; i<run_tree->GetEntries(); i++){
2764 run_tree->GetEntry(i);
2765 cout << i<< "\t | "<<GetRunInfo()->ID<<"\t "<<GetRunInfo()->NEVENTS<< "\t "<<GetRunInfo()->RUNHEADER_TIME<<" <---> "<<GetRunInfo()->RUNTRAILER_TIME<<endl;
2766 run_tree_clone->Fill();
2767 }
2768 cout << "----------------------------------------------------"<<endl;
2769
2770
2771 sel_tree_clone = new TTree("SelectionList","List of selected events ");
2772 sel_tree_clone->Branch("RunEntry",&irun,"runentry/L");
2773 sel_tree_clone->Branch("EventEntry",&irunentry,"eventry/L");
2774
2775
2776 Int_t i=0;
2777 if(TRK1||TRK2||TRKh){
2778 pam_tree_clone[i] = new TTree("Tracker","PAMELA tracker level2 data ");
2779 if(TRK1) {
2780 pam_tree_clone[i]->Branch("TrkLevel1","TrkLevel1", GetPointerTo("TrkLevel1"));
2781 pam_tree_clone[i]->BranchRef();
2782 cout << "Tracker : branch TrkLevel1"<<endl;
2783 // cout << "CreateCloneTrees " << GetTrkLevel1()<<endl;
2784 };
2785 if(TRK2) {
2786 pam_tree_clone[i]->Branch("TrkLevel2", "TrkLevel2",GetPointerTo("TrkLevel2"));
2787 cout << "Tracker : branch TrkLevel2"<<endl;
2788 };
2789 if(TRKh) {
2790 pam_tree_clone[i]->Branch("TrkHough","TrkHough", GetPointerTo("TrkHough"));
2791 cout << "Tracker : branch TrkHough"<<endl;
2792 };
2793 i++;
2794 }
2795
2796 // Calorimeter
2797 if(CAL1||CAL2){
2798 pam_tree_clone[i] = new TTree("Calorimeter","PAMELA calorimeter level2 data ");
2799 if(CAL1) {
2800 pam_tree_clone[i]->Branch("CaloLevel1", "CaloLevel1", GetPointerTo("CaloLevel1"));
2801 cout << "Calorimeter : branch CaloLevel1"<<endl;
2802 };
2803 if(CAL2) {
2804 pam_tree_clone[i]->Branch("CaloLevel2","CaloLevel2", GetPointerTo("CaloLevel2"));
2805 cout << "Calorimeter : branch CaloLevel2"<<endl;
2806 };
2807 i++;
2808 }
2809
2810 // ToF
2811 if(TOF) {
2812 pam_tree_clone[i] = new TTree("ToF","PAMELA ToF level2 data ");
2813 pam_tree_clone[i]->Branch("ToFLevel2","ToFLevel2", GetPointerTo("ToFLevel2"));
2814 cout << "ToF : branch ToFLevel2"<<endl;
2815 i++;
2816 };
2817 // Trigger
2818 if(TRG) {
2819 pam_tree_clone[i] = new TTree("Trigger","PAMELA trigger level2 data ");
2820 pam_tree_clone[i]->Branch("TrigLevel2","TrigLevel2", GetPointerTo("TrigLevel2"));
2821 cout << "Trigger : branch TrigLevel2"<<endl;
2822 i++;
2823 };
2824 // S4
2825 if(S4) {
2826 pam_tree_clone[i] = new TTree("S4","PAMELA S4 level2 data ");
2827 pam_tree_clone[i]->Branch("S4Level2","S4Level2", GetPointerTo("S4Level2"));
2828 cout << "S4 : branch S4Level2"<<endl;
2829 i++;
2830 };
2831 // Neutron Detector
2832 if(ND) {
2833 pam_tree_clone[i] = new TTree("NeutronD","PAMELA neutron detector level2 data ");
2834 pam_tree_clone[i]->Branch("NDLevel2","NDLevel2", GetPointerTo("NDLevel2"));
2835 cout << "NeutronD : branch NDLevel2"<<endl;
2836 i++;
2837 };
2838 // Anticounters
2839 if(AC) {
2840 pam_tree_clone[i] = new TTree("Anticounter","PAMELA anticounter detector level2 data ");
2841 pam_tree_clone[i]->Branch("AcLevel2","AcLevel2", GetPointerTo("AcLevel2"));
2842 cout << "Anticounter : branch AcLevel2"<<endl;
2843 i++;
2844 };
2845 // OrbitalInfo
2846 if(ORB) {
2847 pam_tree_clone[i] = new TTree("OrbitalInfo","PAMELA oribital info ");
2848 pam_tree_clone[i]->Branch("OrbitalInfo","OrbitalInfo", GetPointerTo("OrbitalInfo"));
2849 cout << "OrbitalInfo : branch OrbitalInfo"<<endl;
2850 i++;
2851 };
2852 cout << "+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+" <<endl;
2853
2854 }
2855
2856
2857 /**
2858 * Fill tree (created with CreateCloneTrees)
2859 *
2860 */
2861 //void PamLevel2::FillNewPamTree(TTree *T){
2862 void PamLevel2::FillCloneTrees(){
2863
2864 // cout << "PamLevel2::FillCloneTrees()" << irunentry << endl;
2865
2866 for(Int_t i=0; i<NCLONES; i++){
2867 if(pam_tree_clone[i])pam_tree_clone[i]->Fill();
2868 }
2869 if(sel_tree_clone)sel_tree_clone->Fill();
2870
2871 }
2872
2873
2874 TTree* PamLevel2::GetCloneTree(TString name){
2875
2876 for(Int_t i=0; i<NCLONES; i++){
2877 if(pam_tree_clone[i]){
2878 TString na = pam_tree_clone[i]->GetName();
2879 if(!name.CompareTo(na))return pam_tree_clone[i];
2880 };
2881 }
2882 if(run_tree_clone){
2883 TString na = run_tree_clone->GetName();
2884 if(!name.CompareTo(na))return run_tree_clone;
2885 }
2886 if(sel_tree_clone){
2887 TString na = sel_tree_clone->GetName();
2888 if(!name.CompareTo(na))return sel_tree_clone;
2889 }
2890 return NULL;
2891
2892 }
2893 void PamLevel2::WriteCloneTrees(){
2894 cout << "+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+" <<endl;
2895 cout << "Write clones of PAMELA trees "<<endl;
2896 cout << run_tree_clone->GetName()<<endl;
2897 run_tree_clone->Write();
2898 cout << sel_tree_clone->GetName()<<endl;
2899 sel_tree_clone->Write();
2900 for(Int_t i=0; i<NCLONES; i++){
2901 if(pam_tree_clone[i]){
2902 cout << pam_tree_clone[i]->GetName()<<endl;
2903 pam_tree_clone[i]->Write();
2904 };
2905 }
2906 cout << "+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+" <<endl;
2907
2908 }
2909
2910 /**
2911 * Method to get level2-trees entry, the corresponding run entry and (if required) the level0 entry.
2912 */
2913 //Int_t PamLevel2::GetEntry(Int_t iee){
2914 Int_t PamLevel2::GetEntry(Long64_t iee){
2915
2916 // cout << "-------------------------------------"<<endl;
2917 // cout << "Int_t PamLevel2::GetEntry("<<iee<<")"<<endl;
2918
2919 if(!pam_tree){
2920 cout << " Int_t PamLevel2::GetEntry(Int_t) -- ERROR -- level2 trees not loaded"<<endl;
2921 return 0;
2922 }
2923
2924
2925 //
2926 // This is a sort of bug: if you don't have the run tree you don't want to exit here you want to have loaded the event anyway...
2927 //
2928 // if(!run_tree ){
2929 // cout << " Int_t PamLevel2::GetEntry(Int_t) -- ERROR -- run tree not loeaded"<<endl;
2930 // return 0;
2931 // }
2932
2933 Long64_t ii=0;
2934 //-------------------------------
2935 ii = iee;
2936 if( !pam_tree->GetEntry(ii) ){
2937 cout << " Int_t PamLevel2::GetEntry(Int_t) -- ERROR -- error reading pam tree"<<endl;
2938 return 0;
2939 }
2940 //
2941 // ... that's way I put it here. Notice that nothing change in the code (is backward compatible) since in any case you return with 0.
2942 // in theory one would like to return 1 if run is not loaded but now I don't have the will to add that 2 lines of code and it is not
2943 // a problem if you don't check the return code of getentry.
2944 //
2945 if(!run_tree ){
2946 if ( TRK0 || CAL0 || TOF0 || RUN ) { //forse cosi` va bene per tornare 1?
2947 cout << " Int_t PamLevel2::GetEntry(Int_t) -- ERROR -- run tree not loaded"<<endl;
2948 return 0;
2949 } else {
2950 return 1; //cosi` se non c'e` run esce qua...
2951 }
2952 }
2953
2954 //-------------------------------
2955 ii = iee;
2956 // Bool_t UPDATED = UpdateRunInfo(run_tree,ii);
2957 // Bool_t UPDATED = UpdateRunInfo(ii);
2958 UpdateRunInfo(ii);
2959 if(SELLI==0)irunentry = iee-runfirstentry;
2960 // if(UPDATED && run_tree_clone)run_tree_clone->Fill();
2961
2962 // cout << "PamLevel2::GetEntry("<<iee<<") "<<irun<<" "<<runfirstentry<<" "<<irunentry<<endl;
2963
2964 // cout << "irunentry "<<irunentry << endl;
2965 // cout << "runfirstentry "<<runfirstentry << endl;
2966 // cout << "nevents "<<GetRunInfo()->NEVENTS<< endl;
2967 if( TRK0 || CAL0 || TOF0 ){
2968 if( !GetYodaEntry( ) ){
2969 cout << " Int_t PamLevel2::GetEntry(Int_t) -- ERROR -- error reading level0 tree"<<endl;
2970 return 0;
2971 }
2972 }
2973 return 1;
2974
2975 }
2976
2977 /**
2978 * Method to retrieve the level0 tree (YODA tree) that contains the current event.
2979 * Given the run ID (...), if needed it query the DB and load the proper file.
2980 * @return Pointer to the tree
2981 */
2982
2983
2984 TTree* PamLevel2::GetYodaTree( ){
2985
2986 // cout << "TTree* PamLevel2::GetYodaTree( )"<<endl;
2987 //===================================
2988 // check if iroot has changed
2989 //===================================
2990 if( irun<0 ){
2991 cout << "PamLevel2::GetYodaTree() -- ERROR -- irun = "<<irun<<endl;
2992 // cout << "In order to use this method you have to first load the RunInfo tree "<<endl;
2993 return NULL;
2994 }
2995 Int_t irootnew = run_obj->ID_ROOT_L0;
2996 // cout << "iroot "<<iroot<<endl;
2997 // cout << "irootnew "<<irootnew<<endl;
2998
2999 //===================================
3000 // load the level0 file
3001 // (if not already loaded)
3002 //===================================
3003 if( iroot != irootnew || !l0_tree){
3004 iroot = irootnew;
3005 //===================================
3006 // open the DB connection
3007 // (if not already opened)
3008 //===================================
3009 if(!dbc || (dbc && !dbc->IsConnected())){
3010 cout<<"~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~"<<endl;
3011 cout<<"Connecting to DB"<<endl;
3012 cout<<"HOST "<<host<<endl;
3013 cout<<"USER "<<user<<endl;
3014 cout<<"PSW "<<psw<<endl;
3015 dbc = TSQLServer::Connect(host.Data(),user.Data(),psw.Data());
3016 if( !dbc )return NULL;
3017 if( !dbc->IsConnected() )return NULL;
3018 // cout<<"...done"<<endl;
3019 cout<<"~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~"<<endl;
3020 }else{
3021 // cout<<"DB already connected"<<endl;
3022 }
3023 GL_ROOT glroot = GL_ROOT();
3024 if( glroot.Query_GL_ROOT(iroot,dbc) ){
3025 cout << "TTree* PamLevel2::GetYodaTree( ) -- ERROR -- level0 file iroot = "<<iroot<< " does not exists"<<endl;
3026 return NULL;
3027 };
3028 TString filename = glroot.PATH + glroot.NAME;
3029 if(l0_file){
3030 l0_file->Close();
3031 l0_file->Delete();
3032 }
3033 cout << "Opening LEVEL0 file: "<< filename << endl;
3034 FileStat_t t;
3035 if( gSystem->GetPathInfo(filename.Data(),t) ){
3036 cout << " PamLevel2::GetYodaTree() -- ERROR opening file "<<endl;
3037 return NULL;
3038 }
3039 l0_file = new TFile(filename);
3040 if( !l0_file )return NULL;
3041 l0_tree = (TTree*)l0_file->Get("Physics");
3042 if(!h0_obj)h0_obj = new EventHeader();
3043 l0_tree->SetBranchAddress("Header" ,&h0_obj);
3044 prevshift = 0;
3045 //---------------------------------------------------
3046 // TRACKER:
3047 if(TRK0){
3048 if(!trk0_obj){
3049 trk0_obj = new TrkLevel0();
3050 trk0_obj->Set();
3051 };
3052 l0_tree->SetBranchAddress("Tracker" ,trk0_obj->GetPointerToTrackerEvent());
3053 }
3054 //---------------------------------------------------
3055 // CALORIMETER:
3056 if(CAL0){
3057 if(!calo0_obj){
3058 calo0_obj = new CaloLevel0();
3059 calo0_obj->Set();
3060 };
3061 l0_tree->SetBranchAddress("Calorimeter" ,calo0_obj->GetPointerToCalorimeterEvent());
3062 // cout << "PamLevel2::GetYodaTree() --- level0 calorimeter not implemented "<<endl;
3063 }
3064 //---------------------------------------------------
3065 // TOF:
3066 if(TOF0){
3067 cout << "PamLevel2::GetYodaTree() --- level0 TOF not implemented "<<endl;
3068 }
3069
3070 };
3071
3072 if(!dbc || (dbc && !dbc->IsConnected())){
3073 cout << " TTree* PamLevel2::GetYodaTree( ) -- no DB connected... hai fatto qualche cazzata "<<endl;
3074 }
3075
3076 if ( TRK0 ){
3077 TrkParams::Load(6);
3078 if( !TrkParams::IsLoaded(6) ){
3079 cout << " TTree* PamLevel2::GetYodaTree( ) -- WARNING -- VK-mask not loaded"<<endl;
3080 };
3081 TrkParams::SetCalib(run_obj,dbc);
3082 TrkParams::LoadCalib( );
3083 if( !TrkParams::CalibIsLoaded() ){
3084 cout << " TTree* PamLevel2::GetYodaTree( ) -- WARNING -- Calibration not loaded"<<endl;
3085 };
3086 }
3087
3088 // cout << l0_tree << endl;
3089
3090 return l0_tree;
3091
3092 }
3093
3094 /**
3095 * Method to retrieve the level0 tree (YODA tree) that contains the current event.
3096 */
3097 Int_t PamLevel2::GetYodaEntry(){
3098
3099 // cout << "Int_t PamLevel2::GetYodaEntry()"<<endl;
3100 if(!GetYodaTree())return 0;
3101
3102 // patch
3103 if( irunentry < 0){
3104 // cout << "Int_t PamLevel2::GetYodaEntry() -- ATTENZIONE -- irunentry negativo?!?! "<<(Int_t)irunentry<<endl;
3105 irunentry=0LL;
3106 }
3107 // ---------------------------------
3108 // if file is NOT a preselected file
3109 // ---------------------------------
3110 Long64_t quellagiusta = irunentry + (Long64_t)(run_obj->EV_FROM);
3111 // cout << " irunentry "<<irunentry<<endl;
3112 // cout << " EV_FROM "<<run_obj->EV_FROM<<endl;
3113 // cout << " quellagiusta = irunentry + EV_FROM "<< quellagiusta << endl;
3114
3115 // cout << " irun "<< irun << " irunentry "<< irunentry<<" run_obj->EV_FROM "<<run_obj->EV_FROM <<" quella giusta "<<quellagiusta << endl;
3116 // cout << " iroot "<<iroot<<" run_obj->ID_ROOT_L0 "<<run_obj->ID_ROOT_L0<<endl;
3117 // cout << " time "<< GetOrbitalInfo()->absTime << endl;
3118 // cout << " trk_calib_used "<<run_obj->TRK_CALIB_USED<< endl;
3119
3120 if( !GetOrbitalInfo() ){
3121 cout << "Int_t PamLevel2::GetYodaEntry() -- ERROR -- missing OrbitalInfo "<<endl;
3122 return 0;
3123 }
3124 if( GetOrbitalInfo()->OBT==0 && GetOrbitalInfo()->pkt_num==0 ){
3125 cout << "Int_t PamLevel2::GetYodaEntry() -- ERROR -- level2 event corrupted ?? "<<endl;
3126 return 0;
3127 }
3128
3129 // ---------------------------------------------------------------------
3130 // ATTENTION!!!
3131 // If data are not pre-processed with cleaner, the level0 tree may contain
3132 // spurious nested physics packets.
3133 // The GL_RUN variables EV_FROM, EV_TO and NEVENTS counts also these entries
3134 // while level2 tree DOES NOT!!
3135 // This means that "quellagiusta" in these cases is not correct.
3136 // In order to retrieve the correct level0 event, I implemented a check
3137 // of the OBT and pkt-number. In case of mismatch, the level0 entry number
3138 // is shift forward until when the packets match.
3139 // ---------------------------------------------------------------------
3140 Int_t answer = 0;
3141 Int_t shift =0;
3142 // printf(" siamo qui %i %i \n",shift,prevshift);
3143 Int_t maxshift = 10;
3144 do{
3145 if(shift>0){
3146 cout << " PKTNUM L2 --- "<< GetOrbitalInfo()->pkt_num << " --- L0 --- "<< GetEventHeader()->GetPscuHeader()->GetCounter()<<endl;
3147 cout << " RUN: ID "<< GetRunInfo()->ID << " ID_ROOT_L0 "<<run_obj->ID_ROOT_L0<<" ID_RUN_FRAG "<<GetRunInfo()->ID_RUN_FRAG << " EV_FROM "<<GetRunInfo()->EV_FROM <<endl;
3148 cout << " L2 <--> L0 mismatch ( irun "<<irun<<" irunentry "<<irunentry<<" shift "<<shift<<" prevshift "<<prevshift<<" )"<<endl;
3149 }
3150 answer = l0_tree->GetEntry(quellagiusta+(Long64_t)shift+(Long64_t)prevshift);
3151 shift++;
3152 if( !GetEventHeader() ){
3153 cout << "Int_t PamLevel2::GetYodaEntry() -- ERROR -- missing EventHeader "<<endl;
3154 return 0;
3155 }
3156 // cout << "PKTNUM "<<shift<<" == L2 --- "<< GetOrbitalInfo()->pkt_num << " --- L0 --- "<< GetEventHeader()->GetPscuHeader()->GetCounter()<<endl;
3157 // cout << " L2 --- "<< GetOrbitalInfo()->OBT << " --- L0 --- "<< GetEventHeader()->GetPscuHeader()->GetOrbitalTime()<<endl;
3158 // if( (quellagiusta+shift) == l0_tree->GetEntries()+1 )cout << ">>> end of level0 tree <<<"<<endl;
3159 // cout << " GetOrbitalInfo()->OBT "<< GetOrbitalInfo()->OBT << endl;
3160 // cout << " GetEventHeader()->GetPscuHeader()->GetOrbitalTime() "<< GetEventHeader()->GetPscuHeader()->GetOrbitalTime() << endl;
3161 // cout << " GetOrbitalInfo()->pkt_num "<< GetOrbitalInfo()->pkt_num << endl;
3162 // cout << " GetEventHeader()->GetPscuHeader()->GetCounter() "<< GetEventHeader()->GetPscuHeader()->GetCounter() << endl;
3163 // printf(" IDRUN %u \n",GetRunInfo()->ID);
3164 //
3165 if ( prevshift != 0 && (quellagiusta+(Long64_t)shift) == GetYodaTree()->GetEntries() ){
3166 prevshift = 0;
3167 shift = -1;
3168 };
3169
3170 }while( ( GetOrbitalInfo()->OBT != (UInt_t)(GetEventHeader()->GetPscuHeader()->GetOrbitalTime()) || GetOrbitalInfo()->pkt_num != (UInt_t)(GetEventHeader()->GetPscuHeader()->GetCounter())) && (quellagiusta+(Long64_t)shift) < GetYodaTree()->GetEntries() && shift < maxshift);
3171
3172 if ( (quellagiusta+(Long64_t)shift+(Long64_t)prevshift) > GetYodaTree()->GetEntries() || shift == maxshift ) {
3173 cout << " Big trouble here, no such event in Level0 data! " <<endl;
3174 return 0;
3175 }
3176 // cout << "LA ENTRY GIUSTA E`: "<<quellagiusta<<" (spero...)"<<endl;
3177 // return GetYodaTree()->GetEntry(quellagiusta);
3178 if ( shift > 1 ) prevshift += (shift-1);
3179
3180 return answer;
3181
3182 }
3183
3184

  ViewVC Help
Powered by ViewVC 1.1.23