1 /****************************************************************************
3 * This macro determines the impact parameter of each track in pp event *
4 * with respect to the primaty vertex position reconstructed using all *
5 * the other tracks in the event. This procedure is necessary to avoid *
6 * biases in the impact parameter estimate. *
8 * Output: TNtuple (ntd0) on file ImpactParameters.root *
9 * There is one entry per track. The elements of the ntuple are: *
10 * 1) Event: Event number *
11 * 2) nTraks: Number of tracks for this event *
12 * 3) pt: Transverse momentum for the current track *
13 * 4) d0rphi: transverse impact parameter of the current track *
14 * 5) vx: x coordinate of the vertex *
15 * 6) vy: y coordinate of the vertex *
16 * origin: A.Dainese andrea.dainese@pd.infn.it *
17 ****************************************************************************/
18 #if !defined(__CINT__) || defined(__MAKECINT__)
19 //-- --- standard headers-------------
20 #include <Riostream.h>
21 //--------Root headers ---------------
26 #include <TStopwatch.h>
30 #include <TParticle.h>
32 //----- AliRoot headers ---------------
35 #include "AliKalmanTrack.h"
36 #include "AliITStrackV2.h"
37 #include "AliHeader.h"
38 #include "AliGenEventHeader.h"
39 #include "AliV0vertex.h"
40 #include "AliV0vertexer.h"
41 #include "AliITSVertex.h"
42 #include "AliITSVertexer.h"
43 #include "AliITSVertexerTracks.h"
45 //-------------------------------------
48 const Double_t kBz = 0.4;
50 // this function ckecks file existence
51 Bool_t GetInputFile();
52 void MoveTracks(TTree& itsTree,TObjArray& array);
55 void AliITSImpactParametersPP(Int_t evFirst=0,Int_t evLast=0) {
57 const Char_t *name="AliITSImpactParametersPP";
58 cerr<<'\n'<<name<<"...\n";
59 gBenchmark->Start(name);
61 AliKalmanTrack::SetConvConst(100/0.299792458/kBz);
64 // check existence of input file
65 if(!GetInputFile()) { cerr<<"No tracks file found"<<endl; return; }
67 TString outName("ImpactParameters.root");
70 Double_t v1[3] = {0.,0.,0,};
72 Int_t ev,itsEntries,i;
76 AliITStrackV2 *itstrack = 0;
79 TNtuple *ntd0 = new TNtuple("ntd0","ntd0","Event:nTrks:pt:d0rphi:vx:vy");
81 // create the AliITSVertexerTracks object
82 AliITSVertexerTracks *vertexer1 = new AliITSVertexerTracks;
83 vertexer1->SetMinTracks(3);
84 vertexer1->SetDebug(0);
87 // Open file with ITS tracks
88 TFile* itstrks = TFile::Open("AliITStracksV2.root");
91 // loop on events in file
92 for(ev=evFirst; ev<=evLast; ev++) {
93 printf(" --- Processing event %d ---\n",ev);
94 sprintf(trksName,"TreeT_ITS_%d",ev);
97 TTree *itsTree =(TTree*)itstrks->Get(trksName);
98 if(!itsTree) continue;
99 itsEntries = (Int_t)itsTree->GetEntries();
100 printf("+++\n+++ Number of tracks in ITS: %d\n+++\n\n",itsEntries);
102 TObjArray tArray(itsEntries);
103 MoveTracks(*itsTree,tArray);
105 // count the total number of events
109 for(i=0; i<itsEntries; i++) {
110 itstrack = (AliITStrackV2*)tArray.At(i);
113 pt = 1./TMath::Abs(itstrack->Get1Pt());
115 // primary vertex from other tracks in event
116 Bool_t goodVtx = kFALSE;
117 vertexer1->SetVtxStart(0.,0.);
120 vertexer1->SetSkipTracks(1,skipped);
121 AliITSVertex *vertex1 =
122 (AliITSVertex*)vertexer1->VertexOnTheFly(*itsTree);
123 if(vertex1->GetNContributors()>0) goodVtx = kTRUE;
125 //vertex1->PrintStatus();
126 // impact parameter w.r.t. v1 (in microns)
127 d0rphi = 10000.*itstrack->GetD(v1[0],v1[1]);
131 if (goodVtx) ntd0->Fill(ev,itsEntries,pt,d0rphi,v1[0],v1[1]);
138 } // loop on events in file
141 printf("\n+++\n+++ Total number of events: %d\n+++\n",nTotEv);
147 // store ntuple in a file
148 TFile* outroot = new TFile(outName.Data(),"recreate");
154 gBenchmark->Stop(name);
155 gBenchmark->Show(name);
159 //___________________________________________________________________________
160 Bool_t GetInputFile() {
161 TString itsName("AliITStracksV2.root");
163 if(gSystem->AccessPathName(itsName.Data(),kFileExists)) return kFALSE;
167 //___________________________________________________________________________
168 void MoveTracks(TTree& itsTree,TObjArray& array) {
170 // this function creates two TObjArrays with tracks
173 Int_t entr = (Int_t)itsTree.GetEntries();
175 // trasfer tracks from tree to array
176 for(Int_t i=0; i<entr; i++) {
178 AliITStrackV2 *t = new AliITStrackV2;
179 itsTree.SetBranchAddress("tracks",&t);