]> git.uio.no Git - u/mrichter/AliRoot.git/blob - PWG3/vertexingHF/AliRDHFCutsDplustoKpipi.cxx
Add selection on track filter mask in D meson tasks
[u/mrichter/AliRoot.git] / PWG3 / vertexingHF / AliRDHFCutsDplustoKpipi.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-2010, 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
16 /* $Id$ */
17
18 /////////////////////////////////////////////////////////////
19 //
20 // Class for cuts on AOD reconstructed D+->Kpipi
21 //
22 // Author: R. Bala, bala@to.infn.it
23 //         G. Ortona, ortona@to.infn.it
24 /////////////////////////////////////////////////////////////
25
26 #include <TDatabasePDG.h>
27 #include <Riostream.h>
28
29 #include "AliRDHFCutsDplustoKpipi.h"
30 #include "AliAODRecoDecayHF3Prong.h"
31 #include "AliAODTrack.h"
32 #include "AliESDtrack.h"
33
34
35 ClassImp(AliRDHFCutsDplustoKpipi)
36
37 //--------------------------------------------------------------------------
38 AliRDHFCutsDplustoKpipi::AliRDHFCutsDplustoKpipi(const char* name) : 
39 AliRDHFCuts(name),
40   fUseStrongPid(0),
41   fMaxPtStrongPid(0.),
42   fUseImpParProdCorrCut(kFALSE)
43 {
44   //
45   // Default Constructor
46   //
47   Int_t nvars=14;
48   SetNVars(nvars);
49   TString varNames[14]={"inv. mass [GeV]",
50                         "pTK [GeV/c]",
51                         "pTPi [GeV/c]",
52                         "d0K [cm]   lower limit!",
53                         "d0Pi [cm]  lower limit!",
54                         "dist12 (cm)",
55                         "sigmavert (cm)",
56                         "dist prim-sec (cm)",
57                         "pM=Max{pT1,pT2,pT3} (GeV/c)",
58                         "cosThetaPoint",
59                         "Sum d0^2 (cm^2)",
60                         "dca cut (cm)",
61                         "dec len XY (cm)",
62                         "cosThetaPointXY"};
63   Bool_t isUpperCut[14]={kTRUE,
64                          kFALSE,
65                          kFALSE,
66                          kFALSE,
67                          kFALSE,
68                          kFALSE,
69                          kTRUE,
70                          kFALSE,
71                          kFALSE,
72                          kFALSE,
73                          kFALSE,
74                          kTRUE,
75                          kFALSE,
76                          kFALSE};
77   SetVarNames(nvars,varNames,isUpperCut);
78   Bool_t forOpt[14]={kFALSE,
79                      kFALSE,
80                      kFALSE,
81                      kFALSE,
82                      kFALSE,
83                      kFALSE,
84                      kTRUE,
85                      kTRUE,
86                      kTRUE,
87                      kTRUE,
88                      kTRUE,
89                      kFALSE,
90                      kTRUE,
91                      kTRUE};
92   SetVarsForOpt(7,forOpt);
93   Float_t limits[2]={0,999999999.};
94   SetPtBins(2,limits);
95   if(fPidHF)delete fPidHF;
96   fPidHF=new AliAODPidHF();
97   Double_t plim[2]={0.6,0.8};
98   Double_t nsigma[5]={2.,1.,2.,3.,0.};
99   
100   fPidHF->SetPLimit(plim);
101   fPidHF->SetAsym(kTRUE);
102   fPidHF->SetSigma(nsigma);
103   fPidHF->SetMatch(1);
104   fPidHF->SetTPC(1);
105   fPidHF->SetTOF(1);
106   fPidHF->SetITS(0);
107   fPidHF->SetTRD(0);
108   fPidHF->SetCompat(kTRUE);
109
110   
111 }
112
113
114
115
116
117
118
119
120 //--------------------------------------------------------------------------
121 AliRDHFCutsDplustoKpipi::AliRDHFCutsDplustoKpipi(const AliRDHFCutsDplustoKpipi &source) :
122   AliRDHFCuts(source),
123   fUseStrongPid(source.fUseStrongPid),
124   fMaxPtStrongPid(source.fMaxPtStrongPid),
125   fUseImpParProdCorrCut(source.fUseImpParProdCorrCut)
126 {
127   //
128   // Copy constructor
129   //
130
131 }
132 //--------------------------------------------------------------------------
133 AliRDHFCutsDplustoKpipi &AliRDHFCutsDplustoKpipi::operator=(const AliRDHFCutsDplustoKpipi &source)
134 {
135   //
136   // assignment operator
137   //
138   if(&source == this) return *this;
139
140   AliRDHFCuts::operator=(source);
141
142   return *this;
143 }
144 //
145
146
147 //---------------------------------------------------------------------------
148 void AliRDHFCutsDplustoKpipi::GetCutVarsForOpt(AliAODRecoDecayHF *d,Float_t *vars,Int_t nvars,Int_t *pdgdaughters,AliAODEvent *aod) {
149   // 
150   // Fills in vars the values of the variables 
151   //
152
153
154   if(nvars!=fnVarsForOpt) {
155     printf("AliRDHFCutsDplustoKpipi::GetCutsVarsForOpt: wrong number of variables\n");
156     return;
157   }
158
159   AliAODRecoDecayHF3Prong *dd = (AliAODRecoDecayHF3Prong*)d;
160  
161   //recalculate vertex w/o daughters
162   Bool_t cleanvtx=kFALSE;
163   AliAODVertex *origownvtx=0x0;
164   if(fRemoveDaughtersFromPrimary) {
165     if(dd->GetOwnPrimaryVtx()) origownvtx=new AliAODVertex(*dd->GetOwnPrimaryVtx());
166     cleanvtx=kTRUE;
167     if(!RecalcOwnPrimaryVtx(dd,aod)) {
168       CleanOwnPrimaryVtx(dd,aod,origownvtx);
169       cleanvtx=kFALSE;
170     }
171   }
172
173   Int_t iter=-1;
174   if(fVarsForOpt[0]){
175     iter++;
176     vars[iter]=dd->InvMassDplus();
177   }
178   if(fVarsForOpt[1]){
179     iter++;
180     for(Int_t iprong=0;iprong<3;iprong++){
181       if(TMath::Abs(pdgdaughters[iprong])==321) {
182         vars[iter]=dd->PtProng(iprong);
183       }
184     }
185   }
186   if(fVarsForOpt[2]){
187     iter++;
188     Float_t minPtDau=1000000.0;
189     for(Int_t iprong=0;iprong<3;iprong++){
190       if(TMath::Abs(pdgdaughters[iprong])==211) {
191         if(dd->PtProng(iprong)<minPtDau){
192           minPtDau=dd->PtProng(iprong);
193         }
194       }
195     }
196     vars[iter]=minPtDau;
197   }
198   if(fVarsForOpt[3]){
199     iter++;
200     for(Int_t iprong=0;iprong<3;iprong++){
201       if(TMath::Abs(pdgdaughters[iprong])==321) {
202         vars[iter]=dd->Getd0Prong(iprong);
203       }
204     }
205   }
206   if(fVarsForOpt[4]){
207     iter++;
208     Float_t minImpParDau=1000000.0;
209     for(Int_t iprong=0;iprong<3;iprong++){
210       if(TMath::Abs(pdgdaughters[iprong])==211) {
211         if(dd->Getd0Prong(iprong)<minImpParDau){
212           minImpParDau=dd->Getd0Prong(iprong);
213         }
214       }
215     }
216    vars[iter]=minImpParDau;
217   }
218   if(fVarsForOpt[5]){
219     iter++;
220     Float_t dist12 = dd->GetDist12toPrim();
221     Float_t dist23 = dd->GetDist23toPrim();
222     if(dist12<dist23)vars[iter]=dist12;
223     else vars[iter]=dist23;
224   }
225   if(fVarsForOpt[6]){
226     iter++;
227     vars[iter]=dd->GetSigmaVert();
228   }
229   if(fVarsForOpt[7]){
230     iter++;
231     vars[iter] = dd->DecayLength();
232   }
233   if(fVarsForOpt[8]){
234     iter++;
235     Float_t ptmax=0;
236     for(Int_t i=0;i<3;i++){
237       if(dd->PtProng(i)>ptmax)ptmax=dd->PtProng(i);
238     }
239     vars[iter]=ptmax;
240   }
241   if(fVarsForOpt[9]){
242     iter++;
243     vars[iter]=dd->CosPointingAngle();
244   }
245   if(fVarsForOpt[10]){
246     iter++;
247     vars[iter]=dd->Getd0Prong(0)*dd->Getd0Prong(0)+dd->Getd0Prong(1)*dd->Getd0Prong(1)+dd->Getd0Prong(2)*dd->Getd0Prong(2);
248   }
249   if(fVarsForOpt[11]){
250     iter++;
251     Float_t maxDCA=0;
252     for(Int_t iprong=0;iprong<3;iprong++){
253       if(dd->GetDCA(iprong)<maxDCA){
254         maxDCA=dd->GetDCA(iprong);
255       }
256     }
257     vars[iter]=maxDCA;
258   }
259   if(fVarsForOpt[12]){
260     iter++;
261     vars[iter]=dd->NormalizedDecayLengthXY()*dd->P()/dd->Pt();
262   }
263   if(fVarsForOpt[13]){
264     iter++;
265     vars[iter]=dd->CosPointingAngleXY();
266   }
267
268   if(cleanvtx)CleanOwnPrimaryVtx(dd,aod,origownvtx);
269
270   return;
271 }
272 //---------------------------------------------------------------------------
273 Bool_t AliRDHFCutsDplustoKpipi::IsInFiducialAcceptance(Double_t pt, Double_t y) const
274 {
275   //
276   // Checking if Dplus is in fiducial acceptance region 
277   //
278
279   if(pt > 5.) {
280     // applying cut for pt > 5 GeV
281     AliDebug(2,Form("pt of D+ = %f (> 5), cutting at |y| < 0.8",pt)); 
282     if (TMath::Abs(y) > 0.8) return kFALSE;
283     
284   } else {
285     // appliying smooth cut for pt < 5 GeV
286     Double_t maxFiducialY = -0.2/15*pt*pt+1.9/15*pt+0.5; 
287     Double_t minFiducialY = 0.2/15*pt*pt-1.9/15*pt-0.5;         
288     AliDebug(2,Form("pt of D+ = %f (< 5), cutting  according to the fiducial zone [%f, %f]\n",pt,minFiducialY,maxFiducialY)); 
289     if (y < minFiducialY || y > maxFiducialY) return kFALSE;    
290   }
291
292   return kTRUE;
293 }
294
295 //---------------------------------------------------------------------------
296 Int_t AliRDHFCutsDplustoKpipi::IsSelectedPID(AliAODRecoDecayHF *rd)
297 {
298   //
299   // PID selection, returns 3 if accepted, 0 if not accepted
300   // 
301   if(!fUsePID || !rd) return 3;
302   //if(fUsePID)printf("i am inside the pid \n");
303   Int_t nkaons=0;
304   Int_t nNotKaons=0;
305   Int_t sign= rd->GetCharge(); 
306   for(Int_t daught=0;daught<3;daught++){
307     AliAODTrack *track=(AliAODTrack*)rd->GetDaughter(daught);
308     Int_t isPion=fPidHF->MakeRawPid(track,AliPID::kPion);
309     Int_t isKaon=fPidHF->MakeRawPid(track,AliPID::kKaon);
310     Int_t isProton=fPidHF->MakeRawPid(track,AliPID::kProton);
311     
312     if(isProton>0 &&  isKaon<0  && isPion<0) return 0;
313     if(isKaon>0 && isPion<0) nkaons++;
314     if(isKaon<0) nNotKaons++;  
315     if(sign==track->Charge()){//pions
316       if(isPion<0)return 0;
317       if(rd->Pt()<fMaxPtStrongPid && isPion<=0 && fUseStrongPid>1)return 0;
318     }
319     else{//kaons
320       if(isKaon<0)return 0;
321         if(rd->Pt()<fMaxPtStrongPid && isKaon<=0 && fUseStrongPid>0)return 0;
322     }
323   }
324   
325   if(nkaons>1)return 0;
326   if(nNotKaons==3)return 0;
327   
328   return 3;   
329 }
330
331
332 //---------------------------------------------------------------------------
333 Int_t AliRDHFCutsDplustoKpipi::IsSelected(TObject* obj,Int_t selectionLevel, AliAODEvent* aod) {
334   //
335   // Apply selection, returns 3 if accepted, 0 if not accepted
336   //
337
338
339   fIsSelectedCuts=0;
340   fIsSelectedPID=0;
341
342   if(!fCutsRD){
343     cout<<"Cut matrix not inizialized. Exit..."<<endl;
344     return 0;
345   }
346   //PrintAll();
347   AliAODRecoDecayHF3Prong* d=(AliAODRecoDecayHF3Prong*)obj; 
348
349   
350   if(!d){
351     cout<<"AliAODRecoDecayHF3Prong null"<<endl;
352     return 0;
353   }
354
355   if(fKeepSignalMC) if(IsSignalMC(d,aod,411)) return 3;
356
357   // PID selection
358   Int_t returnvaluePID=3;
359   Int_t returnvalueCuts=3;
360
361   Double_t pt=d->Pt();
362   if(pt<fMinPtCand) return 0;
363   if(pt>fMaxPtCand) return 0;
364
365   if(d->HasBadDaughters()) return 0;
366   
367   // selection on candidate
368   if(selectionLevel==AliRDHFCuts::kAll || 
369      selectionLevel==AliRDHFCuts::kCandidate) {
370     
371     //recalculate vertex w/o daughters
372     AliAODVertex *origownvtx=0x0;
373     if(fRemoveDaughtersFromPrimary) {
374       if(d->GetOwnPrimaryVtx()) origownvtx=new AliAODVertex(*d->GetOwnPrimaryVtx());
375       if(!RecalcOwnPrimaryVtx(d,aod)) {
376         CleanOwnPrimaryVtx(d,aod,origownvtx);
377         return 0;
378       }
379     }
380
381     Int_t ptbin=PtBin(pt);
382     if (ptbin==-1) {
383       CleanOwnPrimaryVtx(d,aod,origownvtx);
384       return 0;
385     }
386     
387     Double_t mDplusPDG = TDatabasePDG::Instance()->GetParticle(411)->Mass();
388     Double_t mDplus=d->InvMassDplus();
389     if(TMath::Abs(mDplus-mDplusPDG)>fCutsRD[GetGlobalIndex(0,ptbin)]) {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}
390
391     //2track cuts
392     if(d->GetDist12toPrim()<fCutsRD[GetGlobalIndex(5,ptbin)]|| d->GetDist23toPrim()<fCutsRD[GetGlobalIndex(5,ptbin)]) {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}
393
394     Double_t sum2=d->Getd0Prong(0)*d->Getd0Prong(0)+d->Getd0Prong(1)*d->Getd0Prong(1)+d->Getd0Prong(2)*d->Getd0Prong(2);
395     if(sum2<fCutsRD[GetGlobalIndex(10,ptbin)]) {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}
396
397     if(fUseImpParProdCorrCut){
398       if(d->Getd0Prong(0)*d->Getd0Prong(1)<0. && d->Getd0Prong(2)*d->Getd0Prong(1)<0.) {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}
399     }
400
401     //sec vert
402     if(d->GetSigmaVert()>fCutsRD[GetGlobalIndex(6,ptbin)]) {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}
403
404     //DCA
405     for(Int_t i=0;i<3;i++) if(d->GetDCA(i)>fCutsRD[GetGlobalIndex(11,ptbin)]) {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}
406
407     if(d->Pt2Prong(1) < fCutsRD[GetGlobalIndex(1,ptbin)]*fCutsRD[GetGlobalIndex(1,ptbin)] || TMath::Abs(d->Getd0Prong(1))<fCutsRD[GetGlobalIndex(3,ptbin)]) {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}//Kaon
408
409     if(d->Pt2Prong(0) < fCutsRD[GetGlobalIndex(2,ptbin)]*fCutsRD[GetGlobalIndex(2,ptbin)] || TMath::Abs(d->Getd0Prong(0))<fCutsRD[GetGlobalIndex(4,ptbin)]) {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}//Pion1
410
411     if(d->Pt2Prong(2) < fCutsRD[GetGlobalIndex(2,ptbin)]*fCutsRD[GetGlobalIndex(2,ptbin)] || TMath::Abs(d->Getd0Prong(2))<fCutsRD[GetGlobalIndex(4,ptbin)]) {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}//Pion2
412     
413     if(d->Pt2Prong(0)<fCutsRD[GetGlobalIndex(8,ptbin)]*fCutsRD[GetGlobalIndex(8,ptbin)] && d->Pt2Prong(1)<fCutsRD[GetGlobalIndex(8,ptbin)]*fCutsRD[GetGlobalIndex(8,ptbin)] && d->Pt2Prong(2)<fCutsRD[GetGlobalIndex(8,ptbin)]*fCutsRD[GetGlobalIndex(8,ptbin)]) {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}
414
415     if(d->DecayLength2()<fCutsRD[GetGlobalIndex(7,ptbin)]*fCutsRD[GetGlobalIndex(7,ptbin)]) {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}
416
417     if(d->CosPointingAngle()< fCutsRD[GetGlobalIndex(9,ptbin)]) {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}
418
419     if(d->NormalizedDecayLengthXY()*d->P()/pt<fCutsRD[GetGlobalIndex(12,ptbin)]){CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}
420
421     if(d->CosPointingAngleXY()<fCutsRD[GetGlobalIndex(13,ptbin)]){CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}
422
423     // unset recalculated primary vertex when not needed any more
424     CleanOwnPrimaryVtx(d,aod,origownvtx);
425     
426     fIsSelectedCuts=returnvalueCuts;
427
428     //if(!returnvalueCuts) return 0; // returnvalueCuts cannot be 0 here
429   }
430   
431   if(selectionLevel==AliRDHFCuts::kAll || 
432      selectionLevel==AliRDHFCuts::kCandidate ||     
433      selectionLevel==AliRDHFCuts::kPID) {
434     returnvaluePID = IsSelectedPID(d);
435     fIsSelectedPID=returnvaluePID;
436   }
437   if(returnvaluePID==0)return 0;
438
439   // selection on daughter tracks 
440   if(selectionLevel==AliRDHFCuts::kAll || 
441      selectionLevel==AliRDHFCuts::kTracks) {
442     if(!AreDaughtersSelected(d)) return 0;
443   }
444   
445  
446
447
448   return 3;
449 }
450
451
452
453
454 //---------------------------------------------------------------------------
455
456
457 void AliRDHFCutsDplustoKpipi::SetStandardCutsPP2010() {
458   //
459   //STANDARD CUTS USED FOR 2010 pp analysis 
460   //                                            
461   
462   SetName("DplustoKpipiCutsStandard");
463   SetTitle("Standard Cuts for D+ analysis");
464   
465   // PILE UP REJECTION
466   SetOptPileup(AliRDHFCuts::kRejectPileupEvent);
467
468   // EVENT CUTS
469   SetMinVtxContr(1);
470
471   AliESDtrackCuts* esdTrackCuts=new AliESDtrackCuts();
472   esdTrackCuts->SetRequireSigmaToVertex(kFALSE);
473   //default
474   esdTrackCuts->SetRequireTPCRefit(kTRUE);
475   esdTrackCuts->SetRequireITSRefit(kTRUE);
476   //esdTrackCuts->SetMinNClustersITS(4); // default is 5
477   esdTrackCuts->SetMinNClustersTPC(70);
478   esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,
479                                          AliESDtrackCuts::kAny); 
480   // default is kBoth, otherwise kAny
481   esdTrackCuts->SetMinDCAToVertexXY(0.);
482   esdTrackCuts->SetPtRange(0.3,1.e10);
483   
484   AddTrackCuts(esdTrackCuts);
485   
486  
487   const Int_t nptbins =13;
488   const Int_t nvars=14;
489   Float_t ptbins[nptbins+1];
490   ptbins[0]=0.;
491   ptbins[1]=1;  
492   ptbins[2]=2.;
493   ptbins[3]=3.;
494   ptbins[4]=4.;
495   ptbins[5]=5.;
496   ptbins[6]=6.;
497   ptbins[7]=8.;
498   ptbins[8]=10.;
499   ptbins[9]=12.;
500   ptbins[10]=14.;
501   ptbins[11]=16.;
502   ptbins[12]=24.;
503   ptbins[13]=99999.;
504       
505     
506   Float_t** anacutsval;
507   anacutsval=new Float_t*[nvars];
508   
509   for(Int_t ic=0;ic<nvars;ic++){anacutsval[ic]=new Float_t[nptbins];}
510
511   //Double_t cutsDplus[12]={0.2,0.4,0.4,0.,0.,0.01,0.06,0.02,0.,0.85,0.,10000000000.};
512   for(Int_t ipt=0;ipt<nptbins;ipt++){
513     anacutsval[0][ipt]=0.2;
514     anacutsval[3][ipt]=0.;
515     anacutsval[4][ipt]=0.;
516     anacutsval[5][ipt]=0.01;
517     anacutsval[11][ipt]=10000000000.;
518     anacutsval[12][ipt]=0.;
519     anacutsval[13][ipt]=0.;
520   }
521
522   anacutsval[1][0]=0.3;
523   anacutsval[1][1]=0.3;
524   anacutsval[1][2]=0.4; 
525   anacutsval[2][0]=0.3;
526   anacutsval[2][1]=0.3;
527   anacutsval[2][2]=0.4;  
528   for(Int_t ipt=3;ipt<nptbins;ipt++){
529     anacutsval[1][ipt]=0.4;
530     anacutsval[2][ipt]=0.4;
531   }
532   
533   anacutsval[6][0]=0.022100;
534   anacutsval[6][1]=0.022100;
535   anacutsval[6][2]=0.034;
536   anacutsval[6][3]=0.020667;
537   anacutsval[6][4]=0.020667;
538   anacutsval[6][5]=0.023333;
539     
540   
541   anacutsval[7][0]=0.08;
542   anacutsval[7][1]=0.08;
543   anacutsval[7][2]=0.09;  
544   anacutsval[7][3]=0.095;
545   anacutsval[7][4]=0.095;
546    
547   anacutsval[8][0]=0.5;
548   anacutsval[8][1]=0.5;
549   anacutsval[8][2]=1.0;
550   anacutsval[8][3]=0.5;
551   anacutsval[8][4]=0.5;
552      
553     
554   anacutsval[9][0]=0.95;
555   anacutsval[9][1]=0.95;
556   anacutsval[9][2]=0.95; 
557   anacutsval[9][3]=0.95; 
558   anacutsval[9][4]= 0.95;
559   anacutsval[9][5]=0.92;
560   anacutsval[9][6]=0.92;
561   anacutsval[9][7]=0.92;
562   anacutsval[9][8]=0.92;
563   anacutsval[9][9]=0.90;
564   anacutsval[9][10]=0.90;
565   anacutsval[9][11]=0.90;
566   anacutsval[9][12]=0.90;
567   
568   anacutsval[10][0]=0.0055;
569   anacutsval[10][1]=0.0055;
570   anacutsval[10][2]= 0.0028;
571   anacutsval[10][3]=0.000883;
572   anacutsval[10][4]=0.000883;
573
574   
575   for(Int_t ipt=5;ipt<nptbins;ipt++){
576     anacutsval[6][ipt]=0.02333;
577     anacutsval[7][ipt]=0.115;
578     anacutsval[8][ipt]=0.5;
579     anacutsval[10][ipt]=0.000883;
580     }   
581
582
583   SetGlobalIndex(nvars,nptbins);
584   SetPtBins(nptbins+1,ptbins);
585   SetCuts(nvars,nptbins,anacutsval);
586   SetUsePID(kTRUE);
587   SetRemoveDaughtersFromPrim(kTRUE);
588   
589   PrintAll();
590
591   for(Int_t iic=0;iic<nvars;iic++){delete [] anacutsval[iic];}
592   delete [] anacutsval;
593   anacutsval=NULL;
594
595   return;
596 }
597
598
599 void AliRDHFCutsDplustoKpipi::SetStandardCutsPbPb2010() {
600   //
601   //STANDARD CUTS USED FOR 2010 Pb Pb analysis.... not optimized yet 
602   //                                              
603   
604   SetName("DplustoKpipiCutsStandard");
605   SetTitle("Standard Cuts for D+ analysis in PbPb2010 run");
606   
607   // PILE UP REJECTION
608   //SetOptPileup(AliRDHFCuts::kRejectPileupEvent);
609
610   // EVENT CUTS
611   SetMinVtxContr(1);
612
613   
614   AliESDtrackCuts* esdTrackCuts=new AliESDtrackCuts();
615   esdTrackCuts->SetRequireSigmaToVertex(kFALSE);
616   //default
617   esdTrackCuts->SetRequireTPCRefit(kTRUE);
618   esdTrackCuts->SetRequireITSRefit(kTRUE);
619   //esdTrackCuts->SetMinNClustersITS(4); // default is 5
620   esdTrackCuts->SetMinNClustersTPC(70);
621   esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,
622                                          AliESDtrackCuts::kAny); 
623   // default is kBoth, otherwise kAny
624   esdTrackCuts->SetMinDCAToVertexXY(0.);
625   esdTrackCuts->SetPtRange(0.8,1.e10);
626     
627   AddTrackCuts(esdTrackCuts);
628     
629   const Int_t nptbins=10;
630   Float_t* ptbins;
631   ptbins=new Float_t[nptbins+1];
632
633   ptbins[0]=0.;
634   ptbins[1]=1.;
635   ptbins[2]=2.;
636   ptbins[3]=3.;
637   ptbins[4]=4.;
638   ptbins[5]=5.;
639   ptbins[6]=6.;
640   ptbins[7]=8.;
641   ptbins[8]=12.;
642   ptbins[9]=16.;
643   ptbins[10]=24.;
644   const Int_t nvars=14;
645     
646   Float_t** anacutsval;
647   anacutsval=new Float_t*[nvars];
648   
649   for(Int_t ic=0;ic<nvars;ic++){anacutsval[ic]=new Float_t[nptbins];}
650   //Double_t cutsDplus[12]={0.2,0.4,0.4,0.,0.,0.01,0.06,0.02,0.,0.85,0.,10000000000.};
651
652   for(Int_t ipt=0;ipt<nptbins;ipt++){
653     anacutsval[0][ipt]=0.2;
654     anacutsval[1][ipt]=0.8;
655     anacutsval[2][ipt]=0.8;
656     anacutsval[3][ipt]=0.;
657     anacutsval[4][ipt]=0.;
658     anacutsval[5][ipt]=0.01;
659     anacutsval[11][ipt]=10000000000.;
660     anacutsval[12][ipt]=0.;
661     anacutsval[13][ipt]=0.;
662   }
663   anacutsval[1][5]=0.9;
664
665   anacutsval[6][0]=0.022100;
666   anacutsval[6][1]=0.022100;
667   anacutsval[6][2]=0.034;
668   anacutsval[6][3]=0.020667;
669   anacutsval[6][4]=0.020667;
670   anacutsval[6][5]=0.023333;
671   
672   anacutsval[7][0]=0.08;
673   anacutsval[7][1]=0.08;
674   anacutsval[7][2]=0.17;  
675   anacutsval[7][3]=0.14;
676   anacutsval[7][4]=0.14;
677   anacutsval[7][5]=0.19;
678   
679   anacutsval[8][0]=0.8;
680   anacutsval[8][1]=0.8;
681   anacutsval[8][2]=1.1;
682   anacutsval[8][3]=0.5;
683   anacutsval[8][4]=0.5;
684   anacutsval[8][5]=0.5;
685   
686   anacutsval[9][0]=0.995;
687   anacutsval[9][1]=0.995;
688   anacutsval[9][2]=0.997; 
689   anacutsval[9][3]=0.998; 
690   anacutsval[9][4]=0.998;
691   anacutsval[9][5]=0.995;
692   
693   anacutsval[10][0]=0.0055;
694   anacutsval[10][1]=0.0055;
695   anacutsval[10][2]= 0.0028;
696   anacutsval[10][3]=0.000883;
697   anacutsval[10][4]=0.000883;
698   anacutsval[10][5]=0.000883;
699
700   anacutsval[12][5]=12.;
701   anacutsval[13][5]=0.998571;
702   anacutsval[12][6]=10.;
703   anacutsval[13][6]=0.997143;
704
705   for(Int_t ipt=6;ipt<nptbins;ipt++){
706     anacutsval[6][ipt]=0.02333;
707     anacutsval[7][ipt]=0.19;
708     anacutsval[8][ipt]=2.0;
709     anacutsval[9][ipt]=0.997;
710     anacutsval[10][ipt]=0.000883;
711   }   
712   anacutsval[7][6]=0.14;
713   anacutsval[9][6]=0.995;
714
715   SetPtBins(nptbins+1,ptbins);
716   SetCuts(nvars,nptbins,anacutsval);
717   SetUsePID(kTRUE);
718   SetMinCentrality(1E-10);
719   SetMaxCentrality(20.);
720   SetUseCentrality(AliRDHFCuts::kCentV0M);
721   SetRemoveDaughtersFromPrim(kFALSE);
722     
723   PrintAll();
724
725   for(Int_t iic=0;iic<nvars;iic++){delete [] anacutsval[iic];}
726   delete [] anacutsval;
727   anacutsval=NULL;
728
729   delete [] ptbins;
730   ptbins=NULL;
731
732   return;
733 }