/[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.54 by mocchiut, Fri Nov 9 09:20:25 2007 UTC revision 1.55 by pam-fi, Wed Nov 14 11:14:28 2007 UTC
# Line 367  void PamTrack::Delete(){ Line 367  void PamTrack::Delete(){
367   * Default Constructor   * Default Constructor
368   */   */
369  PamLevel2::PamLevel2(){  PamLevel2::PamLevel2(){
370    Initialize();      Initialize();
371  };  };
372    
373  /**  /**
# Line 385  PamLevel2::PamLevel2(){ Line 385  PamLevel2::PamLevel2(){
385   * +TRK2+CAL2+CAL1+TOF+TRG+ND+AC+S4+ORB   * +TRK2+CAL2+CAL1+TOF+TRG+ND+AC+S4+ORB
386   */   */
387  PamLevel2::PamLevel2(TString ddir,TString list,TString detlist){  PamLevel2::PamLevel2(TString ddir,TString list,TString detlist){
388    Initialize();      Initialize();
389    GetPamTree(GetListOfLevel2Files(ddir,list),detlist);      GetPamTree(GetListOfLevel2Files(ddir,list),detlist);
390    GetRunTree(GetListOfLevel2Files(ddir,list));      GetRunTree(GetListOfLevel2Files(ddir,list));
391  };  };
392    
393  /**  /**
# Line 398  PamLevel2::PamLevel2(TString ddir,TStrin Line 398  PamLevel2::PamLevel2(TString ddir,TStrin
398   * +TRK2+CAL2+CAL1+TOF+TRG+ND+AC+S4+ORB   * +TRK2+CAL2+CAL1+TOF+TRG+ND+AC+S4+ORB
399   */   */
400  PamLevel2::PamLevel2(TString ddir,TString list){  PamLevel2::PamLevel2(TString ddir,TString list){
401    Initialize();      Initialize();
402    GetPamTree(GetListOfLevel2Files(ddir,list),"");      GetPamTree(GetListOfLevel2Files(ddir,list),"");
403    GetRunTree(GetListOfLevel2Files(ddir,list));      GetRunTree(GetListOfLevel2Files(ddir,list));
404  };  };
405    
406    
# Line 538  void PamLevel2::Delete(){ Line 538  void PamLevel2::Delete(){
538      if(l0_file)l0_file->Close();      if(l0_file)l0_file->Close();
539      //    if(pam_tree)pam_tree->Delete();;      //    if(pam_tree)pam_tree->Delete();;
540    
541    if ( pam_tree ){      if ( pam_tree ){
542      //          //
543      // 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...          // 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...
544      //          //
545      TList *temp = pam_tree->GetListOfFriends();          TList *temp = pam_tree->GetListOfFriends();
546      TList *contents  = new TList; // create chain friend list          TList *contents  = new TList; // create chain friend list
547      contents->SetOwner();          contents->SetOwner();
548      TIter next(temp);          TIter next(temp);
549      TChain *questo = 0;          TChain *questo = 0;
550      while ( (questo = (TChain*) next()) ){          while ( (questo = (TChain*) next()) ){
551        TString name =  questo->GetName();              TString name =  questo->GetName();
552        contents->Add((TChain*)gROOT->FindObject(name.Data()));// add object to the list              contents->Add((TChain*)gROOT->FindObject(name.Data()));// add object to the list
553      };          };
554      //          //
555      // deleting the main chain          // deleting the main chain
556      //          //
557      pam_tree->Delete();          pam_tree->Delete();
558      //          //
559      // deleting the friends...          // deleting the friends...
560      //          //
561      TIter next2(contents);          TIter next2(contents);
562      TChain *questa = 0;          TChain *questa = 0;
563      while ( questa = (TChain*)next2() ){          while ( questa = (TChain*)next2() ){
564        TString name =  questa->GetName();              TString name =  questa->GetName();
565        questa->Delete();                    questa->Delete();      
566        questa=NULL;              questa=NULL;
567            };
568            //
569      };      };
570      //      pam_tree = NULL;
   };  
   pam_tree = NULL;  
571    
572      if(run_tree)run_tree->Delete();;      if(run_tree)run_tree->Delete();;
573      if(sel_tree)sel_tree->Delete();;      if(sel_tree)sel_tree->Delete();;
# Line 621  void PamLevel2::Clear(){ Line 621  void PamLevel2::Clear(){
621  };  };
622    
623  void PamLevel2::Reset(){  void PamLevel2::Reset(){
   //  
   // First of all clear everything  
   //  
   Clear();  
   //  
   // close and reset chains and pointers  
   //      
   if ( pam_tree ){  
624      //      //
625      // 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...      // First of all clear everything
626        //
627        Clear();
628        //
629        // close and reset chains and pointers
630        //    
631        if ( pam_tree ){
632            //
633            // 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...
634            //
635            TList *temp = pam_tree->GetListOfFriends();
636            TList *contents  = new TList; // create chain friend list
637            contents->SetOwner();
638            TIter next(temp);
639            TChain *questo = 0;
640            while ( (questo = (TChain*) next()) ){
641                TString name =  questo->GetName();
642                contents->Add((TChain*)gROOT->FindObject(name.Data()));// add object to the list
643            };
644            //
645            // deleting the main chain
646            //
647            pam_tree->Delete();
648            //
649            // deleting the friends...
650            //
651            TIter next2(contents);
652            TChain *questa = 0;
653            while ( questa = (TChain*)next2() ){
654                TString name =  questa->GetName();
655                questa->Delete();      
656                questa=NULL;
657            };
658            //
659        };
660        pam_tree = NULL;
661        //
662        if(run_tree)run_tree->Delete();;
663        run_tree = NULL;
664        if(sel_tree)sel_tree->Delete();;
665        sel_tree = NULL;
666        //
667        // Close file
668        //
669        if(l0_file)l0_file->Close("R");
670        l0_file = NULL;
671        //
672        h0_obj    = 0;
673        trk0_obj  = 0;
674        calo0_obj  = 0;
675        //
676        trk2_obj  = 0;
677        trk1_obj  = 0;
678        trkh_obj  = 0;
679        calo1_obj = 0;
680        calo2_obj = 0;
681        tof_obj   = 0;
682        trig_obj  = 0;
683        s4_obj    = 0;
684        nd_obj    = 0;
685        ac_obj    = 0;
686        orb_obj   = 0;
687        gp_obj    = 0;
688        //
689        // Reset run pointers
690        //
691        run_obj   = 0;//new GL_RUN();
692        soft_obj   = 0;// Emiliano
693        irun = -1;
694        runfirstentry = 0ULL;
695        runlastentry = 0ULL;    
696        //
697        totdltime[0] = 0LL;
698        totdltime[1] = 0LL;
699      //      //
     TList *temp = pam_tree->GetListOfFriends();  
     TList *contents  = new TList; // create chain friend list  
     contents->SetOwner();  
     TIter next(temp);  
     TChain *questo = 0;  
     while ( (questo = (TChain*) next()) ){  
       TString name =  questo->GetName();  
       contents->Add((TChain*)gROOT->FindObject(name.Data()));// add object to the list  
     };  
     //  
     // deleting the main chain  
     //  
     pam_tree->Delete();  
     //  
     // deleting the friends...  
     //  
     TIter next2(contents);  
     TChain *questa = 0;  
     while ( questa = (TChain*)next2() ){  
       TString name =  questa->GetName();  
       questa->Delete();        
       questa=NULL;  
     };  
     //  
   };  
   pam_tree = NULL;  
   //  
   if(run_tree)run_tree->Delete();;  
   run_tree = NULL;  
   if(sel_tree)sel_tree->Delete();;  
   sel_tree = NULL;  
   //  
   // Close file  
   //  
   if(l0_file)l0_file->Close("R");  
   l0_file = NULL;  
   //  
   h0_obj    = 0;  
   trk0_obj  = 0;  
   calo0_obj  = 0;  
   //  
   trk2_obj  = 0;  
   trk1_obj  = 0;  
   trkh_obj  = 0;  
   calo1_obj = 0;  
   calo2_obj = 0;  
   tof_obj   = 0;  
   trig_obj  = 0;  
   s4_obj    = 0;  
   nd_obj    = 0;  
   ac_obj    = 0;  
   orb_obj   = 0;  
   gp_obj    = 0;  
   //  
   // Reset run pointers  
   //  
   run_obj   = 0;//new GL_RUN();  
   soft_obj   = 0;// Emiliano  
   irun = -1;  
   runfirstentry = 0ULL;  
   runlastentry = 0ULL;      
   //  
   totdltime[0] = 0LL;  
   totdltime[1] = 0LL;  
   //  
700  };  };
701    
702  Bool_t PamLevel2::IsGood(){  Bool_t PamLevel2::IsGood(){
703    Bool_t goodev=true;      Bool_t goodev=true;
704    //  if(trk2_obj && trk2_obj->UnpackError() != 0 ) goodev = false;      //  if(trk2_obj && trk2_obj->UnpackError() != 0 ) goodev = false;
705    if(calo2_obj && calo2_obj->good != 1) goodev = false;      if(calo2_obj && calo2_obj->good != 1) goodev = false;
706    if(tof_obj && tof_obj->unpackError != 0) goodev = false;        if(tof_obj && tof_obj->unpackError != 0) goodev = false;  
707    if(trig_obj && trig_obj->unpackError != 0) goodev = false;      if(trig_obj && trig_obj->unpackError != 0) goodev = false;
708    if(s4_obj && s4_obj->unpackError != 0) goodev = false;        if(s4_obj && s4_obj->unpackError != 0) goodev = false;  
709    if(nd_obj && nd_obj->unpackError != 0) goodev = false;        if(nd_obj && nd_obj->unpackError != 0) goodev = false;  
710    if(ac_obj && ac_obj->unpackError != 255) goodev = false;        if(ac_obj && ac_obj->unpackError != 255) goodev = false;  
711      //  if(orb_obj)        //  if(orb_obj)  
712    return goodev;      return goodev;
713  };  };
714    
715  //--------------------------------------  //--------------------------------------
# Line 795  void *PamLevel2::GetPointerTo(const char Line 795  void *PamLevel2::GetPointerTo(const char
795   * Retrieves the calorimeter track matching the seqno-th tracker stored track.   * Retrieves the calorimeter track matching the seqno-th tracker stored track.
796   * (If seqno = -1 retrieves the self-trigger calorimeter track)   * (If seqno = -1 retrieves the self-trigger calorimeter track)
797   */   */
798   CaloTrkVar *PamLevel2::GetCaloStoredTrack(int seqno){  CaloTrkVar *PamLevel2::GetCaloStoredTrack(int seqno){
799            
800       if( !calo2_obj )return 0;      if( !calo2_obj )return 0;
801    
802       if( calo2_obj->CaloLevel2::ntrk()==0 ){      if( calo2_obj->CaloLevel2::ntrk()==0 ){
803           cout << "PamLevel2::GetCaloStoredTrack(int) : requested tracker SeqNo "<< seqno <<" but no Calorimeter tracks are stored"<<endl;          cout << "PamLevel2::GetCaloStoredTrack(int) : requested tracker SeqNo "<< seqno <<" but no Calorimeter tracks are stored"<<endl;
804           return NULL;          return NULL;
805       };      };
806            
807       CaloTrkVar *c = 0;      CaloTrkVar *c = 0;
808       Int_t it_calo=0;      Int_t it_calo=0;
809            
810       do{      do{
811           c = calo2_obj->CaloLevel2::GetCaloTrkVar(it_calo);          c = calo2_obj->CaloLevel2::GetCaloTrkVar(it_calo);
812           it_calo++;          it_calo++;
813       } while( c && seqno != c->trkseqno && it_calo < calo2_obj->CaloLevel2::ntrk());          } while( c && seqno != c->trkseqno && it_calo < calo2_obj->CaloLevel2::ntrk());    
814            
815       if(!c || seqno != c->trkseqno){      if(!c || seqno != c->trkseqno){
816           c = 0;          c = 0;
817           if(seqno!=-1)cout << "PamLevel2::GetCaloStoredTrack(int) : requested tracker SeqNo "<< seqno <<" does not match Calorimeter stored tracks"<<endl;          if(seqno!=-1)cout << "PamLevel2::GetCaloStoredTrack(int) : requested tracker SeqNo "<< seqno <<" does not match Calorimeter stored tracks"<<endl;
818       };      };
819       return c;      return c;
820            
821   };  };
822  //--------------------------------------  //--------------------------------------
823   //  //
824   //  //
825  //--------------------------------------  //--------------------------------------
826  /**  /**
827    * Retrieves the ToF track matching the seqno-th tracker stored track.   * Retrieves the ToF track matching the seqno-th tracker stored track.
828    * (If seqno = -1 retrieves the tracker-independent tof track)   * (If seqno = -1 retrieves the tracker-independent tof track)
829   */   */
830   ToFTrkVar *PamLevel2::GetToFStoredTrack(int seqno){  ToFTrkVar *PamLevel2::GetToFStoredTrack(int seqno){
831                    
832       if( !tof_obj )return 0;      if( !tof_obj )return 0;
833    
834       if( tof_obj->ToFLevel2::ntrk()==0 ){      if( tof_obj->ToFLevel2::ntrk()==0 ){
835           cout << "PamLevel2::GetToFStoredTrack(int) : requested tracker SeqNo "<< seqno <<" but no ToF tracks are stored"<<endl;          cout << "PamLevel2::GetToFStoredTrack(int) : requested tracker SeqNo "<< seqno <<" but no ToF tracks are stored"<<endl;
836           return NULL;          return NULL;
837       };      };
838            
839       ToFTrkVar *c = 0;      ToFTrkVar *c = 0;
840       Int_t it_tof=0;      Int_t it_tof=0;
841            
842       do{      do{
843           c = tof_obj->ToFLevel2::GetToFTrkVar(it_tof);          c = tof_obj->ToFLevel2::GetToFTrkVar(it_tof);
844           it_tof++;          it_tof++;
845       } while( c && seqno != c->trkseqno && it_tof < tof_obj->ToFLevel2::ntrk());              } while( c && seqno != c->trkseqno && it_tof < tof_obj->ToFLevel2::ntrk());
846            
847       if(!c || seqno != c->trkseqno){      if(!c || seqno != c->trkseqno){
848           c = 0;          c = 0;
849           if(seqno!=-1)cout << "PamLevel2::GetToFStoredTrack(int) : requested tracker SeqNo "<< seqno <<" does not match ToF stored tracks"<<endl;          if(seqno!=-1)cout << "PamLevel2::GetToFStoredTrack(int) : requested tracker SeqNo "<< seqno <<" does not match ToF stored tracks"<<endl;
850       };      };
851       return c;      return c;
852            
853   };  };
854    
855  //--------------------------------------  //--------------------------------------
856   //  //
857   //  //
858  //--------------------------------------  //--------------------------------------
859  /**  /**
860   * Give the pamela track associated to a tracker track, retrieving related calorimeter and tof track information.   * Give the pamela track associated to a tracker track, retrieving related calorimeter and tof track information.
861   */   */
862   PamTrack* PamLevel2::GetPamTrackAlong(TrkTrack* t){  PamTrack* PamLevel2::GetPamTrackAlong(TrkTrack* t){
863            
864       cout <<"PamLevel2::GetPamTrackAlong(TrkTrack* t) **obsolete** "<<endl;      cout <<"PamLevel2::GetPamTrackAlong(TrkTrack* t) **obsolete** "<<endl;
865       cout <<"(if you use it, remember to delete the PamTrack object)"<<endl;      cout <<"(if you use it, remember to delete the PamTrack object)"<<endl;
866    
867       CaloTrkVar *c = 0;      CaloTrkVar *c = 0;
868       ToFTrkVar  *o = 0;      ToFTrkVar  *o = 0;
869            
870       if(CAL2) c = GetCaloStoredTrack(t->GetSeqNo());      if(CAL2) c = GetCaloStoredTrack(t->GetSeqNo());
871       if(TOF) o = GetToFStoredTrack(t->GetSeqNo());      if(TOF) o = GetToFStoredTrack(t->GetSeqNo());
872            
873  //    if(t && c && o)track = new PamTrack(t,c,o);  //    if(t && c && o)track = new PamTrack(t,c,o);
874       PamTrack *track = new PamTrack(t,c,o);      PamTrack *track = new PamTrack(t,c,o);
875            
876       return track;      return track;
877    
878   };  };
879  //--------------------------------------  //--------------------------------------
880  //  //
881  //  //
# Line 911  PamTrack* PamLevel2::GetStoredTrack(Int_ Line 911  PamTrack* PamLevel2::GetStoredTrack(Int_
911  //  //
912    
913  /**  /**
914   * 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.   * Sort physical (tracker) tracks. Left here as backward compatibility method.
915   **/   **/
916  void PamLevel2::SortTracks(TString how){  void PamLevel2::SortTracks(TString how){
917    printf(" WARNING! obsolete, use SortTracks() and SetSortingMethod(TString) instead! \n Setting sorting method to %s \n",how.Data());      printf(" WARNING! obsolete, use SortTracks() and SetSortingMethod(TString) instead! \n Setting sorting method to %s \n",how.Data());
918    howtosort = how;        howtosort = how;  
919    SortTracks();      SortTracks();
920  };  };
921    
922  //  //
923  //--------------------------------------  //--------------------------------------
924  /**  /**
925   * Sort physical (tracker) tracks and stores them in the TRefArray (of TrkTrack objects) which pointer is  PamLevel2::sorted_tracks.   * Sort physical (tracker) tracks.
926   * @param how String to set the sorting cryterium (es: "CAL" or "TRK+CAL+TOF" ecc...).   * @param how String to set the sorting cryterium (es: "CAL" or "TRK+CAL+TOF" ecc...).
927   * Sorting cryteria:   * Sorting cryteria:
928   * TRK: lower chi**2   * TRK: lower chi**2
929   * CAL: lower Y spatial residual on the first calorimeter plane   * CAL: lower Y spatial residual on the first calorimeter plane
930   * TOF: bigger numebr of hit PMTs along the track, on S12 S21 S32 (where paddles are along the Y axis).   * TOF: bigger numebr of hit PMTs along the track, on S12 S21 S32 (where paddles are along the Y axis).
931     * S1: (ask Emiliano)
932     * S2: (ask Emiliano)
933     * S3: (ask Emiliano)
934     * GP: more GP hits
935   * The default sorting cryterium is "TOF+CAL".   * The default sorting cryterium is "TOF+CAL".
936   *   *
937   * 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).
938   */   */
939  void PamLevel2::SortTracks(){  void PamLevel2::SortTracks(){
940    TString how = howtosort;      TString how = howtosort;
941    
942        //  cout <<" PamLevel2::SortTracks(TString how) "<<endl;
943        if( !trk2_obj ){
944            cout << "void PamLevel2::SortTracks():  TrkLevel2 not loaded !!!";
945            return;
946        };
947        //Save current Object count
948        Int_t ObjectNumber = TProcessID::GetObjectCount();
949    
950    //  cout <<" PamLevel2::SortTracks(TString how) "<<endl;      // create TCloneArrays to store tracks and its images
951    if( !trk2_obj ){      if(!tsorted)tsorted = new TClonesArray("PamTrack",trk2_obj->GetNTracks());
952      cout << "void PamLevel2::SortTracks():  TrkLevel2 not loaded !!!";      tsorted->Delete();
953      return;      TClonesArray &ttsorted = *tsorted;
954    };  
955    //    cout << "call SortTracs() "<<endl;      if(!timage)timage = new TClonesArray("PamTrack",trk2_obj->GetNTracks());
956    //Save current Object count      timage->Delete();
957    Int_t ObjectNumber = TProcessID::GetObjectCount();      TClonesArray &ttimage = *timage;
   
   //    cout << "ObjectNumber  "<<ObjectNumber <<endl;  
   
   //     if(!sorted_tracks)sorted_tracks = new TRefArray();  
   //     sorted_tracks->Clear();  
   //    sorted_tracks.Clear();  
   
   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;  
958    
959    
960    
961    //--------------------------------------------------      //--------------------------------------------------
962    // retrieve sorting method      // retrieve sorting method
963    //--------------------------------------------------      //--------------------------------------------------
964    Bool_t use_TRK = how.Contains("TRK", TString::kIgnoreCase);      Bool_t use_TRK = how.Contains("TRK", TString::kIgnoreCase);
965    Bool_t use_CAL = how.Contains("CAL", TString::kIgnoreCase);      Bool_t use_CAL = how.Contains("CAL", TString::kIgnoreCase);
966    Bool_t use_TOF = how.Contains("TOF", TString::kIgnoreCase);      Bool_t use_TOF = how.Contains("TOF", TString::kIgnoreCase);
967    Bool_t use_S1 = how.Contains("S1", TString::kIgnoreCase);      Bool_t use_S1 = how.Contains("S1", TString::kIgnoreCase);
968    Bool_t use_S2 = how.Contains("S2", TString::kIgnoreCase);      Bool_t use_S2 = how.Contains("S2", TString::kIgnoreCase);
969    Bool_t use_S3 = how.Contains("S3", TString::kIgnoreCase);      Bool_t use_S3 = how.Contains("S3", TString::kIgnoreCase);
970        Bool_t use_GP = how.Contains("GP", TString::kIgnoreCase);
971        
972    if ( use_TOF ){      if ( use_TOF ){
973      use_S1 = true;          use_S1 = true;
974      use_S2 = true;          use_S2 = true;
975      use_S3 = true;          use_S3 = true;
976    };      };
977    if( !CAL2 &&  use_CAL) use_CAL = false;      if( !CAL2 &&  use_CAL) use_CAL = false;
978    if( !TOF ){      if( !TOF ){
979      use_TOF = false;          use_TOF = false;
980      use_S1  = false;          use_S1  = false;
981      use_S2  = false;          use_S2  = false;
982      use_S3  = false;          use_S3  = false;
983    }      }
984        if( !GP ){
985    if( !TRK2 ){          use_GP  = false;
986      //  cout << "SortTracks() : without tracker does not work!!! (not yet)" << endl;      }
987      return;  
988    };      if( !TRK2 ){
989            cout << "SortTracks() : without tracker does not work!!! (not yet)" << endl;
990            return;
991        };
992    
993    //   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;
994        
995    
996    //--------------------------------------------------      //--------------------------------------------------
997    // loop over "physical" tracks sorted by the tracker      // loop over "physical" tracks sorted by the tracker
998    //--------------------------------------------------      //--------------------------------------------------
999    for(Int_t i=0; i < trk2_obj->TrkLevel2::GetNTracks(); i++){      for(Int_t i=0; i < trk2_obj->TrkLevel2::GetNTracks(); i++){
1000                    
1001      TrkTrack   *ts = 0;          TrkTrack   *ts = 0;
1002      CaloTrkVar *cs = 0;          CaloTrkVar *cs = 0;
1003      ToFTrkVar  *os = 0;          ToFTrkVar  *os = 0;
1004                    
1005      // get tracker tracks          // get tracker tracks
1006      TrkTrack   *tp = trk2_obj->TrkLevel2::GetTrack(i); //tracker          TrkTrack   *tp = trk2_obj->TrkLevel2::GetTrack(i); //tracker
1007      CaloTrkVar *cp = GetCaloStoredTrack(tp->GetSeqNo());          CaloTrkVar *cp = GetCaloStoredTrack(tp->GetSeqNo());
1008      ToFTrkVar  *op = GetToFStoredTrack(tp->GetSeqNo());          ToFTrkVar  *op = GetToFStoredTrack(tp->GetSeqNo());
1009    
1010      TrkTrack   *ti = 0; //tracker (image)          TrkTrack   *ti = 0; //tracker (image)
1011      CaloTrkVar *ci = 0;          CaloTrkVar *ci = 0;
1012      ToFTrkVar  *oi = 0;          ToFTrkVar  *oi = 0;
1013      //  cout << "trk track n. "<<i << " "<<hex<< tp <<dec<< endl;          //      cout << "trk track n. "<<i << " "<<hex<< tp <<dec<< endl;
1014      // if track has an image, check image selection          // if track has an image, check image selection
1015    
1016      Int_t tp_score = 0;  //main track sorted by the tracker          Int_t tp_score = 0;  //main track sorted by the tracker
1017      Int_t ti_score = 0;  //image track          Int_t ti_score = 0;  //image track
1018      Int_t totp_score = 0;  //main track sorted by the tracker          Int_t totp_score = 0;  //main track sorted by the tracker
1019      Int_t toti_score = 0;  //image track          Int_t toti_score = 0;  //image track
1020    
1021      if(tp->HasImage()){          if(tp->HasImage()){
1022                            
1023        ti = trk2_obj->TrkLevel2::GetTrackImage(i);              //tracker (image)              ti = trk2_obj->TrkLevel2::GetTrackImage(i);              //tracker (image)
1024        ci = GetCaloStoredTrack(ti->GetSeqNo());              ci = GetCaloStoredTrack(ti->GetSeqNo());
1025        oi = GetToFStoredTrack(ti->GetSeqNo());              oi = GetToFStoredTrack(ti->GetSeqNo());
1026                            
1027        //            cout << "its image "<<i << " "<<hex<< ti <<dec<< endl;              //      cout << "its image "<<i << " "<<hex<< ti <<dec<< endl;
1028    
1029        //assign starting scores              //assign starting scores
1030        tp_score = 0;  //main track sorted by the tracker              tp_score = 0;  //main track sorted by the tracker
1031        ti_score = 0;  //image track              ti_score = 0;  //image track
1032                            
1033        // -----------------------------------------------------------------------------------------              // -----------------------------------------------------------------------------------------
1034        // calorimeter check              // *****************************************************************************************
1035        // -----------------------------------------------------------------------------------------              // -----------------------------------------------------------------------------------------
1036        // check the Y spatial residual on the first calorimeter plane              // calorimeter check
1037        // (cut on calorimeter variables from Emiliano)              // -----------------------------------------------------------------------------------------
1038        if( use_CAL && !calo2_obj ){              if( use_CAL && !calo2_obj ){
1039          cout << "void PamLevel2::SortTracks(): howtosort= "<<how<<" but CaloLevel2 not loaded !!!";                  cout << "void PamLevel2::SortTracks(): howtosort= "<<how<<" but CaloLevel2 not loaded !!!";
1040          return;                  return;
       };  
       if( use_CAL && !cp &&  ci ){  
         ti_score++;  
         toti_score++;  
       };  
       if( use_CAL &&  cp && !ci ){  
         tp_score++;  
         totp_score++;  
       };  
       if(  
          use_CAL            &&  
 //       calo2_obj->npcfit[1] > 5     &&   //no. of fit planes on Y view  
 //       calo2_obj->varcfit[1] < 1000. &&  //fit variance on Y view  
          cp && ci &&  
          true){  
   
                   
 //      Float_t resy_p = cp->tbar[0][1] - calo2_obj->cbar[0][1]; if(resy_p < 0)resy_p= - resy_p;  
 //      Float_t resy_i = ci->tbar[0][1] - calo2_obj->cbar[0][1]; if(resy_i < 0)resy_i= - resy_i;  
                   
 //      if(resy_p <= resy_i) tp_score++;  
 //      else                 ti_score++;  
           
   
         if( cp->npresh > ci->npresh){  
           tp_score++;  
           totp_score++;  
         };  
         if( cp->npresh < ci->npresh){  
           ti_score++;  
           toti_score++;  
         };  
   
         //      cout << "CALO "<<tp_score<<ti_score<<endl;  
   
       };  
       // -----------------------------------------------------------------------------------------  
       // TOF check  
       // -----------------------------------------------------------------------------------------  
       // check the number of hit pmts along the track  
       // on S12 S21 and S32, where paddles are parallel to Y axis  
       if( (use_TOF || use_S1 || use_S2 || use_S3 ) && !tof_obj ){  
         cout << "void PamLevel2::SortTracks(): howtosort= "<<how<<" but ToFLevel2 not loaded !!!";  
         return;  
       };  
       //  
       if( (use_TOF || use_S1 || use_S2 || use_S3 ) && !op && oi ){  
         ti_score++;  
         toti_score++;  
       };  
       if( (use_TOF || use_S1 || use_S2 || use_S3 ) && op && !oi ){  
         tp_score++;  
         totp_score++;  
       };  
       if( (use_TOF || use_S1 || use_S2 || use_S3 ) && op && oi ){  
         //  
         Float_t sen = 0.;  
         for (Int_t ih=0; ih < op->npmtadc; ih++){  
           Int_t pl = tof_obj->GetPlaneIndex( (op->pmtadc).At(ih) );  
           if ( pl == 2 || pl == 3 || pl == 4 || pl == 5 ) sen += (op->dedx).At(ih);  
         };  
         for (Int_t ih=0; ih < oi->npmtadc; ih++){  
           Int_t pl = tof_obj->GetPlaneIndex( (oi->pmtadc).At(ih) );  
           if ( pl == 2 || pl == 3 || pl == 4 || pl == 5 ) sen += (oi->dedx).At(ih);  
         };  
         //  
         if (  sen >= sortthr && false){ // temporary disabled NUCLEI special algorithm since the new one should work for every particle (to be checked!)  
           //printf(" IS A NUCLEUS! en = %f \n",sen);  
           //  
           // is a nucleus use a different algorithm  
           //  
           Int_t nz = 6; Float_t zin[6];                                          // << define TOF z-coordinates  
           for(Int_t ip=0; ip<nz; ip++)  
             zin[ip] = tof_obj->ToFLevel2::GetZTOF(tof_obj->ToFLevel2::GetToFPlaneID(ip));     // << read ToF plane z-coordinates  
           Trajectory *tr = new Trajectory(nz,zin);  
           //  
           Int_t nphit_p =0;  
           Int_t nphit_i =0;  
           Float_t enhit_p = 0.;  
           Float_t enhit_i = 0.;  
           //  
           for (Int_t ih=0; ih < op->npmtadc; ih++){  
             Int_t pl = tof_obj->GetPlaneIndex( (op->pmtadc).At(ih) );  
             if(pl == 1 || pl == 2 || pl == 5){  
               nphit_p++;  
               enhit_p += (op->dedx).At(ih);  
             };  
           };  
           //  
           tp->DoTrack2(tr);  
           //  
           if ( fabs(tr->y[0]-oi->ytofpos[0]) < 2. ){  
             for (Int_t ih=0; ih < op->npmtadc; ih++){  
               Int_t pl = tof_obj->GetPlaneIndex( (op->pmtadc).At(ih) );  
               if(pl == 0){  
                 nphit_p++;  
                 enhit_p += (op->dedx).At(ih);  
               };  
1041              };              };
1042            };              if( use_CAL && !cp &&  ci ){
1043            if ( fabs(tr->y[3]-oi->ytofpos[1]) < 2. ){                  ti_score++;
1044              for (Int_t ih=0; ih < op->npmtadc; ih++){                  toti_score++;
               Int_t pl = tof_obj->GetPlaneIndex( (op->pmtadc).At(ih) );  
               if(pl == 3){  
                 nphit_p++;  
                 enhit_p += (op->dedx).At(ih);  
               };  
1045              };              };
1046            };              if( use_CAL &&  cp && !ci ){
1047            if ( fabs(tr->y[4]-oi->ytofpos[2]) < 2. ){                  tp_score++;
1048              for (Int_t ih=0; ih < op->npmtadc; ih++){                  totp_score++;
               Int_t pl = tof_obj->GetPlaneIndex( (op->pmtadc).At(ih) );  
               if(pl == 4){  
                 nphit_p++;  
                 enhit_p += (op->dedx).At(ih);  
               };  
1049              };              };
1050            };              if(
1051                    use_CAL            &&
1052                    cp && ci &&
1053                    true){
1054    
1055                    
1056                    if(
1057                        cp->npresh > ci->npresh &&
1058                        true){
1059                        tp_score++;
1060                        totp_score++;
1061                    };
1062                    if(
1063                        cp->npresh < ci->npresh &&
1064                        true){
1065                        ti_score++;
1066                        toti_score++;
1067                    };
1068    
1069                    //      cout << "CALO "<<tp_score<<ti_score<<endl;
1070    
           for (Int_t ih=0; ih < oi->npmtadc; ih++){  
             Int_t pl = tof_obj->GetPlaneIndex( (oi->pmtadc).At(ih) );  
             if(pl == 1 || pl == 2 || pl == 5){  
               nphit_i++;          
               enhit_i += (op->dedx).At(ih);  
1071              };              };
1072            };              // -----------------------------------------------------------------------------------------
1073            //              // *****************************************************************************************
1074            ti->DoTrack2(tr);              // -----------------------------------------------------------------------------------------
1075            //              // TOF check
1076            if ( fabs(tr->y[0]-oi->ytofpos[0]) < 2. ){              // -----------------------------------------------------------------------------------------
1077              for (Int_t ih=0; ih < oi->npmtadc; ih++){              // check the number of hit pmts along the track
1078                Int_t pl = tof_obj->GetPlaneIndex( (oi->pmtadc).At(ih) );              // on S12 S21 and S32, where paddles are parallel to Y axis
1079                if(pl == 0){              if( (use_TOF || use_S1 || use_S2 || use_S3 ) && !tof_obj ){
1080                  nphit_i++;                  cout << "void PamLevel2::SortTracks(): howtosort= "<<how<<" but ToFLevel2 not loaded !!!";
1081                  enhit_i += (op->dedx).At(ih);                  return;
               };  
1082              };              };
1083            };              //
1084            if ( fabs(tr->y[3]-oi->ytofpos[1]) < 2. ){              if( (use_TOF || use_S1 || use_S2 || use_S3 ) && !op && oi ){
1085              for (Int_t ih=0; ih < oi->npmtadc; ih++){                  ti_score++;
1086                Int_t pl = tof_obj->GetPlaneIndex( (oi->pmtadc).At(ih) );                  toti_score++;
               if(pl == 3){  
                 nphit_i++;  
                 enhit_i += (op->dedx).At(ih);  
               };  
1087              };              };
1088            };              if( (use_TOF || use_S1 || use_S2 || use_S3 ) && op && !oi ){
1089            if ( fabs(tr->y[4]-oi->ytofpos[2]) < 2. ){                  tp_score++;
1090              for (Int_t ih=0; ih < oi->npmtadc; ih++){                  totp_score++;
               Int_t pl = tof_obj->GetPlaneIndex( (oi->pmtadc).At(ih) );  
               if(pl == 4){  
                 nphit_i++;  
                 enhit_i += (op->dedx).At(ih);  
               };  
1091              };              };
1092            };                          if( (use_TOF || use_S1 || use_S2 || use_S3 ) && op && oi ){
1093                    //
1094                    Float_t sen = 0.;
1095                    for (Int_t ih=0; ih < op->npmtadc; ih++){
1096                        Int_t pl = tof_obj->GetPlaneIndex( (op->pmtadc).At(ih) );
1097                        if ( pl == 2 || pl == 3 || pl == 4 || pl == 5 ) sen += (op->dedx).At(ih);
1098                    };
1099                    for (Int_t ih=0; ih < oi->npmtadc; ih++){
1100                        Int_t pl = tof_obj->GetPlaneIndex( (oi->pmtadc).At(ih) );
1101                        if ( pl == 2 || pl == 3 || pl == 4 || pl == 5 ) sen += (oi->dedx).At(ih);
1102                    };
1103                    //
1104                    if (  sen >= sortthr && false){ // temporary disabled NUCLEI special algorithm since the new one should work for every particle (to be checked!)
1105                        //printf(" IS A NUCLEUS! en = %f \n",sen);
1106                        //
1107                        // is a nucleus use a different algorithm
1108                        //
1109                        Int_t nz = 6; Float_t zin[6];                                          // << define TOF z-coordinates
1110                        for(Int_t ip=0; ip<nz; ip++)
1111                            zin[ip] = tof_obj->ToFLevel2::GetZTOF(tof_obj->ToFLevel2::GetToFPlaneID(ip));     // << read ToF plane z-coordinates
1112                        Trajectory *tr = new Trajectory(nz,zin);
1113                        //
1114                        Int_t nphit_p =0;
1115                        Int_t nphit_i =0;
1116                        Float_t enhit_p = 0.;
1117                        Float_t enhit_i = 0.;
1118                        //
1119                        for (Int_t ih=0; ih < op->npmtadc; ih++){
1120                            Int_t pl = tof_obj->GetPlaneIndex( (op->pmtadc).At(ih) );
1121                            if(pl == 1 || pl == 2 || pl == 5){
1122                                nphit_p++;
1123                                enhit_p += (op->dedx).At(ih);
1124                            };
1125                        };
1126                        //
1127                        tp->DoTrack2(tr);
1128                        //
1129                        if ( fabs(tr->y[0]-oi->ytofpos[0]) < 2. ){
1130                            for (Int_t ih=0; ih < op->npmtadc; ih++){
1131                                Int_t pl = tof_obj->GetPlaneIndex( (op->pmtadc).At(ih) );
1132                                if(pl == 0){
1133                                    nphit_p++;
1134                                    enhit_p += (op->dedx).At(ih);
1135                                };
1136                            };
1137                        };
1138                        if ( fabs(tr->y[3]-oi->ytofpos[1]) < 2. ){
1139                            for (Int_t ih=0; ih < op->npmtadc; ih++){
1140                                Int_t pl = tof_obj->GetPlaneIndex( (op->pmtadc).At(ih) );
1141                                if(pl == 3){
1142                                    nphit_p++;
1143                                    enhit_p += (op->dedx).At(ih);
1144                                };
1145                            };
1146                        };
1147                        if ( fabs(tr->y[4]-oi->ytofpos[2]) < 2. ){
1148                            for (Int_t ih=0; ih < op->npmtadc; ih++){
1149                                Int_t pl = tof_obj->GetPlaneIndex( (op->pmtadc).At(ih) );
1150                                if(pl == 4){
1151                                    nphit_p++;
1152                                    enhit_p += (op->dedx).At(ih);
1153                                };
1154                            };
1155                        };
1156    
1157                        for (Int_t ih=0; ih < oi->npmtadc; ih++){
1158                            Int_t pl = tof_obj->GetPlaneIndex( (oi->pmtadc).At(ih) );
1159                            if(pl == 1 || pl == 2 || pl == 5){
1160                                nphit_i++;  
1161                                enhit_i += (op->dedx).At(ih);
1162                            };
1163                        };
1164                        //
1165                        ti->DoTrack2(tr);
1166                        //
1167                        if ( fabs(tr->y[0]-oi->ytofpos[0]) < 2. ){
1168                            for (Int_t ih=0; ih < oi->npmtadc; ih++){
1169                                Int_t pl = tof_obj->GetPlaneIndex( (oi->pmtadc).At(ih) );
1170                                if(pl == 0){
1171                                    nphit_i++;
1172                                    enhit_i += (op->dedx).At(ih);
1173                                };
1174                            };
1175                        };
1176                        if ( fabs(tr->y[3]-oi->ytofpos[1]) < 2. ){
1177                            for (Int_t ih=0; ih < oi->npmtadc; ih++){
1178                                Int_t pl = tof_obj->GetPlaneIndex( (oi->pmtadc).At(ih) );
1179                                if(pl == 3){
1180                                    nphit_i++;
1181                                    enhit_i += (op->dedx).At(ih);
1182                                };
1183                            };
1184                        };
1185                        if ( fabs(tr->y[4]-oi->ytofpos[2]) < 2. ){
1186                            for (Int_t ih=0; ih < oi->npmtadc; ih++){
1187                                Int_t pl = tof_obj->GetPlaneIndex( (oi->pmtadc).At(ih) );
1188                                if(pl == 4){
1189                                    nphit_i++;
1190                                    enhit_i += (op->dedx).At(ih);
1191                                };
1192                            };
1193                        };          
1194    
1195                                    
1196            if(                      if(
1197               (use_TOF || use_S1 || use_S2 || use_S3)            &&                          (use_TOF || use_S1 || use_S2 || use_S3)            &&
1198               (nphit_p+nphit_i) !=0 &&                            (nphit_p+nphit_i) !=0 &&        
1199               true){                          true){
1200                                            
1201              //      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);                          //          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);
1202              //      printf(" score p %i score i %i \n",tp_score,ti_score);                          //          printf(" score p %i score i %i \n",tp_score,ti_score);
1203              //            if( enhit_p > enhit_i ) tp_score++;                          //                if( enhit_p > enhit_i ) tp_score++;
1204              //            if( nphit_p >= nphit_i && enhit_p > enhit_i ) tp_score++;                          //                if( nphit_p >= nphit_i && enhit_p > enhit_i ) tp_score++;
1205              if ( nphit_p > nphit_i ) tp_score++;                          if ( nphit_p > nphit_i ) tp_score++;
1206              if ( nphit_p < nphit_i ) ti_score++;                          if ( nphit_p < nphit_i ) ti_score++;
1207              if ( nphit_p == nphit_i ){                          if ( nphit_p == nphit_i ){
1208                if ( enhit_p > enhit_i ) tp_score++;                              if ( enhit_p > enhit_i ) tp_score++;
1209                else ti_score++;                              else ti_score++;
1210              };                          };
1211              //      printf(" dopo score p %i score i %i \n",tp_score,ti_score);                          //          printf(" dopo score p %i score i %i \n",tp_score,ti_score);
1212            };                      };
1213            delete tr;                      delete tr;
1214            //                      //
1215          } else {                  } else {
1216            //                      // -------------
1217            // NOT a NUCLEUS                      // NOT a NUCLEUS
1218            //                      // -------------
1219            //printf(" NOT a NUCLEUS! en = %f \n",sen);                      //printf(" NOT a NUCLEUS! en = %f \n",sen);
1220    
1221            Int_t nphit_p =0;                      Int_t nphit_p =0;
1222            Int_t nphit_i =0;                      Int_t nphit_i =0;
1223                                    
1224                                    
1225            /*                            cout << "track: npmtadc "<< op->npmtadc << endl;                      /*                          cout << "track: npmtadc "<< op->npmtadc << endl;
1226              cout << "track: npmttdc "<< op->npmttdc << endl;                                                  cout << "track: npmttdc "<< op->npmttdc << endl;
1227              cout << "image: npmtadc "<< oi->npmtadc << endl;                                                  cout << "image: npmtadc "<< oi->npmtadc << endl;
1228              cout << "image: npmttdc "<< oi->npmttdc << endl;*/                                                  cout << "image: npmttdc "<< oi->npmttdc << endl;*/
1229                                    
1230  //        for (Int_t ih=0; ih < op->npmtadc; ih++){  //        for (Int_t ih=0; ih < op->npmtadc; ih++){
1231  //          Int_t pl = tof_obj->GetPlaneIndex( (op->pmtadc).At(ih) );  //          Int_t pl = tof_obj->GetPlaneIndex( (op->pmtadc).At(ih) );
# Line 1236  void PamLevel2::SortTracks(){ Line 1236  void PamLevel2::SortTracks(){
1236  //          Int_t pl = tof_obj->GetPlaneIndex( (oi->pmtadc).At(ih) );  //          Int_t pl = tof_obj->GetPlaneIndex( (oi->pmtadc).At(ih) );
1237  //          if(pl == 1 || pl == 2 || pl == 5)nphit_i++;  //          if(pl == 1 || pl == 2 || pl == 5)nphit_i++;
1238  //        };  //        };
1239            // --- modified to count tdc signals (more efficient?)                      // --- modified to count tdc signals (more efficient?)
1240            // --- and to implement check on tdcflag                      // --- and to implement check on tdcflag
1241            for (Int_t ih=0; ih < op->npmttdc; ih++){                      for (Int_t ih=0; ih < op->npmttdc; ih++){
1242              Int_t pl = tof_obj->GetPlaneIndex( (op->pmttdc).At(ih) );                          Int_t pl = tof_obj->GetPlaneIndex( (op->pmttdc).At(ih) );
1243  //          if( (op->tdcflag).At(ih)==0 && (pl == 1 || pl == 2 || pl == 5) )nphit_p++;  //          if( (op->tdcflag).At(ih)==0 && (pl == 1 || pl == 2 || pl == 5) )nphit_p++;
1244              if ( (use_S1 && ( pl == 0 || pl == 1 )) || (use_S2 && ( pl == 2 || pl == 3 )) || (use_S3 && ( pl == 4 || pl == 5 )) ){                          if ( (use_S1 && ( pl == 0 || pl == 1 )) || (use_S2 && ( pl == 2 || pl == 3 )) || (use_S3 && ( pl == 4 || pl == 5 )) ){
1245                if( (op->tdcflag).At(ih)==0 )nphit_p++;                              if( (op->tdcflag).At(ih)==0 )nphit_p++;
1246              };                          };
1247            };                      };
1248                                    
1249            for (Int_t ih=0; ih < oi->npmttdc; ih++){                      for (Int_t ih=0; ih < oi->npmttdc; ih++){
1250              Int_t pl = tof_obj->GetPlaneIndex( (oi->pmttdc).At(ih) );                          Int_t pl = tof_obj->GetPlaneIndex( (oi->pmttdc).At(ih) );
1251  //          if( (oi->tdcflag).At(ih)==0 && (pl == 1 || pl == 2 || pl == 5) )nphit_i++;    //          if( (oi->tdcflag).At(ih)==0 && (pl == 1 || pl == 2 || pl == 5) )nphit_i++;  
1252              if ( (use_S1 && ( pl == 0 || pl == 1 )) || (use_S2 && ( pl == 2 || pl == 3 )) || (use_S3 && ( pl == 4 || pl == 5 )) ){                          if ( (use_S1 && ( pl == 0 || pl == 1 )) || (use_S2 && ( pl == 2 || pl == 3 )) || (use_S3 && ( pl == 4 || pl == 5 )) ){
1253                if( (oi->tdcflag).At(ih)==0 )nphit_i++;                                if( (oi->tdcflag).At(ih)==0 )nphit_i++;    
1254              };                          };
1255            };                      };
1256                                    
1257            if(                      if(
1258               (nphit_p+nphit_i) !=0 &&                            (nphit_p+nphit_i) !=0 &&        
1259               true){                          true){
1260                                            
1261              if ( nphit_p != nphit_i ){                          if ( nphit_p != nphit_i ){
1262                totp_score += nphit_p;                              totp_score += nphit_p;
1263                toti_score += nphit_i;                              toti_score += nphit_i;
1264                tp_score+=nphit_p;                              tp_score+=nphit_p;
1265                ti_score+=nphit_i;                              ti_score+=nphit_i;
1266                            };
1267                            //          if     ( nphit_p > nphit_i) tp_score+=nphit_p;
1268                            //          else if( nphit_p < nphit_i) ti_score+=nphit_i;
1269                            //          else ;//niente
1270                        };
1271                    };
1272                    //      cout << "TOF "<<tp_score<<ti_score<<endl;
1273              };              };
             //      if     ( nphit_p > nphit_i) tp_score+=nphit_p;  
             //      else if( nphit_p < nphit_i) ti_score+=nphit_i;  
             //      else ;//niente  
           };  
         };  
         //      cout << "TOF "<<tp_score<<ti_score<<endl;  
       };  
1274    
1275    
1276  //      if(tp_score == ti_score) use_TRK = true;  //      if(tp_score == ti_score) use_TRK = true;
1277    
1278    
1279        // -----------------------------------------------------------------------------------------              // -----------------------------------------------------------------------------------------
1280        // tracker check              // *****************************************************************************************
1281        // -----------------------------------------------------------------------------------------              // -----------------------------------------------------------------------------------------
1282        // chi**2 difference is not always large enough to distinguish among              // tracker check
1283        // the real track and its image.              // -----------------------------------------------------------------------------------------
1284        // Tracker check will be applied always when calorimeter and tof information is ambiguous.              // chi**2 difference is not always large enough to distinguish among
1285        // *** MODIFIED ON AUGUST 2007 ***              // the real track and its image.
1286        if(use_TRK){              // Tracker check will be applied always when calorimeter and tof information is ambiguous.
1287                // *** MODIFIED ON AUGUST 2007 ***
1288                if(use_TRK){
1289  //      if(      tp->chi2 > 0 && tp->chi2 < ti->chi2 ) tp_score++ ;  //      if(      tp->chi2 > 0 && tp->chi2 < ti->chi2 ) tp_score++ ;
1290  //      else if( ti->chi2 > 0 && ti->chi2 < tp->chi2 ) ti_score++ ;  //      else if( ti->chi2 > 0 && ti->chi2 < tp->chi2 ) ti_score++ ;
1291    
1292            // CHECK 1 : number of points along X                    // CHECK 1 : number of points along X  
1293          if ( tp->GetNX() >= ti->GetNX() ){                  if ( tp->GetNX() >= ti->GetNX() ){
1294            tp_score++ ;                      tp_score++ ;
1295            totp_score++ ;                      totp_score++ ;
1296          };                  };
1297          if ( tp->GetNX() <= ti->GetNX() ){                  if ( tp->GetNX() <= ti->GetNX() ){
1298            ti_score++ ;                      ti_score++ ;
1299            toti_score++ ;                      toti_score++ ;
1300          };                  };
1301            // CHECK 2 : number of points along Y                    // CHECK 2 : number of points along Y  
1302          if ( tp->GetNY() >= ti->GetNY() ){                  if ( tp->GetNY() >= ti->GetNY() ){
1303            tp_score++ ;                      tp_score++ ;
1304            totp_score++ ;                      totp_score++ ;
1305          };                  };
1306          if ( tp->GetNY() <= ti->GetNY() ){                  if ( tp->GetNY() <= ti->GetNY() ){
1307            ti_score++ ;                      ti_score++ ;
1308            toti_score++ ;                      toti_score++ ;
1309          };                  };
           // CHECK 3 : chi**2 along X    
 //        if(      tp->GetChi2X() > 0 && (tp->GetChi2X() < ti->GetChi2X() || ti->GetChi2X() <=0) ) tp_score++ ;  
 //        if(      ti->GetChi2X() > 0 && (ti->GetChi2X() < tp->GetChi2X() || tp->GetChi2X() <=0) ) ti_score++ ;  
           // CHECK 4 : chi**2 along Y    
 //        if(      tp->GetChi2Y() > 0 && (tp->GetChi2Y() < ti->GetChi2Y() || ti->GetChi2Y() <=0) ) tp_score++ ;  
 //        if(      ti->GetChi2Y() > 0 && (ti->GetChi2Y() < tp->GetChi2Y() || tp->GetChi2Y() <=0) ) ti_score++ ;  
           // CHECK 5 : total chi**2    
 //        if(      tp->chi2 > 0 && (tp->chi2 < ti->chi2 || ti->chi2 <=0) ) tp_score++ ;  
 //        if(      ti->chi2 > 0 && (ti->chi2 < tp->chi2 || tp->chi2 <=0) ) ti_score++ ;  
1310    
1311          //      cout << "TRK "<<tp_score<<ti_score<<endl;                  //      cout << "TRK "<<tp_score<<ti_score<<endl;
1312        };              };
1313                            
1314    
1315                // -----------------------------------------------------------------------------------------
1316                // *****************************************************************************************
1317                // -----------------------------------------------------------------------------------------
1318                // GPamela check
1319                // -----------------------------------------------------------------------------------------
1320                
1321                //---------------------------------------------------
1322                // count the number of GP hits
1323                //---------------------------------------------------
1324                if(use_GP){
1325                    int ngphits_p=0;
1326                    int ngphits_i=0;
1327                    float toll = 0.02; //200 micron
1328                    for(int ih=0; ih<GetGPamela()->Nthspe; ih++){
1329                        int ip = (Int_t) GetGPamela()->Itrpb[ih] - 1;
1330                        if(
1331                            tp &&
1332                            tp->YGood(ip) &&
1333                            fabs(tp->ym[ip]- GetGPamela()->Yavspe[ih])<toll &&
1334                            true) ngphits_p++;
1335                        if( ti &&
1336                            ti->YGood(ip) &&
1337                            fabs(ti->ym[ip]- GetGPamela()->Yavspe[ih])<toll &&
1338                            true) ngphits_i++;
1339                    }
1340                    if(
1341                        ngphits_p > ngphits_i &&
1342                        true){
1343                        tp_score++ ;
1344                        totp_score++ ;                  
1345                    }
1346                    if(
1347                        ngphits_p < ngphits_i &&
1348                        true){
1349                        ti_score++ ;
1350                        toti_score++ ;                  
1351                    }
1352                }
1353    
1354    
1355        // *-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*              // *-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
1356        // the winner is....              // the winner is....
1357        // *-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*              // *-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
1358        if      (tp_score > ti_score) {              if      (tp_score > ti_score) {
1359                    
1360    
1361        }else if (tp_score < ti_score) {              }else if (tp_score < ti_score) {
1362    
1363    
1364          ts = ti;//its image!!                  ts = ti;//its image!!
1365          cs = ci;                  cs = ci;
1366          os = oi;                  os = oi;
1367          Int_t totis = toti_score;                  Int_t totis = toti_score;
1368    
1369          ti = tp;//its image!!                  ti = tp;//its image!!
1370          ci = cp;                  ci = cp;
1371          oi = op;                  oi = op;
1372    
1373          tp = ts;//its image!!                  tp = ts;//its image!!
1374          cp = cs;                  cp = cs;
1375          op = os;                  op = os;
1376    
1377          toti_score = totp_score;                  toti_score = totp_score;
1378          totp_score = totis;                  totp_score = totis;
1379    
1380                                    
1381        }else {              }else {
1382    
1383  //      cout << "Warning - track image ambiguity not solved" << endl;  //      cout << "Warning - track image ambiguity not solved" << endl;
1384    
1385        };              };
1386                            
1387      }else{          }else{
1388        totp_score = 1;              totp_score = 1;
1389        toti_score = 0;              toti_score = 0;
1390                
1391        //            ts = tp;              //      ts = tp;
1392        //            cs = cp;              //      cs = cp;
1393        //            os = op;              //      os = op;
1394      };          };
1395                    
1396      //  cout <<" SortTracks() "<<i<<" -- "<<ts<<endl;          //      cout <<" SortTracks() "<<i<<" -- "<<ts<<endl;
1397      //  sorted_tracks->Add(ts);//save the track in the sorted array          //      sorted_tracks->Add(ts);//save the track in the sorted array
1398      //  sorted_tracks.Add(ts);//save the track in the sorted array          //      sorted_tracks.Add(ts);//save the track in the sorted array
1399      //  sorted_tracks.Add(tp);//save the track in the sorted array          //      sorted_tracks.Add(tp);//save the track in the sorted array
1400      //  cout << "SortTracks:: sorted_tracks->Add(it) "<<i<<" "<<ts<<endl;          //      cout << "SortTracks:: sorted_tracks->Add(it) "<<i<<" "<<ts<<endl;
1401      //  cout<<"o "<<tp<<endl;          //      cout<<"o "<<tp<<endl;
1402      //  cout<<"o "<<cp<<endl;          //      cout<<"o "<<cp<<endl;
1403      //  cout<<"o "<<op<<endl;          //      cout<<"o "<<op<<endl;
1404    
1405      new(ttsorted[i]) PamTrack(tp,cp,op);          new(ttsorted[i]) PamTrack(tp,cp,op);
1406      new(ttimage[i])  PamTrack(ti,ci,oi);          new(ttimage[i])  PamTrack(ti,ci,oi);
1407    
1408      ((PamTrack*)(ttsorted[i]))->SetPScore(totp_score);          ((PamTrack*)(ttsorted[i]))->SetPScore(totp_score);
1409      ((PamTrack*)(ttsorted[i]))->SetIScore(toti_score);          ((PamTrack*)(ttsorted[i]))->SetIScore(toti_score);
1410      ((PamTrack*)(ttimage[i]))->SetPScore(totp_score);          ((PamTrack*)(ttimage[i]))->SetPScore(totp_score);
1411      ((PamTrack*)(ttimage[i]))->SetIScore(toti_score);          ((PamTrack*)(ttimage[i]))->SetIScore(toti_score);
1412    };      };
1413    
1414    if( tsorted->GetEntries() != trk2_obj->GetNTracks() ){      if( tsorted->GetEntries() != trk2_obj->GetNTracks() ){
1415      cout << "void PamLevel2::SortTracks(): 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;
1416      tsorted->Delete(); tsorted=0;          tsorted->Delete(); tsorted=0;
1417      timage->Delete(); timage=0;          timage->Delete(); timage=0;    
1418    }      }
1419    
1420    //Restore Object count      //Restore Object count
1421    //To save space in the table keeping track of all referenced objects      //To save space in the table keeping track of all referenced objects
1422    //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.
1423    TProcessID::SetObjectCount(ObjectNumber);      TProcessID::SetObjectCount(ObjectNumber);
1424            
1425  };  };
1426  //--------------------------------------  //--------------------------------------
# Line 1417  TClonesArray *PamLevel2::GetTracks(){ Line 1449  TClonesArray *PamLevel2::GetTracks(){
1449    
1450      return tsorted;      return tsorted;
1451            
1452   };  };
1453  //--------------------------------------  //--------------------------------------
1454   //  //
1455   //  //
1456  //--------------------------------------  //--------------------------------------
1457  /**  /**
1458   * Retrieves the it-th Pamela "physical" track.   * Retrieves the it-th Pamela "physical" track.
# Line 1552  TTree *PamLevel2::GetPamTree(TFile *f, T Line 1584  TTree *PamLevel2::GetPamTree(TFile *f, T
1584          printf("WARNING: TTree *PamLevel2::GetPamTree(TFile *fl, TString detlist) -- pam_tree already exists!\n ");          printf("WARNING: TTree *PamLevel2::GetPamTree(TFile *fl, TString detlist) -- pam_tree already exists!\n ");
1585          return pam_tree;          return pam_tree;
1586      };      };
1587    //      //
1588    
1589      cout << "TTree *PamLevel2::GetPamTree(TFile *f, TString detlist ) -- obsolte "<<endl;      cout << "TTree *PamLevel2::GetPamTree(TFile *f, TString detlist ) -- obsolte "<<endl;
1590    
# Line 1704  TTree *PamLevel2::GetPamTree(TFile *f, T Line 1736  TTree *PamLevel2::GetPamTree(TFile *f, T
1736  //      cout << "SelectionList: set branch address EventEntry"<<endl;  //      cout << "SelectionList: set branch address EventEntry"<<endl;
1737  //      if(!Trout)Trout=O;  //      if(!Trout)Trout=O;
1738  //      else Trout->AddFriend("SelectionList");  //      else Trout->AddFriend("SelectionList");
1739          cout << " TTree *PamLevel2::GetPamTree(TFile, TString) >>> SelectionList not implemente!!!!"<<endl;          cout << " TTree *PamLevel2::GetPamTree(TFile, TString) >>> SelectionList not implemented!!!!"<<endl;
1740          sel_tree = 0;          sel_tree = 0;
1741      }else{      }else{
1742          cout << "SelectionList  : missing tree"<<endl;          cout << "SelectionList  : missing tree"<<endl;
# Line 1749  TList*  PamLevel2::GetListOfLevel2Files( Line 1781  TList*  PamLevel2::GetListOfLevel2Files(
1781      //    if(ddir=="")ddir = wdir;      //    if(ddir=="")ddir = wdir;
1782  //      TSystemDirectory *datadir = new TSystemDirectory(gSystem->BaseName(ddir),ddir);  //      TSystemDirectory *datadir = new TSystemDirectory(gSystem->BaseName(ddir),ddir);
1783      if ( ddir != ""){      if ( ddir != ""){
1784        cout << "Level2 data directory : "<<  ddir << endl;          cout << "Level2 data directory : "<<  ddir << endl;
1785      } else {      } else {
1786        cout << "Level2 data directory not given as input: trying to evaluate from list or taking working directory " << endl;          cout << "Level2 data directory not given as input: trying to evaluate from list or taking working directory " << endl;
1787      };      };
1788      TList *contents  = new TList; // create output list      TList *contents  = new TList; // create output list
1789      contents->SetOwner();      contents->SetOwner();
# Line 1768  TList*  PamLevel2::GetListOfLevel2Files( Line 1800  TList*  PamLevel2::GetListOfLevel2Files(
1800  //          return 0;  //          return 0;
1801  //      }        //      }      
1802  //      flisttxt = fullpath;  //      flisttxt = fullpath;
1803        if ( !flisttxt.EndsWith(".root") ){          if ( !flisttxt.EndsWith(".root") ){
1804    
1805          flisttxt = gSystem->ConcatFileName(gSystem->DirName(flisttxt),gSystem->BaseName(flisttxt));              flisttxt = gSystem->ConcatFileName(gSystem->DirName(flisttxt),gSystem->BaseName(flisttxt));
1806    
1807          if( !gSystem->ChangeDirectory(ddir) ){              if( !gSystem->ChangeDirectory(ddir) ){
1808              cout << "Cannot change directory : "<<ddir<<endl;                  cout << "Cannot change directory : "<<ddir<<endl;
1809              return 0;                  return 0;
         }  
           
         cout <<"Input file list : " << flisttxt <<endl;  
         ifstream in;  
         in.open(flisttxt, ios::in); //open input file list  
         if(!in.good()){  
             cout <<" ERROR opening the file "<<endl;  
             gSystem->ChangeDirectory(wdir); // back to the working directory  
             return 0;  
         }        
         int line=0;  
         while (1) {  
             TString file;  
             in >> file;  
             if (!in.good()) break;  
             line++;  
 //          cout <<"(1) " << file << endl;  
             if(file.IsNull()){  
                 cout << "-- list interrupted at line "<<line <<endl;  
                 break;  
1810              }              }
1811              if(file.Contains("#"))file = file(0,file.First("#"));          
1812                cout <<"Input file list : " << flisttxt <<endl;
1813                ifstream in;
1814                in.open(flisttxt, ios::in); //open input file list
1815                if(!in.good()){
1816                    cout <<" ERROR opening the file "<<endl;
1817                    gSystem->ChangeDirectory(wdir); // back to the working directory
1818                    return 0;
1819                }  
1820                int line=0;
1821                while (1) {
1822                    TString file;
1823                    in >> file;
1824                    if (!in.good()) break;
1825                    line++;
1826    //          cout <<"(1) " << file << endl;
1827                    if(file.IsNull()){
1828                        cout << "-- list interrupted at line "<<line <<endl;
1829                        break;
1830                    }
1831                    if(file.Contains("#"))file = file(0,file.First("#"));
1832  //          cout <<"(2) " << file << endl;  //          cout <<"(2) " << file << endl;
1833  //          if( gSystem->IsFileInIncludePath(file,&fullpath) ){  //          if( gSystem->IsFileInIncludePath(file,&fullpath) ){
1834  //          if( (fullpath = gSystem->FindFile(ddir,file)) ){  //          if( (fullpath = gSystem->FindFile(ddir,file)) ){
1835              if( file.EndsWith(".root") ){                  if( file.EndsWith(".root") ){
1836                TString filedir;                      TString filedir;
1837                if (ddir != ""){                      if (ddir != ""){
1838                  filedir = ddir; // take the input dir                          filedir = ddir; // take the input dir
1839                } else {                      } else {
1840                  gSystem->ChangeDirectory(wdir); // back to the working directory                          gSystem->ChangeDirectory(wdir); // back to the working directory
1841                  filedir = gSystem->DirName(file); // this will take the path if exist in the list otherwise it will return automatically the working dir                          filedir = gSystem->DirName(file); // this will take the path if exist in the list otherwise it will return automatically the working dir
1842                };                      };
1843                  char *fullpath = gSystem->ConcatFileName(gSystem->DirName(filedir),gSystem->BaseName(file));                      char *fullpath = gSystem->ConcatFileName(gSystem->DirName(filedir),gSystem->BaseName(file));
1844                  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
1845                  delete fullpath;                      delete fullpath;
1846              }                  }
1847  //          }else{  //          }else{
1848  //              if(file.Data()!="")cout << "File: "<<file<<" ---> missing "<< endl;  //              if(file.Data()!="")cout << "File: "<<file<<" ---> missing "<< endl;
1849  //          };  //          };
1850          };              };
1851          in.close();              in.close();
1852        } else {          } else {
1853            if(flisttxt.Contains("#"))flisttxt = flisttxt(0,flisttxt.First("#"));              if(flisttxt.Contains("#"))flisttxt = flisttxt(0,flisttxt.First("#"));
1854            char *fullpath = gSystem->ConcatFileName(gSystem->DirName(flisttxt),gSystem->BaseName(flisttxt));              char *fullpath = gSystem->ConcatFileName(gSystem->DirName(flisttxt),gSystem->BaseName(flisttxt));
1855            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
1856            delete fullpath;              delete fullpath;
1857        };                  };      
1858      }else{      }else{
1859                    
1860          cout << "No input file list given."<<endl;          cout << "No input file list given."<<endl;
1861          cout << "Check for existing root files."<<endl;          cout << "Check for existing root files."<<endl;
1862  //              cout << "Warking directory: "<< gSystem->WorkingDirectory()<< endl;  //              cout << "Warking directory: "<< gSystem->WorkingDirectory()<< endl;
1863          if ( ddir == "" ){          if ( ddir == "" ){
1864            ddir = wdir;              ddir = wdir;
1865            cout << "Level2 data directory : "<<  ddir << endl;              cout << "Level2 data directory : "<<  ddir << endl;
1866          };          };
1867    
1868          TSystemDirectory *datadir = new TSystemDirectory(gSystem->BaseName(ddir),ddir);          TSystemDirectory *datadir = new TSystemDirectory(gSystem->BaseName(ddir),ddir);
# Line 1902  TChain *PamLevel2::GetPamTree(TList *fl, Line 1934  TChain *PamLevel2::GetPamTree(TList *fl,
1934    
1935  //    if( !detlist.IsNull() )SetWhichTrees(detlist);  //    if( !detlist.IsNull() )SetWhichTrees(detlist);
1936  //      //    
1937       TChain *T = 0;          TChain *T = 0;      
1938       TChain *C = 0;      TChain *C = 0;
1939       TChain *O = 0;      TChain *O = 0;
1940       TChain *R = 0;      TChain *R = 0;
1941       TChain *S = 0;      TChain *S = 0;
1942       TChain *N = 0;      TChain *N = 0;
1943       TChain *A = 0;      TChain *A = 0;
1944       TChain *B = 0;      TChain *B = 0;
1945       TChain *G = 0;      TChain *G = 0;
1946    
1947       TChain *L = 0;      TChain *L = 0;
1948            
1949      if(TRK2||TRK1||TRKh) T = new TChain("Tracker");          if(TRK2||TRK1||TRKh) T = new TChain("Tracker");    
1950      if(CAL2||CAL1)       C = new TChain("Calorimeter");      if(CAL2||CAL1)       C = new TChain("Calorimeter");
# Line 2055  TChain *PamLevel2::GetPamTree(TList *fl, Line 2087  TChain *PamLevel2::GetPamTree(TList *fl,
2087      // Selection List      // Selection List
2088      if(L && SELLI==1) {      if(L && SELLI==1) {
2089          cout<<">>> Found selection-list <<<"<<endl;          cout<<">>> Found selection-list <<<"<<endl;
2090    //      L->SetBranchAddress("RunEntry",&irun);
2091          L->SetBranchAddress("RunEntry",&irun);          L->SetBranchAddress("RunEntry",&irun);
2092          cout << "SelectionList: set branch address RunEntry"<<endl;          cout << "SelectionList: set branch address RunEntry"<<endl;
2093          L->SetBranchAddress("EventEntry",&irunentry);          L->SetBranchAddress("EventEntry",&irunentry);
# Line 2179  void PamLevel2::SetBranchAddress(TTree * Line 2212  void PamLevel2::SetBranchAddress(TTree *
2212      };      };
2213      // SelectionList      // SelectionList
2214      if(SELLI==1) {      if(SELLI==1) {
2215    //      t->SetBranchAddress("RunEntry", &irun);
2216          t->SetBranchAddress("RunEntry", &irun);          t->SetBranchAddress("RunEntry", &irun);
2217          cout << "SelectionList: set branch address RunEntry"<<endl;          cout << "SelectionList: set branch address RunEntry"<<endl;
2218          t->SetBranchAddress("EventEntry", &irunentry);          t->SetBranchAddress("EventEntry", &irunentry);
# Line 2205  void PamLevel2::SetBranchAddress(TChain Line 2239  void PamLevel2::SetBranchAddress(TChain
2239  //    GP     = GP & t->GetBranchStatus("h20");  //    GP     = GP & t->GetBranchStatus("h20");
2240    
2241      // Tracker      // Tracker
2242       if(TRK2) {      if(TRK2) {
2243          t->SetBranchAddress("TrkLevel2", GetPointerTo("TrkLevel2"));          t->SetBranchAddress("TrkLevel2", GetPointerTo("TrkLevel2"));
2244          cout << "Tracker      : set branch address TrkLevel2"<<endl;          cout << "Tracker      : set branch address TrkLevel2"<<endl;
2245      };      };
# Line 2213  void PamLevel2::SetBranchAddress(TChain Line 2247  void PamLevel2::SetBranchAddress(TChain
2247          t->SetBranchAddress("TrkLevel1", GetPointerTo("TrkLevel1"));          t->SetBranchAddress("TrkLevel1", GetPointerTo("TrkLevel1"));
2248          cout << "Tracker      : set branch address TrkLevel1"<<endl;          cout << "Tracker      : set branch address TrkLevel1"<<endl;
2249      };      };
2250       if(TRKh) {      if(TRKh) {
2251          t->SetBranchAddress("TrkHough", GetPointerTo("TrkHough"));          t->SetBranchAddress("TrkHough", GetPointerTo("TrkHough"));
2252          cout << "Tracker      : set branch address TrkHough"<<endl;          cout << "Tracker      : set branch address TrkHough"<<endl;
2253       };      };
2254            
2255      // Calorimeter      // Calorimeter
2256      if(CAL2) {      if(CAL2) {
# Line 2268  void PamLevel2::SetBranchAddress(TChain Line 2302  void PamLevel2::SetBranchAddress(TChain
2302      };      };
2303      // SelectionList      // SelectionList
2304      if(SELLI==1) {      if(SELLI==1) {
2305    //      t->SetBranchAddress("RunEntry", &irun);
2306          t->SetBranchAddress("RunEntry", &irun);          t->SetBranchAddress("RunEntry", &irun);
2307          cout << "SelectionList: set branch address RunEntry"<<endl;          cout << "SelectionList: set branch address RunEntry"<<endl;
2308          t->SetBranchAddress("EventEntry", &irunentry);          t->SetBranchAddress("EventEntry", &irunentry);
# Line 2287  void PamLevel2::SetBranchAddress(TChain Line 2322  void PamLevel2::SetBranchAddress(TChain
2322   * @return Pointer to a TChain   * @return Pointer to a TChain
2323   */   */
2324  TChain *PamLevel2::GetRunTree(TList *fl){  TChain *PamLevel2::GetRunTree(TList *fl){
2325    //      //
2326    //      //
2327    //      //
2328    if ( run_tree ){      if ( run_tree ){
2329      printf("WARNING: TChain *PamLevel2::GetRunTree(TList *fl) -- run_tree already exists!\n ");          printf("WARNING: TChain *PamLevel2::GetRunTree(TList *fl) -- run_tree already exists!\n ");
2330      return run_tree;          return run_tree;
2331    };          };  
2332    //      //
2333    TChain *R = new TChain("Run");      TChain *R = new TChain("Run");
2334            
2335      // loop over files and create chains              // loop over files and create chains        
2336      TIter next(fl);      TIter next(fl);
# Line 2334  TChain *PamLevel2::GetRunTree(TList *fl) Line 2369  TChain *PamLevel2::GetRunTree(TList *fl)
2369   * @return Pointer to a TTree   * @return Pointer to a TTree
2370   */   */
2371  TTree *PamLevel2::GetRunTree(TFile *f){  TTree *PamLevel2::GetRunTree(TFile *f){
2372    if ( run_tree ){      if ( run_tree ){
2373      printf("WARNING: TTree *PamLevel2::GetRunTree(TFile *f) -- run_tree already exists!\n ");          printf("WARNING: TTree *PamLevel2::GetRunTree(TFile *f) -- run_tree already exists!\n ");
2374      return run_tree;          return run_tree;
2375    };      };
2376                    
2377    
2378      cout << "TTree *PamLevel2::GetRunTree(TFile *f) -- obsolte "<<endl;      cout << "TTree *PamLevel2::GetRunTree(TFile *f) -- obsolte "<<endl;
# Line 2362  TTree *PamLevel2::GetRunTree(TFile *f){ Line 2397  TTree *PamLevel2::GetRunTree(TFile *f){
2397   * @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
2398   */   */
2399  Bool_t PamLevel2::UpdateRunInfo(TChain *run, Long64_t iev){  Bool_t PamLevel2::UpdateRunInfo(TChain *run, Long64_t iev){
2400    //      //
2401    // check if we have already called once GetEntry, if not call it      // check if we have already called once GetEntry, if not call it
2402    //      //
2403      cout << "Bool_t PamLevel2::UpdateRunInfo(TChain *run, Long64_t iev) --- ATTENZIONE --- NON E` MANTENUTA!!!!!!!.... "<<endl;      cout << "Bool_t PamLevel2::UpdateRunInfo(TChain *run, Long64_t iev) --- ATTENZIONE --- NON E` MANTENUTA!!!!!!!.... "<<endl;
2404      if(!run){      if(!run){
2405            cout << "Bool_t PamLevel2::UpdateRunInfo(TChain *run, ULong64_t iev) -- ERROR -- missing RunInfo tree "<<endl;          cout << "Bool_t PamLevel2::UpdateRunInfo(TChain *run, ULong64_t iev) -- ERROR -- missing RunInfo tree "<<endl;
2406            return(false);                  return(false);  
2407      }      }
2408      if ( run->GetEntries() <= 0 ) return(false);      if ( run->GetEntries() <= 0 ) return(false);
2409    //      //
2410        
2411    //  Int_t oldrun = irun;      //  Int_t oldrun = irun;
2412    Long64_t oldrun = irun;      Long64_t oldrun = irun;
2413    
2414    // --------------------------------------      // --------------------------------------
2415    // if it is a full file (not preselected)      // if it is a full file (not preselected)
2416    // --------------------------------------      // --------------------------------------
2417    if(SELLI==0){      if(SELLI==0){
2418    
2419        //          //
2420        // the absolute time is necessary to relate the event with the run          // the absolute time is necessary to relate the event with the run
2421        //          //
2422        if( !GetOrbitalInfo() && !ISGP){          if( !GetOrbitalInfo() && !ISGP){
2423            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;
2424            return(false);              return(false);
2425        }          }
2426    
2427        ULong64_t abstime = 0;          ULong64_t abstime = 0;
2428        if( GetOrbitalInfo() )abstime=GetOrbitalInfo()->absTime;          if( GetOrbitalInfo() )abstime=GetOrbitalInfo()->absTime;
2429    
2430    
2431        //          //
2432        // the first time the routine is called, set run search from the beginning          // the first time the routine is called, set run search from the beginning
2433        //          //
2434        if ( irun < 0LL ){          if ( irun < 0LL ){
2435            irun = 0LL;              irun = 0LL;
2436            run->GetEntry(irun);              run->GetEntry(irun);
2437            runfirstentry = 0LL;              runfirstentry = 0LL;
2438            runlastentry += (Long64_t)(this->GetRunInfo()->NEVENTS);              runlastentry += (Long64_t)(this->GetRunInfo()->NEVENTS);
2439            if ( (Long64_t)(this->GetRunInfo()->NEVENTS) > 0LL ) runlastentry -= 1LL;              if ( (Long64_t)(this->GetRunInfo()->NEVENTS) > 0LL ) runlastentry -= 1LL;
2440                        
2441            if( ISGP && run->GetEntries()!=1 ){              if( ISGP && run->GetEntries()!=1 ){
2442                cout << "** WARNING ** simulated files are assumed to have 1 single run, not "<< run->GetEntries() << endl;                  cout << "** WARNING ** simulated files are assumed to have 1 single run, not "<< run->GetEntries() << endl;
2443                cout << "** WARNING ** run will not be updated"<<endl;                  cout << "** WARNING ** run will not be updated"<<endl;
2444            }              }
2445    
2446        };                };      
2447        //          //
2448        if(ISGP)abstime = GetRunInfo()->RUNHEADER_TIME; // BARBATRUCCO          if(ISGP)abstime = GetRunInfo()->RUNHEADER_TIME; // BARBATRUCCO
2449        //          //
2450        if ( irun == run->GetEntries()-1LL &&          if ( irun == run->GetEntries()-1LL &&
2451             !(abstime >= GetRunInfo()->RUNHEADER_TIME &&               !(abstime >= GetRunInfo()->RUNHEADER_TIME &&
2452               abstime <= GetRunInfo()->RUNTRAILER_TIME)                 abstime <= GetRunInfo()->RUNTRAILER_TIME)
2453              ){              ){
2454           irun = -1LL;              irun = -1LL;
2455           runfirstentry = 0LL;              runfirstentry = 0LL;
2456           runlastentry = -1LL;              runlastentry = -1LL;
2457         };          };
2458        // modificato il controllo sull'aggiornamento del run, per evitare problemi          // modificato il controllo sull'aggiornamento del run, per evitare problemi
2459        // dovuti agli eventi annidati (NB! NEVENTS conta anche questi!!)          // dovuti agli eventi annidati (NB! NEVENTS conta anche questi!!)
2460        //          //
2461        bool fromfirst = true;          bool fromfirst = true;
2462        //          //
2463        while ( !(abstime >= GetRunInfo()->RUNHEADER_TIME && abstime <= GetRunInfo()->RUNTRAILER_TIME) && irun < run->GetEntries()-1LL ){          while ( !(abstime >= GetRunInfo()->RUNHEADER_TIME && abstime <= GetRunInfo()->RUNTRAILER_TIME) && irun < run->GetEntries()-1LL ){
2464  //      while ( iev > (runfirstentry+(ULong64_t)(this->GetRunInfo()->NEVENTS-1)) && irun < run->GetEntries() ){  //      while ( iev > (runfirstentry+(ULong64_t)(this->GetRunInfo()->NEVENTS-1)) && irun < run->GetEntries() ){
2465            irun++;              irun++;
2466            run->GetEntry(irun);              run->GetEntry(irun);
2467            runfirstentry = runlastentry;              runfirstentry = runlastentry;
2468            if ( (Long64_t)(this->GetRunInfo()->NEVENTS) > 0LL ) runfirstentry += 1LL;              if ( (Long64_t)(this->GetRunInfo()->NEVENTS) > 0LL ) runfirstentry += 1LL;
2469            runlastentry += (Long64_t)(this->GetRunInfo()->NEVENTS);              runlastentry += (Long64_t)(this->GetRunInfo()->NEVENTS);
2470  //        cout << " ))))) UPDATE RUN INFO (((((  @iev "<<iev<<" run "<<GetRunInfo()->ID<<" irun "<<irun<<endl;  //        cout << " ))))) UPDATE RUN INFO (((((  @iev "<<iev<<" run "<<GetRunInfo()->ID<<" irun "<<irun<<endl;
2471  //        cout << "runfirstentry "<<runfirstentry<<endl;  //        cout << "runfirstentry "<<runfirstentry<<endl;
2472            //      printf(" iev %llu %u %llu \n",iev,this->GetRunInfo()->NEVENTS,(ULong64_t)(runfirstentry+(ULong64_t)(this->GetRunInfo()->NEVENTS)));              //    printf(" iev %llu %u %llu \n",iev,this->GetRunInfo()->NEVENTS,(ULong64_t)(runfirstentry+(ULong64_t)(this->GetRunInfo()->NEVENTS)));
2473            //      printf(" abstime %u trailertime %u \n",abstime,GetRunInfo()->RUNTRAILER_TIME);              //    printf(" abstime %u trailertime %u \n",abstime,GetRunInfo()->RUNTRAILER_TIME);
2474            //      printf(" IDRUN %u \n",GetRunInfo()->ID);              //    printf(" IDRUN %u \n",GetRunInfo()->ID);
2475            //              //
2476  //        prevshift = 0;  //        prevshift = 0;
2477            //              //
2478            if ( irun == (Long64_t)(run->GetEntries()-1LL) && fromfirst && !(abstime >= GetRunInfo()->RUNHEADER_TIME && abstime <= GetRunInfo()->RUNTRAILER_TIME)){              if ( irun == (Long64_t)(run->GetEntries()-1LL) && fromfirst && !(abstime >= GetRunInfo()->RUNHEADER_TIME && abstime <= GetRunInfo()->RUNTRAILER_TIME)){
2479                printf(" resetting irun  (it should NOT happen!!!)\n");                  printf(" resetting irun  (it should NOT happen!!!)\n");
2480                fromfirst = false;                  fromfirst = false;
2481                irun = 0;                  irun = 0;
2482                run->GetEntry(irun);                  run->GetEntry(irun);
2483                runfirstentry = 0ULL;                  runfirstentry = 0ULL;
2484                runlastentry += (Long64_t)(this->GetRunInfo()->NEVENTS);                  runlastentry += (Long64_t)(this->GetRunInfo()->NEVENTS);
2485                if ( (Long64_t)(this->GetRunInfo()->NEVENTS) > 0LL ) runlastentry -= 1LL;                  if ( (Long64_t)(this->GetRunInfo()->NEVENTS) > 0LL ) runlastentry -= 1LL;
2486            };              };
2487            //              //
2488        };          };
2489        //          //
2490        if (          if (
2491             !(abstime >= GetRunInfo()->RUNHEADER_TIME &&              !(abstime >= GetRunInfo()->RUNHEADER_TIME &&
2492               abstime <= GetRunInfo()->RUNTRAILER_TIME)                abstime <= GetRunInfo()->RUNTRAILER_TIME)
2493            ) {              ) {
2494          printf(" Something very wrong here: cannot find RUN containing absolute time %llu \n",abstime);              printf(" Something very wrong here: cannot find RUN containing absolute time %llu \n",abstime);
2495            return false;              return false;
2496        }          }
2497        //          //
2498        if ( irun == oldrun || irun >= run->GetEntries() ) return(false);          if ( irun == oldrun || irun >= run->GetEntries() ) return(false);
2499        //          //
2500        //  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);
2501        //          //
2502        prevshift = 0;          prevshift = 0;
2503        cout << " ))))) UPDATE RUN INFO (((((  @iev "<<iev<<" run "<<GetRunInfo()->ID<<" irun "<<irun<<endl;          cout << " ))))) UPDATE RUN INFO (((((  @iev "<<iev<<" run "<<GetRunInfo()->ID<<" irun "<<irun<<endl;
2504  //      cout << "runfirstentry "<<runfirstentry<<endl;  //      cout << "runfirstentry "<<runfirstentry<<endl;
2505        return(true);              return(true);    
2506    };      };
2507    // ----------------------------------------------------      // ----------------------------------------------------
2508    // if it is a preselected file (there is SelectionList)      // if it is a preselected file (there is SelectionList)
2509    // NBNB - the event tree MUST be read first      // NBNB - the event tree MUST be read first
2510    // ----------------------------------------------------      // ----------------------------------------------------
2511    if(SELLI==1){            if(SELLI==1){      
2512        sel_tree->GetEntry(iev);          sel_tree->GetEntry(iev);
2513  //      cout << irun << " "<< irunentry << endl;  //      cout << irun << " "<< irunentry << endl;
2514        if(irun != oldrun){          if(irun != oldrun){
2515            run->GetEntry(irun);              run->GetEntry(irun);
2516            cout << " ))))) UPDATE RUN INFO (((((  @iev "<<iev<<" run "<<GetRunInfo()->ID<<" irun "<<irun<<endl;              cout << " ))))) UPDATE RUN INFO (((((  @iev "<<iev<<" run "<<GetRunInfo()->ID<<" irun "<<irun<<endl;
2517            prevshift = 0;              prevshift = 0;
2518            return true;              return true;
2519        }          }
2520        return false;          return false;
2521    }      }
2522    
2523    return false;      return false;
2524    //      //
2525  };  };
2526    
2527  Bool_t PamLevel2::UpdateRunInfo(Long64_t iev){  Bool_t PamLevel2::UpdateRunInfo(Long64_t iev){
# Line 2579  Bool_t PamLevel2::UpdateRunInfo(Long64_t Line 2614  Bool_t PamLevel2::UpdateRunInfo(Long64_t
2614  //              || !(irunentry < GetRunInfo()->NEVENTS-1-prevshift) // ERRORE!!! fa saltare i run con 1 evento  //              || !(irunentry < GetRunInfo()->NEVENTS-1-prevshift) // ERRORE!!! fa saltare i run con 1 evento
2615                  || !(irunentry <= GetRunInfo()->NEVENTS-1-prevshift)                  || !(irunentry <= GetRunInfo()->NEVENTS-1-prevshift)
2616                  )                  )
2617                  && irun < run_tree->GetEntries() ){              && irun < run_tree->GetEntries() ){
2618    
2619              // - - - - - - - - - - - - -              // - - - - - - - - - - - - -
2620              // irunentry = position of current entry, relative to the run              // irunentry = position of current entry, relative to the run
# Line 2761  Bool_t PamLevel2::UpdateRunInfo(Long64_t Line 2796  Bool_t PamLevel2::UpdateRunInfo(Long64_t
2796   * @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
2797   */   */
2798  Bool_t PamLevel2::UpdateRunInfo(TTree *run, Long64_t iev){  Bool_t PamLevel2::UpdateRunInfo(TTree *run, Long64_t iev){
2799    return(UpdateRunInfo((TChain*)run,iev));      return(UpdateRunInfo((TChain*)run,iev));
2800  };  };
2801    
2802  //--------------------------------------  //--------------------------------------
# Line 2771  Bool_t PamLevel2::UpdateRunInfo(TTree *r Line 2806  Bool_t PamLevel2::UpdateRunInfo(TTree *r
2806  /**  /**
2807   * Set which trees shoul be analysed   * Set which trees shoul be analysed
2808   * @param detlist TString containing the sequence of trees required   * @param detlist TString containing the sequence of trees required
2809  */   */
2810  void PamLevel2::SetWhichTrees(TString detlist){  void PamLevel2::SetWhichTrees(TString detlist){
2811            
2812  //    if(detlist.IsNull() || detlist.Contains("+ALL", TString::kIgnoreCase)){  //    if(detlist.IsNull() || detlist.Contains("+ALL", TString::kIgnoreCase)){
# Line 3488  void PamLevel2::CreateCloneTrees(TFile * Line 3523  void PamLevel2::CreateCloneTrees(TFile *
3523    
3524    
3525      sel_tree_clone = new TTree("SelectionList","List of selected events ");      sel_tree_clone = new TTree("SelectionList","List of selected events ");
3526    //    sel_tree_clone->Branch("RunEntry",&irun,"runentry/L");
3527      sel_tree_clone->Branch("RunEntry",&irun,"runentry/L");      sel_tree_clone->Branch("RunEntry",&irun,"runentry/L");
3528      sel_tree_clone->Branch("EventEntry",&irunentry,"eventry/L");      sel_tree_clone->Branch("EventEntry",&irunentry,"eventry/L");
3529            
# Line 3919  Int_t PamLevel2::GetYodaEntry(){ Line 3955  Int_t PamLevel2::GetYodaEntry(){
3955  //      printf(" IDRUN %u \n",GetRunInfo()->ID);  //      printf(" IDRUN %u \n",GetRunInfo()->ID);
3956  //  //
3957          if ( prevshift != 0 && (quellagiusta+(Long64_t)shift) == GetYodaTree()->GetEntries() ){          if ( prevshift != 0 && (quellagiusta+(Long64_t)shift) == GetYodaTree()->GetEntries() ){
3958            prevshift = 0;              prevshift = 0;
3959            shift = -1;              shift = -1;
3960          };          };
3961    
3962      }while( ( obt != (UInt_t)(GetEventHeader()->GetPscuHeader()->GetOrbitalTime()) ||      }while( ( obt != (UInt_t)(GetEventHeader()->GetPscuHeader()->GetOrbitalTime()) ||

Legend:
Removed from v.1.54  
changed lines
  Added in v.1.55

  ViewVC Help
Powered by ViewVC 1.1.23