/[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.30 by mocchiut, Sun Mar 18 20:25:54 2007 UTC revision 1.43 by pam-fi, Fri Aug 17 16:57:12 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                    
 //     trk2_obj  = TrkLevel2::GetTrkLevel2();  
 //     trk1_obj  = TrkLevel1::GetTrkLevel1();  
 //     trkh_obj  = TrkHough::GetTrkHough();  
 //     calo1_obj = CaloLevel1::GetCaloLevel1();  
 //     calo2_obj = CaloLevel2::GetCaloLevel2();  
 //     tof_obj   = ToFLevel2::GetToFLevel2();  
 //     trig_obj  = TrigLevel2::GetTrigLevel2();  
 //     s4_obj    = S4Level2::GetS4Level2();  
 //     nd_obj    = NDLevel2::GetNDLevel2();  
 //     ac_obj    = AcLevel2::GetAcLevel2();  
 //     orb_obj   = OrbitalInfo::GetOrbitalInfo();  
       
115      h0_obj    = 0;      h0_obj    = 0;
116      trk0_obj  = 0;      trk0_obj  = 0;
117      calo0_obj  = 0;      calo0_obj  = 0;
# Line 119  PamLevel2::PamLevel2(){ Line 130  PamLevel2::PamLevel2(){
130    
131      run_obj   = 0;//new GL_RUN();      run_obj   = 0;//new GL_RUN();
132      soft_obj   = 0;// Emiliano      soft_obj   = 0;// Emiliano
133      irun = -1;      irun = -1LL;
134      runfirstentry = 0ULL;      runfirstentry = 0LL;
135      runlastentry = 0ULL;      runlastentry = 0LL;
136    
137      l0_file = NULL;      l0_file = NULL;
138      l0_tree = NULL;      l0_tree = NULL;
# Line 130  PamLevel2::PamLevel2(){ Line 141  PamLevel2::PamLevel2(){
141    
142      prevshift = 0;      prevshift = 0;
143            
     irun = -1;  
144      run_tree = NULL;      run_tree = NULL;
145      run_tree_clone = NULL;      run_tree_clone = NULL;
146      sel_tree = NULL;      sel_tree = NULL;
147      sel_tree_clone = NULL;      sel_tree_clone = NULL;
148            
149      irunentry = -1;      irunentry = -1LL;
150      pam_tree = NULL;      pam_tree = NULL;
151      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;
152    
# Line 176  PamLevel2::PamLevel2(){ Line 186  PamLevel2::PamLevel2(){
186    
187      SELLI = -1;      SELLI = -1;
188    
189        ISGP = false;
190    
191      tsorted=0;      tsorted=0;
192      timage=0;      timage=0;
193        
194        howtosort = "+CAL+TOF";
195        //howtosort = "+TOF";
196        sortthr = 100.;
197    
198  };  };
199  /**  /**
# Line 219  void PamLevel2::Delete(){ Line 235  void PamLevel2::Delete(){
235      }      }
236            
237      if(l0_file)l0_file->Close();      if(l0_file)l0_file->Close();
238      if(pam_tree)pam_tree->Delete();;      //    if(pam_tree)pam_tree->Delete();;
239    
240      if ( pam_tree ){
241        //
242        // 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...
243        //
244        TList *temp = pam_tree->GetListOfFriends();
245        TList *contents  = new TList; // create chain friend list
246        contents->SetOwner();
247        TIter next(temp);
248        TChain *questo = 0;
249        while ( (questo = (TChain*) next()) ){
250          TString name =  questo->GetName();
251          contents->Add((TChain*)gROOT->FindObject(name.Data()));// add object to the list
252        };
253        //
254        // deleting the main chain
255        //
256        pam_tree->Delete();
257        //
258        // deleting the friends...
259        //
260        TIter next2(contents);
261        TChain *questa = 0;
262        while ( questa = (TChain*)next2() ){
263          TString name =  questa->GetName();
264          questa->Delete();      
265          questa=NULL;
266        };
267        //
268      };
269      pam_tree = NULL;
270    
271      if(run_tree)run_tree->Delete();;      if(run_tree)run_tree->Delete();;
272      if(sel_tree)sel_tree->Delete();;      if(sel_tree)sel_tree->Delete();;
273      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 245  void PamLevel2::Clear(){ Line 293  void PamLevel2::Clear(){
293  //    if(soft_obj) soft_obj->Clear();  //    if(soft_obj) soft_obj->Clear();
294    
295      if(h0_obj)   h0_obj->Clear();            if(h0_obj)   h0_obj->Clear();      
296      if(trk0_obj) trk0_obj->Clear();      //    if(trk0_obj) trk0_obj->Clear();  
     if(calo0_obj) calo0_obj->Clear();    
297      if(trk1_obj) trk1_obj->Clear();          if(trk1_obj) trk1_obj->Clear();    
298      if(trk2_obj) trk2_obj->Clear();      if(trk2_obj) trk2_obj->Clear();
299      if(trkh_obj) trkh_obj->Clear();      if(trkh_obj) trkh_obj->Clear();
300        if(calo0_obj) calo0_obj->Clear();  
301      if(calo1_obj)calo1_obj->Clear();      if(calo1_obj)calo1_obj->Clear();
302      if(calo2_obj)calo2_obj->Clear();      if(calo2_obj)calo2_obj->Clear();
303      if(tof_obj)  tof_obj->Clear();      if(tof_obj)  tof_obj->Clear();
# Line 268  void PamLevel2::Clear(){ Line 316  void PamLevel2::Clear(){
316      if(timage){      if(timage){
317          timage->Delete();          timage->Delete();
318      }      }
       
319  };  };
320    
321    void PamLevel2::Reset(){
322      //
323      // First of all clear everything
324      //
325      Clear();
326      //
327      // close and reset chains and pointers
328      //    
329      if ( pam_tree ){
330        //
331        // 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...
332        //
333        TList *temp = pam_tree->GetListOfFriends();
334        TList *contents  = new TList; // create chain friend list
335        contents->SetOwner();
336        TIter next(temp);
337        TChain *questo = 0;
338        while ( (questo = (TChain*) next()) ){
339          TString name =  questo->GetName();
340          contents->Add((TChain*)gROOT->FindObject(name.Data()));// add object to the list
341        };
342        //
343        // deleting the main chain
344        //
345        pam_tree->Delete();
346        //
347        // deleting the friends...
348        //
349        TIter next2(contents);
350        TChain *questa = 0;
351        while ( questa = (TChain*)next2() ){
352          TString name =  questa->GetName();
353          questa->Delete();      
354          questa=NULL;
355        };
356        //
357      };
358      pam_tree = NULL;
359      //
360      if(run_tree)run_tree->Delete();;
361      run_tree = NULL;
362      if(sel_tree)sel_tree->Delete();;
363      sel_tree = NULL;
364      //
365      // Close file
366      //
367      if(l0_file)l0_file->Close("R");
368      l0_file = NULL;
369      //
370      h0_obj    = 0;
371      trk0_obj  = 0;
372      calo0_obj  = 0;
373      //
374      trk2_obj  = 0;
375      trk1_obj  = 0;
376      trkh_obj  = 0;
377      calo1_obj = 0;
378      calo2_obj = 0;
379      tof_obj   = 0;
380      trig_obj  = 0;
381      s4_obj    = 0;
382      nd_obj    = 0;
383      ac_obj    = 0;
384      orb_obj   = 0;
385      //
386      // Reset run pointers
387      //
388      run_obj   = 0;//new GL_RUN();
389      soft_obj   = 0;// Emiliano
390      irun = -1;
391      runfirstentry = 0ULL;
392      runlastentry = 0ULL;    
393      //
394    };
395    
396    Bool_t PamLevel2::IsGood(){
397      Bool_t goodev=true;
398      //  if(trk2_obj && trk2_obj->UnpackError() != 0 ) goodev = false;
399      if(calo2_obj && calo2_obj->good != 1) goodev = false;
400      if(tof_obj && tof_obj->unpackError != 0) goodev = false;  
401      if(trig_obj && trig_obj->unpackError != 0) goodev = false;
402      if(s4_obj && s4_obj->unpackError != 0) goodev = false;  
403      if(nd_obj && nd_obj->unpackError != 0) goodev = false;  
404      if(ac_obj && ac_obj->unpackError != 255) goodev = false;  
405        //  if(orb_obj)  
406      return goodev;
407    };
408    
409  //--------------------------------------  //--------------------------------------
410  //  //
# Line 465  PamTrack* PamLevel2::GetStoredTrack(Int_ Line 599  PamTrack* PamLevel2::GetStoredTrack(Int_
599  }  }
600  //--------------------------------------  //--------------------------------------
601  //  //
602    
603    /**
604     * 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.
605     **/
606    void PamLevel2::SortTracks(TString how){
607      printf(" WARNING! obsolete, use SortTracks() and SetSortingMethod(TString) instead! \n Setting sorting method to %s \n",how.Data());
608      howtosort = how;  
609      SortTracks();
610    };
611    
612  //  //
613  //--------------------------------------  //--------------------------------------
614  /**  /**
# Line 478  PamTrack* PamLevel2::GetStoredTrack(Int_ Line 622  PamTrack* PamLevel2::GetStoredTrack(Int_
622   *   *
623   * 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).
624   */   */
625  void PamLevel2::SortTracks(TString how){  void PamLevel2::SortTracks(){
626      TString how = howtosort;
627    
628  //    cout <<" PamLevel2::SortTracks(TString how) "<<endl;    //    cout <<" PamLevel2::SortTracks(TString how) "<<endl;
629      if( !trk2_obj ){    if( !trk2_obj ){
630          cout << "void PamLevel2::SortTracks(TString how):  TrkLevel2 not loaded !!!";      cout << "void PamLevel2::SortTracks():  TrkLevel2 not loaded !!!";
631          return;      return;
632      };    };
633  //    cout << "call SortTracs() "<<endl;    //    cout << "call SortTracs() "<<endl;
634    //Save current Object count    //Save current Object count
635      Int_t ObjectNumber = TProcessID::GetObjectCount();    Int_t ObjectNumber = TProcessID::GetObjectCount();
636    
637      //    cout << "ObjectNumber  "<<ObjectNumber <<endl;
638    
639  //    cout << "ObjectNumber  "<<ObjectNumber <<endl;    //     if(!sorted_tracks)sorted_tracks = new TRefArray();
640      //     sorted_tracks->Clear();
641      //    sorted_tracks.Clear();
642    
643  //     if(!sorted_tracks)sorted_tracks = new TRefArray();    if(!tsorted)tsorted = new TClonesArray("PamTrack",trk2_obj->GetNTracks());
644  //     sorted_tracks->Clear();    tsorted->Delete();
645  //    sorted_tracks.Clear();    TClonesArray &ttsorted = *tsorted;
646      if(!timage)timage = new TClonesArray("PamTrack",trk2_obj->GetNTracks());
647      timage->Delete();
648      TClonesArray &ttimage = *timage;
649    
     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;  
650    
651    
652      // loop over the tracks sorted by the tracker    //--------------------------------------------------
653      Bool_t use_TRK = how.Contains("TRK", TString::kIgnoreCase);    // retrieve sorting method
654      Bool_t use_CAL = how.Contains("CAL", TString::kIgnoreCase);    //--------------------------------------------------
655      Bool_t use_TOF = how.Contains("TOF", TString::kIgnoreCase);    Bool_t use_TRK = how.Contains("TRK", TString::kIgnoreCase);
656      Bool_t use_CAL = how.Contains("CAL", TString::kIgnoreCase);
657      Bool_t use_TOF = how.Contains("TOF", TString::kIgnoreCase);
658            
659      if( !CAL2 &&  use_CAL) use_CAL = false;    if( !CAL2 &&  use_CAL) use_CAL = false;
660      if( !TOF &&  use_TOF) use_TOF = false;    if( !TOF &&  use_TOF) use_TOF = false;
661            
662      if( !TRK2 ){    if( !TRK2 ){
663  //      cout << "SortTracks() : without tracker does not work!!! (not yet)" << endl;      //  cout << "SortTracks() : without tracker does not work!!! (not yet)" << endl;
664          return;      return;
665      };    };
666    
667      //   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;
668          
669      for(Int_t i=0; i < trk2_obj->TrkLevel2::GetNTracks(); i++){  
670              //--------------------------------------------------
671          TrkTrack   *ts = 0;    // loop over "physical" tracks sorted by the tracker
672          CaloTrkVar *cs = 0;    //--------------------------------------------------
673          ToFTrkVar  *os = 0;    for(Int_t i=0; i < trk2_obj->TrkLevel2::GetNTracks(); i++){
674                    
675          // get tracker tracks      TrkTrack   *ts = 0;
676          TrkTrack   *tp = trk2_obj->TrkLevel2::GetTrack(i);                    //tracker      CaloTrkVar *cs = 0;
677          CaloTrkVar *cp = GetCaloStoredTrack(tp->GetSeqNo());      ToFTrkVar  *os = 0;
678          ToFTrkVar  *op = GetToFStoredTrack(tp->GetSeqNo());          
679        // get tracker tracks
680          TrkTrack   *ti = 0;              //tracker (image)      TrkTrack   *tp = trk2_obj->TrkLevel2::GetTrack(i); //tracker
681          CaloTrkVar *ci = 0;      CaloTrkVar *cp = GetCaloStoredTrack(tp->GetSeqNo());
682          ToFTrkVar  *oi = 0;      ToFTrkVar  *op = GetToFStoredTrack(tp->GetSeqNo());
683  //      cout << "trk track n. "<<i << " "<<hex<< tp <<dec<< endl;  
684          // if track has an image, check image selection      TrkTrack   *ti = 0; //tracker (image)
685          if(tp->HasImage()){      CaloTrkVar *ci = 0;
686        ToFTrkVar  *oi = 0;
687        //  cout << "trk track n. "<<i << " "<<hex<< tp <<dec<< endl;
688        // if track has an image, check image selection
689        if(tp->HasImage()){
690                            
691              ti = trk2_obj->TrkLevel2::GetTrackImage(i);              //tracker (image)        ti = trk2_obj->TrkLevel2::GetTrackImage(i);              //tracker (image)
692              ci = GetCaloStoredTrack(ti->GetSeqNo());        ci = GetCaloStoredTrack(ti->GetSeqNo());
693              oi = GetToFStoredTrack(ti->GetSeqNo());        oi = GetToFStoredTrack(ti->GetSeqNo());
694                            
695  //          cout << "its image "<<i << " "<<hex<< ti <<dec<< endl;        //            cout << "its image "<<i << " "<<hex<< ti <<dec<< endl;
696    
697              //assign starting scores        //assign starting scores
698              Int_t tp_score = 0;  //main track sorted by the tracker        Int_t tp_score = 0;  //main track sorted by the tracker
699              Int_t ti_score = 0;  //image track        Int_t ti_score = 0;  //image track
700                            
701              // ------------------------        // -----------------------------------------------------------------------------------------
702              // calorimeter check        // calorimeter check
703              // ------------------------        // -----------------------------------------------------------------------------------------
704              // check the Y spatial residual on the first calorimeter plane        // check the Y spatial residual on the first calorimeter plane
705              // (cut on calorimeter variables from Emiliano)        // (cut on calorimeter variables from Emiliano)
706              if( use_CAL && !calo2_obj ){        if( use_CAL && !calo2_obj ){
707                  cout << "void PamLevel2::SortTracks(TString how): how= "<<how<<" but CaloLevel2 not loaded !!!";          cout << "void PamLevel2::SortTracks(): howtosort= "<<how<<" but CaloLevel2 not loaded !!!";
708                  return;          return;
709              };        };
710              if(        if(
711                  use_CAL            &&           use_CAL            &&
712                  calo2_obj->npcfit[1] > 15     &&   //no. of fit planes on Y view  //       calo2_obj->npcfit[1] > 5     &&   //no. of fit planes on Y view
713                  calo2_obj->varcfit[1] < 1000. &&  //fit variance on Y view  //       calo2_obj->varcfit[1] < 1000. &&  //fit variance on Y view
714                  cp && ci &&           cp && ci &&
715                  true){           true){
716    
717                                    
718                  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;
719                  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;
720                                    
721                  if(resy_p <= resy_i) tp_score++;  //      if(resy_p <= resy_i) tp_score++;
722                  else                 ti_score++;  //      else                 ti_score++;
723            
724    
725              if( cp->npresh > ci->npresh) tp_score++;
726              if( cp->npresh < ci->npresh) ti_score++;
727              else ;//niente
728    
729  //              cout << "CALO "<<tp_score<<ti_score<<endl;          //      cout << "CALO "<<tp_score<<ti_score<<endl;
730    
731          };
732          // -----------------------------------------------------------------------------------------
733          // TOF check
734          // -----------------------------------------------------------------------------------------
735          // check the number of hit pmts along the track
736          // on S12 S21 and S32, where paddles are parallel to Y axis
737          if( use_TOF && !tof_obj ){
738            cout << "void PamLevel2::SortTracks(): howtosort= "<<how<<" but ToFLevel2 not loaded !!!";
739            return;
740          };
741          if( use_TOF && op && oi ){
742                    
743            //
744            Float_t sen = 0.;
745            for (Int_t ih=0; ih < op->npmtadc; ih++){
746              Int_t pl = tof_obj->GetPlaneIndex( (op->pmtadc).At(ih) );
747              if ( pl == 2 || pl == 3 || pl == 4 || pl == 5 ) sen += (op->dedx).At(ih);
748            };
749            for (Int_t ih=0; ih < oi->npmtadc; ih++){
750              Int_t pl = tof_obj->GetPlaneIndex( (oi->pmtadc).At(ih) );
751              if ( pl == 2 || pl == 3 || pl == 4 || pl == 5 ) sen += (oi->dedx).At(ih);
752            };
753            //
754            if (  sen >= sortthr ){
755              //printf(" IS A NUCLEUS! en = %f \n",sen);
756              //
757              // is a nucleus use a different algorithm
758              //
759              Int_t nz = 6; Float_t zin[6];                                          // << define TOF z-coordinates
760              for(Int_t ip=0; ip<nz; ip++)
761                zin[ip] = tof_obj->ToFLevel2::GetZTOF(tof_obj->ToFLevel2::GetToFPlaneID(ip));     // << read ToF plane z-coordinates
762              Trajectory *tr = new Trajectory(nz,zin);
763              //
764              Int_t nphit_p =0;
765              Int_t nphit_i =0;
766              Float_t enhit_p = 0.;
767              Float_t enhit_i = 0.;
768              //
769              for (Int_t ih=0; ih < op->npmtadc; ih++){
770                Int_t pl = tof_obj->GetPlaneIndex( (op->pmtadc).At(ih) );
771                if(pl == 1 || pl == 2 || pl == 5){
772                  nphit_p++;
773                  enhit_p += (op->dedx).At(ih);
774              };              };
775              // ------------------------            };
776              // TOF check            //
777              // ------------------------                tp->DoTrack2(tr);
778              // check the number of hit pmts along the track            //
779              // on S12 S21 and S32, where paddles are parallel to Y axis            if ( fabs(tr->y[0]-oi->ytofpos[0]) < 2. ){
780              if( use_TOF && !tof_obj ){              for (Int_t ih=0; ih < op->npmtadc; ih++){
781                  cout << "void PamLevel2::SortTracks(TString how): how= "<<how<<" but ToFLevel2 not loaded !!!";                Int_t pl = tof_obj->GetPlaneIndex( (op->pmtadc).At(ih) );
782                  return;                if(pl == 0){
783                    nphit_p++;
784                    enhit_p += (op->dedx).At(ih);
785                  };
786                };
787              };
788              if ( fabs(tr->y[3]-oi->ytofpos[1]) < 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 == 3){
792                    nphit_p++;
793                    enhit_p += (op->dedx).At(ih);
794                  };
795                };
796              };
797              if ( fabs(tr->y[4]-oi->ytofpos[2]) < 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 == 4){
801                    nphit_p++;
802                    enhit_p += (op->dedx).At(ih);
803                  };
804                };
805              };
806    
807              for (Int_t ih=0; ih < oi->npmtadc; ih++){
808                Int_t pl = tof_obj->GetPlaneIndex( (oi->pmtadc).At(ih) );
809                if(pl == 1 || pl == 2 || pl == 5){
810                  nphit_i++;        
811                  enhit_i += (op->dedx).At(ih);
812                };
813              };
814              //
815              ti->DoTrack2(tr);
816              //
817              if ( fabs(tr->y[0]-oi->ytofpos[0]) < 2. ){
818                for (Int_t ih=0; ih < oi->npmtadc; ih++){
819                  Int_t pl = tof_obj->GetPlaneIndex( (oi->pmtadc).At(ih) );
820                  if(pl == 0){
821                    nphit_i++;
822                    enhit_i += (op->dedx).At(ih);
823                  };
824                };
825              };
826              if ( fabs(tr->y[3]-oi->ytofpos[1]) < 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 == 3){
830                    nphit_i++;
831                    enhit_i += (op->dedx).At(ih);
832                  };
833                };
834              };
835              if ( fabs(tr->y[4]-oi->ytofpos[2]) < 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 == 4){
839                    nphit_i++;
840                    enhit_i += (op->dedx).At(ih);
841                  };
842                };
843              };            
844    
845                    
846              if(
847                 use_TOF            &&
848                 (nphit_p+nphit_i) !=0 &&  
849                 true){
850                        
851                //      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);
852                //      printf(" score p %i score i %i \n",tp_score,ti_score);
853                //            if( enhit_p > enhit_i ) tp_score++;
854                //            if( nphit_p >= nphit_i && enhit_p > enhit_i ) tp_score++;
855                if ( nphit_p > nphit_i ) tp_score++;
856                if ( nphit_p < nphit_i ) ti_score++;
857                if ( nphit_p == nphit_i ){
858                  if ( enhit_p > enhit_i ) tp_score++;
859                  else ti_score++;
860              };              };
861              if( use_TOF && op && oi ){              //      printf(" dopo score p %i score i %i \n",tp_score,ti_score);
862              };
863              delete tr;
864              //
865            } else {
866              //
867              // NOT a NUCLEUS
868              //
869              //printf(" NOT a NUCLEUS! en = %f \n",sen);
870    
871              Int_t nphit_p =0;
872              Int_t nphit_i =0;
873                                    
                 Int_t nphit_p =0;  
                 Int_t nphit_i =0;  
874                                    
875              /*                            cout << "track: npmtadc "<< op->npmtadc << endl;
876                cout << "track: npmttdc "<< op->npmttdc << endl;
877                cout << "image: npmtadc "<< oi->npmtadc << endl;
878                cout << "image: npmttdc "<< oi->npmttdc << endl;*/
879                                    
880  /*                              cout << "track: npmtadc "<< op->npmtadc << endl;  //        for (Int_t ih=0; ih < op->npmtadc; ih++){
881                                  cout << "track: npmttdc "<< op->npmttdc << endl;  //          Int_t pl = tof_obj->GetPlaneIndex( (op->pmtadc).At(ih) );
882                                  cout << "image: npmtadc "<< oi->npmtadc << endl;  //          if(pl == 1 || pl == 2 || pl == 5)nphit_p++;
883                                  cout << "image: npmttdc "<< oi->npmttdc << endl;*/  //        };
884                                    
885                  for (Int_t ih=0; ih < op->npmtadc; ih++){  //        for (Int_t ih=0; ih < oi->npmtadc; ih++){
886                      Int_t pl = tof_obj->GetPlaneIndex( (op->pmtadc).At(ih) );  //          Int_t pl = tof_obj->GetPlaneIndex( (oi->pmtadc).At(ih) );
887                      if(pl == 1 || pl == 2 || pl == 5)nphit_p++;  //          if(pl == 1 || pl == 2 || pl == 5)nphit_i++;
888                  };  //        };
889              // --- modified to count tdc signals (more efficient?)
890              // --- and to implement check on tdcflag
891              for (Int_t ih=0; ih < op->npmttdc; ih++){
892                Int_t pl = tof_obj->GetPlaneIndex( (op->pmttdc).At(ih) );
893    //          if( (op->tdcflag).At(ih)==0 && (pl == 1 || pl == 2 || pl == 5) )nphit_p++;
894                if( (op->tdcflag).At(ih)==0 )nphit_p++;
895              };
896                                    
897                  for (Int_t ih=0; ih < oi->npmtadc; ih++){            for (Int_t ih=0; ih < oi->npmttdc; ih++){
898                      Int_t pl = tof_obj->GetPlaneIndex( (oi->pmtadc).At(ih) );              Int_t pl = tof_obj->GetPlaneIndex( (oi->pmttdc).At(ih) );
899                      if(pl == 1 || pl == 2 || pl == 5)nphit_i++;  //          if( (oi->tdcflag).At(ih)==0 && (pl == 1 || pl == 2 || pl == 5) )nphit_i++;  
900                  };              if( (oi->tdcflag).At(ih)==0 )nphit_i++;    
901              };
902                                    
903                  if(            if(
904                      use_TOF            &&               (nphit_p+nphit_i) !=0 &&  
905                      (nphit_p+nphit_i) !=0 &&                   true){
                     true){  
906                                            
907                      if( nphit_p >= nphit_i) tp_score++;              if     ( nphit_p > nphit_i) tp_score++;
908                      else ti_score++;              else if( nphit_p < nphit_i) ti_score++;
909                  };              else ;//niente
910  //              cout << "TOF "<<tp_score<<ti_score<<endl;            };
911              };          };
912              if(tp_score == ti_score) use_TRK = true;          //      cout << "TOF "<<tp_score<<ti_score<<endl;
913              // ------------------------        };
914              // tracker check  
915              // ------------------------  
916              // chi**2 difference is not always large enough to distinguish among  //      if(tp_score == ti_score) use_TRK = true;
917              // the real track and its image.  
918              // Tracker check will be applied always when calorimeter and tof information is ambiguous.  
919              if(use_TRK){        // -----------------------------------------------------------------------------------------
920                  if(      tp->chi2 > 0 && tp->chi2 < ti->chi2 ) tp_score++ ;        // tracker check
921                  else if( ti->chi2 > 0 && ti->chi2 < tp->chi2 ) ti_score++ ;        // -----------------------------------------------------------------------------------------
922  //              cout << "TRK "<<tp_score<<ti_score<<endl;        // chi**2 difference is not always large enough to distinguish among
923              };        // the real track and its image.
924          // Tracker check will be applied always when calorimeter and tof information is ambiguous.
925          // *** MODIFIED ON AUGUST 2007 ***
926          if(use_TRK){
927    //      if(      tp->chi2 > 0 && tp->chi2 < ti->chi2 ) tp_score++ ;
928    //      else if( ti->chi2 > 0 && ti->chi2 < tp->chi2 ) ti_score++ ;
929    
930              // CHECK 1 : number of points along X  
931              if ( tp->GetNX() >= ti->GetNX() )tp_score++ ;
932              if ( tp->GetNX() <= ti->GetNX() )ti_score++ ;
933              // CHECK 2 : number of points along Y  
934              if ( tp->GetNY() >= ti->GetNY() )tp_score++ ;
935              if ( tp->GetNY() <= ti->GetNY() )ti_score++ ;
936              // CHECK 3 : chi**2 along X  
937    //        if(      tp->GetChi2X() > 0 && (tp->GetChi2X() < ti->GetChi2X() || ti->GetChi2X() <=0) ) tp_score++ ;
938    //        if(      ti->GetChi2X() > 0 && (ti->GetChi2X() < tp->GetChi2X() || tp->GetChi2X() <=0) ) ti_score++ ;
939              // CHECK 4 : chi**2 along Y  
940    //        if(      tp->GetChi2Y() > 0 && (tp->GetChi2Y() < ti->GetChi2Y() || ti->GetChi2Y() <=0) ) tp_score++ ;
941    //        if(      ti->GetChi2Y() > 0 && (ti->GetChi2Y() < tp->GetChi2Y() || tp->GetChi2Y() <=0) ) ti_score++ ;
942              // CHECK 5 : total chi**2  
943    //        if(      tp->chi2 > 0 && (tp->chi2 < ti->chi2 || ti->chi2 <=0) ) tp_score++ ;
944    //        if(      ti->chi2 > 0 && (ti->chi2 < tp->chi2 || tp->chi2 <=0) ) ti_score++ ;
945    
946            //      cout << "TRK "<<tp_score<<ti_score<<endl;
947          };
948                            
949              // ------------------------  
950              // the winner is....  
951              // ------------------------            // *-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
952              if      (tp_score > ti_score) {        // the winner is....
953  //              ts = tp;//the track sorted by the tracker!!        // *-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
954  //              cs = cp;        if      (tp_score > ti_score) {
955  //              os = op;          
956              }else if (tp_score < ti_score) {  
957                  ts = ti;//its image!!        }else if (tp_score < ti_score) {
958                  cs = ci;  
959                  os = oi;  
960            ts = ti;//its image!!
961                  ti = tp;//its image!!          cs = ci;
962                  ci = cp;          os = oi;
963                  oi = op;  
964            ti = tp;//its image!!
965                  tp = ts;//its image!!          ci = cp;
966                  cp = cs;          oi = op;
967                  op = os;  
968            tp = ts;//its image!!
969            cp = cs;
970            op = os;
971    
972                                    
973              }else {        }else {
974  //              ts = tp;  
975  //              cs = cp;  //      cout << "Warning - track image ambiguity not solved" << endl;
976  //              os = op;  
977  //                              cout << "Warning - track image ambiguity not solved" << endl;        };
 //                              cout << ts->GetNtot() << " "<< ts->chi2 << " " << npcfit[1] << " "<< nphit_p << endl;  
             };  
