]>
Commit | Line | Data |
---|---|---|
a144a99d | 1 | |
2 | /************************************************************************** | |
3 | * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * | |
4 | * * | |
5 | * Author: The ALICE Off-line Project. * | |
6 | * Contributors are mentioned in the code where appropriate. * | |
7 | * * | |
8 | * Permission to use, copy, modify and distribute this software and its * | |
9 | * documentation strictly for non-commercial purposes is hereby granted * | |
10 | * without fee, provided that the above copyright notice appears in all * | |
11 | * copies and that both the copyright notice and this permission notice * | |
12 | * appear in the supporting documentation. The authors make no claims * | |
13 | * about the suitability of this software for any purpose. It is * | |
14 | * provided "as is" without express or implied warranty. * | |
15 | **************************************************************************/ | |
16 | ||
17 | /* | |
18 | // Load libraries | |
19 | ||
20 | gSystem->Load("libANALYSIS"); | |
21 | gSystem->Load("libTPCcalib"); | |
22 | ||
23 | ||
24 | .x ~/NimStyle.C | |
25 | gSystem->Load("libANALYSIS"); | |
26 | gSystem->Load("libTPCcalib"); | |
27 | ||
28 | // analyze results | |
29 | ||
30 | TFile f("CalibObjectsTrain2.root"); | |
31 | AliTPCcalibMaterial *calibMaterial = (AliTPCcalibMaterial *)f->Get("alignMaterial"); | |
32 | ||
33 | ||
34 | */ | |
35 | ||
36 | #include "Riostream.h" | |
37 | #include "TChain.h" | |
38 | #include "TTree.h" | |
39 | #include "TH1F.h" | |
40 | #include "TH2F.h" | |
41 | #include "TH3F.h" | |
42 | #include "THnSparse.h" | |
43 | #include "TList.h" | |
44 | #include "TMath.h" | |
45 | #include "TCanvas.h" | |
46 | #include "TFile.h" | |
47 | #include "TF1.h" | |
48 | #include "TVectorD.h" | |
49 | #include "TProfile.h" | |
50 | #include "TGraphErrors.h" | |
51 | #include "TCanvas.h" | |
52 | ||
53 | #include "AliTPCclusterMI.h" | |
54 | #include "AliTPCseed.h" | |
55 | #include "AliESDVertex.h" | |
56 | #include "AliESDEvent.h" | |
57 | #include "AliESDfriend.h" | |
58 | #include "AliESDInputHandler.h" | |
59 | #include "AliAnalysisManager.h" | |
60 | ||
61 | #include "AliTracker.h" | |
62 | #include "AliMagF.h" | |
63 | #include "AliTPCCalROC.h" | |
64 | ||
65 | #include "AliLog.h" | |
66 | ||
67 | #include "AliTPCcalibMaterial.h" | |
68 | ||
69 | #include "TTreeStream.h" | |
70 | #include "AliTPCTracklet.h" | |
71 | #include "TTimeStamp.h" | |
72 | #include "AliTPCcalibDB.h" | |
73 | #include "AliTPCcalibLaser.h" | |
74 | #include "AliDCSSensorArray.h" | |
75 | #include "AliDCSSensor.h" | |
76 | ||
4af75575 | 77 | #include "AliKFParticle.h" |
78 | #include "AliKFVertex.h" | |
79 | ||
80 | #include "AliESDTZERO.h" | |
81 | ||
a144a99d | 82 | ClassImp(AliTPCcalibMaterial) |
83 | ||
84 | AliTPCcalibMaterial::AliTPCcalibMaterial(): | |
85 | AliTPCcalibBase("calibMaterial","calibMaterial"), | |
86 | fHisMaterial(0), | |
87 | fHisMaterialRPhi(0) | |
88 | { | |
89 | ||
90 | } | |
91 | ||
92 | AliTPCcalibMaterial::AliTPCcalibMaterial(const char * name, const char * title): | |
93 | AliTPCcalibBase(name,title), | |
94 | fHisMaterial(0), | |
95 | fHisMaterialRPhi(0) | |
96 | { | |
97 | // | |
98 | // | |
99 | // | |
100 | } | |
101 | ||
102 | AliTPCcalibMaterial::~AliTPCcalibMaterial(){ | |
103 | // | |
104 | // delete histograms | |
105 | // class is owner of all histograms | |
106 | // | |
107 | if (!fHisMaterial) return; | |
108 | delete fHisMaterial; | |
109 | delete fHisMaterialRPhi; | |
110 | fHisMaterial=0; | |
111 | } | |
112 | ||
113 | ||
114 | Long64_t AliTPCcalibMaterial::Merge(TCollection *li) { | |
115 | // | |
116 | // Merge histograms | |
117 | // | |
118 | TIterator* iter = li->MakeIterator(); | |
119 | AliTPCcalibMaterial* cal = 0; | |
120 | ||
121 | while ((cal = (AliTPCcalibMaterial*)iter->Next())) { | |
122 | if (!cal->InheritsFrom(AliTPCcalibMaterial::Class())) { | |
123 | Error("Merge","Attempt to add object of class %s to a %s", cal->ClassName(), this->ClassName()); | |
124 | return -1; | |
125 | } | |
126 | AliTPCcalibMaterial* calib= (AliTPCcalibMaterial*)(cal); | |
a144a99d | 127 | // |
128 | if (!fHisMaterial) fHisMaterial=MakeHisto(); | |
129 | fHisMaterial->Add(calib->fHisMaterial); | |
130 | fHisMaterialRPhi->Add(calib->fHisMaterialRPhi); | |
131 | } | |
132 | return 0; | |
133 | } | |
134 | ||
135 | ||
136 | ||
137 | void AliTPCcalibMaterial::Process(AliESDEvent *event){ | |
138 | // | |
139 | // | |
140 | // | |
141 | const Int_t kMinCl=40; | |
142 | const Float_t kMinRatio=0.7; | |
143 | const Float_t kMaxS=0.05; | |
144 | const Float_t kMinDist=5; | |
145 | const Double_t kStep=1.; | |
146 | if (!event) return; | |
4af75575 | 147 | ProcessPairs(event); |
148 | return; | |
a144a99d | 149 | // TTreeSRedirector * cstream = GetDebugStreamer(); |
150 | // | |
151 | if (!fHisMaterial){ | |
152 | MakeHisto(); | |
153 | } | |
154 | ||
155 | // | |
156 | // fill histogram of track prolongations | |
157 | Float_t dca[2]; | |
158 | Int_t ntracks = event->GetNumberOfTracks(); | |
a144a99d | 159 | for (Int_t itrack=0; itrack<ntracks; itrack++){ |
160 | AliESDtrack *track=event->GetTrack(itrack); | |
161 | if (!track) continue; | |
162 | if (track->GetTPCNcls()<=kMinCl) continue; | |
163 | if ((1.+track->GetTPCNcls())/(1.+track->GetTPCNclsF())<=kMinRatio) continue; | |
164 | if ((1.+track->GetTPCnclsS())/(1.+track->GetTPCNcls())>kMaxS) continue; | |
165 | if (!track->GetInnerParam()) continue; | |
166 | if (track->GetKinkIndex(0)!=0) continue; | |
167 | // | |
168 | track->GetImpactParameters(dca[0],dca[1]); | |
169 | if (TMath::Abs(dca[0])<kMinDist && TMath::Abs(dca[1])<kMinDist) continue; | |
170 | AliExternalTrackParam param(*(track->GetInnerParam())); | |
171 | if (!AliTracker::PropagateTrackTo(¶m,90,0.0005,10,kTRUE)) continue; | |
172 | Double_t x[5]={0,0,0,TMath::Sqrt(TMath::Abs(param.GetP()))*param.GetSign(),TMath::Sqrt(TMath::Abs(track->GetTPCsignal()))}; | |
173 | // | |
174 | // | |
175 | for (Float_t radius=90; radius>0; radius-=kStep){ | |
176 | if (!AliTracker::PropagateTrackTo(¶m,radius,0.0005,kStep*0.5,kTRUE)) break; | |
177 | if (TMath::Abs(param.GetSnp())>0.8) break; | |
178 | param.GetXYZ(x); | |
179 | Double_t weight=1./TMath::Sqrt(1.+param.GetSnp()*param.GetSnp()+param.GetTgl()*param.GetTgl()); | |
180 | fHisMaterial->Fill(x,weight); | |
181 | Double_t r = TMath::Sqrt(x[0]*x[0]+x[1]*x[1]); | |
182 | Double_t phi = TMath::ATan2(x[1],x[0]); | |
183 | x[0]=r; | |
184 | x[1]=phi; | |
185 | fHisMaterialRPhi->Fill(x,weight); | |
186 | } | |
187 | } | |
188 | } | |
189 | ||
190 | THnSparse *AliTPCcalibMaterial::MakeHisto(){ | |
191 | // | |
192 | // Make track prolongation histogram | |
193 | // | |
194 | // | |
195 | // gX gY gz p dEdx | |
196 | Int_t bins[5] = {100, 100, 300, 40, 100}; | |
197 | Double_t xmin[5] = {-100, -100, -300, -2, 5}; | |
198 | Double_t xmax[5] = {100, 100, 300, 2, 33}; | |
199 | TString axisName[5]={ | |
200 | "gx", | |
201 | "gy", | |
202 | "gz", | |
203 | "p", | |
204 | "dedx" | |
205 | }; | |
206 | TString axisTitle[5]={ | |
207 | "x (cm)", | |
208 | "y (cm)", | |
209 | "z (cm)", | |
210 | "p (GeV)", | |
211 | "dedx (a.u)" | |
212 | }; | |
213 | ||
214 | Int_t binsR[5] = {30, 360, 300, 40, 100}; | |
215 | Double_t xminR[5] = { 0, -3.14, -300, -2, 5}; | |
216 | Double_t xmaxR[5] = {30, 3.14, 300, 2, 33}; | |
217 | TString axisNameR[5]={ | |
218 | "r", | |
219 | "rphi", | |
220 | "z", | |
221 | "p", | |
222 | "dedx" | |
223 | }; | |
224 | TString axisTitleR[5]={ | |
225 | "r (cm)", | |
226 | "rphi (cm)", | |
227 | "z (cm)", | |
228 | "p (GeV)", | |
229 | "dedx (a.u)" | |
230 | }; | |
231 | ||
232 | THnSparse *sparse = new THnSparseF("his_Material", "His Material", 5, bins, xmin, xmax); | |
233 | THnSparse *sparseR = new THnSparseF("his_MaterialRPhi", "His Material Rphi", 5, binsR, xminR, xmaxR); | |
234 | for (Int_t iaxis=0; iaxis<5; iaxis++){ | |
235 | sparse->GetAxis(iaxis)->SetName(axisName[iaxis]); | |
236 | sparse->GetAxis(iaxis)->SetTitle(axisTitle[iaxis]); | |
237 | sparseR->GetAxis(iaxis)->SetName(axisNameR[iaxis]); | |
238 | sparseR->GetAxis(iaxis)->SetTitle(axisTitleR[iaxis]); | |
239 | } | |
240 | fHisMaterial=sparse; | |
241 | fHisMaterialRPhi=sparseR; | |
242 | return sparse; | |
243 | } | |
4af75575 | 244 | |
245 | ||
246 | ||
247 | ||
248 | void AliTPCcalibMaterial::ProcessPairs(AliESDEvent *event){ | |
249 | // | |
250 | // Process pairs of tracks to get a material budget map | |
251 | // | |
252 | // | |
253 | if (!event) return; | |
254 | if (event) AliKFParticle::SetField(event->GetMagneticField()); // set mean magnetic field for KF particles | |
255 | // | |
256 | // 1. Calculate total dEdx for all TPC tracks | |
257 | // | |
258 | const Int_t kMinCl=70; | |
259 | const Double_t kEpsilon=0.000001; | |
260 | const Float_t kMinRatio=0.7; | |
261 | const Float_t kMinDist=1.5; | |
262 | const Float_t kMinDistChi2=8; // | |
263 | const Float_t kMaxDistZ=280; // max distanceZ | |
264 | const Float_t kMaxDistR=250; // max distanceR | |
265 | const Double_t kMaxChi2 =36; // maximal chi2 to define the vertex | |
266 | const Double_t kMaxDistVertexSec=3; // maximal distance to secondary vertex | |
267 | ||
268 | if (!event) return; | |
269 | AliESDVertex *vertexSPD = (AliESDVertex *)event->GetPrimaryVertexSPD(); | |
270 | AliESDVertex * spdVertex = (AliESDVertex *)event->GetPrimaryVertexSPD(); | |
271 | AliESDVertex * trackVertex = (AliESDVertex *)event->GetPrimaryVertexTracks(); | |
272 | AliESDVertex * tpcVertex = (AliESDVertex *)event->GetPrimaryVertexTPC(); | |
273 | AliESDTZERO * tzero = (AliESDTZERO *)event->GetESDTZERO() ; | |
274 | // | |
275 | Double_t tpcSignalTotPrim=0; | |
276 | Double_t tpcSignalTotSec=0; | |
277 | Int_t ntracksTPC=0; | |
278 | Int_t nTPCPrim=0; | |
279 | Int_t nTPCSec=0; | |
280 | Int_t ntracks=event->GetNumberOfTracks(); | |
281 | if ( ntracks<=2 ) return; | |
282 | ||
283 | // | |
ac2d7e53 | 284 | Float_t dca[2]={0}; |
285 | Float_t cov[3]={0}; | |
286 | Float_t dca0[2]={0}; | |
287 | Float_t dca1[2]={0}; | |
4af75575 | 288 | // |
289 | //1. Calculate total dEdx for primary and secondary tracks | |
290 | // and count primaries and secondaries | |
291 | Int_t *rejectTrack = new Int_t[ntracks]; | |
292 | ||
293 | for (Int_t itrack=0; itrack<ntracks; itrack++){ | |
294 | AliESDtrack *track=event->GetTrack(itrack); | |
295 | rejectTrack[itrack]=0; | |
296 | if (!track) continue; | |
297 | if (track->GetTPCNcls()<=kMinCl) continue; // skip short tracks | |
298 | ntracksTPC++; | |
299 | if ((1.+track->GetTPCNcls())/(1.+track->GetTPCNclsF())<=kMinRatio) continue; | |
300 | if (!track->GetInnerParam()) continue; // skip not TPC tracks | |
301 | if (track->GetKinkIndex(0)!=0) {rejectTrack[itrack]+=16;continue;} // skip kinks | |
302 | track->GetImpactParameters(dca[0],dca[1]); | |
303 | if (TMath::Abs(dca[0])>kMaxDistR && TMath::Abs(dca[1])>kMaxDistZ) continue; | |
304 | // remove too dip secondaries | |
305 | if (TMath::Abs(dca[0])<kMinDist && TMath::Abs(dca[1])<kMinDist){ | |
306 | tpcSignalTotPrim+=track->GetTPCsignal(); | |
307 | nTPCPrim++; | |
308 | }else{ | |
309 | tpcSignalTotSec+=track->GetTPCsignal(); | |
310 | nTPCSec++; | |
311 | }; | |
312 | if (cov[0]>kEpsilon &&TMath::Abs(dca[0])>kEpsilon &&TMath::Sqrt(dca[0]*dca[0]/(TMath::Abs(cov[0])))<kMinDistChi2) rejectTrack[itrack]+=1; // primary | |
313 | if (cov[0]>kEpsilon &&TMath::Abs(dca[0])>kEpsilon &&TMath::Abs(dca[0])<kMinDist) rejectTrack[itrack]+=1; // primary | |
314 | if (track->GetTPCsignal()<40) rejectTrack[itrack]+=16; | |
315 | // | |
316 | if (CheckLooper(itrack, event)) rejectTrack[itrack]+=2; // looper | |
317 | if (CheckV0(itrack,event)) rejectTrack[itrack]+=4; | |
318 | } | |
319 | ||
320 | ||
321 | // | |
322 | // 2. Find secondary vertices - double loop | |
323 | // | |
324 | for (Int_t itrack0=0; itrack0<ntracks; itrack0++){ | |
325 | AliESDtrack *track0=event->GetTrack(itrack0); | |
326 | if (!track0) continue; | |
327 | if (track0->GetTPCNcls()<=kMinCl) continue; // skip short tracks | |
328 | if ((1.+track0->GetTPCNcls())/(1.+track0->GetTPCNclsF())<=kMinRatio) continue; | |
329 | if (!track0->GetInnerParam()) continue; // skip not TPC tracks | |
330 | if (track0->GetKinkIndex(0)>0) continue; // skip kinks | |
331 | if (rejectTrack[itrack0]) continue; // skip | |
332 | track0->GetImpactParameters(dca[0],dca[1]); | |
333 | track0->GetImpactParameters(dca0[0],dca0[1]); | |
334 | if (TMath::Abs(dca[0])>kMaxDistR && TMath::Abs(dca[1])>kMaxDistZ) continue; | |
335 | // remove too dip secondaries | |
336 | ||
337 | AliKFParticle part0(*track0,211); //assuming pion mass | |
338 | if (track0->Charge()*part0.Q()<0) part0.Q()*=-1; // change sign if opposite | |
339 | // | |
340 | for (Int_t itrack1=itrack0+1; itrack1<ntracks; itrack1++){ | |
341 | AliESDtrack *track1=event->GetTrack(itrack1); | |
342 | if (!track1) continue; | |
343 | if (rejectTrack[itrack1]) continue; // skip | |
344 | if (track1->GetTPCNcls()<=kMinCl) continue; // skip short tracks | |
345 | if ((1.+track1->GetTPCNcls())/(1.+track1->GetTPCNclsF())<=kMinRatio) continue; | |
346 | if (!track1->GetInnerParam()) continue; // skip not TPC tracks | |
347 | if (track1->GetKinkIndex(0)!=0) continue; // skip kinks | |
348 | track1->GetImpactParameters(dca1[0],dca1[1]); | |
349 | track1->GetImpactParameters(dca[0],dca[1]); | |
350 | if (TMath::Abs(dca[0])<kMinDist && TMath::Abs(dca[1])<kMinDist) continue; | |
351 | if (TMath::Abs(dca[0])>kMaxDistR && TMath::Abs(dca[1])>kMaxDistZ) continue; | |
352 | AliKFParticle part1(*track1,211); // assuming pion mass | |
353 | if (track1->Charge()*part1.Q()<0) part1.Q()*=-1; // change sign if opposite | |
354 | ||
355 | // | |
356 | // | |
357 | AliKFVertex vertex; | |
358 | vertex+=part0; | |
359 | vertex+=part1; | |
360 | if ((vertex.GetChi2()/vertex.GetNDF())> kMaxChi2) continue; | |
361 | if (TMath::Abs(vertex.GetX())>kMaxDistR) continue; | |
362 | if (TMath::Abs(vertex.GetY())>kMaxDistR) continue; | |
363 | if (TMath::Abs(vertex.GetZ())>kMaxDistZ) continue; | |
364 | Double_t errX2=vertex.GetErrX(); | |
365 | Double_t errY2=vertex.GetErrY(); | |
366 | Double_t errZ2=vertex.GetErrZ(); | |
367 | // | |
368 | Double_t err3D=TMath::Sqrt(errX2*errX2+errY2*errY2+errZ2*errZ2/10.); | |
369 | Double_t err2D=TMath::Sqrt(errX2*errX2+errY2*errY2); | |
370 | if (err3D>kMaxDistVertexSec) continue; | |
371 | if (err3D*TMath::Sqrt(vertex.GetChi2()+0.00001)>kMaxDistVertexSec) continue; | |
372 | ||
373 | Double_t dvertex=0; | |
374 | dvertex += (vertexSPD->GetX()-vertex.GetX())*(vertexSPD->GetX()-vertex.GetX()); | |
375 | dvertex += (vertexSPD->GetY()-vertex.GetY())*(vertexSPD->GetY()-vertex.GetY()); | |
376 | dvertex += (vertexSPD->GetZ()-vertex.GetZ())*(vertexSPD->GetZ()-vertex.GetZ()); | |
377 | dvertex=TMath::Sqrt(dvertex+0.00000001); | |
378 | if (err3D>0.2*dvertex) continue; | |
379 | if (err3D*TMath::Sqrt(vertex.GetChi2()+0.000001)>0.1*dvertex) continue; | |
380 | Double_t radius = TMath::Sqrt((vertex.GetX()*vertex.GetX()+vertex.GetY()*vertex.GetY())); | |
381 | // | |
382 | AliKFVertex vertex2; | |
383 | vertex2+=part0; | |
384 | vertex2+=part1; | |
385 | // | |
386 | ||
387 | for (Int_t itrack2=0; itrack2<ntracks; itrack2++){ | |
388 | if (itrack2==itrack0) continue; | |
389 | if (itrack2==itrack1) continue; | |
390 | if (rejectTrack[itrack2]) continue; // skip | |
391 | AliESDtrack *track2=event->GetTrack(itrack2); | |
392 | if (!track2) continue; | |
393 | if (track2->GetTPCNcls()<=kMinCl) continue; // skip short tracks | |
394 | if ((1.+track2->GetTPCNcls())/(1.+track2->GetTPCNclsF())<=kMinRatio) continue; | |
395 | if (!track2->GetInnerParam()) continue; // skip not TPC tracks | |
396 | if (track2->GetKinkIndex(0)>0) continue; // skip kinks | |
397 | track2->GetImpactParameters(dca[0],dca[1]); | |
398 | if (TMath::Abs(dca[0])<kMinDist && TMath::Abs(dca[1])<kMinDist) continue; | |
399 | if (TMath::Abs(dca[0])>kMaxDistR && TMath::Abs(dca[1])>kMaxDistZ) continue; | |
400 | if (TMath::Abs(track2->GetD(vertex.GetX(), vertex.GetY(),event->GetMagneticField()))>kMaxDistVertexSec) continue; | |
401 | Double_t vtxx[3]={vertex2.GetX(),vertex2.GetY(),vertex2.GetZ()}; | |
402 | Double_t svtxx[3]={vertex.GetErrX(),vertex.GetErrY(),vertex.GetErrZ()}; | |
403 | AliESDVertex vtx(vtxx,svtxx); | |
404 | AliExternalTrackParam param=*track2; | |
405 | Double_t delta[2]={0,0}; | |
406 | if (!param.PropagateToDCA(&vtx,event->GetMagneticField(),kMaxDistVertexSec,delta)) continue; | |
407 | if (TMath::Abs(delta[0])>kMaxDistVertexSec) continue; | |
408 | if (TMath::Abs(delta[1])>kMaxDistVertexSec) continue; | |
409 | if (TMath::Abs(delta[0])>6.*TMath::Sqrt(param.GetSigmaY2()+vertex2.GetErrY()*vertex2.GetErrY())+0.1) continue; | |
410 | if (TMath::Abs(delta[1])>6.*TMath::Sqrt(param.GetSigmaZ2()+vertex2.GetErrZ()*vertex2.GetErrZ())+0.5) continue; | |
411 | // | |
412 | AliKFParticle part2(param,211); // assuming pion mass | |
413 | if (track2->Charge()*part2.Q()<0) part2.Q()*=-1; // change sign if opposite | |
414 | vertex2+=part2; | |
415 | rejectTrack[itrack0]+=10; // do noit reuse the track | |
416 | rejectTrack[itrack1]+=10; // do not reuse the track | |
417 | rejectTrack[itrack2]+=10; | |
418 | } | |
419 | ||
420 | ||
421 | TTreeSRedirector *pcstream = GetDebugStreamer(); | |
422 | if (pcstream){ | |
423 | // | |
424 | // | |
425 | // | |
426 | Float_t dedx0= track0->GetTPCsignal(); | |
427 | Float_t dedx1= track1->GetTPCsignal(); | |
428 | AliExternalTrackParam * p0= (AliExternalTrackParam *)track0->GetInnerParam(); | |
429 | AliExternalTrackParam * p1= (AliExternalTrackParam *)track1->GetInnerParam(); | |
430 | Double_t errX=vertex2.GetErrX(); | |
431 | Double_t errY=vertex2.GetErrY(); | |
432 | Double_t errZ=vertex2.GetErrZ(); | |
433 | Double_t vx = vertex2.GetX(); | |
434 | Double_t vy = vertex2.GetY(); | |
435 | Double_t vz = vertex2.GetZ(); | |
436 | (*pcstream)<<"mapTPC"<< | |
437 | "run="<<fRun<< // run | |
438 | "event="<<fEvent<< // event number | |
439 | "time="<<fTime<< // timeStamp | |
440 | "trigger="<<fTrigger<< // trigger | |
441 | "mag="<<fMagF<< // magnetic field | |
442 | // | |
443 | "spdV.="<<spdVertex<< // spd vertex | |
444 | "trackV.="<<trackVertex<< // track vertex | |
445 | "tpcV.="<<tpcVertex<< // track vertex | |
446 | "tzero.="<<tzero<< // tzero info | |
447 | // | |
448 | "ntracks="<<ntracks<< | |
449 | "ntracksTPC="<<ntracksTPC<< | |
450 | "nPrim="<<nTPCPrim<< // number of primaries | |
451 | "nSec="<<nTPCSec<< // number of secondaries | |
452 | "sigPrim="<<tpcSignalTotPrim<< // total dEdx in primaries | |
453 | "sigSec="<<tpcSignalTotSec<< // total dEdx in secondaries | |
454 | "dedx0="<<dedx0<< // dedx part 0 | |
455 | "dedx1="<<dedx1<< // dedx part 1 | |
456 | "p0.="<<p0<< // part 0 | |
457 | "p1.="<<p1<< //part 1 | |
458 | "v.="<<&vertex<< // KF vertex | |
459 | "v2.="<<&vertex2<< // KF vertex all tracks | |
460 | "z0="<<dca0[1]<< | |
461 | "z1="<<dca1[1]<< | |
462 | "rphi0="<<dca0[0]<< | |
463 | "rphi1="<<dca1[0]<< | |
464 | "dvertex="<<dvertex<< | |
465 | "radius="<<radius<< | |
466 | "vx="<<vx<< | |
467 | "vy="<<vy<< | |
468 | "vz="<<vz<< | |
469 | "errX="<<errX<< | |
470 | "errY="<<errY<< | |
471 | "errZ="<<errZ<< | |
472 | "err2D="<<err2D<< | |
473 | "err3D="<<err3D<< | |
474 | "\n"; | |
475 | // | |
476 | if (vertex2.GetNDF()>2){ | |
477 | (*pcstream)<<"mapVertex"<< | |
478 | "run="<<fRun<< // run | |
479 | "event="<<fEvent<< // event number | |
480 | "time="<<fTime<< // timeStamp | |
481 | "trigger="<<fTrigger<< // trigger | |
482 | "mag="<<fMagF<< // magnetic field | |
483 | // | |
484 | "spdV.="<<spdVertex<< // spd vertex | |
485 | "trackV.="<<trackVertex<< // track vertex | |
486 | "tpcV.="<<tpcVertex<< // track vertex | |
487 | "tzero.="<<tzero<< // tzero info | |
488 | // | |
489 | // | |
490 | "ntracks="<<ntracks<< | |
491 | "ntracksTPC="<<ntracksTPC<< | |
492 | "nPrim="<<nTPCPrim<< // number of primaries | |
493 | "nSec="<<nTPCSec<< // number of secondaries | |
494 | "sigPrim="<<tpcSignalTotPrim<< // total dEdx in primaries | |
495 | "sigSec="<<tpcSignalTotSec<< // total dEdx in secondaries | |
496 | "dedx0="<<dedx0<< // dedx part 0 | |
497 | "dedx1="<<dedx1<< // dedx part 1 | |
498 | "p0.="<<p0<< // part 0 | |
499 | "p1.="<<p1<< //part 1 | |
500 | "v.="<<&vertex<< // KF vertex | |
501 | "v2.="<<&vertex2<< // KF vertex all tracks | |
502 | "z0="<<dca0[1]<< | |
503 | "z1="<<dca1[1]<< | |
504 | "rphi0="<<dca0[0]<< | |
505 | "rphi1="<<dca1[0]<< | |
506 | "radius="<<radius<< | |
507 | "vx="<<vx<< | |
508 | "vy="<<vy<< | |
509 | "vz="<<vz<< | |
510 | "errX="<<errX<< | |
511 | "errY="<<errY<< | |
512 | "errZ="<<errZ<< | |
513 | "err2D="<<err2D<< | |
514 | "err3D="<<err3D<< | |
515 | "dvertex="<<dvertex<< | |
516 | "\n"; | |
517 | } | |
518 | } | |
519 | } | |
520 | } | |
521 | TTreeSRedirector *pcstream = GetDebugStreamer(); | |
522 | if (pcstream){ | |
523 | (*pcstream)<<"statTPC"<< | |
524 | "run="<<fRun<< // run | |
525 | "time="<<fTime<< // timeStamp | |
526 | "trigger="<<fTrigger<< // trigger | |
527 | "mag="<<fMagF<< // magnetic field | |
528 | "ntracks="<<ntracks<< | |
529 | "ntracksTPC="<<ntracksTPC<< | |
530 | // | |
531 | "nPrim="<<nTPCPrim<< // number of primaries | |
532 | "nSec="<<nTPCSec<< // number of secondaries | |
533 | "sigPrim="<<tpcSignalTotPrim<< // total dEdx in primaries | |
534 | "sigSec="<<tpcSignalTotSec<< // total dEdx in secondaries | |
535 | // | |
536 | "spdV.="<<spdVertex<< // spd vertex | |
537 | "trackV.="<<trackVertex<< // track vertex | |
538 | "tpcV.="<<tpcVertex<< // track vertex | |
539 | "tzero.="<<tzero<< // tzero info | |
540 | "\n"; | |
541 | } | |
542 | delete [] rejectTrack; | |
543 | } | |
544 | ||
545 | ||
546 | Bool_t AliTPCcalibMaterial::CheckLooper(Int_t index, AliESDEvent *event){ | |
547 | // | |
548 | // check if given track is looper candidate | |
549 | // if looper return kTRUE | |
550 | // | |
551 | Int_t ntracks=event->GetNumberOfTracks(); | |
552 | Int_t index1=-1; | |
553 | const Double_t ktglCut=0.03; | |
554 | const Double_t kalphaCut=0.4; | |
555 | static Int_t counter=0; | |
556 | // | |
557 | AliESDtrack * track0 = event->GetTrack(index); | |
558 | AliESDtrack * track1P = 0; | |
559 | for (Int_t itrack1=0; itrack1<ntracks; itrack1++){ | |
560 | if (itrack1==index) continue; | |
561 | AliESDtrack *track1=event->GetTrack(itrack1); | |
562 | if (!track1) continue; | |
563 | if (TMath::Abs(TMath::Abs(track1->GetTgl())-TMath::Abs(track0->GetTgl()))>ktglCut) continue; | |
564 | if (TMath::Abs(TMath::Abs(track1->GetAlpha())-TMath::Abs(track0->GetAlpha()))>kalphaCut) continue; | |
565 | index1=index; | |
566 | track1P=track1; | |
567 | } | |
568 | if (index1>=0){ | |
569 | TTreeSRedirector *pcstream = GetDebugStreamer(); | |
570 | if (pcstream &&counter<100000){ | |
571 | counter++; | |
572 | AliExternalTrackParam p0(*track0); | |
573 | AliExternalTrackParam p1(*track1P); | |
574 | (*pcstream)<<"checkLooper"<< | |
575 | "p0.="<<&p0<< | |
576 | "p1.="<<&p1<< | |
577 | "\n"; | |
578 | } | |
579 | return kTRUE; | |
580 | } | |
581 | return kFALSE; | |
582 | } | |
583 | ||
c57f81bf | 584 | Bool_t AliTPCcalibMaterial::CheckV0(Int_t index, AliESDEvent *event){ |
4af75575 | 585 | // |
586 | // check if given track is V0 candidata | |
587 | // if looper return kTRUE | |
588 | // | |
589 | return kFALSE; | |
590 | Int_t ntracks=event->GetNumberOfTracks(); | |
591 | Int_t index1=-1; | |
592 | const Double_t kSigmaMass=0.001; | |
593 | const Int_t kChi2Cut=10; | |
594 | static Int_t counter=0; | |
595 | // | |
596 | AliESDtrack * track0 = event->GetTrack(index); | |
597 | AliExternalTrackParam pL(*track0); | |
598 | AliKFParticle part0El(*track0, 11); //assuming mass e | |
599 | AliKFParticle part0Pi(*track0, 211); //assuming mass pi0 | |
600 | AliKFParticle part0P(*track0, 2212); //assuming mass proton | |
601 | if (track0->Charge()*part0El.Q()<0) { | |
602 | part0El.Q()*=-1; // change sign if opposite | |
603 | part0Pi.Q()*=-1; // change sign if opposite | |
604 | part0P.Q()*=-1; // change sign if opposite | |
605 | } | |
606 | Bool_t isGamma=0; | |
607 | Bool_t isK0=0; | |
608 | ||
609 | for (Int_t itrack1=0; itrack1<ntracks; itrack1++){ | |
610 | if (itrack1==index) continue; | |
611 | AliESDtrack *track1=event->GetTrack(itrack1); | |
612 | if (!track1) continue; | |
613 | if (track1->Charge()*track0->Charge()>0) continue; | |
614 | AliKFParticle part1El(*track1, 11); //assuming mass e | |
615 | AliKFParticle part1Pi(*track1, 211); //assuming mass e | |
616 | AliKFParticle part1P(*track1, 2212); //assuming mass e | |
617 | if (track1->Charge()*part1El.Q()<0) { | |
618 | part1El.Q()*=-1; // change sign if opposite | |
619 | part1Pi.Q()*=-1; // change sign if opposite | |
620 | part1P.Q()*=-1; // change sign if opposite | |
621 | } | |
622 | // | |
623 | AliKFVertex vertexG; // gamma conversion candidate | |
624 | vertexG+=part0El; | |
625 | vertexG+=part1El; | |
626 | AliKFVertex vertexGC; // gamma conversion candidate | |
627 | vertexGC+=part0El; | |
628 | vertexGC+=part1El; | |
629 | vertexGC.SetMassConstraint(0,kSigmaMass); | |
630 | AliKFVertex vertexK0; // gamma conversion candidate | |
631 | vertexK0+=part0Pi; | |
632 | vertexK0+=part1Pi; | |
633 | AliKFVertex vertexK0C; // gamma conversion candidate | |
634 | vertexK0C+=part0Pi; | |
635 | vertexK0C+=part1Pi; | |
636 | vertexK0C.SetMassConstraint(0.497614,kSigmaMass); | |
637 | if (vertexGC.GetChi2()<kChi2Cut && vertexG.GetMass()<0.06) isGamma=kTRUE; | |
638 | if (vertexK0C.GetChi2()<kChi2Cut&&TMath::Abs(vertexK0.GetMass()-0.5)<0.06) isK0=kTRUE; | |
639 | if (isGamma||isK0) { | |
640 | index1=index; | |
641 | TTreeSRedirector *pcstream = GetDebugStreamer(); | |
642 | if (pcstream&&counter<2000){ | |
643 | counter++; | |
644 | AliExternalTrackParam p0(*track0); | |
645 | AliExternalTrackParam p1(*track1); | |
646 | (*pcstream)<<"checkV0"<< | |
647 | "p0.="<<&p0<< //particle 0 | |
648 | "p1.="<<&p1<< //particle 1 | |
649 | "isGamma="<<isGamma<< // is gamma candidate | |
650 | "isK0s="<<isK0<< // is k0 candidate | |
651 | "vG.="<<&vertexG<< // is gamma candidate | |
652 | "vGC.="<<&vertexGC<< // is gamma candidate | |
653 | "vK.="<<&vertexK0<< // is K0s candidate | |
654 | "vKC.="<<&vertexK0C<< // is K0s candidate | |
655 | "\n"; | |
656 | } | |
657 | break; | |
658 | } | |
659 | } | |
660 | if (index1>0) return kTRUE; | |
661 | return kFALSE; | |
662 | } | |
663 | ||
664 | /* | |
665 | //AliXRDPROOFtoolkit::FilterList("mater.list","* mapTPC",1); | |
666 | AliXRDPROOFtoolkit toolkit; | |
667 | TChain *chain = toolkit.MakeChainRandom("mater.list.Good","mapVertex",0,4000); | |
668 | TChain *chainTPC = toolkit.MakeChainRandom("mater.list.Good","mapTPC",0,50000); | |
669 | ||
670 | ||
671 | TCut cutErr="sqrt(v.fChi2)*err3D<1.0"; | |
672 | TCut cutOccu="sqrt(v.fChi2*(errX^2+errY^2))/min(sqrt(v.fP[0]^2+v.fP[1]^2+v.fP[2]^2/20),20)<0.2&&sqrt((errX^2+errY^2))/min(sqrt(v.fP[0]^2+v.fP[1]^2+v.fP[2]^2/20),20)<0.3&&sqrt(v.fChi2)*err3D<1.5"; | |
673 | // | |
674 | chainTPC->Draw(">>listTPC",cutOccu,"entryList"); | |
675 | TEntryList *elistTPC = (TEntryList*)gDirectory->Get("listTPC"); | |
676 | chainTPC->SetEntryList(elistTPC); | |
677 | ||
678 | // | |
679 | // | |
680 | ||
681 | chainTPC.Draw("v.fP[1]:v.fP[0]","sqrt(v.fP[0]^2+v.fP[1]^2)<100","",10000000); | |
682 | ||
683 | TCut cutITS="abs(v.fP[2])-10<sqrt(v.fP[0]^2+v.fP[1]^2)"; | |
684 | chainTPC.Draw("v.fP[1]:v.fP[0]",cutITS+"(sqrt(v.fP[0]^2+v.fP[1]^2)<50)&&err2D<1&&err3D/dvertex<0.05","",500000); | |
685 | ||
686 | ||
687 | chainTPC.Draw("v.fP[1]:v.fP[0]>>his(400,-100,100,400,-100,100)","err2D<1&&err2D/dvertex<0.2","colz",50000); | |
688 | ||
689 | ||
690 | chainTPC.Draw("radius:vz>>his(400,-100,100,400,0,100)","err2D<1&&err2D/dvertex<0.05","colz",5000); | |
691 | ||
692 | TTree * tree = chain->CloneTree(); | |
693 | TFile f("material.root","recreate"); | |
694 | tree->Write(); | |
695 | f.Close(); | |
696 | ||
697 | ||
698 | ||
699 | TCut cutChi2="v.fChi2<2&&sqrt(errX^2+errY^2)<2."; | |
700 | TFile f("material.root") | |
701 | TTree * tree = (TTree*)f.Get("mapVertex"); | |
702 | ||
703 | tree->Draw("v.fChi2",cutChi2,"",1000) | |
704 | ||
705 | ||
706 | ||
707 | ||
708 | */ | |
709 | ||
710 |