]>
Commit | Line | Data |
---|---|---|
8c555625 | 1 | void AliTPCTestTracking() { |
2 | const char *pname="Param1"; | |
3 | const char *tname="TreeD0_Param1"; | |
4 | ||
5 | // Dynamically link some shared libs | |
6 | if (gClassTable->GetID("AliRun") < 0) { | |
7 | gROOT->LoadMacro("loadlibs.C"); | |
8 | loadlibs(); | |
9 | } else { | |
10 | delete gAlice; | |
11 | gAlice=0; | |
12 | } | |
13 | ||
14 | // Connect the Root Galice file containing Geometry, Kine and Hits | |
15 | TFile *file = (TFile*)gROOT->GetListOfFiles()->FindObject("galice.root"); | |
16 | if (!file) file = new TFile("galice.root"); | |
17 | ||
18 | // Get AliRun object from file or create it if not on file | |
19 | if (!gAlice) { | |
20 | gAlice = (AliRun*)file->Get("gAlice"); | |
21 | if (gAlice) printf("AliRun object found on file\n"); | |
22 | if (!gAlice) gAlice = new AliRun("gAlice","Alice test program"); | |
23 | } | |
24 | ||
25 | gAlice->GetEvent(0); | |
26 | ||
27 | TClonesArray *particles=gAlice->Particles(); | |
28 | int np=particles->GetEntriesFast(); | |
3c0f9266 | 29 | int *good=new int[np]; |
30 | for (int ii=0; ii<np; ii++) good[ii]=0; | |
8c555625 | 31 | |
32 | AliTPC *TPC = (AliTPC*)gAlice->GetDetector("TPC"); | |
33 | int ver=TPC->IsVersion(); | |
34 | cerr<<"TPC version "<<ver<<" has been found !\n"; | |
35 | ||
36 | AliTPCD *digp= (AliTPCD*)file->Get(pname); | |
37 | if (digp!=0) TPC->SetDigParam(digp); | |
38 | else cerr<<"Warning: default TPC parameters will be used !\n"; | |
39 | ||
3c0f9266 | 40 | int nrow_up=TPC->GetDigParam()->GetParam().GetNRowUp(); |
41 | int nrows=TPC->GetDigParam()->GetParam().GetNRowLow()+nrow_up; | |
42 | int zero=TPC->GetDigParam()->GetParam().GetZeroSup(); | |
43 | int gap=int(0.125*nrows); | |
44 | int good_number=int(0.4*nrows); | |
45 | ||
8c555625 | 46 | switch (ver) { |
47 | case 1: | |
48 | cerr<<"Making clusters...\n"; | |
49 | TPC->Hits2Clusters(); | |
50 | TClonesArray *clusters=TPC->Clusters(); | |
51 | if (!clusters) {cerr<<"No clusters found !\n"; return;} | |
52 | int n=clusters->GetEntriesFast(); | |
53 | cerr<<"Number of clusters "<<n<<" \n"; | |
54 | ||
55 | cerr<<"Marking \"good\" tracks... \n"; | |
56 | for (int i=0; i<n; i++) { | |
57 | AliTPCcluster *c=(AliTPCcluster*)clusters->UncheckedAt(i); | |
58 | int lab=c->fTracks[0]; | |
59 | if (lab<0) continue; //noise cluster | |
60 | lab=TMath::Abs(lab); | |
3c0f9266 | 61 | int row=c->fPadRow; |
62 | if (row==nrow_up-1 ) good[lab]|=0x1000; | |
63 | if (row==nrow_up-1-gap) good[lab]|=0x800; | |
64 | good[lab]++; | |
8c555625 | 65 | } |
66 | ||
67 | break; | |
68 | case 2: | |
69 | cerr<<"Looking for clusters...\n"; | |
70 | TPC->Digits2Clusters(); | |
71 | TClonesArray *clusters=TPC->Clusters(); | |
72 | if (!clusters) {cerr<<"No clusters found !\n"; return;} | |
73 | int n=clusters->GetEntriesFast(); | |
74 | cerr<<"Number of clusters "<<n<<" \n"; | |
3c0f9266 | 75 | |
8c555625 | 76 | cerr<<"Marking \"good\" tracks... \n"; |
77 | TTree *TD=gDirectory->Get(tname); | |
78 | TClonesArray *digits=TPC->Digits(); | |
79 | TD->GetBranch("Digits")->SetAddress(&digits); | |
80 | ||
81 | int *count = new int[np]; | |
82 | int i; | |
83 | for (i=0; i<np; i++) count[i]=0; | |
84 | int sectors_by_rows=(int)TD->GetEntries(); | |
85 | for (i=0; i<sectors_by_rows; i++) { | |
86 | if (!TD->GetEvent(i)) continue; | |
3c0f9266 | 87 | int row; |
8c555625 | 88 | int ndigits=digits->GetEntriesFast(); |
89 | int j; | |
90 | for (j=0; j<ndigits; j++) { | |
91 | AliTPCdigit *dig = (AliTPCdigit*)digits->UncheckedAt(j); | |
92 | int idx0=dig->fTracks[0]; | |
93 | int idx1=dig->fTracks[1]; | |
94 | int idx2=dig->fTracks[2]; | |
3c0f9266 | 95 | row=dig->fPadRow; |
96 | if (idx0>=0 && dig->fSignal>=zero) count[idx0]+=1; | |
97 | if (idx1>=0 && dig->fSignal>=zero) count[idx1]+=1; | |
98 | if (idx2>=0 && dig->fSignal>=zero) count[idx2]+=1; | |
8c555625 | 99 | } |
100 | for (j=0; j<np; j++) { | |
8c555625 | 101 | if (count[j]>1) { |
102 | int ks; | |
3c0f9266 | 103 | if (row==nrow_up-1 ) good[j]|=0x1000; |
104 | if (row==nrow_up-1-gap) good[j]|=0x800; | |
105 | good[j]++; | |
8c555625 | 106 | } |
107 | count[j]=0; | |
108 | } | |
109 | } | |
110 | delete[] count; | |
111 | ||
112 | break; | |
113 | default: | |
114 | cerr<<"Invalid TPC version !\n"; | |
115 | return; | |
116 | } | |
117 | ||
118 | cerr<<"Looking for tracks...\n"; | |
119 | TPC->Clusters2Tracks(); | |
120 | int nt=0; | |
121 | TClonesArray *tracks=TPC->Tracks(); | |
122 | if (tracks) nt=tracks->GetEntriesFast(); | |
123 | cerr<<"Number of found tracks "<<nt<<endl; | |
124 | ||
125 | ///////////////////////////////////////////////////////////////////////// | |
126 | cerr<<"Doing comparison...\n"; | |
127 | TH1F *hp=new TH1F("hp","PHI resolution",50,-100.,100.); hp->SetFillColor(4); | |
128 | TH1F *hl=new TH1F("hl","LAMBDA resolution",50,-100,100); hl->SetFillColor(4); | |
129 | TH1F *hpt=new TH1F("hpt","Relative Pt resolution",30,-10.,10.); | |
130 | hpt->SetFillColor(2); | |
131 | TH1F *hd=new TH1F("hd","Impact parameter distribution ",30,0,25); | |
132 | hd->SetFillColor(6); | |
133 | ||
3c0f9266 | 134 | TH1F *hgood=new TH1F("hgood","Good tracks",20,0,10); |
135 | TH1F *hfound=new TH1F("hfound","Found tracks",20,0,10); | |
136 | TH1F *hfake=new TH1F("hfake","Fake tracks",20,0,10); | |
137 | TH1F *hg=new TH1F("hg","",20,0,10); //efficiency for good tracks | |
8c555625 | 138 | hg->SetLineColor(4); hg->SetLineWidth(2); |
3c0f9266 | 139 | TH1F *hf=new TH1F("hf","Efficiency for fake tracks",20,0,10); |
8c555625 | 140 | hf->SetFillColor(1); hf->SetFillStyle(3013); hf->SetLineWidth(2); |
141 | ||
3c0f9266 | 142 | TH1F *he =new TH1F("he","dE/dX for pions with 0.4<p<0.5 GeV/c",50,0.,500.); |
143 | TH2F *hep=new TH2F("hep","dE/dX vs momentum",50,0.,2.,50,0.,500.); | |
144 | ||
8c555625 | 145 | for (int i=0; i<np; i++) { |
3c0f9266 | 146 | TParticle *p = (TParticle*)particles->UncheckedAt(i); |
147 | if (p->GetFirstMother()>=0) continue; //secondary particle | |
148 | if (good[i] < 0x1000+0x800+2+good_number) continue; | |
149 | Double_t ptg=p->Pt(),pxg=p->Px(),pyg=p->Py(),pzg=p->Pz(); | |
8c555625 | 150 | if (ptg<0.100) continue; |
151 | if (fabs(pzg/ptg)>0.999) continue; | |
152 | ||
153 | //cout<<i<<endl; | |
154 | ||
155 | hgood->Fill(ptg); | |
156 | int found=0; | |
157 | for (int j=0; j<nt; j++) { | |
158 | AliTPCtrack *track=(AliTPCtrack*)tracks->UncheckedAt(j); | |
3c0f9266 | 159 | int lab=track->GetLabel(nrows); |
8c555625 | 160 | if (fabs(lab)!=i) continue; |
3c0f9266 | 161 | |
8c555625 | 162 | found=1; |
8c555625 | 163 | |
3c0f9266 | 164 | Double_t xk=76.; |
8c555625 | 165 | track->PropagateTo(xk); |
166 | xk-=0.11; | |
167 | track->PropagateTo(xk,42.7,2.27); //C | |
3c0f9266 | 168 | xk-=26.; |
8c555625 | 169 | track->PropagateTo(xk,36.2,1.98e-3); //C02 |
170 | xk-=0.051; | |
171 | track->PropagateTo(xk,42.7,2.27); //C | |
3c0f9266 | 172 | xk-=25.; |
173 | track->PropagateTo(xk,36.7,1.29e-3);//Air | |
174 | xk-=0.4; // + | |
175 | track->PropagateTo(xk,21.82,2.33);//ITS+beam_pipe+etc (approximately) | |
8c555625 | 176 | |
177 | track->PropagateToVertex(); //comparison should be done at the vertex | |
178 | ||
179 | if (lab==i) hfound->Fill(ptg); | |
180 | else { hfake->Fill(ptg); cerr<<lab<<" fake\n";} | |
181 | Double_t px,py,pz,pt=fabs(track->GetPt());track->GetPxPyPz(px,py,pz); | |
182 | Double_t phig=TMath::ATan(pyg/pxg); | |
183 | Double_t phi =TMath::ATan(py /px ); | |
184 | hp->Fill((phi - phig)*1000.); | |
185 | Double_t lamg=TMath::ATan(pzg/ptg); | |
186 | Double_t lam =TMath::ATan(pz /pt ); | |
187 | hl->Fill((lam - lamg)*1000.); | |
188 | hpt->Fill((pt - ptg)/ptg*100.); | |
3c0f9266 | 189 | Double_t x=track->GetX(), y=track->GetY(), z=track->GetZ(); |
8c555625 | 190 | hd->Fill(sqrt(x*x + y*y + z*z)*10.); |
3c0f9266 | 191 | |
192 | Double_t mom=track->GetP(); | |
193 | Double_t dedx=track->GetdEdX(0.05,0.80)/ | |
194 | digp->GetParam().GetPadPitchLength(); | |
195 | ||
196 | hep->Fill(mom,dedx,1.); | |
197 | if (p->GetPdgCode()==211 || p->GetPdgCode()==-211) | |
198 | if (mom>0.4 && mom<0.5) | |
199 | he->Fill(dedx,1.); | |
200 | ||
8c555625 | 201 | break; |
202 | } | |
203 | if (!found) cerr<<"Track number "<<i<<" was not found !\n"; | |
204 | } | |
3c0f9266 | 205 | |
206 | delete[] good; | |
207 | ||
8c555625 | 208 | Stat_t ngood =hgood->GetEntries(); cerr<<"Good tracks "<<ngood<<endl; |
209 | Stat_t nfound=hfound->GetEntries(); | |
210 | if (ngood!=0) | |
211 | cerr<<"Integral efficiency is about "<<nfound/ngood*100.<<" %\n"; | |
212 | ||
213 | gStyle->SetOptStat(111110); | |
214 | gStyle->SetOptFit(1); | |
215 | ||
216 | TCanvas *c1=new TCanvas("c1","",0,0,700,850); | |
217 | ||
218 | TPad *p1=new TPad("p1","",0,0.3,.5,.6); p1->Draw(); | |
219 | p1->cd(); p1->SetFillColor(42); p1->SetFrameFillColor(10); | |
220 | hp->SetFillColor(4); hp->SetXTitle("(mrad)"); hp->Fit("gaus"); c1->cd(); | |
221 | ||
222 | TPad *p2=new TPad("p2","",0.5,.3,1,.6); p2->Draw(); | |
223 | p2->cd(); p2->SetFillColor(42); p2->SetFrameFillColor(10); | |
224 | hl->SetXTitle("(mrad)"); hl->Fit("gaus"); c1->cd(); | |
225 | ||
226 | TPad *p3=new TPad("p3","",0,0,0.5,0.3); p3->Draw(); | |
227 | p3->cd(); p3->SetFillColor(42); p3->SetFrameFillColor(10); | |
228 | hpt->SetXTitle("(%)"); hpt->Fit("gaus"); c1->cd(); | |
229 | ||
230 | TPad *p4=new TPad("p4","",0.5,0,1,0.3); p4->Draw(); | |
231 | p4->cd(); p4->SetFillColor(42); p4->SetFrameFillColor(10); | |
232 | hd->SetXTitle("(mm)"); hd->Draw(); c1->cd(); | |
233 | ||
234 | TPad *p5=new TPad("p5","",0,0.6,1,1); p5->Draw(); p5->cd(); | |
3c0f9266 | 235 | p5->SetFillColor(41); p5->SetFrameFillColor(10); |
8c555625 | 236 | hfound->Sumw2(); hgood->Sumw2(); hfake->Sumw2(); |
237 | hg->Divide(hfound,hgood,1,1.,"b"); | |
238 | hf->Divide(hfake,hgood,1,1.,"b"); | |
239 | hg->SetMaximum(1.4); | |
240 | hg->SetYTitle("Tracking efficiency"); | |
241 | hg->SetXTitle("Pt (GeV/c)"); | |
242 | hg->Draw(); | |
243 | ||
244 | TLine *line1 = new TLine(0,1.0,2,1.0); line1->SetLineStyle(4); | |
245 | line1->Draw("same"); | |
246 | TLine *line2 = new TLine(0,0.9,2,0.9); line2->SetLineStyle(4); | |
247 | line2->Draw("same"); | |
248 | ||
249 | hf->SetFillColor(1); | |
250 | hf->SetFillStyle(3013); | |
251 | hf->SetLineColor(2); | |
252 | hf->SetLineWidth(2); | |
253 | hf->Draw("histsame"); | |
254 | TText *text = new TText(0.461176,0.248448,"Fake tracks"); | |
255 | text->SetTextSize(0.05); | |
256 | text->Draw(); | |
257 | text = new TText(0.453919,1.11408,"Good tracks"); | |
258 | text->SetTextSize(0.05); | |
259 | text->Draw(); | |
3c0f9266 | 260 | |
261 | ||
262 | ||
263 | TCanvas *c2=new TCanvas("c2","",320,32,530,590); | |
264 | ||
265 | TPad *p6=new TPad("p6","",0.,0.,1.,.5); p6->Draw(); | |
266 | p6->cd(); p6->SetFillColor(42); p6->SetFrameFillColor(10); | |
267 | he->SetFillColor(2); he->SetFillStyle(3005); | |
268 | he->SetXTitle("Arbitrary Units"); | |
269 | he->Fit("gaus"); c2->cd(); | |
270 | ||
271 | TPad *p7=new TPad("p7","",0.,0.5,1.,1.); p7->Draw(); | |
272 | p7->cd(); p7->SetFillColor(42); p7->SetFrameFillColor(10); | |
273 | hep->SetXTitle("p (Gev/c)"); hep->SetYTitle("dE/dX (Arb. Units)"); | |
274 | hep->Draw(); c1->cd(); | |
275 | ||
276 | ||
8c555625 | 277 | } |
278 |