]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG2/KINK/AliResonanceKink.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
fix for Coverity 10651 (B.Hippolyte)
[u/mrichter/AliRoot.git] / PWG2 / KINK / AliResonanceKink.cxx
CommitLineData
10eaad41 1/**************************************************************************
2 * Author: Paraskevi Ganoti, University of Athens (pganoti@phys.uoa.gr) *
3 * Contributors are mentioned in the code where appropriate. *
4 * *
5 * Permission to use, copy, modify and distribute this software and its *
6 * documentation strictly for non-commercial purposes is hereby granted *
7 * without fee, provided that the above copyright notice appears in all *
8 * copies and that both the copyright notice and this permission notice *
9 * appear in the supporting documentation. The authors make no claims *
10 * about the suitability of this software for any purpose. It is *
11 * provided "as is" without express or implied warranty. *
12 **************************************************************************/
13
e3c296cd 14//----------------------------------------------------------------------------------------------------------------
f9afc48d 15// class AliResonanceKink
e3c296cd 16// Example of an analysis task for reconstructing resonances having at least one kaon-kink in their decay
17// products. It provides basic plots as well as plots helping to calculate the corrections.
18// Usage: To analyse a resonance having a kaon in its decay products, one has to modify the integer
19// variables resonancePDG, fdaughter1pdg and fdaughter2pdg accordingly as well as daughter1pdgMass and daughter2pdgMass.
20// Also, depending on the analysis mode (ESD or MC), fAnalysisType in the constructor must also be changed
21//-----------------------------------------------------------------------------------------------------------------
10eaad41 22
10eaad41 23#include "TH2D.h"
24#include "TParticle.h"
92adf4f6 25#include "TDatabasePDG.h"
10eaad41 26#include "TF1.h"
f9afc48d 27#include "TList.h"
894840ad 28
10eaad41 29#include "AliMCEvent.h"
f9afc48d 30#include "AliResonanceKink.h"
10eaad41 31#include "AliESDkink.h"
32#include "AliStack.h"
f9afc48d 33#include "AliESDtrack.h"
34#include "AliESDEvent.h"
92adf4f6 35#include "AliExternalTrackParam.h"
7ccf0419 36#include "AliMCParticle.h"
10eaad41 37
f9afc48d 38ClassImp(AliResonanceKink)
10eaad41 39
40//________________________________________________________________________
f9afc48d 41AliResonanceKink::AliResonanceKink()
e3c296cd 42 : TObject(), fDebug(0), fListOfHistos(0), fOpeningAngle(0), fInvariantMass(0), fInvMassTrue(0), fPhiBothKinks(0), fetabins(0), floweta(0), fupeta(0), fRecPt(0), fRecEta(0), fRecEtaPt(0), fSimPt(0), fSimEta(0), fSimEtaPt(0), fSimPtKink(0), fSimEtaKink(0), fSimEtaPtKink(0), f1(0), f2(0), fAnalysisType(), fvtxz(0), fNbins(0), fLowX(0), fHighX(0), fdaughter1pdg(0), fdaughter2pdg(0), fresonancePDGcode(0), fMaxNSigmaToVertex(0), fMinPtTrackCut(0), fMaxDCAxy(0), fMaxDCAzaxis(0),
92adf4f6 43fMinTPCclusters(0),fMaxChi2PerTPCcluster(0), fMaxCov0(0), fMaxCov2(0), fMaxCov5(0) , fMaxCov9(0), fMaxCov14(0) //, fTPCrefitFlag(kFALSE)
e3c296cd 44, fInvmassPt(0), fInvmassPtTrue(0), fMCInvmassPt(0), fMCInvmassPtTrue(0), fminKinkRadius(0), fmaxKinkRadius(0), fminQt(0), fmaxQt(0), fptbins(0), flowpt(0), fupperpt(0), fmaxAbsEtaCut(0)
10eaad41 45{
46 // Constructor
47}
48
49//________________________________________________________________________
894840ad 50AliResonanceKink::AliResonanceKink(Int_t nbins, Float_t nlowx, Float_t nhighx, Int_t netabins, Float_t nloweta, Float_t nupeta, Int_t nptbins, Float_t nlowpt, Float_t nupperpt, Int_t daughter1, Int_t daughter2, Int_t resonancePDG)
e3c296cd 51 : TObject(), fDebug(0), fListOfHistos(0), fOpeningAngle(0), fInvariantMass(0), fInvMassTrue(0), fPhiBothKinks(0), fetabins(netabins), floweta(nloweta), fupeta(nupeta), fRecPt(0), fRecEta(0), fRecEtaPt(0), fSimPt(0), fSimEta(0), fSimEtaPt(0), fSimPtKink(0), fSimEtaKink(0), fSimEtaPtKink(0), f1(0), f2(0), fAnalysisType(), fvtxz(0), fNbins(nbins), fLowX(nlowx), fHighX(nhighx), fdaughter1pdg(daughter1), fdaughter2pdg(daughter2), fresonancePDGcode(resonancePDG), fMaxNSigmaToVertex(0), fMinPtTrackCut(0),
92adf4f6 52fMaxDCAxy(0), fMaxDCAzaxis(0), fMinTPCclusters(0), fMaxChi2PerTPCcluster(0), fMaxCov0(0), fMaxCov2(0), fMaxCov5(0), fMaxCov9(0), fMaxCov14(0) //, fTPCrefitFlag(kFALSE)
e3c296cd 53, fInvmassPt(0), fInvmassPtTrue(0), fMCInvmassPt(0), fMCInvmassPtTrue(0), fminKinkRadius(0), fmaxKinkRadius(0), fminQt(0), fmaxQt(0), fptbins(nptbins), flowpt(nlowpt), fupperpt(nupperpt), fmaxAbsEtaCut(0)
10eaad41 54{
f9afc48d 55 // Constructor
10eaad41 56
10eaad41 57 fOpeningAngle=new TH1D("fOpeningAngle"," ", 100,-1.0,1.0);
58
f9afc48d 59 fInvariantMass=new TH1D("fInvariantMass"," ",fNbins,fLowX,fHighX);
60 fInvMassTrue=new TH1D("fInvMassTrue"," ",fNbins,fLowX,fHighX);
61 fPhiBothKinks=new TH1D("fPhiBothKinks"," ",fNbins,fLowX,fHighX); // Applicable for phi(1020)
10eaad41 62
894840ad 63 fRecPt=new TH1D("fRecPt"," ", nptbins, nlowpt, nupperpt);
64 fRecEta=new TH1D("fRecEta"," ",netabins, nloweta, nupeta);
65 fRecEtaPt=new TH2D("fRecEtaPt"," ", nptbins, nlowpt, nupperpt, netabins, nloweta, nupeta);
66 fSimPt=new TH1D("fSimPt"," ", nptbins, nlowpt, nupperpt);
67 fSimEta=new TH1D("fSimEta"," ", netabins, nloweta, nupeta);
68 fSimEtaPt=new TH2D("fSimEtaPt"," ", nptbins, nlowpt, nupperpt, netabins, nloweta, nupeta);
69 fSimPtKink=new TH1D("fSimPtKink"," ", nptbins, nlowpt, nupperpt);
70 fSimEtaKink=new TH1D("fSimEtaKink"," ", netabins, nloweta, nupeta);
71 fSimEtaPtKink=new TH2D("fSimEtaPtKink"," ", nptbins, nlowpt, nupperpt, netabins, nloweta, nupeta);
f9afc48d 72
73 f1=new TF1("f1","((atan([0]*[1]*(1.0/(sqrt((x^2)*(1.0-([1]^2))-([0]^2)*([1]^2))))))*180.)/[2]",1.1,10.0);
74 f1->SetParameter(0,0.493677);
75 f1->SetParameter(1,0.9127037);
76 f1->SetParameter(2,TMath::Pi());
77
78 f2=new TF1("f2","((atan([0]*[1]*(1.0/(sqrt((x^2)*(1.0-([1]^2))-([0]^2)*([1]^2))))))*180.)/[2]",0.1,10.0);
79 f2->SetParameter(0,0.13957018);
80 f2->SetParameter(1,0.2731374);
81 f2->SetParameter(2,TMath::Pi());
82
10eaad41 83 fvtxz=new TH1D("fvtxz"," ", 100,-20.0,20.0);
92adf4f6 84
7ccf0419 85 fInvmassPt=new TH2D("fInvmassPt"," ",fNbins,fLowX,fHighX,100,0.0,10.0);
86 fInvmassPtTrue=new TH2D("fInvmassPtTrue"," ",fNbins,fLowX,fHighX,100,0.0,10.0);
87 fMCInvmassPt=new TH2D("fMCInvmassPt"," ",fNbins,fLowX,fHighX,100,0.0,10.0);
88 fMCInvmassPtTrue=new TH2D("fMCInvmassPtTrue"," ",fNbins,fLowX,fHighX,100,0.0,10.0);
10eaad41 89}
90
91//________________________________________________________________________
f9afc48d 92AliResonanceKink:: ~AliResonanceKink()
10eaad41 93{
f9afc48d 94 // Destructor
95 if(fOpeningAngle) delete fOpeningAngle;
96 if(fInvariantMass) delete fInvariantMass;
97 if(fInvMassTrue) delete fInvMassTrue;
98 if(fPhiBothKinks) delete fPhiBothKinks;
99 if(fRecPt) delete fRecPt;
100 if(fRecEta) delete fRecEta;
101 if(fRecEtaPt) delete fRecEtaPt;
102 if(fSimPt) delete fSimPt;
103 if(fSimEta) delete fSimEta;
104 if(fSimEtaPt) delete fSimEtaPt;
105 if(fSimPtKink) delete fSimPtKink;
106 if(fSimEtaKink) delete fSimEtaKink;
107 if(fSimEtaPtKink) delete fSimEtaPtKink;
7ccf0419 108 if(fvtxz) delete fvtxz;
109 if(fInvmassPt) delete fInvmassPt;
110 if(fInvmassPtTrue) delete fInvmassPtTrue;
111 if(fMCInvmassPt) delete fMCInvmassPt;
112 if(fMCInvmassPtTrue) delete fMCInvmassPtTrue;
f9afc48d 113}
114//________________________________________________________________________
115TList* AliResonanceKink::GetHistogramList()
116{
117 // Adding histograms to the list
118 fListOfHistos=new TList();
10eaad41 119
f9afc48d 120 fListOfHistos->Add(fOpeningAngle);
121 fListOfHistos->Add(fInvariantMass);
122 fListOfHistos->Add(fInvMassTrue);
123 fListOfHistos->Add(fPhiBothKinks);
124 fListOfHistos->Add(fRecPt);
125 fListOfHistos->Add(fRecEta);
126 fListOfHistos->Add(fRecEtaPt);
127 fListOfHistos->Add(fSimPt);
128 fListOfHistos->Add(fSimEta);
129 fListOfHistos->Add(fSimEtaPt);
130 fListOfHistos->Add(fSimPtKink);
131 fListOfHistos->Add(fSimEtaKink);
132 fListOfHistos->Add(fSimEtaPtKink);
f9afc48d 133 fListOfHistos->Add(fvtxz);
7ccf0419 134 fListOfHistos->Add(fInvmassPt);
135 fListOfHistos->Add(fInvmassPtTrue);
136 fListOfHistos->Add(fMCInvmassPt);
137 fListOfHistos->Add(fMCInvmassPtTrue);
138
f9afc48d 139 return fListOfHistos;
140}
141
142//________________________________________________________________________
894840ad 143void AliResonanceKink::InitOutputHistograms(Int_t nbins, Float_t nlowx, Float_t nhighx, Int_t netabins, Float_t nloweta, Float_t nupeta, Int_t nptbins, Float_t nlowpt, Float_t nupperpt)
f9afc48d 144{
145 // Initialisation of the output histograms
146 fNbins=nbins;
147 fLowX=nlowx;
148 fHighX=nhighx;
149
150 fOpeningAngle=new TH1D("fOpeningAngle"," ", 100,-1.0,1.0);
151
152 fInvariantMass=new TH1D("fInvariantMass"," ",fNbins,fLowX,fHighX);
153 fInvMassTrue=new TH1D("fInvMassTrue"," ",fNbins,fLowX,fHighX);
154 fPhiBothKinks=new TH1D("fPhiBothKinks"," ",fNbins,fLowX,fHighX); // Applicable for phi(1020)
155
894840ad 156 fRecPt=new TH1D("fRecPt"," ", nptbins, nlowpt, nupperpt);
157 fRecEta=new TH1D("fRecEta"," ", netabins, nloweta, nupeta);
158 fRecEtaPt=new TH2D("fRecEtaPt"," ", nptbins, nlowpt, nupperpt, netabins, nloweta, nupeta);
159 fSimPt=new TH1D("fSimPt"," ", nptbins, nlowpt, nupperpt);
160 fSimEta=new TH1D("fSimEta"," ",netabins, nloweta, nupeta);
161 fSimEtaPt=new TH2D("fSimEtaPt"," ", nptbins, nlowpt, nupperpt, netabins, nloweta, nupeta);
162 fSimPtKink=new TH1D("fSimPtKink"," ", nptbins, nlowpt, nupperpt);
163 fSimEtaKink=new TH1D("fSimEtaKink"," ", netabins, nloweta, nupeta);
164 fSimEtaPtKink=new TH2D("fSimEtaPtKink"," ", nptbins, nlowpt, nupperpt, netabins, nloweta, nupeta);
10eaad41 165
f9afc48d 166 f1=new TF1("f1","((atan([0]*[1]*(1.0/(sqrt((x^2)*(1.0-([1]^2))-([0]^2)*([1]^2))))))*180.)/[2]",1.1,10.0);
167 f1->SetParameter(0,0.493677);
168 f1->SetParameter(1,0.9127037);
169 f1->SetParameter(2,TMath::Pi());
170
171 f2=new TF1("f2","((atan([0]*[1]*(1.0/(sqrt((x^2)*(1.0-([1]^2))-([0]^2)*([1]^2))))))*180.)/[2]",0.1,10.0);
172 f2->SetParameter(0,0.13957018);
173 f2->SetParameter(1,0.2731374);
174 f2->SetParameter(2,TMath::Pi());
f27d6e67 175
f9afc48d 176 fvtxz=new TH1D("fvtxz"," ", 100,-20.0,20.0);
7ccf0419 177
894840ad 178 fInvmassPt=new TH2D("fInvmassPt"," ",fNbins,fLowX,fHighX, nptbins, nlowpt, nupperpt);
179 fInvmassPtTrue=new TH2D("fInvmassPtTrue"," ",fNbins,fLowX,fHighX, nptbins, nlowpt, nupperpt);
180 fMCInvmassPt=new TH2D("fMCInvmassPt"," ",fNbins,fLowX,fHighX,nptbins, nlowpt, nupperpt);
181 fMCInvmassPtTrue=new TH2D("fMCInvmassPtTrue"," ",fNbins,fLowX,fHighX, nptbins, nlowpt, nupperpt);
7ccf0419 182
f9afc48d 183}
be1a7181 184
f9afc48d 185//________________________________________________________________________
186void AliResonanceKink::Analyse(AliESDEvent* esd, AliMCEvent* mcEvent)
187{
188 // Main loop
189 // Called for each event
190 Int_t resonancePDGcode, antiresonancePDGcode;
be1a7181 191 Double_t daughter1pdgMass, daughter2pdgMass;
f9afc48d 192
be1a7181 193 if (fdaughter1pdg==kKPlus) {
f9afc48d 194 resonancePDGcode=fresonancePDGcode;
195 antiresonancePDGcode=-fresonancePDGcode;
be1a7181 196 daughter1pdgMass=TDatabasePDG::Instance()->GetParticle(fdaughter1pdg)->Mass();
197 daughter2pdgMass=TDatabasePDG::Instance()->GetParticle(fdaughter2pdg)->Mass();
10eaad41 198 }
be1a7181 199
200 if (fdaughter1pdg!=kKPlus) {
f9afc48d 201 resonancePDGcode=-fresonancePDGcode;
202 antiresonancePDGcode=fresonancePDGcode;
be1a7181 203 daughter1pdgMass=TDatabasePDG::Instance()->GetParticle(fdaughter2pdg)->Mass();
204 daughter2pdgMass=TDatabasePDG::Instance()->GetParticle(fdaughter1pdg)->Mass();
205 } // to ensure that daughter1pdgMass has always the kaon mass
206
f9afc48d 207 if (fdaughter1pdg==fdaughter2pdg) {
208 resonancePDGcode=fresonancePDGcode;
209 antiresonancePDGcode=fresonancePDGcode;
210 }
f9afc48d 211 if (!esd) {
894840ad 212 Printf("ERROR: esd not available");
10eaad41 213 return;
894840ad 214 }
894840ad 215 if (!mcEvent) {
894840ad 216 return;
217 }
218
10eaad41 219 AliStack* stack=mcEvent->Stack();
220
221 if(fAnalysisType == "MC") {
7ccf0419 222
223 const AliESDVertex* vertex = GetEventVertex(esd);
224 if(!vertex) return;
225 Double_t vtx[3];
226 vertex->GetXYZ(vtx);
227 fvtxz->Fill(vtx[2]);
228
10eaad41 229 for (Int_t iMc = 0; iMc < stack->GetNprimary(); ++iMc)
230 {
231 TParticle* particle = stack->Particle(iMc);
232
233 if (!particle)
234 {
1ac5c71d 235 if (fDebug > 0) Printf("UNEXPECTED: particle with label %d not found in stack (mc loop)", iMc);
10eaad41 236 continue;
237 }
238
7ccf0419 239 if(TMath::Abs(particle->GetPdgCode())==fresonancePDGcode) {
04c3c355 240
241 if(particle->GetNDaughters()>2) continue;
7ccf0419 242
10eaad41 243 Int_t firstD=particle->GetFirstDaughter();
244 Int_t lastD=particle->GetLastDaughter();
7ccf0419 245
246 if((firstD<0)||(firstD>=stack->GetNtrack())) continue;
247 if((lastD<0)||(lastD>=stack->GetNtrack())) continue;
10eaad41 248
eaef9478 249 TParticle *daughterParticle1 = 0;
250 TParticle *daughterParticle2 = 0;
251 AliMCParticle *mcDaughter1 = 0;
252 AliMCParticle *mcDaughter2 = 0;
7ccf0419 253
be1a7181 254 if(fdaughter1pdg==kKPlus) {
7ccf0419 255 daughterParticle1=stack->Particle(firstD);
256 daughterParticle2=stack->Particle(lastD);
257 mcDaughter1= (AliMCParticle*) mcEvent->GetTrack(firstD);
258 mcDaughter2= (AliMCParticle*) mcEvent->GetTrack(lastD);
10eaad41 259 }
79b67623 260 else
be1a7181 261 if(fdaughter2pdg==kKPlus) {
7ccf0419 262 daughterParticle1=stack->Particle(lastD);
263 daughterParticle2=stack->Particle(firstD);
264 mcDaughter1= (AliMCParticle*) mcEvent->GetTrack(lastD);
265 mcDaughter2= (AliMCParticle*) mcEvent->GetTrack(firstD);
266 } //to ensure that the first daughter is always the kaon
79b67623 267
268 if(TMath::Abs(daughterParticle1->GetPdgCode())!=321) continue;
269
be1a7181 270 TParticle * daughters1Daughter=0;
271 TParticle * daughters2Daughter=0;
272 Int_t mcProcessDaughters1Daughter = -999;
273 Int_t mcProcessDaughters2Daughter = -999;
eaef9478 274 AliMCParticle *mcDaughters1Daughter = 0;
7ccf0419 275
276 if(mcDaughter1->Charge()==0) continue;
277 if(mcDaughter2->Charge()==0) continue; //accept resonance decays in two charged daughters
be1a7181 278 Int_t nDecayKaonDaughter=-99;
7ccf0419 279 for(Int_t ia=0; ia<daughterParticle1->GetNDaughters(); ia++) {
280 if(((daughterParticle1->GetFirstDaughter()+ia)>0)&&((daughterParticle1->GetFirstDaughter()+ia)<stack->GetNtrack())) {
281 daughters1Daughter=stack->Particle(daughterParticle1->GetFirstDaughter()+ia);
282 mcProcessDaughters1Daughter=daughters1Daughter->GetUniqueID();
283 if(mcProcessDaughters1Daughter==4) {
284 nDecayKaonDaughter=daughterParticle1->GetFirstDaughter()+ia;
285 break;
286 }
287 }
288 }
289
be1a7181 290 Int_t nProcessDaughter=-99;
7ccf0419 291 for(Int_t ib=0; ib<daughterParticle2->GetNDaughters(); ib++) {
292 if(((daughterParticle2->GetFirstDaughter()+ib)>0)&&((daughterParticle2->GetFirstDaughter()+ib)<stack->GetNtrack())) {
293 daughters2Daughter=stack->Particle(daughterParticle2->GetFirstDaughter()+ib);
294 mcProcessDaughters2Daughter=daughters2Daughter->GetUniqueID();
295 if((mcProcessDaughters2Daughter==4)||(mcProcessDaughters2Daughter==13)) {
296 nProcessDaughter=mcProcessDaughters2Daughter;
297 break;
298 }
299 }
300 }
301
be1a7181 302 Int_t numberOfCharged=0;
7ccf0419 303 if((mcProcessDaughters1Daughter==4)&&(nDecayKaonDaughter>=0)) {
304 for(Int_t ic=nDecayKaonDaughter; ic<=daughterParticle1->GetLastDaughter(); ic++) {
305 if ((ic>=0)&&(ic<stack->GetNtrack())) mcDaughters1Daughter= dynamic_cast<AliMCParticle*>(mcEvent->GetTrack(ic));
306 else continue;
95618715 307 if(mcDaughters1Daughter && mcDaughters1Daughter->Charge()!=0) numberOfCharged=numberOfCharged+1;
7ccf0419 308 }
309 }
310
be1a7181 311 if(numberOfCharged>=2) continue; // leave out kaon decays to more than one charged daughter
e3c296cd 312 if ((particle->Pt()>fMinPtTrackCut)&&(TMath::Abs(particle->Eta())<fmaxAbsEtaCut)) {
313
10eaad41 314 fSimEta->Fill(particle->Eta());
10eaad41 315 fSimEtaPt->Fill(particle->Pt(), particle->Eta());
7ccf0419 316 fSimPt->Fill(particle->Pt());
be1a7181 317 fMCInvmassPtTrue->Fill(particle->GetMass(), particle->Pt());
7ccf0419 318
e3c296cd 319 if((daughterParticle1->Pt()>fMinPtTrackCut)&&(TMath::Abs(daughterParticle1->Eta())<fmaxAbsEtaCut)&&(daughterParticle2->Pt()>fMinPtTrackCut)&&(TMath::Abs(daughterParticle2->Eta())<fmaxAbsEtaCut)) {
894840ad 320 if((mcProcessDaughters1Daughter==4)&&(daughters1Daughter->R()>fminKinkRadius)&&(daughters1Daughter->R()<fmaxKinkRadius)&&( (nProcessDaughter<0)||((daughters2Daughter->R()>fminKinkRadius)&&(nProcessDaughter>0)))) { //below are the findable
7ccf0419 321 fSimPtKink->Fill(particle->Pt());
322 fSimEtaKink->Fill(particle->Eta());
323 fSimEtaPtKink->Fill(particle->Pt(), particle->Eta());
324 fMCInvmassPt->Fill(particle->GetMass(), particle->Pt());
325 }
326 }
327
328 } // for the generated spectra
329
330 } //for the particular resonance
331 } //MC loop
10eaad41 332
894840ad 333 } // end fAnalysisType==MC
10eaad41 334
335 if(fAnalysisType == "ESD") {
f9afc48d 336 const AliESDVertex* vertex = GetEventVertex(esd);
10eaad41 337 if(!vertex) return;
338 Double_t vtx[3];
339 vertex->GetXYZ(vtx);
340 fvtxz->Fill(vtx[2]);
10eaad41 341 Double_t ptrackpos[3], ptrackneg[3];
342
343 TLorentzVector p4pos, anp4pos;
344 TLorentzVector p4neg, anp4neg;
345 TLorentzVector p4comb, anp4comb;
346
f9afc48d 347 for (Int_t iTracks = 0; iTracks < esd->GetNumberOfTracks(); iTracks++) {
348 AliESDtrack* trackpos = esd->GetTrack(iTracks);
10eaad41 349 if (!trackpos) {
1ac5c71d 350 if (fDebug > 0) Printf("ERROR: Could not receive track %d", iTracks);
10eaad41 351 continue;
352 }
353 if (trackpos->GetSign() < 0) continue;
354
92adf4f6 355 AliExternalTrackParam *tpcTrackpos = (AliExternalTrackParam *)trackpos->GetTPCInnerParam();
356 if(!tpcTrackpos) continue;
357 ptrackpos[0]=tpcTrackpos->Px();
358 ptrackpos[1]=tpcTrackpos->Py();
359 ptrackpos[2]=tpcTrackpos->Pz();
10eaad41 360
92adf4f6 361 Bool_t firstLevelAcceptPosTrack=IsAcceptedForKink(esd, vertex, trackpos);
362 if(firstLevelAcceptPosTrack==kFALSE) continue;
10eaad41 363
364 TVector3 posTrackMom(ptrackpos[0],ptrackpos[1],ptrackpos[2]);
92adf4f6 365
10eaad41 366 TParticle * partpos = stack->Particle(TMath::Abs(trackpos->GetLabel()));
367 if (!partpos) continue;
368 Int_t pdgpos = partpos->GetPdgCode();
369 Int_t mumlabelpos=partpos->GetFirstMother();
370 mumlabelpos = TMath::Abs(mumlabelpos);
371 TParticle * mumpos=stack->Particle(mumlabelpos);
372 if (!mumpos) continue;
373 Int_t mumpdgpos = mumpos->GetPdgCode();
374
375 Int_t indexKinkPos=trackpos->GetKinkIndex(0);
be1a7181 376
377 if(indexKinkPos>0) continue;
378
92adf4f6 379 Bool_t posKaonKinkFlag=0;
10eaad41 380
be1a7181 381 if(indexKinkPos<0) {
382 posKaonKinkFlag=IsKink(esd, indexKinkPos, posTrackMom);
383
384 if(posKaonKinkFlag==1) anp4pos.SetVectM(posTrackMom,daughter1pdgMass);
385 if(posKaonKinkFlag==0) continue;
386 }
10eaad41 387
388 if(indexKinkPos==0) {
92adf4f6 389
390 Bool_t secondLevelAcceptPosTrack=IsAcceptedForTrack(esd, vertex, trackpos);
391 if(secondLevelAcceptPosTrack==kFALSE) continue;
392
f9afc48d 393 p4pos.SetVectM(posTrackMom, daughter2pdgMass);
10eaad41 394
395 }
396
f9afc48d 397 for (Int_t j=0; j<esd->GetNumberOfTracks(); j++) {
10eaad41 398 if(iTracks==j) continue;
f9afc48d 399 AliESDtrack* trackneg=esd->GetTrack(j);
10eaad41 400 if (trackneg->GetSign() > 0) continue;
401
92adf4f6 402 AliExternalTrackParam *tpcTrackneg = (AliExternalTrackParam *)trackneg->GetTPCInnerParam();
403 if(!tpcTrackneg) continue;
404 ptrackneg[0]=tpcTrackneg->Px();
405 ptrackneg[1]=tpcTrackneg->Py();
406 ptrackneg[2]=tpcTrackneg->Pz();
10eaad41 407
92adf4f6 408 Bool_t firstLevelAcceptNegTrack=IsAcceptedForKink(esd, vertex, trackneg);
409 if(firstLevelAcceptNegTrack==kFALSE) continue;
10eaad41 410
411 TVector3 negTrackMom(ptrackneg[0],ptrackneg[1],ptrackneg[2]);
10eaad41 412
413 TParticle * partneg = stack->Particle(TMath::Abs(trackneg->GetLabel()));
414 if (!partneg) continue;
415 Int_t pdgneg = partneg->GetPdgCode();
416 Int_t mumlabelneg=partneg->GetFirstMother();
417 mumlabelneg = TMath::Abs(mumlabelneg);
418 TParticle * mumneg=stack->Particle(mumlabelneg);
419 if (!mumneg) continue;
420 Int_t mumpdgneg = mumneg->GetPdgCode();
421
422 Int_t indexKinkNeg=trackneg->GetKinkIndex(0);
be1a7181 423
424 if(indexKinkNeg>0) continue;
425
92adf4f6 426 Bool_t negKaonKinkFlag=0;
10eaad41 427
be1a7181 428 if(indexKinkNeg<0) {
429 negKaonKinkFlag=IsKink(esd, indexKinkNeg, negTrackMom);
430
431 if(negKaonKinkFlag==1) p4neg.SetVectM(negTrackMom,daughter1pdgMass);
432 if(negKaonKinkFlag==0) continue;
433 }
10eaad41 434
435 if(indexKinkNeg==0) {
92adf4f6 436
437 Bool_t secondLevelAcceptNegTrack=IsAcceptedForTrack(esd, vertex, trackneg);
438 if(secondLevelAcceptNegTrack==kFALSE) continue;
439
f9afc48d 440 anp4neg.SetVectM(negTrackMom, daughter2pdgMass);
10eaad41 441
442 }
443
444 Double_t openingAngle=(ptrackpos[0]*ptrackneg[0]+ptrackpos[1]*ptrackneg[1]+ptrackpos[2]*ptrackneg[2])/(posTrackMom.Mag()*negTrackMom.Mag());
445
be1a7181 446 if((posKaonKinkFlag==1)&&(negKaonKinkFlag==1)) {
447 p4comb=anp4pos;
448 p4comb+=p4neg;
e3c296cd 449 if((p4comb.Vect().Pt()>fMinPtTrackCut)&&(TMath::Abs(anp4pos.Vect().Eta())<fmaxAbsEtaCut)&&(TMath::Abs(p4neg.Vect().Eta())<fmaxAbsEtaCut)&&(p4comb.Vect().Eta()<fmaxAbsEtaCut)) {
be1a7181 450 if(openingAngle>0.6) fPhiBothKinks->Fill(p4comb.M());
451 }
10eaad41 452 }
453
454 if(negKaonKinkFlag==1) {
455 p4comb=p4pos;
456 p4comb+=p4neg;
be1a7181 457
894840ad 458 if(p4comb.Vect().Pt()<=fMinPtTrackCut) continue;
be1a7181 459
e3c296cd 460 if((TMath::Abs(p4pos.Vect().Eta())<fmaxAbsEtaCut)&&(TMath::Abs(p4neg.Vect().Eta())<fmaxAbsEtaCut)&&(p4comb.Vect().Eta()<fmaxAbsEtaCut)) {
be1a7181 461 fInvariantMass->Fill(p4comb.M());
462 fInvmassPt->Fill(p4comb.M(), p4comb.Vect().Pt());
463 if ((mumpdgpos==(antiresonancePDGcode))&&(mumpdgneg==(antiresonancePDGcode))&&(mumlabelpos==mumlabelneg)
464 &&(pdgpos==fdaughter2pdg)&&(pdgneg==(-fdaughter1pdg))&&(TMath::Abs(trackpos->GetLabel())>=0)&&(TMath::Abs(trackneg->GetLabel())>=0)&&(mumlabelpos>=0)&&(mumlabelneg>=0)) {
465 fOpeningAngle->Fill(openingAngle);
466 fInvMassTrue->Fill(p4comb.M());
467 fInvmassPtTrue->Fill(p4comb.M(), p4comb.Vect().Pt());
468
10eaad41 469 fRecPt->Fill(p4comb.Vect().Pt());
470 fRecEta->Fill(p4comb.Vect().Eta());
471 fRecEtaPt->Fill(p4comb.Vect().Perp(),p4comb.Vect().Eta());
472
473 }
474
475 }
476
477 }
478
92adf4f6 479 if(posKaonKinkFlag==1) {
10eaad41 480 anp4comb=anp4pos;
481 anp4comb+=anp4neg;
be1a7181 482
894840ad 483 if(anp4comb.Vect().Pt()<=fMinPtTrackCut) continue;
be1a7181 484
e3c296cd 485 if((TMath::Abs(anp4neg.Vect().Eta())<fmaxAbsEtaCut)&&(TMath::Abs(anp4pos.Vect().Eta())<fmaxAbsEtaCut)&&(anp4comb.Vect().Eta()<fmaxAbsEtaCut)) {
be1a7181 486 fInvariantMass->Fill(anp4comb.M());
487 fInvmassPt->Fill(anp4comb.M(), anp4comb.Vect().Pt());
488 if ((mumpdgpos==resonancePDGcode)&&(mumpdgneg==resonancePDGcode)&&(mumlabelpos==mumlabelneg)
489 &&(pdgpos==fdaughter1pdg)&&(pdgneg==(-fdaughter2pdg))&&(TMath::Abs(trackpos->GetLabel())>=0)&&(TMath::Abs(trackneg->GetLabel())>=0)&&(mumlabelpos>=0) &&(mumlabelneg>=0)) {
490 fOpeningAngle->Fill(openingAngle);
491 fInvMassTrue->Fill(anp4comb.M());
492 fInvmassPtTrue->Fill(anp4comb.M(), anp4comb.Vect().Pt());
493
494 fRecPt->Fill(anp4comb.Vect().Pt());
495 fRecEta->Fill(anp4comb.Vect().Eta());
496 fRecEtaPt->Fill(anp4comb.Vect().Pt(), anp4comb.Vect().Eta());
10eaad41 497 }
498
499 }
500
501 }
502
503 } //inner track loop
504
505 } //outer track loop
506
507 } // end fAnalysisType == ESD
894840ad 508
10eaad41 509}
510
894840ad 511//____________________________________________________________________//
512void AliResonanceKink::Analyse(AliESDEvent* esd)
513{
514 // Main loop
515 // Called for each event
516 Double_t daughter1pdgMass, daughter2pdgMass;
517
518 if (fdaughter1pdg==kKPlus) {
519 daughter1pdgMass=TDatabasePDG::Instance()->GetParticle(fdaughter1pdg)->Mass();
520 daughter2pdgMass=TDatabasePDG::Instance()->GetParticle(fdaughter2pdg)->Mass();
521 }
522
523 if (fdaughter1pdg!=kKPlus) {
524 daughter1pdgMass=TDatabasePDG::Instance()->GetParticle(fdaughter2pdg)->Mass();
525 daughter2pdgMass=TDatabasePDG::Instance()->GetParticle(fdaughter1pdg)->Mass();
526 } // to ensure that daughter1pdgMass has always the kaon mass
527
528 if (!esd) {
529 Printf("ERROR: esd not available");
530 return;
531 }
532
533 if(fAnalysisType == "DATA") {
534 const AliESDVertex* vertex = GetEventVertex(esd);
535 if(!vertex) return;
536 Double_t vtx[3];
537 vertex->GetXYZ(vtx);
538 fvtxz->Fill(vtx[2]);
539 Double_t ptrackpos[3], ptrackneg[3];
540
541 TLorentzVector p4pos, anp4pos;
542 TLorentzVector p4neg, anp4neg;
543 TLorentzVector p4comb, anp4comb;
544
545 for (Int_t iTracks = 0; iTracks < esd->GetNumberOfTracks(); iTracks++) {
546 AliESDtrack* trackpos = esd->GetTrack(iTracks);
547 if (!trackpos) {
548 if (fDebug > 0) Printf("ERROR: Could not receive track %d", iTracks);
549 continue;
550 }
551 if (trackpos->GetSign() < 0) continue;
552
553 AliExternalTrackParam *tpcTrackpos = (AliExternalTrackParam *)trackpos->GetTPCInnerParam();
554 if(!tpcTrackpos) continue;
555 ptrackpos[0]=tpcTrackpos->Px();
556 ptrackpos[1]=tpcTrackpos->Py();
557 ptrackpos[2]=tpcTrackpos->Pz();
558
559 Bool_t firstLevelAcceptPosTrack=IsAcceptedForKink(esd, vertex, trackpos);
560 if(firstLevelAcceptPosTrack==kFALSE) continue;
561
562 TVector3 posTrackMom(ptrackpos[0],ptrackpos[1],ptrackpos[2]);
563
564 Int_t indexKinkPos=trackpos->GetKinkIndex(0);
565
566 if(indexKinkPos>0) continue;
567
568 Bool_t posKaonKinkFlag=0;
569
570 if(indexKinkPos<0) {
571 posKaonKinkFlag=IsKink(esd, indexKinkPos, posTrackMom);
572
573 if(posKaonKinkFlag==1) anp4pos.SetVectM(posTrackMom,daughter1pdgMass);
574 if(posKaonKinkFlag==0) continue;
575 }
576
577 if(indexKinkPos==0) {
578
579 Bool_t secondLevelAcceptPosTrack=IsAcceptedForTrack(esd, vertex, trackpos);
580 if(secondLevelAcceptPosTrack==kFALSE) continue;
581
582 p4pos.SetVectM(posTrackMom, daughter2pdgMass);
583
584 }
585
586 for (Int_t j=0; j<esd->GetNumberOfTracks(); j++) {
587 if(iTracks==j) continue;
588 AliESDtrack* trackneg=esd->GetTrack(j);
589 if (trackneg->GetSign() > 0) continue;
590
591 AliExternalTrackParam *tpcTrackneg = (AliExternalTrackParam *)trackneg->GetTPCInnerParam();
592 if(!tpcTrackneg) continue;
593 ptrackneg[0]=tpcTrackneg->Px();
594 ptrackneg[1]=tpcTrackneg->Py();
595 ptrackneg[2]=tpcTrackneg->Pz();
596
597 Bool_t firstLevelAcceptNegTrack=IsAcceptedForKink(esd, vertex, trackneg);
598 if(firstLevelAcceptNegTrack==kFALSE) continue;
599
600 TVector3 negTrackMom(ptrackneg[0],ptrackneg[1],ptrackneg[2]);
601
602 Int_t indexKinkNeg=trackneg->GetKinkIndex(0);
603
604 if(indexKinkNeg>0) continue;
605
606 Bool_t negKaonKinkFlag=0;
607
608 if(indexKinkNeg<0) {
609 negKaonKinkFlag=IsKink(esd, indexKinkNeg, negTrackMom);
610
611 if(negKaonKinkFlag==1) p4neg.SetVectM(negTrackMom,daughter1pdgMass);
612 if(negKaonKinkFlag==0) continue;
613 }
614
615 if(indexKinkNeg==0) {
616
617 Bool_t secondLevelAcceptNegTrack=IsAcceptedForTrack(esd, vertex, trackneg);
618 if(secondLevelAcceptNegTrack==kFALSE) continue;
619
620 anp4neg.SetVectM(negTrackMom, daughter2pdgMass);
621
622 }
623
624 Double_t openingAngle=(ptrackpos[0]*ptrackneg[0]+ptrackpos[1]*ptrackneg[1]+ptrackpos[2]*ptrackneg[2])/(posTrackMom.Mag()*negTrackMom.Mag());
625
626 if((posKaonKinkFlag==1)&&(negKaonKinkFlag==1)) {
627 p4comb=anp4pos;
628 p4comb+=p4neg;
e3c296cd 629 if((p4comb.Vect().Pt()>fMinPtTrackCut)&&(TMath::Abs(anp4pos.Vect().Eta())<fmaxAbsEtaCut)&&(TMath::Abs(p4neg.Vect().Eta())<fmaxAbsEtaCut)&&(p4comb.Vect().Eta()<fmaxAbsEtaCut)) {
894840ad 630 if(openingAngle>0.6) fPhiBothKinks->Fill(p4comb.M());
631 }
632 }
633
634 if(negKaonKinkFlag==1) {
635 p4comb=p4pos;
636 p4comb+=p4neg;
637
638 if(p4comb.Vect().Pt()<=fMinPtTrackCut) continue;
639
e3c296cd 640 if((TMath::Abs(p4pos.Vect().Eta())<fmaxAbsEtaCut)&&(TMath::Abs(p4neg.Vect().Eta())<fmaxAbsEtaCut)&&(p4comb.Vect().Eta()<fmaxAbsEtaCut)) {
894840ad 641 fInvariantMass->Fill(p4comb.M());
642 fInvmassPt->Fill(p4comb.M(), p4comb.Vect().Pt());
e3c296cd 643 fRecPt->Fill(p4comb.Vect().Pt());
644 fRecEta->Fill(p4comb.Vect().Eta());
645 fRecEtaPt->Fill(p4comb.Vect().Perp(),p4comb.Vect().Eta());
894840ad 646 }
647 }
648
649 if(posKaonKinkFlag==1) {
650 anp4comb=anp4pos;
651 anp4comb+=anp4neg;
652
653 if(anp4comb.Vect().Pt()<=fMinPtTrackCut) continue;
654
e3c296cd 655 if((TMath::Abs(anp4neg.Vect().Eta())<fmaxAbsEtaCut)&&(TMath::Abs(anp4pos.Vect().Eta())<fmaxAbsEtaCut)&&(anp4comb.Vect().Eta()<fmaxAbsEtaCut)) {
894840ad 656 fInvariantMass->Fill(anp4comb.M());
657 fInvmassPt->Fill(anp4comb.M(), anp4comb.Vect().Pt());
e3c296cd 658 fRecPt->Fill(anp4comb.Vect().Pt());
659 fRecEta->Fill(anp4comb.Vect().Eta());
660 fRecEtaPt->Fill(anp4comb.Vect().Pt(), anp4comb.Vect().Eta());
894840ad 661 }
662
663 }
664
665 } //inner track loop
666
667 } //outer track loop
668
669 } // end fAnalysisType == DATA
670
671}
10eaad41 672//____________________________________________________________________//
f9afc48d 673Float_t AliResonanceKink::GetSigmaToVertex(AliESDtrack* esdTrack) const {
10eaad41 674 // Calculates the number of sigma to the vertex.
675
676 Float_t b[2];
677 Float_t bRes[2];
678 Float_t bCov[3];
679
92adf4f6 680 esdTrack->GetImpactParametersTPC(b,bCov);
10eaad41 681
682 if (bCov[0]<=0 || bCov[2]<=0) {
10eaad41 683 bCov[0]=0; bCov[2]=0;
684 }
f9afc48d 685
10eaad41 686 bRes[0] = TMath::Sqrt(bCov[0]);
687 bRes[1] = TMath::Sqrt(bCov[2]);
688
689 if (bRes[0] == 0 || bRes[1] ==0) return -1;
690
691 Float_t d = TMath::Sqrt(TMath::Power(b[0]/bRes[0],2) + TMath::Power(b[1]/bRes[1],2));
692
693 if (TMath::Exp(-d * d / 2) < 1e-10) return 1000;
694
695 d = TMath::ErfInverse(1 - TMath::Exp(-d * d / 2)) * TMath::Sqrt(2);
696
697 return d;
698}
699
f9afc48d 700//________________________________________________________________________
701const AliESDVertex* AliResonanceKink::GetEventVertex(const AliESDEvent* esd) const
10eaad41 702{
703 // Get the vertex
704
705 const AliESDVertex* vertex = esd->GetPrimaryVertex();
706
707 if((vertex->GetStatus()==kTRUE)&&(vertex->GetNContributors()>2)) return vertex;
708 else
709 {
710 vertex = esd->GetPrimaryVertexSPD();
711 if((vertex->GetStatus()==kTRUE)&&(vertex->GetNContributors()>2)) return vertex;
712 else
713 return 0;
714 }
715}
716
92adf4f6 717//________________________________________________________________________
718
719 Bool_t AliResonanceKink::IsAcceptedForKink(AliESDEvent *localesd,
720 const AliESDVertex *localvertex, AliESDtrack* localtrack) {
721 // Apply the selections for each kink
722
723 Double_t gPt = 0.0, gPx = 0.0, gPy = 0.0, gPz = 0.0;
724 Double_t dca[2] = {0.0,0.0}, cov[3] = {0.0,0.0,0.0}; //The impact parameters and their covariance.
725 Double_t dca3D = 0.0;
726
727 AliExternalTrackParam *tpcTrack = (AliExternalTrackParam *)localtrack->GetTPCInnerParam();
728 if(!tpcTrack) {
729 gPt = 0.0; gPx = 0.0; gPy = 0.0; gPz = 0.0;
730 dca[0] = -100.; dca[1] = -100.; dca3D = -100.;
731 cov[0] = -100.; cov[1] = -100.; cov[2] = -100.;
732 }
733 else {
734 gPt = tpcTrack->Pt();
735 gPx = tpcTrack->Px();
736 gPy = tpcTrack->Py();
737 gPz = tpcTrack->Pz();
738 tpcTrack->PropagateToDCA(localvertex,
739 localesd->GetMagneticField(),100.,dca,cov);
740 }
741
742 if(GetSigmaToVertex(localtrack) > fMaxNSigmaToVertex) {
1ac5c71d 743 if (fDebug > 1) Printf("IsAcceptedKink: Track rejected because it has a %lf sigmas to vertex TPC (max. requested: %lf)", GetSigmaToVertex(localtrack),fMaxNSigmaToVertex);
92adf4f6 744 return kFALSE;
745 }
746
747 if(TMath::Abs(dca[0]) > fMaxDCAxy) {
1ac5c71d 748 if (fDebug > 1) Printf("IsAcceptedKink: Track rejected because it has a value of dca(xy) (TPC) of %lf (max. requested: %lf)", TMath::Abs(dca[0]), fMaxDCAxy);
92adf4f6 749 return kFALSE;
750 }
751
752 if(TMath::Abs(dca[1]) > fMaxDCAzaxis) {
1ac5c71d 753 if (fDebug > 1) Printf("IsAcceptedKink: Track rejected because it has a value of dca(z) of %lf (max. requested: %lf)", TMath::Abs(dca[1]), fMaxDCAzaxis);
92adf4f6 754 return kFALSE;
755 }
756
be1a7181 757 if(gPt <= fMinPtTrackCut) {
1ac5c71d 758 if (fDebug > 1) Printf("IsAcceptedKink: Track rejected because it has a min value of pt of %lf (min. requested: %lf)", gPt, fMinPtTrackCut);
92adf4f6 759 return kFALSE;
760 }
761
762 return kTRUE;
763}
764
765//________________________________________________________________________
766Bool_t AliResonanceKink::IsAcceptedForTrack(AliESDEvent *localesd, const AliESDVertex *localvertex, AliESDtrack *localtrack) {
767 // Apply the selections for each track
f9afc48d 768
92adf4f6 769 Double_t gPt = 0.0, gPx = 0.0, gPy = 0.0, gPz = 0.0;
770 Double_t dca[2] = {0.0,0.0}, cov[3] = {0.0,0.0,0.0}; //The impact parameters and their covariance.
771 Double_t dca3D = 0.0;
772
773 AliExternalTrackParam *tpcTrack = (AliExternalTrackParam *)localtrack->GetTPCInnerParam();
774 if(!tpcTrack) {
775 gPt = 0.0; gPx = 0.0; gPy = 0.0; gPz = 0.0;
776 dca[0] = -100.; dca[1] = -100.; dca3D = -100.;
777 cov[0] = -100.; cov[1] = -100.; cov[2] = -100.;
778 }
779 else {
780 gPt = tpcTrack->Pt();
781 gPx = tpcTrack->Px();
782 gPy = tpcTrack->Py();
783 gPz = tpcTrack->Pz();
784 tpcTrack->PropagateToDCA(localvertex,
785 localesd->GetMagneticField(),100.,dca,cov);
786 }
787
788 Int_t fcls[200];
789 Int_t nClustersTPC=localtrack->GetTPCclusters(fcls);
790 Float_t chi2perTPCcluster=-1.0;
791 if(nClustersTPC!=0) chi2perTPCcluster=(localtrack->GetTPCchi2())/Float_t(nClustersTPC);
792
793 Double_t extCov[15];
794 localtrack->GetExternalCovariance(extCov);
795
796 if((localtrack->GetStatus() & AliESDtrack::kTPCrefit) == 0) {
1ac5c71d 797 if (fDebug > 1) Printf("IsAccepted: Track rejected because of no refited in TPC");
92adf4f6 798 return kFALSE;
799 }
800
801 if(nClustersTPC < fMinTPCclusters) {
44b5c57e 802 if (fDebug > 1) Printf("IsAccepted: Track rejected because it has a value of nclusters (TPC) of %d (min. requested: %d)", nClustersTPC, fMinTPCclusters);
92adf4f6 803 return kFALSE;
804 }
805
806 if(chi2perTPCcluster > fMaxChi2PerTPCcluster) {
1ac5c71d 807 if (fDebug > 1) Printf("IsAccepted: Track rejected because it has a value of chi2perTPCcluster of %lf (max. requested: %lf)", chi2perTPCcluster, fMaxChi2PerTPCcluster);
92adf4f6 808 return kFALSE;
809 }
810
811 if(extCov[0] > fMaxCov0) {
9ee1c45f 812 if (fDebug > 1) Printf("IsAccepted: Track rejected because it has a value of cov[0] of %lf (max. requested: %lf)", extCov[0], fMaxCov0);
92adf4f6 813 return kFALSE;
814 }
815
816 if(extCov[2] > fMaxCov2) {
9ee1c45f 817 if (fDebug > 1) Printf("IsAccepted: Track rejected because it has a value of cov[2] of %lf (max. requested: %lf)", extCov[2], fMaxCov2);
92adf4f6 818 return kFALSE;
819 }
820
821 if(extCov[5] > fMaxCov5) {
9ee1c45f 822 if (fDebug > 1) Printf("IsAccepted: Track rejected because it has a value of cov[5] of %lf (max. requested: %lf)", extCov[5], fMaxCov5);
92adf4f6 823 return kFALSE;
824 }
825
826 if(extCov[9] > fMaxCov9) {
9ee1c45f 827 if (fDebug > 1) Printf("IsAccepted: Track rejected because it has a value of cov[9] of %lf (max. requested: %lf)", extCov[9], fMaxCov9);
92adf4f6 828 return kFALSE;
829 }
830
831 if(extCov[14] > fMaxCov14) {
9ee1c45f 832 if (fDebug > 1) Printf("IsAccepted: Track rejected because it has a value of cov[14] of %lf (max. requested: %lf)", extCov[14], fMaxCov14);
92adf4f6 833 return kFALSE;
834 }
f9afc48d 835
92adf4f6 836 return kTRUE;
837
838}
839
840//________________________________________________________________________
841Bool_t AliResonanceKink::IsKink(AliESDEvent *localesd, Int_t kinkIndex, TVector3 trackMom)
842{
843 // Test some kinematical criteria for each kink
844
845 AliESDkink *kink=localesd->GetKink(TMath::Abs(kinkIndex)-1);
846 const TVector3 motherMfromKink(kink->GetMotherP());
847 const TVector3 daughterMKink(kink->GetDaughterP());
848 Float_t qt=kink->GetQt();
849
850 Double_t maxDecAngKmu=f1->Eval(motherMfromKink.Mag(),0.,0.,0.);
851 Double_t maxDecAngpimu=f2->Eval(motherMfromKink.Mag(),0.,0.,0.);
852
853 Float_t kinkAngle=TMath::RadToDeg()*kink->GetAngle(2);
854
855 Float_t energyDaughterMu=TMath::Sqrt(daughterMKink.Mag()*daughterMKink.Mag()+0.105658*0.105658);
856 Float_t p1XM= motherMfromKink.Px();
857 Float_t p1YM= motherMfromKink.Py();
858 Float_t p1ZM= motherMfromKink.Pz();
859 Float_t p2XM= daughterMKink.Px();
860 Float_t p2YM= daughterMKink.Py();
861 Float_t p2ZM= daughterMKink.Pz();
862 Float_t p3Daughter=TMath::Sqrt(((p1XM-p2XM)*(p1XM-p2XM))+((p1YM-p2YM)*(p1YM-p2YM))+((p1ZM-p2ZM)*(p1ZM-p2ZM)));
863 Double_t invariantMassKmu= TMath::Sqrt((energyDaughterMu+p3Daughter)*(energyDaughterMu+p3Daughter)-motherMfromKink.Mag()*motherMfromKink.Mag());
864
e3c296cd 865 if((kinkAngle>maxDecAngpimu)&&(qt>fminQt)&&(qt<fmaxQt)&&((kink->GetR()>fminKinkRadius)&&(kink->GetR()<fmaxKinkRadius))&&(TMath::Abs(trackMom.Eta())<fmaxAbsEtaCut)&&(invariantMassKmu<0.6)) {
92adf4f6 866
867 if(trackMom.Mag()<=1.1) {
868 return kTRUE;
869 }
870 else
871 if (kinkAngle<maxDecAngKmu) {
872 return kTRUE;
873 }
874 }
875 return kFALSE;
876}