978                            
979          }else{      }else{
980  //          ts = tp;        //            ts = tp;
981  //          cs = cp;        //            cs = cp;
982  //          os = op;        //            os = op;
983          };      };
984                    
985  //      cout <<" SortTracks() "<<i<<" -- "<<ts<<endl;      //  cout <<" SortTracks() "<<i<<" -- "<<ts<<endl;
986  //      sorted_tracks->Add(ts);//save the track in the sorted array      //  sorted_tracks->Add(ts);//save the track in the sorted array
987  //      sorted_tracks.Add(ts);//save the track in the sorted array      //  sorted_tracks.Add(ts);//save the track in the sorted array
988  //      sorted_tracks.Add(tp);//save the track in the sorted array      //  sorted_tracks.Add(tp);//save the track in the sorted array
989  //      cout << "SortTracks:: sorted_tracks->Add(it) "<<i<<" "<<ts<<endl;      //  cout << "SortTracks:: sorted_tracks->Add(it) "<<i<<" "<<ts<<endl;
990  //      cout<<"o "<<tp<<endl;      //  cout<<"o "<<tp<<endl;
991  //      cout<<"o "<<cp<<endl;      //  cout<<"o "<<cp<<endl;
992  //      cout<<"o "<<op<<endl;      //  cout<<"o "<<op<<endl;
993          new(ttsorted[i]) PamTrack(tp,cp,op);      new(ttsorted[i]) PamTrack(tp,cp,op);
994          new(ttimage[i])  PamTrack(ti,ci,oi);      new(ttimage[i])  PamTrack(ti,ci,oi);
995      };    };
996    
997      if( tsorted->GetEntries() != trk2_obj->GetNTracks() ){    if( tsorted->GetEntries() != trk2_obj->GetNTracks() ){
998          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;
999          tsorted->Delete(); tsorted=0;      tsorted->Delete(); tsorted=0;
1000          timage->Delete(); timage=0;          timage->Delete(); timage=0;
1001      }    }
1002    
1003      //Restore Object count    //Restore Object count
1004      //To save space in the table keeping track of all referenced objects    //To save space in the table keeping track of all referenced objects
1005      //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.
1006      TProcessID::SetObjectCount(ObjectNumber);    TProcessID::SetObjectCount(ObjectNumber);
1007            
1008  };  };
1009  //--------------------------------------  //--------------------------------------
# Line 705  void PamLevel2::SortTracks(TString how){ Line 1027  void PamLevel2::SortTracks(TString how){
1027  TClonesArray *PamLevel2::GetTracks(){  TClonesArray *PamLevel2::GetTracks(){
1028    
1029  //  *-*-*-*-*-*-*-*-*-*-*-*-*  //  *-*-*-*-*-*-*-*-*-*-*-*-*
1030      SortTracks("+CAL+TOF");      SortTracks();
1031  //  *-*-*-*-*-*-*-*-*-*-*-*-*  //  *-*-*-*-*-*-*-*-*-*-*-*-*
1032    
1033      return tsorted;      return tsorted;
# Line 724  PamTrack *PamLevel2::GetTrack(int it){ Line 1046  PamTrack *PamLevel2::GetTrack(int it){
1046    
1047  //    if(!trk2_obj) return 0;  //    if(!trk2_obj) return 0;
1048                    
1049  // //  *-*-*-*-*-*-*-*-*-*-*-*-*  // //  *-*-*-*+-*-*-*-*-*-*-*-*-*
1050  //     SortTracks("+CAL+TOF");  //     SortTracks("+CAL+TOF");
1051  // //  *-*-*-*-*-*-*-*-*-*-*-*-*  // //  *-*-*-*-*-*-*-*-*-*-*-*-*
1052  // //    if(!sorted_tracks)return 0;  // //    if(!sorted_tracks)return 0;
# Line 743  PamTrack *PamLevel2::GetTrack(int it){ Line 1065  PamTrack *PamLevel2::GetTrack(int it){
1065    
1066  //    cout << "PamLevel2::GetTrack(int it) "<<endl;  //    cout << "PamLevel2::GetTrack(int it) "<<endl;
1067  //  *-*-*-*-*-*-*-*-*-*-*-*-*  //  *-*-*-*-*-*-*-*-*-*-*-*-*
1068      SortTracks("+CAL+TOF");      SortTracks();
1069  //  *-*-*-*-*-*-*-*-*-*-*-*-*  //  *-*-*-*-*-*-*-*-*-*-*-*-*
1070      if(!tsorted)return 0;      if(!tsorted)return 0;
1071      if(!tsorted->GetEntries())return 0;      if(!tsorted->GetEntries())return 0;
# Line 800  PamTrack *PamLevel2::GetTrackImage(int i Line 1122  PamTrack *PamLevel2::GetTrackImage(int i
1122    
1123    
1124  //  *-*-*-*-*-*-*-*-*-*-*-*-*  //  *-*-*-*-*-*-*-*-*-*-*-*-*
1125      SortTracks("+CAL+TOF");      SortTracks();
1126  //  *-*-*-*-*-*-*-*-*-*-*-*-*  //  *-*-*-*-*-*-*-*-*-*-*-*-*
1127      if(!timage)return 0;      if(!timage)return 0;
1128      if(!timage->GetEntries())return 0;      if(!timage->GetEntries())return 0;
# Line 837  TTree *PamLevel2::GetPamTree(TFile *f, T Line 1159  TTree *PamLevel2::GetPamTree(TFile *f, T
1159    
1160  //     if( !detlist.IsNull() ) SetWhichTrees(detlist);  //     if( !detlist.IsNull() ) SetWhichTrees(detlist);
1161  //     else                    GetWhichTrees(f);  //     else                    GetWhichTrees(f);
1162          if ( pam_tree ){
1163        printf("WARNING: TTree *PamLevel2::GetPamTree(TFile *fl, TString detlist) -- pam_tree already exists!\n ");
1164        return pam_tree;
1165      };
1166      //
1167    
1168      cout << "TTree *PamLevel2::GetPamTree(TFile *f, TString detlist ) -- obsolte "<<endl;      cout << "TTree *PamLevel2::GetPamTree(TFile *f, TString detlist ) -- obsolte "<<endl;
1169    
1170      SetWhichTrees(detlist);      SetWhichTrees(detlist);
# Line 913  TTree *PamLevel2::GetPamTree(TFile *f, T Line 1240  TTree *PamLevel2::GetPamTree(TFile *f, T
1240      if(R && TRG) {      if(R && TRG) {
1241          R->SetBranchAddress("TrigLevel2", GetPointerTo("TrigLevel2"));          R->SetBranchAddress("TrigLevel2", GetPointerTo("TrigLevel2"));
1242          cout << "Trigger      : set branch address TrigLevel2"<<endl;          cout << "Trigger      : set branch address TrigLevel2"<<endl;
1243          if(!Trout)Trout=O;          if(!Trout)Trout=R;
1244          else Trout->AddFriend(R);          else Trout->AddFriend(R);
1245      }else{      }else{
1246          cout << "Trigger      : missing tree"<<endl;          cout << "Trigger      : missing tree"<<endl;
# Line 923  TTree *PamLevel2::GetPamTree(TFile *f, T Line 1250  TTree *PamLevel2::GetPamTree(TFile *f, T
1250      if(S && S4) {      if(S && S4) {
1251          S->SetBranchAddress("S4Level2", GetPointerTo("S4Level2"));          S->SetBranchAddress("S4Level2", GetPointerTo("S4Level2"));
1252          cout << "S4           : set branch address S4Level2"<<endl;          cout << "S4           : set branch address S4Level2"<<endl;
1253          if(!Trout)Trout=O;          if(!Trout)Trout=S;
1254          else Trout->AddFriend(S);          else Trout->AddFriend(S);
1255      }else{      }else{
1256          cout << "S4           : missing tree"<<endl;          cout << "S4           : missing tree"<<endl;
# Line 933  TTree *PamLevel2::GetPamTree(TFile *f, T Line 1260  TTree *PamLevel2::GetPamTree(TFile *f, T
1260      if(N && ND) {      if(N && ND) {
1261          N->SetBranchAddress("NDLevel2", GetPointerTo("NDLevel2"));          N->SetBranchAddress("NDLevel2", GetPointerTo("NDLevel2"));
1262          cout << "NeutronD     : set branch address NDLevel2"<<endl;          cout << "NeutronD     : set branch address NDLevel2"<<endl;
1263          if(!Trout)Trout=O;          if(!Trout)Trout=N;
1264          else Trout->AddFriend(N);          else Trout->AddFriend(N);
1265      }else{      }else{
1266          cout << "NeutronD     : missing tree"<<endl;          cout << "NeutronD     : missing tree"<<endl;
# Line 943  TTree *PamLevel2::GetPamTree(TFile *f, T Line 1270  TTree *PamLevel2::GetPamTree(TFile *f, T
1270      if(A && AC) {      if(A && AC) {
1271          A->SetBranchAddress("AcLevel2", GetPointerTo("AcLevel2"));          A->SetBranchAddress("AcLevel2", GetPointerTo("AcLevel2"));
1272          cout << "Anticounter  : set branch address AcLevel2"<<endl;          cout << "Anticounter  : set branch address AcLevel2"<<endl;
1273          if(!Trout)Trout=O;          if(!Trout)Trout=A;
1274          else Trout->AddFriend(A);          else Trout->AddFriend(A);
1275      }else{      }else{
1276          cout << "Anticounter  : missing tree"<<endl;          cout << "Anticounter  : missing tree"<<endl;
# Line 953  TTree *PamLevel2::GetPamTree(TFile *f, T Line 1280  TTree *PamLevel2::GetPamTree(TFile *f, T
1280      if(B && ORB) {      if(B && ORB) {
1281          B->SetBranchAddress("OrbitalInfo", GetPointerTo("OrbitalInfo"));          B->SetBranchAddress("OrbitalInfo", GetPointerTo("OrbitalInfo"));
1282          cout << "OrbitalInfo  : set branch address OrbitalInfo"<<endl;          cout << "OrbitalInfo  : set branch address OrbitalInfo"<<endl;
1283          if(!Trout)Trout=O;          if(!Trout)Trout=B;
1284          else Trout->AddFriend(B);          else Trout->AddFriend(B);
1285      }else{      }else{
1286          cout << "OrbitalInfo  : missing tree"<<endl;          cout << "OrbitalInfo  : missing tree"<<endl;
# Line 995  TList*  PamLevel2::GetListOfLevel2Files( Line 1322  TList*  PamLevel2::GetListOfLevel2Files(
1322                    
1323      TString wdir = gSystem->WorkingDirectory();      TString wdir = gSystem->WorkingDirectory();
1324            
1325      if(ddir=="")ddir = wdir;      //    if(ddir=="")ddir = wdir;
1326  //      TSystemDirectory *datadir = new TSystemDirectory(gSystem->BaseName(ddir),ddir);  //      TSystemDirectory *datadir = new TSystemDirectory(gSystem->BaseName(ddir),ddir);
1327      cout << "Level2 data directory : "<<  ddir << endl;      if ( ddir != ""){
1328              cout << "Level2 data directory : "<<  ddir << endl;
1329        } else {
1330          cout << "Level2 data directory not given as input: trying to evaluate from list or taking working directory " << endl;
1331        };
1332      TList *contents  = new TList; // create output list      TList *contents  = new TList; // create output list
1333      contents->SetOwner();      contents->SetOwner();
1334            
# Line 1047  TList*  PamLevel2::GetListOfLevel2Files( Line 1377  TList*  PamLevel2::GetListOfLevel2Files(
1377  //          if( gSystem->IsFileInIncludePath(file,&fullpath) ){  //          if( gSystem->IsFileInIncludePath(file,&fullpath) ){
1378  //          if( (fullpath = gSystem->FindFile(ddir,file)) ){  //          if( (fullpath = gSystem->FindFile(ddir,file)) ){
1379              if( file.EndsWith(".root") ){              if( file.EndsWith(".root") ){
1380                  char *fullpath = gSystem->ConcatFileName(gSystem->DirName(ddir),gSystem->BaseName(file));                TString filedir;
1381                  if (ddir != ""){
1382                    filedir = ddir; // take the input dir
1383                  } else {
1384                    gSystem->ChangeDirectory(wdir); // back to the working directory
1385                    filedir = gSystem->DirName(file); // this will take the path if exist in the list otherwise it will return automatically the working dir
1386                  };
1387                    char *fullpath = gSystem->ConcatFileName(gSystem->DirName(filedir),gSystem->BaseName(file));
1388                  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
1389                  delete fullpath;                  delete fullpath;
1390              }              }
# Line 1057  TList*  PamLevel2::GetListOfLevel2Files( Line 1394  TList*  PamLevel2::GetListOfLevel2Files(
1394          };          };
1395          in.close();          in.close();
1396        } else {        } else {
1397              if(flisttxt.Contains("#"))flisttxt = flisttxt(0,flisttxt.First("#"));
1398            char *fullpath = gSystem->ConcatFileName(gSystem->DirName(flisttxt),gSystem->BaseName(flisttxt));            char *fullpath = gSystem->ConcatFileName(gSystem->DirName(flisttxt),gSystem->BaseName(flisttxt));
1399            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
1400            delete fullpath;            delete fullpath;
# Line 1066  TList*  PamLevel2::GetListOfLevel2Files( Line 1404  TList*  PamLevel2::GetListOfLevel2Files(
1404          cout << "No input file list given."<<endl;          cout << "No input file list given."<<endl;
1405          cout << "Check for existing root files."<<endl;          cout << "Check for existing root files."<<endl;
1406  //              cout << "Warking directory: "<< gSystem->WorkingDirectory()<< endl;  //              cout << "Warking directory: "<< gSystem->WorkingDirectory()<< endl;
1407                    if ( ddir == "" ){
1408              ddir = wdir;
1409              cout << "Level2 data directory : "<<  ddir << endl;
1410            };
1411    
1412          TSystemDirectory *datadir = new TSystemDirectory(gSystem->BaseName(ddir),ddir);          TSystemDirectory *datadir = new TSystemDirectory(gSystem->BaseName(ddir),ddir);
1413          TList *temp = datadir->GetListOfFiles();          TList *temp = datadir->GetListOfFiles();
1414  //              temp->Print();  //              temp->Print();
# Line 1110  TList*  PamLevel2::GetListOfLevel2Files( Line 1452  TList*  PamLevel2::GetListOfLevel2Files(
1452   * @return Pointer to a TChain   * @return Pointer to a TChain
1453   */   */
1454  TChain *PamLevel2::GetPamTree(TList *fl, TString detlist ){  TChain *PamLevel2::GetPamTree(TList *fl, TString detlist ){
1455      //
1456      //
1457      //
1458      if ( pam_tree ){
1459        printf("WARNING: TChain *PamLevel2::GetPamTree(TList *fl, TString detlist) -- pam_tree already exists!\n ");
1460        return pam_tree;
1461      };
1462      //
1463            
1464      TChain *Trout =0;      TChain *Trout =0;
1465    
# Line 1130  TChain *PamLevel2::GetPamTree(TList *fl, Line 1480  TChain *PamLevel2::GetPamTree(TList *fl,
1480      if(ORB)cout<<"ORB ";      if(ORB)cout<<"ORB ";
1481      cout << endl;      cout << endl;
1482    
1483      TChain *T = 0;             TChain *T = 0;    
1484      TChain *C = 0;       TChain *C = 0;
1485      TChain *O = 0;       TChain *O = 0;
1486      TChain *R = 0;       TChain *R = 0;
1487      TChain *S = 0;       TChain *S = 0;
1488      TChain *N = 0;       TChain *N = 0;
1489      TChain *A = 0;       TChain *A = 0;
1490      TChain *B = 0;       TChain *B = 0;
1491    
1492      TChain *L = 0;       TChain *L = 0;
1493            
1494      if(TRK2||TRK1||TRKh) T = new TChain("Tracker");          if(TRK2||TRK1||TRKh) T = new TChain("Tracker");    
1495      if(CAL2||CAL1)       C = new TChain("Calorimeter");      if(CAL2||CAL1)       C = new TChain("Calorimeter");
# Line 1149  TChain *PamLevel2::GetPamTree(TList *fl, Line 1499  TChain *PamLevel2::GetPamTree(TList *fl,
1499      if(ND)               N = new TChain("NeutronD");      if(ND)               N = new TChain("NeutronD");
1500      if(AC)               A = new TChain("Anticounter");      if(AC)               A = new TChain("Anticounter");
1501      if(ORB)              B = new TChain("OrbitalInfo");      if(ORB)              B = new TChain("OrbitalInfo");
   
1502       L = new TChain("SelectionList");       L = new TChain("SelectionList");
1503            
1504      // loop over files and create chains              // loop over files and create chains        
# Line 1169  TChain *PamLevel2::GetPamTree(TList *fl, Line 1518  TChain *PamLevel2::GetPamTree(TList *fl,
1518              if(ORB)              B->Add(name);              if(ORB)              B->Add(name);
1519              if(SELLI==1)         L->Add(name);              if(SELLI==1)         L->Add(name);
1520          };          };
1521      }      };
1522            
1523      cout << "done chain \n";      cout << "done chain \n";
1524        cout << "+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+" <<endl;
1525    
1526  //    UInt_t *found=0;  //    UInt_t *found=0;
1527  // Tracker  // Tracker
# Line 1186  TChain *PamLevel2::GetPamTree(TList *fl, Line 1536  TChain *PamLevel2::GetPamTree(TList *fl,
1536  //      else    T->SetBranchStatus("TrkHough",0,found);  //      else    T->SetBranchStatus("TrkHough",0,found);
1537          if(TRKh)cout << "Tracker      : set branch address TrkHough"<<endl;          if(TRKh)cout << "Tracker      : set branch address TrkHough"<<endl;
1538          if(!Trout)Trout=T;          if(!Trout)Trout=T;
1539          else Trout->AddFriend("Tracker");          else Trout->AddFriend("Tracker");
1540      }else{      }else{
1541          cout << "Tracker      : missing tree"<<endl;          cout << "Tracker      : missing tree"<<endl;
1542      };      };
# Line 1257  TChain *PamLevel2::GetPamTree(TList *fl, Line 1607  TChain *PamLevel2::GetPamTree(TList *fl,
1607      }else{      }else{
1608          cout << "OrbitalInfo  : missing tree"<<endl;          cout << "OrbitalInfo  : missing tree"<<endl;
1609      };      };
1610        cout << "+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+" <<endl;
1611    
1612      // Selection List      // Selection List
1613      if(L && SELLI==1) {      if(L && SELLI==1) {
# Line 1268  TChain *PamLevel2::GetPamTree(TList *fl, Line 1619  TChain *PamLevel2::GetPamTree(TList *fl,
1619          sel_tree = L;          sel_tree = L;
1620  //      if(!Trout)Trout=O;  //      if(!Trout)Trout=O;
1621  //      else Trout->AddFriend("SelectionList");  //      else Trout->AddFriend("SelectionList");
1622            cout << "----------------------------------------------------" <<endl;
1623      }else{      }else{
1624  //      cout << "SelectionList  : missing tree"<<endl;  //      cout << "SelectionList  : missing tree"<<endl;
1625          if(L)L->Delete();          if(L)L->Delete();
# Line 1459  void PamLevel2::SetBranchAddress(TChain Line 1811  void PamLevel2::SetBranchAddress(TChain
1811   * @return Pointer to a TChain   * @return Pointer to a TChain
1812   */   */
1813  TChain *PamLevel2::GetRunTree(TList *fl){  TChain *PamLevel2::GetRunTree(TList *fl){
1814              //
1815      TChain *R = new TChain("Run");    //
1816      //
1817      if ( run_tree ){
1818        printf("WARNING: TChain *PamLevel2::GetRunTree(TList *fl) -- run_tree already exists!\n ");
1819        return run_tree;
1820      };    
1821      //
1822      TChain *R = new TChain("Run");
1823            
1824      // loop over files and create chains              // loop over files and create chains        
1825      TIter next(fl);      TIter next(fl);
# Line 1499  TChain *PamLevel2::GetRunTree(TList *fl) Line 1858  TChain *PamLevel2::GetRunTree(TList *fl)
1858   * @return Pointer to a TTree   * @return Pointer to a TTree
1859   */   */
1860  TTree *PamLevel2::GetRunTree(TFile *f){  TTree *PamLevel2::GetRunTree(TFile *f){
1861      if ( run_tree ){
1862        printf("WARNING: TTree *PamLevel2::GetRunTree(TFile *f) -- run_tree already exists!\n ");
1863        return run_tree;
1864      };
1865            
1866    
1867      cout << "TTree *PamLevel2::GetRunTree(TFile *f) -- obsolte "<<endl;      cout << "TTree *PamLevel2::GetRunTree(TFile *f) -- obsolte "<<endl;
1868    
1869      TTree *R = (TTree*)f->Get("Run");      TTree *T = (TTree*)f->Get("Run");
1870            
1871      if(R){      if(T){
1872          R->SetBranchAddress("RunInfo", GetPointerTo("RunInfo"));          T->SetBranchAddress("RunInfo", GetPointerTo("RunInfo"));
1873          cout << "Run         : set branch address RunInfo"<<endl;          cout << "Run         : set branch address RunInfo"<<endl;
1874          R->SetBranchAddress("SoftInfo", GetPointerTo("SoftInfo")); // Emiliano          T->SetBranchAddress("SoftInfo", GetPointerTo("SoftInfo")); // Emiliano
1875          cout << "Software    : set branch address SoftInfo"<<endl; // Emiliano          cout << "Software    : set branch address SoftInfo"<<endl; // Emiliano
1876      }      }
1877    
1878      run_tree = (TChain*)R;      run_tree = (TChain*)T;
1879    
1880      return R;      return T;
1881            
1882  }  }
1883  /**  /**
# Line 1521  TTree *PamLevel2::GetRunTree(TFile *f){ Line 1885  TTree *PamLevel2::GetRunTree(TFile *f){
1885   * @param run Pointer to the chain/tree which contains run infos   * @param run Pointer to the chain/tree which contains run infos
1886   * @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
1887   */   */
1888  Bool_t PamLevel2::UpdateRunInfo(TChain *run, ULong64_t iev){  Bool_t PamLevel2::UpdateRunInfo(TChain *run, Long64_t iev){
1889    //    //
1890    // check if we have already called once GetEntry, if not call it    // check if we have already called once GetEntry, if not call it
1891    //    //
1892    if ( run->GetEntries() <= 0 ) return(false);      cout << "Bool_t PamLevel2::UpdateRunInfo(TChain *run, Long64_t iev) --- NON FUNZIONA BENISSIMO.... "<<endl;
1893        if(!run){
1894              cout << "Bool_t PamLevel2::UpdateRunInfo(TChain *run, ULong64_t iev) -- ERROR -- missing RunInfo tree "<<endl;
1895              return(false);        
1896        }
1897        if ( run->GetEntries() <= 0 ) return(false);
1898    //    //
1899        
1900    Int_t oldrun = irun;    //  Int_t oldrun = irun;
1901      Long64_t oldrun = irun;
1902    
1903    // --------------------------------------    // --------------------------------------
1904    // if it is a full file (not preselected)    // if it is a full file (not preselected)
1905    // --------------------------------------    // --------------------------------------
1906    if(SELLI==0){    if(SELLI==0){
       if ( irun < 0 ){  
           irun = 0;  
           run->GetEntry(irun);  
           runfirstentry = 0ULL;  
           runlastentry += (ULong64_t)(this->GetRunInfo()->NEVENTS);  
           if ( (ULong64_t)(this->GetRunInfo()->NEVENTS) > 0ULL ) runlastentry -= 1ULL;  
           return(true);  
       };        
1907    
1908        if( !GetOrbitalInfo() ){        //
1909          // the absolute time is necessary to relate the event with the run
1910          //
1911          if( !GetOrbitalInfo() && !ISGP){
1912            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;
1913            return(false);            return(false);
1914        }        }
1915    
1916          ULong64_t abstime = 0;
1917          if( GetOrbitalInfo() )abstime=GetOrbitalInfo()->absTime;
1918    
1919    
1920          //
1921          // the first time the routine is called, set run search from the beginning
1922          //
1923          if ( irun < 0LL ){
1924              irun = 0LL;
1925              run->GetEntry(irun);
1926              runfirstentry = 0LL;
1927              runlastentry += (Long64_t)(this->GetRunInfo()->NEVENTS);
1928              if ( (Long64_t)(this->GetRunInfo()->NEVENTS) > 0LL ) runlastentry -= 1LL;
1929              
1930              if( ISGP && run->GetEntries()!=1 ){
1931                  cout << "** WARNING ** simulated files are assumed to have 1 single run, not "<< run->GetEntries() << endl;
1932                  cout << "** WARNING ** run will not be updated"<<endl;
1933              }
1934    
1935          };      
1936          //
1937          if(ISGP)abstime = GetRunInfo()->RUNHEADER_TIME; // BARBATRUCCO
1938          //
1939          if ( irun == run->GetEntries()-1LL &&
1940               !(abstime >= GetRunInfo()->RUNHEADER_TIME &&
1941                 abstime <= GetRunInfo()->RUNTRAILER_TIME)
1942                ){
1943             irun = -1LL;
1944             runfirstentry = 0LL;
1945             runlastentry = -1LL;
1946           };
1947        // modificato il controllo sull'aggiornamento del run, per evitare problemi        // modificato il controllo sull'aggiornamento del run, per evitare problemi
1948        // dovuti agli eventi annidati (NB! NEVENTS conta anche questi!!)        // dovuti agli eventi annidati (NB! NEVENTS conta anche questi!!)
1949        while ( GetOrbitalInfo()->absTime > GetRunInfo()->RUNTRAILER_TIME && irun < run->GetEntries() ){        //
1950          bool fromfirst = true;
1951          //
1952          while ( !(abstime >= GetRunInfo()->RUNHEADER_TIME && abstime <= GetRunInfo()->RUNTRAILER_TIME) && irun < run->GetEntries()-1LL ){
1953  //      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)));  
1954            irun++;            irun++;
1955            run->GetEntry(irun);            run->GetEntry(irun);
1956            runfirstentry = runlastentry;            runfirstentry = runlastentry;
1957            if ( (ULong64_t)(this->GetRunInfo()->NEVENTS) > 0ULL ) runfirstentry += 1ULL;            if ( (Long64_t)(this->GetRunInfo()->NEVENTS) > 0LL ) runfirstentry += 1LL;
1958            runlastentry += (ULong64_t)(this->GetRunInfo()->NEVENTS);            runlastentry += (Long64_t)(this->GetRunInfo()->NEVENTS);
1959    //        cout << " ))))) UPDATE RUN INFO (((((  @iev "<<iev<<" run "<<GetRunInfo()->ID<<" irun "<<irun<<endl;
1960    //        cout << "runfirstentry "<<runfirstentry<<endl;
1961              //      printf(" iev %llu %u %llu \n",iev,this->GetRunInfo()->NEVENTS,(ULong64_t)(runfirstentry+(ULong64_t)(this->GetRunInfo()->NEVENTS)));
1962              //      printf(" abstime %u trailertime %u \n",abstime,GetRunInfo()->RUNTRAILER_TIME);
1963              //      printf(" IDRUN %u \n",GetRunInfo()->ID);
1964              //
1965    //        prevshift = 0;
1966              //
1967              if ( irun == (Long64_t)(run->GetEntries()-1LL) && fromfirst && !(abstime >= GetRunInfo()->RUNHEADER_TIME && abstime <= GetRunInfo()->RUNTRAILER_TIME)){
1968                  printf(" resetting irun  (it should NOT happen!!!)\n");
1969                  fromfirst = false;
1970                  irun = 0;
1971                  run->GetEntry(irun);
1972                  runfirstentry = 0ULL;
1973                  runlastentry += (Long64_t)(this->GetRunInfo()->NEVENTS);
1974                  if ( (Long64_t)(this->GetRunInfo()->NEVENTS) > 0LL ) runlastentry -= 1LL;
1975              };
1976              //
1977        };        };
1978                //
1979          if (
1980               !(abstime >= GetRunInfo()->RUNHEADER_TIME &&
1981                 abstime <= GetRunInfo()->RUNTRAILER_TIME)
1982              ) {
1983              printf(" Something very wrong here: cannot find RUN containing absolute time %u \n",abstime);
1984              return false;
1985          }
1986        //        //
1987        if ( irun == oldrun || irun >= run->GetEntries() ) return(false);        if ( irun == oldrun || irun >= run->GetEntries() ) return(false);
1988        //        //
1989        //  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);
1990        //        //
1991          prevshift = 0;
1992          cout << " ))))) UPDATE RUN INFO (((((  @iev "<<iev<<" run "<<GetRunInfo()->ID<<" irun "<<irun<<endl;
1993    //      cout << "runfirstentry "<<runfirstentry<<endl;
1994        return(true);            return(true);    
1995    };    };
1996    // ----------------------------------------------------    // ----------------------------------------------------
# Line 1571  Bool_t PamLevel2::UpdateRunInfo(TChain * Line 1999  Bool_t PamLevel2::UpdateRunInfo(TChain *
1999    // ----------------------------------------------------    // ----------------------------------------------------
2000    if(SELLI==1){          if(SELLI==1){      
2001        sel_tree->GetEntry(iev);        sel_tree->GetEntry(iev);
2002    //      cout << irun << " "<< irunentry << endl;
2003        if(irun != oldrun){        if(irun != oldrun){
2004            run->GetEntry(irun);            run->GetEntry(irun);
2005              cout << " ))))) UPDATE RUN INFO (((((  @iev "<<iev<<" run "<<GetRunInfo()->ID<<" irun "<<irun<<endl;
2006              prevshift = 0;
2007            return true;            return true;
2008        }        }
2009        return false;        return false;
# Line 1581  Bool_t PamLevel2::UpdateRunInfo(TChain * Line 2012  Bool_t PamLevel2::UpdateRunInfo(TChain *
2012    return false;    return false;
2013    //    //
2014  };  };
2015    
2016    Bool_t PamLevel2::UpdateRunInfo(Long64_t iev){
2017    
2018        if(!run_tree){
2019            cout << " Bool_t PamLevel2::UpdateRunInfo(ULong64_t iev) -- ERROR -- run tree not loaded"<<endl;
2020            return false;
2021        }
2022        if ( run_tree->GetEntries() <= 0 ) {
2023            cout << " Bool_t PamLevel2::UpdateRunInfo(ULong64_t iev) -- ERROR -- run tree is empty"<<endl;
2024            return(false);
2025        }
2026      
2027        Int_t oldrun = irun; // store current run index
2028    
2029        // +-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-
2030        // if it is a full file (not preselected)
2031        // +-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-
2032        if(SELLI==0){
2033    
2034            ULong64_t abstime = 0;
2035            ULong64_t obt     = 0;
2036            if( GetOrbitalInfo() ){
2037                abstime = GetOrbitalInfo()->absTime;
2038                obt     = GetOrbitalInfo()->OBT;
2039            }
2040            // ---------------------------------------------------------------
2041            // the absolute time is necessary to relate the event with the run
2042            // ---------------------------------------------------------------
2043            if( !GetOrbitalInfo() && !ISGP ){
2044                cout << "Bool_t PamLevel2::UpdateRunInfo(Long64_t iev) -- ERROR -- missing OrbitalInfo "<<endl;
2045                return(false);
2046            }
2047          
2048            // -----------------------------------------------------------------------
2049            // the first time the routine is called, set run search from the beginning
2050            // -----------------------------------------------------------------------
2051            if ( irun < 0 ){
2052                irun = 0LL;
2053                irunentry = 0;
2054                prevshift = 0;
2055                run_tree->GetEntry(irun);
2056    
2057                if( ISGP && run_tree->GetEntries()!=1 ){
2058                    cout << "** WARNING ** simulated files are assumed to have 1 single run, not "<< run_tree->GetEntries() << endl;
2059                    cout << "** WARNING ** run will not be updated"<<endl;
2060                }
2061            };      
2062    
2063            // -----------------------------------------------------------------------
2064            // if it is simulation, assigh abstime by end (temporaneo!!!)
2065            // -----------------------------------------------------------------------
2066            if(ISGP){
2067                abstime = GetRunInfo()->RUNHEADER_TIME; // BARBATRUCCO
2068                obt     = GetRunInfo()->RUNHEADER_OBT; // BARBATRUCCO
2069            }
2070            //
2071            bool fromfirst = true; // first loop over runs
2072    
2073    //      Bool_t hasfrag = false;
2074    //      if( GetRunInfo()->ID_RUN_FRAG!=0 && GetRunInfo()->ID_RUN_FRAG != GetRunInfo()->ID)hasfrag=true;
2075    //      ULong64_t fragid = GetRunInfo()->ID_RUN_FRAG;
2076    
2077            // ------------------------------------------------------
2078            // loop over runs to find the one that contains the event
2079            // ------------------------------------------------------
2080            while (
2081                (
2082                    (
2083                        !(abstime >= GetRunInfo()->RUNHEADER_TIME &&    // check on absolute time (s)
2084                          abstime <= GetRunInfo()->RUNTRAILER_TIME) &&
2085                        !(obt >= GetRunInfo()->RUNHEADER_OBT &&         // additional check on OBT (ms)
2086                          obt <= GetRunInfo()->RUNTRAILER_OBT)
2087                        )
2088                    || GetRunInfo()->NEVENTS==0
2089                    || !(irunentry < GetRunInfo()->NEVENTS-1-prevshift)
2090                    )
2091                    && irun < run_tree->GetEntries() ){
2092    
2093                irun++;
2094                // ------------------------------------
2095                // if the end of run tree is reached...
2096                // ------------------------------------
2097                if( irun == run_tree->GetEntries() ){
2098                    if(!fromfirst){
2099                        // -----------------------------------------------------
2100                        // if it happened already once and the run was not found
2101                        // ---> exit with error
2102                        // -----------------------------------------------------
2103                        cout << "Bool_t PamLevel2::UpdateRunInfo(Long64_t iev) -- ERROR -- event entry #"<<iev<<" does not belong to any run (should not happen)" <<endl;
2104                        return false;
2105                    }
2106                    // -----------------------------------------
2107                    // ...otherwise repeat search from beginning
2108                    // -----------------------------------------
2109                    cout << "Bool_t PamLevel2::UpdateRunInfo(Long64_t iev) -- WARNING -- reached end of run tree. searchin again from beginning " <<endl;
2110                    fromfirst = false;
2111                    irun = 0LL;
2112                    runfirstentry = 0LL;
2113                }
2114                // -------------------------------------------------------------------
2115                // save the index of the first entry of the run, relative to pam_tree,
2116                // and read a new run
2117                // -------------------------------------------------------------------
2118                if(irun>0)runfirstentry += (GetRunInfo()->NEVENTS)-prevshift;
2119                irunentry = 0;
2120                prevshift = 0;
2121                run_tree->GetEntry(irun);      
2122                if(GetRunInfo()->RUNHEADER_OBT>GetRunInfo()->RUNTRAILER_OBT ){
2123                    cout << "Bool_t PamLevel2::UpdateRunInfo(Long64_t iev) -- WARNING -- irun "<<irun<<"  has RUNHEADER_OBT>=RUNTRAILER_OBT " <<endl;
2124                    cout << "                                                            (NB!! in this case some events are assigned to a wrong run)"<<endl;
2125                }
2126    //          if(hasfrag &&  fragid != GetRunInfo()->ID){
2127    //              cout << "... where is the next fragment ??"<<endl;
2128    //          }
2129            };
2130    
2131    
2132            // --------------------------------------
2133            // if there was no need to update the run
2134            // ---> exit with FALSE
2135            // --------------------------------------
2136            if ( irun == oldrun ) return(false);
2137    
2138            // --------------------------------------
2139            // ... otherwise
2140            // ---> exit with TRUE
2141            // --------------------------------------
2142            cout << endl << " ))))) UPDATE RUN INFO (((((  @iev "<<iev<<" run "<<GetRunInfo()->ID<<" irun "<<irun<<endl;
2143            // ----------------------------------------------------
2144            // then check if the run has a fragment
2145            // in this case we have to switch to the next fragment
2146            // when the end of the first fragment is reached
2147            // ----------------------------------------------------
2148            if(
2149                GetRunInfo()->ID_RUN_FRAG != 0 &&
2150    //          GetRunInfo()->ID_RUN_FRAG != GetRunInfo()->ID &&
2151                true ){
2152                cout << "* fragment *"<<endl;              
2153            }
2154            
2155            return(true);    
2156        };
2157        // +-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-
2158        // if it is a preselected file (there is SelectionList)
2159        // +-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-
2160        if(SELLI==1){      
2161            sel_tree->GetEntry(iev);
2162            if(irun != oldrun){
2163                run_tree->GetEntry(irun);
2164                cout << endl << " ))))) UPDATE RUN INFO (((((  @iev "<<iev<<" run "<<GetRunInfo()->ID<<" irun "<<irun<<endl;
2165                prevshift = 0;
2166                return true;
2167            }
2168            return false;
2169        }
2170    
2171        return false;
2172        //
2173    };
2174  /**  /**
2175   * 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)
2176   * @param run Pointer to the chain/tree which contains run infos   * @param run Pointer to the chain/tree which contains run infos
2177   * @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
2178   */   */
2179  Bool_t PamLevel2::UpdateRunInfo(TTree *run, ULong64_t iev){  Bool_t PamLevel2::UpdateRunInfo(TTree *run, Long64_t iev){
2180    return(UpdateRunInfo((TChain*)run,iev));    return(UpdateRunInfo((TChain*)run,iev));
2181  };  };
2182    
# Line 2210  void PamLevel2::CreateCloneTrees0( TChai Line 2800  void PamLevel2::CreateCloneTrees0( TChai
2800   */   */
2801  void PamLevel2::CreateCloneTrees(TFile *ofile){  void PamLevel2::CreateCloneTrees(TFile *ofile){
2802    
2803    
2804    //  if the pointer is null, create a default file
2805    
2806        if(!ofile){
2807            cout << "void PamLevel2::CreateCloneTrees(TFile*) -- WARNING -- Creating file: clone-tree.root "<<endl;
2808            ofile = new TFile("clone-tree.root","recreate");
2809        }
2810    
2811      ofile->cd();      ofile->cd();
2812    
2813      cout << "+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+" <<endl;      cout << "+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+" <<endl;
# Line 2221  void PamLevel2::CreateCloneTrees(TFile * Line 2819  void PamLevel2::CreateCloneTrees(TFile *
2819      cout << "Run          : branch RunInfo"<<endl;      cout << "Run          : branch RunInfo"<<endl;
2820      run_tree_clone->Branch("SoftInfo","SoftInfo",GetPointerTo("SoftInfo"));      run_tree_clone->Branch("SoftInfo","SoftInfo",GetPointerTo("SoftInfo"));
2821      cout << "Run          : branch SoftInfo"<<endl;      cout << "Run          : branch SoftInfo"<<endl;
2822                    // ------------------
2823        // replicate run tree
2824        // ------------------
2825        cout << "----------------------------------------------------"<<endl;
2826        cout << "irun\t | RUN\t NEVENTS\t absolute time"<<endl;
2827        for (Int_t i=0; i<run_tree->GetEntries(); i++){
2828            run_tree->GetEntry(i);
2829            cout << i<< "\t | "<<GetRunInfo()->ID<<"\t "<<GetRunInfo()->NEVENTS<< "\t "<<GetRunInfo()->RUNHEADER_TIME<<" <---> "<<GetRunInfo()->RUNTRAILER_TIME<<endl;
2830            run_tree_clone->Fill();
2831        }
2832        cout << "----------------------------------------------------"<<endl;
2833    
2834    
2835      sel_tree_clone = new TTree("SelectionList","List of selected events ");      sel_tree_clone = new TTree("SelectionList","List of selected events ");
2836      sel_tree_clone->Branch("RunEntry",&irun,"runentry/I");      sel_tree_clone->Branch("RunEntry",&irun,"runentry/L");
2837      sel_tree_clone->Branch("EventEntry",&irunentry,"eventry/I");      sel_tree_clone->Branch("EventEntry",&irunentry,"eventry/L");
2838            
2839    
2840      Int_t i=0;      Int_t i=0;
# Line 2236  void PamLevel2::CreateCloneTrees(TFile * Line 2844  void PamLevel2::CreateCloneTrees(TFile *
2844              pam_tree_clone[i]->Branch("TrkLevel1","TrkLevel1", GetPointerTo("TrkLevel1"));              pam_tree_clone[i]->Branch("TrkLevel1","TrkLevel1", GetPointerTo("TrkLevel1"));
2845              pam_tree_clone[i]->BranchRef();              pam_tree_clone[i]->BranchRef();
2846              cout << "Tracker      : branch TrkLevel1"<<endl;              cout << "Tracker      : branch TrkLevel1"<<endl;
2847              cout << "CreateCloneTrees " << GetTrkLevel1()<<endl;  //          cout << "CreateCloneTrees " << GetTrkLevel1()<<endl;
2848          };          };
2849          if(TRK2) {          if(TRK2) {
2850              pam_tree_clone[i]->Branch("TrkLevel2", "TrkLevel2",GetPointerTo("TrkLevel2"));              pam_tree_clone[i]->Branch("TrkLevel2", "TrkLevel2",GetPointerTo("TrkLevel2"));
# Line 2366  void PamLevel2::WriteCloneTrees(){ Line 2974  void PamLevel2::WriteCloneTrees(){
2974  /**  /**
2975   * 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.
2976   */   */
2977  Int_t PamLevel2::GetEntry(Int_t iee){  //Int_t PamLevel2::GetEntry(Int_t iee){
2978    Int_t PamLevel2::GetEntry(Long64_t iee){
2979            
2980    //     cout << "-------------------------------------"<<endl;
2981    //     cout << "Int_t PamLevel2::GetEntry("<<iee<<")"<<endl;
2982    
2983      if(!pam_tree){      if(!pam_tree){
2984          cout << " Int_t PamLevel2::GetEntry(Int_t) -- ERROR -- level2 trees not loaded"<<endl;          cout << " Int_t PamLevel2::GetEntry(Int_t) -- ERROR -- level2 trees not loaded"<<endl;
2985          return 0;          return 0;
# Line 2382  Int_t PamLevel2::GetEntry(Int_t iee){ Line 2994  Int_t PamLevel2::GetEntry(Int_t iee){
2994      //  return 0;      //  return 0;
2995      //    }      //    }
2996    
2997      Int_t ii=0;      Long64_t ii=0;
2998      //-------------------------------      //-------------------------------
2999      ii = iee;      ii = iee;
3000      if( !pam_tree->GetEntry(ii) ) return 0;      if( !pam_tree->GetEntry(ii) ){      
3001            cout << " Int_t PamLevel2::GetEntry(Int_t) -- ERROR -- error reading pam tree"<<endl;
3002            return 0;
3003        }
3004      //      //
3005      // ... 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.      // ... 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.
3006      // 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      // 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
# Line 2394  Int_t PamLevel2::GetEntry(Int_t iee){ Line 3008  Int_t PamLevel2::GetEntry(Int_t iee){
3008      //      //
3009      if(!run_tree ){      if(!run_tree ){
3010          if ( TRK0 || CAL0 || TOF0 || RUN ) { //forse cosi` va bene per tornare 1?          if ( TRK0 || CAL0 || TOF0 || RUN ) { //forse cosi` va bene per tornare 1?
3011                  cout << " Int_t PamLevel2::GetEntry(Int_t) -- ERROR -- run tree not loaded"<<endl;              cout << " Int_t PamLevel2::GetEntry(Int_t) -- ERROR -- run tree not loaded"<<endl;
3012                  return 0;              return 0;
3013          }  else {          }  else {
3014                  return 1; //cosi` se non c'e` run esce qua...              return 1; //cosi` se non c'e` run esce qua...
3015          }          }
3016      }      }
3017    
3018      //-------------------------------      //-------------------------------
3019      ii = iee;      ii = iee;
3020      Bool_t UPDATED = UpdateRunInfo(run_tree,ii);  //    Bool_t UPDATED = UpdateRunInfo(run_tree,ii);
3021    //    Bool_t UPDATED = UpdateRunInfo(ii);
3022        UpdateRunInfo(ii);
3023      if(SELLI==0)irunentry = iee-runfirstentry;      if(SELLI==0)irunentry = iee-runfirstentry;
3024      if(UPDATED && run_tree_clone)run_tree_clone->Fill();  //    if(UPDATED && run_tree_clone)run_tree_clone->Fill();
3025      
3026  //    cout << "PamLevel2::GetEntry("<<iee<<") "<<irun<<" "<<runfirstentry<<" "<<irunentry<<endl;  //    cout << "PamLevel2::GetEntry("<<iee<<") "<<irun<<" "<<runfirstentry<<" "<<irunentry<<endl;
3027    
3028      if( TRK0 || CAL0 || TOF0 )GetYodaEntry( );  //     cout << "irunentry     "<<irunentry << endl;
3029    //     cout << "runfirstentry "<<runfirstentry << endl;
3030    //     cout << "nevents       "<<GetRunInfo()->NEVENTS<< endl;
3031    
3032    //     if( TRK0 || CAL0 || TOF0 ){
3033    //      if( !GetYodaEntry( ) ){
3034    //          cout << " Int_t PamLevel2::GetEntry(Int_t) -- ERROR -- error reading level0 tree"<<endl;
3035    //          return 0;
3036    //      }
3037    //     }
3038    
3039    
3040      return 1;      return 1;
3041    
3042  }  }
3043    
3044    TrkLevel0    *PamLevel2::GetTrkLevel0(){
3045        if( !TRK0 )return NULL;
3046        if( !GetYodaEntry( ) ){
3047            cout << " Int_t PamLevel2::GetTrkLevel0() -- ERROR -- error reading level0 tree"<<endl;
3048            return 0;
3049        }
3050        return trk0_obj;
3051    };
3052    CaloLevel0    *PamLevel2::GetCaloLevel0(){
3053        if( !CAL0 )return NULL;
3054        if( !GetYodaEntry( ) ){
3055            cout << " Int_t PamLevel2::GetCaloLevel0() -- ERROR -- error reading level0 tree"<<endl;
3056            return 0;
3057        }
3058        return calo0_obj;
3059    };
3060    
3061    
3062  /**  /**
3063   * 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.
3064   * 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 2429  TTree* PamLevel2::GetYodaTree( ){ Line 3073  TTree* PamLevel2::GetYodaTree( ){
3073      // check if iroot has changed      // check if iroot has changed
3074      //===================================      //===================================
3075      if( irun<0 ){      if( irun<0 ){
3076          cout << "PamLevel2::GetYodaTree() -- ERROR "<<endl;          cout << "PamLevel2::GetYodaTree() -- ERROR -- irun = "<<irun<<endl;
3077          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;  
3078          return NULL;          return NULL;
3079      }          }    
3080      Int_t irootnew = run_obj->ID_ROOT_L0;      Int_t irootnew = run_obj->ID_ROOT_L0;
3081  //     cout << "iroot    "<<iroot<<endl;  //      cout << "iroot    "<<iroot<<endl;
3082  //     cout << "irootnew "<<irootnew<<endl;  //      cout << "irootnew "<<irootnew<<endl;
3083    
3084      //===================================      //===================================
3085      // load the level0 file      // load the level0 file
# Line 2494  TTree* PamLevel2::GetYodaTree( ){ Line 3136  TTree* PamLevel2::GetYodaTree( ){
3136              };              };
3137              l0_tree->SetBranchAddress("Tracker" ,trk0_obj->GetPointerToTrackerEvent());              l0_tree->SetBranchAddress("Tracker" ,trk0_obj->GetPointerToTrackerEvent());
3138          }          }
3139          //---------------------------------------------------          //--------------------------------------------------
3140          // CALORIMETER:          // CALORIMETER:
3141          if(CAL0){          if(CAL0){
3142              if(!calo0_obj){              if(!calo0_obj){
# Line 2516  TTree* PamLevel2::GetYodaTree( ){ Line 3158  TTree* PamLevel2::GetYodaTree( ){
3158          cout << " TTree* PamLevel2::GetYodaTree( ) -- no DB connected... hai fatto qualche cazzata "<<endl;          cout << " TTree* PamLevel2::GetYodaTree( ) -- no DB connected... hai fatto qualche cazzata "<<endl;
3159      }      }
3160    
3161      if ( TRK0 ) TrkParams::SetCalib(run_obj,dbc);      if ( TRK0 ){
3162            TrkParams::Load(6);
3163            if( !TrkParams::IsLoaded(6) ){
3164                cout << " TTree* PamLevel2::GetYodaTree( ) -- WARNING -- VK-mask not loaded"<<endl;
3165            };
3166            TrkParams::SetCalib(run_obj,dbc);
3167            TrkParams::LoadCalib( );
3168            if( !TrkParams::CalibIsLoaded() ){
3169                cout << " TTree* PamLevel2::GetYodaTree( ) -- WARNING -- Calibration not loaded"<<endl;
3170            };
3171        }
3172    
3173  //    cout << l0_tree << endl;  //    cout << l0_tree << endl;
3174            
# Line 2534  Int_t PamLevel2::GetYodaEntry(){ Line 3185  Int_t PamLevel2::GetYodaEntry(){
3185      if(!GetYodaTree())return 0;      if(!GetYodaTree())return 0;
3186            
3187      // patch      // patch
3188      if( (Int_t)irunentry < 0){      if( irunentry < 0){
3189          cout << "Int_t PamLevel2::GetYodaEntry() -- ATTENZIONE -- irunentry negativo?!?! "<<(Int_t)irunentry<<endl;  //      cout << "Int_t PamLevel2::GetYodaEntry() -- ATTENZIONE -- irunentry negativo?!?! "<<(Int_t)irunentry<<endl;
3190          irunentry=0;          irunentry=0LL;
3191      }      }
3192      //  ---------------------------------      //  ---------------------------------
3193      //  if file is NOT a preselected file      //  if file is NOT a preselected file
3194      //  ---------------------------------      //  ---------------------------------
3195      UInt_t quellagiusta = irunentry + run_obj->EV_FROM;      Long64_t quellagiusta = irunentry + (Long64_t)(run_obj->EV_FROM);
3196        //     cout << " irunentry "<<irunentry<<endl;
3197    //     cout << " EV_FROM "<<run_obj->EV_FROM<<endl;
3198    //     cout << " quellagiusta = irunentry + EV_FROM "<< quellagiusta << endl;
3199    
3200  //     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;
3201  //     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;
3202  //     cout << " time "<< GetOrbitalInfo()->absTime << endl;  //     cout << " time "<< abstime << endl;
3203  //     cout << " trk_calib_used "<<run_obj->TRK_CALIB_USED<< endl;  //     cout << " trk_calib_used "<<run_obj->TRK_CALIB_USED<< endl;
3204            
3205      if( !GetOrbitalInfo() ){      ULong64_t obt     = 0;
3206        ULong64_t pktn    = 0;
3207        if( GetOrbitalInfo() ){
3208            obt     = GetOrbitalInfo()->OBT;
3209            pktn    = GetOrbitalInfo()->pkt_num;
3210        }
3211    
3212        if( !GetOrbitalInfo() && !ISGP){
3213          cout << "Int_t PamLevel2::GetYodaEntry() -- ERROR -- missing OrbitalInfo "<<endl;          cout << "Int_t PamLevel2::GetYodaEntry() -- ERROR -- missing OrbitalInfo "<<endl;
3214          return 0;          return 0;
3215      }      }
3216        if( obt==0 && pktn==0 && !ISGP ){
3217            cout << "Int_t PamLevel2::GetYodaEntry() -- ERROR -- level2 event corrupted ?? "<<endl;
3218            return 0;
3219        }
3220    
3221      // ---------------------------------------------------------------------      // ---------------------------------------------------------------------
3222      // ATTENTION!!!      // ATTENTION!!!
# Line 2567  Int_t PamLevel2::GetYodaEntry(){ Line 3231  Int_t PamLevel2::GetYodaEntry(){
3231      // ---------------------------------------------------------------------      // ---------------------------------------------------------------------
3232      Int_t answer = 0;      Int_t answer = 0;
3233      Int_t shift =0;      Int_t shift =0;
3234        //    printf(" siamo qui %i %i \n",shift,prevshift);
3235        Int_t maxshift = 10;
3236      do{      do{
3237          if(shift>0){              if(shift>0){    
3238            cout << " level0 <--> level2 mismatch ( irun "<<irun<<" irunentry "<<irunentry<<" shift "<<shift<<" )"<<endl;              cout << " PKTNUM  L2 --- "<< pktn << " --- L0 --- "<< GetEventHeader()->GetPscuHeader()->GetCounter()<<endl;
3239                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;
3240                cout << "         L2 <--> L0 mismatch ( irun "<<irun<<" irunentry "<<irunentry<<" shift "<<shift<<" prevshift "<<prevshift<<" )"<<endl;
3241          }          }
3242          answer = l0_tree->GetEntry(quellagiusta+shift+prevshift);          answer = l0_tree->GetEntry(quellagiusta+(Long64_t)shift+(Long64_t)prevshift);
3243          shift++;          shift++;
3244          if( !GetEventHeader() ){          if( !GetEventHeader() ){
3245              cout << "Int_t PamLevel2::GetYodaEntry() -- ERROR -- missing EventHeader "<<endl;              cout << "Int_t PamLevel2::GetYodaEntry() -- ERROR -- missing EventHeader "<<endl;
3246              return 0;              return 0;
3247          }          }
3248    
3249            if(ISGP){
3250                obt = (UInt_t)(GetEventHeader()->GetPscuHeader()->GetOrbitalTime()); //BARBATRUCCO
3251                pktn = (UInt_t)(GetEventHeader()->GetPscuHeader()->GetCounter()); //BARBATRUCCO
3252            }
3253    
3254    //      cout << "PKTNUM "<<shift<<" ==  L2 --- "<< pktn << " --- L0 --- "<< GetEventHeader()->GetPscuHeader()->GetCounter()<<endl;
3255    //      cout << " L2 --- "<< obt << " --- L0 --- "<< GetEventHeader()->GetPscuHeader()->GetOrbitalTime()<<endl;
3256  //      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;
3257  //      cout << " GetOrbitalInfo()->OBT "<< GetOrbitalInfo()->OBT << endl;  //      cout << " obt "<< obt << endl;
3258  //      cout << " GetEventHeader()->GetPscuHeader()->GetOrbitalTime() "<< GetEventHeader()->GetPscuHeader()->GetOrbitalTime() << endl;  //      cout << " GetEventHeader()->GetPscuHeader()->GetOrbitalTime() "<< GetEventHeader()->GetPscuHeader()->GetOrbitalTime() << endl;
3259  //      cout << " GetOrbitalInfo()->pkt_num "<< GetOrbitalInfo()->pkt_num << endl;  //      cout << " pktn "<< pktn << endl;
3260  //      cout << " GetEventHeader()->GetPscuHeader()->GetCounter() "<< GetEventHeader()->GetPscuHeader()->GetCounter() << endl;  //      cout << " GetEventHeader()->GetPscuHeader()->GetCounter() "<< GetEventHeader()->GetPscuHeader()->GetCounter() << endl;
3261          if ( prevshift != 0 && (quellagiusta+shift) == GetYodaTree()->GetEntries() ){  //      printf(" IDRUN %u \n",GetRunInfo()->ID);
3262    //
3263            if ( prevshift != 0 && (quellagiusta+(Long64_t)shift) == GetYodaTree()->GetEntries() ){
3264            prevshift = 0;            prevshift = 0;
3265            shift = -1;            shift = -1;
3266          };          };
3267    
3268      }while( ( GetOrbitalInfo()->OBT != (UInt_t)(GetEventHeader()->GetPscuHeader()->GetOrbitalTime()) || GetOrbitalInfo()->pkt_num != (UInt_t)(GetEventHeader()->GetPscuHeader()->GetCounter())) && (quellagiusta+shift) < GetYodaTree()->GetEntries() );      }while( ( obt != (UInt_t)(GetEventHeader()->GetPscuHeader()->GetOrbitalTime()) ||
3269                  pktn != (UInt_t)(GetEventHeader()->GetPscuHeader()->GetCounter())       )
3270      if ( (quellagiusta+shift) == GetYodaTree()->GetEntries() ) cout << " Big trouble here, no such event in Level0 data! " <<endl;              && (quellagiusta+(Long64_t)shift) < GetYodaTree()->GetEntries() && shift < maxshift);
3271    
3272        if ( (quellagiusta+(Long64_t)shift+(Long64_t)prevshift) > GetYodaTree()->GetEntries() || shift == maxshift ) {
3273            cout << " Big trouble here, no such event in Level0 data! " <<endl;
3274            return 0;
3275        }
3276  //    cout << "LA ENTRY GIUSTA E`: "<<quellagiusta<<" (spero...)"<<endl;  //    cout << "LA ENTRY GIUSTA E`: "<<quellagiusta<<" (spero...)"<<endl;
3277  //    return GetYodaTree()->GetEntry(quellagiusta);  //    return GetYodaTree()->GetEntry(quellagiusta);
3278      if ( shift > 1 ) prevshift += (shift-1);      if ( shift > 1 ) prevshift += (shift-1);
# Line 2598  Int_t PamLevel2::GetYodaEntry(){ Line 3280  Int_t PamLevel2::GetYodaEntry(){
3280      return answer;      return answer;
3281            
3282  }  }
3283    /**
3284     * \Brief Set DB connection
3285     */
3286    Bool_t PamLevel2::SetDBConnection(){
3287    
3288        cout<<"~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~"<<endl;
3289        cout<<"Connecting to DB"<<endl;
3290        cout<<"HOST "<<host<<endl;
3291        cout<<"USER "<<user<<endl;
3292        cout<<"PSW  "<<psw<<endl;
3293        dbc = TSQLServer::Connect(host.Data(),user.Data(),psw.Data());
3294        if( !dbc )return false;
3295        if( !dbc->IsConnected() )return false;    
3296        cout<<"~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~"<<endl;
3297        return true;
3298    
3299    }
Legend:
Removed from v.1.30  
changed lines
  Added in v.1.43
  ViewVC Help
Powered by ViewVC 1.1.23