[u/mrichter/AliRoot.git] / PWGHF / vertexingHF / macros / makeTFile4CutsD0toKpi.C

#include <Riostream.h>
#include <TFile.h>
#include <AliRDHFCutsD0toKpi.h>
#include <TClonesArray.h>
#include <TParameter.h>


//Use:
//Set hard coded commentet with //set this!!
// root[] .L makeInputD0tasks.C++
// root[] makeInputAliAnalysisTaskSED0Mass()
// root[] makeInputAliAnalysisTaskSESignificanceMaximization()
//similar macros for the other D mesons

//Author: Chiara Bianchin, cbianchi@pd.infn.it


//macro to make a .root file which contains an AliRDHFCutsD0toKpi for AliAnalysisTaskSED0Mass task
void ModifyFromStandardCuts(Int_t system=1 /*0=pp, 1=PbPb, 2=pp 2.76TeV*/){
  AliRDHFCutsD0toKpi* RDHFD0toKpi=new AliRDHFCutsD0toKpi();
  TString info="",cent="";

  if(system==0) {
    RDHFD0toKpi->SetStandardCutsPP2010();
    info+="ppStd";
  }
  if(system==1) {
    RDHFD0toKpi->SetStandardCutsPbPb2010();
    //Change centrality if needed
    /*
    Float_t minc=20,maxc=80;
    RDHFD0toKpi->SetMinCentrality(minc);
    RDHFD0toKpi->SetMaxCentrality(maxc);
    */
    RDHFD0toKpi->SetUseCentrality(AliRDHFCuts::kCentOff); //kCentOff,kCentV0M,kCentTRK,kCentTKL,kCentCL1,kCentInvalid

    cent=Form("%.0f%.0f",RDHFD0toKpi->GetMinCentrality(),RDHFD0toKpi->GetMaxCentrality());
    info+="PbPbStd";

  }
  if(system==2) {
    RDHFD0toKpi->SetStandardCutsPP2011_276TeV();
    info+="pp276Std";
  }

  RDHFD0toKpi->SetName("D0toKpiCuts");
  RDHFD0toKpi->SetTitle("Cuts for D0 analysis");
  

  //here add what you need! Find examples below

  //Trigger mask
  /*
  RDHFD0toKpi->SetTriggerMask(0);
  RDHFD0toKpi->SetTriggerMask(AliVEvent::kEMC1 | AliVEvent::kEMC7);
  RDHFD0toKpi->SetTriggerClass("CEMC");

  info+="EMCTr";
  */

  //Event selection
  RDHFD0toKpi->SetUsePhysicsSelection(kFALSE);
  info+="noPhysSel";

  cout<<"This is the odject I'm going to save:"<<endl;
  cout<<info<<endl;
  RDHFD0toKpi->PrintAll();
  TFile* fout=new TFile(Form("D0toKpiCuts%s%sRecVtx%sPileupRej%s.root", RDHFD0toKpi->GetUseCentrality()==0 ? "" : cent.Data(),RDHFD0toKpi->GetIsPrimaryWithoutDaughters() ? "" : "No",RDHFD0toKpi->GetOptPileUp() ? "" : "No",info.Data()),"recreate");   //set this!! 

  fout->cd();
  RDHFD0toKpi->Write();
  fout->Close();

}

