From 1633a7ac7daa1fe90fd7069b8d383924eaaab355 Mon Sep 17 00:00:00 2001 From: hristov Date: Wed, 30 May 2001 06:04:58 +0000 Subject: [PATCH] Changes made to be consitant with changes in TPC tracking classes (B.Nilsen) --- ITS/AliITS.cxx | 804 +++++++++++++++++++++---------------------- ITS/AliTPCtracknew.C | 370 ++++++++++++++++++++ 2 files changed, 771 insertions(+), 403 deletions(-) create mode 100644 ITS/AliTPCtracknew.C diff --git a/ITS/AliITS.cxx b/ITS/AliITS.cxx index 2fddab339bf..f31ed6e3e60 100644 --- a/ITS/AliITS.cxx +++ b/ITS/AliITS.cxx @@ -15,6 +15,9 @@ /* $Log$ +Revision 1.51 2001/05/16 14:57:15 alibrary +New files for folders and Stack + Revision 1.50 2001/05/11 09:15:21 barbera Corrected to make fast point creation working with PPR geometry @@ -247,7 +250,7 @@ the AliITS class. #include "AliITSRad.h" #include "../TPC/AliTPC.h" #include "../TPC/AliTPCParam.h" - +#include "../TPC/AliTPCtracker.h" ClassImp(AliITS) @@ -717,9 +720,6 @@ void AliITS::SetDefaults() } } - - - //_____________________________________________________________________________ void AliITS::SetDefaultSimulation() { @@ -1290,8 +1290,6 @@ void AliITS::DigitsToRecPoints(Int_t evNumber,Int_t lastentry,Option_t *opt) treeC->Write(hname,TObject::kOverwrite); treeC->Reset(); } - - //____________________________________________________________________________ void AliITS::HitsToFastRecPoints(Int_t evNumber,Int_t bgrev,Int_t size, Option_t *option,Option_t *opt,Text_t *filename) @@ -1307,7 +1305,7 @@ Option_t *option,Option_t *opt,Text_t *filename) //if(ver!=5) return; const char *all = strstr(opt,"All"); - const char *det[3] = {strstr(opt,"SPD"),strstr(opt,"SDD"),strstr(opt,"SSD")}; + const char *det[3] ={strstr(opt,"SPD"),strstr(opt,"SDD"),strstr(opt,"SSD")}; Int_t nmodules; InitModules(size,nmodules); @@ -1335,7 +1333,8 @@ Option_t *option,Option_t *opt,Text_t *filename) AliITSDetType *iDetType=DetType(id); sim = (AliITSsimulation*)iDetType->GetSimulationModel(); if (!sim) { - Error("HitsToFastPoints","The simulation class was not instantiated!"); + Error("HitsToFastPoints", + "The simulation class was not instantiated!"); exit(1); // or SetDefaultSimulation(); } @@ -1368,14 +1367,13 @@ Option_t *option,Option_t *opt,Text_t *filename) delete [] random; } - //________________________________________________________________ -AliITStrack AliITS::Tracking(AliITStrack &track, AliITStrack *reference,TObjArray *fastpoints, Int_t -**vettid, Bool_t flagvert, AliITSRad *rl ) { - -//Origin A. Badala' and G.S. Pappalardo: e-mail Angela.Badala@ct.infn.it, Giuseppe.S.Pappalardo@ct.infn.it +AliITStrack AliITS::Tracking(AliITStrack &track, AliITStrack *reference, + TObjArray *fastpoints, Int_t **vettid, + Bool_t flagvert, AliITSRad *rl ) { +// Origin A. Badala' and G.S. Pappalardo: e-mail Angela.Badala@ct.infn.it, +// Giuseppe.S.Pappalardo@ct.infn.it - TList *list= new TList(); AliITStrack tr(track); @@ -1385,7 +1383,8 @@ AliITStrack AliITS::Tracking(AliITStrack &track, AliITStrack *reference,TObjArr Double_t Pt=(tr).GetPt(); // cout << "\n Pt = " << Pt <<"\n"; //stampa - AliITStracking obj(list, reference, this, fastpoints,TMath::Abs(Pt),vettid, flagvert, rl); + AliITStracking obj(list, reference, this, fastpoints, + TMath::Abs(Pt),vettid, flagvert, rl); list->Delete(); delete list; @@ -1395,17 +1394,18 @@ AliITStrack AliITS::Tracking(AliITStrack &track, AliITStrack *reference,TObjArr for(lay=5; lay>=0; lay--) { TVector VecLabref(3); VecLabref=(*reference).GetLabTrack(lay); - Float_t ClustZ=(*reference).GetZclusterTrack( lay); //aggiunta il 5-3-2001 - for(k=0; k<3; k++){ //{itot++; VecTotLabref(itot)=VecLabref(k);} //cambiata il 5-3-2001 - Int_t lpp=(Int_t)VecLabref(k); - if(lpp>=0) { - TParticle *p=(TParticle*) gAlice->Particle(lpp); - Int_t pcode=p->GetPdgCode(); - if(pcode==11) VecLabref(k)=p->GetFirstMother(); - } - itot++; VecTotLabref(itot)=VecLabref(k); - if(VecLabref(k)==0. && ClustZ == 0.) VecTotLabref(itot) =-3.; } - } + Float_t ClustZ=(*reference).GetZclusterTrack( lay); //aggiunta il 5-3-2001 + for(k=0; k<3; k++){ + //{itot++; VecTotLabref(itot)=VecLabref(k);} //cambiata il 5-3-2001 + Int_t lpp=(Int_t)VecLabref(k); + if(lpp>=0) { + TParticle *p=(TParticle*) gAlice->Particle(lpp); + Int_t pcode=p->GetPdgCode(); + if(pcode==11) VecLabref(k)=p->GetFirstMother(); + } // end if lpp>=0 + itot++; VecTotLabref(itot)=VecLabref(k); + if(VecLabref(k)==0. && ClustZ == 0.) VecTotLabref(itot) =-3.; } + } // end for lay Long_t labref; Int_t freq; (*reference).Search(VecTotLabref, labref, freq); @@ -1418,145 +1418,144 @@ AliITStrack AliITS::Tracking(AliITStrack &track, AliITStrack *reference,TObjArr return *reference; } +//________________________________________________________________------ +void AliITS::DoTracking(Int_t evNumber, Int_t min_t, Int_t max_t, + TFile *file, Bool_t flagvert) { + // ex macro for tracking ITS + printf("begin DoTracking - file %p\n",file); + //const char *pname="75x40_100x60"; + + Int_t imax=200,jmax=450; + AliITSRad *rl = new AliITSRad(imax,jmax); + //cout<<" dopo costruttore AliITSRad\n"; getchar(); -//________________________________________________________________ + struct GoodTrack { + Int_t lab,code; + Float_t px,py,pz,x,y,z,pxg,pyg,pzg,ptg; + Bool_t flag; + }; + gAlice->GetEvent(0); + AliKalmanTrack *kkprov; + kkprov->SetConvConst(100/0.299792458/0.2/gAlice->Field()->Factor()); -void AliITS::DoTracking(Int_t evNumber, Int_t min_t, Int_t max_t, TFile *file, Bool_t flagvert) { + /* //modificato il 26-4-2001 + AliTPC *TPC=(AliTPC*)gAlice->GetDetector("TPC"); + AliTPCParam *digp = (AliTPCParam*)file->Get(pname); + if (digp!=0) TPC->SetParam(digp); + */ + TFile *cf=TFile::Open("AliTPCclusters.root"); + AliTPCParam *digp= (AliTPCParam*)cf->Get("75x40_100x60"); + if (!digp) { cerr<<"TPC parameters have not been found !\n"; getchar();} + + AliTPCtracker *tracker = new AliTPCtracker(digp); //aggiunto il 23-5 + //aggiunto il 23-5 + // Load clusters + tracker->LoadInnerSectors(); + tracker->LoadOuterSectors(); + + GoodTrack gt[15000]; + Int_t ngood=0; + ifstream in("itsgood_tracks"); + + cerr<<"Reading itsgood tracks...\n"; + while (in>>gt[ngood].lab>>gt[ngood].code + >>gt[ngood].px >>gt[ngood].py>>gt[ngood].pz + >>gt[ngood].x >>gt[ngood].y >>gt[ngood].z + >>gt[ngood].pxg >>gt[ngood].pyg >>gt[ngood].pzg + >>gt[ngood].ptg >>gt[ngood].flag) { + ngood++; + cerr<IsOpen()) {cerr<<"Can't open AliTPCtracks.root !\n"; return ;} + TObjArray tracks(200000); + //TTree *tracktree=(TTree*)tf->Get("TreeT"); + TTree *tracktree=(TTree*)tf->Get("TPCf"); //aggiunto il 23-5 + if (!tracktree) {cerr<<"Can't get a tree with TPC tracks !\n";} + TBranch *tbranch=tracktree->GetBranch("tracks"); + Int_t nentr=(Int_t)tracktree->GetEntries(); + Int_t kk; + /* commentato il 26-4-2001 + for (kk=0; kkSetAddress(&iotrack); + tracktree->GetEvent(kk); + tracks.AddLast(iotrack); + } + */ + AliTPCtrack *iotracktpc=0; + for (kk=0; kkSetAddress(&iotracktpc); + tracktree->GetEvent(kk); + tracker->CookLabel(iotracktpc,0.1); //aggiunto 23-5 + tracks.AddLast(iotracktpc); + } // end for kk + delete tracker; + tf->Close(); -// ex macro for tracking ITS + Int_t nt = tracks.GetEntriesFast(); + cerr<<"Number of found tracks "<GetEvent(0); - - // AliKalmanTrack *kkprov; - // kkprov->SetConvConst(100/0.299792458/0.2/gAlice->Field()->Factor()); - - /* //modificato il 26-4-2001 - AliTPC *TPC=(AliTPC*)gAlice->GetDetector("TPC"); - AliTPCParam *digp = (AliTPCParam*)file->Get(pname); - if (digp!=0) TPC->SetParam(digp); - */ - TFile *cf=TFile::Open("AliTPCclusters.root"); - AliTPCParam *digp= (AliTPCParam*)cf->Get("75x40_100x60"); - if (!digp) { cerr<<"TPC parameters have not been found !\n"; getchar();} - - - GoodTrack gt[15000]; - Int_t ngood=0; - ifstream in("itsgood_tracks"); - - cerr<<"Reading itsgood tracks...\n"; - while (in>>gt[ngood].lab>>gt[ngood].code - >>gt[ngood].px >>gt[ngood].py>>gt[ngood].pz - >>gt[ngood].x >>gt[ngood].y >>gt[ngood].z - >>gt[ngood].pxg >>gt[ngood].pyg >>gt[ngood].pzg - >>gt[ngood].ptg >>gt[ngood].flag) { - ngood++; - cerr<IsOpen()) {cerr<<"Can't open AliTPCtracks.root !\n"; return ;} - TObjArray tracks(200000); - TTree *tracktree=(TTree*)tf->Get("TreeT"); - TBranch *tbranch=tracktree->GetBranch("tracks"); - Int_t nentr=(Int_t)tracktree->GetEntries(); - Int_t kk; - /* commentato il 26-4-2001 - for (kk=0; kkSetAddress(&iotrack); - tracktree->GetEvent(kk); - tracks.AddLast(iotrack); - } - */ - AliTPCtrack *iotracktpc=0; - for (kk=0; kkSetAddress(&iotracktpc); - tracktree->GetEvent(kk); - tracks.AddLast(iotracktpc); - } - tf->Close(); - - - Int_t nt = tracks.GetEntriesFast(); - cerr<<"Number of found tracks "<TreeR(); - Int_t nent=(Int_t)TR->GetEntries(); - TClonesArray *recPoints = RecPoints(); - Int_t numbpoints; - Int_t totalpoints=0; - Int_t *np = new Int_t[nent]; - Int_t **vettid = new Int_t* [nent]; - Int_t mod; - for (mod=0; modResetRecPoints(); - //gAlice->TreeR()->GetEvent(mod+1); //first entry in TreeR is empty - gAlice->TreeR()->GetEvent(mod); //first entry in TreeR is empty - numbpoints = recPoints->GetEntries(); - totalpoints+=numbpoints; - np[mod] = numbpoints; - //cout<<" mod = "<TreeR(); + Int_t nent=(Int_t)TR->GetEntries(); + TClonesArray *recPoints = RecPoints(); + Int_t numbpoints; + Int_t totalpoints=0; + Int_t *np = new Int_t[nent]; + Int_t **vettid = new Int_t* [nent]; + Int_t mod; + for (mod=0; modResetRecPoints(); + //gAlice->TreeR()->GetEvent(mod+1); //first entry in TreeR is empty + gAlice->TreeR()->GetEvent(mod); //first entry in TreeR is empty + numbpoints = recPoints->GetEntries(); + totalpoints+=numbpoints; + np[mod] = numbpoints; + //cout<<" mod = "<GetLabel()); - Int_t iii; - for (iii=0;iiiPropagateTo(xk); - xk-=0.11; - track->PropagateTo(xk,42.7,2.27); //C - xk-=2.6; - track->PropagateTo(xk,36.2,1.98e-3); //C02 - xk-=0.051; - track->PropagateTo(xk,42.7,2.27); //C - /////////////////////////////////////////////////// - */ - - ////// new propagation to the end of TPC ////////////// - Double_t xk=77.415; - track->PropagateTo(xk, 28.94, 1.204e-3); //Ne - xk -=0.01; - track->PropagateTo(xk, 44.77, 1.71); //Tedlar - xk -=0.04; - track->PropagateTo(xk, 44.86, 1.45); //Kevlar - xk -=2.0; - track->PropagateTo(xk, 41.28, 0.029); //Nomex - xk-=16; - track->PropagateTo(xk,36.2,1.98e-3); //C02 - xk -=0.01; - track->PropagateTo(xk, 24.01, 2.7); //Al - xk -=0.01; - track->PropagateTo(xk, 44.77, 1.71); //Tedlar - xk -=0.04; - track->PropagateTo(xk, 44.86, 1.45); //Kevlar - xk -=0.5; - track->PropagateTo(xk, 41.28, 0.029); //Nomex - - /////////////////////////////////////////////////////////////// - - /////////////////////////////////////////////////////////////// - AliITStrack trackITS(*track); - AliITStrack result(*track); - AliITStrack primarytrack(*track); + TTree tracktree1("TreeT","Tree with ITS tracks"); + AliITSiotrack *iotrack=0; + tracktree1.Branch("ITStracks","AliITSiotrack",&iotrack,32000,0); + + ofstream out ("AliITSTra.out"); + //ofstream outprova ("AliITSprova.out"); //commentato il 26-4-2001 + + Int_t j; + for (j=min_t; j<=max_t; j++) { + track=(AliTPCtrack*)tracks.UncheckedAt(j); + Int_t flaglab=0; + if (!track) continue; + ////// elimination of not good tracks //////////// + Int_t ilab=TMath::Abs(track->GetLabel()); + Int_t iii; + for (iii=0;iiiPropagateTo(xk); + xk-=0.11; + track->PropagateTo(xk,42.7,2.27); //C + xk-=2.6; + track->PropagateTo(xk,36.2,1.98e-3); //C02 + xk-=0.051; + track->PropagateTo(xk,42.7,2.27); //C + /////////////////////////////////////////////////// + */ + ////// new propagation to the end of TPC ////////////// + Double_t xk=77.415; + track->PropagateTo(xk, 28.94, 1.204e-3); //Ne + xk -=0.01; + track->PropagateTo(xk, 44.77, 1.71); //Tedlar + xk -=0.04; + track->PropagateTo(xk, 44.86, 1.45); //Kevlar + xk -=2.0; + track->PropagateTo(xk, 41.28, 0.029); //Nomex + xk-=16; + track->PropagateTo(xk,36.2,1.98e-3); //C02 + xk -=0.01; + track->PropagateTo(xk, 24.01, 2.7); //Al + xk -=0.01; + track->PropagateTo(xk, 44.77, 1.71); //Tedlar + xk -=0.04; + track->PropagateTo(xk, 44.86, 1.45); //Kevlar + xk -=0.5; + track->PropagateTo(xk, 41.28, 0.029); //Nomex + /////////////////////////////////////////////////////////////// + + /////////////////////////////////////////////////////////////// + AliITStrack trackITS(*track); + AliITStrack result(*track); + AliITStrack primarytrack(*track); -/////////////////////////////////////////////////////////////////////////////////////////////// - TVector Vgeant(3); - Vgeant=result.GetVertex(); - - // Definition of Dv and Zv for vertex constraint - Double_t sigmaDv=0.0050; Double_t sigmaZv=0.010; - //Double_t sigmaDv=0.0015; Double_t sigmaZv=0.0015; + ///////////////////////////////////////////////////////////////// + TVector Vgeant(3); + Vgeant=result.GetVertex(); + + // Definition of Dv and Zv for vertex constraint + Double_t sigmaDv=0.0050; Double_t sigmaZv=0.010; + //Double_t sigmaDv=0.0015; Double_t sigmaZv=0.0015; Double_t uniform= gRandom->Uniform(); Double_t signdv; if(uniform<=0.5) signdv=-1.; - else - signdv=1.; - - Double_t Vr=TMath::Sqrt(Vgeant(0)*Vgeant(0)+ Vgeant(1)*Vgeant(1)); - Double_t Dv=gRandom->Gaus(signdv*Vr,(Float_t)sigmaDv); - Double_t Zv=gRandom->Gaus(Vgeant(2),(Float_t)sigmaZv); - - //cout<<" Dv e Zv = "<GetLabel(); -// cout << " TPC track label = " << lab <<"\n"; // stampa - - -//propagation to vertex - - Double_t rbeam=3.; - - result.Propagation(rbeam); - - Double_t C00,C10,C11,C20,C21,C22,C30,C31,C32,C33,C40,C41,C42,C43,C44; - result.GetCElements(C00,C10,C11,C20,C21,C22,C30,C31,C32,C33,C40,C41,C42,C43,C44); - - Double_t pt=TMath::Abs(result.GetPt()); - Double_t Dr=result.GetD(); - Double_t Z=result.GetZ(); - Double_t tgl=result.GetTgl(); - Double_t C=result.GetC(); - Double_t Cy=C/2.; - Double_t Dz=Z-(tgl/Cy)*TMath::ASin(result.arga(rbeam)); - Dz-=Vgeant(2); - - // cout<<" Dr e dz alla fine = "<duepi) phivertex-=duepi; - if(phivertex<0.) phivertex+=duepi; - Double_t Dtot=TMath::Sqrt(Dr*Dr+Dz*Dz); - -////////////////////////////////////////////////////////////////////////////////////////// - - Int_t idmodule,idpoint; - if(NumofCluster >=5) { // cinque - sei - //if(NumofCluster ==6) { // cinque - sei - - - AliITSiotrack outtrack; - - iotrack=&outtrack; - - iotrack->SetStatePhi(phi); - iotrack->SetStateZ(Z); - iotrack->SetStateD(Dr); - iotrack->SetStateTgl(tgl); - iotrack->SetStateC(C); - Double_t radius=result.Getrtrack(); - iotrack->SetRadius(radius); - Int_t charge; - if(C>0.) charge=-1; else charge=1; - iotrack->SetCharge(charge); - - + else + signdv=1.; - iotrack->SetCovMatrix(C00,C10,C11,C20,C21,C22,C30,C31,C32,C33,C40,C41,C42,C43,C44); - - Double_t px=pt*TMath::Cos(phivertex); - Double_t py=pt*TMath::Sin(phivertex); - Double_t pz=pt*tgl; - - Double_t xtrack=Dr*TMath::Sin(phivertex); - Double_t ytrack=Dr*TMath::Cos(phivertex); - Double_t ztrack=Dz+Vgeant(2); - - - iotrack->SetPx(px); - iotrack->SetPy(py); - iotrack->SetPz(pz); - iotrack->SetX(xtrack); - iotrack->SetY(ytrack); - iotrack->SetZ(ztrack); - iotrack->SetLabel(labITS); - - Int_t il; - for(il=0;il<6; il++){ - iotrack->SetIdPoint(il,result.GetIdPoint(il)); - iotrack->SetIdModule(il,result.GetIdModule(il)); - } - tracktree1.Fill(); - - //cout<<" labITS = "<SetIdPoint(il,idpoint); - iotrack->SetIdModule(il,idmodule); - } - - // cout<<" +++++++++++++ pt e ptg = "<Gaus(signdv*Vr,(Float_t)sigmaDv); + Double_t Zv=gRandom->Gaus(Vgeant(2),(Float_t)sigmaZv); - - } // end if on NumofCluster - //gObjectTable->Print(); // stampa memoria - } // end for (int j=min_t; j<=max_t; j++) - - out.close(); - //outprova.close(); + //cout<<" Dv e Zv = "<GetLabel(); + // cout << " TPC track label = " << lab <<"\n"; // stampa + + //propagation to vertex + + Double_t rbeam=3.; + + result.Propagation(rbeam); + + Double_t C00,C10,C11,C20,C21,C22,C30,C31,C32,C33,C40,C41,C42,C43,C44; + result.GetCElements(C00,C10,C11,C20,C21,C22,C30,C31,C32,C33,C40,C41,C42,C43,C44); + Double_t pt=TMath::Abs(result.GetPt()); + Double_t Dr=result.GetD(); + Double_t Z=result.GetZ(); + Double_t tgl=result.GetTgl(); + Double_t C=result.GetC(); + Double_t Cy=C/2.; + Double_t Dz=Z-(tgl/Cy)*TMath::ASin(result.arga(rbeam)); + Dz-=Vgeant(2); + + // cout<<" Dr e dz alla fine = "<duepi) phivertex-=duepi; + if(phivertex<0.) phivertex+=duepi; + Double_t Dtot=TMath::Sqrt(Dr*Dr+Dz*Dz); + + ///////////////////////////////////////////////////////////////////// + + Int_t idmodule,idpoint; + if(NumofCluster >=5) { // cinque - sei + //if(NumofCluster ==6) { // cinque - sei + + AliITSiotrack outtrack; + + iotrack=&outtrack; + + iotrack->SetStatePhi(phi); + iotrack->SetStateZ(Z); + iotrack->SetStateD(Dr); + iotrack->SetStateTgl(tgl); + iotrack->SetStateC(C); + Double_t radius=result.Getrtrack(); + iotrack->SetRadius(radius); + Int_t charge; + if(C>0.) charge=-1; else charge=1; + iotrack->SetCharge(charge); + + + iotrack->SetCovMatrix(C00,C10,C11,C20,C21,C22,C30,C31,C32,C33,C40,C41,C42,C43,C44); + + Double_t px=pt*TMath::Cos(phivertex); + Double_t py=pt*TMath::Sin(phivertex); + Double_t pz=pt*tgl; + + Double_t xtrack=Dr*TMath::Sin(phivertex); + Double_t ytrack=Dr*TMath::Cos(phivertex); + Double_t ztrack=Dz+Vgeant(2); + + + iotrack->SetPx(px); + iotrack->SetPy(py); + iotrack->SetPz(pz); + iotrack->SetX(xtrack); + iotrack->SetY(ytrack); + iotrack->SetZ(ztrack); + iotrack->SetLabel(labITS); + + Int_t il; + for(il=0;il<6; il++){ + iotrack->SetIdPoint(il,result.GetIdPoint(il)); + iotrack->SetIdModule(il,result.GetIdModule(il)); + } // end for il + tracktree1.Fill(); + + //cout<<" labITS = "<SetIdPoint(il,idpoint); + iotrack->SetIdModule(il,idmodule); + } // end for il + + // cout<<" +++++++++++++ pt e ptg = "<Print(); // stampa memoria + } // end for (int j=min_t; j<=max_t; j++) + + out.close(); + //outprova.close(); - static Bool_t first=kTRUE; - static TFile *tfile; + static Bool_t first=kTRUE; + static TFile *tfile; - if(first) { - tfile=new TFile("itstracks.root","RECREATE"); - //cout<<"I have opened itstracks.root file "<cd(); - tfile->ls(); + if(first) { + tfile=new TFile("itstracks.root","RECREATE"); + //cout<<"I have opened itstracks.root file "<cd(); + tfile->ls(); - char hname[30]; - sprintf(hname,"TreeT%d",evNumber); + char hname[30]; + sprintf(hname,"TreeT%d",evNumber); - tracktree1.Write(hname); + tracktree1.Write(hname); - TTree *fAli=gAlice->TreeK(); - TFile *fileAli=0; + TTree *fAli=gAlice->TreeK(); + TFile *fileAli=0; - if (fAli) fileAli =fAli->GetCurrentFile(); - fileAli->cd(); + if (fAli) fileAli =fAli->GetCurrentFile(); + fileAli->cd(); - //////////////////////////////////////////////////////////////////////////////////////////////// + /////////////////////////////////////////////////////////////////// - printf("delete vectors\n"); - if(np) delete [] np; - if(vettid) delete [] vettid; - + printf("delete vectors\n"); + if(np) delete [] np; + if(vettid) delete [] vettid; } - diff --git a/ITS/AliTPCtracknew.C b/ITS/AliTPCtracknew.C new file mode 100644 index 00000000000..b63909ab3ee --- /dev/null +++ b/ITS/AliTPCtracknew.C @@ -0,0 +1,370 @@ +/**************************************************************************** + * This macro is supposed to do reconstruction in the barrel ALICE trackers * + * (Make sure you have TPC digits and ITS hits before using this macro !!!) * + * It does the following steps (April 12, 2001): * + * 1) TPC cluster finding * + * 2) TPC track finding * + * 3) ITS cluster finding V2 (via fast points !) * + * 4) ITS track finding V2 * + * 5) ITS back track propagation V2 * + * 6) TPC back track propagation * + * (Origin: I.Belikov, CERN, Jouri.Belikov@cern.ch) * + ****************************************************************************/ + +#ifndef __CINT__ + #include "alles.h" + #include "AliMagF.h" + #include "AliTPCtracker.h" + + #include "AliITS.h" + #include "AliITSgeom.h" + #include "AliITSRecPoint.h" + #include "AliITSclusterV2.h" + #include "AliITSsimulationFastPoints.h" + #include "AliITStrackerV2.h" +#endif + +Int_t TPCFindClusters(const Char_t *inname, const Char_t *outname); +Int_t TPCFindTracks(const Char_t *inname, const Char_t *outname); +Int_t TPCSortTracks(const Char_t *inname, const Char_t *outname); +Int_t TPCPropagateBack(const Char_t *inname, const Char_t *outname); + +Int_t ITSFindClusters(const Char_t *inname, const Char_t *outname); +Int_t ITSFindTracks(const Char_t *inname, const Char_t *outname); +Int_t ITSPropagateBack(const Char_t *inname, const Char_t *outname); + + +Int_t AliTPCtracknew() { + const Char_t *TPCdigName="galice.root"; + const Char_t *TPCclsName="AliTPCclusters.root"; + const Char_t *TPCtrkName="AliTPCtracks.root"; + + const Char_t *ITSdigName="galice.root"; + const Char_t *ITSclsName="AliITSclustersV2.root"; + const Char_t *ITStrkName="AliITStracksV2.root"; + + AliKalmanTrack::SetConvConst(100/0.299792458/0.2/gAlice->Field()->Factor()); + // AliKalmanTrack::SetConvConst(100/0.299792458/0.2/2); +// ********** Find TPC clusters *********** // + if (TPCFindClusters(TPCdigName,TPCclsName)) { + cerr<<"Failed to get TPC clusters !\n"; + return 1; + } + +// ********** Find TPC tracks *********** // + if (TPCFindTracks(TPCclsName,TPCtrkName)) { + cerr<<"Failed to get TPC tracks !\n"; + return 1; + } + +// ********** Sort and label TPC tracks *********** // + if (TPCSortTracks(TPCclsName,TPCtrkName)) { + cerr<<"Failed to sort TPC tracks !\n"; + return 1; + } +/* +// ********** Find ITS clusters *********** // + if (ITSFindClusters(ITSdigName,ITSclsName)) { + cerr<<"Failed to get ITS clusters !\n"; + return 1; + } + +// ********** Find ITS tracks *********** // + {TFile *clsFile=TFile::Open(ITSclsName); + if (ITSFindTracks(TPCtrkName,ITStrkName)) { + cerr<<"Failed to get ITS tracks !\n"; + return 1; + } + clsFile->Close();} + +// ********** Back propagation of the ITS tracks *********** // + {TFile *clsFile=TFile::Open(ITSclsName); + if (ITSPropagateBack(ITStrkName,TPCtrkName)) { + cerr<<"Failed to propagate back the ITS tracks !\n"; + return 1; + } + clsFile->Close();} + + +// ********** Back propagation of the TPC tracks *********** // + {TFile *clsFile=TFile::Open(TPCclsName); + if (TPCPropagateBack(TPCtrkName,TPCtrkName)) { + cerr<<"Failed to propagate back the TPC tracks !\n"; + return 1; + } + clsFile->Close();} +*/ + return 0; +} + + +Int_t TPCFindClusters(const Char_t *inname, const Char_t *outname) { + Int_t rc=0; + const Char_t *name="TPCFindClusters"; + cerr<<'\n'<Start(name); + + TFile *out=TFile::Open(outname,"recreate"); + TFile *in =TFile::Open(inname); + + AliTPCParam *param=(AliTPCParam *)in->Get("75x40_100x60"); + if (!param) {cerr<<"TPC parameters have not been found !\n"; return 1;} + AliTPCv2 tpc; + tpc.SetParam(param); + + tpc.Digits2Clusters(out); + + in->Close(); + out->Close(); + gBenchmark->Stop(name); + gBenchmark->Show(name); + + return rc; +} + +Int_t TPCFindTracks(const Char_t *inname, const Char_t *outname) { + Int_t rc=0; + const Char_t *name="TPCFindTracks"; + cerr<<'\n'<Start(name); + TFile *out=TFile::Open(outname,"recreate"); + TFile *in =TFile::Open(inname); + AliTPCParam *param=(AliTPCParam *)in->Get("75x40_100x60"); + if (!param) {cerr<<"TPC parameters have not been found !\n"; return 1;} + + AliTPCtracker *tracker=new AliTPCtracker(param); + rc=tracker->Clusters2Tracks(0,out); + delete tracker; + + in->Close(); + out->Close(); + gBenchmark->Stop(name); + gBenchmark->Show(name); + + return rc; +} + +Int_t TPCSortTracks(const Char_t *inname, const Char_t *outname) { + Int_t rc=0; + const Char_t *name="TPCSortTracks"; + cerr<<'\n'<Start(name); + + TFile *out=TFile::Open(outname); + TFile *in =TFile::Open(inname); + AliTPCParam *param=(AliTPCParam *)in->Get("75x40_100x60"); + if (!param) {cerr<<"TPC parameters have not been found !\n"; return 1;} + + TObjArray tarray(10000); + AliTPCtrack *iotrack=0; + Int_t i; + + out->cd(); + TTree *tracktree=(TTree*)out->Get("TPCf"); + TBranch *tbranch=tracktree->GetBranch("tracks"); + Int_t nentr=(Int_t)tracktree->GetEntries(); + for (i=0; iSetAddress(&iotrack); + tracktree->GetEvent(i); + tarray.AddLast(iotrack); + } + tarray.Sort(); + out->Close(); + + //assign thacks GEANT labels + in->cd(); + AliTPCtracker *tracker = new AliTPCtracker(param); + tracker->LoadInnerSectors(); + tracker->LoadOuterSectors(); + for (i=0; iCookLabel(iotrack,0.1); + } + delete tracker; + in->Close(); + //end of GEANT label assignment + + out=TFile::Open(outname,"recreate"); + tracktree=new TTree("TPCf","Tree with TPC tracks"); + tracktree->Branch("tracks","AliTPCtrack",&iotrack,32000,0); + for (i=0; iFill(); + } + tracktree->Write(); + out->Close(); + + gBenchmark->Stop(name); + gBenchmark->Show(name); + + return rc; +} + +Int_t ITSFindClusters(const Char_t *inname, const Char_t *outname) { + Int_t rc=0; + const Char_t *name="ITSFindClusters"; + cerr<<'\n'<Start(name); + TFile *out=TFile::Open(outname,"recreate"); + TFile *in =TFile::Open(inname,"update"); + + if (!(gAlice=(AliRun*)in->Get("gAlice"))) { + cerr<<"Can't get gAlice !\n"; + return 1; + } + Int_t ev=0; + gAlice->GetEvent(ev); + AliITS *ITS = (AliITS*)gAlice->GetModule("ITS"); + if (!ITS) { cerr<<"Can't get the ITS !\n"; return 1;} + + gAlice->MakeTree("R"); ITS->MakeBranch("R",0); + //////////////// Taken from ITSHitsToFastPoints.C /////////////////////// + AliITSsimulationFastPoints *sim = new AliITSsimulationFastPoints(); + for (Int_t i=0;i<3;i++) { ITS->SetSimulationModel(i,sim); } + Int_t nsignal=25; + Int_t size=-1; + Int_t bgr_ev=Int_t(ev/nsignal); + ITS->HitsToFastRecPoints(ev,bgr_ev,size," ","All"," "); + ////////////////////////////////////////////////////////////////////////// + + gAlice->GetEvent(ev); + + AliITSgeom *geom=ITS->GetITSgeom(); + out->cd(); + geom->Write(); + + TClonesArray *clusters=new TClonesArray("AliITSclusterV2",10000); + TTree *cTree=new TTree("cTree","ITS clusters"); + cTree->Branch("Clusters",&clusters); + + TTree *pTree=gAlice->TreeR(); + if (!pTree) { cerr<<"Can't get TreeR !\n"; return 1; } + TBranch *branch=pTree->GetBranch("ITSRecPoints"); + if (!branch) { cerr<<"Can't get ITSRecPoints branch !\n"; return 1;} + TClonesArray *points=new TClonesArray("AliITSRecPoint",10000); + branch->SetAddress(&points); + + TClonesArray &cl=*clusters; + Int_t nclusters=0; + Int_t nentr=(Int_t)pTree->GetEntries(); + for (Int_t i=0; iGetEvent(i)) {cTree->Fill(); continue;} + Int_t lay,lad,det; geom->GetModuleId(i,lay,lad,det); + Float_t x,y,zshift; geom->GetTrans(lay,lad,det,x,y,zshift); + Double_t rot[9]; geom->GetRotMatrix(lay,lad,det,rot); + Double_t yshift = x*rot[0] + y*rot[1]; + Int_t ndet=(lad-1)*geom->GetNdetectors(lay) + (det-1); + Int_t ncl=points->GetEntriesFast(); + nclusters+=ncl; + for (Int_t j=0; jUncheckedAt(j); + Float_t lp[5]; + lp[0]=-p->GetX()-yshift; if (lay==1) lp[0]=-lp[0]; + lp[1]=p->GetZ()+zshift; + lp[2]=p->GetSigmaX2(); + lp[3]=p->GetSigmaZ2(); + lp[4]=p->GetQ(); + Int_t lab[6]; + lab[0]=p->GetLabel(0);lab[1]=p->GetLabel(1);lab[2]=p->GetLabel(2); + lab[3]=ndet; + + Int_t label=lab[0]; + TParticle *part=(TParticle*)gAlice->Particle(label); + label=-3; + while (part->P() < 0.005) { + Int_t m=part->GetFirstMother(); + if (m<0) {cerr<<"Primary momentum: "<P()<Particle(label); + } + if (lab[1]<0) lab[1]=label; + else if (lab[2]<0) lab[2]=label; + else cerr<<"No empty labels !\n"; + + new(cl[j]) AliITSclusterV2(lab,lp); + } + cTree->Fill(); clusters->Delete(); + points->Delete(); + } + cTree->Write(); + cerr<<"Number of clusters: "<Close(); + out->Close(); + gBenchmark->Stop(name); + gBenchmark->Show(name); + + return rc; +} + +Int_t ITSFindTracks(const Char_t *inname, const Char_t *outname) { + Int_t rc=0; + const Char_t *name="ITSFindTracks"; + cerr<<'\n'<Start(name); + + AliITSgeom *geom=(AliITSgeom*)gFile->Get("AliITSgeom"); + if (!geom) { cerr<<"can't get ITS geometry !\n"; return 1;} + AliITStrackerV2 tracker(geom); + + TFile *out=TFile::Open(outname,"recreate"); + TFile *in =TFile::Open(inname); + rc=tracker.Clusters2Tracks(in,out); + in->Close(); + out->Close(); + + gBenchmark->Stop(name); + gBenchmark->Show(name); + + return rc; +} + +Int_t ITSPropagateBack(const Char_t *inname, const Char_t *outname) { + Int_t rc=0; + const Char_t *name="ITSPropagateBack"; + cerr<<'\n'<Start(name); + + AliITSgeom *geom=(AliITSgeom*)gFile->Get("AliITSgeom"); + if (!geom) { cerr<<"can't get ITS geometry !\n"; return 1;} + AliITStrackerV2 tracker(geom); + + TFile *out=TFile::Open(outname,"update"); + TFile *in =TFile::Open(inname); + rc=tracker.PropagateBack(in,out); + in->Close(); + out->Close(); + + gBenchmark->Stop(name); + gBenchmark->Show(name); + + return rc; +} + +Int_t TPCPropagateBack(const Char_t *inname, const Char_t *outname) { + Int_t rc=0; + const Char_t *name="TPCPropagateBack"; + cerr<<'\n'<Start(name); + + AliTPCParam *param=(AliTPCParam *)gFile->Get("75x40_100x60"); + if (!param) {cerr<<"TPC parameters have not been found !\n"; return 1;} + AliTPCtracker *tracker=new AliTPCtracker(param); + + TFile *out=TFile::Open(outname,"update"); + TFile *in =TFile::Open(inname); + rc=tracker->PropagateBack(in,out); + delete tracker; + in->Close(); + out->Close(); + + gBenchmark->Stop(name); + gBenchmark->Show(name); + + return rc; +} + -- 2.43.0