[u/mrichter/AliRoot.git] / PWGHF / hfe / macros / configs / PbPb / ConfigHFEpbpb.C

TF1* GetEtaCorrection(TString listname){
  
  TString etaMap="$ALICE_ROOT/PWGDQ/dielectron/files/EtaCorrMaps.root";
  
  if (gSystem->AccessPathName(gSystem->ExpandPathName(etaMap.Data()))){
    Error("ConfigHFEpbpb","Eta map not found: %s",etaMap.Data());
    return 0;
  }

  TFile f(etaMap.Data());
  if (!f.IsOpen()) return 0;
  gROOT->cd();
  TList *keys=f.GetListOfKeys();

  for (Int_t i=0; i<keys->GetEntries(); ++i){
    TString kName=keys->At(i)->GetName();
    TPRegexp reg(kName);
    if (reg.MatchB(listname)){
      printf("Using Eta Correction Function: %s\n",kName.Data());
      return (TF1*)f.Get(kName.Data());
    }
  }
  return 0;
}


AliAnalysisTaskHFE* ConfigHFEpbpb(Bool_t useMC=kFALSE,
				  TString appendix,
				  Int_t aodfilter=-1,
				  UChar_t TPCcl=70, UChar_t TPCclPID = 80, 
				  Double_t TPCclRatio = 0.6, Double_t TPCclshared = 1.1,
				  UChar_t ITScl=3,  Double_t ITSchi2perclusters=99999999.,
				  Int_t itspixelcut=AliHFEextraCuts::kFirst,
                                  Double_t dcaxy=1.0, Double_t dcaz=2.0,
				  Double_t TOFs=3.,
				  Double_t IpSig=3., 				  
				  Float_t prodlow=0., Float_t prodhigh=100., 
				  Bool_t beauty=kFALSE,
				  Bool_t kMCQA = kFALSE, 
				  Bool_t kDEStep = kFALSE,
				  Int_t addflag=0, Int_t ptbin=0, 
				  Int_t etacor=0, TString listname="",
				  Int_t nondefaultcentr=0, Float_t* arraycentr=NULL,
				  Double_t* tpcdEdxcut=NULL,Double_t tpcu=3.0){
  //
  // HFE standard task configuration
  //

  AliHFEcuts *hfecuts = new AliHFEcuts(appendix,"HFE cuts pbpb TOF TPC");
  hfecuts->CreateStandardCuts();

  hfecuts->SetMinNClustersTPC(TPCcl);
  hfecuts->SetMinNClustersTPCPID(TPCclPID);
  hfecuts->SetMinRatioTPCclusters(TPCclRatio);
  hfecuts->SetTPCmodes(AliHFEextraCuts::kFound, AliHFEextraCuts::kFoundOverFindable);
  hfecuts->SetFractionOfSharedTPCClusters(TPCclshared);

  hfecuts->SetMinNClustersITS(ITScl);
  hfecuts->SetMaxChi2perClusterITS(ITSchi2perclusters);
  hfecuts->SetCutITSpixel(itspixelcut);
  hfecuts->SetCheckITSLayerStatus(kFALSE);

  hfecuts->SetIPcutParam(0,0,0,IpSig,kTRUE,kTRUE);
  // if(useMC && beauty) hfecuts->SetProductionVertex(prodlow,prodhigh,prodlow,prodhigh);
  if(useMC) hfecuts->SetProductionVertex(prodlow,prodhigh,prodlow,prodhigh);

  hfecuts->SetMaxImpactParam(dcaxy,dcaz);

  // event cuts
  hfecuts->SetUseMixedVertex(kTRUE);
  hfecuts->SetVertexRange(10.);

  // others
  hfecuts->SetTOFPIDStep(kTRUE);
  //hfecuts->SetMaxChi2perClusterITS(36);
  //hfecuts->SetTOFMISMATCHStep(kTRUE);
  //hfecuts->SetTPCPIDCleanUpStep(kTRUE);
  hfecuts->SetQAOn();

  AliAnalysisTaskHFE *task = new AliAnalysisTaskHFE(appendix);
  task->SetHFECuts(hfecuts);
  task->SetPbPbAnalysis(kTRUE);
  task->SetRemovePileUp(kFALSE);
  task->GetPIDQAManager()->SetHighResolutionHistos();
  if(useMC) task->SetHasMCData(kTRUE); // necessary for AOD
  printf("AOD filter %d On/OFF?\n",aodfilter);
  if(aodfilter > 0) {
    printf("ON AOD filter %d\n",aodfilter);
    task->SetUseFilterAOD(kTRUE);
    task->SetFilter(aodfilter);
  }
 
  
  if((nondefaultcentr!=0) && arraycentr) {
    for(Int_t i=0;i<12;i++)
      {
	task->SetPbPbUserCentralityLimit(kTRUE);
	task->SetPbPbUserCentralityArray(i,arraycentr[i]);
	
      }
  }
  
  // Define Variables

  if(ptbin==0)
  {
      Double_t ptbinning[36] = {0., 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1., 1.1, 1.2,
      1.3, 1.4, 1.5, 1.75, 2., 2.25, 2.5, 2.75, 3., 3.5, 4., 4.5, 5.,
      5.5, 6., 7., 8., 10., 12., 14., 16., 18., 20.};
  }

  if(ptbin==1)
  {
      Double_t ptbinning[19] = {0., 0.1, 0.3, 0.5, 0.7, 0.9, 1.1, 1.3, 1.5, 2., 2.5, 3., 4., 5., 6., 8., 12., 16., 20.};
  }


  Double_t etabinning[17] = {-0.8, -0.7, -0.6, -0.5, -0.4, -0.3, -0.2, -0.1, 0., 0.1, 0.2, 
			     0.3, 0.4, 0.5, 0.6, 0.7, 0.8};
  AliHFEvarManager *vm = task->GetVarManager();
  Int_t sizept=(sizeof(ptbinning)/sizeof(double))-1;
  Int_t sizeeta=(sizeof(etabinning)/sizeof(double))-1;

  //  printf("ptbinning: %i \n",sizept);

  //vm->AddVariable("pt");
  //vm->AddVariable("eta");
  vm->AddVariable("pt", sizept, ptbinning);
  // vm->AddVariable("pt", 18, ptbinning);
  //  vm->AddVariable("eta", 16, etabinning);
  vm->AddVariable("eta", sizeeta, etabinning);
  vm->AddVariable("phi");
  vm->AddVariable("charge");
  vm->AddVariable("source");
  vm->AddVariable("centrality");
  
  if(!useMC){
    
    for(Int_t a=0;a<12;a++)
      {
	TF1 *hBackground = new TF1("hadronicBackgroundFunction","TMath::Exp([0]/x + [1])", 0., 20.);
	hBackground->SetParameter(0, -43.87);
	hBackground->SetParameter(1, 2.85);
	task->SetBackGroundFactorsFunction(hBackground,a);
      }
  }
  
  // Define PID
  AliHFEpid *pid = task->GetPID();
  if(useMC) pid->SetHasMCData(kTRUE);
  pid->AddDetector("TOF", 0);
  pid->AddDetector("TPC", 1);
  //pid->ConfigureTPCrejection();

 
  if(!useMC){
    
    // 0-10% 0.16
    // 10-20% 0.29
    // 20-30% 0.38
    // 30-40% 0.44
    Double_t paramsTPCdEdxcut[12] ={0.0, 0.0, 0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0};
    if((nondefaultcentr!=0) && tpcdEdxcut) memcpy(paramsTPCdEdxcut,tpcdEdxcut,sizeof(paramsTPCdEdxcut));

    if(etacor==1&&nondefaultcentr==0)
      {
        Double_t tpcdEdxcutetacor[12]={-0.08, 0.016, 0.12,0.22,0.3,0.3,0.3,0.3,0.3,0.3,0.3,0.3};
	memcpy(paramsTPCdEdxcut,tpcdEdxcutetacor,sizeof(paramsTPCdEdxcut));
      }
    
    char *cutmodel;
    cutmodel="pol0";
    
    for(Int_t a=0;a<11;a++)
    {
     //   cout << a << " " << paramsTPCdEdxcut[a] << endl;
	Double_t tpcparam[1]={paramsTPCdEdxcut[a]};

	pid->ConfigureTPCcentralityCut(a,cutmodel,tpcparam,tpcu);
    }

    if(etacor==1)
      { 
	// Apply eta correction
	AliHFEpidTPC *tpcpid = pid->GetDetPID(AliHFEpid::kTPCpid);
	TF1 *etacorrection = GetEtaCorrection(listname);
	if(etacorrection) tpcpid->SetEtaCorrection(etacorrection);
	
      }

  }
  pid->ConfigureTOF(TOFs);

  AliHFEpidTOF *tofpid = pid->GetDetPID(AliHFEpid::kTOFpid);
  if(TOFs<3.) tofpid->SetTOFnSigmaBand(-3,TOFs);

  // QA
  task->SetQAOn(AliAnalysisTaskHFE::kPIDqa);
  //task->SetFillSignalOnly(kFALSE);    // for DE pluging for MC
  if(kMCQA) task->SetQAOn(AliAnalysisTaskHFE::kMCqa);
  //task->SwitchOnPlugin(AliAnalysisTaskHFE::kIsElecBackGround);
  //task->SwitchOnPlugin(AliAnalysisTaskHFE::kSecVtx);
  if(kDEStep) task->SwitchOnPlugin(AliAnalysisTaskHFE::kDEstep);
  if(useMC && addflag==1) task->SetDebugStreaming();

  task->SelectCollisionCandidates(AliVEvent::kMB | AliVEvent::kCentral | AliVEvent::kSemiCentral);

 
  printf("*************************************\n");
  printf("Configuring task PbPb \n"); 
  //if(isLHC10) printf("Configuring TPC1 Task 2010 :\n");
  //if(isLHC11) printf("Configuring TPC1 Task 2011 :\n");
  task->Print();
  pid->PrintStatus();
  printf("*************************************\n"); 
  return task;
}
Commit	Line	Data
a86a8fda	1	TF1* GetEtaCorrection(TString listname){
a86a8fda	2
eeb1e952	3	TString etaMap="$ALICE_ROOT/PWGDQ/dielectron/files/EtaCorrMaps.root";
6b349323	4
6b349323	5	if (gSystem->AccessPathName(gSystem->ExpandPathName(etaMap.Data()))){
eeb1e952	6	Error("ConfigHFEpbpb","Eta map not found: %s",etaMap.Data());
6b349323	7	return 0;
	8	}
	9
	10	TFile f(etaMap.Data());
	11	if (!f.IsOpen()) return 0;
	12	gROOT->cd();
	13	TList *keys=f.GetListOfKeys();
	14
	15	for (Int_t i=0; i<keys->GetEntries(); ++i){
	16	TString kName=keys->At(i)->GetName();
	17	TPRegexp reg(kName);
a86a8fda	18	if (reg.MatchB(listname)){
6b349323	19	printf("Using Eta Correction Function: %s\n",kName.Data());
	20	return (TF1*)f.Get(kName.Data());
	21	}
	22	}
	23	return 0;
	24	}
	25
	26