void makeInputAliAnalysisTaskSED0Mass(){

  AliRDHFCutsD0toKpi* RDHFD0toKpi=new AliRDHFCutsD0toKpi();
  RDHFD0toKpi->SetName("D0toKpiCuts");
  RDHFD0toKpi->SetTitle("Cuts for D0 analysis");

  AliESDtrackCuts* esdTrackCuts=new AliESDtrackCuts();
  esdTrackCuts->SetRequireSigmaToVertex(kFALSE);
  //default
  esdTrackCuts->SetRequireTPCRefit(kTRUE);
  esdTrackCuts->SetRequireITSRefit(kTRUE);
  esdTrackCuts->SetEtaRange(-0.8,0.8);
  //esdTrackCuts->SetMinNClustersITS(4); // default is 5
  //esdTrackCuts->SetMinNClustersTPC(120);
  esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,
					 AliESDtrackCuts::kAny); 
 // default is kBoth, otherwise kAny
  esdTrackCuts->SetMinDCAToVertexXY(0.);
  esdTrackCuts->SetPtRange(0.8,1.e10);


  RDHFD0toKpi->AddTrackCuts(esdTrackCuts);

  const Int_t nvars=11;

  const Int_t nptbins=13;
  Float_t* ptbins;
  ptbins=new Float_t[nptbins+1];
  ptbins[0]=0.;
  ptbins[1]=0.5;
  ptbins[2]=1.;
  ptbins[3]=2.;
  ptbins[4]=3.;
  ptbins[5]=4.;
  ptbins[6]=5.;
  ptbins[7]=6.;
  ptbins[8]=8.;
  ptbins[9]=12.;
  ptbins[10]=16.;
  ptbins[11]=20.;
  ptbins[12]=24.;
  ptbins[13]=9999.;

  RDHFD0toKpi->SetPtBins(nptbins+1,ptbins);
  

  Float_t** rdcutsvalmine;
  rdcutsvalmine=new Float_t*[nvars];
  for(Int_t iv=0;iv<nvars;iv++){
    rdcutsvalmine[iv]=new Float_t[nptbins];
  }

  //setting my cut values
    //cuts order
    //       printf("    |M-MD0| [GeV]    < %f\n",fD0toKpiCuts[0]);
    //     printf("    dca    [cm]  < %f\n",fD0toKpiCuts[1]);
    //     printf("    cosThetaStar     < %f\n",fD0toKpiCuts[2]);
    //     printf("    pTK     [GeV/c]    > %f\n",fD0toKpiCuts[3]);
    //     printf("    pTpi    [GeV/c]    > %f\n",fD0toKpiCuts[4]);
    //     printf("    |d0K|  [cm]  < %f\n",fD0toKpiCuts[5]);
    //     printf("    |d0pi| [cm]  < %f\n",fD0toKpiCuts[6]);
    //     printf("    d0d0  [cm^2] < %f\n",fD0toKpiCuts[7]);
    //     printf("    cosThetaPoint    > %f\n",fD0toKpiCuts[8]);
    //     printf("    |cosThetaPointXY| < %f\n",fD0toKpiCuts[9]);
    //     printf("    NormDecayLenghtXY    > %f\n",fD0toKpiCuts[10]);


  Double_t arrcuts[11]={0.3,0.03,0.8,0.8,0.8,0.1,0.1,-0.0004,0.9,0.998,5.}; //put the last 2 values at 0. for pp

  //setting my cut values
  //0-0.5 GeV
  for(Int_t ic=0;ic<nvars;ic++) rdcutsvalmine[ic][0]=arrcuts[ic];

  //0.5-1 GeV/c
  for(Int_t ic=0;ic<nvars;ic++) rdcutsvalmine[ic][1]=arrcuts[ic];

  //1-2 GeV 
  arrcuts[1]=0.025;  arrcuts[7]=-0.0003;
  for(Int_t ic=0;ic<nvars;ic++) rdcutsvalmine[ic][2]=arrcuts[ic];

  //2-3 GeV
  arrcuts[7]=-0.00026;
  arrcuts[9]=0.998;
  for(Int_t ic=0;ic<nvars;ic++) rdcutsvalmine[ic][3]=arrcuts[ic];

  //3-4 GeV
  arrcuts[7]=-0.00015; arrcuts[8]=0.85;
  for(Int_t ic=0;ic<nvars;ic++) rdcutsvalmine[ic][4]=arrcuts[ic];

  //4-5 GeV 
  for(Int_t ic=0;ic<nvars;ic++) rdcutsvalmine[ic][5]=arrcuts[ic];

  //5-6 GeV
  arrcuts[7]=-0.0001; arrcuts[8]=0.85;
  for(Int_t ic=0;ic<nvars;ic++) rdcutsvalmine[ic][6]=arrcuts[ic];

  //6-8 GeV
  arrcuts[2]=1.;
  for(Int_t ic=0;ic<nvars;ic++) rdcutsvalmine[ic][7]=arrcuts[ic];

  //8-12 GeV
  arrcuts[8]=0.8;
  for(Int_t ic=0;ic<nvars;ic++) rdcutsvalmine[ic][8]=arrcuts[ic];

  //12-16 GeV
  arrcuts[1]=0.03;
  for(Int_t ic=0;ic<nvars;ic++) rdcutsvalmine[ic][9]=arrcuts[ic];

  //16-20 GeV
  arrcuts[1]=0.035;
  for(Int_t ic=0;ic<nvars;ic++) rdcutsvalmine[ic][10]=arrcuts[ic];

  //20-24 GeV
  for(Int_t ic=0;ic<nvars;ic++) rdcutsvalmine[ic][11]=arrcuts[ic];

  //24-9999 GeV
  for(Int_t ic=0;ic<nvars;ic++) rdcutsvalmine[ic][12]=arrcuts[ic];

  RDHFD0toKpi->SetCuts(nvars,nptbins,rdcutsvalmine);

  Bool_t pidflag=kTRUE;
  RDHFD0toKpi->SetUsePID(pidflag);
  if(pidflag) cout<<"PID is used"<<endl;
  else cout<<"PID is not used"<<endl;

    //pid settings
  AliAODPidHF* pidObj=new AliAODPidHF();
  //pidObj->SetName("pid4D0");
  Int_t mode=1;
  const Int_t nlims=2;
  Double_t plims[nlims]={0.6,0.8}; //TPC limits in momentum [GeV/c]
  Bool_t compat=kTRUE; //effective only for this mode
  Bool_t asym=kTRUE;
  Double_t sigmas[5]={2.,1.,0.,3.,0.}; //to be checked and to be modified with new implementation of setters by Rossella
  pidObj->SetAsym(asym);// if you want to use the asymmetric bands in TPC
  pidObj->SetMatch(mode);
  pidObj->SetPLimit(plims,nlims);
  pidObj->SetSigma(sigmas);
  pidObj->SetCompat(compat);
  pidObj->SetTPC(kTRUE);
  pidObj->SetTOF(kTRUE);
  RDHFD0toKpi->SetPidHF(pidObj);

  RDHFD0toKpi->SetUseDefaultPID(kFALSE); //to use the AliAODPidHF

  //activate pileup rejection (for pp)
  //RDHFD0toKpi->SetOptPileup(AliRDHFCuts::kRejectPileupEvent);

  //Do not recalculate the vertex
  RDHFD0toKpi->SetRemoveDaughtersFromPrim(kFALSE); //activate for pp

  TString cent="";
  //centrality selection (Pb-Pb)
  Float_t minc=20,maxc=80;
  RDHFD0toKpi->SetMinCentrality(minc);
  RDHFD0toKpi->SetMaxCentrality(maxc);
  cent=Form("%.0f%.0f",minc,maxc);
  RDHFD0toKpi->SetUseCentrality(AliRDHFCuts::kCentV0M); //kCentOff,kCentV0M,kCentTRK,kCentTKL,kCentCL1,kCentInvalid

  //temporary
  //RDHFD0toKpi->SetFixRefs();

  cout<<"This is the odject I'm going to save:"<<endl;
  RDHFD0toKpi->PrintAll();
  TFile* fout=new TFile(Form("D0toKpiCuts%s%s%sRecVtx%sPileupRej.root", RDHFD0toKpi->GetUseCentrality()==0 ? "pp" : "PbPb",cent.Data(),RDHFD0toKpi->GetIsPrimaryWithoutDaughters() ? "" : "No",RDHFD0toKpi->GetOptPileUp() ? "" : "No"),"recreate");   //set this!! 

  fout->cd();
  RDHFD0toKpi->Write();
  fout->Close();

}
 
//macro to make a .root file (for significance maximization) which contains an AliRDHFCutsD0toKpi with loose set of cuts  and TParameter with the tighest value of these cuts

