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

Diff of /PamelaLevel2/src/PamLevel2.cpp

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 1.1 by pam-fi, Fri Jun 16 16:43:55 2006 UTC revision 1.38 by mocchiut, Fri Apr 27 13:58:04 2007 UTC
# Line 7  Line 7 
7   * Default constructor   * Default constructor
8   */   */
9  PamTrack::PamTrack(){  PamTrack::PamTrack(){
10      trk_track  = this->TrkTrack::GetTrkTrack();      trk_track  = 0;
11      calo_track = this->CaloTrkVar::GetCaloTrkVar();      calo_track = 0;
12      tof_track  = this->ToFTrkVar::GetToFTrkVar();      tof_track  = 0;
13        candeleteobj = 0;
14  };  };
15  //--------------------------------------  //--------------------------------------
16  //  //
# Line 19  PamTrack::PamTrack(){ Line 20  PamTrack::PamTrack(){
20   * Constructor   * Constructor
21   */   */
22  PamTrack::PamTrack(TrkTrack* t, CaloTrkVar* c, ToFTrkVar* o){  PamTrack::PamTrack(TrkTrack* t, CaloTrkVar* c, ToFTrkVar* o){
23      trk_track  = this->TrkTrack::GetTrkTrack();  
24      calo_track = this->CaloTrkVar::GetCaloTrkVar();      trk_track  = 0;
25      tof_track  = this->ToFTrkVar::GetToFTrkVar();      calo_track = 0;
26      if(t) *trk_track  = *t;      tof_track  = 0;
27      if(c) *calo_track = *c;  //     if(t)trk_track  = new TrkTrack(*t);
28      if(o) *tof_track  = *o;  //     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   * Default Constructor
89   */   */
90  PamLevel2::PamLevel2(){  PamLevel2::PamLevel2(){
91      trk_obj  = this->TrkLevel2::GetTrkLevel2();    Initialize();
92      calo_obj = this->CaloLevel2::GetCaloLevel2();  };
93      tof_obj  = this->ToFLevel2::GetToFLevel2();  
94      trig_obj = this->TrigLevel2::GetTrigLevel2();  /**
95      s4_obj   = this->S4Level2::GetS4Level2();   * Constructor with given dir, list and detectors
96      nd_obj   = this->NDLevel2::GetNDLevel2();   */
97      ac_obj   = this->AcLevel2::GetAcLevel2();  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.   * Retrieves the it-th stored track.
585   * It override TrkLevel2::GetTrack(int it).   * It override TrkLevel2::GetTrack(int it).
586   * @param itrk Track number, ranging from 0 to GetNTracks().   * @param itrk Track number, ranging from 0 to GetNTracks().
# Line 55  PamLevel2::PamLevel2(){ Line 588  PamLevel2::PamLevel2(){
588    
589  PamTrack* PamLevel2::GetStoredTrack(Int_t itrk){  PamTrack* PamLevel2::GetStoredTrack(Int_t itrk){
590            
591      // retrieve itrk-th tracker stored track      cout <<"PamLevel2::GetStoredTrack(Int_t itrk) **to-be-updated** "<<endl;
592      TrkTrack   *t = 0;      cout <<"for the moment, better use separately the methods: TrkLevel2::GetStoredTrack(seqno) CaloLevel2::GetCaloTrkVar(Int_t notrack) ToFLevel2::GetToFTrkVar(Int_t notrack)"<<endl;
593      CaloTrkVar *c = 0;      cout <<"(if you use it, remember to delete the PamTrack object)"<<endl;
594      ToFTrkVar  *o = 0;      PamTrack *track = 0;
595        
596      if( itrk >=0 && itrk < TrkLevel2::ntrk() ){      if( itrk >=0 && itrk < trk2_obj->TrkLevel2::ntrk() ){
597          t = TrkLevel2::GetStoredTrack(itrk);              
598          c = CaloLevel2::GetCaloTrkVar(t->GetSeqNo());          TrkTrack *t = trk2_obj->TrkLevel2::GetStoredTrack(itrk);
599          o = ToFLevel2::GetToFTrkVar(t->GetSeqNo());          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      };      };
     //  retrieve related ToF-track  
606            
     // hence create a "PamTrack" object  
     PamTrack *track = 0;  
     if(t && c && o)track = new PamTrack(t,c,o);  
607      return track;      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 a TObjectArray (of TrkTrack objects).   * Sort physical (tracker) tracks and stores them in the TRefArray (of TrkTrack objects) which pointer is  PamLevel2::sorted_tracks.
626   * The total number of physical tracks is given by GetNTracks() and the it-th physical track can be retrieved by means of the method GetTrack(int it).   * @param how String to set the sorting cryterium (es: "CAL" or "TRK+CAL+TOF" ecc...).
627   * This method overrides TrkLevel2::GetTracks(), where sorting is done by decreasing number of fit points and increasing chi^2.   * Sorting cryteria:
628   * PamLevel2::GetTracks() keeps the same track order given by TrkLevel2::GetTracks(), but checks image selection by using calorimeter and ToF tracking information.   * 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  TClonesArray *PamLevel2::GetTracks(){  void PamLevel2::SortTracks(){
636      TString how = howtosort;
637    
638      // define new array for sorted tracker tracks    //    cout <<" PamLevel2::SortTracks(TString how) "<<endl;
639      TClonesArray *aa = new TClonesArray("TrkTrack");    if( !trk2_obj ){
640      TClonesArray &sorted = *aa;      cout << "void PamLevel2::SortTracks():  TrkLevel2 not loaded !!!";
641        return;
642      // loop over the tracks sorted by the tracker    };
643      for(Int_t i=0; i < TrkLevel2::GetNTracks(); i++){    //    cout << "call SortTracs() "<<endl;
644      //Save current Object count
645          TrkTrack *ts=0;    Int_t ObjectNumber = TProcessID::GetObjectCount();
646    
647          // get tracker tracks    //    cout << "ObjectNumber  "<<ObjectNumber <<endl;
648          TrkTrack   *tp = TrkLevel2::GetTrack(i);                    //tracker  
649          CaloTrkVar *cp = CaloLevel2::GetCaloTrkVar(tp->GetSeqNo()); //calorimeter    //     if(!sorted_tracks)sorted_tracks = new TRefArray();
650          ToFTrkVar  *op = ToFLevel2::GetToFTrkVar(tp->GetSeqNo());   //tof    //     sorted_tracks->Clear();
651      //    sorted_tracks.Clear();
652          // if track has an image, check image selection  
653          if(tp->HasImage()){    if(!tsorted)tsorted = new TClonesArray("PamTrack",trk2_obj->GetNTracks());
654              TrkTrack   *ti = TrkLevel2::GetTrackImage(i);              //tracker (image)    tsorted->Delete();
655              CaloTrkVar *ci = CaloLevel2::GetCaloTrkVar(ti->GetSeqNo());//calorimeter (image)    TClonesArray &ttsorted = *tsorted;
656              ToFTrkVar  *oi = ToFLevel2::GetToFTrkVar(ti->GetSeqNo());  //tof (image)    if(!timage)timage = new TClonesArray("PamTrack",trk2_obj->GetNTracks());
657      timage->Delete();
658              //assign starting scores    TClonesArray &ttimage = *timage;
659              Int_t tp_score = 1;  
660              Int_t ti_score = 0;  
661              // ------------------------    // loop over the tracks sorted by the tracker
662              // calorimeter check    Bool_t use_TRK = how.Contains("TRK", TString::kIgnoreCase);
663              // ------------------------    Bool_t use_CAL = how.Contains("CAL", TString::kIgnoreCase);
664              if(    Bool_t use_TOF = how.Contains("TOF", TString::kIgnoreCase);
665                        npcfit[1] > 3   &&  //no. of fit planes on Y view      
666  //              varcfit[1] < 50.&&  //fit variance on Y view    if( !CAL2 &&  use_CAL) use_CAL = false;
667                        true){    if( !TOF &&  use_TOF) use_TOF = false;
668        
669                  Float_t resy_p = cp->tbar[0][1] - cbar[0][1];if(resy_p < 0)resy_p= - resy_p;    if( !TRK2 ){
670                  Float_t resy_i = ci->tbar[0][1] - cbar[0][1];if(resy_i < 0)resy_i= - resy_i;      //  cout << "SortTracks() : without tracker does not work!!! (not yet)" << endl;
671        return;
672                  if(resy_p <= resy_i) tp_score++;    };
                 else                 ti_score++;  
                       };  
             // ------------------------  
             // TOF check  
             // ------------------------      
   
             // ------------------------  
             // the winner is....  
             // ------------------------      
                       if(tp_score > ti_score) ts = tp;//the track sorted by the tracker!!  
                       else                    ts = ti;//its image!!  
673    
674          }else{    //   cout << "use_CAL "<<use_CAL<<" use_TOF "<<use_TOF<<" use_TRK "<<use_TRK <<endl;
675              ts = tp;      
676          };    for(Int_t i=0; i < trk2_obj->TrkLevel2::GetNTracks(); i++){
677                    
678          new(sorted[i]) TrkTrack(*ts); //save the track in the sorted array      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      return aa;        };
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.   * Retrieves the it-th Pamela "physical" track.
1009   * It override TrkLevel2::GetTrack(int it).   * It override TrkLevel2::GetTrack(int it).
1010   * @param it Track number, ranging from 0 to GetNTracks().   * @param it Track number, ranging from 0 to GetNTracks().
1011   */   */
1012  PamTrack *PamLevel2::GetTrack(int it){  PamTrack *PamLevel2::GetTrack(int it){
1013    
1014      TClonesArray *aa = this->GetTracks();  //    if(!trk2_obj) return 0;
1015      TClonesArray &sorted = *aa;          
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      TrkTrack   *t = 0;  //    cout << "PamLevel2::GetTrack(int it) "<<endl;
1034      CaloTrkVar *c = 0;  //  *-*-*-*-*-*-*-*-*-*-*-*-*
1035      ToFTrkVar  *o = 0;      SortTracks();
1036    //  *-*-*-*-*-*-*-*-*-*-*-*-*
1037        if(!tsorted)return 0;
1038        if(!tsorted->GetEntries())return 0;
1039    
1040      if( it >=0 && it < TrkLevel2::GetNTracks() ){      PamTrack *track = 0;
1041          t = (TrkTrack*)sorted[it];      
1042          c = CaloLevel2::GetCaloTrkVar(t->GetSeqNo());      if( it >=0 && it < trk2_obj->TrkLevel2::GetNTracks()){
1043          o = ToFLevel2::GetToFTrkVar(t->GetSeqNo());  //      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            
     // hence create a "PamTrack" object  
     PamTrack *track = 0;  
     if(t && c && o)track = new PamTrack(t,c,o);  
1052      return track;      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::GetTracks().   * 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().   * @param it Track number, ranging from 0 to GetNTracks().
1063   */   */
1064  PamTrack *PamLevel2::GetTrackImage(int it){  PamTrack *PamLevel2::GetTrackImage(int it){
1065    
1066      TClonesArray *aa = this->GetTracks();  //     if(!trk2_obj) return 0;
1067      TClonesArray &sorted = *aa;  
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    
     TrkTrack   *t = 0;  
     CaloTrkVar *c = 0;  
     ToFTrkVar  *o = 0;  
1090    
1091      if( it >=0 && it < TrkLevel2::GetNTracks() ){  //  *-*-*-*-*-*-*-*-*-*-*-*-*
1092          TrkTrack *temp = (TrkTrack*)sorted[it];      SortTracks();
1093          if( temp->HasImage() ){  //  *-*-*-*-*-*-*-*-*-*-*-*-*
1094              t = TrkLevel2::GetStoredTrack(temp->GetImageSeqNo());      if(!timage)return 0;
1095              c = CaloLevel2::GetCaloTrkVar(temp->GetImageSeqNo());      if(!timage->GetEntries())return 0;
1096              o = ToFLevel2::GetToFTrkVar(temp->GetImageSeqNo());      
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            
     // hence create a "PamTrack" object  
     PamTrack *image = 0;  
     if(t && c && o)image = new PamTrack(t,c,o);  
1112      return image;      return image;
       
1113  }  }
1114    
1115  //--------------------------------------  //--------------------------------------
# Line 209  PamTrack *PamLevel2::GetTrackImage(int i Line 1117  PamTrack *PamLevel2::GetTrackImage(int i
1117  //  //
1118  //--------------------------------------  //--------------------------------------
1119  /**  /**
1120   * Get the Pamela detector trees and make them friends.   * 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::LoadPamTrees(TFile *f){  TTree *PamLevel2::GetPamTree(TFile *f, TString detlist ){
1126        
1127      TTree *Tout =0;  //     if( !detlist.IsNull() ) SetWhichTrees(detlist);
1128        //     else                    GetWhichTrees(f);
1129      // Tracker    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");      TTree *T = (TTree*)f->Get("Tracker");
1165      if(T) {      if(T && (TRK2||TRK1||TRKh)) {
1166          T->SetBranchAddress("TrkLevel2", GetPointerToTrk());          if(TRK2)T->SetBranchAddress("TrkLevel2", GetPointerTo("TrkLevel2"));
1167          cout << "Tracker      : set branch address TrkLevel2"<<endl;  //      else    T->SetBranchStatus("TrkLevel2",0,found);
1168          if(!Tout)Tout=T;          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{      }else{
1178          cout << "Tracker      : missing tree"<<endl;          cout << "Tracker      : missing tree"<<endl;
1179      };      };
1180      // Calorimeter      // Calorimeter
1181      TTree *C = (TTree*)f->Get("Calorimeter");      TTree *C = (TTree*)f->Get("Calorimeter");
1182      if(C) {      if(C && (CAL2||CAL1)) {
1183          C->SetBranchAddress("CaloLevel2", GetPointerToCalo());          if(CAL2)C->SetBranchAddress("CaloLevel2", GetPointerTo("CaloLevel2"));
1184          cout << "Calorimeter  : set branch address CaloLevel2"<<endl;  //      else    C->SetBranchStatus("CaloLevel2",0,found);
1185          if(!Tout)Tout=C;          if(CAL2)cout << "Calorimeter  : set branch address CaloLevel2"<<endl;
1186          else Tout->AddFriend(C);          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{      }else{
1192          cout << "Calorimeter  : missing tree"<<endl;          cout << "Calorimeter  : missing tree"<<endl;
1193      };      };
1194    
1195      // ToF          // ToF    
1196      TTree *O = (TTree*)f->Get("ToF");      TTree *O = (TTree*)f->Get("ToF");
1197      if(O) {      if(O && TOF) {
1198          O->SetBranchAddress("ToFLevel2", GetPointerToToF());          O->SetBranchAddress("ToFLevel2", GetPointerTo("ToFLevel2"));
1199          cout << "ToF          : set branch address ToFLevel2"<<endl;          cout << "ToF          : set branch address ToFLevel2"<<endl;
1200          if(!Tout)Tout=O;          if(!Trout)Trout=O;
1201          else Tout->AddFriend(O);          else Trout->AddFriend(O);
1202      }else{      }else{
1203          cout << "ToF         : missing tree"<<endl;          cout << "ToF         : missing tree"<<endl;
1204      };      };
1205      // Trigger      // Trigger
1206      TTree *R = (TTree*)f->Get("Trigger");      TTree *R = (TTree*)f->Get("Trigger");
1207      if(R) {      if(R && TRG) {
1208          R->SetBranchAddress("TrigLevel2", GetPointerToTrig());          R->SetBranchAddress("TrigLevel2", GetPointerTo("TrigLevel2"));
1209          cout << "Trigger      : set branch address TrigLevel2"<<endl;          cout << "Trigger      : set branch address TrigLevel2"<<endl;
1210          if(!Tout)Tout=O;          if(!Trout)Trout=R;
1211          else Tout->AddFriend(R);          else Trout->AddFriend(R);
1212      }else{      }else{
1213          cout << "Trigger      : missing tree"<<endl;          cout << "Trigger      : missing tree"<<endl;
1214      };      };
1215      // S4      // S4
1216      TTree *S = (TTree*)f->Get("S4");      TTree *S = (TTree*)f->Get("S4");
1217      if(S) {      if(S && S4) {
1218          S->SetBranchAddress("S4Level2", GetPointerToS4());          S->SetBranchAddress("S4Level2", GetPointerTo("S4Level2"));
1219          cout << "S4           : set branch address S4Level2"<<endl;          cout << "S4           : set branch address S4Level2"<<endl;
1220          if(!Tout)Tout=O;          if(!Trout)Trout=S;
1221          else Tout->AddFriend(S);          else Trout->AddFriend(S);
1222      }else{      }else{
1223          cout << "S4           : missing tree"<<endl;          cout << "S4           : missing tree"<<endl;
1224      };      };
1225      // Neutron Detector      // Neutron Detector
1226      TTree *N = (TTree*)f->Get("NeutronD");      TTree *N = (TTree*)f->Get("NeutronD");
1227      if(N) {      if(N && ND) {
1228          N->SetBranchAddress("NDLevel2", GetPointerToND());          N->SetBranchAddress("NDLevel2", GetPointerTo("NDLevel2"));
1229          cout << "NeutronD     : set branch address NDLevel2"<<endl;          cout << "NeutronD     : set branch address NDLevel2"<<endl;
1230          if(!Tout)Tout=O;          if(!Trout)Trout=N;
1231          else Tout->AddFriend(N);          else Trout->AddFriend(N);
1232      }else{      }else{
1233          cout << "NeutronD     : missing tree"<<endl;          cout << "NeutronD     : missing tree"<<endl;
1234      };      };
1235      // Anticounters      // Anticounters
1236      TTree *A = (TTree*)f->Get("Anticounter");      TTree *A = (TTree*)f->Get("Anticounter");
1237      if(A) {      if(A && AC) {
1238          A->SetBranchAddress("AcLevel2", GetPointerToAc());          A->SetBranchAddress("AcLevel2", GetPointerTo("AcLevel2"));
1239          cout << "Anticounter  : set branch address AcLevel2"<<endl;          cout << "Anticounter  : set branch address AcLevel2"<<endl;
1240          if(!Tout)Tout=O;          if(!Trout)Trout=A;
1241          else Tout->AddFriend(A);          else Trout->AddFriend(A);
1242      }else{      }else{
1243          cout << "Anticounter  : missing tree"<<endl;          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 Tout;      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                run_tree->GetEntry(irun);
2003            };      
2004            //
2005            bool fromfirst = true; // first loop over runs
2006            // ------------------------------------------------------
2007            // loop over runs to find the one that contains the event
2008            // ------------------------------------------------------
2009            while ( !(GetOrbitalInfo()->absTime >= GetRunInfo()->RUNHEADER_TIME &&    // check on absolute time (s)
2010                      GetOrbitalInfo()->absTime <= GetRunInfo()->RUNTRAILER_TIME) &&
2011                    !(GetOrbitalInfo()->OBT >= GetRunInfo()->RUNHEADER_OBT &&         // additional check on OBT (ms)
2012                      GetOrbitalInfo()->OBT <= GetRunInfo()->RUNTRAILER_OBT) &&
2013                    irun < run_tree->GetEntries() ){
2014                irun++;
2015                // ------------------------------------
2016                // if the end of run tree is reached...
2017                // ------------------------------------
2018                if( irun == run_tree->GetEntries() ){
2019                    if(!fromfirst){
2020                        // -----------------------------------------------------
2021                        // if it happened already once and the run was not found
2022                        // ---> exit with error
2023                        // -----------------------------------------------------
2024                        cout << "Bool_t PamLevel2::UpdateRunInfo(Long64_t iev) -- ERROR -- event entry #"<<iev<<" does not belong to any run (should not happen)" <<endl;
2025                        return false;
2026                    }
2027                    // -----------------------------------------
2028                    // ...otherwise repeat search from beginning
2029                    // -----------------------------------------
2030                    cout << "Bool_t PamLevel2::UpdateRunInfo(Long64_t iev) -- WARNING -- reached end of run tree. searchin again from beginning " <<endl;
2031                    fromfirst = false;
2032                    irun = 0LL;
2033                    runfirstentry = 0LL;
2034                }
2035                // -------------------------------------------------------------------
2036                // save the index of the first entry of the run, relative to pam_tree,
2037                // and read a new run
2038                // -------------------------------------------------------------------
2039                if(irun>0)runfirstentry += (GetRunInfo()->NEVENTS);
2040                run_tree->GetEntry(irun);      
2041                if(GetRunInfo()->RUNHEADER_OBT>=GetRunInfo()->RUNTRAILER_OBT ){
2042                    cout << "Bool_t PamLevel2::UpdateRunInfo(Long64_t iev) -- WARNING -- irun "<<irun<<"  has RUNHEADER_OBT>=RUNTRAILER_OBT " <<endl;
2043                    cout << "                                                            (NB!! in this case some events are assigned to a wrong run)"<<endl;
2044                }
2045            };
2046            // --------------------------------------
2047            // if there was no need to update the run
2048            // ---> exit with FALSE
2049            // --------------------------------------
2050            if ( irun == oldrun ) return(false);
2051    
2052            // --------------------------------------
2053            // ... otherwise
2054            // ---> exit with TRUE
2055            // --------------------------------------
2056            cout << endl << " ))))) UPDATE RUN INFO (((((  @iev "<<iev<<" run "<<GetRunInfo()->ID<<" irun "<<irun<<endl;
2057            prevshift = 0;
2058            return(true);    
2059        };
2060        // +-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-
2061        // if it is a preselected file (there is SelectionList)
2062        // +-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-
2063        if(SELLI==1){      
2064            sel_tree->GetEntry(iev);
2065            if(irun != oldrun){
2066                run_tree->GetEntry(irun);
2067                cout << endl << " ))))) UPDATE RUN INFO (((((  @iev "<<iev<<" run "<<GetRunInfo()->ID<<" irun "<<irun<<endl;
2068                prevshift = 0;
2069                return true;
2070            }
2071            return false;
2072        }
2073    
2074        return false;
2075        //
2076    };
2077    /**
2078     * Update the runinfo informations (to be used to have Run infos event by event basis)
2079     * @param run Pointer to the chain/tree which contains run infos
2080     * @return true if a new run has been read, false if it is still the same run
2081     */
2082    Bool_t PamLevel2::UpdateRunInfo(TTree *run, Long64_t iev){
2083      return(UpdateRunInfo((TChain*)run,iev));
2084    };
2085    
2086    //--------------------------------------
2087    //
2088    //
2089    //--------------------------------------
2090    /**
2091     * Set which trees shoul be analysed
2092     * @param detlist TString containing the sequence of trees required
2093    */
2094    void PamLevel2::SetWhichTrees(TString detlist){
2095            
2096        if(detlist.IsNull() || detlist.Contains("+ALL", TString::kIgnoreCase)){
2097            CAL0 = false;
2098            CAL1 = true;
2099            CAL2 = true;
2100            TRK2 = true;
2101            TRK1 = false;
2102            TRKh = false;
2103            TRK0 = false;
2104            TRG = true;
2105            TOF = true;
2106            TOF0 = false;
2107            S4  = true;
2108            ND  = true;
2109            AC  = true;
2110            ORB = true;
2111        }else if( detlist.Contains("-ALL", TString::kIgnoreCase) ){
2112            CAL0 = false;
2113            CAL1 = false;
2114            CAL2 = false;
2115            TRK2 = false;
2116            TRK1 = false;
2117            TRKh = false;
2118            TRK0 = false;
2119            TRG = false;
2120            TOF = false;
2121            TOF0 = false;
2122            S4  = false;
2123            ND  = false;
2124            AC  = false;
2125            ORB = false;
2126        };
2127        
2128    //  -------------------------------------------------------------------------
2129        if( detlist.Contains("CAL1", TString::kIgnoreCase) ){
2130            if ( detlist.Contains("-CAL1", TString::kIgnoreCase) )CAL1=false;
2131            if ( detlist.Contains("+CAL1", TString::kIgnoreCase) )CAL1=true;
2132        };
2133    
2134        if( detlist.Contains("CAL0", TString::kIgnoreCase) ){
2135            if ( detlist.Contains("-CAL0", TString::kIgnoreCase) )CAL0=false;
2136            if ( detlist.Contains("+CAL0", TString::kIgnoreCase) )CAL0=true;
2137        };
2138            
2139        if( detlist.Contains("CAL2", TString::kIgnoreCase)){
2140            if ( detlist.Contains("-CAL2", TString::kIgnoreCase) )CAL2=false;
2141            if ( detlist.Contains("+CAL2", TString::kIgnoreCase) )CAL2=true;
2142        };
2143            
2144        if( detlist.Contains("+CAL", TString::kIgnoreCase) && !CAL1 && !CAL2 )CAL2=true;
2145        if( detlist.Contains("-CAL", TString::kIgnoreCase) && CAL1 && CAL2 ){
2146            CAL2=false;
2147            CAL1=false;
2148        }
2149    //  -------------------------------------------------------------------------
2150        if( detlist.Contains("TRK0", TString::kIgnoreCase) ){
2151            if ( detlist.Contains("-TRK0", TString::kIgnoreCase) )TRK0=false;
2152            if ( detlist.Contains("+TRK0", TString::kIgnoreCase) )TRK0=true;
2153        };
2154    
2155        if( detlist.Contains("TRK1", TString::kIgnoreCase) ){
2156            if ( detlist.Contains("-TRK1", TString::kIgnoreCase) )TRK1=false;
2157            if ( detlist.Contains("+TRK1", TString::kIgnoreCase) )TRK1=true;
2158        };
2159            
2160        if( detlist.Contains("TRK2", TString::kIgnoreCase)){
2161            if ( detlist.Contains("-TRK2", TString::kIgnoreCase) )TRK2=false;
2162            if ( detlist.Contains("+TRK2", TString::kIgnoreCase) )TRK2=true;
2163        };
2164    
2165        if( detlist.Contains("TRKh", TString::kIgnoreCase)){
2166            if ( detlist.Contains("-TRKh", TString::kIgnoreCase) )TRKh=false;
2167            if ( detlist.Contains("+TRKh", TString::kIgnoreCase) )TRKh=true;
2168        };
2169            
2170        if( detlist.Contains("+TRK", TString::kIgnoreCase) && !TRK1 && !TRK2 && !TRKh )TRK2=true;
2171        if( detlist.Contains("-TRK", TString::kIgnoreCase) && TRK1 && TRK2 && TRKh){
2172            TRK2=false;
2173            TRK1=false;
2174            TRKh=false;
2175        }
2176    //  -------------------------------------------------------------------------
2177        
2178        if( detlist.Contains("-TRG", TString::kIgnoreCase) )TRG = false;
2179        else if( detlist.Contains("+TRG", TString::kIgnoreCase) )TRG = true;
2180        
2181        if( detlist.Contains("-TOF", TString::kIgnoreCase) )TOF = false;
2182        else if( detlist.Contains("+TOF", TString::kIgnoreCase) )TOF = true;
2183    
2184        if( detlist.Contains("-TOF0", TString::kIgnoreCase) )TOF0 = false;
2185        else if( detlist.Contains("+TOF0", TString::kIgnoreCase) )TOF0 = true;
2186        
2187        if( detlist.Contains("-S4",  TString::kIgnoreCase) )S4  = false;
2188        else if( detlist.Contains("+S4",  TString::kIgnoreCase) )S4  = true;
2189        
2190        if( detlist.Contains("-ND",  TString::kIgnoreCase) )ND  = false;
2191        else if( detlist.Contains("+ND",  TString::kIgnoreCase) )ND  = true;
2192        
2193        if( detlist.Contains("-AC",  TString::kIgnoreCase) )AC  = false;
2194        else if( detlist.Contains("+AC",  TString::kIgnoreCase) )AC  = true;
2195        
2196        if( detlist.Contains("-ORB", TString::kIgnoreCase) )ORB = false;
2197        else if( detlist.Contains("+ORB", TString::kIgnoreCase) )ORB = true;
2198    
2199    //     cout<< "Set detector list --> ";
2200    //     if(TRK1)cout<<"TRK1 ";
2201    //     if(TRK2)cout<<"TRK2 ";
2202    //     if(TRKh)cout<<"TRKH ";
2203    //     if(CAL1)cout<<"CAL1 ";
2204    //     if(CAL2)cout<<"CAL2 ";
2205    //     if(TOF)cout<<"TOF ";
2206    //     if(TRG)cout<<"TRG ";
2207    //     if(AC)cout<<"AC ";
2208    //     if(ND)cout<<"ND ";
2209    //     if(S4)cout<<"S4 ";
2210    //     if(ORB)cout<<"ORB ";
2211    //     cout << endl;
2212        
2213    };
2214    
2215    
2216    /**
2217     * Set tree/branch detector flags from the content of a tree
2218     */
2219    void  PamLevel2::GetWhichTrees(TFile* f){
2220        
2221    
2222    
2223        cout << "void  PamLevel2::GetWhichTrees(TFile* f) --- WARNING!! --- ...potrebbe non funzionare "<<endl;
2224        // -----------
2225        // reset flags
2226        // -----------
2227        CAL1   = false;    
2228        CAL2   = false;    
2229        TRK2   = false;    
2230        TRK1   = false;    
2231        TRKh   = false;    
2232        TRG    = false;    
2233        TOF    = false;    
2234        S4     = false;    
2235        ND     = false;    
2236        AC     = false;    
2237        ORB    = false;    
2238        
2239        RUN    = false;
2240            
2241        cout << "Checking file: "<<f->GetName()<<endl;
2242        if( !f || f->IsZombie() ){
2243            cout << "File: "<< f->GetName() <<" Non valid root file"<< endl;
2244            return;
2245        }
2246    
2247        TList *lk = f->GetListOfKeys();
2248        TIter next(lk);
2249        TKey *key =0;
2250    
2251        Int_t nev = 0;
2252    
2253        while( (key = (TKey*)next()) ){
2254            
2255            if( !strcmp(key->GetName(),"Run"        ) )RUN = true;
2256    
2257            //=========================================================    
2258            if( !strcmp(key->GetName(),"Trigger"    ) ){
2259                TRG = true;
2260                Int_t nevt = ((TTree*)f->Get("Trigger"))->GetEntries();
2261                if( nev && nevt!=nev){
2262                    cout << "File: "<< f->GetName() <<" Trigger tree has "<<nevt<<" events instead of "<<nev<< endl;
2263                    TRG = false;
2264                }else nev=nevt;
2265            }
2266            //=========================================================    
2267            if( !strcmp(key->GetName(),"ToF"        ) ){
2268                TOF = true;
2269                Int_t nevt = ((TTree*)f->Get("ToF"))->GetEntries();
2270                if( nev && nevt!=nev){
2271                    cout << "File: "<< f->GetName() <<"     ToF tree has "<<nevt<<" events instead of "<<nev<< endl;
2272                    TOF = false;
2273                }else nev=nevt;
2274            }
2275            //=========================================================  
2276            if( !strcmp(key->GetName(),"S4"         ) ){
2277                S4 = true;
2278                Int_t nevt = ((TTree*)f->Get("S4"))->GetEntries();
2279                if( nev && nevt!=nev){
2280                    cout << "File: "<< f->GetName() <<"      S4 tree has "<<nevt<<" events instead of "<<nev<< endl;
2281                    S4 = false;
2282                }else nev=nevt;
2283            }
2284            //=========================================================  
2285    
2286            if( !strcmp(key->GetName(),"NeutronD"   ) ){
2287                ND = true;
2288                Int_t nevt = ((TTree*)f->Get("NeutronD"))->GetEntries();
2289                if( nev && nevt!=nev){
2290                    cout << "File: "<< f->GetName() <<"NeutronD tree has "<<nevt<<" events instead of "<<nev<< endl;
2291                    ND =false;
2292                }else nev=nevt;
2293            }      
2294            //=========================================================  
2295            if( !strcmp(key->GetName(),"Anticounter") ){
2296                AC = true;
2297                Int_t nevt = ((TTree*)f->Get("Anticounter"))->GetEntries();
2298                if( nev && nevt!=nev){
2299                    cout << "File: "<< f->GetName() <<" Anticounter tree has "<<nevt<<" events instead of "<<nev<< endl;
2300                    AC =false;
2301                }else nev=nevt;
2302            }
2303            //=========================================================  
2304            if( !strcmp(key->GetName(),"OrbitalInfo") ){
2305                ORB = true;
2306                Int_t nevt = ((TTree*)f->Get("OrbitalInfo"))->GetEntries();
2307                if( nev && nevt!=nev){
2308                    cout << "File: "<< f->GetName() <<" OrbitalInfo tree has "<<nevt<<" events instead of "<<nev<< endl;
2309                    ORB = false;
2310                }else nev=nevt;
2311            }
2312            //=========================================================  
2313            if( !strcmp(key->GetName(),"Tracker"    ) ){
2314                TTree *T = (TTree*)f->Get("Tracker");
2315                for(Int_t i=0; i<T->GetListOfBranches()->GetEntries(); i++){
2316                    TString name = T->GetListOfBranches()->At(i)->GetName();
2317                    if( !name.CompareTo("TrkLevel1") )TRK1=true;
2318                    if( !name.CompareTo("TrkLevel2") )TRK2=true;
2319                    if( !name.CompareTo("TrkHough") )TRKh=true;
2320                };      
2321                Int_t nevt = T->GetEntries();
2322                if( nev && nevt!=nev){
2323                    cout << "File: "<< f->GetName() <<" Tracker tree has "<<nevt<<" events instead of "<<nev<< endl;
2324                    TRK1 = false;
2325                    TRK2 = false;
2326                    TRKh = false;
2327                }else nev=nevt;
2328                T->Delete();
2329            };
2330            //=========================================================  
2331            if( !strcmp(key->GetName(),"Calorimeter"    ) ){
2332                TTree *T = (TTree*)f->Get("Calorimeter");
2333                for(Int_t i=0; i<T->GetListOfBranches()->GetEntries(); i++){
2334                    TString name = T->GetListOfBranches()->At(i)->GetName();
2335                    if( !name.CompareTo("CaloLevel1") )CAL1=true;
2336                    if( !name.CompareTo("CaloLevel2") )CAL2=true;
2337                };    
2338                Int_t nevt = T->GetEntries();
2339                if( nev && nevt!=nev){
2340                    cout << "File: "<< f->GetName() <<"  Calorimeter tree has "<<nevt<<" events instead of "<<nev<< endl;
2341                    CAL1 = false;
2342                    CAL2 = false;
2343                }else nev=nevt;
2344                T->Delete();
2345            };      
2346    
2347        };
2348        
2349        delete lk;
2350    
2351    //     cout<< "Get detector list from input file --> ";
2352    //     if(TRK1)cout<<"TRK1 ";
2353    //     if(TRK2)cout<<"TRK2 ";
2354    //     if(TRKh)cout<<"TRKH ";
2355    //     if(CAL1)cout<<"CAL1 ";
2356    //     if(CAL2)cout<<"CAL2 ";
2357    //     if(TOF)cout<<"TOF ";
2358    //     if(TRG)cout<<"TRG ";
2359    //     if(AC)cout<<"AC ";
2360    //     if(ND)cout<<"ND ";
2361    //     if(S4)cout<<"S4 ";
2362    //     if(ORB)cout<<"ORB ";
2363    //     cout << endl;
2364          
2365        return ;
2366            
2367    };
2368    
2369    
2370    //--------------------------------------
2371    //
2372    //
2373    //--------------------------------------
2374    /**
2375     * Check if a file contains selected Pamela Level2 trees.
2376     * @param name File name
2377     * @return true if the file is ok.
2378     */
2379    Bool_t  PamLevel2::CheckLevel2File(TString name){
2380            
2381        Bool_t CAL1__ok   = false;    
2382        Bool_t CAL2__ok   = false;    
2383        Bool_t TRK2__ok   = false;    
2384        Bool_t TRK1__ok   = false;    
2385        Bool_t TRKh__ok   = false;    
2386        Bool_t TRG__ok    = false;    
2387        Bool_t TOF__ok    = false;    
2388        Bool_t S4__ok     = false;    
2389        Bool_t ND__ok     = false;    
2390        Bool_t AC__ok     = false;    
2391        Bool_t ORB__ok    = false;    
2392        
2393        Bool_t RUN__ok    = false;
2394        
2395        Bool_t SELLI__ok  = false;
2396    
2397        cout << "Checking file: "<<name<<endl;
2398        TFile *f = new TFile(name.Data());
2399        if( !f || f->IsZombie() ){
2400            cout << "File: "<< f->GetName() <<" discarded ---- Non valid root file"<< endl; return false;
2401        }
2402    //    cout << "Get list of keys: "<<f<<endl;
2403        TList *lk = f->GetListOfKeys();
2404    //    lk->Print();
2405        TIter next(lk);
2406        TKey *key =0;
2407    
2408        Int_t nev = 0;
2409    
2410        while( (key = (TKey*)next()) ){
2411            
2412    //      cout << key->GetName() << endl;
2413    //      cout << key->GetName() << ""<<key->GetClassName()<<endl;
2414    //              cout << " Get tree: " << f->Get(key->GetName())<<endl;
2415    //      nev_previous = nev;
2416    //      cout << " n.entries  "<< nev <<endl;
2417    //      if( key->GetClassName()=="TTree" && nev_previous && nev != nev_previous ){
2418    //          nev = ((TTree*)f->Get(key->GetName()))->GetEntries();
2419    //          cout << "File: "<< f->GetName() <<" discarded ---- "<< key->GetName() << " tree: n.entries does not match "<<nev<<" "<<nev_previous<< endl;
2420    //          return false;
2421    //      };
2422    
2423    
2424            if( !strcmp(key->GetName(),"Run"        ) )RUN__ok = true;      
2425    
2426            //=========================================================    
2427            if( !strcmp(key->GetName(),"SelectionList"    ) ){
2428                SELLI__ok = true;
2429                if(SELLI==1){
2430                    Int_t nevt = ((TTree*)f->Get("SelectionList"))->GetEntries();
2431                    if( nev && nevt!=nev){
2432                        cout << "File: "<< f->GetName() <<" discarded ---- SelectionList tree has "<<nevt<<" events instead of "<<nev<< endl;
2433                        return false;
2434                    }
2435                    nev=nevt;
2436                }
2437            }
2438    
2439            //=========================================================    
2440            if( !strcmp(key->GetName(),"Trigger"    ) ){
2441                TRG__ok = true;
2442                if(TRG){
2443                    Int_t nevt = ((TTree*)f->Get("Trigger"))->GetEntries();
2444                    if( nev && nevt!=nev){
2445                        cout << "File: "<< f->GetName() <<" discarded ---- Trigger tree has "<<nevt<<" events instead of "<<nev<< endl;
2446                        return false;
2447                    }
2448                    nev=nevt;
2449                }
2450            }
2451            //=========================================================    
2452            if( !strcmp(key->GetName(),"ToF"        ) ){
2453                TOF__ok = true;
2454                if(TOF){
2455                    Int_t nevt = ((TTree*)f->Get("ToF"))->GetEntries();
2456                    if( nev && nevt!=nev){
2457                        cout << "File: "<< f->GetName() <<" discarded ---- ToF tree has "<<nevt<<" events instead of "<<nev<< endl;
2458                        return false;
2459                    }
2460                    nev=nevt;
2461                }
2462            }
2463            //=========================================================  
2464            if( !strcmp(key->GetName(),"S4"         ) ){
2465                S4__ok = true;
2466                if(S4){
2467                    Int_t nevt = ((TTree*)f->Get("S4"))->GetEntries();
2468                    if( nev && nevt!=nev){
2469                        cout << "File: "<< f->GetName() <<" discarded ---- S4 tree has "<<nevt<<" events instead of "<<nev<< endl;
2470                        return false;
2471                    }
2472                    nev=nevt;
2473                }
2474            }
2475            //=========================================================  
2476    
2477            if( !strcmp(key->GetName(),"NeutronD"   ) ){
2478                ND__ok = true;
2479                if(ND){
2480                    Int_t nevt = ((TTree*)f->Get("NeutronD"))->GetEntries();
2481                    if( nev && nevt!=nev){
2482                        cout << "File: "<< f->GetName() <<" discarded ---- NeutronD tree has "<<nevt<<" events instead of "<<nev<< endl;
2483                        return false;
2484                    }
2485                    nev=nevt;
2486                }
2487            }      
2488            //=========================================================  
2489            if( !strcmp(key->GetName(),"Anticounter") ){
2490                AC__ok = true;
2491                if(AC){
2492                    Int_t nevt = ((TTree*)f->Get("Anticounter"))->GetEntries();
2493                    if( nev && nevt!=nev){
2494                        cout << "File: "<< f->GetName() <<" discarded ---- Anticounter tree has "<<nevt<<" events instead of "<<nev<< endl;
2495                        return false;
2496                    }
2497                    nev=nevt;
2498                }
2499            }
2500            //=========================================================  
2501            if( !strcmp(key->GetName(),"OrbitalInfo") ){
2502                ORB__ok = true;
2503                if(ORB){
2504                    Int_t nevt = ((TTree*)f->Get("OrbitalInfo"))->GetEntries();
2505                    if( nev && nevt!=nev){
2506                        cout << "File: "<< f->GetName() <<" discarded ---- OrbitalInfo tree has "<<nevt<<" events instead of "<<nev<< endl;
2507                        return false;
2508                    }
2509                    nev=nevt;
2510                }
2511            }
2512            //=========================================================  
2513            if( !strcmp(key->GetName(),"Tracker"    ) ){
2514                TTree *T = (TTree*)f->Get("Tracker");
2515                if(TRK1||TRK2||TRKh){
2516                    Int_t nevt = T->GetEntries();
2517                    if( nev && nevt!=nev){
2518                        cout << "File: "<< f->GetName() <<" discarded ---- Tracker tree has "<<nevt<<" events instead of "<<nev<< endl;
2519                        return false;
2520                    }
2521                    nev=nevt;
2522                }
2523                for(Int_t i=0; i<T->GetListOfBranches()->GetEntries(); i++){
2524                    TString name = T->GetListOfBranches()->At(i)->GetName();
2525                    if( !name.CompareTo("TrkLevel1") )TRK1__ok=true;
2526                    if( !name.CompareTo("TrkLevel2") )TRK2__ok=true;
2527                    if( !name.CompareTo("TrkHough") )TRKh__ok=true;
2528                };      
2529                T->Delete();
2530            };
2531            //=========================================================  
2532            if( !strcmp(key->GetName(),"Calorimeter"    ) ){
2533                TTree *T = (TTree*)f->Get("Calorimeter");
2534                if(CAL1||CAL2){
2535                    Int_t nevt = T->GetEntries();
2536                    if( nev && nevt!=nev){
2537                        cout << "File: "<< f->GetName() <<" discarded ---- Calorimeter tree has "<<nevt<<" events instead of "<<nev<< endl;
2538                        return false;
2539                    }
2540                    nev=nevt;
2541                }
2542                for(Int_t i=0; i<T->GetListOfBranches()->GetEntries(); i++){
2543                    TString name = T->GetListOfBranches()->At(i)->GetName();
2544                    if( !name.CompareTo("CaloLevel1") )CAL1__ok=true;
2545                    if( !name.CompareTo("CaloLevel2") )CAL2__ok=true;
2546                };    
2547                T->Delete();
2548            };      
2549    
2550        };
2551    
2552        if( SELLI==-1 )SELLI = (Int_t)SELLI__ok;
2553        if( SELLI==0 && SELLI__ok ){
2554            cout << "File: "<< f->GetName() <<" discarded ---- found SelectionList (it is not a full-event file)" << endl;
2555            return false;  
2556        }
2557        if( SELLI==1 && !SELLI__ok ){
2558            cout << "File: "<< f->GetName() <<" discarded ---- SelectionList missing" << endl;
2559            return false;
2560        }
2561        
2562    //    cout << "SELLI "<<SELLI<<endl;
2563    
2564    //     cout<< "CheckLevel2File(TString): detector list --> ";
2565    //     if(TRK1__ok)cout<<"TRK1 ";
2566    //     if(TRK2__ok)cout<<"TRK2 ";
2567    //     if(TRKh__ok)cout<<"TRKH ";
2568    //     if(CAL1__ok)cout<<"CAL1 ";
2569    //     if(CAL2__ok)cout<<"CAL2 ";
2570    //     if(TOF__ok)cout<<"TOF ";
2571    //     if(TRG__ok)cout<<"TRG ";
2572    //     if(AC__ok)cout<<"AC ";
2573    //     if(ND__ok)cout<<"ND ";
2574    //     if(S4__ok)cout<<"S4 ";
2575    //     if(ORB__ok)cout<<"ORB ";
2576    //     cout << endl;
2577    
2578    
2579        if(TRK2 && TRK1__ok)TRK1=1;
2580    // ----------------------------------------------------------------------------
2581    // NOTA
2582    // se c'e` il level1, lo devo necessarimente leggere.
2583    // infatti (non ho capito perche`) i cluster vengono letti e allocati in memoria
2584    // comunque, ma non vengono disallocati da PamLevel2::Clear()
2585    // ----------------------------------------------------------------------------
2586    
2587    
2588        if(!RUN__ok) {
2589            cout << "File: "<< f->GetName() <<" *WARNING* ---- Missing RunInfo tree"<< endl;
2590    //      return false;
2591        };
2592    
2593        if(CAL1 && !CAL1__ok){
2594            cout << "File: "<< f->GetName() <<" discarded ---- Missing CaloLevel1 branch"<< endl;
2595            return false;
2596        };
2597        if(CAL2 && !CAL2__ok){
2598            cout << "File: "<< f->GetName() <<" discarded ---- Missing CaloLevel2 branch"<< endl;
2599            return false;
2600        };
2601        if(TRK2 && !TRK2__ok){
2602            cout << "File: "<< f->GetName() <<" discarded ---- Missing TrkLevel2 branch"<< endl;
2603            return false;
2604        };
2605        if(TRK1 && !TRK1__ok){
2606            cout << "File: "<< f->GetName() <<" discarded ---- Missing TrkLevel1 branch"<< endl;
2607            return false;
2608        };
2609        if(TRKh && !TRKh__ok){
2610            cout << "File: "<< f->GetName() <<" discarded ---- Missing TrkHough branch"<< endl;
2611            return false;
2612        };
2613        if(ORB && !ORB__ok){
2614            cout << "File: "<< f->GetName() <<" discarded ---- Missing ORB tree"<< endl;
2615            return false;
2616        };
2617        if(AC && !AC__ok){
2618            cout << "File: "<< f->GetName() <<" discarded ---- Missing AC tree"<< endl;
2619            return false;
2620        };
2621        if(S4 && !S4__ok){
2622            cout << "File: "<< f->GetName() <<" discarded ---- Missing S4 tree"<< endl;
2623            return false;
2624        };
2625        if(TOF && !TOF__ok){
2626            cout << "File: "<< f->GetName() <<" discarded ---- Missing ToF tree"<< endl;
2627            return false;
2628        };
2629    
2630        if(ND && !ND__ok){
2631            cout << "File: "<< f->GetName() <<" discarded ---- Missing ND tree"<< endl;
2632            return false;
2633        };
2634        if(TRG && !TRG__ok){
2635            cout << "File: "<< f->GetName() <<" discarded ---- Missing Trigger tree"<< endl;
2636            return false;
2637        };
2638        
2639    
2640    //    lk->Delete();
2641    //    delete lk;
2642        f->Close();
2643    
2644    //     cout<< "CheckLevel2File(TString): detector list --> ";
2645    //     if(TRK1)cout<<"TRK1 ";
2646    //     if(TRK2)cout<<"TRK2 ";
2647    //     if(TRKh)cout<<"TRKH ";
2648    //     if(CAL1)cout<<"CAL1 ";
2649    //     if(CAL2)cout<<"CAL2 ";
2650    //     if(TOF)cout<<"TOF ";
2651    //     if(TRG)cout<<"TRG ";
2652    //     if(AC)cout<<"AC ";
2653    //     if(ND)cout<<"ND ";
2654    //     if(S4)cout<<"S4 ";
2655    //     if(ORB)cout<<"ORB ";
2656    //     cout << endl;
2657          
2658        return true;
2659            
2660    };
2661    
2662    
2663    /**
2664     * Create clone-trees
2665     */
2666    void PamLevel2::CreateCloneTrees0( TChain *fChain, TFile *ofile ){
2667    
2668        cout << "+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+" <<endl;
2669        cout << "Create clones of PAMELA trees "<<endl;
2670        
2671        Int_t i=0;
2672        pam_tree_clone[i] = fChain->GetTree()->CloneTree(0);
2673        TString name =  pam_tree_clone[i]->GetName();
2674        name.Append("_clone");
2675    //    pam_tree_clone[i]->SetName(name.Data());
2676        cout << pam_tree_clone[i]->GetName() <<endl;
2677        i++;
2678    
2679        TList *li = fChain->GetListOfFriends();
2680        TIter next(li);
2681        TFriendElement* T_friend=0;
2682        ofile->cd();
2683        while( (T_friend = (TFriendElement*)next()) ){
2684    //      cout<<T_friend->IsA()->GetName()<<" "<<T_friend->GetName()<<hex << T_friend->GetTree() << dec<<endl;
2685    //      cout<<T_friend->GetTree()->GetName()<< endl;
2686            pam_tree_clone[i] = T_friend->GetTree()->CloneTree(0);
2687            pam_tree_clone[i]->SetAutoSave(1000000);
2688            name =  pam_tree_clone[i]->GetName();
2689            name.Append("_clone");
2690    //      pam_tree_clone[i]->SetName(name.Data());
2691            cout << pam_tree_clone[i]->GetName() << endl;
2692            i++;
2693        }
2694        
2695        delete li;
2696    
2697        cout << "+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+" <<endl;
2698    
2699    }
2700    
2701    /**
2702     * Create clone-trees
2703     */
2704    void PamLevel2::CreateCloneTrees(TFile *ofile){
2705    
2706        ofile->cd();
2707    
2708        cout << "+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+" <<endl;
2709        cout << "Create new PAMELA trees "<<endl;
2710            
2711    
2712        run_tree_clone = new TTree("Run","PAMELA Level2 data from the GL_RUN table ");
2713        run_tree_clone->Branch("RunInfo","GL_RUN",GetPointerTo("RunInfo"));
2714        cout << "Run          : branch RunInfo"<<endl;
2715        run_tree_clone->Branch("SoftInfo","SoftInfo",GetPointerTo("SoftInfo"));
2716        cout << "Run          : branch SoftInfo"<<endl;
2717        // ------------------
2718        // replicate run tree
2719        // ------------------
2720        cout << "----------------------------------------------------"<<endl;
2721        cout << "irun\t | RUN\t NEVENTS\t absolute time"<<endl;
2722        for (Int_t i=0; i<run_tree->GetEntries(); i++){
2723            run_tree->GetEntry(i);
2724            cout << i<< "\t | "<<GetRunInfo()->ID<<"\t "<<GetRunInfo()->NEVENTS<< "\t "<<GetRunInfo()->RUNHEADER_TIME<<" <---> "<<GetRunInfo()->RUNTRAILER_TIME<<endl;
2725            run_tree_clone->Fill();
2726        }
2727        cout << "----------------------------------------------------"<<endl;
2728    
2729    
2730        sel_tree_clone = new TTree("SelectionList","List of selected events ");
2731        sel_tree_clone->Branch("RunEntry",&irun,"runentry/L");
2732        sel_tree_clone->Branch("EventEntry",&irunentry,"eventry/L");
2733        
2734    
2735        Int_t i=0;
2736        if(TRK1||TRK2||TRKh){
2737            pam_tree_clone[i] = new TTree("Tracker","PAMELA tracker level2 data ");
2738            if(TRK1) {
2739                pam_tree_clone[i]->Branch("TrkLevel1","TrkLevel1", GetPointerTo("TrkLevel1"));
2740                pam_tree_clone[i]->BranchRef();
2741                cout << "Tracker      : branch TrkLevel1"<<endl;
2742    //          cout << "CreateCloneTrees " << GetTrkLevel1()<<endl;
2743            };
2744            if(TRK2) {
2745                pam_tree_clone[i]->Branch("TrkLevel2", "TrkLevel2",GetPointerTo("TrkLevel2"));
2746                cout << "Tracker      : branch TrkLevel2"<<endl;        
2747            };
2748            if(TRKh) {
2749                pam_tree_clone[i]->Branch("TrkHough","TrkHough", GetPointerTo("TrkHough"));
2750                cout << "Tracker      : branch TrkHough"<<endl;
2751            };
2752            i++;
2753        }
2754    
2755        // Calorimeter
2756        if(CAL1||CAL2){
2757            pam_tree_clone[i] = new TTree("Calorimeter","PAMELA calorimeter level2 data ");
2758            if(CAL1) {
2759                pam_tree_clone[i]->Branch("CaloLevel1", "CaloLevel1", GetPointerTo("CaloLevel1"));
2760                cout << "Calorimeter  : branch CaloLevel1"<<endl;
2761            };
2762            if(CAL2) {
2763                pam_tree_clone[i]->Branch("CaloLevel2","CaloLevel2", GetPointerTo("CaloLevel2"));
2764                cout << "Calorimeter  : branch CaloLevel2"<<endl;
2765            };
2766            i++;
2767        }    
2768    
2769        // ToF    
2770        if(TOF) {
2771            pam_tree_clone[i] = new TTree("ToF","PAMELA ToF level2 data ");
2772            pam_tree_clone[i]->Branch("ToFLevel2","ToFLevel2", GetPointerTo("ToFLevel2"));
2773            cout << "ToF          : branch ToFLevel2"<<endl;
2774            i++;
2775        };
2776        // Trigger
2777        if(TRG) {
2778            pam_tree_clone[i] = new TTree("Trigger","PAMELA trigger level2 data ");
2779            pam_tree_clone[i]->Branch("TrigLevel2","TrigLevel2", GetPointerTo("TrigLevel2"));
2780            cout << "Trigger      : branch TrigLevel2"<<endl;
2781            i++;
2782        };
2783        // S4
2784        if(S4) {
2785            pam_tree_clone[i] = new TTree("S4","PAMELA S4 level2 data ");  
2786            pam_tree_clone[i]->Branch("S4Level2","S4Level2", GetPointerTo("S4Level2"));
2787            cout << "S4           : branch S4Level2"<<endl;
2788            i++;
2789        };
2790        // Neutron Detector
2791        if(ND) {
2792            pam_tree_clone[i] = new TTree("NeutronD","PAMELA neutron detector level2 data ");      
2793            pam_tree_clone[i]->Branch("NDLevel2","NDLevel2", GetPointerTo("NDLevel2"));
2794            cout << "NeutronD     : branch NDLevel2"<<endl;
2795            i++;
2796        };
2797        // Anticounters
2798        if(AC) {
2799            pam_tree_clone[i] = new TTree("Anticounter","PAMELA anticounter detector level2 data ");        
2800            pam_tree_clone[i]->Branch("AcLevel2","AcLevel2", GetPointerTo("AcLevel2"));
2801            cout << "Anticounter  : branch AcLevel2"<<endl;
2802            i++;
2803        };
2804        // OrbitalInfo
2805        if(ORB) {
2806            pam_tree_clone[i] = new TTree("OrbitalInfo","PAMELA oribital info  ");  
2807            pam_tree_clone[i]->Branch("OrbitalInfo","OrbitalInfo", GetPointerTo("OrbitalInfo"));
2808            cout << "OrbitalInfo  : branch OrbitalInfo"<<endl;
2809            i++;
2810        };
2811        cout << "+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+" <<endl;
2812    
2813    }
2814    
2815    
2816    /**
2817     * Fill tree (created with CreateCloneTrees)
2818     *
2819     */
2820    //void PamLevel2::FillNewPamTree(TTree *T){
2821    void PamLevel2::FillCloneTrees(){
2822        
2823    //    cout << "PamLevel2::FillCloneTrees()" << irunentry << endl;
2824    
2825        for(Int_t i=0; i<NCLONES; i++){
2826            if(pam_tree_clone[i])pam_tree_clone[i]->Fill();
2827        }
2828        if(sel_tree_clone)sel_tree_clone->Fill();    
2829    
2830    }
2831    
2832    
2833    TTree* PamLevel2::GetCloneTree(TString name){
2834    
2835        for(Int_t i=0; i<NCLONES; i++){    
2836            if(pam_tree_clone[i]){
2837                TString na = pam_tree_clone[i]->GetName();
2838                if(!name.CompareTo(na))return pam_tree_clone[i];
2839            };
2840        }
2841        if(run_tree_clone){
2842            TString na = run_tree_clone->GetName();
2843            if(!name.CompareTo(na))return run_tree_clone;
2844        }
2845        if(sel_tree_clone){
2846            TString na = sel_tree_clone->GetName();
2847            if(!name.CompareTo(na))return sel_tree_clone;
2848        }
2849        return NULL;
2850    
2851    }
2852    void PamLevel2::WriteCloneTrees(){
2853        cout << "+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+" <<endl;
2854        cout << "Write clones of PAMELA trees "<<endl;
2855        cout << run_tree_clone->GetName()<<endl;    
2856        run_tree_clone->Write();
2857        cout << sel_tree_clone->GetName()<<endl;    
2858        sel_tree_clone->Write();
2859        for(Int_t i=0; i<NCLONES; i++){    
2860            if(pam_tree_clone[i]){
2861                cout << pam_tree_clone[i]->GetName()<<endl;
2862                pam_tree_clone[i]->Write();
2863            };
2864        }
2865        cout << "+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+" <<endl;
2866    
2867    }
2868    
2869    /**
2870     * Method to get level2-trees entry, the corresponding run entry and (if required) the level0 entry.
2871     */
2872    //Int_t PamLevel2::GetEntry(Int_t iee){
2873    Int_t PamLevel2::GetEntry(Long64_t iee){
2874        
2875        if(!pam_tree){
2876            cout << " Int_t PamLevel2::GetEntry(Int_t) -- ERROR -- level2 trees not loaded"<<endl;
2877            return 0;
2878        }
2879    
2880    
2881        //
2882        // 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...
2883        //
2884        //    if(!run_tree ){
2885        //  cout << " Int_t PamLevel2::GetEntry(Int_t) -- ERROR -- run tree not loeaded"<<endl;
2886        //  return 0;
2887        //    }
2888    
2889        Long64_t ii=0;
2890        //-------------------------------
2891        ii = iee;
2892        if( !pam_tree->GetEntry(ii) ){      
2893            cout << " Int_t PamLevel2::GetEntry(Int_t) -- ERROR -- error reading pam tree"<<endl;
2894            return 0;
2895        }
2896        //
2897        // ... 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.
2898        // 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
2899        // a problem if you don't check the return code of getentry.
2900        //
2901        if(!run_tree ){
2902            if ( TRK0 || CAL0 || TOF0 || RUN ) { //forse cosi` va bene per tornare 1?
2903                cout << " Int_t PamLevel2::GetEntry(Int_t) -- ERROR -- run tree not loaded"<<endl;
2904                return 0;
2905            }  else {
2906                return 1; //cosi` se non c'e` run esce qua...
2907            }
2908        }
2909    
2910        //-------------------------------
2911        ii = iee;
2912    //    Bool_t UPDATED = UpdateRunInfo(run_tree,ii);
2913    //    Bool_t UPDATED = UpdateRunInfo(ii);
2914        UpdateRunInfo(ii);
2915        if(SELLI==0)irunentry = iee-runfirstentry;
2916    //    if(UPDATED && run_tree_clone)run_tree_clone->Fill();
2917      
2918    //    cout << "PamLevel2::GetEntry("<<iee<<") "<<irun<<" "<<runfirstentry<<" "<<irunentry<<endl;
2919    
2920    //    cout << " irunentry "<<irunentry << endl;
2921    
2922        if( TRK0 || CAL0 || TOF0 ){
2923            if( !GetYodaEntry( ) ){
2924                cout << " Int_t PamLevel2::GetEntry(Int_t) -- ERROR -- error reading level0 tree"<<endl;
2925                return 0;
2926            }
2927        }
2928        return 1;
2929    
2930    }
2931    
2932    /**
2933     * Method to retrieve the level0 tree (YODA tree) that contains the current event.
2934     * Given the run ID (...), if needed it query the DB and load the proper file.
2935     * @return Pointer to the tree
2936     */
2937    
2938    
2939    TTree* PamLevel2::GetYodaTree( ){
2940        
2941    //    cout << "TTree* PamLevel2::GetYodaTree( )"<<endl;
2942        //===================================
2943        // check if iroot has changed
2944        //===================================
2945        if( irun<0 ){
2946            cout << "PamLevel2::GetYodaTree() -- ERROR -- irun = "<<irun<<endl;
2947    //      cout << "In order to use this method you have to first load the RunInfo tree "<<endl;
2948            return NULL;
2949        }    
2950        Int_t irootnew = run_obj->ID_ROOT_L0;
2951    //     cout << "iroot    "<<iroot<<endl;
2952    //     cout << "irootnew "<<irootnew<<endl;
2953    
2954        //===================================
2955        // load the level0 file
2956        // (if not already loaded)
2957        //===================================
2958        if( iroot != irootnew || !l0_tree){
2959            iroot = irootnew;
2960            //===================================
2961            // open the DB connection
2962            // (if not already opened)
2963            //===================================
2964            if(!dbc || (dbc && !dbc->IsConnected())){
2965                cout<<"~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~"<<endl;
2966                cout<<"Connecting to DB"<<endl;
2967                cout<<"HOST "<<host<<endl;
2968                cout<<"USER "<<user<<endl;
2969                cout<<"PSW  "<<psw<<endl;
2970                dbc = TSQLServer::Connect(host.Data(),user.Data(),psw.Data());
2971                if( !dbc )return NULL;
2972                if( !dbc->IsConnected() )return NULL;    
2973    //          cout<<"...done"<<endl;
2974                cout<<"~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~"<<endl;
2975            }else{
2976    //          cout<<"DB already connected"<<endl;
2977            }
2978            GL_ROOT glroot = GL_ROOT();
2979            if( glroot.Query_GL_ROOT(iroot,dbc) ){
2980                cout << "TTree* PamLevel2::GetYodaTree( ) -- ERROR -- level0 file iroot = "<<iroot<< " does not exists"<<endl;
2981                return NULL;
2982            };
2983            TString filename = glroot.PATH + glroot.NAME;
2984            if(l0_file){
2985                l0_file->Close();
2986                l0_file->Delete();
2987            }
2988            cout << "Opening LEVEL0 file: "<< filename << endl;
2989            FileStat_t t;
2990            if( gSystem->GetPathInfo(filename.Data(),t) ){
2991                cout << " PamLevel2::GetYodaTree() -- ERROR opening file "<<endl;
2992                return NULL;
2993            }
2994            l0_file = new TFile(filename);
2995            if( !l0_file )return NULL;
2996            l0_tree = (TTree*)l0_file->Get("Physics");
2997            if(!h0_obj)h0_obj = new EventHeader();
2998            l0_tree->SetBranchAddress("Header" ,&h0_obj);
2999            prevshift = 0;
3000            //---------------------------------------------------
3001            // TRACKER:
3002            if(TRK0){
3003                if(!trk0_obj){
3004                    trk0_obj = new TrkLevel0();
3005                    trk0_obj->Set();
3006                };
3007                l0_tree->SetBranchAddress("Tracker" ,trk0_obj->GetPointerToTrackerEvent());
3008            }
3009            //---------------------------------------------------
3010            // CALORIMETER:
3011            if(CAL0){
3012                if(!calo0_obj){
3013                    calo0_obj = new CaloLevel0();
3014                    calo0_obj->Set();
3015                };
3016                l0_tree->SetBranchAddress("Calorimeter" ,calo0_obj->GetPointerToCalorimeterEvent());
3017                //      cout << "PamLevel2::GetYodaTree() --- level0 calorimeter not implemented "<<endl;
3018            }
3019            //---------------------------------------------------
3020            // TOF:
3021            if(TOF0){
3022                cout << "PamLevel2::GetYodaTree() --- level0 TOF not implemented "<<endl;
3023            }
3024    
3025        };
3026    
3027        if(!dbc || (dbc && !dbc->IsConnected())){
3028            cout << " TTree* PamLevel2::GetYodaTree( ) -- no DB connected... hai fatto qualche cazzata "<<endl;
3029        }
3030    
3031        if ( TRK0 ){
3032            TrkParams::Load(6);
3033            if( !TrkParams::IsLoaded(6) ){
3034                cout << " TTree* PamLevel2::GetYodaTree( ) -- WARNING -- VK-mask not loaded"<<endl;
3035            };
3036            TrkParams::SetCalib(run_obj,dbc);
3037            TrkParams::LoadCalib( );
3038            if( !TrkParams::CalibIsLoaded() ){
3039                cout << " TTree* PamLevel2::GetYodaTree( ) -- WARNING -- Calibration not loaded"<<endl;
3040            };
3041        }
3042    
3043    //    cout << l0_tree << endl;
3044        
3045        return l0_tree;
3046    
3047    }
3048    
3049    /**
3050     * Method to retrieve the level0 tree (YODA tree) that contains the current event.
3051     */
3052    Int_t PamLevel2::GetYodaEntry(){
3053    
3054    //    cout << "Int_t PamLevel2::GetYodaEntry()"<<endl;
3055        if(!GetYodaTree())return 0;
3056        
3057        // patch
3058        if( irunentry < 0){
3059    //      cout << "Int_t PamLevel2::GetYodaEntry() -- ATTENZIONE -- irunentry negativo?!?! "<<(Int_t)irunentry<<endl;
3060            irunentry=0LL;
3061        }
3062        //  ---------------------------------
3063        //  if file is NOT a preselected file
3064        //  ---------------------------------
3065        Long64_t quellagiusta = irunentry + (Long64_t)(run_obj->EV_FROM);
3066    //     cout << " irunentry "<<irunentry<<endl;
3067    //     cout << " EV_FROM "<<run_obj->EV_FROM<<endl;
3068    //     cout << " quellagiusta = irunentry + EV_FROM "<< quellagiusta << endl;
3069    
3070    //     cout << " irun "<< irun << " irunentry "<< irunentry<<" run_obj->EV_FROM "<<run_obj->EV_FROM <<" quella giusta "<<quellagiusta << endl;
3071    //     cout << " iroot "<<iroot<<" run_obj->ID_ROOT_L0 "<<run_obj->ID_ROOT_L0<<endl;
3072    //     cout << " time "<< GetOrbitalInfo()->absTime << endl;
3073    //     cout << " trk_calib_used "<<run_obj->TRK_CALIB_USED<< endl;
3074        
3075        if( !GetOrbitalInfo() ){
3076            cout << "Int_t PamLevel2::GetYodaEntry() -- ERROR -- missing OrbitalInfo "<<endl;
3077            return 0;
3078        }
3079        if( GetOrbitalInfo()->OBT==0 && GetOrbitalInfo()->pkt_num==0 ){
3080            cout << "Int_t PamLevel2::GetYodaEntry() -- ERROR -- level2 event corrupted ?? "<<endl;
3081            return 0;
3082        }
3083    
3084        // ---------------------------------------------------------------------
3085        // ATTENTION!!!
3086        // If data are not pre-processed with cleaner, the level0 tree may contain
3087        // spurious nested physics packets.
3088        // The GL_RUN variables EV_FROM, EV_TO and NEVENTS counts also these entries
3089        // while level2 tree DOES NOT!!
3090        // This means that "quellagiusta" in these cases is not correct.
3091        // In order to retrieve the correct level0 event, I implemented a check
3092        // of the OBT and pkt-number. In case of mismatch, the level0 entry number
3093        // is shift forward until when the packets match.
3094        // ---------------------------------------------------------------------
3095        Int_t answer = 0;
3096        Int_t shift =0;
3097        //    printf(" siamo qui %i %i \n",shift,prevshift);
3098        Int_t maxshift = 100;
3099        do{
3100            if(shift>0){    
3101                cout << " PKTNUM  L2 --- "<< GetOrbitalInfo()->pkt_num << " --- L0 --- "<< GetEventHeader()->GetPscuHeader()->GetCounter()<<endl;
3102                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;
3103                cout << "         L2 <--> L0 mismatch ( irun "<<irun<<" irunentry "<<irunentry<<" shift "<<shift<<" prevshift "<<prevshift<<" )"<<endl;
3104            }
3105            answer = l0_tree->GetEntry(quellagiusta+(Long64_t)shift+(Long64_t)prevshift);
3106            shift++;
3107            if( !GetEventHeader() ){
3108                cout << "Int_t PamLevel2::GetYodaEntry() -- ERROR -- missing EventHeader "<<endl;
3109                return 0;
3110            }
3111    //      cout << "PKTNUM "<<shift<<" ==  L2 --- "<< GetOrbitalInfo()->pkt_num << " --- L0 --- "<< GetEventHeader()->GetPscuHeader()->GetCounter()<<endl;
3112    //      cout << " L2 --- "<< GetOrbitalInfo()->OBT << " --- L0 --- "<< GetEventHeader()->GetPscuHeader()->GetOrbitalTime()<<endl;
3113    //      if( (quellagiusta+shift) == l0_tree->GetEntries()+1 )cout << ">>> end of level0 tree <<<"<<endl;
3114    //      cout << " GetOrbitalInfo()->OBT "<< GetOrbitalInfo()->OBT << endl;
3115    //      cout << " GetEventHeader()->GetPscuHeader()->GetOrbitalTime() "<< GetEventHeader()->GetPscuHeader()->GetOrbitalTime() << endl;
3116    //      cout << " GetOrbitalInfo()->pkt_num "<< GetOrbitalInfo()->pkt_num << endl;
3117    //      cout << " GetEventHeader()->GetPscuHeader()->GetCounter() "<< GetEventHeader()->GetPscuHeader()->GetCounter() << endl;
3118    //      printf(" IDRUN %u \n",GetRunInfo()->ID);
3119    //
3120            if ( prevshift != 0 && (quellagiusta+(Long64_t)shift) == GetYodaTree()->GetEntries() ){
3121              prevshift = 0;
3122              shift = -1;
3123            };
3124    
3125        }while( ( GetOrbitalInfo()->OBT != (UInt_t)(GetEventHeader()->GetPscuHeader()->GetOrbitalTime()) || GetOrbitalInfo()->pkt_num != (UInt_t)(GetEventHeader()->GetPscuHeader()->GetCounter())) && (quellagiusta+(Long64_t)shift) < GetYodaTree()->GetEntries() && shift < maxshift);
3126    
3127        if ( (quellagiusta+(Long64_t)shift) == GetYodaTree()->GetEntries() || shift == maxshift ) {
3128            cout << " Big trouble here, no such event in Level0 data! " <<endl;
3129            return 0;
3130        }
3131    //    cout << "LA ENTRY GIUSTA E`: "<<quellagiusta<<" (spero...)"<<endl;
3132    //    return GetYodaTree()->GetEntry(quellagiusta);
3133        if ( shift > 1 ) prevshift += (shift-1);
3134        
3135        return answer;
3136        
3137    }
3138    
3139    

Legend:
Removed from v.1.1  
changed lines
  Added in v.1.38

  ViewVC Help
Powered by ViewVC 1.1.23