a86a8fda	27	AliAnalysisTaskHFE* ConfigHFEpbpb(Bool_t useMC=kFALSE,
	28	TString appendix,
	29	Int_t aodfilter=-1,
d84c68e9	30	UChar_t TPCcl=70, UChar_t TPCclPID = 80,
	31	Double_t TPCclRatio = 0.6, Double_t TPCclshared = 1.1,
	32	UChar_t ITScl=3, Double_t ITSchi2perclusters=99999999.,
a86a8fda	33	Int_t itspixelcut=AliHFEextraCuts::kFirst,
fda4fd95	34	Double_t dcaxy=1.0, Double_t dcaz=2.0,
a86a8fda	35	Double_t TOFs=3.,
	36	Double_t IpSig=3.,
	37	Float_t prodlow=0., Float_t prodhigh=100.,
	38	Bool_t beauty=kFALSE,
	39	Bool_t kMCQA = kFALSE,
	40	Bool_t kDEStep = kFALSE,
	41	Int_t addflag=0, Int_t ptbin=0,
	42	Int_t etacor=0, TString listname="",
fda4fd95	43	Int_t nondefaultcentr=0, Float_t* arraycentr=NULL,
fda4fd95	44	Double_t* tpcdEdxcut=NULL,Double_t tpcu=3.0){
d84c68e9	45	//
	46	// HFE standard task configuration
	47	//
d84c68e9	48
	49	AliHFEcuts *hfecuts = new AliHFEcuts(appendix,"HFE cuts pbpb TOF TPC");
	50	hfecuts->CreateStandardCuts();
	51
	52	hfecuts->SetMinNClustersTPC(TPCcl);
	53	hfecuts->SetMinNClustersTPCPID(TPCclPID);
	54	hfecuts->SetMinRatioTPCclusters(TPCclRatio);
	55	hfecuts->SetTPCmodes(AliHFEextraCuts::kFound, AliHFEextraCuts::kFoundOverFindable);
	56	hfecuts->SetFractionOfSharedTPCClusters(TPCclshared);
	57
	58	hfecuts->SetMinNClustersITS(ITScl);
	59	hfecuts->SetMaxChi2perClusterITS(ITSchi2perclusters);
fda4fd95	60	hfecuts->SetCutITSpixel(itspixelcut);
d84c68e9	61	hfecuts->SetCheckITSLayerStatus(kFALSE);
d84c68e9	62
fda4fd95	63	hfecuts->SetIPcutParam(0,0,0,IpSig,kTRUE,kTRUE);
	64	// if(useMC && beauty) hfecuts->SetProductionVertex(prodlow,prodhigh,prodlow,prodhigh);
	65	if(useMC) hfecuts->SetProductionVertex(prodlow,prodhigh,prodlow,prodhigh);
d84c68e9	66
	67	hfecuts->SetMaxImpactParam(dcaxy,dcaz);
	68
	69	// event cuts
	70	hfecuts->SetUseMixedVertex(kTRUE);
	71	hfecuts->SetVertexRange(10.);
	72
	73	// others
	74	hfecuts->SetTOFPIDStep(kTRUE);
	75	//hfecuts->SetMaxChi2perClusterITS(36);
	76	//hfecuts->SetTOFMISMATCHStep(kTRUE);
	77	//hfecuts->SetTPCPIDCleanUpStep(kTRUE);
	78	hfecuts->SetQAOn();
	79
	80	AliAnalysisTaskHFE *task = new AliAnalysisTaskHFE(appendix);
	81	task->SetHFECuts(hfecuts);
	82	task->SetPbPbAnalysis(kTRUE);
fda4fd95	83	task->SetRemovePileUp(kFALSE);
d84c68e9	84	task->GetPIDQAManager()->SetHighResolutionHistos();
a86a8fda	85	if(useMC) task->SetHasMCData(kTRUE); // necessary for AOD
	86	printf("AOD filter %d On/OFF?\n",aodfilter);
	87	if(aodfilter > 0) {
	88	printf("ON AOD filter %d\n",aodfilter);
	89	task->SetUseFilterAOD(kTRUE);
	90	task->SetFilter(aodfilter);
	91	}
fda4fd95	92
	93
	94	if((nondefaultcentr!=0) && arraycentr) {
	95	for(Int_t i=0;i<12;i++)
	96	{
	97	task->SetPbPbUserCentralityLimit(kTRUE);
	98	task->SetPbPbUserCentralityArray(i,arraycentr[i]);
	99
	100	}
	101	}
d84c68e9	102
d84c68e9	103	// Define Variables
fda4fd95	104
	105	if(ptbin==0)
	106	{
	107	Double_t ptbinning[36] = {0., 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1., 1.1, 1.2,
	108	1.3, 1.4, 1.5, 1.75, 2., 2.25, 2.5, 2.75, 3., 3.5, 4., 4.5, 5.,
	109	5.5, 6., 7., 8., 10., 12., 14., 16., 18., 20.};
	110	}
	111
	112	if(ptbin==1)
	113	{
	114	Double_t ptbinning[19] = {0., 0.1, 0.3, 0.5, 0.7, 0.9, 1.1, 1.3, 1.5, 2., 2.5, 3., 4., 5., 6., 8., 12., 16., 20.};
	115	}
	116
	117
d84c68e9	118	Double_t etabinning[17] = {-0.8, -0.7, -0.6, -0.5, -0.4, -0.3, -0.2, -0.1, 0., 0.1, 0.2,
	119	0.3, 0.4, 0.5, 0.6, 0.7, 0.8};
	120	AliHFEvarManager *vm = task->GetVarManager();
fda4fd95	121	Int_t sizept=(sizeof(ptbinning)/sizeof(double))-1;
	122	Int_t sizeeta=(sizeof(etabinning)/sizeof(double))-1;
	123
	124	// printf("ptbinning: %i \n",sizept);
	125
d84c68e9	126	//vm->AddVariable("pt");
d84c68e9	127	//vm->AddVariable("eta");
fda4fd95	128	vm->AddVariable("pt", sizept, ptbinning);
	129	// vm->AddVariable("pt", 18, ptbinning);
	130	// vm->AddVariable("eta", 16, etabinning);
	131	vm->AddVariable("eta", sizeeta, etabinning);
d84c68e9	132	vm->AddVariable("phi");
	133	vm->AddVariable("charge");
	134	vm->AddVariable("source");
	135	vm->AddVariable("centrality");
	136
	137	if(!useMC){
	138
	139	for(Int_t a=0;a<12;a++)
	140	{
	141	TF1 *hBackground = new TF1("hadronicBackgroundFunction","TMath::Exp([0]/x + [1])", 0., 20.);
	142	hBackground->SetParameter(0, -43.87);
	143	hBackground->SetParameter(1, 2.85);
	144	task->SetBackGroundFactorsFunction(hBackground,a);
	145	}
	146	}
	147
	148	// Define PID
	149	AliHFEpid *pid = task->GetPID();
	150	if(useMC) pid->SetHasMCData(kTRUE);
	151	pid->AddDetector("TOF", 0);
	152	pid->AddDetector("TPC", 1);
	153	//pid->ConfigureTPCrejection();
fda4fd95	154
	155
	156
d84c68e9	157	if(!useMC){
d84c68e9	158
fda4fd95	159	// 0-10% 0.16
	160	// 10-20% 0.29
	161	// 20-30% 0.38
	162	// 30-40% 0.44
	163	Double_t paramsTPCdEdxcut[12] ={0.0, 0.0, 0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0};
	164	if((nondefaultcentr!=0) && tpcdEdxcut) memcpy(paramsTPCdEdxcut,tpcdEdxcut,sizeof(paramsTPCdEdxcut));
6b349323	165
	166	if(etacor==1&&nondefaultcentr==0)
	167	{
	168	Double_t tpcdEdxcutetacor[12]={-0.08, 0.016, 0.12,0.22,0.3,0.3,0.3,0.3,0.3,0.3,0.3,0.3};
	169	memcpy(paramsTPCdEdxcut,tpcdEdxcutetacor,sizeof(paramsTPCdEdxcut));
	170	}
fda4fd95	171
d84c68e9	172	char *cutmodel;
	173	cutmodel="pol0";
	174
d84c68e9	175	for(Int_t a=0;a<11;a++)
fda4fd95	176	{
	177	// cout << a << " " << paramsTPCdEdxcut[a] << endl;
	178	Double_t tpcparam[1]={paramsTPCdEdxcut[a]};
	179
	180	pid->ConfigureTPCcentralityCut(a,cutmodel,tpcparam,tpcu);
	181	}
	182
6b349323	183	if(etacor==1)
	184	{
	185	// Apply eta correction
	186	AliHFEpidTPC *tpcpid = pid->GetDetPID(AliHFEpid::kTPCpid);
a86a8fda	187	TF1 *etacorrection = GetEtaCorrection(listname);
6b349323	188	if(etacorrection) tpcpid->SetEtaCorrection(etacorrection);
	189
	190	}
	191
d84c68e9	192	}
	193	pid->ConfigureTOF(TOFs);
	194
fda4fd95	195	AliHFEpidTOF *tofpid = pid->GetDetPID(AliHFEpid::kTOFpid);
	196	if(TOFs<3.) tofpid->SetTOFnSigmaBand(-3,TOFs);
	197
d84c68e9	198	// QA
	199	task->SetQAOn(AliAnalysisTaskHFE::kPIDqa);
	200	//task->SetFillSignalOnly(kFALSE); // for DE pluging for MC
a86a8fda	201	if(kMCQA) task->SetQAOn(AliAnalysisTaskHFE::kMCqa);
d84c68e9	202	//task->SwitchOnPlugin(AliAnalysisTaskHFE::kIsElecBackGround);
d84c68e9	203	//task->SwitchOnPlugin(AliAnalysisTaskHFE::kSecVtx);
a86a8fda	204	if(kDEStep) task->SwitchOnPlugin(AliAnalysisTaskHFE::kDEstep);
fda4fd95	205	if(useMC && addflag==1) task->SetDebugStreaming();
d84c68e9	206
61fb1609	207	task->SelectCollisionCandidates(AliVEvent::kMB \| AliVEvent::kCentral \| AliVEvent::kSemiCentral);
	208
	209
	210
d84c68e9	211	printf("*************************************\n");
	212	printf("Configuring task PbPb \n");
	213	//if(isLHC10) printf("Configuring TPC1 Task 2010 :\n");
	214	//if(isLHC11) printf("Configuring TPC1 Task 2011 :\n");
	215	task->Print();
	216	pid->PrintStatus();
	217	printf("*************************************\n");
	218	return task;
	219	}