void makeInputAliAnalysisTaskSESignificanceMaximization(){

  AliRDHFCutsD0toKpi* RDHFD0toKpi=new AliRDHFCutsD0toKpi();
  RDHFD0toKpi->SetName("loosercuts");
  RDHFD0toKpi->SetTitle("Cuts for significance maximization");

  AliESDtrackCuts* esdTrackCuts=new AliESDtrackCuts();
  esdTrackCuts->SetRequireSigmaToVertex(kFALSE);
  //default
  esdTrackCuts->SetRequireTPCRefit(kTRUE);
  esdTrackCuts->SetRequireITSRefit(kTRUE);
  //esdTrackCuts->SetMinNClustersITS(4);
  //esdTrackCuts->SetMinNClustersTPC(120);

  esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,AliESDtrackCuts::kAny);
  esdTrackCuts->SetMinDCAToVertexXY(0.);
  esdTrackCuts->SetEtaRange(-0.8,0.8);
  esdTrackCuts->SetPtRange(0.8,1.e10);
  
  RDHFD0toKpi->AddTrackCuts(esdTrackCuts);

  const Int_t nvars=11;

  const Int_t nptbins=14; //change this when adding pt bins!
  Float_t ptbins[nptbins+1];
  ptbins[0]=0.;
  ptbins[1]=0.5;
  ptbins[2]=1.;
  ptbins[3]=2.;
  ptbins[4]=3.;
  ptbins[5]=4.;
  ptbins[6]=5.;
  ptbins[7]=6.;
  ptbins[8]=8.;
  ptbins[9]=10.;
  ptbins[10]=12.;
  ptbins[11]=16.;
  ptbins[12]=20.;
  ptbins[13]=24.;
  ptbins[14]=9999.;

  RDHFD0toKpi->SetPtBins(nptbins+1,ptbins);

  Float_t** rdcutsvalmine;
  rdcutsvalmine=new Float_t*[nvars];
  for(Int_t iv=0;iv<nvars;iv++){
    rdcutsvalmine[iv]=new Float_t[nptbins];
  }

  //setting my cut values
    //cuts order
    //       printf("    |M-MD0| [GeV]    < %f\n",fD0toKpiCuts[0]);
    //     printf("    dca    [cm]  < %f\n",fD0toKpiCuts[1]);
    //     printf("    cosThetaStar     < %f\n",fD0toKpiCuts[2]);
    //     printf("    pTK     [GeV/c]    > %f\n",fD0toKpiCuts[3]);
    //     printf("    pTpi    [GeV/c]    > %f\n",fD0toKpiCuts[4]);
    //     printf("    |d0K|  [cm]  < %f\n",fD0toKpiCuts[5]);
    //     printf("    |d0pi| [cm]  < %f\n",fD0toKpiCuts[6]);
    //     printf("    d0d0  [cm^2] < %f\n",fD0toKpiCuts[7]);
    //     printf("    cosThetaPoint    > %f\n",fD0toKpiCuts[8]);
    //     printf("    |cosThetaPointXY| < %f\n",fD0toKpiCuts[9]);
    //     printf("    NormDecayLenghtXY    > %f\n",fD0toKpiCuts[10]);


    Float_t cutsMatrixD0toKpiStand[nptbins][nvars]={{0.400,400.*1E-4,0.8,0.3,0.3,1000.*1E-4,1000.*1E-4,-40000.*1E-8,0.75,0.,2.},/* pt<0.5*/
						  {0.400,400.*1E-4,0.8,0.3,0.3,1000.*1E-4,1000.*1E-4,-40000.*1E-8,0.75,0.,2.},/* 0.5<pt<1*/
						  {0.400,400.*1E-4,0.8,0.4,0.4,1000.*1E-4,1000.*1E-4,-33000.*1E-8,0.75,0.,2.},/* 1<pt<2 */
						  {0.400,250.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-30000.*1E-8,0.85,0.994,2.},/* 2<pt<3 */
						  {0.400,250.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-26000.*1E-8,0.85,0.994,2.},/* 3<pt<4 */
						  {0.400,250.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-15000.*1E-8,0.85,0.994,2.},/* 4<pt<5 */
						  {0.400,250.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-11000.*1E-8,0.82,0.994,2.},/* 5<pt<6 */
						  {0.400,270.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-10000.*1E-8,0.78,0.994,2.},/* 6<pt<8 */
						  {0.400,300.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-1000.*1E-8,0.7,0.994,2.},/* 8<pt<10 */
						  {0.400,300.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-1000.*1E-8,0.7,0.994,2.},/* 10<pt<12 */
						  {0.400,350.*1E-4,1.0,0.7,0.7,1000.*1E-4,1000.*1E-4,-500.*1E-8,0.7,0.994,2.},/* 12<pt<16 */
						  {0.400,400.*1E-4,1.0,0.7,0.7,1000.*1E-4,1000.*1E-4,-500.*1E-8,0.7,0.994,2.},/* 16<pt<20 */
						  {0.400,400.*1E-4,1.0,0.7,0.7,1000.*1E-4,1000.*1E-4,-500.*1E-8,0.7,0.994,2.},/* 20<pt<24 */
						  {0.400,400.*1E-4,1.0,0.7,0.7,1000.*1E-4,1000.*1E-4,-500.*1E-8,0.7,0.994,2.}};/* pt>24 */

  //CREATE TRANSPOSE MATRIX...REVERSE INDICES as required by AliRDHFCuts
  Float_t **cutsMatrixTransposeStand=new Float_t*[nvars];
  for(Int_t iv=0;iv<nvars;iv++)cutsMatrixTransposeStand[iv]=new Float_t[nptbins];
  for (Int_t ibin=0;ibin<nptbins;ibin++){
    for (Int_t ivar = 0; ivar<nvars; ivar++){
      cutsMatrixTransposeStand[ivar][ibin]=cutsMatrixD0toKpiStand[ibin][ivar];
    }
  }
  RDHFD0toKpi->SetCuts(nvars,nptbins,cutsMatrixTransposeStand);


  Bool_t stdvaropt=kFALSE;
  Int_t dim=4; //set this!!
  Bool_t *boolforopt;
  boolforopt=new Bool_t[nvars];
  if(stdvaropt){
    boolforopt=RDHFD0toKpi->GetVarsForOpt();
  }else{
    TString *names;
    names=new TString[nvars];
    TString answer="";
    Int_t checktrue=0;
    names=RDHFD0toKpi->GetVarNames();
    for(Int_t i=0;i<nvars;i++){
      cout<<names[i]<<" for opt? (y/n)"<<endl;
      cin>>answer;
      if(answer=="y") {
	boolforopt[i]=kTRUE;
	checktrue++;
      }
      else boolforopt[i]=kFALSE;
    }
    if (checktrue!=dim) {
      cout<<"Error! You set "<<checktrue<<" kTRUE instead of "<<dim<<endl;
      return;
    }
    RDHFD0toKpi->SetVarsForOpt(dim,boolforopt);
  }


  Float_t tighterval[dim][nptbins];
  //dca  
  //costhetastar
  //d0d0 <-this 
  //costhetapoint <-this 
  //cosThetaPointXY <-this 
  //NormDecayLength <-this 
  
  //number of steps for each variable is set in the AddTask arguments (default=8)
  //set this!!
  //0-0.5
  tighterval[0][0]=-0.00065;
  tighterval[1][0]=1.;
  tighterval[2][0]=0.4;
  tighterval[3][0]=6.;

  //0.5-1.
  tighterval[0][1]=-0.00065;
  tighterval[1][1]=1.;
  tighterval[2][1]=0.4;
  tighterval[3][1]=6.;

  //1-2
  tighterval[0][2]=-0.00065;
  tighterval[1][2]=1.;
  tighterval[2][2]=0.4;
  tighterval[3][2]=6.;
 
  //2-3
  tighterval[0][3]=-0.0006;
  tighterval[1][3]=1.;
  tighterval[2][3]=1.;
  tighterval[3][3]=6.;

  //3-4
  tighterval[0][4]=-0.00046;
  tighterval[1][4]=1.;
  tighterval[2][4]=1.;
  tighterval[3][4]=6.;
 
  //4-5
  tighterval[0][5]=-0.00045;
  tighterval[1][5]=1.;
  tighterval[2][5]=1.;
  tighterval[3][5]=6.;
 
  //5-6
  tighterval[0][6]=-0.00031;
  tighterval[1][6]=1.;
  tighterval[2][6]=1.;
  tighterval[3][6]=6.;

  //6-8
  tighterval[0][7]=-0.00021;
  tighterval[1][7]=0.98;
  tighterval[2][7]=1.;
  tighterval[3][7]=6.;

  //8-10
  tighterval[0][8]=-0.0001;
  tighterval[1][8]=0.98;
  tighterval[2][8]=1.;
  tighterval[3][8]=6.;

  //10-12
  tighterval[0][9]=-0.0001;
  tighterval[1][9]=0.9;
  tighterval[2][9]=1.;
  tighterval[3][9]=6.;
 
  //12-16
  tighterval[0][10]=-0.00005;
  tighterval[1][10]=0.9;
  tighterval[2][10]=1.;
  tighterval[3][10]=6.;

  //16-20
  tighterval[0][11]=-0.00005;
  tighterval[1][11]=0.9;
  tighterval[2][11]=1.;
  tighterval[3][11]=6.;

  //20-24
  tighterval[0][12]=-0.00005;
  tighterval[1][12]=0.9;
  tighterval[2][12]=1.;
  tighterval[3][12]=6.;

  //>24
  tighterval[0][13]=-0.00005;
  tighterval[1][13]=0.9;
  tighterval[2][13]=1.;
  tighterval[3][13]=6.;


  TString name=""; 
  Int_t arrdim=dim*nptbins;
  cout<<"Will save "<<arrdim<<" TParameter<float>"<<endl;
  TClonesArray max("TParameter<float>",arrdim);
  for(Int_t ival=0;ival<dim;ival++){
    for(Int_t jpt=0;jpt<nptbins;jpt++){
      name=Form("par%dptbin%d",ival,jpt);
      cout<<"Setting "<<name.Data()<<" to "<<tighterval[ival][jpt]<<endl;
      new(max[jpt*dim+ival])TParameter<float>(name.Data(),tighterval[ival][jpt]);
    }
  }

  Bool_t flagPID=kTRUE;
  RDHFD0toKpi->SetUsePID(flagPID);

  RDHFD0toKpi->PrintAll();
  printf("Use PID? %s\n",flagPID ? "yes" : "no");

  //pid settings
  AliAODPidHF* pidObj=new AliAODPidHF();
  //pidObj->SetName("pid4D0");
  Int_t mode=1;
  const Int_t nlims=2;
  Double_t plims[nlims]={0.6,0.8}; //TPC limits in momentum [GeV/c]
  Bool_t compat=kTRUE; //effective only for this mode
  Bool_t asym=kTRUE;
  Double_t sigmas[5]={2.,1.,0.,3.,0.}; //to be checked and to be modified with new implementation of setters by Rossella
  pidObj->SetAsym(asym);// if you want to use the asymmetric bands in TPC
  pidObj->SetMatch(mode);
  pidObj->SetPLimit(plims,nlims);
  pidObj->SetSigma(sigmas);
  pidObj->SetCompat(compat);
  pidObj->SetTPC(kTRUE);
  pidObj->SetTOF(kTRUE);
  RDHFD0toKpi->SetPidHF(pidObj);

  RDHFD0toKpi->SetUseDefaultPID(kFALSE); //to use the AliAODPidHF
  
  //activate pileup rejection (for pp)
  //RDHFD0toKpi->SetOptPileup(AliRDHFCuts::kRejectPileupEvent);

  //Do not recalculate the vertex
  RDHFD0toKpi->SetRemoveDaughtersFromPrim(kFALSE); //activate for pp

  TString cent="";
  //centrality selection (Pb-Pb)
  Float_t minc=20,maxc=80;
  RDHFD0toKpi->SetMinCentrality(minc);
  RDHFD0toKpi->SetMaxCentrality(maxc);
  cent=Form("%.0f%.0f",minc,maxc);
  RDHFD0toKpi->SetUseCentrality(AliRDHFCuts::kCentV0M); //kCentOff,kCentV0M,kCentTRK,kCentTKL,kCentCL1,kCentInvalid

  //temporary
  RDHFD0toKpi->SetFixRefs();

  TFile* fout=new TFile(Form("cuts4SignifMaxim%s%s%sRecVtx%sPileupRej.root", RDHFD0toKpi->GetUseCentrality()==0 ? "pp" : "PbPb",cent.Data(),RDHFD0toKpi->GetIsPrimaryWithoutDaughters() ? "" : "No",RDHFD0toKpi->GetOptPileUp() ? "" : "No"),"recreate");   //set this!! 

  fout->cd();
  RDHFD0toKpi->Write();
  max.Write();
  fout->Close();
 
}
Commit	Line	Data
d52f7b50	1	#include <Riostream.h>
	2	#include <TFile.h>
	3	#include <AliRDHFCutsD0toKpi.h>
	4	#include <TClonesArray.h>
	5	#include <TParameter.h>
