1 /**************************************************************************
2 * Copyright(c) 1998-2010, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
16 /////////////////////////////////////////////////////////////
18 // Class for cuts on AOD reconstructed D+->Kpipi
20 // Author: R. Bala, bala@to.infn.it
21 // G. Ortona, ortona@to.infn.it
22 /////////////////////////////////////////////////////////////
24 #include <TDatabasePDG.h>
25 #include <Riostream.h>
27 #include "AliRDHFCutsDplustoKpipi.h"
28 #include "AliAODRecoDecayHF3Prong.h"
29 #include "AliAODTrack.h"
30 #include "AliESDtrack.h"
33 ClassImp(AliRDHFCutsDplustoKpipi)
35 //--------------------------------------------------------------------------
36 AliRDHFCutsDplustoKpipi::AliRDHFCutsDplustoKpipi(const char* name) :
40 // Default Constructor
44 TString varNames[12]={"inv. mass [GeV]",
47 "d0K [cm] lower limit!",
48 "d0Pi [cm] lower limit!",
52 "pM=Max{pT1,pT2,pT3} (GeV/c)",
56 Bool_t isUpperCut[12]={kTRUE,
68 SetVarNames(nvars,varNames,isUpperCut);
69 Bool_t forOpt[12]={kFALSE,
81 SetVarsForOpt(5,forOpt);
82 Float_t limits[2]={0,999999999.};
84 if(fPidHF)delete fPidHF;
85 fPidHF=new AliAODPidHF();
86 Double_t plim[2]={0.6,0.8};
87 Double_t nsigma[5]={2.,1.,2.,3.,0.};
89 fPidHF->SetPLimit(plim);
90 fPidHF->SetAsym(kTRUE);
91 fPidHF->SetSigma(nsigma);
97 fPidHF->SetCompat(kTRUE);
101 //--------------------------------------------------------------------------
102 AliRDHFCutsDplustoKpipi::AliRDHFCutsDplustoKpipi(const AliRDHFCutsDplustoKpipi &source) :
110 //--------------------------------------------------------------------------
111 AliRDHFCutsDplustoKpipi &AliRDHFCutsDplustoKpipi::operator=(const AliRDHFCutsDplustoKpipi &source)
114 // assignment operator
116 if(&source == this) return *this;
118 AliRDHFCuts::operator=(source);
125 //---------------------------------------------------------------------------
126 void AliRDHFCutsDplustoKpipi::GetCutVarsForOpt(AliAODRecoDecayHF *d,Float_t *vars,Int_t nvars,Int_t *pdgdaughters) {
128 // Fills in vars the values of the variables
132 if(nvars!=fnVarsForOpt) {
133 printf("AliRDHFCutsDplustoKpipi::GetCutsVarsForOpt: wrong number of variables\n");
137 AliAODRecoDecayHF3Prong *dd = (AliAODRecoDecayHF3Prong*)d;
142 vars[iter]=dd->InvMassDplus();
146 for(Int_t iprong=0;iprong<3;iprong++){
147 if(TMath::Abs(pdgdaughters[iprong])==321) {
148 vars[iter]=dd->PtProng(iprong);
154 Float_t minPtDau=1000000.0;
155 for(Int_t iprong=0;iprong<3;iprong++){
156 if(TMath::Abs(pdgdaughters[iprong])==211) {
157 if(dd->PtProng(iprong)<minPtDau){
158 minPtDau=dd->PtProng(iprong);
166 for(Int_t iprong=0;iprong<3;iprong++){
167 if(TMath::Abs(pdgdaughters[iprong])==321) {
168 vars[iter]=dd->Getd0Prong(iprong);
174 Float_t minImpParDau=1000000.0;
175 for(Int_t iprong=0;iprong<3;iprong++){
176 if(TMath::Abs(pdgdaughters[iprong])==211) {
177 if(dd->Getd0Prong(iprong)<minImpParDau){
178 minImpParDau=dd->Getd0Prong(iprong);
182 vars[iter]=minImpParDau;
186 Float_t dist12 = dd->GetDist12toPrim();
187 Float_t dist23 = dd->GetDist23toPrim();
188 if(dist12<dist23)vars[iter]=dist12;
189 else vars[iter]=dist23;
193 vars[iter]=dd->GetSigmaVert();
197 vars[iter] = dd->DecayLength();
202 for(Int_t i=0;i<3;i++){
203 if(dd->PtProng(i)>ptmax)ptmax=dd->PtProng(i);
209 vars[iter]=dd->CosPointingAngle();
213 vars[iter]=dd->Getd0Prong(0)*dd->Getd0Prong(0)+dd->Getd0Prong(1)*dd->Getd0Prong(1)+dd->Getd0Prong(2)*dd->Getd0Prong(2);
218 for(Int_t iprong=0;iprong<3;iprong++){
219 if(dd->GetDCA(iprong)<maxDCA){
220 maxDCA=dd->GetDCA(iprong);
227 //---------------------------------------------------------------------------
228 Bool_t AliRDHFCutsDplustoKpipi::IsInFiducialAcceptance(Double_t pt, Double_t y) const
231 // Checking if Dplus is in fiducial acceptance region
235 // applying cut for pt > 5 GeV
236 AliDebug(2,Form("pt of D+ = %f (> 5), cutting at |y| < 0.8",pt));
237 if (TMath::Abs(y) > 0.8) return kFALSE;
240 // appliying smooth cut for pt < 5 GeV
241 Double_t maxFiducialY = -0.2/15*pt*pt+1.9/15*pt+0.5;
242 Double_t minFiducialY = 0.2/15*pt*pt-1.9/15*pt-0.5;
243 AliDebug(2,Form("pt of D+ = %f (< 5), cutting according to the fiducial zone [%f, %f]\n",pt,minFiducialY,maxFiducialY));
244 if (y < minFiducialY || y > maxFiducialY) return kFALSE;
250 //---------------------------------------------------------------------------
251 Int_t AliRDHFCutsDplustoKpipi::IsSelectedPID(AliAODRecoDecayHF *rd)
254 // PID selection, returns 3 if accepted, 0 if not accepted
256 if(!fUsePID || !rd) return 3;
257 //if(fUsePID)printf("i am inside the pid \n");
260 Int_t sign= rd->GetCharge();
261 for(Int_t daught=0;daught<3;daught++){
262 AliAODTrack *track=(AliAODTrack*)rd->GetDaughter(daught);
263 Int_t isPion=fPidHF->MakeRawPid(track,AliPID::kPion);
264 Int_t isKaon=fPidHF->MakeRawPid(track,AliPID::kKaon);
265 Int_t isProton=fPidHF->MakeRawPid(track,AliPID::kProton);
267 if(isProton>0 && isKaon<0 && isPion<0) return 0;
268 if(isKaon>0 && isPion<0) nkaons++;
269 if(isKaon<0) nNotKaons++;
270 if(sign==track->Charge()){//pions
271 if(isPion<0)return 0;
274 if(isKaon<0)return 0;
280 if(nkaons>1)return 0;
281 if(nNotKaons==3)return 0;
288 //---------------------------------------------------------------------------
289 Int_t AliRDHFCutsDplustoKpipi::IsSelected(TObject* obj,Int_t selectionLevel, AliAODEvent* aod) {
291 // Apply selection, returns 3 if accepted, 0 if not accepted
295 cout<<"Cut matrix not inizialized. Exit..."<<endl;
299 AliAODRecoDecayHF3Prong* d=(AliAODRecoDecayHF3Prong*)obj;
303 cout<<"AliAODRecoDecayHF3Prong null"<<endl;
309 // selection on daughter tracks
310 if(selectionLevel==AliRDHFCuts::kAll ||
311 selectionLevel==AliRDHFCuts::kTracks) {
312 if(!AreDaughtersSelected(d)) return 0;
316 Int_t returnvaluePID=3;
317 Int_t returnvalueCuts=3;
320 //if(selectionLevel==AliRDHFCuts::kAll ||
321 if(selectionLevel==AliRDHFCuts::kCandidate ||
322 selectionLevel==AliRDHFCuts::kPID) {
323 returnvaluePID = IsSelectedPID(d);
325 if(returnvaluePID==0)return 0;
328 // selection on candidate
329 if(selectionLevel==AliRDHFCuts::kAll ||
330 selectionLevel==AliRDHFCuts::kCandidate) {
332 //recalculate vertex w/o daughters
333 AliAODVertex *origownvtx=0x0;
334 AliAODVertex *recvtx=0x0;
335 if(fRemoveDaughtersFromPrimary) {
337 AliError("Can not remove daughters from vertex without AOD event");
340 if(d->GetOwnPrimaryVtx()) origownvtx=new AliAODVertex(*d->GetOwnPrimaryVtx());
341 recvtx=d->RemoveDaughtersFromPrimaryVtx(aod);
343 AliDebug(2,"Removal of daughter tracks failed");
344 //recvtx=d->GetPrimaryVtx();
351 //set recalculed primary vertex
352 d->SetOwnPrimaryVtx(recvtx);
353 delete recvtx; recvtx=NULL;
358 Int_t ptbin=PtBin(pt);
361 d->SetOwnPrimaryVtx(origownvtx);
365 else d->UnsetOwnPrimaryVtx();
369 Double_t mDplusPDG = TDatabasePDG::Instance()->GetParticle(411)->Mass();
370 Double_t mDplus=d->InvMassDplus();
371 if(TMath::Abs(mDplus-mDplusPDG)>fCutsRD[GetGlobalIndex(0,ptbin)])returnvalueCuts=0;
372 // if(d->PtProng(1) < fCutsRD[GetGlobalIndex(3,ptbin)] || d->PtProng(0) < fCutsRD[GetGlobalIndex(4,ptbin)]) okD0 = 0;
373 if(TMath::Abs(d->PtProng(1)) < fCutsRD[GetGlobalIndex(1,ptbin)] || TMath::Abs(d->Getd0Prong(1))<fCutsRD[GetGlobalIndex(3,ptbin)])returnvalueCuts=0;//Kaon
374 if(TMath::Abs(d->PtProng(0)) < fCutsRD[GetGlobalIndex(2,ptbin)] || TMath::Abs(d->Getd0Prong(0))<fCutsRD[GetGlobalIndex(4,ptbin)])returnvalueCuts=0;//Pion1
375 if(TMath::Abs(d->PtProng(2)) < fCutsRD[GetGlobalIndex(2,ptbin)] || TMath::Abs(d->Getd0Prong(2))<fCutsRD[GetGlobalIndex(4,ptbin)])returnvalueCuts=0;//Pion2
380 if(d->GetDist12toPrim()<fCutsRD[GetGlobalIndex(5,ptbin)]|| d->GetDist23toPrim()<fCutsRD[GetGlobalIndex(5,ptbin)])returnvalueCuts=0;
381 if(d->Getd0Prong(0)*d->Getd0Prong(1)<0. && d->Getd0Prong(2)*d->Getd0Prong(1)<0.)returnvalueCuts=0;
384 if(d->GetSigmaVert()>fCutsRD[GetGlobalIndex(6,ptbin)])returnvalueCuts=0;
386 if(d->DecayLength()<fCutsRD[GetGlobalIndex(7,ptbin)])returnvalueCuts=0;
388 if(TMath::Abs(d->PtProng(0))<fCutsRD[GetGlobalIndex(8,ptbin)] && TMath::Abs(d->PtProng(1))<fCutsRD[GetGlobalIndex(8,ptbin)] && TMath::Abs(d->PtProng(2))<fCutsRD[GetGlobalIndex(8,ptbin)])returnvalueCuts=0;
389 if(d->CosPointingAngle()< fCutsRD[GetGlobalIndex(9,ptbin)])returnvalueCuts=0;
390 Double_t sum2=d->Getd0Prong(0)*d->Getd0Prong(0)+d->Getd0Prong(1)*d->Getd0Prong(1)+d->Getd0Prong(2)*d->Getd0Prong(2);
391 if(sum2<fCutsRD[GetGlobalIndex(10,ptbin)])returnvalueCuts=0;
394 for(Int_t i=0;i<3;i++) if(d->GetDCA(i)>fCutsRD[GetGlobalIndex(11,ptbin)]) returnvalueCuts=0;
395 // unset recalculated primary vertex when not needed any more
397 d->SetOwnPrimaryVtx(origownvtx);
400 } else if(fRemoveDaughtersFromPrimary) {
401 d->UnsetOwnPrimaryVtx();
402 AliDebug(3,"delete new vertex\n");
405 return returnvalueCuts;
410 //---------------------------------------------------------------------------