Adding Kink directory and tasks from Evi (pganoti@phys.uoa.gr) and Martha (mspyrop...
[u/mrichter/AliRoot.git] / PWG2 / KINK / AliResonanceKinkLikeSign.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
14//----------------------------------------------------------------------------------------------------------------
15// class AliResonanceKinkLikeSign
16// Example of an analysis task for producing a like-sign background for resonances having at least one
17// kaon-kink in their decay products.
18// Background is calculated from a negative kaon kink and a negative track.
19//-----------------------------------------------------------------------------------------------------------------
20
21#include "TChain.h"
22#include "TTree.h"
23#include "TVector3.h"
24#include "TF1.h"
25#include "TH1D.h"
26#include "TLorentzVector.h"
27#include "AliAnalysisTaskSE.h"
28#include "AliAnalysisManager.h"
29
30#include "AliESDInputHandler.h"
31
32#include "AliResonanceKinkLikeSign.h"
33#include "AliESDkink.h"
34
35ClassImp(AliResonanceKinkLikeSign)
36
37//________________________________________________________________________
38AliResonanceKinkLikeSign::AliResonanceKinkLikeSign()
39 : AliAnalysisTaskSE(), fESD(0), fListOfHistos(0), f1(0), f2(0), fNegKaonLikeSign(0)
40{
41 // Constructor
42}
43
44//________________________________________________________________________
45AliResonanceKinkLikeSign::AliResonanceKinkLikeSign(const char *name)
46 : AliAnalysisTaskSE(name), fESD(0), fListOfHistos(0), f1(0), f2(0), fNegKaonLikeSign(0)
47
48{
49 // Constructor
50
51 // Define input and output slots here
52 // Input slot #0 works with a TChain
53 DefineInput(0, TChain::Class());
54 DefineOutput(1, TList::Class());
55}
56
57//________________________________________________________________________
58void AliResonanceKinkLikeSign::ConnectInputData(Option_t *)
59{
60 // Connect ESD or AOD here
61 // Called once
62 AliPID *pidObject;
63 pidObject->Init();
64
65 TTree* tree = dynamic_cast<TTree*> (GetInputData(0));
66 if (!tree) {
67 Printf("ERROR: Could not read chain from input slot 0");
68 } else {
69 // Disable all branches, we want to process only MC
70 tree->SetBranchStatus("*",kTRUE );
71 tree->SetBranchStatus("*Calo*", kFALSE);
72
73 AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
74
75 if (!esdH) {
76 Printf("ERROR: Could not get ESDInputHandler");
77 } else
78 fESD = esdH->GetEvent();
79 }
80}
81
82//________________________________________________________________________
83void AliResonanceKinkLikeSign::CreateOutputObjects()
84{
85 // Create histograms
86 // Called once
87
88 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);
89 f1->SetParameter(0,0.493677);
90 f1->SetParameter(1,0.9127037);
91 f1->SetParameter(2,TMath::Pi());
92
93 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);
94 f2->SetParameter(0,0.13957018);
95 f2->SetParameter(1,0.2731374);
96 f2->SetParameter(2,TMath::Pi());
97
98 //OpenFile(1); //uncomment for proof analysis
99
100 // for K*0(892)
101 fNegKaonLikeSign=new TH1D("fNegKaonLikeSign"," ", 60,0.6,1.2);
102
103 // for phi(1020)
104 // fNegKaonLikeSign=new TH1D("fNegKaonLikeSign"," ", 70,0.99,1.088);
105
106
107 fListOfHistos=new TList();
108 fListOfHistos->Add(fNegKaonLikeSign);
109
110}
111
112//________________________________________________________________________
113void AliResonanceKinkLikeSign::Exec(Option_t *)
114{
115 // Main loop
116 // Called for each event
117
118 Float_t daughter1Mass=AliPID::ParticleMass(AliPID::kKaon);
119 Float_t daughter2Mass=AliPID::ParticleMass(AliPID::kPion);
120
121 if (!fESD) {
122 Printf("ERROR: fESD not available");
123 return;
124 }
125
126 const AliESDVertex* vertex = GetEventVertex(fESD);
127 if(!vertex) return;
128
129 Double_t ptrackpos[3], ptrackneg[3];
130
131 TLorentzVector p4pos;
132 TLorentzVector p4neg;
133 TLorentzVector p4comb;
134
135 for (Int_t iTracks = 0; iTracks < fESD->GetNumberOfTracks(); iTracks++) {
136 AliESDtrack* trackpos = fESD->GetTrack(iTracks);
137 if (!trackpos) {
138 Printf("ERROR: Could not receive track %d", iTracks);
139 continue;
140 }
141
142 if (trackpos->GetSign() > 0) continue;
143
144 trackpos->GetPxPyPz(ptrackpos);
145
146 Float_t nSigmaToVertex = GetSigmaToVertex(trackpos);
147
148 Float_t bpos[2];
149 Float_t bCovpos[3];
150 trackpos->GetImpactParameters(bpos,bCovpos);
151
152 if (bCovpos[0]<=0 || bCovpos[2]<=0) {
153 Printf("Estimated b resolution lower or equal zero!");
154 bCovpos[0]=0; bCovpos[2]=0;
155 }
156
157 Float_t dcaToVertexXYpos = bpos[0];
158 Float_t dcaToVertexZpos = bpos[1];
159
160 if(nSigmaToVertex>=4) continue;
161 if((dcaToVertexXYpos>3.0)||(dcaToVertexZpos>3.0)) continue;
162
163 TVector3 posTrackMom(ptrackpos[0],ptrackpos[1],ptrackpos[2]);
164
165 if(posTrackMom.Perp()<=0.25) continue;
166
167 //Uncomment the following block if the Like Sign is made of K- kink + negative track
168
169 //Int_t indexKink=trackpos->GetKinkIndex(0);
170 //Int_t kaonKinkFlag=0;
171// if(indexKink<0){
172//
173// AliESDkink *kink=fESD->GetKink(TMath::Abs(IndexKink)-1);
174// const TVector3 motherMfromKink(kink->GetMotherP());
175// const TVector3 daughterMKink(kink->GetDaughterP());
176// Float_t qt=kink->GetQt();
177//
178// Double_t maxDecAngKmu=f1->Eval(motherMfromKink.Mag(),0.,0.,0.);
179// Double_t maxDecAngpimu=f2->Eval(motherMfromKink.Mag(),0.,0.,0.);
180//
181// Float_t kinkAngle=TMath::RadToDeg()*kink->GetAngle(2);
182//
183// Float_t energyDaughterMu=TMath::Sqrt(daughterMKink.Mag()*daughterMKink.Mag()+0.105658*0.105658);
184// Float_t p1XM= motherMfromKink.Px();
185// Float_t p1YM= motherMfromKink.Py();
186// Float_t p1ZM= motherMfromKink.Pz();
187// Float_t p2XM= daughterMKink.Px();
188// Float_t p2YM= daughterMKink.Py();
189// Float_t p2ZM= daughterMKink.Pz();
190// Float_t p3Daughter=TMath::Sqrt(((p1XM-p2XM)*(p1XM-p2XM))+((p1YM-p2YM)*(p1YM-p2YM))+((p1ZM-p2ZM)*(p1ZM-p2ZM)));
191// Double_t invariantMassKmu= TMath::Sqrt((energyDaughterMu+p3Daughter)*(energyDaughterMu+p3Daughter)-motherMfromKink.Mag()*motherMfromKink.Mag());
192//
193// if((kinkAngle>maxDecAngpimu)&&(qt>0.05)&&(qt<0.25)&&((kink->GetR()>110.)&&(kink->GetR()<230.))&&(TMath::Abs(posTrackMom.Eta())<1.1)&&(invariantMassKmu<0.6)) {
194//
195// if(posTrackMom.Mag()<=1.1) {
196// kaonKinkFlag=1;
197// }
198// else
199// if (kinkAngle<maxDecAngKmu) {
200// kaonKinkFlag=1;
201// }
202// }
203//
204// } //End Kink Information
205//
206// if(kaonKinkFlag==0) continue;
207// if(kaonKinkFlag==1) p4pos.SetVectM(posTrackMom,daughter1Mass);
208
209 // Comment the following statements till the "for" if the Like Sign of K- kink + negative track is needed
210
211 UInt_t status=trackpos->GetStatus();
212 if((status&AliESDtrack::kTPCrefit)==0) continue;
213 if(trackpos->GetTPCclusters(0)<50) continue;
214 if((trackpos->GetTPCchi2()/trackpos->GetTPCclusters(0))>3.5) continue;
215 Double_t extCovPos[15];
216 trackpos->GetExternalCovariance(extCovPos);
217 if(extCovPos[0]>2) continue;
218 if(extCovPos[2]>2) continue;
219 if(extCovPos[5]>0.5) continue;
220 if(extCovPos[9]>0.5) continue;
221 if(extCovPos[14]>2) continue;
222
223 p4pos.SetVectM(posTrackMom,daughter1Mass);
224
225 for (Int_t j=0; j<fESD->GetNumberOfTracks(); j++) {
226 if(iTracks==j) continue;
227 AliESDtrack* trackneg=fESD->GetTrack(j);
228 if (trackneg->GetSign() > 0) continue;
229
230 trackneg->GetPxPyPz(ptrackneg);
231 Float_t negSigmaToVertex = GetSigmaToVertex(trackneg);
232
233 Float_t bneg[2];
234 Float_t bCovneg[3];
235 trackneg->GetImpactParameters(bneg,bCovneg);
236 if (bCovneg[0]<=0 || bCovneg[2]<=0) {
237 Printf("Estimated b resolution lower or equal zero!");
238 bCovneg[0]=0; bCovneg[2]=0;
239 }
240
241 Float_t dcaToVertexXYneg = bneg[0];
242 Float_t dcaToVertexZneg = bneg[1];
243
244 if(negSigmaToVertex>=4) continue;
245 if((dcaToVertexXYneg>3.0)||(dcaToVertexZneg>3.0)) continue;
246
247 TVector3 negTrackMom(ptrackneg[0],ptrackneg[1],ptrackneg[2]);
248
249 if(negTrackMom.Perp()<=0.25) continue;
250
251 UInt_t statusneg=trackneg->GetStatus();
252
253 if((statusneg&AliESDtrack::kTPCrefit)==0) continue;
254
255 if(trackneg->GetTPCclusters(0)<50) continue;
256 if((trackneg->GetTPCchi2()/trackneg->GetTPCclusters(0))>3.5) continue;
257 Double_t extCovneg[15];
258 trackneg->GetExternalCovariance(extCovneg);
259 if(extCovneg[0]>2) continue;
260 if(extCovneg[2]>2) continue;
261 if(extCovneg[5]>0.5) continue;
262 if(extCovneg[9]>0.5) continue;
263 if(extCovneg[14]>2) continue;
264
265 p4neg.SetVectM(negTrackMom, daughter2Mass);
266
267 p4comb=p4pos;
268 p4comb+=p4neg;
269 fNegKaonLikeSign->Fill(p4comb.M());
270
271 } //inner track loop
272
273 } //outer track loop
274
275 PostData(1, fListOfHistos);
276}
277
278//________________________________________________________________________
279void AliResonanceKinkLikeSign::Terminate(Option_t *)
280{
281
282}
283
284//____________________________________________________________________//
285
286Float_t AliResonanceKinkLikeSign::GetSigmaToVertex(AliESDtrack* esdTrack) const {
287 // Calculates the number of sigma to the vertex.
288
289 Float_t b[2];
290 Float_t bRes[2];
291 Float_t bCov[3];
292
293 esdTrack->GetImpactParameters(b,bCov);
294
295 if (bCov[0]<=0 || bCov[2]<=0) {
296 //AliDebug(1, "Estimated b resolution lower or equal zero!");
297 bCov[0]=0; bCov[2]=0;
298 }
299 bRes[0] = TMath::Sqrt(bCov[0]);
300 bRes[1] = TMath::Sqrt(bCov[2]);
301
302 if (bRes[0] == 0 || bRes[1] ==0) return -1;
303
304 Float_t d = TMath::Sqrt(TMath::Power(b[0]/bRes[0],2) + TMath::Power(b[1]/bRes[1],2));
305
306 if (TMath::Exp(-d * d / 2) < 1e-10) return 1000;
307
308 d = TMath::ErfInverse(1 - TMath::Exp(-d * d / 2)) * TMath::Sqrt(2);
309
310 return d;
311}
312
313//____________________________________________________________________//
314
315const AliESDVertex* AliResonanceKinkLikeSign::GetEventVertex(const AliESDEvent* esd) const
316
317{
318 // Get the vertex
319
320 const AliESDVertex* vertex = esd->GetPrimaryVertex();
321
322 if((vertex->GetStatus()==kTRUE)&&(vertex->GetNContributors()>2)) return vertex;
323 else
324 {
325 vertex = esd->GetPrimaryVertexSPD();
326 if((vertex->GetStatus()==kTRUE)&&(vertex->GetNContributors()>2)) return vertex;
327 else
328 return 0;
329 }
330}
331