cbedddce	6
cbedddce	7
	8	//Use:
	9	//Set hard coded commentet with //set this!!
3b79985d	10	// root[] .L makeInputD0tasks.C++
	11	// root[] makeInputAliAnalysisTaskSED0Mass()
	12	// root[] makeInputAliAnalysisTaskSESignificanceMaximization()
cbedddce	13	//similar macros for the other D mesons
	14
	15	//Author: Chiara Bianchin, cbianchi@pd.infn.it
	16
3b79985d	17
3b79985d	18	//macro to make a .root file which contains an AliRDHFCutsD0toKpi for AliAnalysisTaskSED0Mass task
6fff1d5d	19	void ModifyFromStandardCuts(Int_t system=1 /0=pp, 1=PbPb, 2=pp 2.76TeV/){
	20	AliRDHFCutsD0toKpi* RDHFD0toKpi=new AliRDHFCutsD0toKpi();
	21	TString info="",cent="";
	22
	23	if(system==0) {
	24	RDHFD0toKpi->SetStandardCutsPP2010();
	25	info+="ppStd";
	26	}
	27	if(system==1) {
	28	RDHFD0toKpi->SetStandardCutsPbPb2010();
	29	//Change centrality if needed
	30	/*
	31	Float_t minc=20,maxc=80;
	32	RDHFD0toKpi->SetMinCentrality(minc);
	33	RDHFD0toKpi->SetMaxCentrality(maxc);
	34	*/
	35	RDHFD0toKpi->SetUseCentrality(AliRDHFCuts::kCentOff); //kCentOff,kCentV0M,kCentTRK,kCentTKL,kCentCL1,kCentInvalid
	36
	37	cent=Form("%.0f%.0f",RDHFD0toKpi->GetMinCentrality(),RDHFD0toKpi->GetMaxCentrality());
	38	info+="PbPbStd";
	39
	40	}
	41	if(system==2) {
	42	RDHFD0toKpi->SetStandardCutsPP2011_276TeV();
	43	info+="pp276Std";
	44	}
	45
	46	RDHFD0toKpi->SetName("D0toKpiCuts");
	47	RDHFD0toKpi->SetTitle("Cuts for D0 analysis");
	48
	49
	50	//here add what you need! Find examples below
	51
	52	//Trigger mask
	53	/*
	54	RDHFD0toKpi->SetTriggerMask(0);
	55	RDHFD0toKpi->SetTriggerMask(AliVEvent::kEMC1 \| AliVEvent::kEMC7);
	56	RDHFD0toKpi->SetTriggerClass("CEMC");
	57
	58	info+="EMCTr";
	59	*/
	60
	61	//Event selection
	62	RDHFD0toKpi->SetUsePhysicsSelection(kFALSE);
	63	info+="noPhysSel";
	64
	65	cout<<"This is the odject I'm going to save:"<<endl;
	66	cout<<info<<endl;
	67	RDHFD0toKpi->PrintAll();
	68	TFile* fout=new TFile(Form("D0toKpiCuts%s%sRecVtx%sPileupRej%s.root", RDHFD0toKpi->GetUseCentrality()==0 ? "" : cent.Data(),RDHFD0toKpi->GetIsPrimaryWithoutDaughters() ? "" : "No",RDHFD0toKpi->GetOptPileUp() ? "" : "No",info.Data()),"recreate"); //set this!!
	69
	70	fout->cd();
	71	RDHFD0toKpi->Write();
	72	fout->Close();
	73
	74	}
