/[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.27 by pam-fi, Fri Mar 16 20:34:08 2007 UTC revision 1.41 by pam-fi, Thu Aug 2 13:59:47 2007 UTC
# Line 85  void PamTrack::Delete(){ Line 85  void PamTrack::Delete(){
85  //  //
86  //--------------------------------------  //--------------------------------------
87  /**  /**
88   * Constructor   * Default Constructor
89   */   */
90  PamLevel2::PamLevel2(){  PamLevel2::PamLevel2(){
91      Initialize();
92    };
93    
94    /**
95     * Constructor with given dir, list and detectors
96     */
97    PamLevel2::PamLevel2(TString ddir,TString list,TString detlist){
98      Initialize();
99      GetPamTree(GetListOfLevel2Files(ddir,list),detlist);
100      GetRunTree(GetListOfLevel2Files(ddir,list));
101    };
102    
103    /**
104     * Constructor with given dir and list
105     */
106    PamLevel2::PamLevel2(TString ddir,TString list){
107      Initialize();
108      GetPamTree(GetListOfLevel2Files(ddir,list),"");
109      GetRunTree(GetListOfLevel2Files(ddir,list));
110    };
111    
112    
113    void PamLevel2::Initialize(){
114                    
115  //     trk2_obj  = TrkLevel2::GetTrkLevel2();  //     trk2_obj  = TrkLevel2::GetTrkLevel2();
116  //     trk1_obj  = TrkLevel1::GetTrkLevel1();  //     trk1_obj  = TrkLevel1::GetTrkLevel1();
# Line 103  PamLevel2::PamLevel2(){ Line 126  PamLevel2::PamLevel2(){
126            
127      h0_obj    = 0;      h0_obj    = 0;
128      trk0_obj  = 0;      trk0_obj  = 0;
129        calo0_obj  = 0;
130    
131      trk2_obj  = 0;      trk2_obj  = 0;
132      trk1_obj  = 0;      trk1_obj  = 0;
# Line 118  PamLevel2::PamLevel2(){ Line 142  PamLevel2::PamLevel2(){
142    
143      run_obj   = 0;//new GL_RUN();      run_obj   = 0;//new GL_RUN();
144      soft_obj   = 0;// Emiliano      soft_obj   = 0;// Emiliano
145      irun = -1;      irun = -1LL;
146      runfirstentry = 0ULL;      runfirstentry = 0LL;
147      runlastentry = 0ULL;      runlastentry = 0LL;
148    
149      l0_file = NULL;      l0_file = NULL;
150      l0_tree = NULL;      l0_tree = NULL;
151      iroot   = -1;      iroot   = -1;
152      dbc     = 0;      dbc     = 0;
153    
154        prevshift = 0;
155            
     irun = -1;  
156      run_tree = NULL;      run_tree = NULL;
157      run_tree_clone = NULL;      run_tree_clone = NULL;
158      sel_tree = NULL;      sel_tree = NULL;
159      sel_tree_clone = NULL;      sel_tree_clone = NULL;
160            
161      irunentry = -1;      irunentry = -1LL;
162      pam_tree = NULL;      pam_tree = NULL;
163      for(Int_t i=0; i<NCLONES; i++ )pam_tree_clone[i]=NULL;      for(Int_t i=0; i<NCLONES; i++ )pam_tree_clone[i]=NULL;
164    
# Line 175  PamLevel2::PamLevel2(){ Line 200  PamLevel2::PamLevel2(){
200    
201      tsorted=0;      tsorted=0;
202      timage=0;      timage=0;
203        
204        howtosort = "+CAL+TOF";
205        //howtosort = "+TOF";
206        sortthr = 100.;
207    
208  };  };
209  /**  /**
# Line 188  void PamLevel2::Delete(){ Line 217  void PamLevel2::Delete(){
217  //    cout << "void PamLevel2::Clear()"<<endl;  //    cout << "void PamLevel2::Clear()"<<endl;
218      if(h0_obj)   delete h0_obj;      if(h0_obj)   delete h0_obj;
219      if(trk0_obj) delete trk0_obj;            if(trk0_obj) delete trk0_obj;      
220        if(calo0_obj) delete calo0_obj;    
221      if(trk1_obj) delete trk1_obj;            if(trk1_obj) delete trk1_obj;      
222      if(trk2_obj) delete trk2_obj;      if(trk2_obj) delete trk2_obj;
223      if(trkh_obj) delete trkh_obj;      if(trkh_obj) delete trkh_obj;
# Line 215  void PamLevel2::Delete(){ Line 245  void PamLevel2::Delete(){
245      }      }
246            
247      if(l0_file)l0_file->Close();      if(l0_file)l0_file->Close();
248      if(pam_tree)pam_tree->Delete();;      //    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();;      if(run_tree)run_tree->Delete();;
282      if(sel_tree)sel_tree->Delete();;      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();;      for(Int_t i=0; i<NCLONES; i++ )if(pam_tree_clone[i])pam_tree_clone[i]->Delete();;
# Line 241  void PamLevel2::Clear(){ Line 303  void PamLevel2::Clear(){
303  //    if(soft_obj) soft_obj->Clear();  //    if(soft_obj) soft_obj->Clear();
304    
305      if(h0_obj)   h0_obj->Clear();            if(h0_obj)   h0_obj->Clear();      
306      if(trk0_obj) trk0_obj->Clear();      //    if(trk0_obj) trk0_obj->Clear();  
307      if(trk1_obj) trk1_obj->Clear();          if(trk1_obj) trk1_obj->Clear();    
308      if(trk2_obj) trk2_obj->Clear();      if(trk2_obj) trk2_obj->Clear();
309      if(trkh_obj) trkh_obj->Clear();      if(trkh_obj) trkh_obj->Clear();
310        if(calo0_obj) calo0_obj->Clear();  
311      if(calo1_obj)calo1_obj->Clear();      if(calo1_obj)calo1_obj->Clear();
312      if(calo2_obj)calo2_obj->Clear();      if(calo2_obj)calo2_obj->Clear();
313      if(tof_obj)  tof_obj->Clear();      if(tof_obj)  tof_obj->Clear();
# Line 263  void PamLevel2::Clear(){ Line 326  void PamLevel2::Clear(){
326      if(timage){      if(timage){
327          timage->Delete();          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  //  //
# Line 460  PamTrack* PamLevel2::GetStoredTrack(Int_ Line 609  PamTrack* PamLevel2::GetStoredTrack(Int_
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  /**  /**
# Line 473  PamTrack* PamLevel2::GetStoredTrack(Int_ Line 632  PamTrack* PamLevel2::GetStoredTrack(Int_
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).   * The total number of physical tracks is always given by GetNTracks() and the it-th physical track can be retrieved by means of the methods GetTrack(int it) and GetTrack(int it, TString how).
634   */   */
635  void PamLevel2::SortTracks(TString how){  void PamLevel2::SortTracks(){
636      TString how = howtosort;
637    
638  //    cout <<" PamLevel2::SortTracks(TString how) "<<endl;    //    cout <<" PamLevel2::SortTracks(TString how) "<<endl;
639      if( !trk2_obj ){    if( !trk2_obj ){
640          cout << "void PamLevel2::SortTracks(TString how):  TrkLevel2 not loaded !!!";      cout << "void PamLevel2::SortTracks():  TrkLevel2 not loaded !!!";
641          return;      return;
642      };    };
643  //    cout << "call SortTracs() "<<endl;    //    cout << "call SortTracs() "<<endl;
644    //Save current Object count    //Save current Object count
645      Int_t ObjectNumber = TProcessID::GetObjectCount();    Int_t ObjectNumber = TProcessID::GetObjectCount();
646    
647      //    cout << "ObjectNumber  "<<ObjectNumber <<endl;
648    
649  //    cout << "ObjectNumber  "<<ObjectNumber <<endl;    //     if(!sorted_tracks)sorted_tracks = new TRefArray();
650      //     sorted_tracks->Clear();
651      //    sorted_tracks.Clear();
652    
653  //     if(!sorted_tracks)sorted_tracks = new TRefArray();    if(!tsorted)tsorted = new TClonesArray("PamTrack",trk2_obj->GetNTracks());
654  //     sorted_tracks->Clear();    tsorted->Delete();
655  //    sorted_tracks.Clear();    TClonesArray &ttsorted = *tsorted;
656      if(!timage)timage = new TClonesArray("PamTrack",trk2_obj->GetNTracks());
657      timage->Delete();
658      TClonesArray &ttimage = *timage;
659    
     if(!tsorted)tsorted = new TClonesArray("PamTrack",trk2_obj->GetNTracks());  
     tsorted->Delete();  
     TClonesArray &ttsorted = *tsorted;  
     if(!timage)timage = new TClonesArray("PamTrack",trk2_obj->GetNTracks());  
     timage->Delete();  
     TClonesArray &ttimage = *timage;  
660    
661    
662      // loop over the tracks sorted by the tracker    //--------------------------------------------------
663      Bool_t use_TRK = how.Contains("TRK", TString::kIgnoreCase);    // retrieve sorting method
664      Bool_t use_CAL = how.Contains("CAL", TString::kIgnoreCase);    //--------------------------------------------------
665      Bool_t use_TOF = how.Contains("TOF", TString::kIgnoreCase);    Bool_t use_TRK = how.Contains("TRK", TString::kIgnoreCase);
666      Bool_t use_CAL = how.Contains("CAL", TString::kIgnoreCase);
667      Bool_t use_TOF = how.Contains("TOF", TString::kIgnoreCase);
668            
669      if( !CAL2 &&  use_CAL) use_CAL = false;    if( !CAL2 &&  use_CAL) use_CAL = false;
670      if( !TOF &&  use_TOF) use_TOF = false;    if( !TOF &&  use_TOF) use_TOF = false;
671            
672      if( !TRK2 ){    if( !TRK2 ){
673  //      cout << "SortTracks() : without tracker does not work!!! (not yet)" << endl;      //  cout << "SortTracks() : without tracker does not work!!! (not yet)" << endl;
674          return;      return;
675      };    };
676    
677      //   cout << "use_CAL "<<use_CAL<<" use_TOF "<<use_TOF<<" use_TRK "<<use_TRK <<endl;    //   cout << "use_CAL "<<use_CAL<<" use_TOF "<<use_TOF<<" use_TRK "<<use_TRK <<endl;
678          
679      for(Int_t i=0; i < trk2_obj->TrkLevel2::GetNTracks(); i++){  
680              //--------------------------------------------------
681          TrkTrack   *ts = 0;    // loop over "physical" tracks sorted by the tracker
682          CaloTrkVar *cs = 0;    //--------------------------------------------------
683          ToFTrkVar  *os = 0;    for(Int_t i=0; i < trk2_obj->TrkLevel2::GetNTracks(); i++){
684                    
685          // get tracker tracks      TrkTrack   *ts = 0;
686          TrkTrack   *tp = trk2_obj->TrkLevel2::GetTrack(i);                    //tracker      CaloTrkVar *cs = 0;
687          CaloTrkVar *cp = GetCaloStoredTrack(tp->GetSeqNo());      ToFTrkVar  *os = 0;
688          ToFTrkVar  *op = GetToFStoredTrack(tp->GetSeqNo());          
689        // get tracker tracks
690          TrkTrack   *ti = 0;              //tracker (image)      TrkTrack   *tp = trk2_obj->TrkLevel2::GetTrack(i); //tracker
691          CaloTrkVar *ci = 0;      CaloTrkVar *cp = GetCaloStoredTrack(tp->GetSeqNo());
692          ToFTrkVar  *oi = 0;      ToFTrkVar  *op = GetToFStoredTrack(tp->GetSeqNo());
693  //      cout << "trk track n. "<<i << " "<<hex<< tp <<dec<< endl;  
694          // if track has an image, check image selection      TrkTrack   *ti = 0; //tracker (image)
695          if(tp->HasImage()){      CaloTrkVar *ci = 0;
696        ToFTrkVar  *oi = 0;
697        //  cout << "trk track n. "<<i << " "<<hex<< tp <<dec<< endl;
698        // if track has an image, check image selection
699        if(tp->HasImage()){
700                            
701              ti = trk2_obj->TrkLevel2::GetTrackImage(i);              //tracker (image)        ti = trk2_obj->TrkLevel2::GetTrackImage(i);              //tracker (image)
702              ci = GetCaloStoredTrack(ti->GetSeqNo());        ci = GetCaloStoredTrack(ti->GetSeqNo());
703              oi = GetToFStoredTrack(ti->GetSeqNo());        oi = GetToFStoredTrack(ti->GetSeqNo());
704                            
705  //          cout << "its image "<<i << " "<<hex<< ti <<dec<< endl;        //            cout << "its image "<<i << " "<<hex<< ti <<dec<< endl;
706    
707              //assign starting scores        //assign starting scores
708              Int_t tp_score = 0;  //main track sorted by the tracker        Int_t tp_score = 0;  //main track sorted by the tracker
709              Int_t ti_score = 0;  //image track        Int_t ti_score = 0;  //image track
710                            
711              // ------------------------        // -----------------------------------------------------------------------------------------
712              // calorimeter check        // calorimeter check
713              // ------------------------        // -----------------------------------------------------------------------------------------
714              // check the Y spatial residual on the first calorimeter plane        // check the Y spatial residual on the first calorimeter plane
715              // (cut on calorimeter variables from Emiliano)        // (cut on calorimeter variables from Emiliano)
716              if( use_CAL && !calo2_obj ){        if( use_CAL && !calo2_obj ){
717                  cout << "void PamLevel2::SortTracks(TString how): how= "<<how<<" but CaloLevel2 not loaded !!!";          cout << "void PamLevel2::SortTracks(): howtosort= "<<how<<" but CaloLevel2 not loaded !!!";
718                  return;          return;
719              };        };
720              if(        if(
721                  use_CAL            &&           use_CAL            &&
722                  calo2_obj->npcfit[1] > 15     &&   //no. of fit planes on Y view  //       calo2_obj->npcfit[1] > 5     &&   //no. of fit planes on Y view
723                  calo2_obj->varcfit[1] < 1000. &&  //fit variance on Y view  //       calo2_obj->varcfit[1] < 1000. &&  //fit variance on Y view
724                  cp && ci &&           cp && ci &&
725                  true){           true){
726    
727                                    
728                  Float_t resy_p = cp->tbar[0][1] - calo2_obj->cbar[0][1]; if(resy_p < 0)resy_p= - resy_p;  //      Float_t resy_p = cp->tbar[0][1] - calo2_obj->cbar[0][1]; if(resy_p < 0)resy_p= - resy_p;
729                  Float_t resy_i = ci->tbar[0][1] - calo2_obj->cbar[0][1]; if(resy_i < 0)resy_i= - resy_i;  //      Float_t resy_i = ci->tbar[0][1] - calo2_obj->cbar[0][1]; if(resy_i < 0)resy_i= - resy_i;
730                                    
731                  if(resy_p <= resy_i) tp_score++;  //      if(resy_p <= resy_i) tp_score++;
732                  else                 ti_score++;  //      else                 ti_score++;
733            
734    
735              if( cp->npresh >= ci->npresh) tp_score++;
736              else ti_score++;
737    
738  //              cout << "CALO "<<tp_score<<ti_score<<endl;          //      cout << "CALO "<<tp_score<<ti_score<<endl;
739    
740          };
741          // -----------------------------------------------------------------------------------------
742          // TOF check
743          // -----------------------------------------------------------------------------------------
744          // check the number of hit pmts along the track
745          // on S12 S21 and S32, where paddles are parallel to Y axis
746          if( use_TOF && !tof_obj ){
747            cout << "void PamLevel2::SortTracks(): howtosort= "<<how<<" but ToFLevel2 not loaded !!!";
748            return;
749          };
750          if( use_TOF && op && oi ){
751                    
752            //
753            Float_t sen = 0.;
754            for (Int_t ih=0; ih < op->npmtadc; ih++){
755              Int_t pl = tof_obj->GetPlaneIndex( (op->pmtadc).At(ih) );
756              if ( pl == 2 || pl == 3 || pl == 4 || pl == 5 ) sen += (op->dedx).At(ih);
757            };
758            for (Int_t ih=0; ih < oi->npmtadc; ih++){
759              Int_t pl = tof_obj->GetPlaneIndex( (oi->pmtadc).At(ih) );
760              if ( pl == 2 || pl == 3 || pl == 4 || pl == 5 ) sen += (oi->dedx).At(ih);
761            };
762            //
763            if (  sen >= sortthr ){
764              //printf(" IS A NUCLEUS! en = %f \n",sen);
765              //
766              // is a nucleus use a different algorithm
767              //
768              Int_t nz = 6; Float_t zin[6];                                          // << define TOF z-coordinates
769              for(Int_t ip=0; ip<nz; ip++)
770                zin[ip] = tof_obj->ToFLevel2::GetZTOF(tof_obj->ToFLevel2::GetToFPlaneID(ip));     // << read ToF plane z-coordinates
771              Trajectory *tr = new Trajectory(nz,zin);
772              //
773              Int_t nphit_p =0;
774              Int_t nphit_i =0;
775              Float_t enhit_p = 0.;
776              Float_t enhit_i = 0.;
777              //
778              for (Int_t ih=0; ih < op->npmtadc; ih++){
779                Int_t pl = tof_obj->GetPlaneIndex( (op->pmtadc).At(ih) );
780                if(pl == 1 || pl == 2 || pl == 5){
781                  nphit_p++;
782                  enhit_p += (op->dedx).At(ih);
783                };
784              };
785              //
786              tp->DoTrack2(tr);
787              //
788              if ( fabs(tr->y[0]-oi->ytofpos[0]) < 2. ){
789                for (Int_t ih=0; ih < op->npmtadc; ih++){
790                  Int_t pl = tof_obj->GetPlaneIndex( (op->pmtadc).At(ih) );
791                  if(pl == 0){
792                    nphit_p++;
793                    enhit_p += (op->dedx).At(ih);
794                  };
795                };
796              };
797              if ( fabs(tr->y[3]-oi->ytofpos[1]) < 2. ){
798                for (Int_t ih=0; ih < op->npmtadc; ih++){
799                  Int_t pl = tof_obj->GetPlaneIndex( (op->pmtadc).At(ih) );
800                  if(pl == 3){
801                    nphit_p++;
802                    enhit_p += (op->dedx).At(ih);
803                  };
804              };              };
805              // ------------------------            };
806              // TOF check            if ( fabs(tr->y[4]-oi->ytofpos[2]) < 2. ){
807              // ------------------------                  for (Int_t ih=0; ih < op->npmtadc; ih++){
808              // check the number of hit pmts along the track                Int_t pl = tof_obj->GetPlaneIndex( (op->pmtadc).At(ih) );
809              // on S12 S21 and S32, where paddles are parallel to Y axis                if(pl == 4){
810              if( use_TOF && !tof_obj ){                  nphit_p++;
811                  cout << "void PamLevel2::SortTracks(TString how): how= "<<how<<" but ToFLevel2 not loaded !!!";                  enhit_p += (op->dedx).At(ih);
812                  return;                };
813              };              };
814              if( use_TOF && op && oi ){            };
815    
816              for (Int_t ih=0; ih < oi->npmtadc; ih++){
817                Int_t pl = tof_obj->GetPlaneIndex( (oi->pmtadc).At(ih) );
818                if(pl == 1 || pl == 2 || pl == 5){
819                  nphit_i++;        
820                  enhit_i += (op->dedx).At(ih);
821                };
822              };
823              //
824              ti->DoTrack2(tr);
825              //
826              if ( fabs(tr->y[0]-oi->ytofpos[0]) < 2. ){
827                for (Int_t ih=0; ih < oi->npmtadc; ih++){
828                  Int_t pl = tof_obj->GetPlaneIndex( (oi->pmtadc).At(ih) );
829                  if(pl == 0){
830                    nphit_i++;
831                    enhit_i += (op->dedx).At(ih);
832                  };
833                };
834              };
835              if ( fabs(tr->y[3]-oi->ytofpos[1]) < 2. ){
836                for (Int_t ih=0; ih < oi->npmtadc; ih++){
837                  Int_t pl = tof_obj->GetPlaneIndex( (oi->pmtadc).At(ih) );
838                  if(pl == 3){
839                    nphit_i++;
840                    enhit_i += (op->dedx).At(ih);
841                  };
842                };
843              };
844              if ( fabs(tr->y[4]-oi->ytofpos[2]) < 2. ){
845                for (Int_t ih=0; ih < oi->npmtadc; ih++){
846                  Int_t pl = tof_obj->GetPlaneIndex( (oi->pmtadc).At(ih) );
847                  if(pl == 4){
848                    nphit_i++;
849                    enhit_i += (op->dedx).At(ih);
850                  };
851                };
852              };            
853    
854                                    
855                  Int_t nphit_p =0;            if(
856                  Int_t nphit_i =0;               use_TOF            &&
857                 (nphit_p+nphit_i) !=0 &&  
858                 true){
859                        
860                //      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);
861                //      printf(" score p %i score i %i \n",tp_score,ti_score);
862                //            if( enhit_p > enhit_i ) tp_score++;
863                //            if( nphit_p >= nphit_i && enhit_p > enhit_i ) tp_score++;
864                if ( nphit_p > nphit_i ) tp_score++;
865                if ( nphit_p < nphit_i ) ti_score++;
866                if ( nphit_p == nphit_i ){
867                  if ( enhit_p > enhit_i ) tp_score++;
868                  else ti_score++;
869                };
870                //      printf(" dopo score p %i score i %i \n",tp_score,ti_score);
871              };
872              delete tr;
873              //
874            } else {
875              //
876              // NOT a NUCLEUS
877              //
878              //printf(" NOT a NUCLEUS! en = %f \n",sen);
879    
880              Int_t nphit_p =0;
881              Int_t nphit_i =0;
882                                    
883                                    
884  /*                              cout << "track: npmtadc "<< op->npmtadc << endl;            /*                            cout << "track: npmtadc "<< op->npmtadc << endl;
885                                  cout << "track: npmttdc "<< op->npmttdc << endl;              cout << "track: npmttdc "<< op->npmttdc << endl;
886                                  cout << "image: npmtadc "<< oi->npmtadc << endl;              cout << "image: npmtadc "<< oi->npmtadc << endl;
887                                  cout << "image: npmttdc "<< oi->npmttdc << endl;*/              cout << "image: npmttdc "<< oi->npmttdc << endl;*/
888                                    
889                  for (Int_t ih=0; ih < op->npmtadc; ih++){  //        for (Int_t ih=0; ih < op->npmtadc; ih++){
890                      Int_t pl = tof_obj->GetPlaneIndex( (op->pmtadc).At(ih) );  //          Int_t pl = tof_obj->GetPlaneIndex( (op->pmtadc).At(ih) );
891                      if(pl == 1 || pl == 2 || pl == 5)nphit_p++;  //          if(pl == 1 || pl == 2 || pl == 5)nphit_p++;
892                  };  //        };
893                                    
894                  for (Int_t ih=0; ih < oi->npmtadc; ih++){  //        for (Int_t ih=0; ih < oi->npmtadc; ih++){
895                      Int_t pl = tof_obj->GetPlaneIndex( (oi->pmtadc).At(ih) );  //          Int_t pl = tof_obj->GetPlaneIndex( (oi->pmtadc).At(ih) );
896                      if(pl == 1 || pl == 2 || pl == 5)nphit_i++;  //          if(pl == 1 || pl == 2 || pl == 5)nphit_i++;
897                  };  //        };
898              // --- modified to count tdc signals (more efficient?)
899              // --- and to implement check on tdcflag
900              for (Int_t ih=0; ih < op->npmttdc; ih++){
901                Int_t pl = tof_obj->GetPlaneIndex( (op->pmttdc).At(ih) );
902    //          if( (op->tdcflag).At(ih)==0 && (pl == 1 || pl == 2 || pl == 5) )nphit_p++;
903                if( (op->tdcflag).At(ih)==0 )nphit_p++;
904              };
905                                    
906                  if(            for (Int_t ih=0; ih < oi->npmttdc; ih++){
907                      use_TOF            &&              Int_t pl = tof_obj->GetPlaneIndex( (oi->pmttdc).At(ih) );
908                      (nphit_p+nphit_i) !=0 &&      //          if( (oi->tdcflag).At(ih)==0 && (pl == 1 || pl == 2 || pl == 5) )nphit_i++;  
909                      true){              if( (oi->tdcflag).At(ih)==0 )nphit_i++;    
910              };
911                    
912              if(
913                 (nphit_p+nphit_i) !=0 &&  
914                 true){
915                                            
916                      if( nphit_p >= nphit_i) tp_score++;              if( nphit_p >= nphit_i) tp_score++;
917                      else ti_score++;              else ti_score++;
918                  };            };
919  //              cout << "TOF "<<tp_score<<ti_score<<endl;          };
920              };          //      cout << "TOF "<<tp_score<<ti_score<<endl;
921              if(tp_score == ti_score) use_TRK = true;        };
922              // ------------------------  
923              // tracker check  
924              // ------------------------  //      if(tp_score == ti_score) use_TRK = true;
925              // chi**2 difference is not always large enough to distinguish among  
926              // the real track and its image.  
927              // Tracker check will be applied always when calorimeter and tof information is ambiguous.        // -----------------------------------------------------------------------------------------
928              if(use_TRK){        // tracker check
929                  if(      tp->chi2 > 0 && tp->chi2 < ti->chi2 ) tp_score++ ;        // -----------------------------------------------------------------------------------------
930                  else if( ti->chi2 > 0 && ti->chi2 < tp->chi2 ) ti_score++ ;        // chi**2 difference is not always large enough to distinguish among
931  //              cout << "TRK "<<tp_score<<ti_score<<endl;        // the real track and its image.
932              };        // Tracker check will be applied always when calorimeter and tof information is ambiguous.
933          // *** MODIFIED ON AUGUST 2007 ***
934          if(use_TRK){
935    //      if(      tp->chi2 > 0 && tp->chi2 < ti->chi2 ) tp_score++ ;
936    //      else if( ti->chi2 > 0 && ti->chi2 < tp->chi2 ) ti_score++ ;
937    
938              // CHECK 1 : number of points along X  
939              if ( tp->GetNX() >= ti->GetNX() )tp_score++ ;
940              if ( tp->GetNX() <= ti->GetNX() )ti_score++ ;
941              // CHECK 2 : number of points along Y  
942              if ( tp->GetNY() >= ti->GetNY() )tp_score++ ;
943              if ( tp->GetNY() <= ti->GetNY() )ti_score++ ;
944              // CHECK 3 : chi**2 along X  
945    //        if(      tp->GetChi2X() > 0 && (tp->GetChi2X() < ti->GetChi2X() || ti->GetChi2X() <=0) ) tp_score++ ;
946    //        if(      ti->GetChi2X() > 0 && (ti->GetChi2X() < tp->GetChi2X() || tp->GetChi2X() <=0) ) ti_score++ ;
947              // CHECK 4 : chi**2 along Y  
948    //        if(      tp->GetChi2Y() > 0 && (tp->GetChi2Y() < ti->GetChi2Y() || ti->GetChi2Y() <=0) ) tp_score++ ;
949    //        if(      ti->GetChi2Y() > 0 && (ti->GetChi2Y() < tp->GetChi2Y() || tp->GetChi2Y() <=0) ) ti_score++ ;
950              // CHECK 5 : total chi**2  
951    //        if(      tp->chi2 > 0 && (tp->chi2 < ti->chi2 || ti->chi2 <=0) ) tp_score++ ;
952    //        if(      ti->chi2 > 0 && (ti->chi2 < tp->chi2 || tp->chi2 <=0) ) ti_score++ ;
953    
954            //      cout << "TRK "<<tp_score<<ti_score<<endl;
955          };
956                            
957              // ------------------------  
958              // the winner is....  
959              // ------------------------            // *-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
960              if      (tp_score > ti_score) {        // the winner is....
961  //              ts = tp;//the track sorted by the tracker!!        // *-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
962  //              cs = cp;        if      (tp_score > ti_score) {
963  //              os = op;          
964              }else if (tp_score < ti_score) {  
965                  ts = ti;//its image!!        }else if (tp_score < ti_score) {
966                  cs = ci;  
967                  os = oi;  
968            ts = ti;//its image!!
969                  ti = tp;//its image!!          cs = ci;
970                  ci = cp;          os = oi;
971                  oi = op;  
972            ti = tp;//its image!!
973                  tp = ts;//its image!!          ci = cp;
974                  cp = cs;          oi = op;
975                  op = os;  
976            tp = ts;//its image!!
977            cp = cs;
978            op = os;
979    
980                                    
981              }else {        }else {
982  //              ts = tp;  
983  //              cs = cp;  //      cout << "Warning - track image ambiguity not solved" << endl;
984  //              os = op;  
985  //                              cout << "Warning - track image ambiguity not solved" << endl;        };
 //                              cout << ts->GetNtot() << " "<< ts->chi2 << " " << npcfit[1] << " "<< nphit_p << endl;  
             };  
986                            
987          }else{      }else{
988  //          ts = tp;        //            ts = tp;
989  //          cs = cp;        //            cs = cp;
990  //          os = op;        //            os = op;
991          };      };
992                    
993  //      cout <<" SortTracks() "<<i<<" -- "<<ts<<endl;      //  cout <<" SortTracks() "<<i<<" -- "<<ts<<endl;
994  //      sorted_tracks->Add(ts);//save the track in the sorted array      //  sorted_tracks->Add(ts);//save the track in the sorted array
995  //      sorted_tracks.Add(ts);//save the track in the sorted array      //  sorted_tracks.Add(ts);//save the track in the sorted array
996  //      sorted_tracks.Add(tp);//save the track in the sorted array      //  sorted_tracks.Add(tp);//save the track in the sorted array
997  //      cout << "SortTracks:: sorted_tracks->Add(it) "<<i<<" "<<ts<<endl;      //  cout << "SortTracks:: sorted_tracks->Add(it) "<<i<<" "<<ts<<endl;
998  //      cout<<"o "<<tp<<endl;      //  cout<<"o "<<tp<<endl;
999  //      cout<<"o "<<cp<<endl;      //  cout<<"o "<<cp<<endl;
1000  //      cout<<"o "<<op<<endl;      //  cout<<"o "<<op<<endl;
1001          new(ttsorted[i]) PamTrack(tp,cp,op);      new(ttsorted[i]) PamTrack(tp,cp,op);
1002          new(ttimage[i])  PamTrack(ti,ci,oi);      new(ttimage[i])  PamTrack(ti,ci,oi);
1003      };    };
1004    
1005      if( tsorted->GetEntries() != trk2_obj->GetNTracks() ){    if( tsorted->GetEntries() != trk2_obj->GetNTracks() ){
1006          cout << "void PamLevel2::SortTracks(TString how): tsorted->GetEntries() "<<tsorted->GetEntries()<<" != trk2_obj->GetNTracks() = "<<trk2_obj->GetNTracks() <<endl;      cout << "void PamLevel2::SortTracks(): tsorted->GetEntries() "<<tsorted->GetEntries()<<" != trk2_obj->GetNTracks() = "<<trk2_obj->GetNTracks() <<endl;
1007          tsorted->Delete(); tsorted=0;      tsorted->Delete(); tsorted=0;
1008          timage->Delete(); timage=0;          timage->Delete(); timage=0;
1009      }    }
1010    
1011      //Restore Object count    //Restore Object count
1012      //To save space in the table keeping track of all referenced objects    //To save space in the table keeping track of all referenced objects
1013      //We reset the object count to what it was at the beginning of the event.    //We reset the object count to what it was at the beginning of the event.
1014      TProcessID::SetObjectCount(ObjectNumber);    TProcessID::SetObjectCount(ObjectNumber);
1015            
1016  };  };
1017  //--------------------------------------  //--------------------------------------
# Line 700  void PamLevel2::SortTracks(TString how){ Line 1035  void PamLevel2::SortTracks(TString how){
1035  TClonesArray *PamLevel2::GetTracks(){  TClonesArray *PamLevel2::GetTracks(){
1036    
1037  //  *-*-*-*-*-*-*-*-*-*-*-*-*  //  *-*-*-*-*-*-*-*-*-*-*-*-*
1038      SortTracks("+CAL+TOF");      SortTracks();
1039  //  *-*-*-*-*-*-*-*-*-*-*-*-*  //  *-*-*-*-*-*-*-*-*-*-*-*-*
1040    
1041      return tsorted;      return tsorted;
# Line 719  PamTrack *PamLevel2::GetTrack(int it){ Line 1054  PamTrack *PamLevel2::GetTrack(int it){
1054    
1055  //    if(!trk2_obj) return 0;  //    if(!trk2_obj) return 0;
1056                    
1057  // //  *-*-*-*-*-*-*-*-*-*-*-*-*  // //  *-*-*-*+-*-*-*-*-*-*-*-*-*
1058  //     SortTracks("+CAL+TOF");  //     SortTracks("+CAL+TOF");
1059  // //  *-*-*-*-*-*-*-*-*-*-*-*-*  // //  *-*-*-*-*-*-*-*-*-*-*-*-*
1060  // //    if(!sorted_tracks)return 0;  // //    if(!sorted_tracks)return 0;
# Line 738  PamTrack *PamLevel2::GetTrack(int it){ Line 1073  PamTrack *PamLevel2::GetTrack(int it){
1073    
1074  //    cout << "PamLevel2::GetTrack(int it) "<<endl;  //    cout << "PamLevel2::GetTrack(int it) "<<endl;
1075  //  *-*-*-*-*-*-*-*-*-*-*-*-*  //  *-*-*-*-*-*-*-*-*-*-*-*-*
1076      SortTracks("+CAL+TOF");      SortTracks();
1077  //  *-*-*-*-*-*-*-*-*-*-*-*-*  //  *-*-*-*-*-*-*-*-*-*-*-*-*
1078      if(!tsorted)return 0;      if(!tsorted)return 0;
1079      if(!tsorted->GetEntries())return 0;      if(!tsorted->GetEntries())return 0;
# Line 795  PamTrack *PamLevel2::GetTrackImage(int i Line 1130  PamTrack *PamLevel2::GetTrackImage(int i
1130    
1131    
1132  //  *-*-*-*-*-*-*-*-*-*-*-*-*  //  *-*-*-*-*-*-*-*-*-*-*-*-*
1133      SortTracks("+CAL+TOF");      SortTracks();
1134  //  *-*-*-*-*-*-*-*-*-*-*-*-*  //  *-*-*-*-*-*-*-*-*-*-*-*-*
1135      if(!timage)return 0;      if(!timage)return 0;
1136      if(!timage->GetEntries())return 0;      if(!timage->GetEntries())return 0;
# Line 832  TTree *PamLevel2::GetPamTree(TFile *f, T Line 1167  TTree *PamLevel2::GetPamTree(TFile *f, T
1167    
1168  //     if( !detlist.IsNull() ) SetWhichTrees(detlist);  //     if( !detlist.IsNull() ) SetWhichTrees(detlist);
1169  //     else                    GetWhichTrees(f);  //     else                    GetWhichTrees(f);
1170          if ( pam_tree ){
1171        printf("WARNING: TTree *PamLevel2::GetPamTree(TFile *fl, TString detlist) -- pam_tree already exists!\n ");
1172        return pam_tree;
1173      };
1174      //
1175    
1176      cout << "TTree *PamLevel2::GetPamTree(TFile *f, TString detlist ) -- obsolte "<<endl;      cout << "TTree *PamLevel2::GetPamTree(TFile *f, TString detlist ) -- obsolte "<<endl;
1177    
1178      SetWhichTrees(detlist);      SetWhichTrees(detlist);
# Line 908  TTree *PamLevel2::GetPamTree(TFile *f, T Line 1248  TTree *PamLevel2::GetPamTree(TFile *f, T
1248      if(R && TRG) {      if(R && TRG) {
1249          R->SetBranchAddress("TrigLevel2", GetPointerTo("TrigLevel2"));          R->SetBranchAddress("TrigLevel2", GetPointerTo("TrigLevel2"));
1250          cout << "Trigger      : set branch address TrigLevel2"<<endl;          cout << "Trigger      : set branch address TrigLevel2"<<endl;
1251          if(!Trout)Trout=O;          if(!Trout)Trout=R;
1252          else Trout->AddFriend(R);          else Trout->AddFriend(R);
1253      }else{      }else{
1254          cout << "Trigger      : missing tree"<<endl;          cout << "Trigger      : missing tree"<<endl;
# Line 918  TTree *PamLevel2::GetPamTree(TFile *f, T Line 1258  TTree *PamLevel2::GetPamTree(TFile *f, T
1258      if(S && S4) {      if(S && S4) {
1259          S->SetBranchAddress("S4Level2", GetPointerTo("S4Level2"));          S->SetBranchAddress("S4Level2", GetPointerTo("S4Level2"));
1260          cout << "S4           : set branch address S4Level2"<<endl;          cout << "S4           : set branch address S4Level2"<<endl;
1261          if(!Trout)Trout=O;          if(!Trout)Trout=S;
1262          else Trout->AddFriend(S);          else Trout->AddFriend(S);
1263      }else{      }else{
1264          cout << "S4           : missing tree"<<endl;          cout << "S4           : missing tree"<<endl;
# Line 928  TTree *PamLevel2::GetPamTree(TFile *f, T Line 1268  TTree *PamLevel2::GetPamTree(TFile *f, T
1268      if(N && ND) {      if(N && ND) {
1269          N->SetBranchAddress("NDLevel2", GetPointerTo("NDLevel2"));          N->SetBranchAddress("NDLevel2", GetPointerTo("NDLevel2"));
1270          cout << "NeutronD     : set branch address NDLevel2"<<endl;          cout << "NeutronD     : set branch address NDLevel2"<<endl;
1271          if(!Trout)Trout=O;          if(!Trout)Trout=N;
1272          else Trout->AddFriend(N);          else Trout->AddFriend(N);
1273      }else{      }else{
1274          cout << "NeutronD     : missing tree"<<endl;          cout << "NeutronD     : missing tree"<<endl;
# Line 938  TTree *PamLevel2::GetPamTree(TFile *f, T Line 1278  TTree *PamLevel2::GetPamTree(TFile *f, T
1278      if(A && AC) {      if(A && AC) {
1279          A->SetBranchAddress("AcLevel2", GetPointerTo("AcLevel2"));          A->SetBranchAddress("AcLevel2", GetPointerTo("AcLevel2"));
1280          cout << "Anticounter  : set branch address AcLevel2"<<endl;          cout << "Anticounter  : set branch address AcLevel2"<<endl;
1281          if(!Trout)Trout=O;          if(!Trout)Trout=A;
1282          else Trout->AddFriend(A);          else Trout->AddFriend(A);
1283      }else{      }else{
1284          cout << "Anticounter  : missing tree"<<endl;          cout << "Anticounter  : missing tree"<<endl;
# Line 948  TTree *PamLevel2::GetPamTree(TFile *f, T Line 1288  TTree *PamLevel2::GetPamTree(TFile *f, T
1288      if(B && ORB) {      if(B && ORB) {
1289          B->SetBranchAddress("OrbitalInfo", GetPointerTo("OrbitalInfo"));          B->SetBranchAddress("OrbitalInfo", GetPointerTo("OrbitalInfo"));
1290          cout << "OrbitalInfo  : set branch address OrbitalInfo"<<endl;          cout << "OrbitalInfo  : set branch address OrbitalInfo"<<endl;
1291          if(!Trout)Trout=O;          if(!Trout)Trout=B;
1292          else Trout->AddFriend(B);          else Trout->AddFriend(B);
1293      }else{      }else{
1294          cout << "OrbitalInfo  : missing tree"<<endl;          cout << "OrbitalInfo  : missing tree"<<endl;
# Line 969  TTree *PamLevel2::GetPamTree(TFile *f, T Line 1309  TTree *PamLevel2::GetPamTree(TFile *f, T
1309      };      };
1310            
1311      cout<<endl<<" Number of entries: "<<Trout->GetEntries()<<endl<<endl;      cout<<endl<<" Number of entries: "<<Trout->GetEntries()<<endl<<endl;
1312    
1313        pam_tree = (TChain*)Trout;
1314            
1315      return Trout;      return Trout;
1316            
# Line 988  TList*  PamLevel2::GetListOfLevel2Files( Line 1330  TList*  PamLevel2::GetListOfLevel2Files(
1330                    
1331      TString wdir = gSystem->WorkingDirectory();      TString wdir = gSystem->WorkingDirectory();
1332            
1333      if(ddir=="")ddir = wdir;      //    if(ddir=="")ddir = wdir;
1334  //      TSystemDirectory *datadir = new TSystemDirectory(gSystem->BaseName(ddir),ddir);  //      TSystemDirectory *datadir = new TSystemDirectory(gSystem->BaseName(ddir),ddir);
1335      cout << "Level2 data directory : "<<  ddir << endl;      if ( ddir != ""){
1336              cout << "Level2 data directory : "<<  ddir << endl;
1337        } else {
1338          cout << "Level2 data directory not given as input: trying to evaluate from list or taking working directory " << endl;
1339        };
1340      TList *contents  = new TList; // create output list      TList *contents  = new TList; // create output list
1341      contents->SetOwner();      contents->SetOwner();
1342            
# Line 1007  TList*  PamLevel2::GetListOfLevel2Files( Line 1352  TList*  PamLevel2::GetListOfLevel2Files(
1352  //          return 0;  //          return 0;
1353  //      }        //      }      
1354  //      flisttxt = fullpath;  //      flisttxt = fullpath;
1355          if ( !flisttxt.EndsWith(".root") ){
1356    
1357          flisttxt = gSystem->ConcatFileName(gSystem->DirName(flisttxt),gSystem->BaseName(flisttxt));          flisttxt = gSystem->ConcatFileName(gSystem->DirName(flisttxt),gSystem->BaseName(flisttxt));
1358    
# Line 1039  TList*  PamLevel2::GetListOfLevel2Files( Line 1385  TList*  PamLevel2::GetListOfLevel2Files(
1385  //          if( gSystem->IsFileInIncludePath(file,&fullpath) ){  //          if( gSystem->IsFileInIncludePath(file,&fullpath) ){
1386  //          if( (fullpath = gSystem->FindFile(ddir,file)) ){  //          if( (fullpath = gSystem->FindFile(ddir,file)) ){
1387              if( file.EndsWith(".root") ){              if( file.EndsWith(".root") ){
1388                  char *fullpath = gSystem->ConcatFileName(gSystem->DirName(ddir),gSystem->BaseName(file));                TString filedir;
1389                  if (ddir != ""){
1390                    filedir = ddir; // take the input dir
1391                  } else {
1392                    gSystem->ChangeDirectory(wdir); // back to the working directory
1393                    filedir = gSystem->DirName(file); // this will take the path if exist in the list otherwise it will return automatically the working dir
1394                  };
1395                    char *fullpath = gSystem->ConcatFileName(gSystem->DirName(filedir),gSystem->BaseName(file));
1396                  contents->Add(new TSystemFile(fullpath,gSystem->DirName(fullpath)));// add file to the list                  contents->Add(new TSystemFile(fullpath,gSystem->DirName(fullpath)));// add file to the list
1397                  delete fullpath;                  delete fullpath;
1398              }              }
# Line 1048  TList*  PamLevel2::GetListOfLevel2Files( Line 1401  TList*  PamLevel2::GetListOfLevel2Files(
1401  //          };  //          };
1402          };          };
1403          in.close();          in.close();
1404                  } else {
1405              if(flisttxt.Contains("#"))flisttxt = flisttxt(0,flisttxt.First("#"));
1406              char *fullpath = gSystem->ConcatFileName(gSystem->DirName(flisttxt),gSystem->BaseName(flisttxt));
1407              contents->Add(new TSystemFile(fullpath,gSystem->DirName(fullpath)));// add file to the list
1408              delete fullpath;
1409          };        
1410      }else{      }else{
1411                    
1412          cout << "No input file list given."<<endl;          cout << "No input file list given."<<endl;
1413          cout << "Check for existing root files."<<endl;          cout << "Check for existing root files."<<endl;
1414  //              cout << "Warking directory: "<< gSystem->WorkingDirectory()<< endl;  //              cout << "Warking directory: "<< gSystem->WorkingDirectory()<< endl;
1415                    if ( ddir == "" ){
1416              ddir = wdir;
1417              cout << "Level2 data directory : "<<  ddir << endl;
1418            };
1419    
1420          TSystemDirectory *datadir = new TSystemDirectory(gSystem->BaseName(ddir),ddir);          TSystemDirectory *datadir = new TSystemDirectory(gSystem->BaseName(ddir),ddir);
1421          TList *temp = datadir->GetListOfFiles();          TList *temp = datadir->GetListOfFiles();
1422  //              temp->Print();  //              temp->Print();
# Line 1098  TList*  PamLevel2::GetListOfLevel2Files( Line 1460  TList*  PamLevel2::GetListOfLevel2Files(
1460   * @return Pointer to a TChain   * @return Pointer to a TChain
1461   */   */
1462  TChain *PamLevel2::GetPamTree(TList *fl, TString detlist ){  TChain *PamLevel2::GetPamTree(TList *fl, TString detlist ){
1463      //
1464      //
1465      //
1466      if ( pam_tree ){
1467        printf("WARNING: TChain *PamLevel2::GetPamTree(TList *fl, TString detlist) -- pam_tree already exists!\n ");
1468        return pam_tree;
1469      };
1470      //
1471            
1472      TChain *Trout =0;      TChain *Trout =0;
1473    
# Line 1118  TChain *PamLevel2::GetPamTree(TList *fl, Line 1488  TChain *PamLevel2::GetPamTree(TList *fl,
1488      if(ORB)cout<<"ORB ";      if(ORB)cout<<"ORB ";
1489      cout << endl;      cout << endl;
1490    
1491      TChain *T = 0;             TChain *T = 0;    
1492      TChain *C = 0;       TChain *C = 0;
1493      TChain *O = 0;       TChain *O = 0;
1494      TChain *R = 0;       TChain *R = 0;
1495      TChain *S = 0;       TChain *S = 0;
1496      TChain *N = 0;       TChain *N = 0;
1497      TChain *A = 0;       TChain *A = 0;
1498      TChain *B = 0;       TChain *B = 0;
1499    
1500      TChain *L = 0;       TChain *L = 0;
1501            
1502      if(TRK2||TRK1||TRKh) T = new TChain("Tracker");          if(TRK2||TRK1||TRKh) T = new TChain("Tracker");    
1503      if(CAL2||CAL1)       C = new TChain("Calorimeter");      if(CAL2||CAL1)       C = new TChain("Calorimeter");
# Line 1137  TChain *PamLevel2::GetPamTree(TList *fl, Line 1507  TChain *PamLevel2::GetPamTree(TList *fl,
1507      if(ND)               N = new TChain("NeutronD");      if(ND)               N = new TChain("NeutronD");
1508      if(AC)               A = new TChain("Anticounter");      if(AC)               A = new TChain("Anticounter");
1509      if(ORB)              B = new TChain("OrbitalInfo");      if(ORB)              B = new TChain("OrbitalInfo");
   
1510       L = new TChain("SelectionList");       L = new TChain("SelectionList");
1511            
1512      // loop over files and create chains              // loop over files and create chains        
# Line 1157  TChain *PamLevel2::GetPamTree(TList *fl, Line 1526  TChain *PamLevel2::GetPamTree(TList *fl,
1526              if(ORB)              B->Add(name);              if(ORB)              B->Add(name);
1527              if(SELLI==1)         L->Add(name);              if(SELLI==1)         L->Add(name);
1528          };          };
1529      }      };
1530            
1531      cout << "done chain \n";      cout << "done chain \n";
1532        cout << "+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+" <<endl;
1533    
1534  //    UInt_t *found=0;  //    UInt_t *found=0;
1535  // Tracker  // Tracker
# Line 1174  TChain *PamLevel2::GetPamTree(TList *fl, Line 1544  TChain *PamLevel2::GetPamTree(TList *fl,
1544  //      else    T->SetBranchStatus("TrkHough",0,found);  //      else    T->SetBranchStatus("TrkHough",0,found);
1545          if(TRKh)cout << "Tracker      : set branch address TrkHough"<<endl;          if(TRKh)cout << "Tracker      : set branch address TrkHough"<<endl;
1546          if(!Trout)Trout=T;          if(!Trout)Trout=T;
1547          else Trout->AddFriend("Tracker");          else Trout->AddFriend("Tracker");
1548      }else{      }else{
1549          cout << "Tracker      : missing tree"<<endl;          cout << "Tracker      : missing tree"<<endl;
1550      };      };
# Line 1245  TChain *PamLevel2::GetPamTree(TList *fl, Line 1615  TChain *PamLevel2::GetPamTree(TList *fl,
1615      }else{      }else{
1616          cout << "OrbitalInfo  : missing tree"<<endl;          cout << "OrbitalInfo  : missing tree"<<endl;
1617      };      };
1618        cout << "+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+" <<endl;
1619    
1620      // Selection List      // Selection List
1621      if(L && SELLI==1) {      if(L && SELLI==1) {
# Line 1256  TChain *PamLevel2::GetPamTree(TList *fl, Line 1627  TChain *PamLevel2::GetPamTree(TList *fl,
1627          sel_tree = L;          sel_tree = L;
1628  //      if(!Trout)Trout=O;  //      if(!Trout)Trout=O;
1629  //      else Trout->AddFriend("SelectionList");  //      else Trout->AddFriend("SelectionList");
1630            cout << "----------------------------------------------------" <<endl;
1631      }else{      }else{
1632  //      cout << "SelectionList  : missing tree"<<endl;  //      cout << "SelectionList  : missing tree"<<endl;
1633          if(L)L->Delete();          if(L)L->Delete();
# Line 1447  void PamLevel2::SetBranchAddress(TChain Line 1819  void PamLevel2::SetBranchAddress(TChain
1819   * @return Pointer to a TChain   * @return Pointer to a TChain
1820   */   */
1821  TChain *PamLevel2::GetRunTree(TList *fl){  TChain *PamLevel2::GetRunTree(TList *fl){
1822              //
1823      TChain *R = new TChain("Run");    //
1824      //
1825      if ( run_tree ){
1826        printf("WARNING: TChain *PamLevel2::GetRunTree(TList *fl) -- run_tree already exists!\n ");
1827        return run_tree;
1828      };    
1829      //
1830      TChain *R = new TChain("Run");
1831            
1832      // loop over files and create chains              // loop over files and create chains        
1833      TIter next(fl);      TIter next(fl);
# Line 1487  TChain *PamLevel2::GetRunTree(TList *fl) Line 1866  TChain *PamLevel2::GetRunTree(TList *fl)
1866   * @return Pointer to a TTree   * @return Pointer to a TTree
1867   */   */
1868  TTree *PamLevel2::GetRunTree(TFile *f){  TTree *PamLevel2::GetRunTree(TFile *f){
1869      if ( run_tree ){
1870        printf("WARNING: TTree *PamLevel2::GetRunTree(TFile *f) -- run_tree already exists!\n ");
1871        return run_tree;
1872      };
1873            
1874    
1875      cout << "TTree *PamLevel2::GetRunTree(TFile *f) -- obsolte "<<endl;      cout << "TTree *PamLevel2::GetRunTree(TFile *f) -- obsolte "<<endl;
1876    
1877      TTree *R = (TTree*)f->Get("Run");      TTree *T = (TTree*)f->Get("Run");
1878            
1879      if(R){      if(T){
1880          R->SetBranchAddress("RunInfo", GetPointerTo("RunInfo"));          T->SetBranchAddress("RunInfo", GetPointerTo("RunInfo"));
1881          cout << "Run         : set branch address RunInfo"<<endl;          cout << "Run         : set branch address RunInfo"<<endl;
1882          R->SetBranchAddress("SoftInfo", GetPointerTo("SoftInfo")); // Emiliano          T->SetBranchAddress("SoftInfo", GetPointerTo("SoftInfo")); // Emiliano
1883          cout << "Software    : set branch address SoftInfo"<<endl; // Emiliano          cout << "Software    : set branch address SoftInfo"<<endl; // Emiliano
1884      }      }
1885      return R;  
1886        run_tree = (TChain*)T;
1887    
1888        return T;
1889            
1890  }  }
1891  /**  /**
# Line 1506  TTree *PamLevel2::GetRunTree(TFile *f){ Line 1893  TTree *PamLevel2::GetRunTree(TFile *f){
1893   * @param run Pointer to the chain/tree which contains run infos   * @param run Pointer to the chain/tree which contains run infos
1894   * @return true if a new run has been read, false if it is still the same run   * @return true if a new run has been read, false if it is still the same run
1895   */   */
1896  Bool_t PamLevel2::UpdateRunInfo(TChain *run, ULong64_t iev){  Bool_t PamLevel2::UpdateRunInfo(TChain *run, Long64_t iev){
1897    //    //
1898    // check if we have already called once GetEntry, if not call it    // check if we have already called once GetEntry, if not call it
1899    //    //
1900    if ( run->GetEntries() <= 0 ) return(false);      cout << "Bool_t PamLevel2::UpdateRunInfo(TChain *run, Long64_t iev) --- NON FUNZIONA BENISSIMO.... "<<endl;
1901        if(!run){
1902              cout << "Bool_t PamLevel2::UpdateRunInfo(TChain *run, ULong64_t iev) -- ERROR -- missing RunInfo tree "<<endl;
1903              return(false);        
1904        }
1905        if ( run->GetEntries() <= 0 ) return(false);
1906    //    //
1907        
1908    Int_t oldrun = irun;    //  Int_t oldrun = irun;
1909      Long64_t oldrun = irun;
1910    
1911    // --------------------------------------    // --------------------------------------
1912    // if it is a full file (not preselected)    // if it is a full file (not preselected)
1913    // --------------------------------------    // --------------------------------------
1914    if(SELLI==0){    if(SELLI==0){
       if ( irun < 0 ){  
           irun = 0;  
           run->GetEntry(irun);  
           runfirstentry = 0ULL;  
           runlastentry += (ULong64_t)(this->GetRunInfo()->NEVENTS) - 1ULL;  
           return(true);  
       };        
1915    
1916          //
1917          // the absolute time is necessary to relate the event with the run
1918          //
1919        if( !GetOrbitalInfo() ){        if( !GetOrbitalInfo() ){
1920            cout << "Bool_t PamLevel2::UpdateRunInfo(TChain *run, ULong64_t iev) -- ERROR -- missing OrbitalInfo "<<endl;            cout << "Bool_t PamLevel2::UpdateRunInfo(TChain *run, ULong64_t iev) -- ERROR -- missing OrbitalInfo "<<endl;
1921            return(false);            return(false);
1922        }        }
1923    
1924          //
1925          // the first time the routine is called, set run search from the beginning
1926          //
1927          if ( irun < 0LL ){
1928              irun = 0LL;
1929              run->GetEntry(irun);
1930              runfirstentry = 0LL;
1931              runlastentry += (Long64_t)(this->GetRunInfo()->NEVENTS);
1932              if ( (Long64_t)(this->GetRunInfo()->NEVENTS) > 0LL ) runlastentry -= 1LL;
1933          };      
1934          //
1935           if ( irun == run->GetEntries()-1LL &&
1936               !(GetOrbitalInfo()->absTime >= GetRunInfo()->RUNHEADER_TIME &&
1937                 GetOrbitalInfo()->absTime <= GetRunInfo()->RUNTRAILER_TIME)
1938                ){
1939             irun = -1LL;
1940             runfirstentry = 0LL;
1941             runlastentry = -1LL;
1942           };
1943        // modificato il controllo sull'aggiornamento del run, per evitare problemi        // modificato il controllo sull'aggiornamento del run, per evitare problemi
1944        // dovuti agli eventi annidati (NB! NEVENTS conta anche questi!!)        // dovuti agli eventi annidati (NB! NEVENTS conta anche questi!!)
1945        while ( GetOrbitalInfo()->absTime > GetRunInfo()->RUNTRAILER_TIME && irun < run->GetEntries() ){        //
1946          bool fromfirst = true;
1947          //
1948          while ( !(GetOrbitalInfo()->absTime >= GetRunInfo()->RUNHEADER_TIME && GetOrbitalInfo()->absTime <= GetRunInfo()->RUNTRAILER_TIME) && irun < run->GetEntries()-1LL ){
1949  //      while ( iev > (runfirstentry+(ULong64_t)(this->GetRunInfo()->NEVENTS-1)) && irun < run->GetEntries() ){  //      while ( iev > (runfirstentry+(ULong64_t)(this->GetRunInfo()->NEVENTS-1)) && irun < run->GetEntries() ){
           //    printf(" iev %llu %u %llu \n",iev,this->GetRunInfo()->NEVENTS,(ULong64_t)(runfirstentry+(ULong64_t)(this->GetRunInfo()->NEVENTS)));  
1950            irun++;            irun++;
1951            run->GetEntry(irun);            run->GetEntry(irun);
1952            runfirstentry = runlastentry+1ULL;            runfirstentry = runlastentry;
1953            runlastentry += (ULong64_t)(this->GetRunInfo()->NEVENTS);            if ( (Long64_t)(this->GetRunInfo()->NEVENTS) > 0LL ) runfirstentry += 1LL;
1954              runlastentry += (Long64_t)(this->GetRunInfo()->NEVENTS);
1955    //        cout << " ))))) UPDATE RUN INFO (((((  @iev "<<iev<<" run "<<GetRunInfo()->ID<<" irun "<<irun<<endl;
1956    //        cout << "runfirstentry "<<runfirstentry<<endl;
1957              //      printf(" iev %llu %u %llu \n",iev,this->GetRunInfo()->NEVENTS,(ULong64_t)(runfirstentry+(ULong64_t)(this->GetRunInfo()->NEVENTS)));
1958              //      printf(" abstime %u trailertime %u \n",GetOrbitalInfo()->absTime,GetRunInfo()->RUNTRAILER_TIME);
1959              //      printf(" IDRUN %u \n",GetRunInfo()->ID);
1960              //
1961    //        prevshift = 0;
1962              //
1963              if ( irun == (Long64_t)(run->GetEntries()-1LL) && fromfirst && !(GetOrbitalInfo()->absTime >= GetRunInfo()->RUNHEADER_TIME && GetOrbitalInfo()->absTime <= GetRunInfo()->RUNTRAILER_TIME)){
1964                  printf(" resetting irun  (it should NOT happen!!!)\n");
1965                  fromfirst = false;
1966                  irun = 0;
1967                  run->GetEntry(irun);
1968                  runfirstentry = 0ULL;
1969                  runlastentry += (Long64_t)(this->GetRunInfo()->NEVENTS);
1970                  if ( (Long64_t)(this->GetRunInfo()->NEVENTS) > 0LL ) runlastentry -= 1LL;
1971              };
1972              //
1973        };        };
         
1974        //        //
1975        if        if (
1976   ( irun == oldrun || irun >= run->GetEntries() ) return(false);             !(GetOrbitalInfo()->absTime >= GetRunInfo()->RUNHEADER_TIME &&
1977                 GetOrbitalInfo()->absTime <= GetRunInfo()->RUNTRAILER_TIME)
1978              ) {
1979              printf(" Something very wrong here: cannot find RUN containing absolute time %u \n",GetOrbitalInfo()->absTime);
1980              return false;
1981          }
1982          //
1983          if ( irun == oldrun || irun >= run->GetEntries() ) return(false);
1984        //        //
1985        //  printf(" iev %llu irun %i nevents %u 1st %llu last %llu \n",iev,irun,this->GetRunInfo()->NEVENTS,(ULong64_t)runfirstentry,(ULong64_t)runlastentry);        //  printf(" iev %llu irun %i nevents %u 1st %llu last %llu \n",iev,irun,this->GetRunInfo()->NEVENTS,(ULong64_t)runfirstentry,(ULong64_t)runlastentry);
1986        //        //
1987          prevshift = 0;
1988          cout << " ))))) UPDATE RUN INFO (((((  @iev "<<iev<<" run "<<GetRunInfo()->ID<<" irun "<<irun<<endl;
1989    //      cout << "runfirstentry "<<runfirstentry<<endl;
1990        return(true);            return(true);    
1991    };    };
1992    // ----------------------------------------------------    // ----------------------------------------------------
# Line 1555  Bool_t PamLevel2::UpdateRunInfo(TChain * Line 1995  Bool_t PamLevel2::UpdateRunInfo(TChain *
1995    // ----------------------------------------------------    // ----------------------------------------------------
1996    if(SELLI==1){          if(SELLI==1){      
1997        sel_tree->GetEntry(iev);        sel_tree->GetEntry(iev);
1998    //      cout << irun << " "<< irunentry << endl;
1999        if(irun != oldrun){        if(irun != oldrun){
2000            run->GetEntry(irun);            run->GetEntry(irun);
2001              cout << " ))))) UPDATE RUN INFO (((((  @iev "<<iev<<" run "<<GetRunInfo()->ID<<" irun "<<irun<<endl;
2002              prevshift = 0;
2003            return true;            return true;
2004        }        }
2005        return false;        return false;
# Line 1565  Bool_t PamLevel2::UpdateRunInfo(TChain * Line 2008  Bool_t PamLevel2::UpdateRunInfo(TChain *
2008    return false;    return false;
2009    //    //
2010  };  };
2011    
2012    Bool_t PamLevel2::UpdateRunInfo(Long64_t iev){
2013    
2014        if(!run_tree){
2015            cout << " Bool_t PamLevel2::UpdateRunInfo(ULong64_t iev) -- ERROR -- run tree not loaded"<<endl;
2016            return false;
2017        }
2018        if ( run_tree->GetEntries() <= 0 ) {
2019            cout << " Bool_t PamLevel2::UpdateRunInfo(ULong64_t iev) -- ERROR -- run tree is empty"<<endl;
2020            return(false);
2021        }
2022      
2023        Int_t oldrun = irun; // store current run index
2024    
2025        // +-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-
2026        // if it is a full file (not preselected)
2027        // +-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-
2028        if(SELLI==0){
2029    
2030            // ---------------------------------------------------------------
2031            // the absolute time is necessary to relate the event with the run
2032            // ---------------------------------------------------------------
2033            if( !GetOrbitalInfo() ){
2034                cout << "Bool_t PamLevel2::UpdateRunInfo(Long64_t iev) -- ERROR -- missing OrbitalInfo "<<endl;
2035                return(false);
2036            }
2037          
2038            // -----------------------------------------------------------------------
2039            // the first time the routine is called, set run search from the beginning
2040            // -----------------------------------------------------------------------
2041            if ( irun < 0 ){
2042                irun = 0LL;
2043                irunentry = 0;
2044                prevshift = 0;
2045                run_tree->GetEntry(irun);
2046            };      
2047            //
2048            bool fromfirst = true; // first loop over runs
2049    
2050    //      Bool_t hasfrag = false;
2051    //      if( GetRunInfo()->ID_RUN_FRAG!=0 && GetRunInfo()->ID_RUN_FRAG != GetRunInfo()->ID)hasfrag=true;
2052    //      ULong64_t fragid = GetRunInfo()->ID_RUN_FRAG;
2053    
2054            // ------------------------------------------------------
2055            // loop over runs to find the one that contains the event
2056            // ------------------------------------------------------
2057            while (
2058                (
2059                    (
2060                        !(GetOrbitalInfo()->absTime >= GetRunInfo()->RUNHEADER_TIME &&    // check on absolute time (s)
2061                          GetOrbitalInfo()->absTime <= GetRunInfo()->RUNTRAILER_TIME) &&
2062                        !(GetOrbitalInfo()->OBT >= GetRunInfo()->RUNHEADER_OBT &&         // additional check on OBT (ms)
2063                          GetOrbitalInfo()->OBT <= GetRunInfo()->RUNTRAILER_OBT)
2064                        )
2065                    || GetRunInfo()->NEVENTS==0
2066                    || !(irunentry < GetRunInfo()->NEVENTS-1-prevshift)
2067                    )
2068                    && irun < run_tree->GetEntries() ){
2069    
2070                irun++;
2071                // ------------------------------------
2072                // if the end of run tree is reached...
2073                // ------------------------------------
2074                if( irun == run_tree->GetEntries() ){
2075                    if(!fromfirst){
2076                        // -----------------------------------------------------
2077                        // if it happened already once and the run was not found
2078                        // ---> exit with error
2079                        // -----------------------------------------------------
2080                        cout << "Bool_t PamLevel2::UpdateRunInfo(Long64_t iev) -- ERROR -- event entry #"<<iev<<" does not belong to any run (should not happen)" <<endl;
2081                        return false;
2082                    }
2083                    // -----------------------------------------
2084                    // ...otherwise repeat search from beginning
2085                    // -----------------------------------------
2086                    cout << "Bool_t PamLevel2::UpdateRunInfo(Long64_t iev) -- WARNING -- reached end of run tree. searchin again from beginning " <<endl;
2087                    fromfirst = false;
2088                    irun = 0LL;
2089                    runfirstentry = 0LL;
2090                }
2091                // -------------------------------------------------------------------
2092                // save the index of the first entry of the run, relative to pam_tree,
2093                // and read a new run
2094                // -------------------------------------------------------------------
2095                if(irun>0)runfirstentry += (GetRunInfo()->NEVENTS)-prevshift;
2096                irunentry = 0;
2097                prevshift = 0;
2098                run_tree->GetEntry(irun);      
2099                if(GetRunInfo()->RUNHEADER_OBT>GetRunInfo()->RUNTRAILER_OBT ){
2100                    cout << "Bool_t PamLevel2::UpdateRunInfo(Long64_t iev) -- WARNING -- irun "<<irun<<"  has RUNHEADER_OBT>=RUNTRAILER_OBT " <<endl;
2101                    cout << "                                                            (NB!! in this case some events are assigned to a wrong run)"<<endl;
2102                }
2103    //          if(hasfrag &&  fragid != GetRunInfo()->ID){
2104    //              cout << "... where is the next fragment ??"<<endl;
2105    //          }
2106            };
2107    
2108    
2109            // --------------------------------------
2110            // if there was no need to update the run
2111            // ---> exit with FALSE
2112            // --------------------------------------
2113            if ( irun == oldrun ) return(false);
2114    
2115            // --------------------------------------
2116            // ... otherwise
2117            // ---> exit with TRUE
2118            // --------------------------------------
2119            cout << endl << " ))))) UPDATE RUN INFO (((((  @iev "<<iev<<" run "<<GetRunInfo()->ID<<" irun "<<irun<<endl;
2120            // ----------------------------------------------------
2121            // then check if the run has a fragment
2122            // in this case we have to switch to the next fragment
2123            // when the end of the first fragment is reached
2124            // ----------------------------------------------------
2125            if(
2126                GetRunInfo()->ID_RUN_FRAG != 0 &&
2127    //          GetRunInfo()->ID_RUN_FRAG != GetRunInfo()->ID &&
2128                true ){
2129                cout << "* fragment *"<<endl;              
2130            }
2131            
2132            return(true);    
2133        };
2134        // +-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-
2135        // if it is a preselected file (there is SelectionList)
2136        // +-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-
2137        if(SELLI==1){      
2138            sel_tree->GetEntry(iev);
2139            if(irun != oldrun){
2140                run_tree->GetEntry(irun);
2141                cout << endl << " ))))) UPDATE RUN INFO (((((  @iev "<<iev<<" run "<<GetRunInfo()->ID<<" irun "<<irun<<endl;
2142                prevshift = 0;
2143                return true;
2144            }
2145            return false;
2146        }
2147    
2148        return false;
2149        //
2150    };
2151  /**  /**
2152   * Update the runinfo informations (to be used to have Run infos event by event basis)   * Update the runinfo informations (to be used to have Run infos event by event basis)
2153   * @param run Pointer to the chain/tree which contains run infos   * @param run Pointer to the chain/tree which contains run infos
2154   * @return true if a new run has been read, false if it is still the same run   * @return true if a new run has been read, false if it is still the same run
2155   */   */
2156  Bool_t PamLevel2::UpdateRunInfo(TTree *run, ULong64_t iev){  Bool_t PamLevel2::UpdateRunInfo(TTree *run, Long64_t iev){
2157    return(UpdateRunInfo((TChain*)run,iev));    return(UpdateRunInfo((TChain*)run,iev));
2158  };  };
2159    
# Line 2194  void PamLevel2::CreateCloneTrees0( TChai Line 2777  void PamLevel2::CreateCloneTrees0( TChai
2777   */   */
2778  void PamLevel2::CreateCloneTrees(TFile *ofile){  void PamLevel2::CreateCloneTrees(TFile *ofile){
2779    
2780    
2781    //  if the pointer is null, create a default file
2782    
2783        if(!ofile){
2784            cout << "void PamLevel2::CreateCloneTrees(TFile*) -- WARNING -- Creating file: clone-tree.root "<<endl;
2785            ofile = new TFile("clone-tree.root","recreate");
2786        }
2787    
2788      ofile->cd();      ofile->cd();
2789    
2790      cout << "+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+" <<endl;      cout << "+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+" <<endl;
# Line 2205  void PamLevel2::CreateCloneTrees(TFile * Line 2796  void PamLevel2::CreateCloneTrees(TFile *
2796      cout << "Run          : branch RunInfo"<<endl;      cout << "Run          : branch RunInfo"<<endl;
2797      run_tree_clone->Branch("SoftInfo","SoftInfo",GetPointerTo("SoftInfo"));      run_tree_clone->Branch("SoftInfo","SoftInfo",GetPointerTo("SoftInfo"));
2798      cout << "Run          : branch SoftInfo"<<endl;      cout << "Run          : branch SoftInfo"<<endl;
2799                    // ------------------
2800        // replicate run tree
2801        // ------------------
2802        cout << "----------------------------------------------------"<<endl;
2803        cout << "irun\t | RUN\t NEVENTS\t absolute time"<<endl;
2804        for (Int_t i=0; i<run_tree->GetEntries(); i++){
2805            run_tree->GetEntry(i);
2806            cout << i<< "\t | "<<GetRunInfo()->ID<<"\t "<<GetRunInfo()->NEVENTS<< "\t "<<GetRunInfo()->RUNHEADER_TIME<<" <---> "<<GetRunInfo()->RUNTRAILER_TIME<<endl;
2807            run_tree_clone->Fill();
2808        }
2809        cout << "----------------------------------------------------"<<endl;
2810    
2811    
2812      sel_tree_clone = new TTree("SelectionList","List of selected events ");      sel_tree_clone = new TTree("SelectionList","List of selected events ");
2813      sel_tree_clone->Branch("RunEntry",&irun,"runentry/I");      sel_tree_clone->Branch("RunEntry",&irun,"runentry/L");
2814      sel_tree_clone->Branch("EventEntry",&irunentry,"eventry/I");      sel_tree_clone->Branch("EventEntry",&irunentry,"eventry/L");
2815            
2816    
2817      Int_t i=0;      Int_t i=0;
# Line 2220  void PamLevel2::CreateCloneTrees(TFile * Line 2821  void PamLevel2::CreateCloneTrees(TFile *
2821              pam_tree_clone[i]->Branch("TrkLevel1","TrkLevel1", GetPointerTo("TrkLevel1"));              pam_tree_clone[i]->Branch("TrkLevel1","TrkLevel1", GetPointerTo("TrkLevel1"));
2822              pam_tree_clone[i]->BranchRef();              pam_tree_clone[i]->BranchRef();
2823              cout << "Tracker      : branch TrkLevel1"<<endl;              cout << "Tracker      : branch TrkLevel1"<<endl;
2824              cout << "CreateCloneTrees " << GetTrkLevel1()<<endl;  //          cout << "CreateCloneTrees " << GetTrkLevel1()<<endl;
2825          };          };
2826          if(TRK2) {          if(TRK2) {
2827              pam_tree_clone[i]->Branch("TrkLevel2", "TrkLevel2",GetPointerTo("TrkLevel2"));              pam_tree_clone[i]->Branch("TrkLevel2", "TrkLevel2",GetPointerTo("TrkLevel2"));
# Line 2350  void PamLevel2::WriteCloneTrees(){ Line 2951  void PamLevel2::WriteCloneTrees(){
2951  /**  /**
2952   * Method to get level2-trees entry, the corresponding run entry and (if required) the level0 entry.   * Method to get level2-trees entry, the corresponding run entry and (if required) the level0 entry.
2953   */   */
2954  Int_t PamLevel2::GetEntry(Int_t iee){  //Int_t PamLevel2::GetEntry(Int_t iee){
2955    Int_t PamLevel2::GetEntry(Long64_t iee){
2956            
2957    //     cout << "-------------------------------------"<<endl;
2958    //     cout << "Int_t PamLevel2::GetEntry("<<iee<<")"<<endl;
2959    
2960      if(!pam_tree){      if(!pam_tree){
2961          cout << " Int_t PamLevel2::GetEntry(Int_t) -- ERROR -- level2 trees not loeaded"<<endl;          cout << " Int_t PamLevel2::GetEntry(Int_t) -- ERROR -- level2 trees not loaded"<<endl;
         return 0;  
     }  
     if(!run_tree ){  
         cout << " Int_t PamLevel2::GetEntry(Int_t) -- ERROR -- run tree not loeaded"<<endl;  
2962          return 0;          return 0;
2963      }      }
2964    
2965      Int_t ii=0;  
2966        //
2967        // 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...
2968        //
2969        //    if(!run_tree ){
2970        //  cout << " Int_t PamLevel2::GetEntry(Int_t) -- ERROR -- run tree not loeaded"<<endl;
2971        //  return 0;
2972        //    }
2973    
2974        Long64_t ii=0;
2975      //-------------------------------      //-------------------------------
2976      ii = iee;      ii = iee;
2977      if( !pam_tree->GetEntry(ii) ) return 0;      if( !pam_tree->GetEntry(ii) ){      
2978            cout << " Int_t PamLevel2::GetEntry(Int_t) -- ERROR -- error reading pam tree"<<endl;
2979            return 0;
2980        }
2981        //
2982        // ... 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.
2983        // 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
2984        // a problem if you don't check the return code of getentry.
2985        //
2986        if(!run_tree ){
2987            if ( TRK0 || CAL0 || TOF0 || RUN ) { //forse cosi` va bene per tornare 1?
2988                cout << " Int_t PamLevel2::GetEntry(Int_t) -- ERROR -- run tree not loaded"<<endl;
2989                return 0;
2990            }  else {
2991                return 1; //cosi` se non c'e` run esce qua...
2992            }
2993        }
2994    
2995      //-------------------------------      //-------------------------------
2996      ii = iee;      ii = iee;
2997      Bool_t UPDATED = UpdateRunInfo(run_tree,ii);  //    Bool_t UPDATED = UpdateRunInfo(run_tree,ii);
2998    //    Bool_t UPDATED = UpdateRunInfo(ii);
2999        UpdateRunInfo(ii);
3000      if(SELLI==0)irunentry = iee-runfirstentry;      if(SELLI==0)irunentry = iee-runfirstentry;
3001      if(UPDATED && run_tree_clone)run_tree_clone->Fill();  //    if(UPDATED && run_tree_clone)run_tree_clone->Fill();
3002      
3003  //    cout << "PamLevel2::GetEntry("<<iee<<") "<<irun<<" "<<runfirstentry<<" "<<irunentry<<endl;  //    cout << "PamLevel2::GetEntry("<<iee<<") "<<irun<<" "<<runfirstentry<<" "<<irunentry<<endl;
3004    
3005      if( TRK0 || CAL0 || TOF0 )GetYodaEntry( );  //     cout << "irunentry     "<<irunentry << endl;
3006    //     cout << "runfirstentry "<<runfirstentry << endl;
3007    //     cout << "nevents       "<<GetRunInfo()->NEVENTS<< endl;
3008    
3009    //     if( TRK0 || CAL0 || TOF0 ){
3010    //      if( !GetYodaEntry( ) ){
3011    //          cout << " Int_t PamLevel2::GetEntry(Int_t) -- ERROR -- error reading level0 tree"<<endl;
3012    //          return 0;
3013    //      }
3014    //     }
3015    
3016    
3017      return 1;      return 1;
3018    
3019  }  }
3020    
3021    TrkLevel0    *PamLevel2::GetTrkLevel0(){
3022        if( !TRK0 )return NULL;
3023        if( !GetYodaEntry( ) ){
3024            cout << " Int_t PamLevel2::GetTrkLevel0() -- ERROR -- error reading level0 tree"<<endl;
3025            return 0;
3026        }
3027        return trk0_obj;
3028    };
3029    CaloLevel0    *PamLevel2::GetCaloLevel0(){
3030        if( !CAL0 )return NULL;
3031        if( !GetYodaEntry( ) ){
3032            cout << " Int_t PamLevel2::GetCaloLevel0() -- ERROR -- error reading level0 tree"<<endl;
3033            return 0;
3034        }
3035        return calo0_obj;
3036    };
3037    
3038    
3039  /**  /**
3040   * Method to retrieve the level0 tree (YODA tree) that contains the current event.   * Method to retrieve the level0 tree (YODA tree) that contains the current event.
3041   * Given the run ID (...), if needed it query the DB and load the proper file.   * Given the run ID (...), if needed it query the DB and load the proper file.
# Line 2393  TTree* PamLevel2::GetYodaTree( ){ Line 3050  TTree* PamLevel2::GetYodaTree( ){
3050      // check if iroot has changed      // check if iroot has changed
3051      //===================================      //===================================
3052      if( irun<0 ){      if( irun<0 ){
3053          cout << "PamLevel2::GetYodaTree() -- ERROR "<<endl;          cout << "PamLevel2::GetYodaTree() -- ERROR -- irun = "<<irun<<endl;
3054          cout << "In order to use this method you have to load the RunInfo tree "<<endl;  //      cout << "In order to use this method you have to first load the RunInfo tree "<<endl;
         cout << "with the method TChain* PamLevel2::GetRunTree(TString, TString)"<<endl;  
 //      cout << " - read the event with PamLevel2::GetEntry(Int_t)"<<endl;  
3055          return NULL;          return NULL;
3056      }          }    
3057      Int_t irootnew = run_obj->ID_ROOT_L0;      Int_t irootnew = run_obj->ID_ROOT_L0;
3058  //     cout << "iroot    "<<iroot<<endl;  //      cout << "iroot    "<<iroot<<endl;
3059  //     cout << "irootnew "<<irootnew<<endl;  //      cout << "irootnew "<<irootnew<<endl;
3060    
3061      //===================================      //===================================
3062      // load the level0 file      // load the level0 file
# Line 2448  TTree* PamLevel2::GetYodaTree( ){ Line 3103  TTree* PamLevel2::GetYodaTree( ){
3103          l0_tree = (TTree*)l0_file->Get("Physics");          l0_tree = (TTree*)l0_file->Get("Physics");
3104          if(!h0_obj)h0_obj = new EventHeader();          if(!h0_obj)h0_obj = new EventHeader();
3105          l0_tree->SetBranchAddress("Header" ,&h0_obj);          l0_tree->SetBranchAddress("Header" ,&h0_obj);
3106            prevshift = 0;
3107          //---------------------------------------------------          //---------------------------------------------------
3108          // TRACKER:          // TRACKER:
3109          if(TRK0){          if(TRK0){
# Line 2457  TTree* PamLevel2::GetYodaTree( ){ Line 3113  TTree* PamLevel2::GetYodaTree( ){
3113              };              };
3114              l0_tree->SetBranchAddress("Tracker" ,trk0_obj->GetPointerToTrackerEvent());              l0_tree->SetBranchAddress("Tracker" ,trk0_obj->GetPointerToTrackerEvent());
3115          }          }
3116          //---------------------------------------------------          //--------------------------------------------------
3117          // CALORIMETER:          // CALORIMETER:
3118          if(CAL0){          if(CAL0){
3119              cout << "PamLevel2::GetYodaTree() --- level0 calorimeter not implemented "<<endl;              if(!calo0_obj){
3120                    calo0_obj = new CaloLevel0();
3121                    calo0_obj->Set();
3122                };
3123                l0_tree->SetBranchAddress("Calorimeter" ,calo0_obj->GetPointerToCalorimeterEvent());
3124                //      cout << "PamLevel2::GetYodaTree() --- level0 calorimeter not implemented "<<endl;
3125          }          }
3126          //---------------------------------------------------          //---------------------------------------------------
3127          // TOF:          // TOF:
# Line 2474  TTree* PamLevel2::GetYodaTree( ){ Line 3135  TTree* PamLevel2::GetYodaTree( ){
3135          cout << " TTree* PamLevel2::GetYodaTree( ) -- no DB connected... hai fatto qualche cazzata "<<endl;          cout << " TTree* PamLevel2::GetYodaTree( ) -- no DB connected... hai fatto qualche cazzata "<<endl;
3136      }      }
3137    
3138      TrkParams::SetCalib(run_obj,dbc);      if ( TRK0 ){
3139            TrkParams::Load(6);
3140            if( !TrkParams::IsLoaded(6) ){
3141                cout << " TTree* PamLevel2::GetYodaTree( ) -- WARNING -- VK-mask not loaded"<<endl;
3142            };
3143            TrkParams::SetCalib(run_obj,dbc);
3144            TrkParams::LoadCalib( );
3145            if( !TrkParams::CalibIsLoaded() ){
3146                cout << " TTree* PamLevel2::GetYodaTree( ) -- WARNING -- Calibration not loaded"<<endl;
3147            };
3148        }
3149    
3150  //    cout << l0_tree << endl;  //    cout << l0_tree << endl;
3151            
# Line 2492  Int_t PamLevel2::GetYodaEntry(){ Line 3162  Int_t PamLevel2::GetYodaEntry(){
3162      if(!GetYodaTree())return 0;      if(!GetYodaTree())return 0;
3163            
3164      // patch      // patch
3165      if( (Int_t)irunentry < 0){      if( irunentry < 0){
3166          cout << "Int_t PamLevel2::GetYodaEntry() -- ATTENZIONE -- irunentry negativo?!?! "<<(Int_t)irunentry<<endl;  //      cout << "Int_t PamLevel2::GetYodaEntry() -- ATTENZIONE -- irunentry negativo?!?! "<<(Int_t)irunentry<<endl;
3167          irunentry=0;          irunentry=0LL;
3168      }      }
3169      //  ---------------------------------      //  ---------------------------------
3170      //  if file is NOT a preselected file      //  if file is NOT a preselected file
3171      //  ---------------------------------      //  ---------------------------------
3172      UInt_t quellagiusta = irunentry + run_obj->EV_FROM;      Long64_t quellagiusta = irunentry + (Long64_t)(run_obj->EV_FROM);
3173        //     cout << " irunentry "<<irunentry<<endl;
3174    //     cout << " EV_FROM "<<run_obj->EV_FROM<<endl;
3175    //     cout << " quellagiusta = irunentry + EV_FROM "<< quellagiusta << endl;
3176    
3177  //     cout << " irun "<< irun << " irunentry "<< irunentry<<" run_obj->EV_FROM "<<run_obj->EV_FROM <<" quella giusta "<<quellagiusta << endl;  //     cout << " irun "<< irun << " irunentry "<< irunentry<<" run_obj->EV_FROM "<<run_obj->EV_FROM <<" quella giusta "<<quellagiusta << endl;
3178  //     cout << " iroot "<<iroot<<" run_obj->ID_ROOT_L0 "<<run_obj->ID_ROOT_L0<<endl;  //     cout << " iroot "<<iroot<<" run_obj->ID_ROOT_L0 "<<run_obj->ID_ROOT_L0<<endl;
# Line 2511  Int_t PamLevel2::GetYodaEntry(){ Line 3183  Int_t PamLevel2::GetYodaEntry(){
3183          cout << "Int_t PamLevel2::GetYodaEntry() -- ERROR -- missing OrbitalInfo "<<endl;          cout << "Int_t PamLevel2::GetYodaEntry() -- ERROR -- missing OrbitalInfo "<<endl;
3184          return 0;          return 0;
3185      }      }
3186        if( GetOrbitalInfo()->OBT==0 && GetOrbitalInfo()->pkt_num==0 ){
3187            cout << "Int_t PamLevel2::GetYodaEntry() -- ERROR -- level2 event corrupted ?? "<<endl;
3188            return 0;
3189        }
3190    
3191      // ---------------------------------------------------------------------      // ---------------------------------------------------------------------
3192      // ATTENTION!!!      // ATTENTION!!!
# Line 2525  Int_t PamLevel2::GetYodaEntry(){ Line 3201  Int_t PamLevel2::GetYodaEntry(){
3201      // ---------------------------------------------------------------------      // ---------------------------------------------------------------------
3202      Int_t answer = 0;      Int_t answer = 0;
3203      Int_t shift =0;      Int_t shift =0;
3204        //    printf(" siamo qui %i %i \n",shift,prevshift);
3205        Int_t maxshift = 10;
3206      do{      do{
3207          if(shift>0){          if(shift>0){    
3208              cout << " level0 <--> level2 mismatch ( irun "<<irun<<" irunentry "<<irunentry<<" shift "<<shift<<" )"<<endl;              cout << " PKTNUM  L2 --- "<< GetOrbitalInfo()->pkt_num << " --- L0 --- "<< GetEventHeader()->GetPscuHeader()->GetCounter()<<endl;
3209                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;
3210                cout << "         L2 <--> L0 mismatch ( irun "<<irun<<" irunentry "<<irunentry<<" shift "<<shift<<" prevshift "<<prevshift<<" )"<<endl;
3211          }          }
3212          answer = l0_tree->GetEntry(quellagiusta+shift);          answer = l0_tree->GetEntry(quellagiusta+(Long64_t)shift+(Long64_t)prevshift);
3213          shift++;          shift++;
3214          if( !GetEventHeader() ){          if( !GetEventHeader() ){
3215              cout << "Int_t PamLevel2::GetYodaEntry() -- ERROR -- missing EventHeader "<<endl;              cout << "Int_t PamLevel2::GetYodaEntry() -- ERROR -- missing EventHeader "<<endl;
3216              return 0;              return 0;
3217          }          }
3218    //      cout << "PKTNUM "<<shift<<" ==  L2 --- "<< GetOrbitalInfo()->pkt_num << " --- L0 --- "<< GetEventHeader()->GetPscuHeader()->GetCounter()<<endl;
3219    //      cout << " L2 --- "<< GetOrbitalInfo()->OBT << " --- L0 --- "<< GetEventHeader()->GetPscuHeader()->GetOrbitalTime()<<endl;
3220  //      if( (quellagiusta+shift) == l0_tree->GetEntries()+1 )cout << ">>> end of level0 tree <<<"<<endl;  //      if( (quellagiusta+shift) == l0_tree->GetEntries()+1 )cout << ">>> end of level0 tree <<<"<<endl;
3221  //      cout << " GetOrbitalInfo()->OBT "<< GetOrbitalInfo()->OBT << endl;  //      cout << " GetOrbitalInfo()->OBT "<< GetOrbitalInfo()->OBT << endl;
3222  //      cout << " GetEventHeader()->GetPscuHeader()->GetOrbitalTime() "<< GetEventHeader()->GetPscuHeader()->GetOrbitalTime() << endl;  //      cout << " GetEventHeader()->GetPscuHeader()->GetOrbitalTime() "<< GetEventHeader()->GetPscuHeader()->GetOrbitalTime() << endl;
3223  //      cout << " GetOrbitalInfo()->pkt_num "<< GetOrbitalInfo()->pkt_num << endl;  //      cout << " GetOrbitalInfo()->pkt_num "<< GetOrbitalInfo()->pkt_num << endl;
3224  //      cout << " GetEventHeader()->GetPscuHeader()->GetCounter() "<< GetEventHeader()->GetPscuHeader()->GetCounter() << endl;  //      cout << " GetEventHeader()->GetPscuHeader()->GetCounter() "<< GetEventHeader()->GetPscuHeader()->GetCounter() << endl;
3225    //      printf(" IDRUN %u \n",GetRunInfo()->ID);
3226    //
3227            if ( prevshift != 0 && (quellagiusta+(Long64_t)shift) == GetYodaTree()->GetEntries() ){
3228              prevshift = 0;
3229              shift = -1;
3230            };
3231    
3232      }while( ( GetOrbitalInfo()->OBT != (UInt_t)(GetEventHeader()->GetPscuHeader()->GetOrbitalTime()) || GetOrbitalInfo()->pkt_num != (UInt_t)(GetEventHeader()->GetPscuHeader()->GetCounter())) && (quellagiusta+shift) < GetYodaTree()->GetEntries() );      }while( ( GetOrbitalInfo()->OBT != (UInt_t)(GetEventHeader()->GetPscuHeader()->GetOrbitalTime()) || GetOrbitalInfo()->pkt_num != (UInt_t)(GetEventHeader()->GetPscuHeader()->GetCounter())) && (quellagiusta+(Long64_t)shift) < GetYodaTree()->GetEntries() && shift < maxshift);
3233    
3234        if ( (quellagiusta+(Long64_t)shift+(Long64_t)prevshift) > GetYodaTree()->GetEntries() || shift == maxshift ) {
3235            cout << " Big trouble here, no such event in Level0 data! " <<endl;
3236            return 0;
3237        }
3238  //    cout << "LA ENTRY GIUSTA E`: "<<quellagiusta<<" (spero...)"<<endl;  //    cout << "LA ENTRY GIUSTA E`: "<<quellagiusta<<" (spero...)"<<endl;
3239  //    return GetYodaTree()->GetEntry(quellagiusta);  //    return GetYodaTree()->GetEntry(quellagiusta);
3240        if ( shift > 1 ) prevshift += (shift-1);
3241            
3242      return answer;      return answer;
3243            
3244  }  }
3245    /**
3246     * \Brief Set DB connection
3247     */
3248    Bool_t PamLevel2::SetDBConnection(){
3249    
3250        cout<<"~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~"<<endl;
3251        cout<<"Connecting to DB"<<endl;
3252        cout<<"HOST "<<host<<endl;
3253        cout<<"USER "<<user<<endl;
3254        cout<<"PSW  "<<psw<<endl;
3255        dbc = TSQLServer::Connect(host.Data(),user.Data(),psw.Data());
3256        if( !dbc )return false;
3257        if( !dbc->IsConnected() )return false;    
3258        cout<<"~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~"<<endl;
3259        return true;
3260    
3261    }
Legend:
Removed from v.1.27  
changed lines
  Added in v.1.41
  ViewVC Help
Powered by ViewVC 1.1.23