/[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.3 by pam-fi, Mon Aug 7 08:24:49 2006 UTC
# Line 24  PamTrack::PamTrack(TrkTrack* t, CaloTrkV Line 24  PamTrack::PamTrack(TrkTrack* t, CaloTrkV
24      tof_track  = this->ToFTrkVar::GetToFTrkVar();      tof_track  = this->ToFTrkVar::GetToFTrkVar();
25      if(t) *trk_track  = *t;      if(t) *trk_track  = *t;
26      if(c) *calo_track = *c;      if(c) *calo_track = *c;
27      if(o) *tof_track  = *o;      if(o) *tof_track  = *o;    
28  };  };
29    
30  //--------------------------------------  //--------------------------------------
# Line 35  PamTrack::PamTrack(TrkTrack* t, CaloTrkV Line 35  PamTrack::PamTrack(TrkTrack* t, CaloTrkV
35   * Constructor   * Constructor
36   */   */
37  PamLevel2::PamLevel2(){  PamLevel2::PamLevel2(){
38            
39      trk_obj  = this->TrkLevel2::GetTrkLevel2();      trk_obj  = this->TrkLevel2::GetTrkLevel2();
40      calo_obj = this->CaloLevel2::GetCaloLevel2();      calo_obj = this->CaloLevel2::GetCaloLevel2();
41      tof_obj  = this->ToFLevel2::GetToFLevel2();      tof_obj  = this->ToFLevel2::GetToFLevel2();
# Line 42  PamLevel2::PamLevel2(){ Line 43  PamLevel2::PamLevel2(){
43      s4_obj   = this->S4Level2::GetS4Level2();      s4_obj   = this->S4Level2::GetS4Level2();
44      nd_obj   = this->NDLevel2::GetNDLevel2();      nd_obj   = this->NDLevel2::GetNDLevel2();
45      ac_obj   = this->AcLevel2::GetAcLevel2();      ac_obj   = this->AcLevel2::GetAcLevel2();
46            orb_obj  = this->OrbitalInfo::GetOrbitalInfo();
47    
48            sorted_tracks = new TRefArray();
49    
50            CAL = true;
51            TRK = true;
52            TRG = true;
53            TOF = true;
54            S4  = true;
55            ND  = true;
56            AC  = true;
57            ORB = true;
58    
59    };
60    /**
61     * Destructor
62     */
63    PamLevel2::~PamLevel2(){
64                    
65            TrkLevel2::Clear();
66            CaloLevel2::Clear();
67            ToFLevel2::Clear();
68            TrigLevel2::Clear();
69            S4Level2::Clear();
70            NDLevel2::Clear();
71            AcLevel2::Clear();
72            OrbitalInfo::Clear();
73            
74            delete sorted_tracks;
75            
76    };
77    /**
78     * Clear the event
79     */
80    void PamLevel2::Clear(){
81            
82            
83            sorted_tracks->Delete();
84            
85  };  };
86    
87    
88    //--------------------------------------
89    //
90    //
91    //--------------------------------------
92    /**
93     * Retrieves the calorimeter track matching the seqno-th tracker stored track.
94     * (If seqno = -1 retrieves the self-trigger calorimeter track)
95     */
96     CaloTrkVar *PamLevel2::GetCaloStoredTrack(int seqno){
97            
98             CaloTrkVar *c = 0;
99             Int_t it_calo=0;
100            
101             do{
102                     c = CaloLevel2::GetCaloTrkVar(it_calo);
103                     it_calo++;
104             } while( seqno != c->trkseqno && it_calo < CaloLevel2::ntrk());        
105            
106             if(seqno != c->trkseqno){
107                     c = 0;
108                     if(seqno!=-1)cout << "PamLevel2::GetCaloStoredTrack(int) : requested tracker SeqNo "<< seqno <<" does not match Calorimeter stored tracks"<<endl;
109             };
110             return c;
111    
112     };
113    //--------------------------------------
114     //
115     //
116    //--------------------------------------
117    /**
118      * Retrieves the ToF track matching the seqno-th tracker stored track.
119      * (If seqno = -1 retrieves the tracker-independent tof track)
120     */
121     ToFTrkVar *PamLevel2::GetToFStoredTrack(int seqno){
122            
123             ToFTrkVar *c = 0;
124             Int_t it_tof=0;
125            
126             do{
127                     c = ToFLevel2::GetToFTrkVar(it_tof);
128                     it_tof++;
129             } while( seqno != c->trkseqno && it_tof < ToFLevel2::ntrk());  
130            
131             if(seqno != c->trkseqno){
132                     c = 0;
133                     if(seqno!=-1)cout << "PamLevel2::GetToFStoredTrack(int) : requested tracker SeqNo "<< seqno <<" does not match ToF stored tracks"<<endl;
134             };
135             return c;
136    
137     };
138    
139    //--------------------------------------
140     //
141     //
142    //--------------------------------------
143    /**
144     * Give the pamela track associated to a tracker track, retrieving related calorimeter and tof track information.
145     */
146     PamTrack* PamLevel2::GetPamTrackAlong(TrkTrack* t){
147            
148             CaloTrkVar *c = 0;
149             ToFTrkVar  *o = 0;
150            
151             c = GetCaloStoredTrack(t->GetSeqNo());
152             o = GetToFStoredTrack(t->GetSeqNo());
153            
154    //    if(t && c && o)track = new PamTrack(t,c,o);
155             PamTrack *track = new PamTrack(t,c,o);
156            
157             return track;
158    
159     };
160  //--------------------------------------  //--------------------------------------
161  //  //
162  //  //
# Line 55  PamLevel2::PamLevel2(){ Line 169  PamLevel2::PamLevel2(){
169    
170  PamTrack* PamLevel2::GetStoredTrack(Int_t itrk){  PamTrack* PamLevel2::GetStoredTrack(Int_t itrk){
171            
172      // retrieve itrk-th tracker stored track          PamTrack *track = 0;
173      TrkTrack   *t = 0;  
174      CaloTrkVar *c = 0;          if( itrk >=0 && itrk < TrkLevel2::ntrk() ){
175      ToFTrkVar  *o = 0;                  
176                    TrkTrack *t = TrkLevel2::GetStoredTrack(itrk);
177      if( itrk >=0 && itrk < TrkLevel2::ntrk() ){                  track = GetPamTrackAlong(t);
178          t = TrkLevel2::GetStoredTrack(itrk);              
179          c = CaloLevel2::GetCaloTrkVar(t->GetSeqNo());          }else{
180          o = ToFLevel2::GetToFTrkVar(t->GetSeqNo());                  cout << "PamLevel2::GetStoredTrack(int) : tracker track SeqNo "<< itrk <<" does not exist (GetNTracks() = "<<TrkLevel2::GetNTracks()<<")"<<endl;
181      };          };
182      //  retrieve related ToF-track  
183                return track;
184      // hence create a "PamTrack" object          
     PamTrack *track = 0;  
     if(t && c && o)track = new PamTrack(t,c,o);  
     return track;  
       
185  }  }
186  //--------------------------------------  //--------------------------------------
187  //  //
188  //  //
189  //--------------------------------------  //--------------------------------------
190  /**  /**
191   * 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.
192   * 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...).
193   * This method overrides TrkLevel2::GetTracks(), where sorting is done by decreasing number of fit points and increasing chi^2.   * Sorting cryteria:
194   * 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
195     * CAL: lower Y spatial residual on the first calorimeter plane
196     * TOF: bigger numebr of hit PMTs along the track, on S12 S21 S32 (where paddles are along the Y axis).
197     * The default sorting cryterium is "TOF+CAL".
198     *
199     * 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).
200   */   */
201  TClonesArray *PamLevel2::GetTracks(){  void PamLevel2::SortTracks(TString how){
   
     // define new array for sorted tracker tracks  
     TClonesArray *aa = new TClonesArray("TrkTrack");  
     TClonesArray &sorted = *aa;  
202    
203            sorted_tracks->Delete(); //temporaneo
204      // loop over the tracks sorted by the tracker      // loop over the tracks sorted by the tracker
205      for(Int_t i=0; i < TrkLevel2::GetNTracks(); i++){          Bool_t use_TRK = how.Contains("TRK", TString::kIgnoreCase);
206            Bool_t use_CAL = how.Contains("CAL", TString::kIgnoreCase);
207            Bool_t use_TOF = how.Contains("TOF", TString::kIgnoreCase);
208            
209            
210            for(Int_t i=0; i < TrkLevel2::GetNTracks(); i++){
211    
212          TrkTrack *ts=0;                  TrkTrack *ts = 0;
213    
214          // get tracker tracks          // get tracker tracks
215          TrkTrack   *tp = TrkLevel2::GetTrack(i);                    //tracker                  TrkTrack   *tp = TrkLevel2::GetTrack(i);                    //tracker
216          CaloTrkVar *cp = CaloLevel2::GetCaloTrkVar(tp->GetSeqNo()); //calorimeter                  CaloTrkVar *cp = GetCaloStoredTrack(tp->GetSeqNo());
217          ToFTrkVar  *op = ToFLevel2::GetToFTrkVar(tp->GetSeqNo());   //tof                  ToFTrkVar  *op = GetToFStoredTrack(tp->GetSeqNo());
218    
219          // if track has an image, check image selection          // if track has an image, check image selection
220          if(tp->HasImage()){                  if(tp->HasImage()){
221              TrkTrack   *ti = TrkLevel2::GetTrackImage(i);              //tracker (image)                          
222              CaloTrkVar *ci = CaloLevel2::GetCaloTrkVar(ti->GetSeqNo());//calorimeter (image)                          TrkTrack   *ti = TrkLevel2::GetTrackImage(i);              //tracker (image)
223              ToFTrkVar  *oi = ToFLevel2::GetToFTrkVar(ti->GetSeqNo());  //tof (image)                          CaloTrkVar *ci = GetCaloStoredTrack(ti->GetSeqNo());
224                            ToFTrkVar  *oi = GetToFStoredTrack(ti->GetSeqNo());
225    
226              //assign starting scores              //assign starting scores
227              Int_t tp_score = 1;                          Int_t tp_score = 0;  //main track sorted by the tracker
228              Int_t ti_score = 0;                          Int_t ti_score = 0;  //image track
229                                                    
230              // ------------------------              // ------------------------
231              // calorimeter check              // calorimeter check
232              // ------------------------              // ------------------------
233              if(                  // check the Y spatial residual on the first calorimeter plane
234                        npcfit[1] > 3   &&  //no. of fit planes on Y view                  // (cut on calorimeter variables from Emiliano)
235  //              varcfit[1] < 50.&&  //fit variance on Y view                          if(
236                        true){                                  use_CAL            &&
237                                    npcfit[1] > 15     &&   //no. of fit planes on Y view
238                  Float_t resy_p = cp->tbar[0][1] - cbar[0][1];if(resy_p < 0)resy_p= - resy_p;                                  varcfit[1] < 1000. &&  //fit variance on Y view
239                  Float_t resy_i = ci->tbar[0][1] - cbar[0][1];if(resy_i < 0)resy_i= - resy_i;                                  true){
240    
241                  if(resy_p <= resy_i) tp_score++;                                  Float_t resy_p = cp->tbar[0][1] - cbar[0][1]; if(resy_p < 0)resy_p= - resy_p;
242                  else                 ti_score++;                                  Float_t resy_i = ci->tbar[0][1] - cbar[0][1]; if(resy_i < 0)resy_i= - resy_i;
243                        };  
244                                    if(resy_p <= resy_i) tp_score++;
245                                    else                 ti_score++;
246                            };
247              // ------------------------              // ------------------------
248              // TOF check              // TOF check
249              // ------------------------                  // ------------------------    
250                    // check the number of hit pmts along the track
251                    // on S12 S21 and S32, where paddles are parallel to Y axis
252                            Int_t nphit_p =0;
253                            Int_t nphit_i =0;
254                                    
255                            for (Int_t ih=0; ih < op->npmtadc; ih++){
256                                    Int_t pl = GetPlaneIndex( (op->pmtadc).At(ih) );
257                                    if(pl == 1 || pl == 2 || pl == 5)nphit_p++;
258                            };
259                                                    
260                            for (Int_t ih=0; ih < oi->npmtadc; ih++){
261                                    Int_t pl = GetPlaneIndex( (oi->pmtadc).At(ih) );
262                                    if(pl == 1 || pl == 2 || pl == 5)nphit_i++;    
263                            };
264                            
265                            if(
266                                    use_TOF            &&
267                                    (nphit_p+nphit_i) !=0 &&        
268                                    true){
269    
270                                    if( nphit_p >= nphit_i) tp_score++;
271                                    else                                              ti_score++;
272                            };
273                
274                            if(tp_score == ti_score) use_TRK = true;
275                    // ------------------------
276                // tracker check
277                // ------------------------
278                    // chi**2 difference is not always large enough to distinguish among
279                    // the real track and its image.
280                    // Tracker check will be applied always when calorimeter and tof information is ambiguous.
281                            if(use_TRK){
282                                    if(      tp->chi2 > 0 && tp->chi2 < ti->chi2 ) tp_score++ ;
283                                    else if( ti->chi2 > 0 && ti->chi2 < tp->chi2 ) ti_score++ ;
284                            };
285    
286              // ------------------------              // ------------------------
287              // the winner is....              // the winner is....
288              // ------------------------                  // ------------------------    
289                        if(tp_score > ti_score) ts = tp;//the track sorted by the tracker!!                          if      (tp_score > ti_score) ts = tp;//the track sorted by the tracker!!
290                        else                    ts = ti;//its image!!                          else if (tp_score < ti_score) ts = ti;//its image!!
291                            else {
292          }else{                                  ts = tp;
293              ts = tp;  //                              cout << "Warning - track image ambiguity not solved" << endl;
294    //                              cout << ts->GetNtot() << " "<< ts->chi2 << " " << npcfit[1] << " "<< nphit_p << endl;
295                            };
296    
297                    }else{
298                            ts = tp;
299                    };
300                            
301                    sorted_tracks->Add(ts);//save the track in the sorted array
302                    
303          };          };
           
         new(sorted[i]) TrkTrack(*ts); //save the track in the sorted array  
     };  
304    
     return aa;  
305  };  };
306  //--------------------------------------  //--------------------------------------
307  //  //
308  //  //
309  //--------------------------------------  //--------------------------------------
310  /**  /**
311     * This method overrides TrkLevel2::GetTracks(), where sorting is done by decreasing number of fit points and increasing chi^2.
312     * PamLevel2::GetTracks() keeps the same track order given by TrkLevel2::GetTracks(), but checks image selection by using calorimeter and ToF tracking information.
313     */
314     TRefArray *PamLevel2::GetTracks(){
315            SortTracks("+CAL+TOF");
316            return sorted_tracks;
317     };
318    //--------------------------------------
319     //
320     //
321    //--------------------------------------
322    /**
323   * Retrieves the it-th Pamela "physical" track.   * Retrieves the it-th Pamela "physical" track.
324   * It override TrkLevel2::GetTrack(int it).   * It override TrkLevel2::GetTrack(int it).
325   * @param it Track number, ranging from 0 to GetNTracks().   * @param it Track number, ranging from 0 to GetNTracks().
326   */   */
327  PamTrack *PamLevel2::GetTrack(int it){  PamTrack *PamLevel2::GetTrack(int it){
328            
329    //  *-*-*-*-*-*-*
330            SortTracks("+CAL+TOF");
331    //  *-*-*-*-*-*-*
332    
333      TClonesArray *aa = this->GetTracks();          PamTrack *track = 0;
334      TClonesArray &sorted = *aa;          
335            if( it >=0 && it < TrkLevel2::GetNTracks() ){
336      TrkTrack   *t = 0;                  TrkTrack   *t = (TrkTrack*)sorted_tracks->At(it);
337      CaloTrkVar *c = 0;                  track = GetPamTrackAlong(t);
338      ToFTrkVar  *o = 0;          }else{
339                    cout << "PamLevel2::GetTrack(int) : tracker track SeqNo "<< it <<" does not exist (GetNTracks() = "<<TrkLevel2::GetNTracks()<<")"<<endl;
340      if( it >=0 && it < TrkLevel2::GetNTracks() ){          };
341          t = (TrkTrack*)sorted[it];          
342          c = CaloLevel2::GetCaloTrkVar(t->GetSeqNo());          return track;
         o = ToFLevel2::GetToFTrkVar(t->GetSeqNo());  
     };  
       
     // hence create a "PamTrack" object  
     PamTrack *track = 0;  
     if(t && c && o)track = new PamTrack(t,c,o);  
     return track;  
343    
344  };  };
345    
346  //--------------------------------------  //--------------------------------------
347  //  //
348  //  //
349  //--------------------------------------  //--------------------------------------
350  /**  /**
351   * 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().
352   * @param it Track number, ranging from 0 to GetNTracks().   * @param it Track number, ranging from 0 to GetNTracks().
353   */   */
354  PamTrack *PamLevel2::GetTrackImage(int it){  PamTrack *PamLevel2::GetTrackImage(int it){
355    
356      TClonesArray *aa = this->GetTracks();          SortTracks("+CAL+TOF");
357      TClonesArray &sorted = *aa;          
358            PamTrack *image = 0;
359    
360      TrkTrack   *t = 0;          if( it >=0 && it < TrkLevel2::GetNTracks() ){
361      CaloTrkVar *c = 0;                  TrkTrack *temp = (TrkTrack*)sorted_tracks->At(it);
362      ToFTrkVar  *o = 0;                  if( temp->HasImage() ){
363                            TrkTrack   *t = TrkLevel2::GetStoredTrack(temp->GetImageSeqNo());
364      if( it >=0 && it < TrkLevel2::GetNTracks() ){                          image = GetPamTrackAlong(t);
365          TrkTrack *temp = (TrkTrack*)sorted[it];                  }else{
366          if( temp->HasImage() ){                          cout <<"PamLevel2::GetTrackImage(int) : Track SeqNo "<<it<<" does not have image"<<endl;
367              t = TrkLevel2::GetStoredTrack(temp->GetImageSeqNo());                  };
368              c = CaloLevel2::GetCaloTrkVar(temp->GetImageSeqNo());          }else{
369              o = ToFLevel2::GetToFTrkVar(temp->GetImageSeqNo());                  cout << "PamLevel2::GetTrackImage(int) : Tracker track SeqNo "<< it <<" does not exist (GetNTracks() = "<<TrkLevel2::GetNTracks()<<")"<<endl;
370          };          };
     };  
       
     // hence create a "PamTrack" object  
     PamTrack *image = 0;  
     if(t && c && o)image = new PamTrack(t,c,o);  
     return image;  
371            
372            return image;
373  }  }
374    
375  //--------------------------------------  //--------------------------------------
# Line 209  PamTrack *PamLevel2::GetTrackImage(int i Line 377  PamTrack *PamLevel2::GetTrackImage(int i
377  //  //
378  //--------------------------------------  //--------------------------------------
379  /**  /**
380   * Get the Pamela detector trees and make them friends.   * Get the Pamela detector trees in a single file and make them friends.
381     * @param f TFile pointer
382     * @param detlist String to select trees to be included
383     * @return Pointer to a TTree
384   */   */
385  TTree *PamLevel2::LoadPamTrees(TFile *f){  TTree *PamLevel2::LoadPamTrees(TFile *f, TString detlist="+ALL"){
386            
387      TTree *Tout =0;      TTree *Tout =0;
388                
389            SetWhichTrees(detlist);
390            
391      // Tracker      // Tracker
392      TTree *T = (TTree*)f->Get("Tracker");      TTree *T = (TTree*)f->Get("Tracker");
393      if(T) {      if(T && TRK) {
394          T->SetBranchAddress("TrkLevel2", GetPointerToTrk());          T->SetBranchAddress("TrkLevel2", GetPointerToTrk());
395          cout << "Tracker      : set branch address TrkLevel2"<<endl;          cout << "Tracker      : set branch address TrkLevel2"<<endl;
396          if(!Tout)Tout=T;          if(!Tout)Tout=T;
# Line 226  TTree *PamLevel2::LoadPamTrees(TFile *f) Line 399  TTree *PamLevel2::LoadPamTrees(TFile *f)
399      };      };
400      // Calorimeter      // Calorimeter
401      TTree *C = (TTree*)f->Get("Calorimeter");      TTree *C = (TTree*)f->Get("Calorimeter");
402      if(C) {      if(C && CAL) {
403          C->SetBranchAddress("CaloLevel2", GetPointerToCalo());          C->SetBranchAddress("CaloLevel2", GetPointerToCalo());
404          cout << "Calorimeter  : set branch address CaloLevel2"<<endl;          cout << "Calorimeter  : set branch address CaloLevel2"<<endl;
405          if(!Tout)Tout=C;          if(!Tout)Tout=C;
# Line 236  TTree *PamLevel2::LoadPamTrees(TFile *f) Line 409  TTree *PamLevel2::LoadPamTrees(TFile *f)
409      };      };
410      // ToF          // ToF    
411      TTree *O = (TTree*)f->Get("ToF");      TTree *O = (TTree*)f->Get("ToF");
412      if(O) {      if(O && TOF) {
413          O->SetBranchAddress("ToFLevel2", GetPointerToToF());          O->SetBranchAddress("ToFLevel2", GetPointerToToF());
414          cout << "ToF          : set branch address ToFLevel2"<<endl;          cout << "ToF          : set branch address ToFLevel2"<<endl;
415          if(!Tout)Tout=O;          if(!Tout)Tout=O;
# Line 246  TTree *PamLevel2::LoadPamTrees(TFile *f) Line 419  TTree *PamLevel2::LoadPamTrees(TFile *f)
419      };      };
420      // Trigger      // Trigger
421      TTree *R = (TTree*)f->Get("Trigger");      TTree *R = (TTree*)f->Get("Trigger");
422      if(R) {      if(R && TRG) {
423          R->SetBranchAddress("TrigLevel2", GetPointerToTrig());          R->SetBranchAddress("TrigLevel2", GetPointerToTrig());
424          cout << "Trigger      : set branch address TrigLevel2"<<endl;          cout << "Trigger      : set branch address TrigLevel2"<<endl;
425          if(!Tout)Tout=O;          if(!Tout)Tout=O;
# Line 256  TTree *PamLevel2::LoadPamTrees(TFile *f) Line 429  TTree *PamLevel2::LoadPamTrees(TFile *f)
429      };      };
430      // S4      // S4
431      TTree *S = (TTree*)f->Get("S4");      TTree *S = (TTree*)f->Get("S4");
432      if(S) {      if(S && S4) {
433          S->SetBranchAddress("S4Level2", GetPointerToS4());          S->SetBranchAddress("S4Level2", GetPointerToS4());
434          cout << "S4           : set branch address S4Level2"<<endl;          cout << "S4           : set branch address S4Level2"<<endl;
435          if(!Tout)Tout=O;          if(!Tout)Tout=O;
# Line 266  TTree *PamLevel2::LoadPamTrees(TFile *f) Line 439  TTree *PamLevel2::LoadPamTrees(TFile *f)
439      };      };
440      // Neutron Detector      // Neutron Detector
441      TTree *N = (TTree*)f->Get("NeutronD");      TTree *N = (TTree*)f->Get("NeutronD");
442      if(N) {      if(N && ND) {
443          N->SetBranchAddress("NDLevel2", GetPointerToND());          N->SetBranchAddress("NDLevel2", GetPointerToND());
444          cout << "NeutronD     : set branch address NDLevel2"<<endl;          cout << "NeutronD     : set branch address NDLevel2"<<endl;
445          if(!Tout)Tout=O;          if(!Tout)Tout=O;
# Line 276  TTree *PamLevel2::LoadPamTrees(TFile *f) Line 449  TTree *PamLevel2::LoadPamTrees(TFile *f)
449      };      };
450      // Anticounters      // Anticounters
451      TTree *A = (TTree*)f->Get("Anticounter");      TTree *A = (TTree*)f->Get("Anticounter");
452      if(A) {      if(A && AC) {
453          A->SetBranchAddress("AcLevel2", GetPointerToAc());          A->SetBranchAddress("AcLevel2", GetPointerToAc());
454          cout << "Anticounter  : set branch address AcLevel2"<<endl;          cout << "Anticounter  : set branch address AcLevel2"<<endl;
455          if(!Tout)Tout=O;          if(!Tout)Tout=O;
# Line 284  TTree *PamLevel2::LoadPamTrees(TFile *f) Line 457  TTree *PamLevel2::LoadPamTrees(TFile *f)
457      }else{      }else{
458          cout << "Anticounter  : missing tree"<<endl;          cout << "Anticounter  : missing tree"<<endl;
459      };      };
460        // Orbital Info
461            TTree *B = (TTree*)f->Get("OrbitalInfo");
462            if(B && ORB) {
463                    B->SetBranchAddress("OrbitalInfo", GetPointerToOrb());
464                    cout << "OrbitalInfo  : set branch address OrbitalInfo"<<endl;
465                    if(!Tout)Tout=O;
466                    else Tout->AddFriend(B);
467            }else{
468                    cout << "OrbitalInfo  : missing tree"<<endl;
469            };
470            
471      return Tout;      return Tout;
472            
473  }  }
474    //--------------------------------------
475    //
476    //
477    //--------------------------------------
478    /**
479     * Get all the Pamela detector trees in a single file and make them friends.
480     * @param f TFile pointer
481     * @return Pointer to a TTree
482     */
483    TTree *PamLevel2::LoadPamTrees(TFile *f){
484        
485            TTree *Tout =0;
486            
487            SetWhichTrees("+ALL");
488            
489        // Tracker
490            TTree *T = (TTree*)f->Get("Tracker");
491            if(T && TRK) {
492                    T->SetBranchAddress("TrkLevel2", GetPointerToTrk());
493                    cout << "Tracker      : set branch address TrkLevel2"<<endl;
494                    if(!Tout)Tout=T;
495            }else{
496                    cout << "Tracker      : missing tree"<<endl;
497            };
498        // Calorimeter
499            TTree *C = (TTree*)f->Get("Calorimeter");
500            if(C && CAL) {
501                    C->SetBranchAddress("CaloLevel2", GetPointerToCalo());
502                    cout << "Calorimeter  : set branch address CaloLevel2"<<endl;
503                    if(!Tout)Tout=C;
504                    else Tout->AddFriend(C);
505            }else{
506                    cout << "Calorimeter  : missing tree"<<endl;
507            };
508        // ToF    
509            TTree *O = (TTree*)f->Get("ToF");
510            if(O && TOF) {
511                    O->SetBranchAddress("ToFLevel2", GetPointerToToF());
512                    cout << "ToF          : set branch address ToFLevel2"<<endl;
513                    if(!Tout)Tout=O;
514                    else Tout->AddFriend(O);
515            }else{
516                    cout << "ToF         : missing tree"<<endl;
517            };
518        // Trigger
519            TTree *R = (TTree*)f->Get("Trigger");
520            if(R && TRG) {
521                    R->SetBranchAddress("TrigLevel2", GetPointerToTrig());
522                    cout << "Trigger      : set branch address TrigLevel2"<<endl;
523                    if(!Tout)Tout=O;
524                    else Tout->AddFriend(R);
525            }else{
526                    cout << "Trigger      : missing tree"<<endl;
527            };
528        // S4
529            TTree *S = (TTree*)f->Get("S4");
530            if(S && S4) {
531                    S->SetBranchAddress("S4Level2", GetPointerToS4());
532                    cout << "S4           : set branch address S4Level2"<<endl;
533                    if(!Tout)Tout=O;
534                    else Tout->AddFriend(S);
535            }else{
536                    cout << "S4           : missing tree"<<endl;
537            };
538        // Neutron Detector
539            TTree *N = (TTree*)f->Get("NeutronD");
540            if(N && ND) {
541                    N->SetBranchAddress("NDLevel2", GetPointerToND());
542                    cout << "NeutronD     : set branch address NDLevel2"<<endl;
543                    if(!Tout)Tout=O;
544                    else Tout->AddFriend(N);
545            }else{
546                    cout << "NeutronD     : missing tree"<<endl;
547            };
548        // Anticounters
549            TTree *A = (TTree*)f->Get("Anticounter");
550            if(A && AC) {
551                    A->SetBranchAddress("AcLevel2", GetPointerToAc());
552                    cout << "Anticounter  : set branch address AcLevel2"<<endl;
553                    if(!Tout)Tout=O;
554                    else Tout->AddFriend(A);
555            }else{
556                    cout << "Anticounter  : missing tree"<<endl;
557            };
558        // Orbital Info
559            TTree *B = (TTree*)f->Get("OrbitalInfo");
560            if(B && ORB) {
561                    B->SetBranchAddress("OrbitalInfo", GetPointerToOrb());
562                    cout << "OrbitalInfo  : set branch address OrbitalInfo"<<endl;
563                    if(!Tout)Tout=O;
564                    else Tout->AddFriend(B);
565            }else{
566                    cout << "OrbitalInfo  : missing tree"<<endl;
567            };
568        
569            return Tout;
570        
571    }
572    //--------------------------------------
573    //
574    //
575    //--------------------------------------
576    /**
577     * Get list of Level2 files.
578     * @param ddir Level2 data directory.
579     * @param flisttxt Name of txt file containing file list.
580     * @return Pointer to a TList of TSystemFiles
581     * If no input file list is given , all the Level2 files inside the directory are processed.
582     */
583    TList*  PamLevel2::GetListOfLevel2Files(TString ddir, TString flisttxt = ""){
584            
585            TString wdir = gSystem->WorkingDirectory();    
586            
587            if(ddir=="")ddir = wdir;
588            TSystemDirectory *datadir = new TSystemDirectory(gSystem->BaseName(ddir),ddir);
589            cout << "Level2 data directory : "<< endl <<  ddir << endl;
590                    
591            TList *contents  = new TList; // create output list
592            contents->SetOwner();
593                    
594            // if no input file list is given:      
595            if ( flisttxt != "" ){
596                    
597                    if( !gSystem->ChangeDirectory(ddir) )return 0;
598            
599                    flisttxt = gSystem->ConcatFileName(wdir,gSystem->BaseName(flisttxt));
600                    
601                    cout <<"Input file list : "<< endl << flisttxt <<endl;
602                    ifstream in;
603                    in.open(flisttxt, ios::in);
604                    while (1) {
605                            TString file;
606                            in >> file;
607                            if (!in.good()) break;
608                            char *fullpath;
609                            if( gSystem->IsFileInIncludePath(file,&fullpath) ){
610    //                              contents->Add(new TSystemDirectory(fullpath,ddir)); // add file to the list
611    //                              contents->Add(new TSystemFile(fullpath,ddir)); // add file to the list
612                                    contents->Add(new TSystemFile(fullpath,gSystem->DirName(fullpath)));// add file to the list
613                            }else{
614                                    cout << "warning! --- File "<<file<<" does not exists"<< endl;
615                            };
616                    };
617                    in.close();
618                    
619            }else{
620            
621                    cout << "No input file list given."<<endl;
622                    cout << "Check for existing root files."<<endl;
623    //              cout << "Warking directory: "<< gSystem->WorkingDirectory()<< endl;
624                    
625                    TList *temp = datadir->GetListOfFiles();
626    //              temp->Print();
627    //              cout << "*************************************" << endl;
628                    
629                    TIter next(temp);
630                    TSystemFile *questo = 0;
631                    
632                    
633                    while ( (questo = (TSystemFile*) next()) ) {
634                            TString name =  questo-> GetName();
635                            if( name.EndsWith(".root") ){
636                                    char *fullpath;
637                                    gSystem->IsFileInIncludePath(name,&fullpath);
638                                    contents->Add(new TSystemFile(fullpath,gSystem->DirName(fullpath)));
639                            };
640                    }
641                    delete temp;
642            
643            };
644            gSystem->ChangeDirectory(wdir);
645            cout << endl << "Selected files:" << endl;
646            contents->Print();
647            cout << endl;
648    //      cout << "Warking directory: "<< gSystem->WorkingDirectory()<< endl;
649            return contents;
650    };
651    //--------------------------------------
652    //
653    //
654    //--------------------------------------
655    /**
656     * Get the Pamela detector chains from a list of files and make them friends.
657     * @param fl Pointer to a TList of TSystemFiles
658     * @param detlist String to select trees to be included
659     * @return Pointer to a TChain
660     */
661    TChain *PamLevel2::LoadPamTrees(TList *fl, TString detlist="+ALL"){
662    
663            TChain *Tout=0;
664            
665            SetWhichTrees(detlist);
666            
667            TChain *T = 0;  
668            TChain *C = 0;
669            TChain *O = 0;
670            TChain *R = 0;
671            TChain *S = 0;
672            TChain *N = 0;
673            TChain *A = 0;
674            TChain *B = 0;
675            
676            if(TRK) T = new TChain("Tracker");      
677            if(CAL) C = new TChain("Calorimeter");
678            if(TOF) O = new TChain("ToF");
679            if(TRG) R = new TChain("Trigger");
680            if(S4)  S = new TChain("S4");
681            if(ND)  N = new TChain("NeutronD");
682            if(AC)  A = new TChain("Anticounter");
683            if(ORB) B = new TChain("OrbitalInfo");
684            
685            // loop over files and create chains    
686            TIter next(fl);
687            TSystemFile *questo = 0;
688            while ( (questo = (TSystemFile*) next()) ) {
689                    TString name =  questo->GetName();
690    //              cout << "File: "<< name << endl;
691                    if( CheckLevel2File(name) ){
692                            if(TRK) T->Add(name);
693                            if(CAL) C->Add(name);
694                            if(TOF) O->Add(name);
695                            if(TRG) R->Add(name);
696                            if(S4)  S->Add(name);
697                            if(ND)  N->Add(name);
698                            if(AC)  A->Add(name);
699                            if(ORB) B->Add(name);
700                    };
701            }
702            
703        // Tracker
704            if(TRK) {
705                    T->SetBranchAddress("TrkLevel2", GetPointerToTrk());
706                    cout << "Tracker      : set branch address TrkLevel2"<<endl;
707                    if(!Tout)Tout=T;
708            };
709        
710            // Calorimeter
711            if(CAL) {
712                    C->SetBranchAddress("CaloLevel2", GetPointerToCalo());
713                    cout << "Calorimeter  : set branch address CaloLevel2"<<endl;
714                    if(!Tout)Tout=C;
715                    else Tout->AddFriend("Calorimeter");
716            };
717            
718        // ToF    
719            if(TOF) {
720                    O->SetBranchAddress("ToFLevel2", GetPointerToToF());
721                    cout << "ToF          : set branch address ToFLevel2"<<endl;
722                    if(!Tout)Tout=O;
723                    else Tout->AddFriend("ToF");
724            };
725        // Trigger
726            if(TRG) {
727                    R->SetBranchAddress("TrigLevel2", GetPointerToTrig());
728                    cout << "Trigger      : set branch address TrigLevel2"<<endl;
729                    if(!Tout)Tout=O;
730                    else Tout->AddFriend("Trigger");
731            };
732        // S4
733            if(S4) {
734                    S->SetBranchAddress("S4Level2", GetPointerToS4());
735                    cout << "S4           : set branch address S4Level2"<<endl;
736                    if(!Tout)Tout=O;
737                    else Tout->AddFriend("S4");
738            };
739        // Neutron Detector
740            if(ND) {
741                    N->SetBranchAddress("NDLevel2", GetPointerToND());
742                    cout << "NeutronD     : set branch address NDLevel2"<<endl;
743                    if(!Tout)Tout=O;
744                    else Tout->AddFriend("NeutronD");
745            };
746        // Anticounters
747            if(AC) {
748                    A->SetBranchAddress("AcLevel2", GetPointerToAc());
749                    cout << "Anticounter  : set branch address AcLevel2"<<endl;
750                    if(!Tout)Tout=O;
751                    else Tout->AddFriend("Anticounter");
752            };
753        // OrbitalInfo
754            if(ORB) {
755                    B->SetBranchAddress("OrbitalInfo", GetPointerToOrb());
756                    cout << "OrbitalInfo  : set branch address OrbitalInfo"<<endl;
757                    if(!Tout)Tout=O;
758                    else Tout->AddFriend("OrbitalInfo");
759            };
760        
761            return Tout;
762        
763    }
764    //--------------------------------------
765    //
766    //
767    //--------------------------------------
768    /**
769     * Set which trees should be analysed
770     * @param detlist TString containing the sequence of trees required
771    */
772    void PamLevel2::SetWhichTrees(TString detlist){
773            
774            if(detlist.Contains("+ALL", TString::kIgnoreCase)){
775                    CAL = true;
776                    TRK = true;
777                    TRG = true;
778                    TOF = true;
779                    S4  = true;
780                    ND  = true;
781                    AC  = true;
782                    ORB = true;
783            }else if( detlist.Contains("-ALL", TString::kIgnoreCase) ){
784                    CAL = false;
785                    TRK = false;
786                    TRG = false;
787                    TOF = false;
788                    S4  = false;
789                    ND  = false;
790                    AC  = false;
791                    ORB = false;
792            };
793            
794            if( detlist.Contains("-CAL", TString::kIgnoreCase) )CAL = false;
795            else if( detlist.Contains("+CAL", TString::kIgnoreCase) )CAL = true;
796            
797            if( detlist.Contains("-TRK", TString::kIgnoreCase) )TRK = false;
798            else if( detlist.Contains("+TRK", TString::kIgnoreCase) )TRK = true;
799            
800            if( detlist.Contains("-TRG", TString::kIgnoreCase) )TRG = false;
801            else if( detlist.Contains("+TRG", TString::kIgnoreCase) )TRG = true;
802            
803            if( detlist.Contains("-TOF", TString::kIgnoreCase) )TOF = false;
804            else if( detlist.Contains("+TOF", TString::kIgnoreCase) )TOF = true;
805            
806            if( detlist.Contains("-S4",  TString::kIgnoreCase) )S4  = false;
807            else if( detlist.Contains("+S4",  TString::kIgnoreCase) )S4  = true;
808            
809            if( detlist.Contains("-ND",  TString::kIgnoreCase) )ND  = false;
810            else if( detlist.Contains("+ND",  TString::kIgnoreCase) )ND  = true;
811            
812            if( detlist.Contains("-AC",  TString::kIgnoreCase) )AC  = false;
813            else if( detlist.Contains("+AC",  TString::kIgnoreCase) )AC  = true;
814            
815            if( detlist.Contains("-ORB", TString::kIgnoreCase) )ORB = false;
816            else if( detlist.Contains("+ORB", TString::kIgnoreCase) )ORB = true;
817            
818    };
819    //--------------------------------------
820    //
821    //
822    //--------------------------------------
823    /**
824     * Check if a file contains selected Pamela Level2 trees.
825     * @param name File name
826     * @return Pointer to a TChain
827     */
828    Bool_t  PamLevel2::CheckLevel2File(TString name){
829            
830            Bool_t CAL__ok = false;
831            Bool_t TRK__ok = false;
832            Bool_t TRG__ok = false;
833            Bool_t TOF__ok = false;
834            Bool_t S4__ok  = false;
835            Bool_t ND__ok  = false;
836            Bool_t AC__ok  = false;
837            Bool_t ORB__ok = false;
838            
839            Bool_t RUN__ok = false;
840            
841            TFile *f = new TFile(name.Data());
842            TList *lk = f->GetListOfKeys();
843    //      lk->Print();
844            TIter next(lk);
845            TKey *key =0;
846            while( (key = (TKey*)next()) ){
847    //              cout << key->GetName() << endl;
848                    if( !strcmp(key->GetName(),"Calorimeter") )CAL__ok = true;
849                    if( !strcmp(key->GetName(),"Tracker"    ) )TRK__ok = true;
850                    if( !strcmp(key->GetName(),"Trigger"    ) )TRG__ok = true;
851                    if( !strcmp(key->GetName(),"ToF"        ) )TOF__ok = true;
852                    if( !strcmp(key->GetName(),"S4"         ) )S4__ok = true;
853                    if( !strcmp(key->GetName(),"NeutronD"   ) )ND__ok = true;
854                    if( !strcmp(key->GetName(),"Anticounter") )AC__ok = true;
855                    if( !strcmp(key->GetName(),"OrbitalInfo") )ORB__ok = true;
856                    if( !strcmp(key->GetName(),"Run"        ) )RUN__ok = true;
857            };
858                    
859            lk->Delete();
860            f->Close();
861            
862            Bool_t FLAG = true;
863            if(!RUN__ok) {
864                    cout << "File: "<< f->GetName() <<" discarded ---- Missing RunInfo tree"<< endl;
865                    FLAG = false;
866            };
867            if(CAL && !CAL__ok){
868                    cout << "File: "<< f->GetName() <<" discarded ---- Missing Calorimeter tree"<< endl;
869                    FLAG = false;
870            };
871            if(TRK && !TRK__ok){
872                    cout << "File: "<< f->GetName() <<" discarded ---- Missing Tracker tree"<< endl;
873                    FLAG = false;
874            };
875            if(TRG && !TRG__ok){
876                    cout << "File: "<< f->GetName() <<" discarded ---- Missing Trigger tree"<< endl;
877                    FLAG = false;
878            };
879            if(TOF && !TOF__ok){
880                    cout << "File: "<< f->GetName() <<" discarded ---- Missing ToF tree"<< endl;
881                    FLAG = false;
882            };
883            if(S4 && !S4__ok){
884                    cout << "File: "<< f->GetName() <<" discarded ---- Missing S4 tree"<< endl;
885                    FLAG = false;
886            };
887            if(ND && !ND__ok){
888                    cout << "File: "<< f->GetName() <<" discarded ---- Missing ND tree"<< endl;
889                    FLAG = false;
890            };
891            if(AC && !AC__ok){
892                    cout << "File: "<< f->GetName() <<" discarded ---- Missing AC tree"<< endl;
893                    FLAG = false;
894            };
895            if(ORB && !ORB__ok){
896                    cout << "File: "<< f->GetName() <<" discarded ---- Missing ORB tree"<< endl;
897                    FLAG = false;
898            };
899            
900            return FLAG;
901            
902    };
903    
904    
Legend:
Removed from v.1.1  
changed lines
  Added in v.1.3
  ViewVC Help
Powered by ViewVC 1.1.23