3b79985d	75
	76	void makeInputAliAnalysisTaskSED0Mass(){
	77
	78	AliRDHFCutsD0toKpi* RDHFD0toKpi=new AliRDHFCutsD0toKpi();
	79	RDHFD0toKpi->SetName("D0toKpiCuts");
	80	RDHFD0toKpi->SetTitle("Cuts for D0 analysis");
	81
	82	AliESDtrackCuts* esdTrackCuts=new AliESDtrackCuts();
	83	esdTrackCuts->SetRequireSigmaToVertex(kFALSE);
	84	//default
	85	esdTrackCuts->SetRequireTPCRefit(kTRUE);
	86	esdTrackCuts->SetRequireITSRefit(kTRUE);
44fb0179	87	esdTrackCuts->SetEtaRange(-0.8,0.8);
0d4398a2	88	//esdTrackCuts->SetMinNClustersITS(4); // default is 5
44fb0179	89	//esdTrackCuts->SetMinNClustersTPC(120);
3b79985d	90	esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,
	91	AliESDtrackCuts::kAny);
	92	// default is kBoth, otherwise kAny
	93	esdTrackCuts->SetMinDCAToVertexXY(0.);
0beb4580	94	esdTrackCuts->SetPtRange(0.8,1.e10);
3b79985d	95
	96
	97	RDHFD0toKpi->AddTrackCuts(esdTrackCuts);
	98
29068278	99	const Int_t nvars=11;
3b79985d	100
0beb4580	101	const Int_t nptbins=13;
3b79985d	102	Float_t* ptbins;
	103	ptbins=new Float_t[nptbins+1];
	104	ptbins[0]=0.;
0beb4580	105	ptbins[1]=0.5;
	106	ptbins[2]=1.;
	107	ptbins[3]=2.;
	108	ptbins[4]=3.;
	109	ptbins[5]=4.;
	110	ptbins[6]=5.;
	111	ptbins[7]=6.;
	112	ptbins[8]=8.;
	113	ptbins[9]=12.;
	114	ptbins[10]=16.;
	115	ptbins[11]=20.;
	116	ptbins[12]=24.;
	117	ptbins[13]=9999.;
	118
3b79985d	119	RDHFD0toKpi->SetPtBins(nptbins+1,ptbins);
	120
	121
	122	Float_t** rdcutsvalmine;
	123	rdcutsvalmine=new Float_t*[nvars];
	124	for(Int_t iv=0;iv<nvars;iv++){
	125	rdcutsvalmine[iv]=new Float_t[nptbins];
	126	}
	127
	128	//setting my cut values
	129	//cuts order
	130	// printf(" \|M-MD0\| [GeV] < %f\n",fD0toKpiCuts[0]);
	131	// printf(" dca [cm] < %f\n",fD0toKpiCuts[1]);
	132	// printf(" cosThetaStar < %f\n",fD0toKpiCuts[2]);
	133	// printf(" pTK [GeV/c] > %f\n",fD0toKpiCuts[3]);
	134	// printf(" pTpi [GeV/c] > %f\n",fD0toKpiCuts[4]);
	135	// printf(" \|d0K\| [cm] < %f\n",fD0toKpiCuts[5]);
	136	// printf(" \|d0pi\| [cm] < %f\n",fD0toKpiCuts[6]);
	137	// printf(" d0d0 [cm^2] < %f\n",fD0toKpiCuts[7]);
	138	// printf(" cosThetaPoint > %f\n",fD0toKpiCuts[8]);
29068278	139	// printf(" \|cosThetaPointXY\| < %f\n",fD0toKpiCuts[9]);
29068278	140	// printf(" NormDecayLenghtXY > %f\n",fD0toKpiCuts[10]);
3b79985d	141
0beb4580	142
29068278	143	Double_t arrcuts[11]={0.3,0.03,0.8,0.8,0.8,0.1,0.1,-0.0004,0.9,0.998,5.}; //put the last 2 values at 0. for pp
0beb4580	144
3b79985d	145	//setting my cut values
0beb4580	146	//0-0.5 GeV
	147	for(Int_t ic=0;ic<nvars;ic++) rdcutsvalmine[ic][0]=arrcuts[ic];
	148
	149	//0.5-1 GeV/c
	150	for(Int_t ic=0;ic<nvars;ic++) rdcutsvalmine[ic][1]=arrcuts[ic];
6b3e3c78	151
0beb4580	152	//1-2 GeV
	153	arrcuts[1]=0.025; arrcuts[7]=-0.0003;
	154	for(Int_t ic=0;ic<nvars;ic++) rdcutsvalmine[ic][2]=arrcuts[ic];
6b3e3c78	155
6b3e3c78	156	//2-3 GeV
0beb4580	157	arrcuts[7]=-0.00026;
29068278	158	arrcuts[9]=0.998;
0beb4580	159	for(Int_t ic=0;ic<nvars;ic++) rdcutsvalmine[ic][3]=arrcuts[ic];
	160
	161	//3-4 GeV
	162	arrcuts[7]=-0.00015; arrcuts[8]=0.85;
	163	for(Int_t ic=0;ic<nvars;ic++) rdcutsvalmine[ic][4]=arrcuts[ic];
	164
	165	//4-5 GeV
	166	for(Int_t ic=0;ic<nvars;ic++) rdcutsvalmine[ic][5]=arrcuts[ic];
	167
	168	//5-6 GeV
	169	arrcuts[7]=-0.0001; arrcuts[8]=0.85;
	170	for(Int_t ic=0;ic<nvars;ic++) rdcutsvalmine[ic][6]=arrcuts[ic];
	171
	172	//6-8 GeV
	173	arrcuts[2]=1.;
	174	for(Int_t ic=0;ic<nvars;ic++) rdcutsvalmine[ic][7]=arrcuts[ic];
	175
	176	//8-12 GeV
	177	arrcuts[8]=0.8;
	178	for(Int_t ic=0;ic<nvars;ic++) rdcutsvalmine[ic][8]=arrcuts[ic];
	179
	180	//12-16 GeV
	181	arrcuts[1]=0.03;
	182	for(Int_t ic=0;ic<nvars;ic++) rdcutsvalmine[ic][9]=arrcuts[ic];
	183
	184	//16-20 GeV
	185	arrcuts[1]=0.035;
	186	for(Int_t ic=0;ic<nvars;ic++) rdcutsvalmine[ic][10]=arrcuts[ic];
	187
	188	//20-24 GeV
	189	for(Int_t ic=0;ic<nvars;ic++) rdcutsvalmine[ic][11]=arrcuts[ic];
	190
	191	//24-9999 GeV
	192	for(Int_t ic=0;ic<nvars;ic++) rdcutsvalmine[ic][12]=arrcuts[ic];
