]>
Commit | Line | Data |
---|---|---|
70c95f95 | 1 | /************************************************************************** |
2 | * Copyright(c) 2006-2008, ALICE Experiment at CERN, All rights reserved. * | |
3 | * * | |
4 | * Author: The ALICE Off-line Project. * | |
5 | * Contributors are mentioned in the code where appropriate. * | |
6 | * * | |
7 | * Permission to use, copy, modify and distribute this software and its * | |
8 | * documentation strictly for non-commercial purposes is hereby granted * | |
9 | * without fee, provided that the above copyright notice appears in all * | |
10 | * copies and that both the copyright notice and this permission notice * | |
11 | * appear in the supporting documentation. The authors make no claims * | |
12 | * about the suitability of this software for any purpose. It is * | |
13 | * provided "as is" without express or implied warranty. * | |
14 | **************************************************************************/ | |
15 | #include <AliESDVertex.h> | |
16 | #include <AliITSVertexer3D.h> | |
17 | #include <AliStrLine.h> | |
18 | #include <AliVertexerTracks.h> | |
19 | #include <Riostream.h> | |
20 | #include <TH3F.h> | |
21 | #include <TTree.h> | |
22 | #include<TClonesArray.h> | |
cd706e57 | 23 | #include "AliLog.h" |
70c95f95 | 24 | #include "AliRunLoader.h" |
25 | #include "AliITSLoader.h" | |
70c95f95 | 26 | #include "AliITSDetTypeRec.h" |
27 | #include "AliITSRecPoint.h" | |
32e63e47 | 28 | #include "AliITSgeomTGeo.h" |
70c95f95 | 29 | ///////////////////////////////////////////////////////////////// |
30 | // this class implements a method to determine | |
31 | // the 3 coordinates of the primary vertex | |
32 | // for p-p collisions | |
33 | // It can be used successfully with Pb-Pb collisions | |
34 | //////////////////////////////////////////////////////////////// | |
35 | ||
36 | ClassImp(AliITSVertexer3D) | |
37 | ||
38 | //______________________________________________________________________ | |
39 | AliITSVertexer3D::AliITSVertexer3D():AliITSVertexer(), | |
40 | fLines(), | |
41 | fVert3D(), | |
42 | fCoarseDiffPhiCut(0.), | |
05d1294c | 43 | fCoarseMaxRCut(0.), |
70c95f95 | 44 | fMaxRCut(0.), |
45 | fZCutDiamond(0.), | |
05d1294c | 46 | fMaxZCut(0.), |
70c95f95 | 47 | fDCAcut(0.), |
27167524 | 48 | fDiffPhiMax(0.) |
49 | { | |
70c95f95 | 50 | // Default constructor |
51 | SetCoarseDiffPhiCut(); | |
05d1294c | 52 | SetCoarseMaxRCut(); |
70c95f95 | 53 | SetMaxRCut(); |
54 | SetZCutDiamond(); | |
05d1294c | 55 | SetMaxZCut(); |
70c95f95 | 56 | SetDCAcut(); |
57 | SetDiffPhiMax(); | |
70c95f95 | 58 | } |
59 | ||
60 | //______________________________________________________________________ | |
61 | AliITSVertexer3D::AliITSVertexer3D(TString fn): AliITSVertexer(fn), | |
62 | fLines(), | |
63 | fVert3D(), | |
05d1294c | 64 | fCoarseDiffPhiCut(0.), |
65 | fCoarseMaxRCut(0.), | |
70c95f95 | 66 | fMaxRCut(0.), |
67 | fZCutDiamond(0.), | |
05d1294c | 68 | fMaxZCut(0.), |
70c95f95 | 69 | fDCAcut(0.), |
27167524 | 70 | fDiffPhiMax(0.) |
71 | { | |
70c95f95 | 72 | // Standard constructor |
73 | fLines = new TClonesArray("AliStrLine",1000); | |
74 | SetCoarseDiffPhiCut(); | |
05d1294c | 75 | SetCoarseMaxRCut(); |
70c95f95 | 76 | SetMaxRCut(); |
77 | SetZCutDiamond(); | |
05d1294c | 78 | SetMaxZCut(); |
70c95f95 | 79 | SetDCAcut(); |
80 | SetDiffPhiMax(); | |
81 | } | |
82 | ||
83 | //______________________________________________________________________ | |
84 | AliITSVertexer3D::~AliITSVertexer3D() { | |
85 | // Destructor | |
27167524 | 86 | if(fLines){ |
70c95f95 | 87 | fLines->Delete(); |
88 | delete fLines; | |
89 | } | |
90 | } | |
91 | ||
27167524 | 92 | //______________________________________________________________________ |
93 | void AliITSVertexer3D::ResetVert3D(){ | |
94 | // | |
95 | fVert3D.SetXv(0.); | |
96 | fVert3D.SetYv(0.); | |
97 | fVert3D.SetZv(0.); | |
98 | fVert3D.SetDispersion(0.); | |
99 | fVert3D.SetNContributors(0); | |
100 | } | |
70c95f95 | 101 | //______________________________________________________________________ |
102 | AliESDVertex* AliITSVertexer3D::FindVertexForCurrentEvent(Int_t evnumber){ | |
103 | // Defines the AliESDVertex for the current event | |
27167524 | 104 | ResetVert3D(); |
cd706e57 | 105 | AliDebug(1,Form("FindVertexForCurrentEvent - 3D - PROCESSING EVENT %d",evnumber)); |
70c95f95 | 106 | if(fLines)fLines->Clear(); |
107 | ||
108 | Int_t nolines = FindTracklets(evnumber,0); | |
109 | fCurrentVertex = 0; | |
110 | if(nolines<2)return fCurrentVertex; | |
05d1294c | 111 | Int_t rc=Prepare3DVertex(0); |
27167524 | 112 | if(rc==0) fVert3D=AliVertexerTracks::TrackletVertexFinder(fLines,0); |
113 | /* uncomment to debug | |
114 | printf("Vertex found in first iteration:\n"); | |
115 | fVert3D.Print(); | |
116 | printf("Start second iteration\n"); | |
117 | end of debug lines */ | |
118 | if(fVert3D.GetNContributors()>0){ | |
70c95f95 | 119 | if(fLines) fLines->Delete(); |
120 | nolines = FindTracklets(evnumber,1); | |
d35ad08f | 121 | if(nolines>=2){ |
05d1294c | 122 | rc=Prepare3DVertex(1); |
27167524 | 123 | if(rc==0) fVert3D=AliVertexerTracks::TrackletVertexFinder(fLines,0); |
d35ad08f | 124 | } |
70c95f95 | 125 | } |
27167524 | 126 | /* uncomment to debug |
127 | printf("Vertex found in second iteration:\n"); | |
128 | fVert3D.Print(); | |
129 | end of debug lines */ | |
130 | ||
131 | Float_t vRadius=TMath::Sqrt(fVert3D.GetXv()*fVert3D.GetXv()+fVert3D.GetYv()*fVert3D.GetYv()); | |
132 | if(vRadius<GetPipeRadius() && fVert3D.GetNContributors()>0){ | |
70c95f95 | 133 | fCurrentVertex = new AliESDVertex(); |
134 | fCurrentVertex->SetTitle("vertexer: 3D"); | |
135 | fCurrentVertex->SetName("Vertex"); | |
27167524 | 136 | fCurrentVertex->SetXv(fVert3D.GetXv()); |
137 | fCurrentVertex->SetYv(fVert3D.GetYv()); | |
138 | fCurrentVertex->SetZv(fVert3D.GetZv()); | |
139 | fCurrentVertex->SetDispersion(fVert3D.GetDispersion()); | |
140 | fCurrentVertex->SetNContributors(fVert3D.GetNContributors()); | |
70c95f95 | 141 | } |
142 | FindMultiplicity(evnumber); | |
143 | return fCurrentVertex; | |
144 | } | |
145 | ||
146 | //______________________________________________________________________ | |
147 | Int_t AliITSVertexer3D::FindTracklets(Int_t evnumber, Int_t optCuts){ | |
148 | // All the possible combinations between recpoints on layer 1and 2 are | |
149 | // considered. Straight lines (=tracklets)are formed. | |
150 | // The tracklets are processed in Prepare3DVertex | |
151 | AliRunLoader *rl =AliRunLoader::GetRunLoader(); | |
152 | AliITSLoader* itsLoader = (AliITSLoader*)rl->GetLoader("ITSLoader"); | |
32e63e47 | 153 | // AliITSgeom* geom = itsLoader->GetITSgeom(); |
70c95f95 | 154 | itsLoader->LoadRecPoints(); |
155 | rl->GetEvent(evnumber); | |
156 | ||
157 | AliITSDetTypeRec detTypeRec; | |
158 | ||
159 | TTree *tR = itsLoader->TreeR(); | |
160 | detTypeRec.SetTreeAddressR(tR); | |
161 | TClonesArray *itsRec = 0; | |
162 | // lc and gc are local and global coordinates for layer 1 | |
163 | Float_t lc[3]; for(Int_t ii=0; ii<3; ii++) lc[ii]=0.; | |
164 | Float_t gc[3]; for(Int_t ii=0; ii<3; ii++) gc[ii]=0.; | |
165 | // lc2 and gc2 are local and global coordinates for layer 2 | |
166 | Float_t lc2[3]; for(Int_t ii=0; ii<3; ii++) lc2[ii]=0.; | |
167 | Float_t gc2[3]; for(Int_t ii=0; ii<3; ii++) gc2[ii]=0.; | |
168 | ||
169 | itsRec = detTypeRec.RecPoints(); | |
170 | TBranch *branch; | |
171 | branch = tR->GetBranch("ITSRecPoints"); | |
172 | ||
173 | // Set values for cuts | |
174 | Float_t xbeam=0., ybeam=0.; | |
05d1294c | 175 | Float_t zvert=0.; |
70c95f95 | 176 | Float_t deltaPhi=fCoarseDiffPhiCut; |
05d1294c | 177 | Float_t deltaR=fCoarseMaxRCut; |
178 | Float_t dZmax=fZCutDiamond; | |
70c95f95 | 179 | if(optCuts){ |
27167524 | 180 | xbeam=fVert3D.GetXv(); |
181 | ybeam=fVert3D.GetYv(); | |
182 | zvert=fVert3D.GetZv(); | |
70c95f95 | 183 | deltaPhi = fDiffPhiMax; |
05d1294c | 184 | deltaR=fMaxRCut; |
185 | dZmax=fMaxZCut; | |
70c95f95 | 186 | } |
70c95f95 | 187 | Int_t nrpL1 = 0; // number of rec points on layer 1 |
188 | Int_t nrpL2 = 0; // number of rec points on layer 2 | |
189 | ||
190 | // By default irstL1=0 and lastL1=79 | |
32e63e47 | 191 | Int_t irstL1 = AliITSgeomTGeo::GetModuleIndex(1,1,1); |
192 | Int_t lastL1 = AliITSgeomTGeo::GetModuleIndex(2,1,1)-1; | |
70c95f95 | 193 | for(Int_t module= irstL1; module<=lastL1;module++){ // count number of recopints on layer 1 |
194 | branch->GetEvent(module); | |
195 | nrpL1+= itsRec->GetEntries(); | |
196 | detTypeRec.ResetRecPoints(); | |
197 | } | |
198 | //By default irstL2=80 and lastL2=239 | |
32e63e47 | 199 | Int_t irstL2 = AliITSgeomTGeo::GetModuleIndex(2,1,1); |
200 | Int_t lastL2 = AliITSgeomTGeo::GetModuleIndex(3,1,1)-1; | |
70c95f95 | 201 | for(Int_t module= irstL2; module<=lastL2;module++){ // count number of recopints on layer 2 |
202 | branch->GetEvent(module); | |
203 | nrpL2+= itsRec->GetEntries(); | |
204 | detTypeRec.ResetRecPoints(); | |
205 | } | |
206 | if(nrpL1 == 0 || nrpL2 == 0){ | |
207 | itsLoader->UnloadRecPoints(); | |
208 | return -1; | |
209 | } | |
cd706e57 | 210 | AliDebug(1,Form("RecPoints on Layer 1,2 = %d, %d\n",nrpL1,nrpL2)); |
70c95f95 | 211 | |
212 | Double_t a[3]={xbeam,ybeam,0.}; | |
213 | Double_t b[3]={xbeam,ybeam,10.}; | |
214 | AliStrLine zeta(a,b,kTRUE); | |
215 | ||
216 | Int_t nolines = 0; | |
217 | // Loop on modules of layer 1 | |
218 | for(Int_t modul1= irstL1; modul1<=lastL1;modul1++){ // Loop on modules of layer 1 | |
27167524 | 219 | UShort_t ladder=int(modul1/4)+1; // ladders are numbered starting from 1 |
70c95f95 | 220 | branch->GetEvent(modul1); |
221 | Int_t nrecp1 = itsRec->GetEntries(); | |
222 | TClonesArray *prpl1 = new TClonesArray("AliITSRecPoint",nrecp1); | |
223 | prpl1->SetOwner(); | |
224 | TClonesArray &rpl1 = *prpl1; | |
225 | for(Int_t j=0;j<nrecp1;j++){ | |
226 | AliITSRecPoint *recp = (AliITSRecPoint*)itsRec->At(j); | |
227 | new(rpl1[j])AliITSRecPoint(*recp); | |
228 | } | |
229 | detTypeRec.ResetRecPoints(); | |
230 | for(Int_t j=0;j<nrecp1;j++){ | |
231 | AliITSRecPoint *recp = (AliITSRecPoint*)prpl1->At(j); | |
232 | // Local coordinates of this recpoint | |
32e63e47 | 233 | /* |
70c95f95 | 234 | lc[0]=recp->GetDetLocalX(); |
235 | lc[2]=recp->GetDetLocalZ(); | |
236 | geom->LtoG(modul1,lc,gc); // global coordinates | |
32e63e47 | 237 | */ |
238 | recp->GetGlobalXYZ(gc); | |
70c95f95 | 239 | Double_t phi1 = TMath::ATan2(gc[1]-ybeam,gc[0]-xbeam); |
240 | if(phi1<0)phi1=2*TMath::Pi()+phi1; | |
27167524 | 241 | for(Int_t ladl2=0 ; ladl2<fLadOnLay2*2+1;ladl2++){ |
242 | for(Int_t k=0;k<4;k++){ | |
243 | Int_t ladmod=fLadders[ladder-1]+ladl2; | |
32e63e47 | 244 | if(ladmod>AliITSgeomTGeo::GetNLadders(2)) ladmod=ladmod-AliITSgeomTGeo::GetNLadders(2); |
245 | Int_t modul2=AliITSgeomTGeo::GetModuleIndex(2,ladmod,k+1); | |
27167524 | 246 | branch->GetEvent(modul2); |
247 | Int_t nrecp2 = itsRec->GetEntries(); | |
248 | for(Int_t j2=0;j2<nrecp2;j2++){ | |
249 | recp = (AliITSRecPoint*)itsRec->At(j2); | |
32e63e47 | 250 | /* |
27167524 | 251 | lc2[0]=recp->GetDetLocalX(); |
252 | lc2[2]=recp->GetDetLocalZ(); | |
253 | geom->LtoG(modul2,lc2,gc2); | |
32e63e47 | 254 | */ |
255 | recp->GetGlobalXYZ(gc2); | |
27167524 | 256 | Double_t phi2 = TMath::ATan2(gc2[1]-ybeam,gc2[0]-xbeam); |
257 | if(phi2<0)phi2=2*TMath::Pi()+phi2; | |
258 | Double_t diff = TMath::Abs(phi2-phi1); | |
259 | if(diff>TMath::Pi())diff=2.*TMath::Pi()-diff; | |
260 | if(diff>deltaPhi)continue; | |
261 | AliStrLine line(gc,gc2,kTRUE); | |
262 | Double_t cp[3]; | |
263 | Int_t retcode = line.Cross(&zeta,cp); | |
264 | if(retcode<0)continue; | |
265 | Double_t dca = line.GetDCA(&zeta); | |
266 | if(dca<0.) continue; | |
267 | if(dca>deltaR)continue; | |
268 | Double_t deltaZ=cp[2]-zvert; | |
269 | if(TMath::Abs(deltaZ)>dZmax)continue; | |
270 | MakeTracklet(gc,gc2,nolines); | |
271 | } | |
272 | detTypeRec.ResetRecPoints(); | |
70c95f95 | 273 | } |
70c95f95 | 274 | } |
275 | } | |
276 | delete prpl1; | |
277 | } | |
278 | if(nolines == 0)return -2; | |
279 | return nolines; | |
280 | } | |
281 | ||
282 | //______________________________________________________________________ | |
283 | void AliITSVertexer3D::FindVertices(){ | |
284 | // computes the vertices of the events in the range FirstEvent - LastEvent | |
285 | AliRunLoader *rl = AliRunLoader::GetRunLoader(); | |
286 | AliITSLoader* itsLoader = (AliITSLoader*) rl->GetLoader("ITSLoader"); | |
287 | itsLoader->ReloadRecPoints(); | |
288 | for(Int_t i=fFirstEvent;i<=fLastEvent;i++){ | |
289 | rl->GetEvent(i); | |
290 | FindVertexForCurrentEvent(i); | |
291 | if(fCurrentVertex){ | |
292 | WriteCurrentVertex(); | |
293 | } | |
294 | } | |
295 | } | |
296 | ||
297 | ||
70c95f95 | 298 | //______________________________________________________________________ |
05d1294c | 299 | Int_t AliITSVertexer3D::Prepare3DVertex(Int_t optCuts){ |
70c95f95 | 300 | // Finds the 3D vertex information using tracklets |
301 | Int_t retcode = -1; | |
302 | ||
05d1294c | 303 | Float_t xbeam=0.; |
304 | Float_t ybeam=0.; | |
305 | Float_t zvert=0.; | |
306 | Float_t deltaR=fCoarseMaxRCut; | |
307 | Float_t dZmax=fZCutDiamond; | |
308 | if(optCuts){ | |
27167524 | 309 | xbeam=fVert3D.GetXv(); |
310 | ybeam=fVert3D.GetYv(); | |
311 | zvert=fVert3D.GetZv(); | |
05d1294c | 312 | deltaR=fMaxRCut; |
313 | dZmax=fMaxZCut; | |
314 | } | |
315 | ||
70c95f95 | 316 | Int_t nbr=50; |
05d1294c | 317 | Float_t rl=-fCoarseMaxRCut; |
318 | Float_t rh=fCoarseMaxRCut; | |
70c95f95 | 319 | Int_t nbz=100; |
320 | Float_t zl=-fZCutDiamond; | |
321 | Float_t zh=fZCutDiamond; | |
322 | Float_t binsizer=(rh-rl)/nbr; | |
323 | Float_t binsizez=(zh-zl)/nbz; | |
324 | TH3F *h3d = new TH3F("h3d","xyz distribution",nbr,rl,rh,nbr,rl,rh,nbz,zl,zh); | |
325 | ||
70c95f95 | 326 | // cleanup of the TCLonesArray of tracklets (i.e. fakes are removed) |
327 | Int_t *validate = new Int_t [fLines->GetEntriesFast()]; | |
328 | for(Int_t i=0; i<fLines->GetEntriesFast();i++)validate[i]=0; | |
70c95f95 | 329 | for(Int_t i=0; i<fLines->GetEntriesFast()-1;i++){ |
05d1294c | 330 | if(validate[i]==1)continue; |
70c95f95 | 331 | AliStrLine *l1 = (AliStrLine*)fLines->At(i); |
332 | for(Int_t j=i+1;j<fLines->GetEntriesFast();j++){ | |
333 | AliStrLine *l2 = (AliStrLine*)fLines->At(j); | |
05d1294c | 334 | Float_t dca=l1->GetDCA(l2); |
335 | if(dca > fDCAcut || dca<0.00001) continue; | |
70c95f95 | 336 | Double_t point[3]; |
337 | Int_t retc = l1->Cross(l2,point); | |
338 | if(retc<0)continue; | |
339 | Double_t rad=TMath::Sqrt(point[0]*point[0]+point[1]*point[1]); | |
05d1294c | 340 | if(rad>fCoarseMaxRCut)continue; |
341 | Double_t deltaX=point[0]-xbeam; | |
342 | Double_t deltaY=point[1]-ybeam; | |
343 | Double_t deltaZ=point[2]-zvert; | |
344 | Double_t raddist=TMath::Sqrt(deltaX*deltaX+deltaY*deltaY); | |
345 | if(TMath::Abs(deltaZ)>dZmax)continue; | |
346 | if(raddist>deltaR)continue; | |
27167524 | 347 | validate[i]=1; |
348 | validate[j]=1; | |
70c95f95 | 349 | h3d->Fill(point[0],point[1],point[2]); |
350 | } | |
351 | } | |
352 | ||
353 | ||
354 | ||
70c95f95 | 355 | Int_t numbtracklets=0; |
356 | for(Int_t i=0; i<fLines->GetEntriesFast();i++)if(validate[i]>=1)numbtracklets++; | |
357 | if(numbtracklets<2){delete [] validate; delete h3d; return retcode; } | |
358 | ||
359 | for(Int_t i=0; i<fLines->GetEntriesFast();i++){ | |
360 | if(validate[i]<1)fLines->RemoveAt(i); | |
361 | } | |
362 | fLines->Compress(); | |
cd706e57 | 363 | AliDebug(1,Form("Number of tracklets (after compress)%d ",fLines->GetEntriesFast())); |
70c95f95 | 364 | delete [] validate; |
365 | ||
05d1294c | 366 | |
367 | // finds peak in histo | |
70c95f95 | 368 | TAxis *xax = h3d->GetXaxis(); |
369 | TAxis *yax = h3d->GetYaxis(); | |
370 | TAxis *zax = h3d->GetZaxis(); | |
70c95f95 | 371 | Double_t peak[3]={0.,0.,0.}; |
372 | Float_t contref = 0.; | |
373 | for(Int_t i=xax->GetFirst();i<=xax->GetLast();i++){ | |
374 | Float_t xval = xax->GetBinCenter(i); | |
375 | for(Int_t j=yax->GetFirst();j<=yax->GetLast();j++){ | |
376 | Float_t yval = yax->GetBinCenter(j); | |
377 | for(Int_t k=zax->GetFirst();k<=zax->GetLast();k++){ | |
378 | Float_t bc = h3d->GetBinContent(i,j,k); | |
379 | Float_t zval = zax->GetBinCenter(k); | |
380 | if(bc>contref){ | |
381 | contref = bc; | |
382 | peak[2] = zval; | |
383 | peak[1] = yval; | |
384 | peak[0] = xval; | |
385 | } | |
386 | } | |
387 | } | |
388 | } | |
389 | delete h3d; | |
390 | ||
391 | // Second selection loop | |
392 | Float_t bs=(binsizer+binsizez)/2.; | |
05d1294c | 393 | for(Int_t i=0; i<fLines->GetEntriesFast();i++){ |
70c95f95 | 394 | AliStrLine *l1 = (AliStrLine*)fLines->At(i); |
395 | if(l1->GetDistFromPoint(peak)>2.5*bs)fLines->RemoveAt(i); | |
396 | } | |
397 | fLines->Compress(); | |
cd706e57 | 398 | AliDebug(1,Form("Number of tracklets (after 2nd compression) %d",fLines->GetEntriesFast())); |
70c95f95 | 399 | |
70c95f95 | 400 | if(fLines->GetEntriesFast()>1){ |
27167524 | 401 | // find a first candidate for the primary vertex |
402 | fVert3D=AliVertexerTracks::TrackletVertexFinder(fLines,0); | |
70c95f95 | 403 | // make a further selection on tracklets based on this first candidate |
27167524 | 404 | fVert3D.GetXYZ(peak); |
405 | AliDebug(1,Form("FIRST V candidate: %f ; %f ; %f",peak[0],peak[1],peak[2])); | |
05d1294c | 406 | for(Int_t i=0; i<fLines->GetEntriesFast();i++){ |
70c95f95 | 407 | AliStrLine *l1 = (AliStrLine*)fLines->At(i); |
408 | if(l1->GetDistFromPoint(peak)> fDCAcut)fLines->RemoveAt(i); | |
409 | } | |
410 | fLines->Compress(); | |
cd706e57 | 411 | AliDebug(1,Form("Number of tracklets (after 3rd compression) %d",fLines->GetEntriesFast())); |
27167524 | 412 | if(fLines->GetEntriesFast()>1) retcode=0; // this new tracklet selection is used |
413 | else retcode =1; // the previous tracklet selection will be used | |
70c95f95 | 414 | } |
415 | else { | |
416 | retcode = 0; | |
417 | } | |
418 | return retcode; | |
419 | } | |
420 | ||
421 | //______________________________________________________________________ | |
422 | void AliITSVertexer3D::MakeTracklet(Double_t *pA, Double_t *pB, Int_t &nolines) { | |
423 | // makes a tracklet | |
424 | TClonesArray &lines = *fLines; | |
425 | if(nolines == 0){ | |
426 | if(fLines->GetEntriesFast()>0)fLines->Clear(); | |
427 | } | |
428 | if(fLines->GetEntriesFast()==fLines->GetSize()){ | |
429 | Int_t newsize=(Int_t) 1.5*fLines->GetEntriesFast(); | |
430 | fLines->Expand(newsize); | |
431 | } | |
432 | ||
433 | new(lines[nolines++])AliStrLine(pA,pB,kTRUE); | |
434 | } | |
435 | ||
436 | //______________________________________________________________________ | |
437 | void AliITSVertexer3D::MakeTracklet(Float_t *pA, Float_t *pB, Int_t &nolines) {// Makes a tracklet | |
438 | // | |
439 | Double_t a[3],b[3]; | |
440 | for(Int_t i=0;i<3;i++){ | |
441 | a[i] = pA[i]; | |
442 | b[i] = pB[i]; | |
443 | } | |
444 | MakeTracklet(a,b,nolines); | |
445 | } | |
446 | ||
447 | //________________________________________________________ | |
448 | void AliITSVertexer3D::PrintStatus() const { | |
449 | // Print current status | |
450 | cout <<"=======================================================\n"; | |
451 | cout << "Loose cut on Delta Phi "<<fCoarseDiffPhiCut<<endl; | |
05d1294c | 452 | cout << "Cut on tracklet DCA to Z axis "<<fCoarseMaxRCut<<endl; |
453 | cout << "Cut on tracklet DCA to beam axis "<<fMaxRCut<<endl; | |
70c95f95 | 454 | cout << "Cut on diamond (Z) "<<fZCutDiamond<<endl; |
455 | cout << "Cut on DCA - tracklet to tracklet and to vertex "<<fDCAcut<<endl; | |
456 | cout <<" Max Phi difference: "<<fDiffPhiMax<<endl; | |
457 | ||
458 | ||
459 | cout <<"=======================================================\n"; | |
460 | } |