]>
Commit | Line | Data |
---|---|---|
973a91f8 | 1 | ///////////////////////////////////////////////////////////////////////////// |
2 | // // | |
3 | // AliFemtoShareQualityPairCut - a pair cut which checks for some pair // | |
4 | // qualities that attempt to identify slit/doubly reconstructed tracks // | |
5 | // // | |
6 | ///////////////////////////////////////////////////////////////////////////// | |
7 | /*************************************************************************** | |
8 | * | |
9 | * $Id: AliFemtoShareQualityPairCut.cxx 50722 2011-07-21 15:18:38Z akisiel $ | |
10 | * | |
11 | * Author: Adam Kisiel, Ohio State, kisiel@mps.ohio-state.edu | |
12 | *************************************************************************** | |
13 | * | |
14 | * Description: part of STAR HBT Framework: AliFemtoMaker package | |
15 | * a cut to remove "shared" and "split" pairs | |
16 | * | |
17 | *************************************************************************** | |
18 | * | |
19 | * | |
20 | **************************************************************************/ | |
21 | ||
22 | #include "AliFemtoV0TrackPairCut.h" | |
23 | #include <string> | |
24 | #include <cstdio> | |
25 | ||
26 | #ifdef __ROOT__ | |
27 | ClassImp(AliFemtoV0TrackPairCut) | |
28 | #endif | |
29 | ||
30 | //__________________ | |
31 | AliFemtoV0TrackPairCut::AliFemtoV0TrackPairCut(): | |
32 | fNPairsPassed(0), | |
33 | fNPairsFailed(0), | |
34 | fV0Max(1.0), | |
35 | fShareQualityMax(1.0), | |
36 | fShareFractionMax(1.0), | |
37 | fRemoveSameLabel(0), | |
1ae130d9 | 38 | fTrackTPCOnly(0), |
39 | fDataType(kAOD), | |
40 | fDTPCMin(0), | |
ce7b3d98 | 41 | fDTPCExitMin(0), |
42 | fKstarCut(0), | |
43 | fFirstParticleType(kLambda), | |
6691ea4f | 44 | fSecondParticleType(kProton), |
45 | fMinAvgSepTrackPos(0), | |
46 | fMinAvgSepTrackNeg(0) | |
973a91f8 | 47 | { |
48 | // Default constructor | |
49 | // Nothing to do | |
50 | } | |
51 | //__________________ | |
52 | AliFemtoV0TrackPairCut::~AliFemtoV0TrackPairCut(){ | |
53 | /* no-op */ | |
54 | } | |
1ae130d9 | 55 | AliFemtoV0TrackPairCut& AliFemtoV0TrackPairCut::operator=(const AliFemtoV0TrackPairCut& cut) |
578b62b9 | 56 | { |
57 | if (this != &cut) { | |
58 | ||
59 | AliFemtoPairCut::operator=(cut); | |
60 | fNPairsPassed = 0; | |
61 | fNPairsFailed =0; | |
62 | fV0Max = 1.0; | |
63 | fShareQualityMax = 1.0; | |
64 | fShareFractionMax = 1.0; | |
65 | fRemoveSameLabel = 0; | |
1ae130d9 | 66 | fTrackTPCOnly = 0; |
67 | fDataType = kAOD; | |
68 | fDTPCMin = 0; | |
69 | fDTPCExitMin = 0; | |
6691ea4f | 70 | fMinAvgSepTrackPos = 0; |
71 | fMinAvgSepTrackNeg = 0; | |
578b62b9 | 72 | } |
73 | ||
74 | return *this; | |
75 | } | |
973a91f8 | 76 | //__________________ |
77 | bool AliFemtoV0TrackPairCut::Pass(const AliFemtoPair* pair){ | |
78 | // Check for pairs that are possibly shared/double reconstruction | |
79 | ||
80 | bool temp = true; | |
81 | //Track1 - V0 | |
82 | //Track2 - track | |
83 | ||
973a91f8 | 84 | if(!(pair->Track1()->V0() && pair->Track2()->Track())) |
85 | { | |
86 | return false; | |
87 | } | |
88 | if(fTrackTPCOnly) | |
89 | { | |
90 | if( (-(pair->Track1()->V0()->IdNeg()+1)) ==pair->Track2()->TrackId() || (-(pair->Track1()->V0()->IdPos()+1)) ==pair->Track2()->TrackId()) | |
91 | { | |
92 | return false; | |
93 | } | |
94 | } | |
95 | else | |
96 | { | |
97 | if(pair->Track1()->V0()->IdNeg()==pair->Track2()->TrackId() || pair->Track1()->V0()->IdPos()==pair->Track2()->TrackId()) | |
98 | { | |
99 | return false; | |
100 | } | |
101 | } | |
102 | ||
1ae130d9 | 103 | bool tempTPCEntrancePos = true; |
104 | bool tempTPCEntranceNeg = true; | |
105 | bool tempTPCExitPos = true; | |
106 | bool tempTPCExitNeg = true; | |
107 | if(fDataType==kESD || fDataType==kAOD) | |
108 | { | |
109 | double distx = pair->Track1()->V0()->NominalTpcEntrancePointPos().x() - pair->Track2()->Track()->NominalTpcEntrancePoint().x(); | |
110 | double disty = pair->Track1()->V0()->NominalTpcEntrancePointPos().y() - pair->Track2()->Track()->NominalTpcEntrancePoint().y(); | |
111 | double distz = pair->Track1()->V0()->NominalTpcEntrancePointPos().z() - pair->Track2()->Track()->NominalTpcEntrancePoint().z(); | |
112 | double distPos = sqrt(distx*distx + disty*disty + distz*distz); | |
113 | ||
114 | distx = pair->Track1()->V0()->NominalTpcEntrancePointNeg().x() - pair->Track2()->Track()->NominalTpcEntrancePoint().x(); | |
115 | disty = pair->Track1()->V0()->NominalTpcEntrancePointNeg().y() - pair->Track2()->Track()->NominalTpcEntrancePoint().y(); | |
116 | distz = pair->Track1()->V0()->NominalTpcEntrancePointNeg().z() - pair->Track2()->Track()->NominalTpcEntrancePoint().z(); | |
117 | double distNeg = sqrt(distx*distx + disty*disty + distz*distz); | |
118 | ||
119 | double distExitX = pair->Track1()->V0()->NominalTpcExitPointPos().x() - pair->Track2()->Track()->NominalTpcExitPoint().x(); | |
120 | double distExitY = pair->Track1()->V0()->NominalTpcExitPointPos().y() - pair->Track2()->Track()->NominalTpcExitPoint().y(); | |
121 | double distExitZ = pair->Track1()->V0()->NominalTpcExitPointPos().z() - pair->Track2()->Track()->NominalTpcExitPoint().z(); | |
122 | double distExitPos = sqrt(distExitX*distExitX + distExitY*distExitY + distExitZ*distExitZ); | |
123 | ||
124 | distExitX = pair->Track1()->V0()->NominalTpcExitPointNeg().x() - pair->Track2()->Track()->NominalTpcExitPoint().x(); | |
125 | distExitY = pair->Track1()->V0()->NominalTpcExitPointNeg().y() - pair->Track2()->Track()->NominalTpcExitPoint().y(); | |
126 | distExitZ = pair->Track1()->V0()->NominalTpcExitPointNeg().z() - pair->Track2()->Track()->NominalTpcExitPoint().z(); | |
127 | double distExitNeg = sqrt(distExitX*distExitX + distExitY*distExitY + distExitZ*distExitZ); | |
128 | ||
129 | tempTPCEntrancePos = distPos > fDTPCMin; | |
130 | tempTPCEntranceNeg = distNeg > fDTPCMin; | |
131 | tempTPCExitPos = distExitPos > fDTPCExitMin; | |
132 | tempTPCExitNeg = distExitNeg > fDTPCExitMin; | |
133 | } | |
134 | ||
135 | if(!tempTPCEntrancePos || !tempTPCEntranceNeg || !tempTPCExitPos || !tempTPCExitNeg) return false; | |
973a91f8 | 136 | |
137 | //reject merged trakcs in TPC | |
138 | ||
139 | ||
140 | //temp = dist > fDTPCMin; | |
141 | //koniec kopii | |
142 | //if(!temp) | |
143 | //return false; | |
144 | ||
145 | //kopia z AliFemtoShareQualityPairCut.cxx | |
146 | Int_t nh = 0; | |
147 | Int_t an = 0; | |
148 | Int_t ns = 0; | |
149 | ||
150 | ||
151 | ||
152 | if ((fShareFractionMax < 1.0) && ( fShareQualityMax < 1.0)) { | |
153 | for (unsigned int imap=0; imap<pair->Track1()->V0()->TPCclustersPos().GetNbits(); imap++) { | |
154 | // If both have clusters in the same row | |
155 | if (pair->Track1()->V0()->TPCclustersPos().TestBitNumber(imap) && | |
156 | pair->Track2()->Track()->TPCclusters().TestBitNumber(imap)) { | |
157 | // Do they share it ? | |
158 | if (pair->Track1()->V0()->TPCsharingPos().TestBitNumber(imap) && | |
159 | pair->Track2()->Track()->TPCsharing().TestBitNumber(imap)) | |
160 | { | |
161 | // cout << "A shared cluster !!!" << endl; | |
162 | // cout << "imap idx1 idx2 " | |
163 | // << imap << " " | |
164 | // << tP1idx[imap] << " " << tP2idx[imap] << endl; | |
165 | an++; | |
166 | nh+=2; | |
167 | ns+=2; | |
168 | } | |
169 | ||
170 | // Different hits on the same padrow | |
171 | else { | |
172 | an--; | |
173 | nh+=2; | |
174 | } | |
175 | } | |
176 | else if (pair->Track1()->V0()->TPCclustersPos().TestBitNumber(imap) || | |
177 | pair->Track2()->Track()->TPCclusters().TestBitNumber(imap)) { | |
178 | // One track has a hit, the other does not | |
179 | an++; | |
180 | nh++; | |
181 | } | |
182 | } | |
183 | ||
184 | Float_t hsmval = 0.0; | |
185 | Float_t hsfval = 0.0; | |
186 | ||
187 | if (nh >0) { | |
188 | hsmval = an*1.0/nh; | |
189 | hsfval = ns*1.0/nh; | |
190 | } | |
191 | // if (hsmval > -0.4) { | |
192 | // cout << "Pair quality: " << hsmval << " " << an << " " << nh << " " | |
193 | // << (pair->Track1()->Track()) << " " | |
194 | // << (pair->Track2()->Track()) << endl; | |
195 | // cout << "Bits: " << pair->Track1()->Track()->TPCclusters().GetNbits() << endl; | |
196 | // } | |
197 | // if (hsfval > 0.0) { | |
198 | // cout << "Pair sharity: " << hsfval << " " << ns << " " << nh << " " << hsmval << " " << an << " " << nh << endl; | |
199 | // } | |
200 | ||
201 | temp = (hsmval < fShareQualityMax) && (hsfval < fShareFractionMax); | |
202 | if(!temp) return false; | |
203 | ||
204 | nh = 0; | |
205 | an = 0; | |
206 | ns = 0; | |
207 | ||
208 | for (unsigned int imap=0; imap<pair->Track1()->V0()->TPCclustersNeg().GetNbits(); imap++) { | |
209 | // If both have clusters in the same row | |
210 | if (pair->Track1()->V0()->TPCclustersNeg().TestBitNumber(imap) && | |
211 | pair->Track2()->Track()->TPCclusters().TestBitNumber(imap)) { | |
212 | // Do they share it ? | |
213 | if (pair->Track1()->V0()->TPCsharingNeg().TestBitNumber(imap) && | |
214 | pair->Track2()->Track()->TPCsharing().TestBitNumber(imap)) | |
215 | { | |
216 | // cout << "A shared cluster !!!" << endl; | |
217 | // cout << "imap idx1 idx2 " | |
218 | // << imap << " " | |
219 | // << tP1idx[imap] << " " << tP2idx[imap] << endl; | |
220 | an++; | |
221 | nh+=2; | |
222 | ns+=2; | |
223 | } | |
224 | ||
225 | // Different hits on the same padrow | |
226 | else { | |
227 | an--; | |
228 | nh+=2; | |
229 | } | |
230 | } | |
231 | else if (pair->Track1()->V0()->TPCclustersNeg().TestBitNumber(imap) || | |
232 | pair->Track2()->Track()->TPCclusters().TestBitNumber(imap)) { | |
233 | // One track has a hit, the other does not | |
234 | an++; | |
235 | nh++; | |
236 | } | |
237 | } | |
238 | ||
239 | hsmval = 0.0; | |
240 | hsfval = 0.0; | |
241 | ||
242 | if (nh >0) { | |
243 | hsmval = an*1.0/nh; | |
244 | hsfval = ns*1.0/nh; | |
245 | } | |
246 | // if (hsmval > -0.4) { | |
247 | // cout << "Pair quality: " << hsmval << " " << an << " " << nh << " " | |
248 | // << (pair->Track1()->Track()) << " " | |
249 | // << (pair->Track2()->Track()) << endl; | |
250 | // cout << "Bits: " << pair->Track1()->Track()->TPCclusters().GetNbits() << endl; | |
251 | // } | |
252 | // if (hsfval > 0.0) { | |
253 | // cout << "Pair sharity: " << hsfval << " " << ns << " " << nh << " " << hsmval << " " << an << " " << nh << endl; | |
254 | // } | |
255 | ||
256 | temp = (hsmval < fShareQualityMax) && (hsfval < fShareFractionMax); | |
257 | ||
258 | ||
259 | } | |
260 | else | |
261 | temp = true; | |
262 | //koniec kopii | |
263 | ||
264 | ||
6691ea4f | 265 | //avg sep pair cut |
266 | double avgSep=0; | |
267 | AliFemtoThreeVector first, second, tmp; | |
268 | for(int i=0; i<8 ;i++) | |
269 | { | |
270 | tmp = pair->Track1()->V0()->NominalTpcPointPos(i); | |
271 | //cout<<"X pos: "<<tmp.x()<<endl; | |
272 | first.SetX((double)(tmp.x())); | |
273 | first.SetY((double)tmp.y()); | |
274 | first.SetZ((double)tmp.z()); | |
275 | ||
276 | tmp = pair->Track2()->Track()->NominalTpcPoint(i); | |
277 | second.SetX((double)tmp.x()); | |
278 | second.SetY((double)tmp.y()); | |
279 | second.SetZ((double)tmp.z()); | |
280 | ||
281 | avgSep += TMath::Sqrt(((double)first.x()-(double)second.x())*((double)first.x()-(double)second.x())+((double)first.y()-(double)second.y())*((double)first.y()-second.y())+((double)first.z()-(double)second.z())*((double)first.z()-(double)second.z())); | |
282 | } | |
283 | avgSep /= 8; | |
284 | ||
285 | if(avgSep<fMinAvgSepTrackPos) return false; | |
286 | ||
287 | avgSep = 0; | |
288 | ||
289 | for(int i=0; i<8 ;i++) | |
290 | { | |
291 | tmp = pair->Track1()->V0()->NominalTpcPointNeg(i); | |
292 | //cout<<"X pos: "<<tmp.x()<<endl; | |
293 | first.SetX((double)(tmp.x())); | |
294 | first.SetY((double)tmp.y()); | |
295 | first.SetZ((double)tmp.z()); | |
296 | ||
297 | tmp = pair->Track2()->Track()->NominalTpcPoint(i); | |
298 | second.SetX((double)tmp.x()); | |
299 | second.SetY((double)tmp.y()); | |
300 | second.SetZ((double)tmp.z()); | |
301 | ||
302 | avgSep += TMath::Sqrt(((double)first.x()-(double)second.x())*((double)first.x()-(double)second.x())+((double)first.y()-(double)second.y())*((double)first.y()-second.y())+((double)first.z()-(double)second.z())*((double)first.z()-(double)second.z())); | |
303 | } | |
304 | avgSep /= 8; | |
305 | ||
306 | if(avgSep<fMinAvgSepTrackNeg) return false; | |
307 | ||
308 | ||
309 | ||
310 | ||
ce7b3d98 | 311 | //Qinv cut (we are trying to get rid of antiresidual correlation between primary protons) |
312 | if(fKstarCut > 0) | |
313 | { | |
314 | ||
315 | //2 daughters of first V0 | |
316 | double mom1PosX = pair->Track1()->V0()->MomPosX(); | |
317 | double mom1PosY = pair->Track1()->V0()->MomPosY(); | |
318 | double mom1PosZ = pair->Track1()->V0()->MomPosZ(); | |
319 | double mom1NegX = pair->Track1()->V0()->MomNegX(); | |
320 | double mom1NegY = pair->Track1()->V0()->MomNegY(); | |
321 | double mom1NegZ = pair->Track1()->V0()->MomNegZ(); | |
322 | ||
323 | ||
324 | //double PionMass = 0.13956995; | |
325 | //double KaonMass = 0.493677; | |
326 | double ProtonMass = 0.938272; | |
327 | //double LambdaMass = 1.115683; | |
328 | ||
329 | AliFemtoLorentzVector fFourMomentum1; // Particle momentum | |
330 | AliFemtoLorentzVector fFourMomentum2; // Particle momentum | |
331 | ||
332 | AliFemtoThreeVector temp1; | |
333 | double ener1=0; | |
334 | if(fFirstParticleType == 0) //lambda | |
335 | { | |
336 | if(fSecondParticleType == 2) //proton | |
337 | { | |
338 | //temp1 = ::sqrt(mom1PosX*mom1PosX+mom1PosY*mom1PosY+mom1PosZ*mom1PosZ); | |
339 | temp1.SetX(mom1PosX); | |
340 | temp1.SetY(mom1PosY); | |
341 | temp1.SetZ(mom1PosZ); | |
342 | ener1 = ::sqrt(temp1.Mag2()+ProtonMass*ProtonMass); | |
343 | } | |
344 | } | |
345 | else if(fFirstParticleType == 1) //antilambda | |
346 | { | |
347 | if(fSecondParticleType == 3) //antiproton | |
348 | { | |
349 | //temp1 = ::sqrt(mom1NegX*mom1NegX+mom1NegY*mom1NegY+mom1NegZ*mom1NegZ); | |
350 | temp1.SetX(mom1NegX); | |
351 | temp1.SetY(mom1NegY); | |
352 | temp1.SetZ(mom1NegZ); | |
353 | ener1 = ::sqrt(temp1.Mag2()+ProtonMass*ProtonMass); | |
354 | } | |
355 | } | |
356 | fFourMomentum1.SetVect(temp1); | |
357 | fFourMomentum1.SetE(ener1); | |
358 | ||
359 | //AliFemtoLorentzVector fFourMomentum2; // Particle momentum | |
360 | AliFemtoThreeVector temp2; | |
361 | double ener2=0; | |
362 | ||
363 | //2 daughters of second V0 | |
364 | temp2 = pair->Track2()->Track()->P(); | |
365 | ||
366 | ||
367 | if(fSecondParticleType == 2 || fSecondParticleType == 3) //proton | |
368 | { | |
369 | ener2 = ::sqrt(temp2.Mag2()+ProtonMass*ProtonMass); | |
370 | } | |
371 | ||
372 | fFourMomentum2.SetVect(temp2); | |
373 | fFourMomentum2.SetE(ener2); | |
973a91f8 | 374 | |
ce7b3d98 | 375 | //QInv calculation |
376 | AliFemtoLorentzVector tDiff = (fFourMomentum1-fFourMomentum2); | |
973a91f8 | 377 | |
ce7b3d98 | 378 | double tQinv = fabs(-1.*tDiff.m()); // note - qInv() will be negative for identical pairs... |
379 | ||
380 | //cout<<"tQinv/2: "<<tQinv/2<<endl; | |
381 | if(tQinv/2 < fKstarCut) | |
382 | { | |
973a91f8 | 383 | |
ce7b3d98 | 384 | temp = false; |
385 | } | |
386 | } | |
973a91f8 | 387 | return temp; |
388 | } | |
389 | //__________________ | |
390 | AliFemtoString AliFemtoV0TrackPairCut::Report(){ | |
391 | // Prepare the report from the execution | |
392 | string stemp = "AliFemtoV0 Pair Cut - remove shared and split pairs\n"; char ctemp[100]; | |
393 | snprintf(ctemp , 100, "Number of pairs which passed:\t%ld Number which failed:\t%ld\n",fNPairsPassed,fNPairsFailed); | |
394 | stemp += ctemp; | |
395 | AliFemtoString returnThis = stemp; | |
396 | return returnThis;} | |
397 | //__________________ | |
398 | ||
399 | void AliFemtoV0TrackPairCut::SetV0Max(Double_t aV0Max) { | |
400 | fV0Max = aV0Max; | |
401 | } | |
402 | ||
403 | Double_t AliFemtoV0TrackPairCut::GetAliFemtoV0Max() const { | |
404 | return fV0Max; | |
405 | } | |
406 | ||
407 | ||
408 | TList *AliFemtoV0TrackPairCut::ListSettings() | |
409 | { | |
410 | // return a list of settings in a writable form | |
411 | TList *tListSetttings = new TList(); | |
412 | char buf[200]; | |
413 | snprintf(buf, 200, "AliFemtoV0TrackPairCut.sharequalitymax=%f", fV0Max); | |
414 | snprintf(buf, 200, "AliFemtoV0TrackPairCut.sharefractionmax=%f", fShareFractionMax); | |
415 | tListSetttings->AddLast(new TObjString(buf)); | |
416 | ||
417 | return tListSetttings; | |
418 | } | |
419 | ||
420 | void AliFemtoV0TrackPairCut::SetRemoveSameLabel(Bool_t aRemove) | |
421 | { | |
422 | fRemoveSameLabel = aRemove; | |
423 | } | |
424 | ||
425 | void AliFemtoV0TrackPairCut::SetTPCOnly(Bool_t tpconly) | |
426 | { | |
427 | fTrackTPCOnly = tpconly; | |
428 | } | |
429 | ||
430 | void AliFemtoV0TrackPairCut::SetShareQualityMax(Double_t aShareQualityMax) { | |
431 | fShareQualityMax = aShareQualityMax; | |
432 | } | |
433 | ||
434 | void AliFemtoV0TrackPairCut::SetShareFractionMax(Double_t aShareFractionMax) { | |
435 | fShareFractionMax = aShareFractionMax; | |
436 | } | |
1ae130d9 | 437 | |
438 | void AliFemtoV0TrackPairCut::SetDataType(AliFemtoDataType type) | |
439 | { | |
440 | fDataType = type; | |
441 | } | |
442 | ||
443 | void AliFemtoV0TrackPairCut::SetTPCEntranceSepMinimum(double dtpc) | |
444 | { | |
445 | fDTPCMin = dtpc; | |
446 | } | |
447 | ||
448 | void AliFemtoV0TrackPairCut::SetTPCExitSepMinimum(double dtpc) | |
449 | { | |
450 | fDTPCExitMin = dtpc; | |
451 | } | |
ce7b3d98 | 452 | |
453 | void AliFemtoV0TrackPairCut::SetKstarCut(double kstar, AliFemtoParticleType firstParticle, AliFemtoParticleType secondParticle) | |
454 | { | |
455 | fKstarCut = kstar; | |
456 | fFirstParticleType = firstParticle; //for kstar - first particle type (V0 type) | |
457 | fSecondParticleType = secondParticle; | |
458 | } | |
6691ea4f | 459 | |
460 | void AliFemtoV0TrackPairCut::SetMinAvgSeparation(int type, double minSep) | |
461 | { | |
462 | if(type == 0) //Track-Pos | |
463 | fMinAvgSepTrackPos = minSep; | |
464 | else if(type == 1) //Track-Neg | |
465 | fMinAvgSepTrackNeg = minSep; | |
466 | } |