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

Contents of /PamelaLevel2/src/PamLevel2.cpp

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1.33 - (show annotations) (download)
Wed Mar 21 09:07:32 2007 UTC (17 years, 11 months ago) by pam-fi
Branch: MAIN
Changes since 1.32: +148 -40 lines
...

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

  ViewVC Help
Powered by ViewVC 1.1.23