]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGLF/totEt/macros/hadEt/GetCorrections.C
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Transition PWG4/JetTasks -> PWGJE and PWG4/totET -> PWGLF/totET
[u/mrichter/AliRoot.git] / PWGLF / totEt / macros / hadEt / GetCorrections.C
CommitLineData
f427cbed 1//Create by Christine Nattrass, Rebecca Scott, Irakli Martashvili
2//University of Tennessee at Knoxville
3//This macro takes an input file created by AliAnalysisHadEtMonteCarlo and creates an AliAnalysisHadEtCorrections for the determination of the corrected Et
4// #include "TFile.h"
5// #include <iostream>
6// #include "AliAnalysisHadEtCorrections.h"
7
8// #include <iostream>
9// #include <TROOT.h>
10// #include <TSystem.h>
11// #include "TStopwatch.h"
12
48a07264 13Float_t CorrNeutral(float ptcut, char *prodname, char *shortprodname, bool ispp = true, bool forSim = true, bool TPC, bool hadronic = false, float etacut = 0.7);
14TH1D *GetHistoCorrNeutral(float cut, char *name, bool ispp, bool forSim, int mycase, bool eta, int color, int marker, bool hadronic = false);
f427cbed 15
48a07264 16Float_t CorrPtCut(float ptcut, char *prodname = "Enter Production Name", char *shortprodname = "EnterProductionName", bool ispp = true, bool forSim = true, int mycase = 0);
17TH1D *GetHistoCorrPtCut(float ptcut = 0.15, char *name, bool ispp = true, bool forSim = true, int mycase);
f427cbed 18
48a07264 19TH1D *GetHistoCorrNotID(float etacut,char *name, bool TPC, bool eta, bool ispp = true, bool forSim = true);
20TH1D *CorrNotID(float etacut,char *name, char *prodname, char *shortprodname, bool TPC, bool ispp = true, bool forSim = true);
21Float_t CorrNotIDConst(float ptcut, float etacut,char *name, char *prodname, char *shortprodname, bool TPC, bool ispp, bool forSim);
f427cbed 22
48a07264 23TH1D *GetHistoNoID(float etacut, char *name, bool eta, bool TPC, bool ispp, bool forSim);
24TH1D *CorrNoID(float etacut,char *name, char *prodname, char *shortprodname, bool ispp, bool forSim);
25Float_t CorrNoIDConst(float etacut, float ptcut,char *name, char *prodname, char *shortprodname, bool ispp, bool forSim);
b64de20c 26
27TH1D* bayneseffdiv(TH1D* numerator, TH1D* denominator,Char_t* name);
48a07264 28TH1D *GetHistoEfficiency(float cut, char *name, int mycase, int color, int marker,bool TPC, bool ITS, int cb = -1, int cblast = -1);
29void CorrEfficiencyPlots(bool TPC, char *prodname, char *shortprodname);
30
31TH1D *GetHistoCorrBkgd(float etacut,char *name, bool TPC,bool ispp,bool forSim);
32void CorrBkgdPlots(char *prodname, char *shortprodname, bool TPC,bool ispp,bool forSim);
33
34//Some variables that we'll use multiple times. We'll declare them here since they don't seem to delete right in the functions
35char prefix[100];
36char histoname[100];
37char epsname[100];
38char pngname[100];
39TFile *file = NULL;//initiated in main function
40const char *mynameTPC = "TPC";
41const char *mynameITS = "ITS";
42const char *mynameTPCITS = "TPCITS";
43const char *detectorEMCAL = "EMCAL";
44const char *detectorPHOS = "PHOS";
45const char *reweightedNo = "";
46const char *reweightedYes = "Reweighted";
b64de20c 47
48//===========================================================================================
f427cbed 49
48a07264 50void GetCorrections(char *prodname = "Enter Production Name", char *shortprodname = "EnterProductionName", bool ispp = true, bool forSim = true, bool TPC = true, char *infilename="Et.ESD.new.sim.merged.root", int dataset = 2009){
f427cbed 51 TStopwatch timer;
52 timer.Start();
53 gSystem->Load("libTree.so");
54 gSystem->Load("libGeom.so");
55 gSystem->Load("libVMC.so");
56 gSystem->Load("libXMLIO.so");
57
58 gSystem->Load("libSTEERBase.so");
59 gSystem->Load("libESD.so");
60 gSystem->Load("libAOD.so");
61
62 gSystem->Load("libANALYSIS");
63 gSystem->Load("libANALYSISalice");
64
65 gSystem->AddIncludePath("-I$ALICE_ROOT/include");
0e866ddc 66 gROOT->ProcessLine(".L AliAnalysisEtCuts.cxx+g");
f427cbed 67 gROOT->ProcessLine(".L AliAnalysisHadEtCorrections.cxx+g");
48a07264 68 file = new TFile(infilename);
f427cbed 69
a02dfa56 70 char outfilename[200];
f43fc416 71 char *sim = "ForData";
72 if(forSim) sim = "ForSimulations";
73 char *system = "PbPb";
74 if(ispp) system = "pp";
0add53f6 75 sprintf(outfilename,"rootFiles/corrections/corrections.%s.%s.%s.root",shortprodname,system,sim);
f427cbed 76 TFile *outfile = new TFile(outfilename,"RECREATE");
77 AliAnalysisHadEtCorrections *hadCorrectionEMCAL = new AliAnalysisHadEtCorrections();
78 hadCorrectionEMCAL->SetName("hadCorrectionEMCAL");
b64de20c 79 float etacut = 0.7;
80 hadCorrectionEMCAL->SetEtaCut(etacut);
48a07264 81 hadCorrectionEMCAL->IsData(!forSim);
82 hadCorrectionEMCAL->IsEMCal(kTRUE);
83 hadCorrectionEMCAL->SetProduction(shortprodname);
84 hadCorrectionEMCAL->SetProductionDescription(prodname);
85 hadCorrectionEMCAL->SetDataSet(dataset);
b64de20c 86 //float etacut = hadCorrectionEMCAL->GetEtaCut();
87 //cout<<"eta cut is "<<etacut<<endl;
f427cbed 88 hadCorrectionEMCAL->SetAcceptanceCorrectionFull(1.0);
48a07264 89 cout<<"Warning: Acceptance corrections will have to be updated to include real acceptance maps of the EMCAL"<<endl;
f427cbed 90 hadCorrectionEMCAL->SetAcceptanceCorrectionPHOS(360.0/60.0);
91 hadCorrectionEMCAL->SetAcceptanceCorrectionEMCAL(360.0/60.0);
b64de20c 92
f427cbed 93 float ptcut = 0.1;
48a07264 94 float neutralCorr = CorrNeutral(ptcut,prodname,shortprodname,ispp,forSim,TPC,false,etacut);
f427cbed 95 hadCorrectionEMCAL->SetNeutralCorrection(neutralCorr);
48a07264 96 //Using error from data, see analysis note for details
97 if(ispp){
98 hadCorrectionEMCAL->SetNeutralCorrectionLowBound(neutralCorr*(1.0-0.013));
99 hadCorrectionEMCAL->SetNeutralCorrectionHighBound(neutralCorr*(1.0+0.013));
100 }
101 else{
102 hadCorrectionEMCAL->SetNeutralCorrectionLowBound(neutralCorr*(1.0-0.049));
103 hadCorrectionEMCAL->SetNeutralCorrectionHighBound(neutralCorr*(1.0+0.049));
104 }
105
106 float hadronicCorr = CorrNeutral(ptcut,prodname,shortprodname,ispp,forSim,TPC,true,etacut);
b64de20c 107 hadCorrectionEMCAL->SetNotHadronicCorrection(hadronicCorr);
48a07264 108 if(ispp){
109 hadCorrectionEMCAL->SetNotHadronicCorrectionLowBound(hadronicCorr*(1.0-0.008));
110 hadCorrectionEMCAL->SetNotHadronicCorrectionHighBound(hadronicCorr*(1.0+0.008));
111 }
112 else{
113 hadCorrectionEMCAL->SetNotHadronicCorrectionLowBound(hadronicCorr*(1.0-0.023));
114 hadCorrectionEMCAL->SetNotHadronicCorrectionHighBound(hadronicCorr*(1.0+0.023));
115 }
116
117 float ptcutITS = CorrPtCut(0.1,prodname,shortprodname,ispp,forSim);
f427cbed 118 hadCorrectionEMCAL->SetpTCutCorrectionITS(ptcutITS);
48a07264 119 float ptcutTPC = CorrPtCut(0.15,prodname,shortprodname,ispp,forSim);
f427cbed 120 hadCorrectionEMCAL->SetpTCutCorrectionTPC(ptcutTPC);
48a07264 121 float ptcutITSLow = CorrPtCut(0.1,prodname,shortprodname,ispp,forSim,-1);
122 float ptcutTPCLow = CorrPtCut(0.15,prodname,shortprodname,ispp,forSim,-1);
123 hadCorrectionEMCAL->SetpTCutCorrectionITSLowBound(ptcutITSLow);
124 hadCorrectionEMCAL->SetpTCutCorrectionTPCLowBound(ptcutTPCLow);
125 float ptcutITSHigh = CorrPtCut(0.1,prodname,shortprodname,ispp,forSim,1);
126 float ptcutTPCHigh = CorrPtCut(0.15,prodname,shortprodname,ispp,forSim,1);
127 hadCorrectionEMCAL->SetpTCutCorrectionITSHighBound(ptcutITSHigh);
128 hadCorrectionEMCAL->SetpTCutCorrectionTPCHighBound(ptcutTPCHigh);
129
130 TH1D *NotIDTPC = CorrNotID(etacut,"CorrNotIDEMCALTPC",prodname,shortprodname,true,ispp,forSim);
131 TH1D *NotIDITS = CorrNotID(etacut,"CorrNotIDEMCALITS",prodname,shortprodname,false,ispp,forSim);
b64de20c 132 hadCorrectionEMCAL->SetNotIDCorrectionTPC(NotIDTPC);
133 hadCorrectionEMCAL->SetNotIDCorrectionITS(NotIDITS);
134
48a07264 135 Float_t NotIDConstTPC = CorrNotIDConst(0.15,etacut,"CorrNotIDEMCALTPC2",prodname,shortprodname,true,ispp,forSim);
136 Float_t NotIDConstITS = CorrNotIDConst(0.10,etacut,"CorrNotIDEMCALTPC2",prodname,shortprodname,true,ispp,forSim);
0e866ddc 137 hadCorrectionEMCAL->SetNotIDConstCorrectionTPC(1.0/NotIDConstTPC);
138 hadCorrectionEMCAL->SetNotIDConstCorrectionITS(1.0/NotIDConstITS);
48a07264 139 if(ispp){
0add53f6 140 hadCorrectionEMCAL->SetNotIDConstCorrectionTPCLowBound(1./NotIDConstTPC*(1.0-0.010));
141 hadCorrectionEMCAL->SetNotIDConstCorrectionITSLowBound(1./NotIDConstITS*(1.0-0.010));
142 hadCorrectionEMCAL->SetNotIDConstCorrectionTPCHighBound(1./NotIDConstTPC*(1.0+0.010));
143 hadCorrectionEMCAL->SetNotIDConstCorrectionITSHighBound(1./NotIDConstITS*(1.0+0.010));
48a07264 144 }
145 else{
0add53f6 146 hadCorrectionEMCAL->SetNotIDConstCorrectionTPCLowBound(1./NotIDConstTPC*(1.0-0.022));
147 hadCorrectionEMCAL->SetNotIDConstCorrectionITSLowBound(1./NotIDConstITS*(1.0-0.022));
148 hadCorrectionEMCAL->SetNotIDConstCorrectionTPCHighBound(1./NotIDConstTPC*(1.0+0.022));
149 hadCorrectionEMCAL->SetNotIDConstCorrectionITSHighBound(1./NotIDConstITS*(1.0+0.022));
48a07264 150 }
151
152 TH1D *NoID = CorrNoID(etacut,"CorrNoIDEMCAL",prodname,shortprodname,ispp,forSim);
b64de20c 153 hadCorrectionEMCAL->SetNotIDCorrectionNoPID(NoID);
154
48a07264 155 Float_t NoIDTPC = CorrNoIDConst(etacut,0.15,"CorrNoIDEMCAL2",prodname,shortprodname,ispp,forSim);
156 Float_t NoIDITS = CorrNoIDConst(etacut,0.1,"CorrNoIDEMCAL2",prodname,shortprodname,ispp,forSim);
0e866ddc 157 hadCorrectionEMCAL->SetNotIDConstCorrectionTPCNoID(1./NoIDTPC);
158 hadCorrectionEMCAL->SetNotIDConstCorrectionITSNoID(1./NoIDITS);
48a07264 159 hadCorrectionEMCAL->SetNotIDConstCorrectionTPCNoIDLowBound(1./NoIDTPC*.98);
160 hadCorrectionEMCAL->SetNotIDConstCorrectionITSNoIDLowBound(1./NoIDITS*.98);
161 hadCorrectionEMCAL->SetNotIDConstCorrectionTPCNoIDHighBound(1./NoIDTPC*1.02);
162 hadCorrectionEMCAL->SetNotIDConstCorrectionITSNoIDHighBound(1./NoIDITS*1.02);
0e866ddc 163
48a07264 164 TH1D *efficiencyPionTPC = GetHistoEfficiency(etacut,"hEfficiencyPionTPC",1,1,20,true,true);
a02dfa56 165 hadCorrectionEMCAL->SetEfficiencyPionTPC(efficiencyPionTPC);
48a07264 166 TH1D *efficiencyKaonTPC = GetHistoEfficiency(etacut,"hEfficiencyKaonTPC",2,1,20,true,true);
a02dfa56 167 hadCorrectionEMCAL->SetEfficiencyKaonTPC(efficiencyKaonTPC);
48a07264 168 TH1D *efficiencyProtonTPC = GetHistoEfficiency(etacut,"hEfficiencyProtonTPC",3,1,20,true,true);
a02dfa56 169 hadCorrectionEMCAL->SetEfficiencyProtonTPC(efficiencyProtonTPC);
48a07264 170 TH1D *efficiencyHadronTPC = GetHistoEfficiency(etacut,"hEfficiencyHadronTPC",0,1,20,true,true);
a02dfa56 171 hadCorrectionEMCAL->SetEfficiencyHadronTPC(efficiencyHadronTPC);
48a07264 172 TH1D *efficiencyPionITS = GetHistoEfficiency(etacut,"hEfficiencyPionITS",1,1,20,false,true);
a02dfa56 173 hadCorrectionEMCAL->SetEfficiencyPionITS(efficiencyPionITS);
48a07264 174 TH1D *efficiencyKaonITS = GetHistoEfficiency(etacut,"hEfficiencyKaonITS",2,1,20,false,true);
a02dfa56 175 hadCorrectionEMCAL->SetEfficiencyKaonITS(efficiencyKaonITS);
48a07264 176 TH1D *efficiencyProtonITS = GetHistoEfficiency(etacut,"hEfficiencyProtonITS",3,1,20,false,true);
a02dfa56 177 hadCorrectionEMCAL->SetEfficiencyProtonITS(efficiencyProtonITS);
48a07264 178 TH1D *efficiencyHadronITS = GetHistoEfficiency(etacut,"hEfficiencyHadronITS",0,1,20,false,true);
a02dfa56 179 hadCorrectionEMCAL->SetEfficiencyHadronITS(efficiencyHadronITS);
180
48a07264 181 if(!ispp){
182 TH1D *efficiencyPionTPCCB0 = GetHistoEfficiency(etacut,"hEfficiencyPionTPCCB0",1,1,20,true,true,0,4);
183 for(int i=0;i<=4;i++) hadCorrectionEMCAL->SetEfficiencyPionTPC((TH1D*)efficiencyPionTPCCB0->Clone(Form("Test%i",i)),i);
184 TH1D *efficiencyPionTPCCB5 = GetHistoEfficiency(etacut,"hEfficiencyPionTPCCB5",1,1,20,true,true,5,9);
185 for(int i=5;i<=9;i++) hadCorrectionEMCAL->SetEfficiencyPionTPC((TH1D*)efficiencyPionTPCCB5->Clone(Form("Test%i",i)),i);
186 TH1D *efficiencyPionTPCCB10 = GetHistoEfficiency(etacut,"hEfficiencyPionTPCCB10",1,1,20,true,true,10,15);
187 for(int i=10;i<=19;i++) hadCorrectionEMCAL->SetEfficiencyPionTPC((TH1D*)efficiencyPionTPCCB10->Clone(Form("Test%i",i)),i);
188
189 TH1D *efficiencyKaonTPCCB0 = GetHistoEfficiency(etacut,"hEfficiencyKaonTPCCB0",2,1,20,true,true,0,4);
190 for(int i=0;i<=4;i++) hadCorrectionEMCAL->SetEfficiencyKaonTPC((TH1D*)efficiencyKaonTPCCB0->Clone(Form("Test%i",i)),i);
191 TH1D *efficiencyKaonTPCCB5 = GetHistoEfficiency(etacut,"hEfficiencyKaonTPCCB5",2,1,20,true,true,5,9);
192 for(int i=5;i<=9;i++) hadCorrectionEMCAL->SetEfficiencyKaonTPC((TH1D*)efficiencyKaonTPCCB5->Clone(Form("Test%i",i)),i);
193 TH1D *efficiencyKaonTPCCB10 = GetHistoEfficiency(etacut,"hEfficiencyKaonTPCCB10",2,1,20,true,true,10,15);
194 for(int i=10;i<=19;i++) hadCorrectionEMCAL->SetEfficiencyKaonTPC((TH1D*)efficiencyKaonTPCCB10->Clone(Form("Test%i",i)),i);//Kaon
195
196 TH1D *efficiencyProtonTPCCB0 = GetHistoEfficiency(etacut,"hEfficiencyProtonTPCCB0",3,1,20,true,true,0,4);
197 for(int i=0;i<=4;i++) hadCorrectionEMCAL->SetEfficiencyProtonTPC((TH1D*)efficiencyProtonTPCCB0->Clone(Form("Test%i",i)),i);
198 TH1D *efficiencyProtonTPCCB5 = GetHistoEfficiency(etacut,"hEfficiencyProtonTPCCB5",3,1,20,true,true,5,9);
199 for(int i=5;i<=9;i++) hadCorrectionEMCAL->SetEfficiencyProtonTPC((TH1D*)efficiencyProtonTPCCB5->Clone(Form("Test%i",i)),i);
200 TH1D *efficiencyProtonTPCCB10 = GetHistoEfficiency(etacut,"hEfficiencyProtonTPCCB10",3,1,20,true,true,10,15);
201 for(int i=10;i<=19;i++) hadCorrectionEMCAL->SetEfficiencyProtonTPC((TH1D*)efficiencyProtonTPCCB10->Clone(Form("Test%i",i)),i);//Proton
202
203 TH1D *efficiencyHadronTPCCB0 = GetHistoEfficiency(etacut,"hEfficiencyHadronTPCCB0",0,1,20,true,true,0,4);
204 for(int i=0;i<=4;i++) hadCorrectionEMCAL->SetEfficiencyHadronTPC((TH1D*)efficiencyHadronTPCCB0->Clone(Form("Test%i",i)),i);
205 TH1D *efficiencyHadronTPCCB5 = GetHistoEfficiency(etacut,"hEfficiencyHadronTPCCB5",0,1,20,true,true,5,9);
206 for(int i=5;i<=9;i++) hadCorrectionEMCAL->SetEfficiencyHadronTPC((TH1D*)efficiencyHadronTPCCB5->Clone(Form("Test%i",i)),i);
207 TH1D *efficiencyHadronTPCCB10 = GetHistoEfficiency(etacut,"hEfficiencyHadronTPCCB10",0,1,20,true,true,10,15);
208 for(int i=10;i<=19;i++) hadCorrectionEMCAL->SetEfficiencyHadronTPC((TH1D*)efficiencyHadronTPCCB10->Clone(Form("Test%i",i)),i);//Hadron
209
210
211 TH1D *efficiencyPionITSCB0 = GetHistoEfficiency(etacut,"hEfficiencyPionITSCB0",1,1,20,false,true,0,4);
212 for(int i=0;i<=4;i++) hadCorrectionEMCAL->SetEfficiencyPionITS((TH1D*)efficiencyPionITSCB0->Clone(Form("Test%i",i)),i);
213 TH1D *efficiencyPionITSCB5 = GetHistoEfficiency(etacut,"hEfficiencyPionITSCB5",1,1,20,false,true,5,9);
214 for(int i=5;i<=9;i++) hadCorrectionEMCAL->SetEfficiencyPionITS((TH1D*)efficiencyPionITSCB5->Clone(Form("Test%i",i)),i);
215 TH1D *efficiencyPionITSCB10 = GetHistoEfficiency(etacut,"hEfficiencyPionITSCB10",1,1,20,false,true,10,15);
216 for(int i=10;i<=19;i++) hadCorrectionEMCAL->SetEfficiencyPionITS((TH1D*)efficiencyPionITSCB10->Clone(Form("Test%i",i)),i);//Pion
217
218 TH1D *efficiencyKaonITSCB0 = GetHistoEfficiency(etacut,"hEfficiencyKaonITSCB0",2,1,20,false,true,0,4);
219 for(int i=0;i<=4;i++) hadCorrectionEMCAL->SetEfficiencyKaonITS((TH1D*)efficiencyKaonITSCB0->Clone(Form("Test%i",i)),i);
220 TH1D *efficiencyKaonITSCB5 = GetHistoEfficiency(etacut,"hEfficiencyKaonITSCB5",2,1,20,false,true,5,9);
221 for(int i=5;i<=9;i++) hadCorrectionEMCAL->SetEfficiencyKaonITS((TH1D*)efficiencyKaonITSCB5->Clone(Form("Test%i",i)),i);
222 TH1D *efficiencyKaonITSCB10 = GetHistoEfficiency(etacut,"hEfficiencyKaonITSCB10",2,1,20,false,true,10,15);
223 for(int i=10;i<=19;i++) hadCorrectionEMCAL->SetEfficiencyKaonITS((TH1D*)efficiencyKaonITSCB10->Clone(Form("Test%i",i)),i);//Kaon
224
225 TH1D *efficiencyProtonITSCB0 = GetHistoEfficiency(etacut,"hEfficiencyProtonITSCB0",3,1,20,false,true,0,4);
226 for(int i=0;i<=4;i++) hadCorrectionEMCAL->SetEfficiencyProtonITS((TH1D*)efficiencyProtonITSCB0->Clone(Form("Test%i",i)),i);
227 TH1D *efficiencyProtonITSCB5 = GetHistoEfficiency(etacut,"hEfficiencyProtonITSCB5",3,1,20,false,true,5,9);
228 for(int i=5;i<=9;i++) hadCorrectionEMCAL->SetEfficiencyProtonITS((TH1D*)efficiencyProtonITSCB5->Clone(Form("Test%i",i)),i);
229 TH1D *efficiencyProtonITSCB10 = GetHistoEfficiency(etacut,"hEfficiencyProtonITSCB10",3,1,20,false,true,10,15);
230 for(int i=10;i<=19;i++) hadCorrectionEMCAL->SetEfficiencyProtonITS((TH1D*)efficiencyProtonITSCB10->Clone(Form("Test%i",i)),i);//Proton
231
232 TH1D *efficiencyHadronITSCB0 = GetHistoEfficiency(etacut,"hEfficiencyHadronITSCB0",0,1,20,false,true,0,4);
233 for(int i=0;i<=4;i++) hadCorrectionEMCAL->SetEfficiencyHadronITS((TH1D*)efficiencyHadronITSCB0->Clone(Form("Test%i",i)),i);
234 TH1D *efficiencyHadronITSCB5 = GetHistoEfficiency(etacut,"hEfficiencyHadronITSCB5",0,1,20,false,true,5,9);
235 for(int i=5;i<=9;i++) hadCorrectionEMCAL->SetEfficiencyHadronITS((TH1D*)efficiencyHadronITSCB5->Clone(Form("Test%i",i)),i);
236 TH1D *efficiencyHadronITSCB10 = GetHistoEfficiency(etacut,"hEfficiencyHadronITSCB10",0,1,20,false,true,10,15);
237 for(int i=10;i<=19;i++) hadCorrectionEMCAL->SetEfficiencyHadronITS((TH1D*)efficiencyHadronITSCB10->Clone(Form("Test%i",i)),i);//Hadron
238 }//EMCAL
239
240 //CorrEfficiencyPlots(true,prodname,shortprodname);
a02dfa56 241 //CorrEfficiencyPlots(false,prodname,shortprodname,infilename);
b64de20c 242
48a07264 243
244 //hadCorrectionEMCAL->GetEfficiencyHadronTPC()->Draw();
245 TH1D *backgroundTPC = GetHistoCorrBkgd(etacut,"hBackgroundTPC",true,ispp,forSim);
246 TH1D *backgroundITS = GetHistoCorrBkgd(etacut,"hBackgroundITS",false,ispp,forSim);
b64de20c 247 hadCorrectionEMCAL->SetBackgroundCorrectionTPC(backgroundTPC);
248 hadCorrectionEMCAL->SetBackgroundCorrectionITS(backgroundITS);
48a07264 249 //CorrBkgdPlots(prodname,shortprodname,true,ispp,forSim);
250 //CorrBkgdPlots(prodname,shortprodname,false,ispp,forSim);
251
252 hadCorrectionEMCAL->Report();
f427cbed 253
f427cbed 254 outfile->cd();
255 hadCorrectionEMCAL->Write();
256 outfile->Write();
a02dfa56 257 delete hadCorrectionEMCAL;
258
259 AliAnalysisHadEtCorrections *hadCorrectionPHOS = new AliAnalysisHadEtCorrections();
260 hadCorrectionPHOS->SetName("hadCorrectionPHOS");
261 float etacut = 0.12;
262 hadCorrectionPHOS->SetEtaCut(etacut);
48a07264 263 hadCorrectionPHOS->IsData(!forSim);
264 hadCorrectionPHOS->IsEMCal(kTRUE);
265 hadCorrectionPHOS->SetProduction(shortprodname);
266 hadCorrectionPHOS->SetProductionDescription(prodname);
267 hadCorrectionPHOS->SetDataSet(dataset);
a02dfa56 268 //float etacut = hadCorrectionPHOS->GetEtaCut();
269 //cout<<"eta cut is "<<etacut<<endl;
a02dfa56 270 hadCorrectionPHOS->SetAcceptanceCorrectionFull(1.0);
48a07264 271 cout<<"Warning: Acceptance corrections will have to be updated to include real acceptance maps of the PHOS"<<endl;
a02dfa56 272 hadCorrectionPHOS->SetAcceptanceCorrectionPHOS(360.0/60.0);
273 hadCorrectionPHOS->SetAcceptanceCorrectionEMCAL(360.0/60.0);
274
275 float ptcut = 0.1;
48a07264 276 float neutralCorr = CorrNeutral(ptcut,prodname,shortprodname,ispp,forSim,TPC,false,etacut);
a02dfa56 277 hadCorrectionPHOS->SetNeutralCorrection(neutralCorr);
48a07264 278 //Using error from data, see analysis note for details
279 if(ispp){
0add53f6 280 hadCorrectionPHOS->SetNeutralCorrectionLowBound(neutralCorr*(1.0-.013));
281 hadCorrectionPHOS->SetNeutralCorrectionHighBound(neutralCorr*(1.0+.013));
48a07264 282 }
283 else{
0add53f6 284 hadCorrectionPHOS->SetNeutralCorrectionLowBound(neutralCorr*(1.0-0.049));
285 hadCorrectionPHOS->SetNeutralCorrectionHighBound(neutralCorr*(1.0+0.049));
48a07264 286 }
287
288 float hadronicCorr = CorrNeutral(ptcut,prodname,shortprodname,ispp,forSim,TPC,true,etacut);
a02dfa56 289 hadCorrectionPHOS->SetNotHadronicCorrection(hadronicCorr);
48a07264 290 if(ispp){
0add53f6 291 hadCorrectionPHOS->SetNotHadronicCorrectionLowBound(neutralCorr*(1.0-0.008));
292 hadCorrectionPHOS->SetNotHadronicCorrectionHighBound(neutralCorr*(1.0+0.008));
48a07264 293 }
294 else{
0add53f6 295 hadCorrectionPHOS->SetNotHadronicCorrectionLowBound(neutralCorr*(1.0-0.023));
296 hadCorrectionPHOS->SetNotHadronicCorrectionHighBound(neutralCorr*(1.0+0.023));
48a07264 297 }
298
299 float ptcutITS = CorrPtCut(0.1,prodname,shortprodname,ispp,forSim);
a02dfa56 300 hadCorrectionPHOS->SetpTCutCorrectionITS(ptcutITS);
48a07264 301 float ptcutTPC = CorrPtCut(0.15,prodname,shortprodname,ispp,forSim);
a02dfa56 302 hadCorrectionPHOS->SetpTCutCorrectionTPC(ptcutTPC);
48a07264 303
304 float ptcutITSLow = CorrPtCut(0.1,prodname,shortprodname,ispp,forSim,-1);
305 float ptcutTPCLow = CorrPtCut(0.15,prodname,shortprodname,ispp,forSim,-1);
306 hadCorrectionPHOS->SetpTCutCorrectionITSLowBound(ptcutITSLow);
307 hadCorrectionPHOS->SetpTCutCorrectionTPCLowBound(ptcutTPCLow);
308 float ptcutITSHigh = CorrPtCut(0.1,prodname,shortprodname,ispp,forSim,1);
309 float ptcutTPCHigh = CorrPtCut(0.15,prodname,shortprodname,ispp,forSim,1);
310 hadCorrectionPHOS->SetpTCutCorrectionITSHighBound(ptcutITSHigh);
311 hadCorrectionPHOS->SetpTCutCorrectionTPCHighBound(ptcutTPCHigh);
312
313 TH1D *NotIDTPC = CorrNotID(etacut,"CorrNotIDPHOSTPC",prodname,shortprodname,true,ispp,forSim);
314 TH1D *NotIDITS = CorrNotID(etacut,"CorrNotIDPHOSITS",prodname,shortprodname,false,ispp,forSim);
a02dfa56 315 hadCorrectionPHOS->SetNotIDCorrectionTPC(NotIDTPC);
316 hadCorrectionPHOS->SetNotIDCorrectionITS(NotIDITS);
317
48a07264 318 Float_t NotIDConstTPC = CorrNotIDConst(0.15,etacut,"CorrNotIDPHOSTPC2",prodname,shortprodname,true,ispp,forSim);
319 Float_t NotIDConstITS = CorrNotIDConst(0.10,etacut,"CorrNotIDPHOSTPC2",prodname,shortprodname,true,ispp,forSim);
0e866ddc 320 hadCorrectionPHOS->SetNotIDConstCorrectionTPC(1./NotIDConstTPC);
321 hadCorrectionPHOS->SetNotIDConstCorrectionITS(1./NotIDConstITS);
48a07264 322 if(ispp){
0add53f6 323 hadCorrectionPHOS->SetNotIDConstCorrectionTPCLowBound(1./NotIDConstTPC*(1.0-0.010));
324 hadCorrectionPHOS->SetNotIDConstCorrectionITSLowBound(1./NotIDConstITS*(1.0-0.010));
325 hadCorrectionPHOS->SetNotIDConstCorrectionTPCHighBound(1./NotIDConstTPC*(1.0+0.010));
326 hadCorrectionPHOS->SetNotIDConstCorrectionITSHighBound(1./NotIDConstITS*(1.0+0.010));
48a07264 327 }
328 else{
0add53f6 329 hadCorrectionPHOS->SetNotIDConstCorrectionTPCLowBound(1./NotIDConstTPC*(1.0-0.022));
330 hadCorrectionPHOS->SetNotIDConstCorrectionITSLowBound(1./NotIDConstITS*(1.0-0.022));
331 hadCorrectionPHOS->SetNotIDConstCorrectionTPCHighBound(1./NotIDConstTPC*(1.0+0.022));
332 hadCorrectionPHOS->SetNotIDConstCorrectionITSHighBound(1./NotIDConstITS*(1.0+0.022));
48a07264 333 }
334
335
336 TH1D *NoID = CorrNoID(etacut,"CorrNoIDPHOS",prodname,shortprodname,ispp,forSim);
a02dfa56 337 hadCorrectionPHOS->SetNotIDCorrectionNoPID(NoID);
338
0e866ddc 339
48a07264 340 Float_t NoIDTPC = CorrNoIDConst(etacut,0.15,"CorrNoIDPHOS2",prodname,shortprodname,ispp,forSim);
341 Float_t NoIDITS = CorrNoIDConst(etacut,0.1,"CorrNoIDPHOS2",prodname,shortprodname,ispp,forSim);
0e866ddc 342 hadCorrectionPHOS->SetNotIDConstCorrectionTPCNoID(1./NoIDTPC);
343 hadCorrectionPHOS->SetNotIDConstCorrectionITSNoID(1./NoIDITS);
48a07264 344 hadCorrectionPHOS->SetNotIDConstCorrectionTPCNoIDLowBound(1./NoIDTPC*.98);
345 hadCorrectionPHOS->SetNotIDConstCorrectionITSNoIDLowBound(1./NoIDITS*.98);
346 hadCorrectionPHOS->SetNotIDConstCorrectionTPCNoIDHighBound(1./NoIDTPC*1.02);
347 hadCorrectionPHOS->SetNotIDConstCorrectionITSNoIDHighBound(1./NoIDITS*1.02);
348
349 TH1D *efficiencyPionTPC = GetHistoEfficiency(etacut,"hEfficiencyPionTPC",1,1,20,true,true);
350 TH1D *efficiencyKaonTPC = GetHistoEfficiency(etacut,"hEfficiencyKaonTPC",2,1,20,true,true);
351 TH1D *efficiencyProtonTPC = GetHistoEfficiency(etacut,"hEfficiencyProtonTPC",3,1,20,true,true);
352 TH1D *efficiencyHadronTPC = GetHistoEfficiency(etacut,"hEfficiencyHadronTPC",0,1,20,true,true);
353 TH1D *efficiencyPionITS = GetHistoEfficiency(etacut,"hEfficiencyPionITS",1,1,20,false,true);
354 TH1D *efficiencyKaonITS = GetHistoEfficiency(etacut,"hEfficiencyKaonITS",2,1,20,false,true);
355 TH1D *efficiencyProtonITS = GetHistoEfficiency(etacut,"hEfficiencyProtonITS",3,1,20,false,true);
356 TH1D *efficiencyHadronITS = GetHistoEfficiency(etacut,"hEfficiencyHadronITS",0,1,20,false,true);
357 //CorrEfficiencyPlots(true,prodname,shortprodname);
358 //CorrEfficiencyPlots(false,prodname,shortprodname);
a02dfa56 359 hadCorrectionPHOS->SetEfficiencyPionTPC(efficiencyPionTPC);
360 hadCorrectionPHOS->SetEfficiencyKaonTPC(efficiencyKaonTPC);
361 hadCorrectionPHOS->SetEfficiencyProtonTPC(efficiencyProtonTPC);
362 hadCorrectionPHOS->SetEfficiencyHadronTPC(efficiencyHadronTPC);
363 hadCorrectionPHOS->SetEfficiencyPionITS(efficiencyPionITS);
364 hadCorrectionPHOS->SetEfficiencyKaonITS(efficiencyKaonITS);
365 hadCorrectionPHOS->SetEfficiencyProtonITS(efficiencyProtonITS);
366 hadCorrectionPHOS->SetEfficiencyHadronITS(efficiencyHadronITS);
367
48a07264 368
369 if(!ispp){
370 TH1D *efficiencyPionTPCCB0 = GetHistoEfficiency(etacut,"hEfficiencyPionTPCCB0",1,1,20,true,true,0,4);
371 for(int i=0;i<=4;i++) hadCorrectionPHOS->SetEfficiencyPionTPC((TH1D*)efficiencyPionTPCCB0->Clone(Form("Test%i",i)),i);
372 TH1D *efficiencyPionTPCCB5 = GetHistoEfficiency(etacut,"hEfficiencyPionTPCCB5",1,1,20,true,true,5,9);
373 for(int i=5;i<=9;i++) hadCorrectionPHOS->SetEfficiencyPionTPC((TH1D*)efficiencyPionTPCCB5->Clone(Form("Test%i",i)),i);
374 TH1D *efficiencyPionTPCCB10 = GetHistoEfficiency(etacut,"hEfficiencyPionTPCCB10",1,1,20,true,true,10,15);
375 for(int i=10;i<=19;i++) hadCorrectionPHOS->SetEfficiencyPionTPC((TH1D*)efficiencyPionTPCCB10->Clone(Form("Test%i",i)),i);
376
377 TH1D *efficiencyKaonTPCCB0 = GetHistoEfficiency(etacut,"hEfficiencyKaonTPCCB0",2,1,20,true,true,0,4);
378 for(int i=0;i<=4;i++) hadCorrectionPHOS->SetEfficiencyKaonTPC((TH1D*)efficiencyKaonTPCCB0->Clone(Form("Test%i",i)),i);
379 TH1D *efficiencyKaonTPCCB5 = GetHistoEfficiency(etacut,"hEfficiencyKaonTPCCB5",2,1,20,true,true,5,9);
380 for(int i=5;i<=9;i++) hadCorrectionPHOS->SetEfficiencyKaonTPC((TH1D*)efficiencyKaonTPCCB5->Clone(Form("Test%i",i)),i);
381 TH1D *efficiencyKaonTPCCB10 = GetHistoEfficiency(etacut,"hEfficiencyKaonTPCCB10",2,1,20,true,true,10,15);
382 for(int i=10;i<=19;i++) hadCorrectionPHOS->SetEfficiencyKaonTPC((TH1D*)efficiencyKaonTPCCB10->Clone(Form("Test%i",i)),i);//Kaon
383
384 TH1D *efficiencyProtonTPCCB0 = GetHistoEfficiency(etacut,"hEfficiencyProtonTPCCB0",3,1,20,true,true,0,4);
385 for(int i=0;i<=4;i++) hadCorrectionPHOS->SetEfficiencyProtonTPC((TH1D*)efficiencyProtonTPCCB0->Clone(Form("Test%i",i)),i);
386 TH1D *efficiencyProtonTPCCB5 = GetHistoEfficiency(etacut,"hEfficiencyProtonTPCCB5",3,1,20,true,true,5,9);
387 for(int i=5;i<=9;i++) hadCorrectionPHOS->SetEfficiencyProtonTPC((TH1D*)efficiencyProtonTPCCB5->Clone(Form("Test%i",i)),i);
388 TH1D *efficiencyProtonTPCCB10 = GetHistoEfficiency(etacut,"hEfficiencyProtonTPCCB10",3,1,20,true,true,10,15);
389 for(int i=10;i<=19;i++) hadCorrectionPHOS->SetEfficiencyProtonTPC((TH1D*)efficiencyProtonTPCCB10->Clone(Form("Test%i",i)),i);//Proton
390
391 TH1D *efficiencyHadronTPCCB0 = GetHistoEfficiency(etacut,"hEfficiencyHadronTPCCB0",0,1,20,true,true,0,4);
392 for(int i=0;i<=4;i++) hadCorrectionPHOS->SetEfficiencyHadronTPC((TH1D*)efficiencyHadronTPCCB0->Clone(Form("Test%i",i)),i);
393 TH1D *efficiencyHadronTPCCB5 = GetHistoEfficiency(etacut,"hEfficiencyHadronTPCCB5",0,1,20,true,true,5,9);
394 for(int i=5;i<=9;i++) hadCorrectionPHOS->SetEfficiencyHadronTPC((TH1D*)efficiencyHadronTPCCB5->Clone(Form("Test%i",i)),i);
395 TH1D *efficiencyHadronTPCCB10 = GetHistoEfficiency(etacut,"hEfficiencyHadronTPCCB10",0,1,20,true,true,10,15);
396 for(int i=10;i<=19;i++) hadCorrectionPHOS->SetEfficiencyHadronTPC((TH1D*)efficiencyHadronTPCCB10->Clone(Form("Test%i",i)),i);//Hadron
397
398
399 TH1D *efficiencyPionITSCB0 = GetHistoEfficiency(etacut,"hEfficiencyPionITSCB0",1,1,20,false,true,0,4);
400 for(int i=0;i<=4;i++) hadCorrectionPHOS->SetEfficiencyPionITS((TH1D*)efficiencyPionITSCB0->Clone(Form("Test%i",i)),i);
401 TH1D *efficiencyPionITSCB5 = GetHistoEfficiency(etacut,"hEfficiencyPionITSCB5",1,1,20,false,true,5,9);
402 for(int i=5;i<=9;i++) hadCorrectionPHOS->SetEfficiencyPionITS((TH1D*)efficiencyPionITSCB5->Clone(Form("Test%i",i)),i);
403 TH1D *efficiencyPionITSCB10 = GetHistoEfficiency(etacut,"hEfficiencyPionITSCB10",1,1,20,false,true,10,15);
404 for(int i=10;i<=19;i++) hadCorrectionPHOS->SetEfficiencyPionITS((TH1D*)efficiencyPionITSCB10->Clone(Form("Test%i",i)),i);//Pion
405
406 TH1D *efficiencyKaonITSCB0 = GetHistoEfficiency(etacut,"hEfficiencyKaonITSCB0",2,1,20,false,true,0,4);
407 for(int i=0;i<=4;i++) hadCorrectionPHOS->SetEfficiencyKaonITS((TH1D*)efficiencyKaonITSCB0->Clone(Form("Test%i",i)),i);
408 TH1D *efficiencyKaonITSCB5 = GetHistoEfficiency(etacut,"hEfficiencyKaonITSCB5",2,1,20,false,true,5,9);
409 for(int i=5;i<=9;i++) hadCorrectionPHOS->SetEfficiencyKaonITS((TH1D*)efficiencyKaonITSCB5->Clone(Form("Test%i",i)),i);
410 TH1D *efficiencyKaonITSCB10 = GetHistoEfficiency(etacut,"hEfficiencyKaonITSCB10",2,1,20,false,true,10,15);
411 for(int i=10;i<=19;i++) hadCorrectionPHOS->SetEfficiencyKaonITS((TH1D*)efficiencyKaonITSCB10->Clone(Form("Test%i",i)),i);//Kaon
412
413 TH1D *efficiencyProtonITSCB0 = GetHistoEfficiency(etacut,"hEfficiencyProtonITSCB0",3,1,20,false,true,0,4);
414 for(int i=0;i<=4;i++) hadCorrectionPHOS->SetEfficiencyProtonITS((TH1D*)efficiencyProtonITSCB0->Clone(Form("Test%i",i)),i);
415 TH1D *efficiencyProtonITSCB5 = GetHistoEfficiency(etacut,"hEfficiencyProtonITSCB5",3,1,20,false,true,5,9);
416 for(int i=5;i<=9;i++) hadCorrectionPHOS->SetEfficiencyProtonITS((TH1D*)efficiencyProtonITSCB5->Clone(Form("Test%i",i)),i);
417 TH1D *efficiencyProtonITSCB10 = GetHistoEfficiency(etacut,"hEfficiencyProtonITSCB10",3,1,20,false,true,10,15);
418 for(int i=10;i<=19;i++) hadCorrectionPHOS->SetEfficiencyProtonITS((TH1D*)efficiencyProtonITSCB10->Clone(Form("Test%i",i)),i);//Proton
419
420 TH1D *efficiencyHadronITSCB0 = GetHistoEfficiency(etacut,"hEfficiencyHadronITSCB0",0,1,20,false,true,0,4);
421 for(int i=0;i<=4;i++) hadCorrectionPHOS->SetEfficiencyHadronITS((TH1D*)efficiencyHadronITSCB0->Clone(Form("Test%i",i)),i);
422 TH1D *efficiencyHadronITSCB5 = GetHistoEfficiency(etacut,"hEfficiencyHadronITSCB5",0,1,20,false,true,5,9);
423 for(int i=5;i<=9;i++) hadCorrectionPHOS->SetEfficiencyHadronITS((TH1D*)efficiencyHadronITSCB5->Clone(Form("Test%i",i)),i);
424 TH1D *efficiencyHadronITSCB10 = GetHistoEfficiency(etacut,"hEfficiencyHadronITSCB10",0,1,20,false,true,10,15);
425 for(int i=10;i<=19;i++) hadCorrectionPHOS->SetEfficiencyHadronITS((TH1D*)efficiencyHadronITSCB10->Clone(Form("Test%i",i)),i);//Hadron
426 }//EMCAL
427
428 TH1D *backgroundTPC = GetHistoCorrBkgd(etacut,"hBackgroundTPC",true,ispp,forSim);
429 TH1D *backgroundITS = GetHistoCorrBkgd(etacut,"hBackgroundITS",false,ispp,forSim);
a02dfa56 430 hadCorrectionPHOS->SetBackgroundCorrectionTPC(backgroundTPC);
431 hadCorrectionPHOS->SetBackgroundCorrectionITS(backgroundITS);
48a07264 432 CorrBkgdPlots(prodname,shortprodname,true,ispp,forSim);
433 CorrBkgdPlots(prodname,shortprodname,false,ispp,forSim);
a02dfa56 434
48a07264 435 hadCorrectionPHOS->Report();
a02dfa56 436 //Write the output
437 outfile->cd();
438 hadCorrectionPHOS->Write();
439 outfile->Write();
f427cbed 440 outfile->Close();
441
f43fc416 442
f427cbed 443 timer.Stop();
444 timer.Print();
445}
446
b64de20c 447//==================================CorrNeutral==============================================
48a07264 448Float_t CorrNeutral(float ptcut, char *prodname, char *shortprodname, bool ispp, bool forSim, bool TPC, bool hadronic, float etacut){
449 if(!forSim){//for data we have evaluated the neutral correction from ALICE data
450 if(hadronic){//for tot et from had et
451 if(ispp){
452 return 1.0/0.571;
453 }
454 else{
455 return 1.0/0.549;
456 }
457 }
458 else{//for had et only
459 if(ispp){
460 return 1.0/0.736;
461 }
462 else{
463 return 1.0/0.689;
464 }
465 }
466 }
f427cbed 467 gStyle->SetOptTitle(0);
468 gStyle->SetOptStat(0);
469 gStyle->SetOptFit(0);
470 TCanvas *c = new TCanvas("c","c",400,400);
471 c->SetTopMargin(0.0);
472 c->SetRightMargin(0.0);
473 c->SetBorderSize(0);
474 c->SetFillColor(0);
475 c->SetFillColor(0);
476 c->SetBorderMode(0);
477 c->SetFrameFillColor(0);
478 c->SetFrameBorderMode(0);
479
480 TPad *ptpad = c->cd(1);
481 ptpad->SetTopMargin(0.04);
482 ptpad->SetRightMargin(0.04);
483 ptpad->SetLeftMargin(0.149288);
484 ptpad->SetBorderSize(0);
485 ptpad->SetFillColor(0);
486 ptpad->SetFillColor(0);
487 ptpad->SetBorderMode(0);
488 ptpad->SetFrameFillColor(0);
489 ptpad->SetFrameBorderMode(0);
490
491 int phosmarker = 20;
492
f427cbed 493 sprintf(prefix,"%s%2.1f",shortprodname,ptcut);
494
f427cbed 495 sprintf(histoname,"%stotal",histoname);
496 int colortotal = 1;
b64de20c 497 int casetotal = 4;
498 if(hadronic) casetotal = 8;
48a07264 499 TH1D *total = GetHistoCorrNeutral(ptcut,histoname,ispp,forSim,casetotal,false,colortotal,phosmarker,hadronic);
f427cbed 500
501 int colorallneutral = 2;
48a07264 502 TH1D *allneutral = GetHistoCorrNeutral(ptcut,"allneutral",ispp,forSim,3,false,colorallneutral,phosmarker,hadronic);
f427cbed 503
504 int colorchargedsecondary = TColor::kViolet-3;
48a07264 505 TH1D *chargedsecondary = GetHistoCorrNeutral(ptcut,"chargedsecondary",ispp,forSim,2,false,colorchargedsecondary,phosmarker,hadronic);
f427cbed 506
507 int colorneutralUndet = 4;
48a07264 508 TH1D *neutralUndet = GetHistoCorrNeutral(ptcut,"neutralUndet",ispp,forSim,1,false,colorneutralUndet,phosmarker,hadronic);
f427cbed 509
510 int colorv0 = TColor::kGreen+2;
48a07264 511 TH1D *v0 = GetHistoCorrNeutral(ptcut,"v0",ispp,forSim,0,false,colorv0,phosmarker,hadronic);
b64de20c 512
513 int colorem = TColor::kCyan;
48a07264 514 TH1D *em = GetHistoCorrNeutral(ptcut,"em",ispp,forSim,9,false,colorem,phosmarker,hadronic);
f427cbed 515
516 TF1 *func = new TF1("func","[0]",-.7,.7);
517 func->SetParameter(0,0.2);
0e866ddc 518 total->Fit(func,"","",-etacut,etacut);
f427cbed 519
520 //total->SetAxisRange(0.0,4);
521 total->GetXaxis()->SetLabelSize(0.05);
522 total->GetYaxis()->SetLabelSize(0.045);
523 total->GetXaxis()->SetTitleSize(0.05);
524 total->GetYaxis()->SetTitleSize(0.06);
b64de20c 525 if(hadronic){
526 total->SetMaximum(0.6);
527 total->SetMinimum(0.0);
528 }
529 else{
530 total->SetMaximum(0.3);
531 total->SetMinimum(0.0);
532 }
f427cbed 533 total->Draw();
534 allneutral->Draw("same");
535 chargedsecondary->Draw("same");
536 neutralUndet->Draw("same");
537 v0->Draw("same");
b64de20c 538 if(hadronic) em->Draw("same");
f427cbed 539
540 TLatex *tex = new TLatex(0.161478,1.0835,prodname);
541 tex->SetTextSize(0.0537634);
542 tex->Draw();
543 TLegend *leg2 = new TLegend(0.518321,0.746873,0.774812,0.955343);
544 leg2->AddEntry(total,"Total");
545 leg2->AddEntry(allneutral,"#Lambda,#bar{#Lambda},K^{0}_{S},K^{0}_{L},n,#bar{n}");
546 leg2->AddEntry(neutralUndet,"K^{0}_{L},n,#bar{n}");
b64de20c 547 if(hadronic) leg2->AddEntry(em,"e^{#pm},#gamma,#eta,#pi^{0},#omega");
f427cbed 548 leg2->SetFillStyle(0);
549 leg2->SetFillColor(0);
550 leg2->SetBorderSize(0);
551 leg2->SetTextSize(0.0548607);
552 leg2->Draw();
b64de20c 553 if(hadronic){
554 sprintf(epsname,"pics/%s/fhadronic.eps",shortprodname);
555 sprintf(pngname,"pics/%s/fhadronic.png",shortprodname);
556 }
557 else{
558 sprintf(epsname,"pics/%s/fneutral.eps",shortprodname);
559 sprintf(pngname,"pics/%s/fneutral.png",shortprodname);
560 }
f427cbed 561 c->SaveAs(epsname);
562 c->SaveAs(pngname);
563
a02dfa56 564 delete total;
565 delete allneutral;
566 delete chargedsecondary;
567 delete neutralUndet;
568 delete v0;
569 delete em;
570 delete c;
f427cbed 571 float corr = func->GetParameter(0);
a02dfa56 572 delete func;
48a07264 573 delete tex;
574 delete leg2;
f427cbed 575 return 1.0/(1.0-corr);
576
577}
48a07264 578TH1D *GetHistoCorrNeutral(float cut, char *name, bool ispp, bool forSim, int mycase, bool eta, int color, int marker, bool hadronic){
579 file->cd();
f43fc416 580 char *reweightname = "";
581 if(!forSim) reweightname = "Reweighted";
f427cbed 582 TH2F *numeratorParent;
583 switch(mycase){
584 case 0:
f43fc416 585 numeratorParent= (TH2F*)((TH2F*) out2->FindObject(Form("EtSimulatedLambda%s",reweightname)))->Clone("v0");
586 numeratorParent->Add((TH2F*) out2->FindObject(Form("EtSimulatedAntiLambda%s",reweightname)));
587 numeratorParent->Add((TH2F*) out2->FindObject(Form("EtSimulatedK0S%s",reweightname)));
f427cbed 588 break;
589 case 1:
f43fc416 590 numeratorParent = (TH2F*)((TH2F*) out2->FindObject(Form("EtSimulatedK0L%s",reweightname)))->Clone("Knnbar");
f427cbed 591 numeratorParent->Add((TH2F*) out2->FindObject("EtSimulatedNeutron"));
592 numeratorParent->Add((TH2F*) out2->FindObject("EtSimulatedAntiNeutron"));
593 break;
594 case 2:
595 numeratorParent = (TH2F*)((TH2F*) out2->FindObject("EtSimulatedOmega"))->Clone("ch2ndary");
596 numeratorParent->Add((TH2F*) out2->FindObject("EtSimulatedAntiOmega"));
597 numeratorParent->Add((TH2F*) out2->FindObject("EtSimulatedXi"));
598 numeratorParent->Add((TH2F*) out2->FindObject("EtSimulatedAntiXi"));
599 numeratorParent->Add((TH2F*) out2->FindObject("EtSimulatedXi0"));
600 numeratorParent->Add((TH2F*) out2->FindObject("EtSimulatedAntiXi0"));
601 break;
602 case 3:
f43fc416 603 numeratorParent= (TH2F*)((TH2F*) out2->FindObject(Form("EtSimulatedLambda%s",reweightname)))->Clone("allneutral");
604 numeratorParent->Add((TH2F*) out2->FindObject(Form("EtSimulatedAntiLambda%s",reweightname)));
605 numeratorParent->Add((TH2F*) out2->FindObject(Form("EtSimulatedK0S%s",reweightname)));
606 numeratorParent->Add((TH2F*) out2->FindObject(Form("EtSimulatedK0L%s",reweightname)));
f427cbed 607 numeratorParent->Add((TH2F*) out2->FindObject("EtSimulatedNeutron"));
608 numeratorParent->Add((TH2F*) out2->FindObject("EtSimulatedAntiNeutron"));
609 break;
610 case 4:
f43fc416 611 numeratorParent= (TH2F*)((TH2F*) out2->FindObject(Form("EtSimulatedLambda%s",reweightname)))->Clone("allneutral");
612 numeratorParent->Add((TH2F*) out2->FindObject(Form("EtSimulatedAntiLambda%s",reweightname)));
613 numeratorParent->Add((TH2F*) out2->FindObject(Form("EtSimulatedK0S%s",reweightname)));
614 numeratorParent->Add((TH2F*) out2->FindObject(Form("EtSimulatedK0L%s",reweightname)));
f427cbed 615 numeratorParent->Add((TH2F*) out2->FindObject("EtSimulatedNeutron"));
616 numeratorParent->Add((TH2F*) out2->FindObject("EtSimulatedAntiNeutron"));
617 numeratorParent->Add((TH2F*) out2->FindObject("EtSimulatedOmega"));
618 numeratorParent->Add((TH2F*) out2->FindObject("EtSimulatedAntiOmega"));
619 numeratorParent->Add((TH2F*) out2->FindObject("EtSimulatedXi"));
620 numeratorParent->Add((TH2F*) out2->FindObject("EtSimulatedAntiXi"));
621 numeratorParent->Add((TH2F*) out2->FindObject("EtSimulatedXi0"));
622 numeratorParent->Add((TH2F*) out2->FindObject("EtSimulatedAntiXi0"));
623 break;
624 case 5:
625 numeratorParent= (TH2F*)((TH2F*) out2->FindObject("EtSimulatedXi"))->Clone("allxi");
626 numeratorParent->Add((TH2F*) out2->FindObject("EtSimulatedAntiXi"));
627 numeratorParent->Add((TH2F*) out2->FindObject("EtSimulatedXi0"));
628 numeratorParent->Add((TH2F*) out2->FindObject("EtSimulatedAntiXi0"));
629 break;
630 case 6:
631 numeratorParent = (TH2F*)((TH2F*) out2->FindObject("EtSimulatedOmega"))->Clone("allomega");
632 numeratorParent->Add((TH2F*) out2->FindObject("EtSimulatedAntiOmega"));
633 break;
634 case 7:
635 numeratorParent = (TH2F*)((TH2F*) out2->FindObject("EtSimulatedSigma"))->Clone("allsigma");
636 numeratorParent->Add((TH2F*) out2->FindObject("EtSimulatedAntiSigma"));
637 break;
b64de20c 638 case 8:
f43fc416 639 numeratorParent= (TH2F*)((TH2F*) out2->FindObject(Form("EtSimulatedLambda%s",reweightname)))->Clone("allneutral");
640 numeratorParent->Add((TH2F*) out2->FindObject(Form("EtSimulatedAntiLambda%s",reweightname)));
641 numeratorParent->Add((TH2F*) out2->FindObject(Form("EtSimulatedK0S%s",reweightname)));
642 numeratorParent->Add((TH2F*) out2->FindObject(Form("EtSimulatedK0L%s",reweightname)));
b64de20c 643 numeratorParent->Add((TH2F*) out2->FindObject("EtSimulatedNeutron"));
644 numeratorParent->Add((TH2F*) out2->FindObject("EtSimulatedAntiNeutron"));
645 numeratorParent->Add((TH2F*) out2->FindObject("EtSimulatedOmega"));
646 numeratorParent->Add((TH2F*) out2->FindObject("EtSimulatedAntiOmega"));
647 numeratorParent->Add((TH2F*) out2->FindObject("EtSimulatedXi"));
648 numeratorParent->Add((TH2F*) out2->FindObject("EtSimulatedAntiXi"));
649 numeratorParent->Add((TH2F*) out2->FindObject("EtSimulatedXi0"));
650 numeratorParent->Add((TH2F*) out2->FindObject("EtSimulatedAntiXi0"));
651 numeratorParent->Add((TH2F*) out2->FindObject("EtSimulatedGamma"));
652 numeratorParent->Add((TH2F*) out2->FindObject("EtSimulatedEta"));
653 numeratorParent->Add((TH2F*) out2->FindObject("EtSimulatedPi0"));
654 numeratorParent->Add((TH2F*) out2->FindObject("EtSimulatedOmega0"));
655 numeratorParent->Add((TH2F*) out2->FindObject("EtSimulatedEPlus"));
656 numeratorParent->Add((TH2F*) out2->FindObject("EtSimulatedEMinus"));
657 break;
658 case 9:
659 numeratorParent= (TH2F*)((TH2F*) out2->FindObject("EtSimulatedGamma"))->Clone("allem");
660 numeratorParent->Add((TH2F*) out2->FindObject("EtSimulatedEta"));
661 numeratorParent->Add((TH2F*) out2->FindObject("EtSimulatedPi0"));
662 numeratorParent->Add((TH2F*) out2->FindObject("EtSimulatedOmega0"));
663 numeratorParent->Add((TH2F*) out2->FindObject("EtSimulatedEPlus"));
664 numeratorParent->Add((TH2F*) out2->FindObject("EtSimulatedEMinus"));
665 break;
f427cbed 666 }
667
668 TH2F *allhad;
f43fc416 669 //allhad=(TH2F*) ((TH2F*) out2->FindObject("EtSimulatedAllHadron"))->Clone("id");
670 allhad=(TH2F*) ((TH2F*) out2->FindObject("EtSimulatedPiPlus"))->Clone("id");
671 allhad->Add((TH2F*) out2->FindObject("EtSimulatedPiMinus"));
672 allhad->Add((TH2F*) out2->FindObject("EtSimulatedKMinus"));
673 allhad->Add((TH2F*) out2->FindObject("EtSimulatedKPlus"));
674 allhad->Add((TH2F*) out2->FindObject("EtSimulatedProton"));
675 allhad->Add((TH2F*) out2->FindObject("EtSimulatedAntiProton"));
676 allhad->Add((TH2F*) out2->FindObject(Form("EtSimulatedLambda%s",reweightname)));
677 allhad->Add((TH2F*) out2->FindObject(Form("EtSimulatedAntiLambda%s",reweightname)));
678 allhad->Add((TH2F*) out2->FindObject(Form("EtSimulatedK0S%s",reweightname)));
679 allhad->Add((TH2F*) out2->FindObject(Form("EtSimulatedK0L%s",reweightname)));
680 allhad->Add((TH2F*) out2->FindObject("EtSimulatedNeutron"));
681 allhad->Add((TH2F*) out2->FindObject("EtSimulatedAntiNeutron"));
682 allhad->Add((TH2F*) out2->FindObject("EtSimulatedEPlus"));
683 allhad->Add((TH2F*) out2->FindObject("EtSimulatedEMinus"));
684 allhad->Add((TH2F*) out2->FindObject("EtSimulatedOmega"));
685 allhad->Add((TH2F*) out2->FindObject("EtSimulatedAntiOmega"));
686 allhad->Add((TH2F*) out2->FindObject("EtSimulatedXi"));
687 allhad->Add((TH2F*) out2->FindObject("EtSimulatedAntiXi"));
688 allhad->Add((TH2F*) out2->FindObject("EtSimulatedSigma"));
689 allhad->Add((TH2F*) out2->FindObject("EtSimulatedAntiSigma"));
690 allhad->Add((TH2F*) out2->FindObject("EtSimulatedXi0"));
691 allhad->Add((TH2F*) out2->FindObject("EtSimulatedAntiXi0"));
692
b64de20c 693 if(hadronic){//if we are getting the correction for the hadronic only case...
694 allhad->Add((TH2F*) out2->FindObject("EtSimulatedGamma"));
695 allhad->Add((TH2F*) out2->FindObject("EtSimulatedEta"));
696 allhad->Add((TH2F*) out2->FindObject("EtSimulatedPi0"));
697 allhad->Add((TH2F*) out2->FindObject("EtSimulatedOmega0"));
698 allhad->Add((TH2F*) out2->FindObject("EtSimulatedEPlus"));
699 allhad->Add((TH2F*) out2->FindObject("EtSimulatedEMinus"));
700 }
f427cbed 701
f427cbed 702 TH1D *denominator;
703 TH1D *numerator;
704 if(eta){
705 int lowbin = numeratorParent->GetYaxis()->FindBin(-cut+.001);//make sure we don't accv0entally get the wrong bin
706 int highbin = numeratorParent->GetYaxis()->FindBin(cut-.001);
a02dfa56 707 //cout<<"Projecting from "<<numeratorParent->GetYaxis()->GetBinLowEdge(lowbin)<<" to "<<numeratorParent->GetYaxis()->GetBinLowEdge(highbin+1)<<endl;
f427cbed 708 denominator = allhad->ProjectionX("name",lowbin,highbin);
709 numerator = numeratorParent->ProjectionX("numerator",lowbin,highbin);
710 }
711 else{
712 int lowbin = allhad->GetXaxis()->FindBin(cut);//make sure we don't accidentally get the wrong bin
713 int highbin = allhad->GetXaxis()->GetNbins();
a02dfa56 714 //cout<<"Projecting from "<<allhad->GetXaxis()->GetBinLowEdge(lowbin)<<" to "<<allhad->GetXaxis()->GetBinLowEdge(highbin+1)<<endl;
f427cbed 715 numerator = numeratorParent->ProjectionY("name",lowbin,highbin);
716 denominator = allhad->ProjectionY("denominator",lowbin,highbin);
717 }
718 numerator->Divide(denominator);
0e866ddc 719 //numerator->Rebin(2);
720 //numerator->Scale(0.5);
f427cbed 721 numerator->SetYTitle("E_{T}^{had,sample}/E_{T}^{had,total}");
722 numerator->GetYaxis()->SetTitleOffset(1.2);
723 numerator->SetMarkerColor(color);
724 numerator->SetLineColor(color);
725 numerator->SetMarkerStyle(marker);
a02dfa56 726 delete denominator;
727 delete numeratorParent;
728 delete allhad;
729 //file->Close();
48a07264 730 numerator->SetName(name);
f427cbed 731 return numerator;
732
733}
734
b64de20c 735//===============================CorrPtCut=========================================
48a07264 736TH1D *GetHistoCorrPtCut(float ptcut, char *name, bool ispp, bool forSim, int mycase){
737 file->cd();
f427cbed 738 TH2F *allhad = ((TH2F*) out2->FindObject("EtSimulatedAllHadron"))->Clone("allhad");
48a07264 739 TH2F *ptlow = ((TH2F*) out2->FindObject("EtSimulatedChargedHadronAssumingNoPt"))->Clone("ptlow");
740 TH2F *pthigh;
741 if(ptcut>0.14){//TPC cut off
742 (TH2F*)pthigh =(TH2F*) ((TH2F*) out2->FindObject("EtSimulatedChargedHadronAssumingPtTPCCut"))->Clone("pthigh");
743 }
744 else{
745 (TH2F*)pthigh =(TH2F*) ((TH2F*) out2->FindObject("EtSimulatedChargedHadronAssumingPtITSCut"))->Clone("pthigh");
746 }
f427cbed 747
748 int lowbin = allhad->GetXaxis()->FindBin(0.0);//make sure we don't accidentally get the wrong bin
749 int highbin = allhad->GetXaxis()->FindBin(ptcut);
750 int nbins = allhad->GetXaxis()->GetNbins();
a02dfa56 751 //cout<<"Projecting from "<<allhad->GetXaxis()->GetBinLowEdge(lowbin)<<" to "<<allhad->GetXaxis()->GetBinLowEdge(highbin+1)<<endl;
752 //cout<<"Projecting from "<<allhad->GetXaxis()->GetBinLowEdge(lowbin)<<" to "<<allhad->GetXaxis()->GetBinLowEdge(nbins)<<endl;
f427cbed 753
a02dfa56 754 //allhad->Sumw2();
48a07264 755 TH1D *numerator;
756 TH1D *denominator;
757 switch(mycase){
758 case -1:
759 numerator = ptlow->ProjectionY("nameLow",lowbin,highbin);
760 denominator = allhad->ProjectionY("denominatorLow",highbin,nbins);
761 denominator->Add(ptlow);
762 break;
763 case 1:
764 numerator = pthigh->ProjectionY("nameHigh",lowbin,highbin);
765 denominator = allhad->ProjectionY("denominatorHigh",highbin,nbins);
766 denominator->Add(pthigh);
767 break;
768 default:
769 numerator = allhad->ProjectionY("name",lowbin,highbin);
770 denominator = allhad->ProjectionY("denominator",lowbin,nbins);
771 }
f427cbed 772 numerator->Divide(denominator);
773 numerator->SetYTitle("E_{T}^{had, p_{T}<cut-off}/E_{T}^{had, all p_{T}}");
774 numerator->GetYaxis()->SetTitleOffset(1.);
775 numerator->GetYaxis()->SetTitleSize(0.08);
776 numerator->GetYaxis()->SetLabelSize(0.05);
777 numerator->GetXaxis()->SetTitleSize(0.08);
778 numerator->GetXaxis()->SetLabelSize(0.05);
779 numerator->GetXaxis()->SetTitleOffset(.6);
780 //numerator->Rebin(2);
781 //numerator->Scale(0.5);
782 //numerator->Draw("e");
a02dfa56 783 delete allhad;
784 delete denominator;
48a07264 785 delete ptlow;
786 delete pthigh;
787 numerator->SetName(name);
f427cbed 788 return numerator;
789
790}
791
48a07264 792Float_t CorrPtCut(float ptcut, char *prodname, char *shortprodname, bool ispp, bool forSim, int mycase){
f427cbed 793
794 gStyle->SetOptTitle(0);
795 gStyle->SetOptStat(0);
796 gStyle->SetOptFit(0);
797 TCanvas *c = new TCanvas("c","c",500,400);
798 c->SetTopMargin(0.04);
799 c->SetRightMargin(0.04);
800 c->SetLeftMargin(0.181452);
801 c->SetBottomMargin(0.134409);
802 c->SetBorderSize(0);
803 c->SetFillColor(0);
804 c->SetFillColor(0);
805 c->SetBorderMode(0);
806 c->SetFrameFillColor(0);
807 c->SetFrameBorderMode(0);
808
809
810
48a07264 811 TH1D *High = GetHistoCorrPtCut(0.15-.001,"High",ispp,forSim,mycase);
812 TH1D *Low = GetHistoCorrPtCut(0.1-.001,"Low",ispp,forSim,mycase);
813 TH1D *Lowest = GetHistoCorrPtCut(0.05-.001,"Lowest",ispp,forSim,mycase);
f427cbed 814
815 TF1 *func = new TF1("func","[0]",-.7,.7);
816 func->SetParameter(0,0.2);
817 if(ptcut<.125){//its cuts
818 Low->Fit(func);
819 }
820 else{//tpc cuts
821 High->Fit(func);
822 }
823
824 High->SetMaximum(0.04);
825 High->SetMinimum(0.0);
826 High->SetMarkerColor(2);
827 Low->SetMarkerColor(4);
828 High->SetLineColor(2);
829 Low->SetLineColor(4);
830 High->SetMinimum(0.0);
831 High->SetMarkerStyle(20);
832 Low->SetMarkerStyle(21);
833 Lowest->SetMarkerStyle(22);
834 High->Draw();
835 Low->Draw("same");
836 Lowest->Draw("same");
837
838
839
840 TLatex *tex = new TLatex(-0.723444,0.0373593,prodname);
841 tex->SetTextSize(0.0537634);
842 tex->Draw();
843 TLegend *leg = new TLegend(0.217742,0.66129,0.477823,0.873656);
844 leg->AddEntry(High,"p_{T} cut-off = 0.15 GeV/c");
845 leg->AddEntry(Low,"p_{T} cut-off = 0.1 GeV/c");
846 leg->AddEntry(Lowest,"p_{T} cut-off = 0.05 GeV/c");
847 leg->SetFillStyle(0);
848 leg->SetFillColor(0);
849 leg->SetBorderSize(0);
850 leg->SetTextSize(0.0537634);
851 leg->Draw();
852
853 if(ptcut<.125){//its cuts
854 c->SaveAs(Form("pics/%s/fptcutITS.eps",shortprodname));
855 c->SaveAs(Form("pics/%s/fptcutITS.png",shortprodname));
856 }
857 else{//tpc cuts
858 c->SaveAs(Form("pics/%s/fptcutTPC.eps",shortprodname));
859 c->SaveAs(Form("pics/%s/fptcutTPC.png",shortprodname));
860 }
861
862 float corr = func->GetParameter(0);
a02dfa56 863 //cout<<"Pt cut correction: "<<1.0/(1.0-corr)<<endl;
864 delete High;
865 delete Low;
866 delete Lowest;
867 delete func;
868 delete c;
48a07264 869 delete tex;
870 delete leg;
f427cbed 871 return 1.0/(1.0-corr);
872}
b64de20c 873
874
875
876//==================================CorrNotID=================================================
48a07264 877TH1D *GetHistoCorrNotID(float etacut,char *name, bool TPC, bool eta, bool ispp, bool forSim){
878 file->cd();
879 char *myname = mynameITS;
880 if(TPC) myname = mynameTPC;
b64de20c 881 TH2F *notid = ((TH2F*) out2->FindObject(Form("EtReconstructed%sUnidentifiedAssumingPion",myname)))->Clone("notid");
882 TH2F *nNotid = ((TH2F*) out2->FindObject(Form("EtNReconstructed%sUnidentified",myname)))->Clone("nNotid");
0e866ddc 883 if(!eta){
48a07264 884 //cout<<"Correction determined for all charged hadrons"<<endl;
0e866ddc 885 notid->Add((TH2F*) out2->FindObject(Form("EtReconstructed%sIdentifiedPiPlus",myname)));
886 notid->Add((TH2F*) out2->FindObject(Form("EtReconstructed%sIdentifiedPiMinus",myname)));
887 notid->Add((TH2F*) out2->FindObject(Form("EtReconstructed%sIdentifiedKPlus",myname)));
888 notid->Add((TH2F*) out2->FindObject(Form("EtReconstructed%sIdentifiedKMinus",myname)));
889 notid->Add((TH2F*) out2->FindObject(Form("EtReconstructed%sIdentifiedProton",myname)));
890 notid->Add((TH2F*) out2->FindObject(Form("EtReconstructed%sIdentifiedAntiProton",myname)));
891 nNotid->Add((TH2F*) out2->FindObject(Form("EtNReconstructed%sPiPlus",myname)));
892 nNotid->Add((TH2F*) out2->FindObject(Form("EtNReconstructed%sPiMinus",myname)));
893 nNotid->Add((TH2F*) out2->FindObject(Form("EtNReconstructed%sKPlus",myname)));
894 nNotid->Add((TH2F*) out2->FindObject(Form("EtNReconstructed%sKMinus",myname)));
895 nNotid->Add((TH2F*) out2->FindObject(Form("EtNReconstructed%sProton",myname)));
896 nNotid->Add((TH2F*) out2->FindObject(Form("EtNReconstructed%sAntiProton",myname)));
897 }
b64de20c 898
b64de20c 899 TH2F *id = ((TH2F*) out2->FindObject(Form("EtReconstructed%sUnidentified",myname)))->Clone("id");
0e866ddc 900 if(!eta){
901 id->Add((TH2F*) out2->FindObject(Form("EtReconstructed%sIdentifiedPiPlus",myname)));
902 id->Add((TH2F*) out2->FindObject(Form("EtReconstructed%sIdentifiedPiMinus",myname)));
903 id->Add((TH2F*) out2->FindObject(Form("EtReconstructed%sIdentifiedKPlus",myname)));
904 id->Add((TH2F*) out2->FindObject(Form("EtReconstructed%sIdentifiedKMinus",myname)));
905 id->Add((TH2F*) out2->FindObject(Form("EtReconstructed%sIdentifiedProton",myname)));
906 id->Add((TH2F*) out2->FindObject(Form("EtReconstructed%sIdentifiedAntiProton",myname)));
907 }
b64de20c 908
0e866ddc 909 TH1D *nNotidProj;
910 TH1D *denominator;
911 TH1D *numerator;
912 if(eta){
913 int lowbin = notid->GetYaxis()->FindBin(-etacut+.001);//make sure we don't accv0entally get the wrong bin
914 int highbin = notid->GetYaxis()->FindBin(etacut-.001);
48a07264 915 //cout<<"Projecting from "<<notid->GetYaxis()->GetBinLowEdge(lowbin)<<" to "<<notid->GetYaxis()->GetBinLowEdge(highbin+1)<<endl;
0e866ddc 916 denominator = id->ProjectionX("name",lowbin,highbin);
917 numerator = notid->ProjectionX("numerator",lowbin,highbin);
918 nNotidProj = nNotid->ProjectionX("nNotidProj",lowbin,highbin);
919 }
920 else{
48a07264 921 //cout<<"Getting eta dependence"<<endl;
0e866ddc 922 int lowbin = id->GetXaxis()->FindBin(etacut);//make sure we don't accidentally get the wrong bin
923 int highbin;
924 if(etacut<0.15){//then we actually have ITS standalone tracks and we only want this to run from 0.1 to 0.15 because this will be used only for ITS standalone tracks
925 highbin = id->GetXaxis()->FindBin(0.15);
926 }
927 else{
928 highbin = id->GetXaxis()->GetNbins();
929 }
48a07264 930 //cout<<"Projecting from "<<id->GetXaxis()->GetBinLowEdge(lowbin)<<" to "<<id->GetXaxis()->GetBinLowEdge(highbin+1)<<endl;
0e866ddc 931 numerator = notid->ProjectionY("name",lowbin,highbin);
932 denominator = id->ProjectionY("denominator",lowbin,highbin);
933 nNotidProj = nNotid->ProjectionY("nNotidProj",lowbin,highbin);
934 }
b64de20c 935 TH1D *result = numerator;
936 if(!denominator){
937 cerr<<"Uh-oh! Can't find denominator!!";
938 return numerator;
939 }
940 else{result->Divide(denominator);}
941 if(result->GetNbinsX() != nNotidProj->GetNbinsX()){
942 cerr<<"Uh-oh! Can't rescale errors! "<<result->GetNbinsX()<<"!="<<nNotidProj->GetNbinsX()<<endl;
943 return result;
944 }
b64de20c 945 if(!nNotidProj){
946 cerr<<"Uh-oh! Can't find histogram!!";
947 return numerator;
948 }
0e866ddc 949 //fixing the errors
950 for(int i=1;i<=result->GetNbinsX();i++){
951 Float_t value = result->GetBinContent(i);
952 Float_t valueerr = result->GetBinError(i);
953 Float_t n = nNotidProj->GetBinContent(i);
954 Float_t err;
955 if(n<=0) err = 0.0;
956 else{err= value/TMath::Power(n,0.5);}
957 result->SetBinError(i,err);
958 //cout<<"Was "<<valueerr<<", setting to "<<err<<endl;
b64de20c 959 }
960 result->SetYTitle("Ratio of E_{T}^{assuming pion}/E_{T}^{real}");
961 result->GetYaxis()->SetTitleOffset(1.2);
a02dfa56 962 delete denominator;
963 delete nNotidProj;
48a07264 964 delete notid;
965 delete nNotid;
966 delete id;
967 result->SetName(name);
b64de20c 968 return result;
969
970}
971
48a07264 972TH1D *CorrNotID(float etacut,char *name, char *prodname, char *shortprodname, bool TPC, bool ispp, bool forSim){
b64de20c 973 gStyle->SetOptTitle(0);
974 gStyle->SetOptStat(0);
975 gStyle->SetOptFit(0);
976 TCanvas *c = new TCanvas("c","c",500,400);
977 c->SetTopMargin(0.04);
978 c->SetRightMargin(0.04);
979 c->SetBorderSize(0);
980 c->SetFillColor(0);
981 c->SetFillColor(0);
982 c->SetBorderMode(0);
983 c->SetFrameFillColor(0);
984 c->SetFrameBorderMode(0);
985
48a07264 986 TH1D *PHOS = GetHistoCorrNotID(etacut,name,TPC,true,ispp,forSim);
b64de20c 987 PHOS->SetMarkerColor(2);
988 PHOS->SetLineColor(2);
b64de20c 989 PHOS->SetAxisRange(0.0,4);
990 if(TPC){
991 PHOS->SetMaximum(1.1);
992 PHOS->SetMinimum(0.85);
993 }
994 else{
b64de20c 995 PHOS->SetAxisRange(0.0,0.5);
996 }
a02dfa56 997 PHOS->SetMarkerStyle(20);
b64de20c 998 PHOS->Draw();
999 TLatex *tex = new TLatex(0.161478,1.0835,prodname);
1000 tex->SetTextSize(0.0537634);
1001 tex->Draw();
48a07264 1002 char *detector = detectorEMCAL;
1003 if(etacut<0.2) detector = detectorPHOS;
b64de20c 1004 if(TPC){
1005 sprintf(epsname,"pics/%s/fnotidTPC%s.eps",shortprodname,detector);
1006 sprintf(pngname,"pics/%s/fnotidTPC%s.png",shortprodname,detector);
1007 }
1008 else{
1009 sprintf(epsname,"pics/%s/fnotidITS%s.eps",shortprodname,detector);
1010 sprintf(pngname,"pics/%s/fnotidITS%s.png",shortprodname,detector);
1011 }
1012
1013 c->SaveAs(epsname);
1014 c->SaveAs(pngname);
a02dfa56 1015 delete c;
48a07264 1016 delete tex;
1017 PHOS->SetName(name);
b64de20c 1018 return PHOS;
1019}
1020
48a07264 1021Float_t CorrNotIDConst(float ptcut, float etacut,char *name, char *prodname, char *shortprodname, bool TPC, bool ispp, bool forSim){
1022 if(!forSim){
1023 if(ispp){
1024 return 0.996;
1025 }
1026 else{
1027 return 0.976;
1028 }
1029 }
0e866ddc 1030 gStyle->SetOptTitle(0);
1031 gStyle->SetOptStat(0);
1032 gStyle->SetOptFit(0);
1033 TCanvas *c = new TCanvas("c","c",500,400);
1034 c->SetTopMargin(0.04);
1035 c->SetRightMargin(0.04);
1036 c->SetBorderSize(0);
1037 c->SetFillColor(0);
1038 c->SetFillColor(0);
1039 c->SetBorderMode(0);
1040 c->SetFrameFillColor(0);
1041 c->SetFrameBorderMode(0);
1042
48a07264 1043 TH1D *PHOS = GetHistoCorrNotID(ptcut,name,TPC,false,ispp,forSim);
0e866ddc 1044 PHOS->SetMarkerColor(2);
1045 PHOS->SetLineColor(2);
1046 PHOS->SetMaximum(1.01);
1047 PHOS->SetMinimum(0.98);
1048 TF1 *func = new TF1("func","[0]",-etacut,etacut);
1049 PHOS->Fit(func,"","",-etacut,etacut);
1050 PHOS->SetMarkerStyle(20);
1051 PHOS->Draw();
1052 TLatex *tex = new TLatex(0.161478,1.0835,prodname);
1053 tex->SetTextSize(0.0537634);
1054 tex->Draw();
48a07264 1055 char *detector = detectorEMCAL;
1056 if(etacut<0.2) detector = detectorPHOS;
0e866ddc 1057 if(TPC){
1058 sprintf(epsname,"pics/%s/fnotidConstTPC%s.eps",shortprodname,detector);
1059 sprintf(pngname,"pics/%s/fnotidConstTPC%s.png",shortprodname,detector);
1060 }
1061 else{
1062 sprintf(epsname,"pics/%s/fnotidConstITS%s.eps",shortprodname,detector);
1063 sprintf(pngname,"pics/%s/fnotidConstITS%s.png",shortprodname,detector);
1064 }
1065
1066 c->SaveAs(epsname);
1067 c->SaveAs(pngname);
1068 delete c;
48a07264 1069 delete PHOS;
1070 float value = func->GetParameter(0);
1071 delete func;
1072 delete tex;
1073 return value;
0e866ddc 1074}
1075
b64de20c 1076//==================================CorrNoID=================================================
48a07264 1077TH1D *GetHistoNoID(float etacut, char *name, bool eta, bool TPC, bool ispp, bool forSim){
1078 file->cd();
1079 char *myname = mynameITS;
1080 if(TPC) myname = mynameTPC;
0e866ddc 1081 TH2F *notid = ((TH2F*) out2->FindObject(Form("EtReconstructed%sChargedHadronAssumingPion",myname)))->Clone("notid");
1082 TH2F *nNotid = ((TH2F*) out2->FindObject(Form("EtNReconstructed%sChargedHadron",myname)))->Clone("nNotid");
b64de20c 1083
0e866ddc 1084 TH2F *id = ((TH2F*) out2->FindObject(Form("EtReconstructed%sChargedHadron",myname)))->Clone("id");
b64de20c 1085 int lowbin = id->GetYaxis()->FindBin(-etacut+.001);//make sure we don't accidentally get the wrong bin
1086 int highbin = id->GetYaxis()->FindBin(etacut-.001);
b64de20c 1087
0e866ddc 1088 TH1D *nNotidProj;
1089 TH1D *denominator;
1090 TH1D *numerator;
1091 if(eta){
1092 int lowbin = notid->GetYaxis()->FindBin(-etacut+.001);//make sure we don't accv0entally get the wrong bin
1093 int highbin = notid->GetYaxis()->FindBin(etacut-.001);
48a07264 1094 //cout<<"Projecting from "<<notid->GetYaxis()->GetBinLowEdge(lowbin)<<" to "<<notid->GetYaxis()->GetBinLowEdge(highbin+1)<<endl;
0e866ddc 1095 denominator = id->ProjectionX("name",lowbin,highbin);
1096 numerator = notid->ProjectionX("numerator",lowbin,highbin);
1097 nNotidProj = nNotid->ProjectionX("nNotidProj",lowbin,highbin);
1098 }
1099 else{
48a07264 1100 //cout<<"Getting eta dependence"<<endl;
0e866ddc 1101 int lowbin = id->GetXaxis()->FindBin(etacut);//make sure we don't accidentally get the wrong bin
1102 int highbin;
1103 if(etacut<0.15){//then we actually have ITS standalone tracks and we only want this to run from 0.1 to 0.15 because this will be used only for ITS standalone tracks
1104 highbin = id->GetXaxis()->FindBin(0.15);
1105 }
1106 else{
1107 highbin = id->GetXaxis()->GetNbins();
1108 }
48a07264 1109 //cout<<"Projecting from "<<id->GetXaxis()->GetBinLowEdge(lowbin)<<" to "<<id->GetXaxis()->GetBinLowEdge(highbin+1)<<endl;
0e866ddc 1110 numerator = notid->ProjectionY("name",lowbin,highbin);
1111 denominator = id->ProjectionY("denominator",lowbin,highbin);
1112 nNotidProj = nNotid->ProjectionY("nNotidProj",lowbin,highbin);
1113 }
a02dfa56 1114 if(denominator) numerator->Divide(denominator);
b64de20c 1115
1116 if(numerator->GetNbinsX() != nNotidProj->GetNbinsX()){
1117 cerr<<"Uh-oh! Can't rescale errors! "<<numerator->GetNbinsX()<<"!="<<nNotidProj->GetNbinsX()<<endl;
1118 return numerator;
1119 }
1120 //fixing the errors
1121 for(int i=1;i<=numerator->GetNbinsX();i++){
1122 Float_t value = numerator->GetBinContent(i);
1123 Float_t valueerr = numerator->GetBinError(i);
1124 Float_t n = nNotidProj->GetBinContent(i);
a02dfa56 1125 Float_t err=0.;
1126 if(n>0.0){
1127 err = value/TMath::Power(n,0.5);
1128 }
1129 numerator->SetBinError(i,err);;
b64de20c 1130 }
1131 numerator->SetYTitle("Ratio of E_{T}^{assuming pion}/E_{T}^{real}");
1132 numerator->GetYaxis()->SetTitleOffset(1.2);
a02dfa56 1133 delete denominator;
1134 delete nNotidProj;
48a07264 1135 delete notid;
1136 delete nNotid;
1137 delete id;
1138 numerator->SetName(name);
b64de20c 1139 return numerator;
1140
1141}
1142
48a07264 1143TH1D *CorrNoID(float etacut,char *name, char *prodname, char *shortprodname, bool ispp, bool forSim){
b64de20c 1144 gStyle->SetOptTitle(0);
1145 gStyle->SetOptStat(0);
1146 gStyle->SetOptFit(0);
1147 TCanvas *c = new TCanvas("c","c",500,400);
1148 c->SetTopMargin(0.04);
1149 c->SetRightMargin(0.04);
1150 c->SetBorderSize(0);
1151 c->SetFillColor(0);
1152 c->SetFillColor(0);
1153 c->SetBorderMode(0);
1154 c->SetFrameFillColor(0);
1155 c->SetFrameBorderMode(0);
1156
48a07264 1157 TH1D *PHOS = GetHistoNoID(etacut,name,true,true,ispp,forSim);
b64de20c 1158 PHOS->SetMarkerColor(2);
1159 PHOS->SetLineColor(2);
1160 PHOS->SetAxisRange(0.0,4);
1161 PHOS->SetMaximum(1.1);
1162 PHOS->SetMinimum(0.85);
1163 PHOS->SetMarkerStyle(20);;
1164 PHOS->Draw();
1165 TLatex *tex = new TLatex(0.161478,1.0835,prodname);
1166 tex->SetTextSize(0.0537634);
1167 tex->Draw();
1168
1169
48a07264 1170 char *detector = detectorEMCAL;
1171 if(etacut<0.2) detector = detectorPHOS;
b64de20c 1172 sprintf(epsname,"pics/%s/fnoid%s.eps",shortprodname,detector);
1173 sprintf(pngname,"pics/%s/fnoid%s.png",shortprodname,detector);
1174
1175 c->SaveAs(epsname);
1176 c->SaveAs(pngname);
a02dfa56 1177 delete c;
48a07264 1178 delete tex;
1179 PHOS->SetName(name);
b64de20c 1180 return PHOS;
1181
0e866ddc 1182}
1183
48a07264 1184Float_t CorrNoIDConst(float etacut, float ptcut,char *name, char *prodname, char *shortprodname, bool ispp, bool forSim){
0e866ddc 1185 gStyle->SetOptTitle(0);
1186 gStyle->SetOptStat(0);
1187 gStyle->SetOptFit(0);
1188 TCanvas *c = new TCanvas("c","c",500,400);
1189 c->SetTopMargin(0.04);
1190 c->SetRightMargin(0.04);
1191 c->SetBorderSize(0);
1192 c->SetFillColor(0);
1193 c->SetFillColor(0);
1194 c->SetBorderMode(0);
1195 c->SetFrameFillColor(0);
1196 c->SetFrameBorderMode(0);
1197
1198 bool TPC = true;
1199 if(ptcut<.15) TPC = false;
48a07264 1200 TH1D *PHOS = GetHistoNoID(ptcut,name,false,TPC,ispp,forSim);
0e866ddc 1201 TF1 *func = new TF1("func","[0]",-etacut,etacut);
1202 PHOS->Fit(func,"","",-etacut,etacut);
1203 PHOS->SetMarkerColor(2);
1204 PHOS->SetLineColor(2);
1205 PHOS->SetMaximum(1.1);
1206 PHOS->SetMinimum(0.85);
1207 PHOS->SetMarkerStyle(20);;
1208 PHOS->Draw();
1209 TLatex *tex = new TLatex(0.161478,1.0835,prodname);
1210 tex->SetTextSize(0.0537634);
1211 tex->Draw();
1212
1213
48a07264 1214 char *detector = detectorEMCAL;
1215 if(etacut<0.2) detector = detectorPHOS;
0e866ddc 1216 if(TPC){
1217 sprintf(epsname,"pics/%s/fnoid%sTPC.eps",shortprodname,detector);
1218 sprintf(pngname,"pics/%s/fnoid%sTPC.png",shortprodname,detector);
1219 }
1220 else{
1221 sprintf(epsname,"pics/%s/fnoid%sITS.eps",shortprodname,detector);
1222 sprintf(pngname,"pics/%s/fnoid%sITS.png",shortprodname,detector);
1223 }
1224
1225 c->SaveAs(epsname);
1226 c->SaveAs(pngname);
1227 delete c;
48a07264 1228 delete PHOS;
1229 float value = func->GetParameter(0);
1230 delete func;
1231 delete tex;
1232 return value;
0e866ddc 1233
b64de20c 1234}
1235//==================================Efficiency=================================================
1236TH1D* bayneseffdiv(TH1D* numerator, TH1D* denominator,Char_t* name)
1237{
1238 if(!numerator){
1239 cerr<<"Error: numerator does not exist!"<<endl;
1240 return NULL;
1241 }
1242 if(!denominator){
1243 cerr<<"Error: denominator does not exist!"<<endl;
1244 return NULL;
1245 }
1246 TH1D* result = (TH1D*)numerator->Clone(name);
1247 Int_t nbins = numerator->GetNbinsX();
1248 for (Int_t ibin=0; ibin<= nbins+1; ++ibin) {
1249 Double_t numeratorVal = numerator->GetBinContent(ibin);
1250 Double_t denominatorVal = denominator->GetBinContent(ibin);
1251 // Check if the errors are right or the thing is scaled
1252 Double_t numeratorValErr = numerator->GetBinError(ibin);
1253 if (!(numeratorValErr==0. || numeratorVal ==0.) ) {
1254 Double_t rescale = numeratorValErr*numeratorValErr/numeratorVal;
1255 numeratorVal /= rescale;
1256 }
1257 Double_t denominatorValErr = denominator->GetBinError(ibin);
1258 if (!(denominatorValErr==0. || denominatorVal==0. )) {
1259 Double_t rescale = denominatorValErr*denominatorValErr/denominatorVal;
1260 denominatorVal /= rescale;
1261 }
1262 Double_t quotient = 0.;
1263 if (denominatorVal!=0.) {
1264 quotient = numeratorVal/denominatorVal;
1265 }
1266 Double_t quotientErr=0;
1267 if(numeratorVal>0.0 && denominatorVal>0.0 && (denominatorVal+2.0)>0.0 && (denominatorVal+3.0) >0.0){
1268 double argument = (numeratorVal+1.0)/(denominatorVal+2.0)*
1269 ((numeratorVal+2.0)/(denominatorVal+3.0)-(numeratorVal+1.0)/(denominatorVal+2.0));
1270 double test = TMath::Power(TMath::Abs(argument),0.5);
1271 quotientErr = TMath::Sqrt( TMath::Abs(
1272 (numeratorVal+1.0)/(denominatorVal+2.0)*
1273 ((numeratorVal+2.0)/(denominatorVal+3.0)-(numeratorVal+1.0)/(denominatorVal+2.0))));
1274 }
1275 result->SetBinContent(ibin,quotient);
1276 result->SetBinError(ibin,quotientErr);
1277 //cout<<"Setting bin "<<ibin<<" to "<<quotient<<" "<<numeratorVal<<"/"<<denominatorVal<<endl;
1278 }
1279 return result;
1280}
1281
1282
1283
48a07264 1284TH1D *GetHistoEfficiency(float cut, char *name, int mycase, int color, int marker,bool TPC,bool ITS, int cb, int cblast){
b64de20c 1285 bool eta = true;
48a07264 1286 file->cd();
1287 char *myname = mynameITS;
1288 if(TPC&&!ITS) myname = mynameTPC;
1289 if(TPC&&ITS) myname = mynameTPCITS;
b64de20c 1290 TH2F *numeratorParent;
1291 switch(mycase){
1292 case 0:
48a07264 1293 if(cblast != -1){//add more centrality bins
1294 for(int i=cb;i<=cblast;i++){
1295 if(i==cb) numeratorParent = (TH2F*)((TH2F*) out2->FindObject(Form("EtNReconstructed%s%sCB%i",myname,"PiPlus",i)))->Clone("RecoHadron");
1296 else{numeratorParent->Add((TH2F*) out2->FindObject(Form("EtNReconstructed%s%sCB%i",myname,"PiPlus",i)));}
1297 numeratorParent->Add((TH2F*) out2->FindObject(Form("EtNReconstructed%s%sCB%i",myname,"PiMinus",i)));
1298 numeratorParent->Add((TH2F*) out2->FindObject(Form("EtNReconstructed%s%sCB%i",myname,"KMinus",i)));
1299 numeratorParent->Add((TH2F*) out2->FindObject(Form("EtNReconstructed%s%sCB%i",myname,"KPlus",i)));
1300 numeratorParent->Add((TH2F*) out2->FindObject(Form("EtNReconstructed%s%sCB%i",myname,"Proton",i)));
1301 numeratorParent->Add((TH2F*) out2->FindObject(Form("EtNReconstructed%s%sCB%i",myname,"AntiProton",i)));
1302 }
1303 }
1304 else{
1305 numeratorParent = (TH2F*)((TH2F*) out2->FindObject(Form("EtNReconstructed%s%s",myname,"PiPlus")))->Clone("RecoHadron");
1306 numeratorParent->Add((TH2F*) out2->FindObject(Form("EtNReconstructed%s%s",myname,"PiMinus")));
1307 numeratorParent->Add((TH2F*) out2->FindObject(Form("EtNReconstructed%s%s",myname,"KMinus")));
1308 numeratorParent->Add((TH2F*) out2->FindObject(Form("EtNReconstructed%s%s",myname,"KPlus")));
1309 numeratorParent->Add((TH2F*) out2->FindObject(Form("EtNReconstructed%s%s",myname,"Proton")));
1310 numeratorParent->Add((TH2F*) out2->FindObject(Form("EtNReconstructed%s%s",myname,"AntiProton")));
1311 }
b64de20c 1312 break;
1313 case 1://pion
48a07264 1314 if(cblast != -1){//add more centrality bins
1315 for(int i=cb;i<=cblast;i++){
1316 if(i==cb) numeratorParent = (TH2F*)((TH2F*) out2->FindObject(Form("EtNReconstructed%s%sCB%i",myname,"PiPlus",i)))->Clone("RecoPion");
1317 else{numeratorParent->Add((TH2F*) out2->FindObject(Form("EtNReconstructed%s%sCB%i",myname,"PiPlus",i)));}
1318 numeratorParent->Add((TH2F*) out2->FindObject(Form("EtNReconstructed%s%sCB%i",myname,"PiMinus",i)));
1319 }
1320 }
1321 else{
1322 numeratorParent = (TH2F*)((TH2F*) out2->FindObject(Form("EtNReconstructed%s%s",myname,"PiPlus")))->Clone("RecoPion");
1323 numeratorParent->Add((TH2F*) out2->FindObject(Form("EtNReconstructed%s%s",myname,"PiMinus")));
1324 }
b64de20c 1325 break;
1326 case 2://kaon
48a07264 1327 if(cblast != -1){//add more centrality bins
1328 for(int i=cb;i<=cblast;i++){
1329 if(i==cb) numeratorParent = (TH2F*)((TH2F*) out2->FindObject(Form("EtNReconstructed%s%sCB%i",myname,"KPlus",i)))->Clone("RecoKaon");
1330 else{numeratorParent->Add((TH2F*) out2->FindObject(Form("EtNReconstructed%s%sCB%i",myname,"KPlus",i)));}
1331 numeratorParent->Add((TH2F*) out2->FindObject(Form("EtNReconstructed%s%sCB%i",myname,"KMinus",i)));
1332 }
1333 }
1334 else{
1335 numeratorParent = (TH2F*)((TH2F*) out2->FindObject(Form("EtNReconstructed%s%s",myname,"KPlus")))->Clone("RecoKaon");
1336 numeratorParent->Add((TH2F*) out2->FindObject(Form("EtNReconstructed%s%s",myname,"KMinus")));
1337 }
b64de20c 1338 break;
1339 case 3://proton
48a07264 1340 if(cblast != -1){//add more centrality bins
1341 for(int i=cb;i<=cblast;i++){
1342 if(i==cb) numeratorParent = (TH2F*)((TH2F*) out2->FindObject(Form("EtNReconstructed%s%sCB%i",myname,"Proton",i)))->Clone("RecoProton");
1343 else{numeratorParent->Add((TH2F*) out2->FindObject(Form("EtNReconstructed%s%sCB%i",myname,"Proton",i)));}
1344 numeratorParent->Add((TH2F*) out2->FindObject(Form("EtNReconstructed%s%sCB%i",myname,"AntiProton",i)));
1345 }
1346 }
1347 else{
1348 numeratorParent = (TH2F*)((TH2F*) out2->FindObject(Form("EtNReconstructed%s%s",myname,"Proton")))->Clone("RecoProton");
1349 numeratorParent->Add((TH2F*) out2->FindObject(Form("EtNReconstructed%s%s",myname,"AntiProton")));
1350 }
b64de20c 1351 break;
1352 case 4://electron
48a07264 1353 numeratorParent = (TH2F*)((TH2F*) out2->FindObject(Form("EtNReconstructed%s%s%s",myname,"EPlus",cbname)))->Clone("RecoElectron");
1354 numeratorParent->Add((TH2F*) out2->FindObject(Form("EtNReconstructed%s%s%s",myname,"EMinus",cbname)));
b64de20c 1355 break;
1356 }
1357 TH2F *denominatorParent;
1358 switch(mycase){
1359 case 0:
48a07264 1360 if(cblast != -1){//add more centrality bins
1361 denominatorParent = (TH2F*)((TH2F*) out2->FindObject(Form("EtNSimulatedChargedHadronCB%i",cb)))->Clone("RecoHadron");
1362 for(int i=cb+1;i<=cblast;i++){
1363 denominatorParent->Add((TH2F*) out2->FindObject(Form("EtNSimulatedChargedHadronCB%i",i)));
1364 }
1365 }
1366 else{
1367 denominatorParent = (TH2F*)((TH2F*) out2->FindObject("EtNSimulatedChargedHadron"))->Clone("RecoHadron");
1368 }
b64de20c 1369 break;
1370 case 1://pion
48a07264 1371 if(cblast != -1){//add more centrality bins
1372 denominatorParent = (TH2F*)((TH2F*) out2->FindObject(Form("EtNSimulatedPiPlusCB%i",cb)))->Clone("RecoPion");
1373 denominatorParent->Add((TH2F*) out2->FindObject(Form("EtNSimulatedPiMinusCB%i",cb)));
1374 for(int i=cb+1;i<=cblast;i++){
1375 denominatorParent->Add((TH2F*) out2->FindObject(Form("EtNSimulatedPiPlusCB%i",i)));
1376 denominatorParent->Add((TH2F*) out2->FindObject(Form("EtNSimulatedPiMinusCB%i",i)));
1377 }
1378 }
1379 else{
1380 denominatorParent = (TH2F*)((TH2F*) out2->FindObject("EtNSimulatedPiPlus"))->Clone("RecoPion");
1381 denominatorParent->Add((TH2F*) out2->FindObject("EtNSimulatedPiMinus"));
1382 }
b64de20c 1383 break;
1384 case 2://kaon
48a07264 1385 if(cblast != -1){//add more centrality bins
1386 denominatorParent = (TH2F*)((TH2F*) out2->FindObject(Form("EtNSimulatedKPlusCB%i",cb)))->Clone("RecoKaon");
1387 denominatorParent->Add((TH2F*) out2->FindObject(Form("EtNSimulatedKMinusCB%i",cb)));
1388 for(int i=cb+1;i<=cblast;i++){
1389 denominatorParent->Add((TH2F*) out2->FindObject(Form("EtNSimulatedKPlusCB%i",i)));
1390 denominatorParent->Add((TH2F*) out2->FindObject(Form("EtNSimulatedKMinusCB%i",i)));
1391 }
1392 }
1393 else{
1394 denominatorParent = (TH2F*)((TH2F*) out2->FindObject("EtNSimulatedKPlus"))->Clone("RecoKaon");
1395 denominatorParent->Add((TH2F*) out2->FindObject("EtNSimulatedKMinus"));
1396 }
b64de20c 1397 break;
1398 case 3://proton
48a07264 1399 if(cblast != -1){//add more centrality bins
1400 for(int i=cb;i<=cblast;i++){
1401 if(cb==i)denominatorParent = (TH2F*)((TH2F*) out2->FindObject(Form("EtNSimulatedProtonCB%i",i)))->Clone("RecoProton");
1402 else{denominatorParent->Add((TH2F*) out2->FindObject(Form("EtNSimulatedProtonCB%i",i)));}
1403 denominatorParent->Add((TH2F*) out2->FindObject(Form("EtNSimulatedAntiProtonCB%i",i)));
1404 }
1405 }
1406 else{
1407 denominatorParent = (TH2F*)((TH2F*) out2->FindObject("EtNSimulatedProton"))->Clone("RecoProton");
1408 denominatorParent->Add((TH2F*) out2->FindObject("EtNSimulatedAntiProton"));
1409 }
b64de20c 1410 break;
1411 case 4://electron
1412 denominatorParent = (TH2F*)((TH2F*) out2->FindObject("EtNSimulatedEPlus"))->Clone("RecoElectron");
1413 denominatorParent->Add((TH2F*) out2->FindObject("EtNSimulatedEMinus"));
1414 break;
1415 }
b64de20c 1416 TH1D *denominator;
1417 TH1D *numerator;
1418 if(eta){
1419 int lowbin = numeratorParent->GetYaxis()->FindBin(-cut+.001);//make sure we don't accv0entally get the wrong bin
1420 int highbin = numeratorParent->GetYaxis()->FindBin(cut-.001);
a02dfa56 1421 //cout<<"Projecting from "<<numeratorParent->GetYaxis()->GetBinLowEdge(lowbin)<<" to "<<numeratorParent->GetYaxis()->GetBinLowEdge(highbin+1)<<endl;
b64de20c 1422 denominator = denominatorParent->ProjectionX(Form("garbage%s",name),lowbin,highbin);
1423 numerator = numeratorParent->ProjectionX(name,lowbin,highbin);
1424 }
1425 else{
1426 int lowbin = denominatorParent->GetXaxis()->FindBin(cut);//make sure we don't accidentally get the wrong bin
1427 int highbin = denominatorParent->GetXaxis()->GetNbins();
a02dfa56 1428 //cout<<"Projecting from "<<denominatorParent->GetXaxis()->GetBinLowEdge(lowbin)<<" to "<<denominatorParent->GetXaxis()->GetBinLowEdge(highbin+1)<<endl;
b64de20c 1429 numerator = numeratorParent->ProjectionY(name,lowbin,highbin);
1430 denominator = denominatorParent->ProjectionY(Form("denominator%s",name),lowbin,highbin);
1431 }
1432 delete numeratorParent;
1433 delete denominatorParent;
1434 //numerator->Divide(denominator);
1435 TH1D *result = bayneseffdiv((TH1D*) numerator,(TH1D*)denominator,name);
1436 //result->Rebin(2);
1437 //result->Scale(0.5);
1438 result->SetYTitle("Efficiency");
1439 result->GetYaxis()->SetTitleOffset(0.8);
1440 result->GetXaxis()->SetTitleOffset(0.8);
1441 result->GetYaxis()->SetLabelSize(0.05);
1442 result->GetXaxis()->SetLabelSize(0.05);
1443 result->GetYaxis()->SetTitleSize(0.05);
1444 result->GetXaxis()->SetTitleSize(0.05);
1445 result->SetMarkerColor(color);
1446 result->SetLineColor(color);
1447 result->SetMarkerStyle(marker);
a02dfa56 1448 //result->SetName(name);
b64de20c 1449 //result->Draw("e");
a02dfa56 1450 delete denominator;
1451 delete numerator;
48a07264 1452 delete numeratorParent;
1453 delete denominatorParent;
1454 result->SetName(name);
b64de20c 1455 return result;
1456
1457}
1458
48a07264 1459void CorrEfficiencyPlots(bool TPC, char *prodname, char *shortprodname){
f43fc416 1460 bool ITS = true;
b64de20c 1461 gStyle->SetOptTitle(0);
1462 gStyle->SetOptStat(0);
1463 gStyle->SetOptFit(0);
1464 TCanvas *c = new TCanvas("c","c",600,400);
1465 c->SetTopMargin(0.02);
1466 c->SetRightMargin(0.02);
1467 c->SetBorderSize(0);
1468 c->SetFillColor(0);
1469 c->SetFillColor(0);
1470 c->SetBorderMode(0);
1471 c->SetFrameFillColor(0);
1472 c->SetFrameBorderMode(0);
1473 //c->SetLogx();
1474
1475 int colortotal = 1;
1476 int colorpi = 2;
1477 int colork = 3;
1478 int colorp = 4;
1479 int phosmarker = 20;
1480 int emcalmarker = 24;
1481 float ptcut1 = 0.05;
1482 float ptcut2 = 0.1;
48a07264 1483 TH1D *PHOStotal = GetHistoEfficiency(0.12,"PHOStotal",0,colortotal,phosmarker,TPC,ITS);
1484 TH1D *PHOSpi = GetHistoEfficiency(0.12,"PHOSpi",1,colorpi,phosmarker,TPC,ITS);
1485 TH1D *PHOSp = GetHistoEfficiency(0.12,"PHOSp",2,colork,phosmarker,TPC,ITS);
1486 TH1D *PHOSk = GetHistoEfficiency(0.12,"PHOSk",3,colorp,phosmarker,TPC,ITS);
b64de20c 1487 if(!TPC){PHOStotal->GetXaxis()->SetRange(PHOStotal->GetXaxis()->FindBin(0.05),PHOStotal->GetXaxis()->FindBin(1.0));}
1488 else{PHOStotal->GetXaxis()->SetRange(PHOStotal->GetXaxis()->FindBin(0.15),PHOStotal->GetXaxis()->FindBin(3.0));}
1489 PHOStotal->SetMinimum(0.0);
1490 PHOStotal->SetMaximum(1.0);
1491 PHOStotal->Draw();
1492 PHOSpi->Draw("same");
1493 PHOSp->Draw("same");
1494 PHOSk->Draw("same");
48a07264 1495 TH1D *EMCALtotal = GetHistoEfficiency(0.7,"EMCALtotal",0,colortotal,emcalmarker,TPC,ITS);
1496 TH1D *EMCALpi = GetHistoEfficiency(0.7,"EMCALpi",1,colorpi,emcalmarker,TPC,ITS);
1497 TH1D *EMCALp = GetHistoEfficiency(0.7,"EMCALp",2,colork,emcalmarker,TPC,ITS);
1498 TH1D *EMCALk = GetHistoEfficiency(0.7,"EMCALk",3,colorp,emcalmarker,TPC,ITS);
b64de20c 1499 EMCALtotal->Draw("same");
1500 EMCALpi->Draw("same");
1501 EMCALp->Draw("same");
1502 EMCALk->Draw("same");
1503
1504
1505 TLegend *leg = new TLegend(0.22651,0.247312,0.370805,0.438172);
1506 leg->AddEntry(PHOStotal,"#pi,K,p");
1507 leg->AddEntry(PHOSpi,"#pi^{#pm}");
1508 leg->AddEntry(PHOSk,"K^{#pm}");
1509 leg->AddEntry(PHOSp,"p,#bar{p}");
1510 leg->SetFillStyle(0);
1511 leg->SetFillColor(0);
1512 leg->SetBorderSize(0);
1513 leg->SetTextSize(0.06);
1514 leg->Draw();
1515
1516 TLine *line = new TLine(0.2,0.0,0.2,1.0);
1517 line->Draw();
1518 line->SetLineWidth(3.0);
1519 //line->SetLineColor(TColor::kYellow);
1520 line->SetLineStyle(2);
1521 TLatex *tex = new TLatex(0.497269,0.0513196,prodname);
1522 tex->SetTextSize(0.0537634);
1523 tex->Draw();
1524 TLatex *tex3 = new TLatex(1.16186,0.28348,"Closed symbols |#eta|<0.12 (PHOS)");
1525 tex3->SetTextSize(0.0537634);
1526 tex3->Draw();
1527 TLatex *tex4 = new TLatex(1.16186,0.213221,"Open symbols |#eta|<0.70 (EMCal)");
1528 tex4->SetTextSize(0.0537634);
1529 tex4->Draw();
1530 TLatex *tex2 = new TLatex(0.241937,0.448436,"Likely TPC cut-off 200 MeV/c");
1531 tex2->SetTextSize(0.0537634);
1532 tex2->Draw();
b64de20c 1533 if(TPC){
f43fc416 1534 if(ITS){
1535 sprintf(epsname,"pics/%s/CorrEfficiencyITSTPC.eps",shortprodname);
1536 sprintf(pngname,"pics/%s/CorrEfficiencyITSTPC.png",shortprodname);
1537 }
1538 else{
1539 sprintf(epsname,"pics/%s/CorrEfficiencyTPC.eps",shortprodname);
1540 sprintf(pngname,"pics/%s/CorrEfficiencyTPC.png",shortprodname);
1541 }
b64de20c 1542 }
1543 else{
1544 sprintf(epsname,"pics/%s/CorrEfficiencyITS.eps",shortprodname);
1545 sprintf(pngname,"pics/%s/CorrEfficiencyITS.png",shortprodname);
1546 }
a02dfa56 1547 delete PHOStotal;
1548 delete PHOSpi;
1549 delete PHOSp;
1550 delete PHOSk;
1551 delete EMCALtotal;
1552 delete EMCALpi;
1553 delete EMCALp;
1554 delete EMCALk;
1555 delete leg;
b64de20c 1556 c->SaveAs(epsname);
1557 c->SaveAs(pngname);
a02dfa56 1558 delete c;
48a07264 1559 delete line;
1560 delete tex;
1561 delete tex3;
b64de20c 1562}
1563
1564//==================================CorrBkgd=================================================
48a07264 1565TH1D *GetHistoCorrBkgd(float etacut,char *name, bool TPC,bool ispp,bool forSim){
1566 file->cd();
1567 char *reweightname = reweightedNo;
1568 if(!forSim) reweightname = reweightedYes;
1569 char *myname = mynameITS;
1570 if(TPC) myname = mynameTPC;
b64de20c 1571 TH2F *signal = ((TH2F*) out2->FindObject(Form("EtReconstructed%sIdentifiedPiPlus",myname)))->Clone("signal");
1572 signal->Add((TH2F*) out2->FindObject(Form("EtReconstructed%sIdentifiedPiMinus",myname)));
1573 signal->Add((TH2F*) out2->FindObject(Form("EtReconstructed%sIdentifiedKMinus",myname)));
1574 signal->Add((TH2F*) out2->FindObject(Form("EtReconstructed%sIdentifiedKPlus",myname)));
1575 signal->Add((TH2F*) out2->FindObject(Form("EtReconstructed%sIdentifiedProton",myname)));
1576 signal->Add((TH2F*) out2->FindObject(Form("EtReconstructed%sIdentifiedAntiProton",myname)));
1577 signal->Add((TH2F*) out2->FindObject(Form("EtReconstructed%sUnidentifiedAssumingPion",myname)));
1578
1579 //Et of all unidentified hadrons (plus hadrons identified as pions) calculated assuming their true mass
1580 TH2F *bkgd = ((TH2F*) out2->FindObject(Form("EtReconstructed%sMisidentifiedElectrons",myname)))->Clone("bkgd");
f43fc416 1581 bkgd->Add((TH2F*) out2->FindObject(Form("EtReconstructed%sLambdaDaughters%s",myname,reweightname)));
1582 bkgd->Add((TH2F*) out2->FindObject(Form("EtReconstructed%sAntiLambdaDaughters%s",myname,reweightname)));
1583 bkgd->Add((TH2F*) out2->FindObject(Form("EtReconstructed%sK0SDaughters%s",myname,reweightname)));
b64de20c 1584 bkgd->Add((TH2F*) out2->FindObject(Form("EtReconstructed%sXiDaughters",myname)));
1585 bkgd->Add((TH2F*) out2->FindObject(Form("EtReconstructed%sAntiXiDaughters",myname)));
1586 bkgd->Add((TH2F*) out2->FindObject(Form("EtReconstructed%sOmegaDaughters",myname)));
1587 bkgd->Add((TH2F*) out2->FindObject(Form("EtReconstructed%sAntiOmegaDaughters",myname)));
1588 int lowbin = bkgd->GetYaxis()->FindBin(-etacut+.001);//make sure we don't accidentally get the wrong bin
1589 int highbin = bkgd->GetYaxis()->FindBin(etacut-.001);
a02dfa56 1590 //cout<<"Projecting from "<<bkgd->GetYaxis()->GetBinLowEdge(lowbin)<<" to "<<bkgd->GetYaxis()->GetBinLowEdge(highbin+1)<<endl;
b64de20c 1591
1592
1593 TH1D *denominator = signal->ProjectionX(name,lowbin,highbin);
1594 TH1D *numerator = bkgd->ProjectionX("numerator",lowbin,highbin);
1595 numerator->Divide(denominator);
1596 numerator->SetYTitle("Ratio of E_{T}^{background}/E_{T}^{real}");
1597 numerator->GetYaxis()->SetTitleOffset(1.2);
a02dfa56 1598 delete signal;
1599 delete bkgd;
1600 delete denominator;
48a07264 1601 numerator->SetName(name);
b64de20c 1602 return numerator;
1603
1604}
1605
48a07264 1606void CorrBkgdPlots(char *prodname, char *shortprodname, bool TPC,bool ispp,bool forSim){
b64de20c 1607 gStyle->SetOptTitle(0);
1608 gStyle->SetOptStat(0);
1609 gStyle->SetOptFit(0);
1610 TCanvas *c = new TCanvas("c","c",500,400);
1611 c->SetTopMargin(0.04);
1612 c->SetRightMargin(0.04);
1613 c->SetBorderSize(0);
1614 c->SetFillColor(0);
1615 c->SetFillColor(0);
1616 c->SetBorderMode(0);
1617 c->SetFrameFillColor(0);
1618 c->SetFrameBorderMode(0);
1619
48a07264 1620 TH1D *PHOS = GetHistoCorrBkgd(0.12,"PHOS2",TPC,ispp,forSim);
1621 TH1D *EMCAL = GetHistoCorrBkgd(0.7,"EMCAL2",TPC,ispp,forSim);
b64de20c 1622 PHOS->SetMarkerColor(2);
1623 EMCAL->SetMarkerColor(4);
1624 PHOS->SetLineColor(2);
1625 EMCAL->SetLineColor(4);
1626 //EMCAL->SetLineWidth(2);
1627 //PHOS->SetAxisRange(0.0,4);
1628 PHOS->SetMaximum(0.2);
1629 PHOS->SetMinimum(0.0);
1630 PHOS->SetMarkerStyle(20);;
1631 EMCAL->SetMarkerStyle(21);
1632 // TF1 *funcEMCAL = new TF1("funcEMCAL","[0]+0.0*x",0.05,4);
1633// funcEMCAL->SetParameter(0,0.95);
1634// funcEMCAL->SetParLimits(0,0.9,1.1);
1635 //EMCAL->Fit(funcEMCAL);//,"","",0.05,3.0);
1636// TF1 *funcPHOS = new TF1("funcPHOS","[0]+0.0*x",0.05,4);
1637// funcPHOS->SetParameter(0,1.0);
1638 //PHOS->Fit(funcPHOS);
1639 if(TPC) PHOS->GetXaxis()->SetRange(PHOS->GetXaxis()->FindBin(0.0),PHOS->GetXaxis()->FindBin(4.));
1640 else{ PHOS->GetXaxis()->SetRange(PHOS->GetXaxis()->FindBin(0.0),PHOS->GetXaxis()->FindBin(1.));}
1641 PHOS->Draw();
1642 EMCAL->Draw("same");
1643 TLatex *tex = new TLatex(0.161478,1.0835,prodname);
1644 tex->SetTextSize(0.0537634);
1645 tex->Draw();
1646 TLegend *leg = new TLegend(0.145161,0.604839,0.40121,0.860215);
1647 leg->AddEntry(PHOS,"|#eta|<0.12");
1648 leg->AddEntry(EMCAL,"|#eta|<0.70");
1649 leg->SetFillStyle(0);
1650 leg->SetFillColor(0);
1651 leg->SetBorderSize(0);
1652 leg->Draw();
1653 if(TPC){
1654 c->SaveAs(Form("pics/%s/bkgdTPC.eps",shortprodname));
1655 c->SaveAs(Form("pics/%s/bkgdTPC.png",shortprodname));
1656 }
1657 else{
1658 c->SaveAs(Form("pics/%s/bkgdITS.eps",shortprodname));
1659 c->SaveAs(Form("pics/%s/bkgdITS.png",shortprodname));
1660 }
a02dfa56 1661 delete c;
1662 delete PHOS;
1663 delete EMCAL;
48a07264 1664 delete tex;
1665 delete leg;
b64de20c 1666
1667}