3b79985d	193
3b79985d	194	RDHFD0toKpi->SetCuts(nvars,nptbins,rdcutsvalmine);
6b3e3c78	195
	196	Bool_t pidflag=kTRUE;
	197	RDHFD0toKpi->SetUsePID(pidflag);
	198	if(pidflag) cout<<"PID is used"<<endl;
	199	else cout<<"PID is not used"<<endl;
	200
d52f7b50	201	//pid settings
6b3e3c78	202	AliAODPidHF* pidObj=new AliAODPidHF();
	203	//pidObj->SetName("pid4D0");
	204	Int_t mode=1;
	205	const Int_t nlims=2;
	206	Double_t plims[nlims]={0.6,0.8}; //TPC limits in momentum [GeV/c]
	207	Bool_t compat=kTRUE; //effective only for this mode
	208	Bool_t asym=kTRUE;
	209	Double_t sigmas[5]={2.,1.,0.,3.,0.}; //to be checked and to be modified with new implementation of setters by Rossella
	210	pidObj->SetAsym(asym);// if you want to use the asymmetric bands in TPC
	211	pidObj->SetMatch(mode);
	212	pidObj->SetPLimit(plims,nlims);
	213	pidObj->SetSigma(sigmas);
	214	pidObj->SetCompat(compat);
	215	pidObj->SetTPC(kTRUE);
	216	pidObj->SetTOF(kTRUE);
	217	RDHFD0toKpi->SetPidHF(pidObj);
	218
	219	RDHFD0toKpi->SetUseDefaultPID(kFALSE); //to use the AliAODPidHF
	220
0beb4580	221	//activate pileup rejection (for pp)
	222	//RDHFD0toKpi->SetOptPileup(AliRDHFCuts::kRejectPileupEvent);
	223
	224	//Do not recalculate the vertex
	225	RDHFD0toKpi->SetRemoveDaughtersFromPrim(kFALSE); //activate for pp
	226
	227	TString cent="";
	228	//centrality selection (Pb-Pb)
	229	Float_t minc=20,maxc=80;
	230	RDHFD0toKpi->SetMinCentrality(minc);
	231	RDHFD0toKpi->SetMaxCentrality(maxc);
	232	cent=Form("%.0f%.0f",minc,maxc);
	233	RDHFD0toKpi->SetUseCentrality(AliRDHFCuts::kCentV0M); //kCentOff,kCentV0M,kCentTRK,kCentTKL,kCentCL1,kCentInvalid
	234
	235	//temporary
29068278	236	//RDHFD0toKpi->SetFixRefs();
d52f7b50	237
3b79985d	238	cout<<"This is the odject I'm going to save:"<<endl;
3b79985d	239	RDHFD0toKpi->PrintAll();
0beb4580	240	TFile* fout=new TFile(Form("D0toKpiCuts%s%s%sRecVtx%sPileupRej.root", RDHFD0toKpi->GetUseCentrality()==0 ? "pp" : "PbPb",cent.Data(),RDHFD0toKpi->GetIsPrimaryWithoutDaughters() ? "" : "No",RDHFD0toKpi->GetOptPileUp() ? "" : "No"),"recreate"); //set this!!
0beb4580	241
3b79985d	242	fout->cd();
	243	RDHFD0toKpi->Write();
	244	fout->Close();
	245
	246	}
29068278	247
3b79985d	248	//macro to make a .root file (for significance maximization) which contains an AliRDHFCutsD0toKpi with loose set of cuts and TParameter with the tighest value of these cuts
	249
	250	void makeInputAliAnalysisTaskSESignificanceMaximization(){
cbedddce	251
	252	AliRDHFCutsD0toKpi* RDHFD0toKpi=new AliRDHFCutsD0toKpi();
	253	RDHFD0toKpi->SetName("loosercuts");
	254	RDHFD0toKpi->SetTitle("Cuts for significance maximization");
	255
	256	AliESDtrackCuts* esdTrackCuts=new AliESDtrackCuts();
	257	esdTrackCuts->SetRequireSigmaToVertex(kFALSE);
	258	//default
	259	esdTrackCuts->SetRequireTPCRefit(kTRUE);
	260	esdTrackCuts->SetRequireITSRefit(kTRUE);
0d4398a2	261	//esdTrackCuts->SetMinNClustersITS(4);
44fb0179	262	//esdTrackCuts->SetMinNClustersTPC(120);
44fb0179	263
cbedddce	264	esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,AliESDtrackCuts::kAny);
cbedddce	265	esdTrackCuts->SetMinDCAToVertexXY(0.);
44fb0179	266	esdTrackCuts->SetEtaRange(-0.8,0.8);
0beb4580	267	esdTrackCuts->SetPtRange(0.8,1.e10);
cbedddce	268
	269	RDHFD0toKpi->AddTrackCuts(esdTrackCuts);
	270
29068278	271	const Int_t nvars=11;
cbedddce	272
29068278	273	const Int_t nptbins=14; //change this when adding pt bins!
cbedddce	274	Float_t ptbins[nptbins+1];
cbedddce	275	ptbins[0]=0.;
6b3e3c78	276	ptbins[1]=0.5;
	277	ptbins[2]=1.;
	278	ptbins[3]=2.;
	279	ptbins[4]=3.;
	280	ptbins[5]=4.;
	281	ptbins[6]=5.;
	282	ptbins[7]=6.;
	283	ptbins[8]=8.;
29068278	284	ptbins[9]=10.;
	285	ptbins[10]=12.;
	286	ptbins[11]=16.;
	287	ptbins[12]=20.;
	288	ptbins[13]=24.;
	289	ptbins[14]=9999.;
