--- PamelaLevel2/src/PamLevel2.cpp 2007/03/21 09:07:32 1.33 +++ PamelaLevel2/src/PamLevel2.cpp 2007/03/28 15:29:14 1.37 @@ -85,9 +85,32 @@ // //-------------------------------------- /** - * Constructor + * Default Constructor */ PamLevel2::PamLevel2(){ + Initialize(); +}; + +/** + * Constructor with given dir, list and detectors + */ +PamLevel2::PamLevel2(TString ddir,TString list,TString detlist){ + Initialize(); + GetPamTree(GetListOfLevel2Files(ddir,list),detlist); + GetRunTree(GetListOfLevel2Files(ddir,list)); +}; + +/** + * Constructor with given dir and list + */ +PamLevel2::PamLevel2(TString ddir,TString list){ + Initialize(); + GetPamTree(GetListOfLevel2Files(ddir,list),""); + GetRunTree(GetListOfLevel2Files(ddir,list)); +}; + + +void PamLevel2::Initialize(){ // trk2_obj = TrkLevel2::GetTrkLevel2(); // trk1_obj = TrkLevel1::GetTrkLevel1(); @@ -177,6 +200,10 @@ tsorted=0; timage=0; + + howtosort = "+CAL+TOF"; + //howtosort = "+TOF"; + sortthr = 100.; }; /** @@ -218,7 +245,39 @@ } if(l0_file)l0_file->Close(); - if(pam_tree)pam_tree->Delete();; + // if(pam_tree)pam_tree->Delete();; + + if ( pam_tree ){ + // + // 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... + // + 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();; if(sel_tree)sel_tree->Delete();; for(Int_t i=0; iDelete();; @@ -247,7 +306,7 @@ // if(trk0_obj) trk0_obj->Clear(); if(calo0_obj) calo0_obj->Clear(); if(trk1_obj) trk1_obj->Clear(); - if(trk2_obj) trk2_obj->Clear(); + if(trk2_obj) trk2_obj->Clear(); if(trkh_obj) trkh_obj->Clear(); if(calo1_obj)calo1_obj->Clear(); if(calo2_obj)calo2_obj->Clear(); @@ -267,7 +326,81 @@ if(timage){ timage->Delete(); } - +}; + +void PamLevel2::Reset(){ + // + // First of all clear everything + // + Clear(); + // + // close and reset chains and pointers + // + if ( pam_tree ){ + // + // 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... + // + 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; + // + // Reset run pointers + // + run_obj = 0;//new GL_RUN(); + soft_obj = 0;// Emiliano + irun = -1; + runfirstentry = 0ULL; + runlastentry = 0ULL; + // }; @@ -464,6 +597,16 @@ } //-------------------------------------- // + +/** + * 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. + **/ +void PamLevel2::SortTracks(TString how){ + printf(" WARNING! obsolete, use SortTracks() and SetSortingMethod(TString) instead! \n Setting sorting method to %s \n",how.Data()); + howtosort = how; + SortTracks(); +}; + // //-------------------------------------- /** @@ -477,210 +620,345 @@ * * 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). */ -void PamLevel2::SortTracks(TString how){ +void PamLevel2::SortTracks(){ + TString how = howtosort; -// cout <<" PamLevel2::SortTracks(TString how) "<Clear(); -// sorted_tracks.Clear(); + // 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; + 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; - // loop over the tracks sorted by the tracker - Bool_t use_TRK = how.Contains("TRK", TString::kIgnoreCase); - Bool_t use_CAL = how.Contains("CAL", TString::kIgnoreCase); - Bool_t use_TOF = how.Contains("TOF", TString::kIgnoreCase); + // loop over the tracks sorted by the tracker + Bool_t use_TRK = how.Contains("TRK", TString::kIgnoreCase); + Bool_t use_CAL = how.Contains("CAL", TString::kIgnoreCase); + Bool_t use_TOF = how.Contains("TOF", TString::kIgnoreCase); - if( !CAL2 && use_CAL) use_CAL = false; - if( !TOF && use_TOF) use_TOF = false; + if( !CAL2 && use_CAL) use_CAL = false; + if( !TOF && use_TOF) use_TOF = false; - if( !TRK2 ){ -// cout << "SortTracks() : without tracker does not work!!! (not yet)" << endl; - return; - }; + if( !TRK2 ){ + // cout << "SortTracks() : without tracker does not work!!! (not yet)" << endl; + return; + }; - // cout << "use_CAL "<TrkLevel2::GetNTracks(); i++){ + for(Int_t i=0; i < trk2_obj->TrkLevel2::GetNTracks(); i++){ - TrkTrack *ts = 0; - CaloTrkVar *cs = 0; - ToFTrkVar *os = 0; - - // get tracker tracks - TrkTrack *tp = trk2_obj->TrkLevel2::GetTrack(i); //tracker - CaloTrkVar *cp = GetCaloStoredTrack(tp->GetSeqNo()); - ToFTrkVar *op = GetToFStoredTrack(tp->GetSeqNo()); - - TrkTrack *ti = 0; //tracker (image) - CaloTrkVar *ci = 0; - ToFTrkVar *oi = 0; -// cout << "trk track n. "<HasImage()){ + TrkTrack *ts = 0; + CaloTrkVar *cs = 0; + ToFTrkVar *os = 0; + + // get tracker tracks + TrkTrack *tp = trk2_obj->TrkLevel2::GetTrack(i); //tracker + CaloTrkVar *cp = GetCaloStoredTrack(tp->GetSeqNo()); + ToFTrkVar *op = GetToFStoredTrack(tp->GetSeqNo()); + + TrkTrack *ti = 0; //tracker (image) + CaloTrkVar *ci = 0; + ToFTrkVar *oi = 0; + // cout << "trk track n. "<HasImage()){ - ti = trk2_obj->TrkLevel2::GetTrackImage(i); //tracker (image) - ci = GetCaloStoredTrack(ti->GetSeqNo()); - oi = GetToFStoredTrack(ti->GetSeqNo()); + ti = trk2_obj->TrkLevel2::GetTrackImage(i); //tracker (image) + ci = GetCaloStoredTrack(ti->GetSeqNo()); + oi = GetToFStoredTrack(ti->GetSeqNo()); -// cout << "its image "<npcfit[1] > 15 && //no. of fit planes on Y view - calo2_obj->varcfit[1] < 1000. && //fit variance on Y view - cp && ci && - true){ + // ------------------------ + // calorimeter check + // ------------------------ + // check the Y spatial residual on the first calorimeter plane + // (cut on calorimeter variables from Emiliano) + if( use_CAL && !calo2_obj ){ + cout << "void PamLevel2::SortTracks(): howtosort= "<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; + 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(resy_p <= resy_i) tp_score++; + else ti_score++; + + // cout << "CALO "<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 ){ + //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; ipToFLevel2::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); + }; + }; + }; + if ( fabs(tr->y[3]-oi->ytofpos[1]) < 2. ){ + for (Int_t ih=0; ih < op->npmtadc; ih++){ + Int_t pl = tof_obj->GetPlaneIndex( (op->pmtadc).At(ih) ); + if(pl == 3){ + nphit_p++; + enhit_p += (op->dedx).At(ih); + }; + }; + }; + if ( fabs(tr->y[4]-oi->ytofpos[2]) < 2. ){ + for (Int_t ih=0; ih < op->npmtadc; ih++){ + Int_t pl = tof_obj->GetPlaneIndex( (op->pmtadc).At(ih) ); + if(pl == 4){ + nphit_p++; + enhit_p += (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 == 1 || pl == 2 || pl == 5){ + nphit_i++; + enhit_i += (op->dedx).At(ih); }; - // ------------------------ - // 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 && !tof_obj ){ - cout << "void PamLevel2::SortTracks(TString how): how= "<DoTrack2(tr); + // + if ( fabs(tr->y[0]-oi->ytofpos[0]) < 2. ){ + for (Int_t ih=0; ih < oi->npmtadc; ih++){ + Int_t pl = tof_obj->GetPlaneIndex( (oi->pmtadc).At(ih) ); + if(pl == 0){ + nphit_i++; + enhit_i += (op->dedx).At(ih); + }; + }; + }; + if ( fabs(tr->y[3]-oi->ytofpos[1]) < 2. ){ + for (Int_t ih=0; ih < oi->npmtadc; ih++){ + Int_t pl = tof_obj->GetPlaneIndex( (oi->pmtadc).At(ih) ); + if(pl == 3){ + nphit_i++; + enhit_i += (op->dedx).At(ih); + }; }; - if( use_TOF && op && oi ){ + }; + if ( fabs(tr->y[4]-oi->ytofpos[2]) < 2. ){ + for (Int_t ih=0; ih < oi->npmtadc; ih++){ + Int_t pl = tof_obj->GetPlaneIndex( (oi->pmtadc).At(ih) ); + if(pl == 4){ + nphit_i++; + enhit_i += (op->dedx).At(ih); + }; + }; + }; + - Int_t nphit_p =0; - Int_t nphit_i =0; + if( + use_TOF && + (nphit_p+nphit_i) !=0 && + true){ + + // 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(" score p %i score i %i \n",tp_score,ti_score); + // if( enhit_p > enhit_i ) tp_score++; + // if( nphit_p >= nphit_i && enhit_p > enhit_i ) tp_score++; + if ( nphit_p > nphit_i ) tp_score++; + if ( nphit_p < nphit_i ) ti_score++; + if ( nphit_p == nphit_i ){ + if ( enhit_p > enhit_i ) tp_score++; + else ti_score++; + }; + // printf(" dopo score p %i score i %i \n",tp_score,ti_score); + }; + delete tr; + // + } else { + // + // NOT a NUCLEUS + // + //printf(" NOT a NUCLEUS! en = %f \n",sen); + + Int_t nphit_p =0; + Int_t nphit_i =0; -/* cout << "track: npmtadc "<< op->npmtadc << endl; - cout << "track: npmttdc "<< op->npmttdc << endl; - cout << "image: npmtadc "<< oi->npmtadc << endl; - cout << "image: npmttdc "<< oi->npmttdc << endl;*/ + /* cout << "track: npmtadc "<< op->npmtadc << endl; + cout << "track: npmttdc "<< op->npmttdc << endl; + cout << "image: npmtadc "<< oi->npmtadc << endl; + cout << "image: npmttdc "<< oi->npmttdc << endl;*/ - 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++; - }; + 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++; + }; - 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++; - }; + 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++; + }; - if( - use_TOF && - (nphit_p+nphit_i) !=0 && - true){ + if( + use_TOF && + (nphit_p+nphit_i) !=0 && + true){ - if( nphit_p >= nphit_i) tp_score++; - else ti_score++; - }; -// cout << "TOF "<chi2 > 0 && tp->chi2 < ti->chi2 ) tp_score++ ; - else if( ti->chi2 > 0 && ti->chi2 < tp->chi2 ) ti_score++ ; -// cout << "TRK "<= nphit_i) tp_score++; + else ti_score++; + }; + }; + // cout << "TOF "<chi2 > 0 && tp->chi2 < ti->chi2 ) tp_score++ ; + else if( ti->chi2 > 0 && ti->chi2 < tp->chi2 ) ti_score++ ; + // cout << "TRK "< ti_score) { -// ts = tp;//the track sorted by the tracker!! -// cs = cp; -// os = op; - }else if (tp_score < ti_score) { - ts = ti;//its image!! - cs = ci; - os = oi; - - ti = tp;//its image!! - ci = cp; - oi = op; - - tp = ts;//its image!! - cp = cs; - op = os; + // ------------------------ + // the winner is.... + // ------------------------ + if (tp_score > ti_score) { + // ts = tp;//the track sorted by the tracker!! + // cs = cp; + // os = op; + // cout << "It's primary!" << endl; + + }else if (tp_score < ti_score) { + + // cout << "It's image!" << endl; + + ts = ti;//its image!! + cs = ci; + os = oi; + + ti = tp;//its image!! + ci = cp; + oi = op; + + tp = ts;//its image!! + cp = cs; + op = os; - }else { -// ts = tp; -// cs = cp; -// os = op; -// cout << "Warning - track image ambiguity not solved" << endl; -// cout << ts->GetNtot() << " "<< ts->chi2 << " " << npcfit[1] << " "<< nphit_p << endl; - }; + }else { + // ts = tp; + // cs = cp; + // os = op; + cout << "Warning - track image ambiguity not solved" << endl; + // cout << ts->GetNtot() << " "<< ts->chi2 << " " << npcfit[1] << " "<< nphit_p << endl; + }; - }else{ -// ts = tp; -// cs = cp; -// os = op; - }; - -// cout <<" SortTracks() "<Add(ts);//save the track in the sorted array -// sorted_tracks.Add(ts);//save the track in the sorted array -// sorted_tracks.Add(tp);//save the track in the sorted array -// cout << "SortTracks:: sorted_tracks->Add(it) "<GetEntries() != trk2_obj->GetNTracks() ){ - cout << "void PamLevel2::SortTracks(TString how): tsorted->GetEntries() "<GetEntries()<<" != trk2_obj->GetNTracks() = "<GetNTracks() <Delete(); tsorted=0; - timage->Delete(); timage=0; - } - - //Restore Object count - //To save space in the table keeping track of all referenced objects - //We reset the object count to what it was at the beginning of the event. - TProcessID::SetObjectCount(ObjectNumber); + }else{ + // ts = tp; + // cs = cp; + // os = op; + }; + + // cout <<" SortTracks() "<Add(ts);//save the track in the sorted array + // sorted_tracks.Add(ts);//save the track in the sorted array + // sorted_tracks.Add(tp);//save the track in the sorted array + // cout << "SortTracks:: sorted_tracks->Add(it) "<GetEntries() != trk2_obj->GetNTracks() ){ + cout << "void PamLevel2::SortTracks(): tsorted->GetEntries() "<GetEntries()<<" != trk2_obj->GetNTracks() = "<GetNTracks() <Delete(); tsorted=0; + timage->Delete(); timage=0; + } + + //Restore Object count + //To save space in the table keeping track of all referenced objects + //We reset the object count to what it was at the beginning of the event. + TProcessID::SetObjectCount(ObjectNumber); }; //-------------------------------------- @@ -704,7 +982,7 @@ TClonesArray *PamLevel2::GetTracks(){ // *-*-*-*-*-*-*-*-*-*-*-*-* - SortTracks("+CAL+TOF"); + SortTracks(); // *-*-*-*-*-*-*-*-*-*-*-*-* return tsorted; @@ -723,7 +1001,7 @@ // if(!trk2_obj) return 0; -// // *-*-*-*-*-*-*-*-*-*-*-*-* +// // *-*-*-*+-*-*-*-*-*-*-*-*-* // SortTracks("+CAL+TOF"); // // *-*-*-*-*-*-*-*-*-*-*-*-* // // if(!sorted_tracks)return 0; @@ -742,7 +1020,7 @@ // cout << "PamLevel2::GetTrack(int it) "<GetEntries())return 0; @@ -799,7 +1077,7 @@ // *-*-*-*-*-*-*-*-*-*-*-*-* - SortTracks("+CAL+TOF"); + SortTracks(); // *-*-*-*-*-*-*-*-*-*-*-*-* if(!timage)return 0; if(!timage->GetEntries())return 0; @@ -836,7 +1114,12 @@ // if( !detlist.IsNull() ) SetWhichTrees(detlist); // else GetWhichTrees(f); - + if ( pam_tree ){ + printf("WARNING: TTree *PamLevel2::GetPamTree(TFile *fl, TString detlist) -- pam_tree already exists!\n "); + return pam_tree; + }; + // + cout << "TTree *PamLevel2::GetPamTree(TFile *f, TString detlist ) -- obsolte "<SetBranchAddress("TrigLevel2", GetPointerTo("TrigLevel2")); cout << "Trigger : set branch address TrigLevel2"<AddFriend(R); }else{ cout << "Trigger : missing tree"<SetBranchAddress("S4Level2", GetPointerTo("S4Level2")); cout << "S4 : set branch address S4Level2"<AddFriend(S); }else{ cout << "S4 : missing tree"<SetBranchAddress("NDLevel2", GetPointerTo("NDLevel2")); cout << "NeutronD : set branch address NDLevel2"<AddFriend(N); }else{ cout << "NeutronD : missing tree"<SetBranchAddress("AcLevel2", GetPointerTo("AcLevel2")); cout << "Anticounter : set branch address AcLevel2"<AddFriend(A); }else{ cout << "Anticounter : missing tree"<SetBranchAddress("OrbitalInfo", GetPointerTo("OrbitalInfo")); cout << "OrbitalInfo : set branch address OrbitalInfo"<AddFriend(B); }else{ cout << "OrbitalInfo : missing tree"<WorkingDirectory(); - if(ddir=="")ddir = wdir; + // if(ddir=="")ddir = wdir; // TSystemDirectory *datadir = new TSystemDirectory(gSystem->BaseName(ddir),ddir); - cout << "Level2 data directory : "<< ddir << endl; - + if ( ddir != ""){ + cout << "Level2 data directory : "<< ddir << endl; + } else { + cout << "Level2 data directory not given as input: trying to evaluate from list or taking working directory " << endl; + }; TList *contents = new TList; // create output list contents->SetOwner(); @@ -1046,7 +1332,14 @@ // if( gSystem->IsFileInIncludePath(file,&fullpath) ){ // if( (fullpath = gSystem->FindFile(ddir,file)) ){ if( file.EndsWith(".root") ){ - char *fullpath = gSystem->ConcatFileName(gSystem->DirName(ddir),gSystem->BaseName(file)); + TString filedir; + if (ddir != ""){ + filedir = ddir; // take the input dir + } else { + gSystem->ChangeDirectory(wdir); // back to the working directory + filedir = gSystem->DirName(file); // this will take the path if exist in the list otherwise it will return automatically the working dir + }; + char *fullpath = gSystem->ConcatFileName(gSystem->DirName(filedir),gSystem->BaseName(file)); contents->Add(new TSystemFile(fullpath,gSystem->DirName(fullpath)));// add file to the list delete fullpath; } @@ -1056,6 +1349,7 @@ }; in.close(); } else { + if(flisttxt.Contains("#"))flisttxt = flisttxt(0,flisttxt.First("#")); char *fullpath = gSystem->ConcatFileName(gSystem->DirName(flisttxt),gSystem->BaseName(flisttxt)); contents->Add(new TSystemFile(fullpath,gSystem->DirName(fullpath)));// add file to the list delete fullpath; @@ -1065,7 +1359,11 @@ cout << "No input file list given."<WorkingDirectory()<< endl; - + if ( ddir == "" ){ + ddir = wdir; + cout << "Level2 data directory : "<< ddir << endl; + }; + TSystemDirectory *datadir = new TSystemDirectory(gSystem->BaseName(ddir),ddir); TList *temp = datadir->GetListOfFiles(); // temp->Print(); @@ -1109,6 +1407,14 @@ * @return Pointer to a TChain */ TChain *PamLevel2::GetPamTree(TList *fl, TString detlist ){ + // + // + // + if ( pam_tree ){ + printf("WARNING: TChain *PamLevel2::GetPamTree(TList *fl, TString detlist) -- pam_tree already exists!\n "); + return pam_tree; + }; + // TChain *Trout =0; @@ -1129,16 +1435,16 @@ if(ORB)cout<<"ORB "; cout << endl; - TChain *T = 0; - TChain *C = 0; - TChain *O = 0; - TChain *R = 0; - TChain *S = 0; - TChain *N = 0; - TChain *A = 0; - TChain *B = 0; + TChain *T = 0; + TChain *C = 0; + TChain *O = 0; + TChain *R = 0; + TChain *S = 0; + TChain *N = 0; + TChain *A = 0; + TChain *B = 0; - TChain *L = 0; + TChain *L = 0; if(TRK2||TRK1||TRKh) T = new TChain("Tracker"); if(CAL2||CAL1) C = new TChain("Calorimeter"); @@ -1148,7 +1454,6 @@ if(ND) N = new TChain("NeutronD"); if(AC) A = new TChain("Anticounter"); if(ORB) B = new TChain("OrbitalInfo"); - L = new TChain("SelectionList"); // loop over files and create chains @@ -1168,9 +1473,10 @@ if(ORB) B->Add(name); if(SELLI==1) L->Add(name); }; - } + }; cout << "done chain \n"; + cout << "+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+" <SetBranchStatus("TrkHough",0,found); if(TRKh)cout << "Tracker : set branch address TrkHough"<AddFriend("Tracker"); + else Trout->AddFriend("Tracker"); }else{ cout << "Tracker : missing tree"<AddFriend("SelectionList"); + cout << "----------------------------------------------------" <Delete(); @@ -1458,8 +1766,15 @@ * @return Pointer to a TChain */ TChain *PamLevel2::GetRunTree(TList *fl){ - - TChain *R = new TChain("Run"); + // + // + // + if ( run_tree ){ + printf("WARNING: TChain *PamLevel2::GetRunTree(TList *fl) -- run_tree already exists!\n "); + return run_tree; + }; + // + TChain *R = new TChain("Run"); // loop over files and create chains TIter next(fl); @@ -1498,21 +1813,26 @@ * @return Pointer to a TTree */ TTree *PamLevel2::GetRunTree(TFile *f){ + if ( run_tree ){ + printf("WARNING: TTree *PamLevel2::GetRunTree(TFile *f) -- run_tree already exists!\n "); + return run_tree; + }; + cout << "TTree *PamLevel2::GetRunTree(TFile *f) -- obsolte "<Get("Run"); + TTree *T = (TTree*)f->Get("Run"); - if(R){ - R->SetBranchAddress("RunInfo", GetPointerTo("RunInfo")); + if(T){ + T->SetBranchAddress("RunInfo", GetPointerTo("RunInfo")); cout << "Run : set branch address RunInfo"<SetBranchAddress("SoftInfo", GetPointerTo("SoftInfo")); // Emiliano + T->SetBranchAddress("SoftInfo", GetPointerTo("SoftInfo")); // Emiliano cout << "Software : set branch address SoftInfo"<GetEntries() <= 0 ) return(false); + cout << "Bool_t PamLevel2::UpdateRunInfo(TChain *run, Long64_t iev) --- NON FUNZIONA BENISSIMO.... "<GetEntries() <= 0 ) return(false); // -// Int_t oldrun = irun; + // Int_t oldrun = irun; Long64_t oldrun = irun; + // -------------------------------------- // if it is a full file (not preselected) // -------------------------------------- if(SELLI==0){ + // + // the absolute time is necessary to relate the event with the run + // if( !GetOrbitalInfo() ){ cout << "Bool_t PamLevel2::UpdateRunInfo(TChain *run, ULong64_t iev) -- ERROR -- missing OrbitalInfo "<GetEntry(irun); runfirstentry = 0LL; runlastentry += (Long64_t)(this->GetRunInfo()->NEVENTS); if ( (Long64_t)(this->GetRunInfo()->NEVENTS) > 0LL ) runlastentry -= 1LL; -// cout << " ))))) UPDATE RUN INFO ((((( @iev "<ID<<" irun "<GetEntries()-1LL && + !(GetOrbitalInfo()->absTime >= GetRunInfo()->RUNHEADER_TIME && + GetOrbitalInfo()->absTime <= GetRunInfo()->RUNTRAILER_TIME) + ){ + irun = -1LL; + runfirstentry = 0LL; + runlastentry = -1LL; + }; // modificato il controllo sull'aggiornamento del run, per evitare problemi // dovuti agli eventi annidati (NB! NEVENTS conta anche questi!!) // bool fromfirst = true; // - while ( !(GetOrbitalInfo()->absTime >= GetRunInfo()->RUNHEADER_TIME && GetOrbitalInfo()->absTime <= GetRunInfo()->RUNTRAILER_TIME) && irun < run->GetEntries()-1 ){ + while ( !(GetOrbitalInfo()->absTime >= GetRunInfo()->RUNHEADER_TIME && GetOrbitalInfo()->absTime <= GetRunInfo()->RUNTRAILER_TIME) && irun < run->GetEntries()-1LL ){ // while ( iev > (runfirstentry+(ULong64_t)(this->GetRunInfo()->NEVENTS-1)) && irun < run->GetEntries() ){ irun++; run->GetEntry(irun); @@ -1570,7 +1907,7 @@ // // prevshift = 0; // - if ( irun == run->GetEntries()-1 && fromfirst && !(GetOrbitalInfo()->absTime >= GetRunInfo()->RUNHEADER_TIME && GetOrbitalInfo()->absTime <= GetRunInfo()->RUNTRAILER_TIME)){ + if ( irun == (Long64_t)(run->GetEntries()-1LL) && fromfirst && !(GetOrbitalInfo()->absTime >= GetRunInfo()->RUNHEADER_TIME && GetOrbitalInfo()->absTime <= GetRunInfo()->RUNTRAILER_TIME)){ printf(" resetting irun (it should NOT happen!!!)\n"); fromfirst = false; irun = 0; @@ -1582,7 +1919,7 @@ // }; // - if ( !fromfirst && + if ( !(GetOrbitalInfo()->absTime >= GetRunInfo()->RUNHEADER_TIME && GetOrbitalInfo()->absTime <= GetRunInfo()->RUNTRAILER_TIME) ) { @@ -1596,7 +1933,7 @@ // prevshift = 0; cout << " ))))) UPDATE RUN INFO ((((( @iev "<ID<<" irun "<GetEntry(iev); +// cout << irun << " "<< irunentry << endl; if(irun != oldrun){ run->GetEntry(irun); cout << " ))))) UPDATE RUN INFO ((((( @iev "<ID<<" irun "<GetEntries() <= 0 ) { -// cout << " Bool_t PamLevel2::UpdateRunInfo(ULong64_t iev) -- ERROR -- run tree is empty"<GetEntries() <= 0 ) { + cout << " Bool_t PamLevel2::UpdateRunInfo(ULong64_t iev) -- ERROR -- run tree is empty"<GetEntry(irun); -// runfirstentry = iev; -// if( iev != 0){ -// cout << "Bool_t PamLevel2::UpdateRunInfo(TChain *run, ULong64_t iev) -- WARNING -- "<absTime >= GetRunInfo()->RUNHEADER_TIME && GetOrbitalInfo()->absTime <= GetRunInfo()->RUNTRAILER_TIME) && irun < run->GetEntries() ){ -// irun++; -// if( irun == run->GetEntries() ){ -// if(!fromfirst){ -// cout << "Bool_t PamLevel2::UpdateRunInfo(TChain *run, ULong64_t iev) -- ERROR -- event entry #"<GetEntry(irun); -// runfirstentry = iev; -// }; -// // -// if ( irun == oldrun ) return(false); -// // -// cout << " ))))) UPDATE RUN INFO ((((( @iev "<ID<<" irun "<GetEntry(iev); -// if(irun != oldrun){ -// run->GetEntry(irun); -// cout << " ))))) UPDATE RUN INFO ((((( @iev "<ID<<" irun "<GetEntry(irun); + }; + // + bool fromfirst = true; // first loop over runs + // ------------------------------------------------------ + // loop over runs to find the one that contains the event + // ------------------------------------------------------ + while ( !(GetOrbitalInfo()->absTime >= GetRunInfo()->RUNHEADER_TIME && // check on absolute time (s) + GetOrbitalInfo()->absTime <= GetRunInfo()->RUNTRAILER_TIME) && + !(GetOrbitalInfo()->OBT >= GetRunInfo()->RUNHEADER_OBT && // additional check on OBT (ms) + GetOrbitalInfo()->OBT <= GetRunInfo()->RUNTRAILER_OBT) && + irun < run_tree->GetEntries() ){ + irun++; + // ------------------------------------ + // if the end of run tree is reached... + // ------------------------------------ + if( irun == run_tree->GetEntries() ){ + if(!fromfirst){ + // ----------------------------------------------------- + // if it happened already once and the run was not found + // ---> exit with error + // ----------------------------------------------------- + cout << "Bool_t PamLevel2::UpdateRunInfo(Long64_t iev) -- ERROR -- event entry #"<0)runfirstentry += (GetRunInfo()->NEVENTS); + run_tree->GetEntry(irun); + if(GetRunInfo()->RUNHEADER_OBT>=GetRunInfo()->RUNTRAILER_OBT ){ + cout << "Bool_t PamLevel2::UpdateRunInfo(Long64_t iev) -- WARNING -- irun "<=RUNTRAILER_OBT " < exit with FALSE + // -------------------------------------- + if ( irun == oldrun ) return(false); + + // -------------------------------------- + // ... otherwise + // ---> exit with TRUE + // -------------------------------------- + cout << endl << " ))))) UPDATE RUN INFO ((((( @iev "<ID<<" irun "<GetEntry(iev); + if(irun != oldrun){ + run_tree->GetEntry(irun); + cout << endl << " ))))) UPDATE RUN INFO ((((( @iev "<ID<<" irun "<Branch("SoftInfo","SoftInfo",GetPointerTo("SoftInfo")); cout << "Run : branch SoftInfo"<GetEntries(); i++){ + run_tree->GetEntry(i); + cout << i<< "\t | "<ID<<"\t "<NEVENTS<< "\t "<RUNHEADER_TIME<<" <---> "<RUNTRAILER_TIME<Fill(); + } + cout << "----------------------------------------------------"<Branch("RunEntry",&irun,"runentry/I"); - sel_tree_clone->Branch("EventEntry",&irunentry,"eventry/I"); + sel_tree_clone->Branch("RunEntry",&irun,"runentry/L"); + sel_tree_clone->Branch("EventEntry",&irunentry,"eventry/L"); Int_t i=0; @@ -2349,7 +2727,7 @@ pam_tree_clone[i]->Branch("TrkLevel1","TrkLevel1", GetPointerTo("TrkLevel1")); pam_tree_clone[i]->BranchRef(); cout << "Tracker : branch TrkLevel1"<Branch("TrkLevel2", "TrkLevel2",GetPointerTo("TrkLevel2")); @@ -2510,24 +2888,31 @@ // if(!run_tree ){ if ( TRK0 || CAL0 || TOF0 || RUN ) { //forse cosi` va bene per tornare 1? - cout << " Int_t PamLevel2::GetEntry(Int_t) -- ERROR -- run tree not loaded"<Fill(); - +// if(UPDATED && run_tree_clone)run_tree_clone->Fill(); + // cout << "PamLevel2::GetEntry("<ID_ROOT_L0; @@ -2661,14 +3044,16 @@ // patch if( irunentry < 0){ - cout << "Int_t PamLevel2::GetYodaEntry() -- ATTENZIONE -- irunentry negativo?!?! "<<(Int_t)irunentry<EV_FROM); - +// cout << " irunentry "<EV_FROM<EV_FROM "<EV_FROM <<" quella giusta "<ID_ROOT_L0 "<ID_ROOT_L0<OBT==0 && GetOrbitalInfo()->pkt_num==0 ){ + cout << "Int_t PamLevel2::GetYodaEntry() -- ERROR -- level2 event corrupted ?? "<0){ - cout << " level0 <--> level2 mismatch ( irun "<pkt_num << " --- L0 --- "<< GetEventHeader()->GetPscuHeader()->GetCounter()<ID << " ID_ROOT_L0 "<ID_ROOT_L0<<" ID_RUN_FRAG "<ID_RUN_FRAG << " EV_FROM "<EV_FROM < L0 mismatch ( irun "<GetEntry(quellagiusta+(Long64_t)shift+(Long64_t)prevshift); shift++; @@ -2704,7 +3096,7 @@ cout << "Int_t PamLevel2::GetYodaEntry() -- ERROR -- missing EventHeader "<pkt_num << " --- L0 --- "<< GetEventHeader()->GetPscuHeader()->GetCounter()<pkt_num << " --- L0 --- "<< GetEventHeader()->GetPscuHeader()->GetCounter()<OBT << " --- L0 --- "<< GetEventHeader()->GetPscuHeader()->GetOrbitalTime()<GetEntries()+1 )cout << ">>> end of level0 tree <<<"<OBT "<< GetOrbitalInfo()->OBT << endl; @@ -2718,10 +3110,12 @@ shift = -1; }; - }while( ( GetOrbitalInfo()->OBT != (UInt_t)(GetEventHeader()->GetPscuHeader()->GetOrbitalTime()) || GetOrbitalInfo()->pkt_num != (UInt_t)(GetEventHeader()->GetPscuHeader()->GetCounter())) && (quellagiusta+(Long64_t)shift) < GetYodaTree()->GetEntries() ); - - if ( (quellagiusta+(Long64_t)shift) == GetYodaTree()->GetEntries() ) cout << " Big trouble here, no such event in Level0 data! " <OBT != (UInt_t)(GetEventHeader()->GetPscuHeader()->GetOrbitalTime()) || GetOrbitalInfo()->pkt_num != (UInt_t)(GetEventHeader()->GetPscuHeader()->GetCounter())) && (quellagiusta+(Long64_t)shift) < GetYodaTree()->GetEntries() && shift < maxshift); + if ( (quellagiusta+(Long64_t)shift) == GetYodaTree()->GetEntries() || shift == maxshift ) { + cout << " Big trouble here, no such event in Level0 data! " <GetEntry(quellagiusta); if ( shift > 1 ) prevshift += (shift-1);