0beb4580	290
cbedddce	291	RDHFD0toKpi->SetPtBins(nptbins+1,ptbins);
cbedddce	292
	293	Float_t** rdcutsvalmine;
	294	rdcutsvalmine=new Float_t*[nvars];
	295	for(Int_t iv=0;iv<nvars;iv++){
	296	rdcutsvalmine[iv]=new Float_t[nptbins];
	297	}
	298
	299	//setting my cut values
	300	//cuts order
	301	// printf(" \|M-MD0\| [GeV] < %f\n",fD0toKpiCuts[0]);
	302	// printf(" dca [cm] < %f\n",fD0toKpiCuts[1]);
	303	// printf(" cosThetaStar < %f\n",fD0toKpiCuts[2]);
	304	// printf(" pTK [GeV/c] > %f\n",fD0toKpiCuts[3]);
	305	// printf(" pTpi [GeV/c] > %f\n",fD0toKpiCuts[4]);
	306	// printf(" \|d0K\| [cm] < %f\n",fD0toKpiCuts[5]);
	307	// printf(" \|d0pi\| [cm] < %f\n",fD0toKpiCuts[6]);
	308	// printf(" d0d0 [cm^2] < %f\n",fD0toKpiCuts[7]);
	309	// printf(" cosThetaPoint > %f\n",fD0toKpiCuts[8]);
29068278	310	// printf(" \|cosThetaPointXY\| < %f\n",fD0toKpiCuts[9]);
	311	// printf(" NormDecayLenghtXY > %f\n",fD0toKpiCuts[10]);
	312
	313
	314	Float_t cutsMatrixD0toKpiStand[nptbins][nvars]={{0.400,400.1E-4,0.8,0.3,0.3,1000.1E-4,1000.1E-4,-40000.1E-8,0.75,0.,2.},/* pt<0.5*/
	315	{0.400,400.1E-4,0.8,0.3,0.3,1000.1E-4,1000.1E-4,-40000.1E-8,0.75,0.,2.},/* 0.5<pt<1*/
	316	{0.400,400.1E-4,0.8,0.4,0.4,1000.1E-4,1000.1E-4,-33000.1E-8,0.75,0.,2.},/* 1<pt<2 */
	317	{0.400,250.1E-4,0.8,0.7,0.7,1000.1E-4,1000.1E-4,-30000.1E-8,0.85,0.994,2.},/* 2<pt<3 */
	318	{0.400,250.1E-4,0.8,0.7,0.7,1000.1E-4,1000.1E-4,-26000.1E-8,0.85,0.994,2.},/* 3<pt<4 */
	319	{0.400,250.1E-4,0.8,0.7,0.7,1000.1E-4,1000.1E-4,-15000.1E-8,0.85,0.994,2.},/* 4<pt<5 */
	320	{0.400,250.1E-4,0.8,0.7,0.7,1000.1E-4,1000.1E-4,-11000.1E-8,0.82,0.994,2.},/* 5<pt<6 */
	321	{0.400,270.1E-4,0.8,0.7,0.7,1000.1E-4,1000.1E-4,-10000.1E-8,0.78,0.994,2.},/* 6<pt<8 */
	322	{0.400,300.1E-4,0.8,0.7,0.7,1000.1E-4,1000.1E-4,-1000.1E-8,0.7,0.994,2.},/* 8<pt<10 */
	323	{0.400,300.1E-4,0.8,0.7,0.7,1000.1E-4,1000.1E-4,-1000.1E-8,0.7,0.994,2.},/* 10<pt<12 */
	324	{0.400,350.1E-4,1.0,0.7,0.7,1000.1E-4,1000.1E-4,-500.1E-8,0.7,0.994,2.},/* 12<pt<16 */
	325	{0.400,400.1E-4,1.0,0.7,0.7,1000.1E-4,1000.1E-4,-500.1E-8,0.7,0.994,2.},/* 16<pt<20 */
	326	{0.400,400.1E-4,1.0,0.7,0.7,1000.1E-4,1000.1E-4,-500.1E-8,0.7,0.994,2.},/* 20<pt<24 */
	327	{0.400,400.1E-4,1.0,0.7,0.7,1000.1E-4,1000.1E-4,-500.1E-8,0.7,0.994,2.}};/* pt>24 */
d52f7b50	328
6b3e3c78	329	//CREATE TRANSPOSE MATRIX...REVERSE INDICES as required by AliRDHFCuts
d52f7b50	330	Float_t *cutsMatrixTransposeStand=new Float_t[nvars];
d52f7b50	331	for(Int_t iv=0;iv<nvars;iv++)cutsMatrixTransposeStand[iv]=new Float_t[nptbins];
6b3e3c78	332	for (Int_t ibin=0;ibin<nptbins;ibin++){
6b3e3c78	333	for (Int_t ivar = 0; ivar<nvars; ivar++){
d52f7b50	334	cutsMatrixTransposeStand[ivar][ibin]=cutsMatrixD0toKpiStand[ibin][ivar];
6b3e3c78	335	}
6b3e3c78	336	}
d52f7b50	337	RDHFD0toKpi->SetCuts(nvars,nptbins,cutsMatrixTransposeStand);
cbedddce	338
cbedddce	339
29068278	340	Bool_t stdvaropt=kFALSE;
29068278	341	Int_t dim=4; //set this!!
cbedddce	342	Bool_t *boolforopt;
cbedddce	343	boolforopt=new Bool_t[nvars];
29068278	344	if(stdvaropt){
	345	boolforopt=RDHFD0toKpi->GetVarsForOpt();
	346	}else{
	347	TString *names;
	348	names=new TString[nvars];
	349	TString answer="";
	350	Int_t checktrue=0;
	351	names=RDHFD0toKpi->GetVarNames();
	352	for(Int_t i=0;i<nvars;i++){
	353	cout<<names[i]<<" for opt? (y/n)"<<endl;
	354	cin>>answer;
	355	if(answer=="y") {
	356	boolforopt[i]=kTRUE;
	357	checktrue++;
cbedddce	358	}
29068278	359	else boolforopt[i]=kFALSE;
cbedddce	360	}
29068278	361	if (checktrue!=dim) {
	362	cout<<"Error! You set "<<checktrue<<" kTRUE instead of "<<dim<<endl;
	363	return;
	364	}
	365	RDHFD0toKpi->SetVarsForOpt(dim,boolforopt);
cbedddce	366	}
	367
	368
	369	Float_t tighterval[dim][nptbins];
44fb0179	370	//dca
cbedddce	371	//costhetastar
	372	//d0d0 <-this
	373	//costhetapoint <-this
29068278	374	//cosThetaPointXY <-this
29068278	375	//NormDecayLength <-this
cbedddce	376
44fb0179	377	//number of steps for each variable is set in the AddTask arguments (default=8)
cbedddce	378	//set this!!
29068278	379	//0-0.5
44fb0179	380	tighterval[0][0]=-0.00065;
44fb0179	381	tighterval[1][0]=1.;
29068278	382	tighterval[2][0]=0.4;
29068278	383	tighterval[3][0]=6.;
cbedddce	384
29068278	385	//0.5-1.
44fb0179	386	tighterval[0][1]=-0.00065;
44fb0179	387	tighterval[1][1]=1.;
29068278	388	tighterval[2][1]=0.4;
29068278	389	tighterval[3][1]=6.;
44fb0179	390
29068278	391	//1-2
44fb0179	392	tighterval[0][2]=-0.00065;
44fb0179	393	tighterval[1][2]=1.;
29068278	394	tighterval[2][2]=0.4;
29068278	395	tighterval[3][2]=6.;
cbedddce	396
29068278	397	//2-3
29068278	398	tighterval[0][3]=-0.0006;
44fb0179	399	tighterval[1][3]=1.;
29068278	400	tighterval[2][3]=1.;
29068278	401	tighterval[3][3]=6.;
cbedddce	402
29068278	403	//3-4
29068278	404	tighterval[0][4]=-0.00046;
44fb0179	405	tighterval[1][4]=1.;
29068278	406	tighterval[2][4]=1.;
29068278	407	tighterval[3][4]=6.;
44fb0179	408
29068278	409	//4-5
29068278	410	tighterval[0][5]=-0.00045;
44fb0179	411	tighterval[1][5]=1.;
29068278	412	tighterval[2][5]=1.;
29068278	413	tighterval[3][5]=6.;
44fb0179	414
29068278	415	//5-6
29068278	416	tighterval[0][6]=-0.00031;
44fb0179	417	tighterval[1][6]=1.;
29068278	418	tighterval[2][6]=1.;
	419	tighterval[3][6]=6.;
	420
	421	//6-8
	422	tighterval[0][7]=-0.00021;
	423	tighterval[1][7]=0.98;
	424	tighterval[2][7]=1.;
	425	tighterval[3][7]=6.;
	426
	427	//8-10
	428	tighterval[0][8]=-0.0001;
	429	tighterval[1][8]=0.98;
	430	tighterval[2][8]=1.;
	431	tighterval[3][8]=6.;
	432
	433	//10-12
	434	tighterval[0][9]=-0.0001;
	435	tighterval[1][9]=0.9;
	436	tighterval[2][9]=1.;
	437	tighterval[3][9]=6.;
44fb0179	438
29068278	439	//12-16
	440	tighterval[0][10]=-0.00005;
	441	tighterval[1][10]=0.9;
	442	tighterval[2][10]=1.;
	443	tighterval[3][10]=6.;
	444
	445	//16-20
	446	tighterval[0][11]=-0.00005;
	447	tighterval[1][11]=0.9;
	448	tighterval[2][11]=1.;
	449	tighterval[3][11]=6.;
	450
	451	//20-24
	452	tighterval[0][12]=-0.00005;
	453	tighterval[1][12]=0.9;
	454	tighterval[2][12]=1.;
	455	tighterval[3][12]=6.;
	456
	457	//>24
	458	tighterval[0][13]=-0.00005;
	459	tighterval[1][13]=0.9;
	460	tighterval[2][13]=1.;
	461	tighterval[3][13]=6.;
44fb0179	462
44fb0179	463
cbedddce	464	TString name="";
	465	Int_t arrdim=dim*nptbins;
	466	cout<<"Will save "<<arrdim<<" TParameter<float>"<<endl;
	467	TClonesArray max("TParameter<float>",arrdim);
3b79985d	468	for(Int_t ival=0;ival<dim;ival++){
	469	for(Int_t jpt=0;jpt<nptbins;jpt++){
	470	name=Form("par%dptbin%d",ival,jpt);
	471	cout<<"Setting "<<name.Data()<<" to "<<tighterval[ival][jpt]<<endl;
	472	new(max[jpt*dim+ival])TParameter<float>(name.Data(),tighterval[ival][jpt]);
cbedddce	473	}
cbedddce	474	}
0beb4580	475
6b3e3c78	476	Bool_t flagPID=kTRUE;
6b3e3c78	477	RDHFD0toKpi->SetUsePID(flagPID);
0beb4580	478
d52f7b50	479	RDHFD0toKpi->PrintAll();
6b3e3c78	480	printf("Use PID? %s\n",flagPID ? "yes" : "no");
	481
	482	//pid settings
	483	AliAODPidHF* pidObj=new AliAODPidHF();
	484	//pidObj->SetName("pid4D0");
	485	Int_t mode=1;
	486	const Int_t nlims=2;
	487	Double_t plims[nlims]={0.6,0.8}; //TPC limits in momentum [GeV/c]
	488	Bool_t compat=kTRUE; //effective only for this mode
	489	Bool_t asym=kTRUE;
	490	Double_t sigmas[5]={2.,1.,0.,3.,0.}; //to be checked and to be modified with new implementation of setters by Rossella
	491	pidObj->SetAsym(asym);// if you want to use the asymmetric bands in TPC
	492	pidObj->SetMatch(mode);
	493	pidObj->SetPLimit(plims,nlims);
	494	pidObj->SetSigma(sigmas);
	495	pidObj->SetCompat(compat);
	496	pidObj->SetTPC(kTRUE);
	497	pidObj->SetTOF(kTRUE);
	498	RDHFD0toKpi->SetPidHF(pidObj);
	499
6b3e3c78	500	RDHFD0toKpi->SetUseDefaultPID(kFALSE); //to use the AliAODPidHF
d52f7b50	501
0beb4580	502	//activate pileup rejection (for pp)
	503	//RDHFD0toKpi->SetOptPileup(AliRDHFCuts::kRejectPileupEvent);
	504
	505	//Do not recalculate the vertex
	506	RDHFD0toKpi->SetRemoveDaughtersFromPrim(kFALSE); //activate for pp
	507
	508	TString cent="";
	509	//centrality selection (Pb-Pb)
	510	Float_t minc=20,maxc=80;
	511	RDHFD0toKpi->SetMinCentrality(minc);
	512	RDHFD0toKpi->SetMaxCentrality(maxc);
	513	cent=Form("%.0f%.0f",minc,maxc);
	514	RDHFD0toKpi->SetUseCentrality(AliRDHFCuts::kCentV0M); //kCentOff,kCentV0M,kCentTRK,kCentTKL,kCentCL1,kCentInvalid
d52f7b50	515
0beb4580	516	//temporary
0beb4580	517	RDHFD0toKpi->SetFixRefs();
6b3e3c78	518
0beb4580	519	TFile* fout=new TFile(Form("cuts4SignifMaxim%s%s%sRecVtx%sPileupRej.root", RDHFD0toKpi->GetUseCentrality()==0 ? "pp" : "PbPb",cent.Data(),RDHFD0toKpi->GetIsPrimaryWithoutDaughters() ? "" : "No",RDHFD0toKpi->GetOptPileUp() ? "" : "No"),"recreate"); //set this!!
44fb0179	520
cbedddce	521	fout->cd();
	522	RDHFD0toKpi->Write();
	523	max.Write();
	524	fout->Close();
	525
	526	}
3b79985d	527