]>
Commit | Line | Data |
---|---|---|
5f2ca1e6 | 1 | #if !defined(__CINT__) || defined(__MAKECINT__) |
4c0d7bc7 | 2 | #include <iostream> |
3 | #include "TH1F.h" | |
4 | #include "TGraphErrors.h" | |
5 | #include "AliBWFunc.h" | |
6 | #include "AliBWTools.h" | |
7 | #include "TF1.h" | |
8 | #include "TFile.h" | |
9 | #include "TDatabasePDG.h" | |
10 | #include "TROOT.h" | |
11 | #include "TCanvas.h" | |
12 | #include "TFolder.h" | |
13 | #include "TStyle.h" | |
14 | #include "AliLatexTable.h" | |
15 | #include "TLegend.h" | |
16 | #include "TVirtualFitter.h" | |
17 | #include "TMath.h" | |
18 | #include "TH2F.h" | |
19 | #include "TSystem.h" | |
20 | #include "TLine.h" | |
21 | #include "TLatex.h" | |
22 | #include "TMath.h" | |
23 | ||
24 | #include "TASImage.h" | |
25 | #include "TPaveText.h" | |
4c0d7bc7 | 26 | #include "StarPPSpectra.C" |
27 | #include "GetE735Ratios.C" | |
3472a634 | 28 | #include "TString.h" |
3725266c | 29 | #include "TObjString.h" |
77233d47 | 30 | #include "TGraphAsymmErrors.h" |
31 | ||
5f2ca1e6 | 32 | #endif |
4c0d7bc7 | 33 | |
3725266c | 34 | |
4c0d7bc7 | 35 | using namespace std; |
36 | ||
37 | ||
38 | // A bunch of useful enums and constants | |
39 | enum {kPion=0,kKaon,kProton,kNPart}; | |
ee08b77d | 40 | enum {kTPC=0,kTOF,kITS,kITSTPC,kK0,kKinks,kCombTOFTPC,kCombAll,kNHist};// "k0" listed here as a kind of PID method... |
4c0d7bc7 | 41 | const Int_t kNDet = kITS+2; |
42 | enum {kPos=0,kNeg,kNCharge}; | |
43 | enum {kPhojet=0,kPyTuneAtlasCSC, kPyTuneCMS6D6T, kPyTunePerugia0, kNTunes} ; | |
44 | enum {kFitLevi=0, kFitUA1, kFitPowerLaw, | |
45 | kFitPhojet, kFitAtlasCSC, kFitCMS6D6T, kFitPerugia0, | |
46 | kNFit}; | |
77233d47 | 47 | enum {kDoFits=0, kDoRatios, kDoSuperposition, kDoDrawWithModels, kDoCompareToStar, kDoDrawSyst, kDoCompareToAllCharged, kDoHelp}; |
058137a6 | 48 | enum {kStatErrors = 0, kSystErrors, kStatSystErrors}; // which errors do we put in the histo that we fit? stat,syst or stat+syst? |
4c0d7bc7 | 49 | |
50 | // flags, labels and names | |
51 | const char * partFlag[] = {"Pion", "Kaon", "Proton"}; | |
d8e6677d | 52 | const char * detFlag[] = {"TPC", "TOF", "ITS standalone", "ITS-TPC", "K0", "Kinks", "Combined TOF + TPC", "Combined TOF + TPC + ITS"}; |
4c0d7bc7 | 53 | const char * chargeFlag[] = {"Pos", "Neg"}; |
54 | const char * chargeLabel[] = {"Positive", "Negative"}; | |
55 | const char * partLabel[kNPart][kNCharge] = {{"#pi^{+}", "#pi^{-}"}, | |
56 | // {"K^{+} (#times 2)", "K^{-} (#times 2)"}, | |
57 | {"K^{+}", "K^{-}"}, | |
58 | {"p" , "#bar{p}"}}; | |
5f2ca1e6 | 59 | const char * partLatex[kNPart][kNCharge] = {{"$\\pi^{+}$", "$\\pi^{-}$"}, |
60 | // {"K^{+} (#times 2)", "K^{-} (#times 2)"}, | |
61 | {"$K^{+}$", "$K^{-}$"}, | |
62 | {"$p$" , "$\\bar{p}$"}}; | |
4c0d7bc7 | 63 | const char * mcTuneName[] = {"Phojet", |
64 | "Pythia - CSC 306", | |
65 | "Pythia - D6T 109", | |
66 | "Pythia - Perugia0 - 320", }; | |
67 | const char * funcName[] = { "Levi", "UA1", "PowerLaw", "Phojet", "AtlasCSC", "CMS6D6T", "Perugia0"}; | |
68 | ||
69 | // Style | |
70 | //const Int_t marker[] = {25,24,28,20,21}; // set for highlithining marek | |
d8e6677d | 71 | //const Int_t marker[] = {20,24,25,28,21}; // standard set |
72 | const Int_t marker[] = {24,25,28,27,21}; // No full symbols | |
4c0d7bc7 | 73 | const Int_t color [] = {1,2,4}; |
74 | ||
75 | const Int_t mcLineColor[] = {kGreen,kRed,kBlue,kBlack}; | |
76 | const Int_t mcLineStyle[] = {1,2,3,4}; | |
77 | ||
78 | ||
79 | // Template needed to combine different detectors | |
ee08b77d | 80 | const Float_t templBins[] = {0.05,0.1,0.12,0.14,0.16,0.18,0.20,0.25,0.30,0.35,0.4,0.45,0.5,0.55,0.6,0.65,0.7,0.75,0.8,0.85,0.9,0.95,1.0,1.1,1.2,1.3,1.4,1.5,1.6,1.7,1.8,2,2.2,2.4,2.6}; |
81 | Int_t nbinsTempl=34; | |
4c0d7bc7 | 82 | |
83 | TH1F * htemplate = new TH1F("htemplate", "htemplate",nbinsTempl, templBins ); | |
84 | ||
85 | // Globals | |
86 | TH1F * hSpectra[kNHist][kNPart][kNCharge]; | |
87 | TH1F * hSpectraMC[kNTunes][kNPart][kNCharge]; | |
3725266c | 88 | TH1F * hSystError[kNHist][kNPart][kNCharge]; |
4c0d7bc7 | 89 | Double_t mass[kNPart]; |
90 | ||
91 | // Functions: | |
92 | // Loading | |
93 | void LoadSpectra() ; | |
94 | void LoadMC() ; | |
95 | ||
96 | // Additional tasks (may be uncommented below) | |
3725266c | 97 | void LoadLibs(); |
4c0d7bc7 | 98 | void DrawStar(Int_t icharge); |
99 | void GetITSResiduals(); | |
100 | void DrawWithModels() ; | |
101 | void DrawAllAndKaons(); | |
102 | void DrawWithJacek(); | |
103 | void DrawRatioToStar(); | |
104 | void DrawRatios(); | |
3725266c | 105 | void DrawWithSyst(); |
5f2ca1e6 | 106 | void FitCombined(); |
3725266c | 107 | void PrintCanvas(TCanvas* c,const TString formats) ; |
4d4a268c | 108 | void RebinHistosForRatios() ; |
c4fe358b | 109 | void Help(); |
4c0d7bc7 | 110 | |
111 | // External stuff | |
5f2ca1e6 | 112 | void ALICEWorkInProgress(TCanvas *c,TString today="11/05/2010", TString label = "ALICE performance"){ |
113 | ||
55056a39 | 114 | TPad *myPadLogo = new TPad("myPadLogo", "Pad for ALICE Logo",0.72,0.72,0.89,0.89); |
115 | //TPad *myPadLogo = new TPad("myPadLogo", "Pad for ALICE Logo",0.72,0.62,0.85,0.75); | |
5f2ca1e6 | 116 | myPadLogo->SetFillColor(0); |
117 | myPadLogo->SetBorderMode(0); | |
118 | myPadLogo->SetBorderSize(2); | |
119 | myPadLogo->SetFrameBorderMode(0); | |
120 | myPadLogo->SetLeftMargin(0.0); | |
121 | myPadLogo->SetTopMargin(0.0); | |
122 | myPadLogo->SetBottomMargin(0.0); | |
123 | myPadLogo->SetRightMargin(0.0); | |
124 | myPadLogo->SetFillStyle(0); | |
125 | myPadLogo->Draw(); | |
126 | myPadLogo->cd(); | |
127 | TASImage *myAliceLogo = new TASImage("alice_logo.png"); | |
128 | myAliceLogo->Draw(); | |
55056a39 | 129 | c->cd(); |
130 | TPaveText* t1=new TPaveText(0.418103, 0.837798, 0.656609, 0.888393,"NDC"); | |
5f2ca1e6 | 131 | t1->SetFillStyle(0); |
132 | t1->SetBorderSize(0); | |
133 | t1->AddText(0.,0.,label); | |
134 | t1->SetTextColor(kRed); | |
135 | t1->SetTextFont(42); | |
55056a39 | 136 | t1->SetTextSize(0.04); |
5f2ca1e6 | 137 | t1->Draw(); |
55056a39 | 138 | TPaveText* t2=new TPaveText(0.418103, 0.80, 0.656609, 0.84,"NDC"); |
5f2ca1e6 | 139 | t2->SetFillStyle(0); |
140 | t2->SetBorderSize(0); | |
55056a39 | 141 | t2->AddText(0.,0.,"pp at #sqrt{s} = 900 GeV (2009 data)"); |
5f2ca1e6 | 142 | t2->SetTextColor(kRed); |
55056a39 | 143 | t2->SetTextFont(42); |
144 | t2->SetTextSize(0.027); | |
5f2ca1e6 | 145 | t2->Draw(); |
55056a39 | 146 | TPaveText* t3=new TPaveText(0.418103, 0.76, 0.656609, 0.80,"NDC"); |
147 | t3->SetFillStyle(0); | |
148 | t3->SetBorderSize(0); | |
149 | t3->AddText(0.,0.,"Statistical and systematic errors"); | |
150 | t3->SetTextColor(kRed); | |
151 | t3->SetTextFont(42); | |
152 | t3->SetTextSize(0.027); | |
153 | t3->Draw(); | |
154 | // TPaveText* t2=new TPaveText(0.65,0.65,0.89,0.7,"NDC"); | |
155 | // t2->SetFillStyle(0); | |
156 | // t2->SetBorderSize(0); | |
157 | // t2->SetTextColor(kRed); | |
158 | // t2->SetTextFont(52); | |
159 | // t2->AddText(0.,0.,today.Data()); | |
160 | // t2->Draw(); | |
5f2ca1e6 | 161 | } |
4c0d7bc7 | 162 | |
163 | // Used to tag plots | |
164 | TDatime dt; | |
165 | TString today = ""; | |
166 | ||
167 | ||
168 | ||
169 | // Switches | |
170 | Bool_t convertToMT = 0; | |
77233d47 | 171 | Bool_t sumCharge = 0; |
78033eeb | 172 | Int_t whatToFit = kStatErrors; |
4c0d7bc7 | 173 | Bool_t scaleKaons = kFALSE; |
0e6f7392 | 174 | Bool_t drawStar = kFALSE; // Overlay star when doing fits |
4c0d7bc7 | 175 | Bool_t correctSecondaries = 1; |
176 | Bool_t correctGeantFlukaXS = 1; | |
5f2ca1e6 | 177 | Int_t iCombInStudy = kCombAll; //kCombTOFTPC |
c4fe358b | 178 | Int_t fitFuncID = kFitLevi; |
5f2ca1e6 | 179 | Bool_t showMC=kTRUE; |
180 | Bool_t showE735=kTRUE; | |
77233d47 | 181 | Bool_t useSecCorrFromDCA=kTRUE; |
2210b002 | 182 | Bool_t flagPreliminary=kFALSE; // Add "preliminary" flag and logo to the plots |
f7bf9679 | 183 | Bool_t doPrint = 1; |
184 | Bool_t applyEtaCutCorrection=kTRUE; // Fixes a bug in the ITSsa analysis | |
2210b002 | 185 | TString printFormats = "C,eps,root"; // format in which canvases will be printed, if PrintCanvas is called (not all prints are based on printcanvas at the moment). This is a comma separated list. |
22abbb27 | 186 | |
4c0d7bc7 | 187 | |
eea641de | 188 | void CombineSpectra(Int_t analysisType=kDoHelp, Int_t locfitFuncID = kFitLevi) { //kDoSuperposition;//kDoDrawWithModels;// kDoFits; //kDoRatios; |
4c0d7bc7 | 189 | |
190 | // This macro is used to combine the 900 GeV spectra from 2009 | |
191 | ||
c4fe358b | 192 | fitFuncID=locfitFuncID; |
193 | ||
4c0d7bc7 | 194 | // Load stuff |
3725266c | 195 | LoadLibs(); |
4c0d7bc7 | 196 | |
eea641de | 197 | // Print Help and quit |
198 | if (analysisType == kDoHelp) { | |
199 | Help(); | |
200 | return; | |
201 | } | |
202 | ||
203 | ||
4c0d7bc7 | 204 | // Set date |
205 | today = today + long(dt.GetDay()) +"/" + long(dt.GetMonth()) +"/"+ long(dt.GetYear()); | |
206 | ||
207 | ||
208 | // Set Masses | |
209 | mass[kPion] = TDatabasePDG::Instance()->GetParticle("pi+")->Mass(); | |
210 | mass[kKaon] = TDatabasePDG::Instance()->GetParticle("K+")->Mass(); | |
211 | mass[kProton] = TDatabasePDG::Instance()->GetParticle("proton")->Mass(); | |
212 | ||
213 | // Load histos | |
214 | LoadSpectra(); | |
215 | LoadMC(); | |
216 | ||
217 | // Additional tasks | |
218 | //DrawStar(icharge); | |
219 | // GetITSResiduals(); | |
5f2ca1e6 | 220 | if(analysisType==kDoSuperposition) DrawAllAndKaons(); |
221 | else if(analysisType==kDoDrawWithModels) DrawWithModels() ; | |
77233d47 | 222 | else if(analysisType==kDoCompareToAllCharged) DrawWithJacek(); |
c4fe358b | 223 | else if(analysisType==kDoCompareToStar) DrawRatioToStar(); |
5f2ca1e6 | 224 | else if(analysisType==kDoRatios) DrawRatios(); |
3725266c | 225 | else if(analysisType==kDoDrawSyst) DrawWithSyst(); |
5f2ca1e6 | 226 | else if(analysisType==kDoFits) FitCombined(); |
b7018bc7 | 227 | return; |
5f2ca1e6 | 228 | } |
4c0d7bc7 | 229 | |
230 | ||
5f2ca1e6 | 231 | void FitCombined() { |
eea641de | 232 | gStyle->SetOptTitle(0); |
233 | gStyle->SetOptStat(0); | |
234 | ||
4c0d7bc7 | 235 | // Draw combined & Fit |
236 | AliBWFunc * fm = new AliBWFunc; | |
237 | fm->SetVarType(AliBWFunc::kdNdpt); | |
238 | if (convertToMT) fm->SetVarType(AliBWFunc::kOneOverMtdNdmt); | |
239 | ||
240 | // table to print results | |
2f7f947e | 241 | AliLatexTable table(10,"c|ccccccccc"); |
4c0d7bc7 | 242 | if (fitFuncID == kFitLevi) { |
55056a39 | 243 | table.InsertCustomRow("Part & Yield & Yield (FIT) & T Slope & n & $\\chi^2$/NDF & Min X & Extrap (%) & $\\langle p_{t} \\rangle (Fit)$ & $\\langle p_{t} \\rangle$ (FIT+DATA)\\\\"); |
4c0d7bc7 | 244 | } else if (fitFuncID == kFitPowerLaw) { |
55056a39 | 245 | table.InsertCustomRow("Part & Yield & Norm & n & pt0 & $\\chi^2$/NDF & Min X & Frac Above & $\\langle p_{t} \\rangle (Fit)$ & $\\langle p_{t} \\rangle$ (FIT+DATA) \\\\"); |
4c0d7bc7 | 246 | } else { |
55056a39 | 247 | table.InsertCustomRow("Part & Yield & Par0 & Par2 & Par1 & $\\chi^2$/NDF & Min X & Frac Above & $\\langle p_{t} \\rangle (Fit)$ & $\\langle p_{t} \\rangle$ (FIT+DATA) \\\\"); |
4c0d7bc7 | 248 | |
249 | } | |
2f7f947e | 250 | table.InsertHline(); |
49ad829f | 251 | AliLatexTable tempTable(3,"c|cc"); |
252 | tempTable.InsertCustomRow("Part & Yield & $\\langle p_{t} \\rangle$ \\\\"); | |
253 | tempTable.InsertHline(); | |
4c0d7bc7 | 254 | |
2210b002 | 255 | TH1F* hRatiosToFit[kNPart][kNCharge]; // Ratio data/fit function, stat error by default |
256 | TH1F* hRatiosToFitSyst[kNPart][kNCharge]; // Ratio data/fit, stat + syst | |
058137a6 | 257 | // Fit all |
258 | Int_t chargeLoop = sumCharge ? 1 : 2; | |
57a944cd | 259 | Bool_t divideRatiosCanvas=kTRUE; |
058137a6 | 260 | for(Int_t icharge = 0; icharge < chargeLoop; icharge++){ |
5f2ca1e6 | 261 | |
262 | TCanvas * c2 = new TCanvas(TString("cCombined")+chargeFlag[icharge]+"_"+funcName[fitFuncID], TString("cCombined")+chargeFlag[icharge],700,700); | |
263 | c2->SetTickx(); | |
264 | c2->SetTicky(); | |
265 | c2->SetLeftMargin(0.14); | |
57a944cd | 266 | TCanvas * c2r = 0; |
267 | if(divideRatiosCanvas) c2r = new TCanvas(TString("cCombinedRatio")+chargeFlag[icharge]+"_"+funcName[fitFuncID], TString("cCombinedRatio")+chargeFlag[icharge],1200,500); | |
268 | else c2r = new TCanvas(TString("cCombinedRatio")+chargeFlag[icharge]+"_"+funcName[fitFuncID], TString("cCombinedRatio")+chargeFlag[icharge],700,700); | |
5f2ca1e6 | 269 | c2->cd(); |
270 | gPad->SetLogy(); | |
271 | TH2F * hempty = new TH2F(TString("hempty")+long(icharge),"hempty",100,0.,2.9, 100, 0.0005,5); | |
4c0d7bc7 | 272 | hempty->SetXTitle("p_{t} (GeV/c)"); |
273 | hempty->SetYTitle("1/N_{ev} d^{2}N / dydp_{t} (GeV/c)^{-1}"); | |
5f2ca1e6 | 274 | hempty->GetYaxis()->SetTitleOffset(1.35); |
275 | hempty->GetXaxis()->SetTitleOffset(1.1); | |
4c0d7bc7 | 276 | hempty->Draw(); |
5f2ca1e6 | 277 | c2r->cd(); |
57a944cd | 278 | if(divideRatiosCanvas) { |
279 | c2r->Divide(3,1); | |
280 | c2r->cd(1); | |
281 | } | |
5f2ca1e6 | 282 | gPad->SetGridy(); |
283 | TH2F * hemptyR = new TH2F(TString("hemptyR")+long(icharge),"hemptyR",100,0.,2.9, 100, 0.5,1.5); | |
284 | hemptyR->SetXTitle("p_{t} (GeV/c)"); | |
285 | hemptyR->SetYTitle("Data/Fit"); | |
286 | hemptyR->Draw(); | |
22abbb27 | 287 | |
57a944cd | 288 | |
289 | TLegend * l = new TLegend(0.176724, 0.181548, 0.477011, 0.462798,chargeLabel[icharge]); | |
4c0d7bc7 | 290 | l->SetFillColor(kWhite); |
5f2ca1e6 | 291 | l->SetTextSize(0.035); |
55056a39 | 292 | |
293 | ||
294 | ||
295 | TPaveText* tf=new TPaveText(0.321839, 0.175595, 0.685345, 0.299107,"NDC"); | |
5f2ca1e6 | 296 | if(fitFuncID == kFitLevi){ |
297 | tf->AddText("#frac{dN}{dp_{t}} #propto p_{t} #left(1+#frac{#sqrt{m^{2}+p_{t}^{2}} -m}{nT} #right)^{-n}"); | |
298 | // tf->SetNDC(); | |
299 | tf->SetTextFont(12); | |
55056a39 | 300 | tf->SetTextSize(0.032); |
5f2ca1e6 | 301 | } |
57a944cd | 302 | |
4c0d7bc7 | 303 | for(Int_t ipart = 0; ipart < kNPart; ipart++){ |
5f2ca1e6 | 304 | printf(" ----- Fit %s %s ------\n",partFlag[ipart],chargeFlag[icharge]); |
4c0d7bc7 | 305 | Float_t fitmin = 0; |
306 | Float_t fitmax = 3; | |
307 | ||
308 | // Get functions | |
309 | TF1 * func = 0; | |
3472a634 | 310 | Int_t normPar = -1; |
4c0d7bc7 | 311 | if(fitFuncID == kFitLevi) { |
3472a634 | 312 | normPar = 0; // The levi is normalized by parameter 0 |
4c0d7bc7 | 313 | if (ipart == kPion) |
314 | func = fm->GetLevi(mass[ipart], 0.12, 7, 1.5); | |
315 | if (ipart == kKaon) | |
316 | func = fm->GetLevi(mass[ipart], 0.17, 7, 0.17); | |
317 | if (ipart == kProton) | |
5f2ca1e6 | 318 | func = fm->GetLevi(mass[ipart], 0.15, 8.5, 0.09); |
319 | } | |
4c0d7bc7 | 320 | else if(fitFuncID == kFitUA1) func = fm->GetUA1(mass[ipart],0.2,1.25,1.25,0.2,1.5); |
321 | else if(fitFuncID == kFitPowerLaw) { | |
22abbb27 | 322 | // fitmin=0; |
323 | // fitmax=1; | |
4c0d7bc7 | 324 | if (ipart == kPion) |
325 | func = fm->GetPowerLaw(1.0, 8.6, 7); | |
326 | if (ipart == kKaon) | |
327 | func = fm->GetPowerLaw(3.0, 12, 2.6); | |
328 | if (ipart == kProton) | |
329 | func = fm->GetPowerLaw(24, 72, 0.8); | |
330 | } | |
331 | else if(fitFuncID == kFitPhojet) func = fm->GetHistoFunc(hSpectraMC[kPhojet] [ipart][icharge]); | |
332 | else if(fitFuncID == kFitAtlasCSC) func = fm->GetHistoFunc(hSpectraMC[kPyTuneAtlasCSC][ipart][icharge]); | |
333 | else if(fitFuncID == kFitPerugia0) func = fm->GetHistoFunc(hSpectraMC[kPyTunePerugia0][ipart][icharge]); | |
334 | else if(fitFuncID == kFitCMS6D6T) func = fm->GetHistoFunc(hSpectraMC[kPyTuneCMS6D6T] [ipart][icharge]); | |
335 | else { | |
336 | cout << "Unknown fit ID " << fitFuncID << endl; | |
337 | return; | |
338 | } | |
339 | ||
340 | if(fitFuncID >= kFitPhojet){ | |
341 | fitmin = 0.0; | |
342 | fitmax = 1.0; | |
343 | } | |
344 | ||
058137a6 | 345 | // Temp: fit histo with sist errors |
22abbb27 | 346 | TH1F * hsyst = new TH1F(*htemplate); |
347 | TH1F * hsyststat = 0; | |
55056a39 | 348 | hsyst->SetFillColor(kGray); |
3725266c | 349 | AliBWTools::GetValueAndError(hsyst,hSpectra[iCombInStudy][ipart][icharge],hSystError[iCombInStudy][ipart][icharge],kTRUE); |
22abbb27 | 350 | hsyststat= new TH1F(*hsyst); |
351 | AliBWTools::GetHistoCombinedErrors(hsyststat,hSpectra[iCombInStudy][ipart][icharge]); // combine syst and stat | |
55056a39 | 352 | hsyststat->SetFillColor(kGray); |
3725266c | 353 | |
058137a6 | 354 | TH1F * hToFit = 0; |
355 | if (whatToFit == kStatErrors) hToFit = hSpectra[iCombInStudy][ipart][icharge]; // Shorthand | |
058137a6 | 356 | if (whatToFit == kStatSystErrors) { |
22abbb27 | 357 | hToFit = hsyststat;// Shorthand |
058137a6 | 358 | } |
359 | if (whatToFit == kSystErrors) hToFit = hsyst; | |
57b5143c | 360 | if (sumCharge) { |
361 | if(whatToFit == kStatErrors) | |
362 | hToFit->Add(hSpectra[iCombInStudy][ipart][1]); | |
363 | else if (whatToFit == kStatSystErrors) { | |
364 | TH1F * hsyst2 = new TH1F(*htemplate); | |
365 | TH1F * hsyststat2 = 0; | |
366 | hsyst2->SetFillColor(kGray); | |
367 | AliBWTools::GetValueAndError(hsyst2,hSpectra[iCombInStudy][ipart][1],hSystError[iCombInStudy][ipart][1],kTRUE); | |
368 | hsyststat2= new TH1F(*hsyst2); | |
369 | AliBWTools::GetHistoCombinedErrors(hsyststat2,hSpectra[iCombInStudy][ipart][1]); // combine syst and stat | |
370 | hToFit->Add(hsyststat2); | |
371 | } | |
372 | } | |
373 | ||
058137a6 | 374 | |
3725266c | 375 | |
376 | if(!AliBWTools::Fit(hToFit,func,fitmin,fitmax)) { | |
4c0d7bc7 | 377 | cout << " FIT ERROR " << endl; |
378 | return; | |
379 | } | |
3725266c | 380 | cout << "DRAWING" << endl; |
22abbb27 | 381 | c2->cd(); |
382 | hsyststat->Draw("same,e2"); | |
383 | // TH1F * hsyststat2 = new TH1F(*hsyststat); | |
384 | // hsyststat2->SetFillColor(kWhite); | |
385 | // hsyststat2->SetLineColor(kBlack); | |
386 | // hsyststat2->Draw("same,e2"); | |
058137a6 | 387 | hToFit->Draw("same"); |
3725266c | 388 | TF1* fitfunc=(TF1*)hToFit->GetListOfFunctions()->At(0); |
5f2ca1e6 | 389 | fitfunc->Draw("same"); |
390 | fitfunc->SetRange(0,4); | |
391 | fitfunc->SetLineColor(hSpectra[iCombInStudy][ipart][icharge]->GetLineColor()); | |
c4fe358b | 392 | if(drawStar) DrawStar(icharge); |
058137a6 | 393 | hRatiosToFit[ipart][icharge]=(TH1F*)hToFit->Clone(Form("hRatio%s%s",chargeFlag[icharge],partFlag[icharge])); // Ratio data/fit |
2210b002 | 394 | hRatiosToFitSyst[ipart][icharge]=(TH1F*)hsyststat->Clone(Form("hRatioSyst%s%s",chargeFlag[icharge],partFlag[icharge])); // Ratio data/fit |
57a944cd | 395 | // Syst error on ratio with no fill |
396 | hRatiosToFitSyst[ipart][icharge]->SetLineColor(hToFit->GetLineColor()); | |
397 | hRatiosToFitSyst[ipart][icharge]->SetFillStyle(0); | |
398 | ||
399 | ||
4d4a268c | 400 | // Compute ratio data/function integrating the function in the bin width |
058137a6 | 401 | for(Int_t iBin=1; iBin<hToFit->GetNbinsX(); iBin++){ |
402 | Double_t lowLim=hToFit->GetBinLowEdge(iBin); | |
403 | Double_t highLim=hToFit->GetBinLowEdge(iBin+1); | |
5f2ca1e6 | 404 | Double_t contFunc=fitfunc->Integral(lowLim,highLim)/(highLim-lowLim); |
058137a6 | 405 | Double_t ratio=hToFit->GetBinContent(iBin)/contFunc; |
406 | Double_t eratio=hToFit->GetBinError(iBin)/contFunc; | |
2210b002 | 407 | Double_t eratioSyst=hsyststat->GetBinError(iBin)/contFunc; |
408 | hRatiosToFitSyst[ipart][icharge]->SetBinContent(iBin,ratio); | |
409 | hRatiosToFitSyst[ipart][icharge]->SetBinError(iBin,eratioSyst); | |
5f2ca1e6 | 410 | hRatiosToFit[ipart][icharge]->SetBinContent(iBin,ratio); |
411 | hRatiosToFit[ipart][icharge]->SetBinError(iBin,eratio); | |
412 | } | |
5f2ca1e6 | 413 | c2->Update(); |
058137a6 | 414 | l->AddEntry(hToFit, |
4c0d7bc7 | 415 | scaleKaons && ipart == kKaon ? |
416 | (TString(partLabel[ipart][icharge])+" #times 2").Data() | |
417 | : partLabel[ipart][icharge]); | |
4c0d7bc7 | 418 | |
419 | // populate table | |
420 | // Float_t yield = func->Integral(0.45,1.05); | |
421 | // Float_t yieldE = func->IntegralError(0.45,1.05); | |
422 | ||
423 | Float_t yield = func->Integral(0.,100); | |
424 | //Float_t yieldE = func->IntegralError(0.,100); | |
425 | ||
426 | Double_t yieldTools, yieldETools; | |
427 | Double_t partialYields[3],partialYieldsErrors[3]; | |
058137a6 | 428 | AliBWTools::GetYield(hToFit, func, yieldTools, yieldETools, |
4c0d7bc7 | 429 | 0, 100, partialYields,partialYieldsErrors); |
430 | Double_t tslope = func->GetParameter(2); | |
431 | Double_t tslopeE = func->GetParError(2); | |
432 | ||
5f2ca1e6 | 433 | table.SetNextCol(partLatex[ipart][icharge]); |
4c0d7bc7 | 434 | //table.SetNextCol(yield,yieldE,-4); |
435 | table.SetNextCol(yieldTools, yieldETools,-4); | |
436 | table.SetNextCol(func->GetParameter(0)); | |
437 | table.SetNextCol(tslope,tslopeE,-4); | |
55056a39 | 438 | table.SetNextCol(func->GetParameter(1),func->GetParError(1),-2); |
4c0d7bc7 | 439 | table.SetNextCol(Form("%2.2f/%d",func->GetChisquare(),func->GetNDF())); |
058137a6 | 440 | Float_t lowestPoint = AliBWTools::GetLowestNotEmptyBinEdge(hToFit); |
4c0d7bc7 | 441 | //Float_t lowestPoint = AliBWTools::GetLowestNotEmptyBinEdge(hSpectra[kITS][ipart][icharge]); |
442 | Float_t yieldAbove = func->Integral(lowestPoint,100); | |
443 | table.SetNextCol(lowestPoint,-2); | |
55056a39 | 444 | table.SetNextCol((1 - yieldAbove/yield)*100,-1); |
445 | Float_t mean, meane; // Mean using only fit function | |
446 | Double_t meanDF, meanDFe; // Mean from data + fit | |
4c0d7bc7 | 447 | Float_t mean2, mean2e; |
55056a39 | 448 | if (fitFuncID < kFitPhojet) { |
449 | AliBWTools::GetMean (func, mean, meane , 0.,100., normPar); | |
450 | AliBWTools::GetMeanSquare(func, mean2, mean2e, 0.,100., normPar); | |
451 | } | |
ffcaf7be | 452 | // AliBWTools::GetMeanDataAndExtrapolation(hToFit, func, meanDF, meanDFe, 0.,100.); |
55056a39 | 453 | AliBWTools::GetMeanDataAndExtrapolation(hToFit, func, meanDF, meanDFe, 0.,100.); |
55056a39 | 454 | table.SetNextCol(mean, meane ,-4); |
455 | table.SetNextCol(meanDF, meanDFe,-4); | |
4c0d7bc7 | 456 | |
457 | // fMean2->IntegralError(0,100)/func->Integral(0,100),-7); | |
458 | table.InsertRow(); | |
459 | ||
460 | ||
461 | /// TEMP TABLE | |
49ad829f | 462 | tempTable.SetNextCol(partLatex[ipart][icharge]); |
4c0d7bc7 | 463 | tempTable.SetNextCol(yieldTools, yieldETools, -4); |
55056a39 | 464 | tempTable.SetNextCol(meanDF, meanDFe ,-4); |
49ad829f | 465 | // tempTable.SetNextCol(partialYields[1], partialYieldsErrors[1], -4); |
466 | // tempTable.SetNextCol(yieldAbove/yield,-2); | |
4c0d7bc7 | 467 | tempTable.InsertRow(); |
5f2ca1e6 | 468 | c2r->cd(); |
57a944cd | 469 | if(divideRatiosCanvas) { |
470 | c2r->cd(ipart+1); | |
471 | gPad->SetGridy(); | |
472 | TH2F * heloc = new TH2F(*hemptyR); | |
473 | heloc->Draw(); | |
474 | } | |
2210b002 | 475 | hRatiosToFitSyst[ipart][icharge]->Draw("e2same"); |
5f2ca1e6 | 476 | hRatiosToFit[ipart][icharge]->Draw("esame"); |
477 | ||
4c0d7bc7 | 478 | } |
5f2ca1e6 | 479 | c2->cd(); |
2210b002 | 480 | if (flagPreliminary) ALICEWorkInProgress(c2,"","ALICE Preliminary"); |
5f2ca1e6 | 481 | l->Draw(); |
2210b002 | 482 | if (flagPreliminary) tf->Draw(); |
5f2ca1e6 | 483 | c2r->cd(); |
2210b002 | 484 | if(flagPreliminary) ALICEWorkInProgress(c2r,"","ALICE Preliminary"); |
4c0d7bc7 | 485 | l->Draw(); |
0e019ec6 | 486 | PrintCanvas(c2,printFormats); |
487 | PrintCanvas(c2r,printFormats); | |
4c0d7bc7 | 488 | |
5f2ca1e6 | 489 | |
4c0d7bc7 | 490 | } |
491 | ||
492 | ||
57b5143c | 493 | table.PrintTable("ASCII"); |
494 | cout << "" << endl; | |
495 | table.PrintTable(""); | |
4c0d7bc7 | 496 | |
497 | cout << "" << endl; | |
498 | tempTable.PrintTable("ASCII"); | |
499 | ||
5f2ca1e6 | 500 | |
501 | ||
4c0d7bc7 | 502 | } |
503 | ||
504 | void LoadSpectra() { | |
505 | ||
3725266c | 506 | // Loads spectra for different detectors and corresponding systematic errors. |
507 | ||
4c0d7bc7 | 508 | TFile * f=0; |
509 | ||
3725266c | 510 | // Systematic errors: initialize histos |
511 | gROOT->cd(); | |
512 | for(Int_t ipart = 0; ipart < kNPart; ipart++){ | |
513 | for(Int_t idet = 0; idet < kNHist; idet++){ | |
514 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ | |
515 | hSystError[idet][ipart][icharge] = new TH1F (TString("hSyst_")+detFlag[idet]+partFlag[ipart]+chargeFlag[icharge], | |
516 | TString("hSyst_")+detFlag[idet]+partFlag[ipart]+chargeFlag[icharge], | |
517 | nbinsTempl,templBins); | |
518 | hSystError[idet][ipart][icharge]->SetMarkerColor (color[ipart] ); | |
519 | hSystError[idet][ipart][icharge]->SetLineColor (color[ipart] ); | |
520 | ||
521 | } | |
522 | ||
523 | } | |
524 | ||
525 | } | |
526 | ||
4c0d7bc7 | 527 | // Load |
528 | ||
529 | ||
530 | // TOF | |
5f2ca1e6 | 531 | // Load Efficiencies |
4d4a268c | 532 | f = new TFile("./Files/effTOFhistosSmooth.root"); |
533 | TFile * ftrk = new TFile("./Files/effhistos.root"); | |
5f2ca1e6 | 534 | TH1D * hEffTrackTOF[kNPart][kNCharge]; |
535 | TH1D * hEffMatchTOF[kNPart][kNCharge]; | |
4d4a268c | 536 | hEffTrackTOF[kPion] [kPos] = (TH1D*) ftrk->Get("hpitrk_pos"); |
537 | hEffTrackTOF[kKaon] [kPos] = (TH1D*) ftrk->Get("hkatrk_pos"); | |
538 | hEffTrackTOF[kProton][kPos] = (TH1D*) ftrk->Get("hprtrk_pos"); | |
5f2ca1e6 | 539 | hEffMatchTOF[kPion] [kPos] = (TH1D*) f->Get("hpieff_pos"); |
540 | hEffMatchTOF[kKaon] [kPos] = (TH1D*) f->Get("hkaeff_pos"); | |
541 | hEffMatchTOF[kProton][kPos] = (TH1D*) f->Get("hpreff_pos"); | |
4d4a268c | 542 | hEffTrackTOF[kPion] [kNeg] = (TH1D*) ftrk->Get("hpitrk_neg"); |
543 | hEffTrackTOF[kKaon] [kNeg] = (TH1D*) ftrk->Get("hkatrk_neg"); | |
544 | hEffTrackTOF[kProton][kNeg] = (TH1D*) ftrk->Get("hprtrk_neg"); | |
5f2ca1e6 | 545 | hEffMatchTOF[kPion] [kNeg] = (TH1D*) f->Get("hpieff_neg"); |
546 | hEffMatchTOF[kKaon] [kNeg] = (TH1D*) f->Get("hkaeff_neg"); | |
547 | hEffMatchTOF[kProton][kNeg] = (TH1D*) f->Get("hpreff_neg"); | |
548 | ||
4c0d7bc7 | 549 | // f = new TFile("./Files/spectra-pos-y.root"); |
5f2ca1e6 | 550 | f = new TFile("./Files/spectraRaw-pos-y.root"); |
4c0d7bc7 | 551 | hSpectra[kTOF][kPion] [kPos]= (TH1F*) f->Get("hpi"); |
552 | hSpectra[kTOF][kProton][kPos]= (TH1F*) f->Get("hpr"); | |
553 | hSpectra[kTOF][kKaon] [kPos]= (TH1F*) f->Get("hka"); | |
554 | hSpectra[kTOF][kPion] [kPos]->SetName("hpiPos"); | |
555 | hSpectra[kTOF][kProton][kPos]->SetName("hprPos"); | |
556 | hSpectra[kTOF][kKaon] [kPos]->SetName("hkaPos"); | |
557 | //f = new TFile("./Files/spectra-neg-y.root"); | |
5f2ca1e6 | 558 | f = new TFile("./Files/spectraRaw-neg-y.root"); |
4c0d7bc7 | 559 | hSpectra[kTOF][kPion] [kNeg]= (TH1F*) f->Get("hpi"); |
560 | hSpectra[kTOF][kProton][kNeg]= (TH1F*) f->Get("hpr"); | |
561 | hSpectra[kTOF][kKaon] [kNeg]= (TH1F*) f->Get("hka"); | |
562 | hSpectra[kTOF][kPion] [kNeg]->SetName("hpiNeg"); | |
563 | hSpectra[kTOF][kProton][kNeg]->SetName("hprNeg"); | |
564 | hSpectra[kTOF][kKaon] [kNeg]->SetName("hkaNeg"); | |
565 | ||
5f2ca1e6 | 566 | // Divide for efficiency |
567 | hSpectra[kTOF][kPion] [kPos]->Divide(hEffTrackTOF[kPion] [kPos]); | |
568 | hSpectra[kTOF][kKaon] [kPos]->Divide(hEffTrackTOF[kKaon] [kPos]); | |
569 | hSpectra[kTOF][kProton][kPos]->Divide(hEffTrackTOF[kProton][kPos]); | |
570 | hSpectra[kTOF][kPion] [kPos]->Divide(hEffMatchTOF[kPion] [kPos]); | |
571 | hSpectra[kTOF][kKaon] [kPos]->Divide(hEffMatchTOF[kKaon] [kPos]); | |
572 | hSpectra[kTOF][kProton][kPos]->Divide(hEffMatchTOF[kProton][kPos]); | |
573 | hSpectra[kTOF][kPion] [kNeg]->Divide(hEffTrackTOF[kPion] [kNeg]); | |
574 | hSpectra[kTOF][kKaon] [kNeg]->Divide(hEffTrackTOF[kKaon] [kNeg]); | |
575 | hSpectra[kTOF][kProton][kNeg]->Divide(hEffTrackTOF[kProton][kNeg]); | |
576 | hSpectra[kTOF][kPion] [kNeg]->Divide(hEffMatchTOF[kPion] [kNeg]); | |
577 | hSpectra[kTOF][kKaon] [kNeg]->Divide(hEffMatchTOF[kKaon] [kNeg]); | |
578 | hSpectra[kTOF][kProton][kNeg]->Divide(hEffMatchTOF[kProton][kNeg]); | |
579 | ||
580 | ||
581 | // Clean UP TOF spectra, removing unwanted points | |
4c0d7bc7 | 582 | cout << "Cleaning Up TOF spectra" << endl; |
583 | Int_t nbin = hSpectra[kTOF][kKaon][kPos]->GetNbinsX(); | |
584 | for(Int_t ibin = 1; ibin <= nbin; ibin++){ | |
585 | Float_t pt = hSpectra[kTOF][kKaon][kPos]->GetBinCenter(ibin); | |
586 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ | |
587 | if (pt > 2.35) { | |
588 | hSpectra[kTOF][kKaon][icharge]->SetBinContent(ibin,0); | |
589 | hSpectra[kTOF][kKaon][icharge]->SetBinError (ibin,0); | |
590 | hSpectra[kTOF][kProton][icharge]->SetBinContent(ibin,0); | |
591 | hSpectra[kTOF][kProton][icharge]->SetBinError (ibin,0); | |
4d4a268c | 592 | } |
593 | if (pt < 0.7) { | |
594 | hSpectra[kTOF][kProton][icharge]->SetBinContent(ibin,0); | |
595 | hSpectra[kTOF][kProton][icharge]->SetBinError (ibin,0); | |
596 | } | |
4c0d7bc7 | 597 | } |
598 | } | |
599 | // cout << "Scaling TOF to TPC" << endl; | |
600 | // // Scale protons to TPC | |
601 | // hSpectra[kTOF][kProton][kPos]->Scale(1./1.05); | |
602 | // // Scale pbar so that pbar/p is compatible with Panos | |
603 | // hSpectra[kTOF][kProton][kNeg]->Scale(1./1.1); | |
604 | ||
3725266c | 605 | // TOF: systematics |
606 | // Load TOF systematic errors: | |
607 | f = new TFile ("./Files/systMatchingPos.root"); | |
608 | AliBWTools::AddHisto(hSystError[kTOF][kPion][kPos] ,(TH1*)gDirectory->Get("hErrPionMatch") ); | |
609 | AliBWTools::AddHisto(hSystError[kTOF][kProton][kPos],(TH1*)gDirectory->Get("hErrProtonMatch")); | |
610 | AliBWTools::AddHisto(hSystError[kTOF][kKaon][kPos] ,(TH1*)gDirectory->Get("hErrKaonMatch") ); | |
611 | f = new TFile ("./Files/systMatchingNeg.root"); | |
612 | AliBWTools::AddHisto(hSystError[kTOF][kPion] [kNeg],(TH1*)gDirectory->Get("hErrPionMatch")); | |
613 | AliBWTools::AddHisto(hSystError[kTOF][kProton][kNeg],(TH1*)gDirectory->Get("hErrProtonMatch")); | |
614 | AliBWTools::AddHisto(hSystError[kTOF][kKaon] [kNeg],(TH1*)gDirectory->Get("hErrKaonMatch")); | |
615 | f = new TFile ("./Files/systPIDall.root"); | |
616 | AliBWTools::AddHisto(hSystError[kTOF][kPion] [kPos],(TH1*)gDirectory->Get("hpiCorr")); | |
617 | AliBWTools::AddHisto(hSystError[kTOF][kProton][kPos],(TH1*)gDirectory->Get("hpCorr")); | |
618 | AliBWTools::AddHisto(hSystError[kTOF][kKaon] [kPos],(TH1*)gDirectory->Get("hkCorr")); | |
619 | AliBWTools::AddHisto(hSystError[kTOF][kPion] [kNeg],(TH1*)gDirectory->Get("hpiCorr")); | |
620 | AliBWTools::AddHisto(hSystError[kTOF][kProton][kNeg],(TH1*)gDirectory->Get("hpCorr")); | |
621 | AliBWTools::AddHisto(hSystError[kTOF][kKaon] [kNeg],(TH1*)gDirectory->Get("hkCorr")); | |
4c0d7bc7 | 622 | |
623 | ||
5f2ca1e6 | 624 | // ITS SA |
2210b002 | 625 | f = new TFile("./Files/ITSsaSpectraCorr-28Sep10.root"); |
4c0d7bc7 | 626 | hSpectra[kITS][kPion] [kPos]= (TH1F*) f->Get("hSpectraPos0"); |
627 | hSpectra[kITS][kKaon] [kPos]= (TH1F*) f->Get("hSpectraPos1"); | |
628 | hSpectra[kITS][kProton][kPos]= (TH1F*) f->Get("hSpectraPos2"); | |
629 | hSpectra[kITS][kPion] [kNeg]= (TH1F*) f->Get("hSpectraNeg0"); | |
630 | hSpectra[kITS][kKaon] [kNeg]= (TH1F*) f->Get("hSpectraNeg1"); | |
631 | hSpectra[kITS][kProton][kNeg]= (TH1F*) f->Get("hSpectraNeg2"); | |
632 | ||
4d4a268c | 633 | // Remove unwanted bins |
4c0d7bc7 | 634 | for(Int_t ipart = 0; ipart < kNPart; ipart++){ |
635 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ | |
2210b002 | 636 | Int_t nbinITS = hSpectra[kITS][ipart][icharge]->GetNbinsX(); |
637 | for(Int_t ibin = 1; ibin <= nbinITS; ibin++){ | |
4c0d7bc7 | 638 | if(hSpectra[kITS][ipart][icharge]->GetBinContent(ibin) < 0 ) { |
639 | hSpectra[kITS][ipart][icharge]->SetBinContent(ibin,0); | |
640 | hSpectra[kITS][ipart][icharge]->SetBinError (ibin,0); | |
641 | } | |
f7bf9679 | 642 | // if(ipart == kProton && (ibin==9 || ibin==10)){ // FIXME |
643 | if(ipart == kProton && (ibin==9)){ | |
5f2ca1e6 | 644 | printf("Kill bin %d (%f - %f GeV/c)for ITS protons\n",ibin,hSpectra[kITS][ipart][icharge]->GetBinLowEdge(ibin),hSpectra[kITS][ipart][icharge]->GetBinLowEdge(ibin)+hSpectra[kITS][ipart][icharge]->GetBinWidth(ibin)); |
645 | hSpectra[kITS][ipart][icharge]->SetBinContent(ibin,0); | |
646 | hSpectra[kITS][ipart][icharge]->SetBinError (ibin,0); | |
647 | } | |
f7bf9679 | 648 | // FIXME: uncomment? |
649 | // if(ipart == kKaon && ibin==7){ | |
650 | // printf("Kill bin %d (%f - %f GeV/c)for ITS kaons\n",ibin,hSpectra[kITS][ipart][icharge]->GetBinLowEdge(ibin),hSpectra[kITS][ipart][icharge]->GetBinLowEdge(ibin)+hSpectra[kITS][ipart][icharge]->GetBinWidth(ibin)); | |
651 | // hSpectra[kITS][ipart][icharge]->SetBinContent(ibin,0); | |
652 | // hSpectra[kITS][ipart][icharge]->SetBinError (ibin,0); | |
653 | // } | |
57a944cd | 654 | |
4c0d7bc7 | 655 | // if ((ipart == kKaon && ibin >= 12) || (ipart == kProton && ibin >= 20)) { |
656 | // hSpectra[kITS][ipart][icharge]->SetBinContent(ibin,0); | |
657 | // hSpectra[kITS][ipart][icharge]->SetBinError (ibin,0); | |
658 | // } | |
659 | } | |
660 | ||
661 | } | |
662 | } | |
663 | ||
eea641de | 664 | if(useSecCorrFromDCA){ |
77233d47 | 665 | TFile* fseccorr = new TFile("./Files/CorrFac-SecProtons3Meth-ITSsa.root"); |
666 | TH1F* hcorrp=(TH1F*)fseccorr->Get("hSecPCorrFromDCAexpo"); | |
667 | TH1F* hcorrpbar=(TH1F*)fseccorr->Get("hSecPbarCorrFromDCAexpo"); | |
eea641de | 668 | hSpectra[kITS][kProton][kPos]->Multiply(hcorrp); |
669 | hSpectra[kITS][kProton][kNeg]->Multiply(hcorrpbar); | |
670 | fseccorr->Close(); | |
671 | } | |
4c0d7bc7 | 672 | |
f7bf9679 | 673 | // Correct for the eta bug (see mails around 30/01/2011) |
674 | if(applyEtaCutCorrection) { | |
675 | TH1F * hCorrectionITSsaEtaCut[kNPart][kNCharge] = {0}; | |
676 | f= TFile::Open("./Files/ITSsaEtaCutCorr.root"); | |
677 | hCorrectionITSsaEtaCut[kPion ][kPos] = (TH1F*) f->Get("ITSsaPos0"); | |
678 | hCorrectionITSsaEtaCut[kKaon ][kPos] = (TH1F*) f->Get("ITSsaPos1"); | |
679 | hCorrectionITSsaEtaCut[kProton][kPos] = (TH1F*) f->Get("ITSsaPos2"); | |
680 | hCorrectionITSsaEtaCut[kPion ][kNeg] = (TH1F*) f->Get("ITSsaNeg0"); | |
681 | hCorrectionITSsaEtaCut[kKaon ][kNeg] = (TH1F*) f->Get("ITSsaNeg1"); | |
682 | hCorrectionITSsaEtaCut[kProton][kNeg] = (TH1F*) f->Get("ITSsaNeg2"); | |
683 | for(Int_t ipart = 0; ipart < kNPart; ipart++){ | |
684 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ | |
685 | Int_t nbin = hSpectra[kITS][ipart][icharge]->GetNbinsX(); | |
686 | for(Int_t ibin = 0; ibin < nbin; ibin++){ | |
687 | Int_t content = hSpectra[kITS][ipart][icharge]->GetBinContent(ibin) * hCorrectionITSsaEtaCut[ipart][icharge]->GetBinContent(ibin); | |
688 | Int_t error = hSpectra[kITS][ipart][icharge]->GetBinError(ibin) * hCorrectionITSsaEtaCut[ipart][icharge]->GetBinContent(ibin); | |
689 | hSpectra[kITS][ipart][icharge]->SetBinContent(ibin, content); | |
690 | hSpectra[kITS][ipart][icharge]->SetBinError (ibin, error); | |
691 | } | |
692 | } | |
693 | } | |
694 | } | |
695 | ||
3725266c | 696 | // Load ITS sa systematics, only pt dependent ones |
4d4a268c | 697 | // f = TFile::Open("./Files/ITSsa-systematics_20100930.root"); |
698 | f = TFile::Open("./Files/ITSsa-systematics-20101014.root"); | |
3725266c | 699 | for(Int_t ipart = 0; ipart < kNPart; ipart++){ |
700 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ | |
701 | AliBWTools::AddHisto(hSystError[kITS][ipart][icharge], (TH1*) gDirectory->Get(Form("hSystTot%s%s",chargeFlag[icharge],partFlag[ipart]))); // Using total error | |
702 | } | |
703 | } | |
704 | ||
4c0d7bc7 | 705 | // ITS + TPC (Marek) |
0e019ec6 | 706 | //f = TFile::Open("./Files/SpectraCorrectedITSBeforeProtons20100720.root"); |
57b5143c | 707 | // f = TFile::Open("./Files/SpectraCorrectedITSBeforeProtons14092010new.root"); |
708 | f = TFile::Open("./Files/SpectraCorrectedITSBeforeProtonsstat29102010.root"); | |
4c0d7bc7 | 709 | TList * list = (TList*) gDirectory->Get("output"); |
710 | hSpectra[kITSTPC][kPion] [kPos]= (TH1F*) list->FindObject("Pions"); | |
711 | hSpectra[kITSTPC][kKaon] [kPos]= (TH1F*) list->FindObject("Kaons"); | |
712 | hSpectra[kITSTPC][kProton][kPos]= (TH1F*) list->FindObject("Protons"); | |
713 | hSpectra[kITSTPC][kPion] [kNeg]= (TH1F*) list->FindObject("AntiPions"); | |
714 | hSpectra[kITSTPC][kKaon] [kNeg]= (TH1F*) list->FindObject("AntiKaons"); | |
715 | hSpectra[kITSTPC][kProton][kNeg]= (TH1F*) list->FindObject("AntiProtons"); | |
57b5143c | 716 | // Systematics |
717 | f = TFile::Open("./Files/SpectraCorrectedITSBeforeProtonssys29102010.root"); | |
718 | list = (TList*) gDirectory->Get("output"); | |
719 | hSystError[kITSTPC][kPion] [kPos]= (TH1F*) AliBWTools::GetRelativeError((TH1*) list->FindObject("Pions")); | |
720 | hSystError[kITSTPC][kKaon] [kPos]= (TH1F*) AliBWTools::GetRelativeError((TH1*) list->FindObject("Kaons")); | |
721 | hSystError[kITSTPC][kProton][kPos]= (TH1F*) AliBWTools::GetRelativeError((TH1*) list->FindObject("Protons")); | |
722 | hSystError[kITSTPC][kPion] [kNeg]= (TH1F*) AliBWTools::GetRelativeError((TH1*) list->FindObject("AntiPions")); | |
723 | hSystError[kITSTPC][kKaon] [kNeg]= (TH1F*) AliBWTools::GetRelativeError((TH1*) list->FindObject("AntiKaons")); | |
724 | hSystError[kITSTPC][kProton][kNeg]= (TH1F*) AliBWTools::GetRelativeError((TH1*) list->FindObject("AntiProtons")); | |
725 | ||
726 | ||
4c0d7bc7 | 727 | |
728 | // TPC | |
4c0d7bc7 | 729 | f = new TFile("./Files/protonSpectra_20100615.root"); |
730 | hSpectra[kTPC][kProton][kPos]= AliBWTools::GetHistoFromGraph((TGraphErrors*)f->Get("protonPosClassic"),htemplate); | |
731 | hSpectra[kTPC][kProton][kNeg]= AliBWTools::GetHistoFromGraph((TGraphErrors*)f->Get("protonNegClassic"),htemplate); | |
732 | f = new TFile("./Files/pionSpectra_20100615.root"); | |
733 | hSpectra[kTPC][kPion][kPos]= AliBWTools::GetHistoFromGraph((TGraphErrors*)f->Get("pionPosClassic"),htemplate); | |
734 | hSpectra[kTPC][kPion][kNeg]= AliBWTools::GetHistoFromGraph((TGraphErrors*)f->Get("pionNegClassic"),htemplate); | |
5b103cfc | 735 | // f = new TFile("./Files/kaonSpectra_20100615.root"); |
736 | f = new TFile("./Files/kaonsTPCcorr_2010_08_31.root"); | |
4c0d7bc7 | 737 | hSpectra[kTPC][kKaon][kPos]= AliBWTools::GetHistoFromGraph((TGraphErrors*)f->Get("kaonPosClassic"),htemplate); |
738 | hSpectra[kTPC][kKaon][kNeg]= AliBWTools::GetHistoFromGraph((TGraphErrors*)f->Get("kaonNegClassic"),htemplate); | |
739 | ||
740 | // Clean UP TPC spectra, removing unwanted points | |
741 | cout << "Cleaning Up TPC spectra" << endl; | |
742 | nbin = hSpectra[kTPC][kKaon][kPos]->GetNbinsX(); | |
743 | for(Int_t ibin = 0; ibin < nbin; ibin++){ | |
744 | Float_t pt = hSpectra[kTPC][kKaon][kPos]->GetBinCenter(ibin); | |
57a944cd | 745 | if (pt > 0.45 || pt<0.25) { |
4c0d7bc7 | 746 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ |
57a944cd | 747 | hSpectra[kTPC][kKaon][icharge]->SetBinContent(ibin,0); |
748 | hSpectra[kTPC][kKaon][icharge]->SetBinError (ibin,0); | |
4c0d7bc7 | 749 | } |
750 | } | |
57a944cd | 751 | // if (pt < 0.25) { |
752 | // hSpectra[kTPC][kKaon][kNeg]->SetBinContent(ibin,0); | |
753 | // hSpectra[kTPC][kKaon][kNeg]->SetBinError (ibin,0); | |
754 | // } | |
4c0d7bc7 | 755 | if (pt < 0.45) { |
756 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ | |
757 | hSpectra[kTPC][kProton][icharge]->SetBinContent(ibin,0); | |
758 | hSpectra[kTPC][kProton][icharge]->SetBinError (ibin,0); | |
759 | } | |
760 | } | |
761 | } | |
762 | ||
3725266c | 763 | // Load TPC systematics |
764 | cout << "WARNING: TPC SYST FOR NEGATIVES TO BE CORRECTED" << endl; | |
765 | f = TFile ::Open("./Files/pionsSystSum.root"); | |
766 | AliBWTools::AddHisto(hSystError[kTPC][kPion][kPos],(TH1*) gDirectory->Get("pionsSyst")); | |
767 | AliBWTools::AddHisto(hSystError[kTPC][kPion][kNeg],(TH1*) gDirectory->Get("pionsSyst")); | |
768 | f = TFile ::Open("./Files/kaonsSystSum.root"); | |
769 | AliBWTools::AddHisto(hSystError[kTPC][kKaon][kPos],(TH1*) gDirectory->Get("kaonsSyst")); | |
770 | AliBWTools::AddHisto(hSystError[kTPC][kKaon][kNeg],(TH1*) gDirectory->Get("kaonsSyst")); | |
771 | f = TFile ::Open("./Files/ProtonsSystSum.root"); | |
772 | AliBWTools::AddHisto(hSystError[kTPC][kProton][kPos],(TH1*) gDirectory->Get("ProtonsSyst")); | |
773 | AliBWTools::AddHisto(hSystError[kTPC][kProton][kNeg],(TH1*) gDirectory->Get("ProtonsSyst")); | |
774 | ||
4c0d7bc7 | 775 | // K0s |
49ad829f | 776 | f = new TFile ("./Files/PtSpectraCorrectedK0sOff_20100803.root"); |
4c0d7bc7 | 777 | // hSpectra[kK0][kKaon][kPos] = (TH1F*) AliBWTools::GetdNdPtFromOneOverPt((TH1*) gDirectory->Get("hSpectraOff")); |
778 | hSpectra[kK0][kKaon][kPos] = (TH1F*) gDirectory->Get("hSpectraOff"); | |
779 | // hSpectra[kK0][kKaon][kPos]->Scale(2*TMath::Pi()); | |
780 | // hSpectra[kK0][kKaon][kPos]->Scale(1./272463); | |
781 | hSpectra[kK0][kKaon][kNeg] = hSpectra[kK0][kKaon][kPos]; | |
782 | ||
49ad829f | 783 | // Kinks |
4c0d7bc7 | 784 | // f = new TFile ("./Files/PtAllKaonKinkRap6Apr24.root"); |
49ad829f | 785 | // f = new TFile ("./Files/PtKaonKinkJune13AllPN_20100615.root"); |
22abbb27 | 786 | // f = new TFile ("./Files/KaonKinkJun16PhySel2N.root"); |
787 | f = new TFile ("./Files/KaonKinkJun16PhySel2NRebPass6.root"); | |
4c0d7bc7 | 788 | hSpectra[kKinks][kKaon][kPos] = (TH1F*)gDirectory->Get("fptallKPA"); |
789 | hSpectra[kKinks][kKaon][kNeg] = (TH1F*)gDirectory->Get("fptallKNA"); | |
49ad829f | 790 | // hSpectra[kKinks][kKaon][kPos]->Scale(0.5/0.7); // different rapidity range for kinks |
791 | // hSpectra[kKinks][kKaon][kNeg]->Scale(0.5/0.7); // different rapidity range for kinks | |
792 | // hSpectra[kKinks][kKaon][kPos]->Scale(276004./263345.); // different N of events | |
793 | // hSpectra[kKinks][kKaon][kNeg]->Scale(276004./263345.); // different N of events | |
22abbb27 | 794 | // This normalization if for KaonKinkJun16PhySel2N |
795 | // hSpectra[kKinks][kKaon][kPos]->Scale(1./303214); // different N of events | |
796 | // hSpectra[kKinks][kKaon][kNeg]->Scale(1./303214); // different N of events | |
797 | // This normalization is for KaonKinkJun16PhySel2NRebPass6 | |
798 | hSpectra[kKinks][kKaon][kPos]->Scale(1./260432.26); // different N of events | |
799 | hSpectra[kKinks][kKaon][kNeg]->Scale(1./260432.26); // different N of events | |
4c0d7bc7 | 800 | |
5f2ca1e6 | 801 | // Apply correction factors |
4c0d7bc7 | 802 | // Secondaries for protons |
5f2ca1e6 | 803 | |
0e019ec6 | 804 | //f = new TFile ("./Files/corrFactorProtons_20100615.root"); |
2210b002 | 805 | f = new TFile ("./Files/corrFactorProtons_2010_09_24.root"); |
4c0d7bc7 | 806 | TH1F * hCorrSecondaries = (TH1F*) gDirectory->Get("corrFactorProtons"); |
807 | if(correctSecondaries) { | |
808 | cout << "CORRECTING SECONDARIES" << endl; | |
809 | ||
810 | for(Int_t idet = 0; idet <= kTOF; idet++){ // TPC and TOF only | |
811 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ | |
812 | Int_t ipart = kProton; | |
813 | TH1* h = hSpectra[idet][ipart][icharge]; // lighter notation below | |
814 | if (h){ | |
815 | Int_t nbins = h->GetNbinsX(); | |
816 | for(Int_t ibin = 0; ibin < nbins; ibin++){ | |
817 | // Float_t pt = convertToMT ? TMath::Sqrt(h->GetBinCenter(ibin)*h->GetBinCenter(ibin)-mass[kProton]*mass[kProton]) : h->GetBinCenter(ibin); | |
818 | Float_t pt = h->GetBinCenter(ibin); | |
819 | if (icharge == kNeg) pt = -pt; | |
820 | Int_t binCorrection = hCorrSecondaries->FindBin(pt); | |
821 | Float_t correction = hCorrSecondaries->GetBinContent(binCorrection); | |
822 | // cout << pt << " " << correction << endl; | |
823 | ||
c4fe358b | 824 | if (correction) {// If the bin is empty this is a 0 |
4c0d7bc7 | 825 | h->SetBinContent(ibin,h->GetBinContent(ibin)/correction); |
826 | h->SetBinError (ibin,h->GetBinError (ibin)/correction); | |
827 | } else if (h->GetBinContent(ibin) > 0) { // If we are skipping a non-empty bin, we notify the user | |
828 | cout << "Not correcting bin "<<ibin << " for protons secondaries, " << detFlag[idet] << " " << chargeFlag[icharge] << endl; | |
829 | cout << " Bin content: " << h->GetBinContent(ibin) << endl; | |
830 | ||
831 | } | |
832 | } | |
833 | } | |
834 | } | |
835 | } | |
836 | } | |
837 | // geant/fluka absorption | |
838 | if(correctGeantFlukaXS) { | |
839 | cout << "CORRECTING GEANT3/FLUKA" << endl; | |
eea641de | 840 | // common file for Kaons |
841 | TFile *fFlukakaons = TFile::Open("./Files/correctionForCrossSection.321.root"); | |
842 | TH1D * hCorrFlukakaon[kNCharge]; | |
843 | hCorrFlukakaon[kPos] = (TH1D*)fFlukakaons->Get("gHistCorrectionForCrossSectionParticles"); | |
844 | hCorrFlukakaon[kNeg] = (TH1D*)fFlukakaons->Get("gHistCorrectionForCrossSectionAntiParticles"); | |
845 | ||
846 | for(Int_t idet = 0; idet < kNDet; idet++){ | |
847 | if( idet != kITS) continue; // comment to use fluka for kaons for all dets | |
848 | if (idet == kTOF) continue; // TOF already corrected | |
849 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ | |
850 | Int_t ipart = kKaon; | |
2f7f947e | 851 | TH1 * h = hSpectra[idet][ipart][icharge]; // only ITS sa |
eea641de | 852 | if (h){ |
853 | Int_t nbins = h->GetNbinsX(); | |
854 | Int_t nbinsy=hCorrFlukakaon[icharge]->GetNbinsY(); | |
855 | for(Int_t ibin = 0; ibin < nbins; ibin++){ | |
856 | Float_t pt = h->GetBinCenter(ibin); | |
857 | Float_t minPtCorrection = hCorrFlukakaon[icharge]->GetXaxis()->GetBinLowEdge(1); | |
858 | Float_t maxPtCorrection = hCorrFlukakaon[icharge]->GetXaxis()->GetBinLowEdge(nbinsy+1); | |
859 | if (pt < minPtCorrection) pt = minPtCorrection+0.0001; | |
860 | if (pt > maxPtCorrection) pt = maxPtCorrection; | |
861 | Float_t correction = hCorrFlukakaon[icharge]->GetBinContent(hCorrFlukakaon[icharge]->GetXaxis()->FindBin(pt)); | |
862 | if (correction != 0) {// If the bin is empty this is a 0 | |
863 | h->SetBinContent(ibin,h->GetBinContent(ibin)*correction); | |
864 | h->SetBinError (ibin,h->GetBinError (ibin)*correction); | |
865 | } else if (h->GetBinContent(ibin) > 0) { // If we are skipping a non-empty bin, we notify the user | |
866 | cout << "Fluka/GEANT: Not correcting bin "<<ibin << " for protons secondaries, ITS, " << chargeFlag[icharge] << endl; | |
867 | cout << " Bin content: " << h->GetBinContent(ibin) << endl; | |
868 | } | |
869 | } | |
870 | } | |
871 | } | |
872 | } | |
873 | ||
2f7f947e | 874 | // PROTONS |
875 | ||
eea641de | 876 | // ITS specific file for protons/antiprotons |
5f2ca1e6 | 877 | TFile* fITS = new TFile ("./Files/correctionForCrossSectionITS_20100719.root"); |
878 | TH2D * hCorrFlukaITS[kNCharge]; | |
879 | hCorrFlukaITS[kPos] = (TH2D*)fITS->Get("gHistCorrectionForCrossSectionProtons"); | |
880 | hCorrFlukaITS[kNeg] = (TH2D*)fITS->Get("gHistCorrectionForCrossSectionAntiProtons"); | |
881 | ||
882 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ | |
883 | Int_t ipart = kProton; | |
884 | TH1 * h = hSpectra[kITS][ipart][icharge]; // only ITS sa | |
885 | if (h){ | |
886 | Int_t nbins = h->GetNbinsX(); | |
887 | Int_t nbinsy=hCorrFlukaITS[icharge]->GetNbinsY(); | |
888 | for(Int_t ibin = 0; ibin < nbins; ibin++){ | |
889 | Float_t pt = h->GetBinCenter(ibin); | |
890 | Float_t minPtCorrection = hCorrFlukaITS[icharge]->GetYaxis()->GetBinLowEdge(1); | |
891 | Float_t maxPtCorrection = hCorrFlukaITS[icharge]->GetYaxis()->GetBinLowEdge(nbinsy+1); | |
892 | if (pt < minPtCorrection) pt = minPtCorrection+0.0001; | |
893 | if (pt > maxPtCorrection) pt = maxPtCorrection; | |
894 | Float_t correction = hCorrFlukaITS[icharge]->GetBinContent(1,hCorrFlukaITS[icharge]->GetYaxis()->FindBin(pt)); | |
895 | if (correction != 0) {// If the bin is empty this is a 0 | |
896 | h->SetBinContent(ibin,h->GetBinContent(ibin)*correction); | |
897 | h->SetBinError (ibin,h->GetBinError (ibin)*correction); | |
898 | } else if (h->GetBinContent(ibin) > 0) { // If we are skipping a non-empty bin, we notify the user | |
899 | cout << "Fluka/GEANT: Not correcting bin "<<ibin << " for protons secondaries, ITS, " << chargeFlag[icharge] << endl; | |
900 | cout << " Bin content: " << h->GetBinContent(ibin) << endl; | |
901 | } | |
902 | } | |
903 | } | |
904 | } | |
905 | ||
4c0d7bc7 | 906 | |
77233d47 | 907 | |
908 | //f = new TFile ("./Files/correctionForCrossSection_20100615.root"); | |
909 | f = new TFile ("./Files/correctionForCrossSection_20100920.root"); | |
4c0d7bc7 | 910 | TH2D * hCorrFluka[kNCharge]; |
911 | hCorrFluka[kPos] = (TH2D*)gDirectory->Get("gHistCorrectionForCrossSectionProtons"); | |
912 | hCorrFluka[kNeg] = (TH2D*)gDirectory->Get("gHistCorrectionForCrossSectionAntiProtons"); | |
913 | for(Int_t idet = 0; idet < kNDet; idet++){ | |
914 | if (idet == kITS) continue; // skip ITS sa | |
915 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ | |
916 | Int_t ipart = kProton; | |
917 | TH1 * h = hSpectra[idet][ipart][icharge]; // lighter notation below | |
918 | if (h){ | |
919 | Int_t nbins = h->GetNbinsX(); | |
920 | for(Int_t ibin = 0; ibin < nbins; ibin++){ | |
921 | // Float_t pt = convertToMT ? | |
922 | // TMath::Sqrt(h->GetBinCenter(ibin)*h->GetBinCenter(ibin)-mass[kProton]*mass[kProton]) : | |
923 | // h->GetBinCenter(ibin); | |
924 | Float_t pt = h->GetBinCenter(ibin); | |
925 | Float_t minPtCorrection = hCorrFluka[icharge]->GetYaxis()->GetBinLowEdge(1); | |
2210b002 | 926 | // hCorrFluka[icharge]->Draw(); |
927 | TH1D * htmpFluka=hCorrFluka[icharge]->ProjectionY(); | |
928 | Float_t maxPtCorrection = AliBWTools::GetHighestNotEmptyBinEdge(htmpFluka);//->GetYaxis()->GetBinLowEdge(hCorrFluka[icharge]->GetNbinsY()+1); | |
4c0d7bc7 | 929 | if (pt < minPtCorrection) pt = minPtCorrection+0.0001; |
2210b002 | 930 | if (pt > maxPtCorrection) pt = maxPtCorrection-0.0001; |
4c0d7bc7 | 931 | Float_t correction = hCorrFluka[icharge]->GetBinContent(1,hCorrFluka[icharge]->GetYaxis()->FindBin(pt)); |
77233d47 | 932 | |
5f2ca1e6 | 933 | // already in the efficiency correction (F. Noferini) |
2f7f947e | 934 | if (idet == kTOF) { |
49ad829f | 935 | if (icharge == kNeg) correction = 1; // antiprotons already corrected in efficiency |
936 | // Scale correction for the different material budget. Recipe by Francesco Noferini | |
2210b002 | 937 | else { |
938 | correction = TMath::Power(correction,0.07162/0.03471); | |
939 | } | |
49ad829f | 940 | } |
4c0d7bc7 | 941 | if (correction != 0) {// If the bin is empty this is a 0 |
77233d47 | 942 | |
4c0d7bc7 | 943 | h->SetBinContent(ibin,h->GetBinContent(ibin)*correction); |
944 | h->SetBinError (ibin,h->GetBinError (ibin)*correction); | |
4c0d7bc7 | 945 | } else if (h->GetBinContent(ibin) > 0) { // If we are skipping a non-empty bin, we notify the user |
946 | cout << "Fluka/GEANT: Not correcting bin "<<ibin << " for protons secondaries, " << detFlag[idet] << " " << chargeFlag[icharge] << endl; | |
947 | cout << " Bin content: " << h->GetBinContent(ibin) << endl; | |
948 | } | |
949 | ||
950 | } | |
951 | ||
952 | } | |
953 | } | |
954 | } | |
955 | } | |
956 | ||
957 | ||
958 | // Set style and scale | |
959 | for(Int_t idet = 0; idet < kNDet; idet++){ | |
960 | for(Int_t ipart = 0; ipart < kNPart; ipart++){ | |
961 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ | |
962 | if (hSpectra[idet][ipart][icharge]){ | |
963 | hSpectra[idet][ipart][icharge]->SetStats(0); // disable stats | |
964 | hSpectra[idet][ipart][icharge]->SetMarkerColor (color[ipart] ); | |
965 | hSpectra[idet][ipart][icharge]->SetLineColor (color[ipart] ); | |
966 | hSpectra[idet][ipart][icharge]->SetMarkerStyle (marker[ipart]); | |
967 | hSpectra[idet][ipart][icharge]->SetXTitle("p_{t} (GeV/c)"); | |
968 | hSpectra[idet][ipart][icharge]->SetYTitle("1/N_{ev} dN/dp_{t} (GeV/c)^{-1}"); | |
969 | if (convertToMT) { | |
970 | TH1F * htmp = (TH1F*) AliBWTools::GetOneOverPtdNdPt(hSpectra[idet][ipart][icharge]); | |
971 | hSpectra[idet][ipart][icharge] = (TH1F*)AliBWTools::GetdNdmtFromdNdpt(htmp,mass[ipart]); | |
972 | hSpectra[idet][ipart][icharge]->SetXTitle("m_{t} (GeV)"); | |
973 | hSpectra[idet][ipart][icharge]->SetYTitle("1/N_{ev} 1/m_{t} dN/dp_{t} (GeV)^{-1}"); | |
974 | } | |
975 | // if (idet == kTOF || idet == kTPC) { | |
976 | // hSpectra[idet][ipart][icharge]->Scale(1./236763); | |
977 | // } | |
978 | // if (idet == kITS){ | |
979 | // hSpectra[idet][ipart][icharge]->Scale(202945./236763); | |
980 | // } | |
981 | if (scaleKaons && ipart == kKaon) { | |
982 | hSpectra[idet][ipart][icharge]->Scale(2.); | |
983 | } | |
984 | } else { | |
985 | printf("Cannot load %s,%s,%s\n",detFlag[idet],partFlag[ipart],chargeFlag[icharge]); | |
986 | } | |
987 | } | |
988 | } | |
989 | } | |
990 | ||
ee08b77d | 991 | |
992 | // Create fake weights for the mean; To be update once we have syst errors | |
2f7f947e | 993 | // Using syste from table in paper. It would be better to have this as a function of pt. |
994 | TH1F * hWeights[3][kNPart]; | |
995 | const Double_t kWeights[3][kNPart] = | |
3725266c | 996 | // {{4, 3, 10.2}, // TPC |
997 | // {4.1,8.8,7.0 }, //TOF | |
998 | // {4.5,5.6,7.0 }}; // ITS | |
999 | {{0.1,0.1,0.1}, // TPC | |
1000 | {0.1,0.1,0.1}, //TOF | |
1001 | {0.1,0.1,0.1}}; // ITS | |
2f7f947e | 1002 | for(Int_t ipart = 0; ipart <= kNPart ; ipart++){ |
1003 | for(Int_t idet = 0; idet <= kITS ; idet++){ | |
1004 | hWeights[idet][ipart] = (TH1F*) hSpectra[idet][ipart][kPos]->Clone(); | |
2210b002 | 1005 | Int_t nbinW = hWeights[idet][ipart]->GetNbinsX(); |
1006 | for(Int_t ibin = 1; ibin <= nbinW; ibin++){ | |
2f7f947e | 1007 | hWeights[idet][ipart]->SetBinError(ibin, kWeights[idet][ipart]); |
1008 | } | |
1009 | } | |
ee08b77d | 1010 | } |
49ad829f | 1011 | const Double_t scaleDet[] = {1.00,1.00,1.00}; // During the weekly meeting 19/08/2010 we decided not to apply any scaling. |
1012 | // const Double_t scaleDet[] = {0.98,1,1.00}; // Scaling factor for the different detectors. Estimated from ratios, it has an estimated uncertainty of ~2% | |
2f7f947e | 1013 | // const Double_t scaleDet[] = {0.88,1,0.88}; // Scaling factor for the different detectors. Estimated from ratios, it has an estimated uncertainty of ~2% |
ee08b77d | 1014 | |
1015 | ||
4c0d7bc7 | 1016 | // Combine detectors |
1017 | for(Int_t ipart = 0; ipart < kNPart; ipart++){ | |
1018 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ | |
1019 | TH1F * htemplLocal = htemplate; // If we are converting to 1/mt dNdmt we need to convert the template as well... | |
1020 | ||
1021 | if (convertToMT) { | |
1022 | TH1F * htmp = (TH1F*) AliBWTools::GetOneOverPtdNdPt(htemplate); | |
1023 | htemplLocal = (TH1F*)AliBWTools::GetdNdmtFromdNdpt(htmp,mass[ipart]); | |
1024 | ||
1025 | } | |
1026 | hSpectra[kCombTOFTPC][ipart][icharge] = AliBWTools::CombineHistos(hSpectra[kTOF][ipart][icharge], | |
1027 | hSpectra[kTPC][ipart][icharge], | |
1028 | htemplLocal,1.);; | |
3725266c | 1029 | |
ee08b77d | 1030 | hSpectra[kCombAll][ipart][icharge] = |
3725266c | 1031 | AliBWTools::Combine3HistosWithErrors(hSpectra[kITS][ipart][icharge], // Histos to combine |
ee08b77d | 1032 | hSpectra[kTPC][ipart][icharge], |
1033 | hSpectra[kTOF][ipart][icharge], | |
3725266c | 1034 | hSystError[kITS][ipart][icharge], // Errors (weights) used for the average |
1035 | hSystError[kTPC][ipart][icharge], | |
1036 | hSystError[kTOF][ipart][icharge], | |
1037 | // hWeights[kITS][ipart], | |
1038 | // hWeights[kTPC][ipart], | |
1039 | // hWeights[kTOF][ipart], | |
1040 | htemplLocal,1, // template, take statistical error from TPC in overlap | |
1041 | scaleDet[kITS], // relative scaling | |
2f7f947e | 1042 | scaleDet[kTPC], |
3725266c | 1043 | scaleDet[kTOF], |
1044 | (TH1**)&hSystError[kCombAll][ipart][icharge], // put combined syst error here | |
1045 | 1 // weights histos contain error in bin content | |
2f7f947e | 1046 | ); |
4c0d7bc7 | 1047 | // if (convertToMT) { |
1048 | // TH1F * htmp = (TH1F*) AliBWTools::GetOneOverPtdNdPt(hSpectra[kCombTOFTPC][ipart][icharge]); | |
1049 | // hSpectra[kCombTOFTPC][ipart][icharge] = (TH1F*)AliBWTools::GetdNdmtFromdNdpt(htmp,mass[ipart]); | |
1050 | // hSpectra[kCombTOFTPC][ipart][icharge]->SetXTitle("m_{t} (GeV)"); | |
1051 | // hSpectra[kCombTOFTPC][ipart][icharge]->SetYTitle("1/N_{ev} 1/m_{t} dN/dp_{t} (GeV)^{-1}"); | |
1052 | // } | |
1053 | } | |
1054 | } | |
1055 | ||
1056 | ||
1057 | // Scale for the number of inelastic collisions and correct for | |
1058 | // efficiency losses due to physics selection: | |
1059 | ||
1060 | Double_t effPhysSel[kNPart]; | |
1061 | effPhysSel[kPion] = 1.012; | |
1062 | effPhysSel[kKaon] = 1.013; | |
1063 | effPhysSel[kProton] = 1.014; | |
1064 | ||
1065 | ||
1066 | for(Int_t idet = 0; idet < kNHist; idet++){ | |
1067 | for(Int_t ipart = 0; ipart < kNPart; ipart++){ | |
1068 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ | |
1069 | if(hSpectra[idet][ipart][icharge]) { | |
1070 | // cout << "Scaling!" << endl; | |
49ad829f | 1071 | if(idet != kKinks && idet != kK0){ // Kinks use a different run list, k0s normalized by Helene |
57b5143c | 1072 | if(idet!=kITSTPC) hSpectra[idet][ipart][icharge]->Scale(1.*effPhysSel[ipart]/278366.15); // Scale PhysSel efficiency (computed by Alexander) |
1073 | else hSpectra[idet][ipart][icharge]->Scale(1./278366.15); // Scale PhysSel efficiency (computed by Alexander) | |
5f2ca1e6 | 1074 | } |
4c0d7bc7 | 1075 | } |
1076 | } | |
1077 | } | |
1078 | } | |
1079 | ||
1080 | ||
1081 | } | |
1082 | ||
1083 | void LoadMC() { | |
1084 | ||
1085 | TFile * f = 0; | |
1086 | const char * evClass= "INEL"; | |
1087 | const char * files[] = {"./Files/dndeta_Phojet_10M_900GeV.root", | |
1088 | "./Files/dndeta_AtlasCSC306_10M_900GeV.root", | |
1089 | "./Files/dndeta_CMSD6T109_10M_900GeV.root", | |
1090 | "./Files/dndeta_Perugia0320_10M_900GeV.root", }; | |
1091 | ||
1092 | // Phojet | |
1093 | for(Int_t itune = 0; itune < kNTunes; itune++){ | |
1094 | f = new TFile(files[itune]); | |
1095 | for(Int_t ipart = 0; ipart < kNPart; ipart++){ | |
1096 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ | |
1097 | hSpectraMC[itune][ipart][icharge] = (TH1F*) f->Get(Form("fHistPtID_%s_%s%s",evClass,partFlag[ipart],icharge==0 ? "Pos" : "Neg")); | |
1098 | } | |
1099 | } | |
1100 | } | |
1101 | ||
1102 | ||
1103 | // Set style | |
1104 | for(Int_t itune = 0; itune < kNTunes; itune++){ | |
1105 | for(Int_t ipart = 0; ipart < kNPart; ipart++){ | |
1106 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ | |
1107 | if (hSpectraMC[itune][ipart][icharge]){ | |
1108 | hSpectraMC[itune][ipart][icharge]->SetName(TString(hSpectraMC[itune][ipart][icharge]->GetName())+mcTuneName[itune]); | |
1109 | hSpectraMC[itune][ipart][icharge]->SetMarkerColor (mcLineColor[itune] ); | |
1110 | hSpectraMC[itune][ipart][icharge]->SetLineColor (mcLineColor[itune] ); | |
1111 | hSpectraMC[itune][ipart][icharge]->SetLineStyle (mcLineStyle[itune] ); | |
1112 | hSpectraMC[itune][ipart][icharge]->SetMarkerStyle (1); | |
1113 | if (convertToMT) { | |
1114 | TH1F * htmp = (TH1F*)AliBWTools::GetOneOverPtdNdPt(hSpectraMC[itune][ipart][icharge]); | |
1115 | hSpectraMC[itune][ipart][icharge] = (TH1F*)AliBWTools::GetdNdmtFromdNdpt(htmp,mass[ipart]); | |
1116 | hSpectraMC[itune][ipart][icharge]->SetXTitle("m_{t} (GeV)"); | |
1117 | hSpectraMC[itune][ipart][icharge]->SetYTitle("1/N_{ev} 1/m_{t} dN/dp_{t} (GeV)^{-1}"); | |
1118 | } | |
1119 | ||
1120 | } else { | |
1121 | printf("Cannot load MC # %d,%s,%s\n",itune,partFlag[ipart],chargeFlag[icharge]); | |
1122 | } | |
1123 | } | |
1124 | } | |
1125 | } | |
1126 | ||
1127 | } | |
1128 | ||
1129 | void DrawStar(Int_t icharge) { | |
1130 | // cout << icharge << endl; | |
1131 | ||
1132 | // gROOT->LoadMacro("StarPPSpectra.C"); | |
1133 | TGraphErrors ** gStar = StarPPSpectra(); | |
1134 | ||
1135 | for(Int_t istar = 0; istar < 6; istar++){ | |
1136 | gStar[istar]->SetMarkerStyle(kOpenStar); | |
1137 | if (icharge==kPos && (istar%2) ) gStar[istar]->Draw("P"); | |
1138 | else if (icharge==kNeg && !(istar%2) ) gStar[istar]->Draw("P"); | |
1139 | else cout << "Skipping Star " << istar << endl; | |
1140 | } | |
1141 | } | |
1142 | ||
1143 | void GetITSResiduals() { | |
1144 | ||
1145 | ||
1146 | for(Int_t ipart = 0; ipart < kNPart; ipart++){ | |
1147 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ | |
1148 | // cout << "1 " << ipart << " " << icharge << endl; | |
1149 | ||
1150 | // gSystem->ProcessEvents(); | |
1151 | // gSystem->Sleep(1000); | |
1152 | TF1* f = (TF1*) hSpectra[kCombTOFTPC][ipart][icharge]->GetListOfFunctions()->At(0); | |
1153 | TH1F * hres1, *hres2; | |
1154 | AliBWTools::GetResiduals(hSpectra[kITS][ipart][icharge], f, &hres1, &hres2); | |
1155 | ||
1156 | TCanvas * c1 = new TCanvas(TString(partFlag[ipart])+"_"+chargeFlag[icharge],TString(partFlag[ipart])+"_"+chargeFlag[icharge]); | |
1157 | c1->SetLogy(); | |
1158 | hSpectra[kCombTOFTPC][ipart][icharge]->Draw(); | |
1159 | hSpectra[kITS][ipart][icharge]->SetMarkerStyle(24); | |
1160 | hSpectra[kCombTOFTPC][ipart][icharge]->SetMarkerStyle(20); | |
1161 | hSpectra[kITS][ipart][icharge]->Draw("same"); | |
1162 | hSpectra[kCombTOFTPC][ipart][icharge]->GetListOfFunctions()->At(0)->Draw("same"); | |
1163 | TLegend* l = new TLegend( 0.182886, 0.192308, 0.505034,0.384615, TString(partLabel[ipart][icharge])+" "+chargeFlag[icharge]); | |
1164 | l->AddEntry(hSpectra[kCombTOFTPC][ipart][icharge], "TOF+TPC"); | |
1165 | l->AddEntry(hSpectra[kITS][ipart][icharge], "ITS"); | |
1166 | l->AddEntry(f, "Fit to TOF+TPC"); | |
1167 | l->Draw(); | |
1168 | ||
1169 | TCanvas * c2 = new TCanvas(TString(partFlag[ipart])+"_"+chargeFlag[icharge]+"_res", | |
1170 | TString(partFlag[ipart])+"_"+chargeFlag[icharge]+"_res"); | |
1171 | c2->SetGridy(); | |
1172 | hres2->SetMinimum(-1); | |
1173 | hres2->SetMaximum(1); | |
1174 | hres2->Draw(); | |
1175 | hres2->GetYaxis()->SetTitleOffset(1.2); | |
1176 | Float_t x = AliBWTools::GetLowestNotEmptyBinEdge(hSpectra[kCombTOFTPC][ipart][icharge]); | |
1177 | TLine * line = new TLine(x,-1,x,1); | |
1178 | line->SetLineStyle(kDashed); | |
1179 | line->Draw("same"); | |
1180 | ||
1181 | if (doPrint) { | |
1182 | c1->Update(); | |
1183 | c2->Update(); | |
1184 | gSystem->ProcessEvents(); | |
1185 | c1->Print(TString(c1->GetName()) + ".png"); | |
1186 | c2->Print(TString(c2->GetName()) + ".png"); | |
1187 | } | |
1188 | } | |
1189 | } | |
1190 | } | |
1191 | ||
1192 | void DrawWithModels() { | |
1193 | ||
4d4a268c | 1194 | Int_t chargeLoop = sumCharge ? 1 : 2; |
1195 | for(Int_t icharge = 0; icharge < chargeLoop; icharge++){ | |
4c0d7bc7 | 1196 | |
1197 | // Draw with models | |
1198 | for(Int_t ipart = 0; ipart < kNPart; ipart++){ | |
1199 | // Multipad canvas | |
4d4a268c | 1200 | TString chargeFlagLocal = chargeFlag[icharge]; |
1201 | if(sumCharge) chargeFlagLocal += chargeFlag[1]; | |
1202 | TCanvas * c1 = new TCanvas(TString("cSpectra")+partFlag[ipart]+chargeFlagLocal,TString("cSpectra")+partFlag[ipart]+chargeFlagLocal,700,700); | |
1203 | TPad *p1 = new TPad(TString("p1")+partFlag[ipart]+chargeFlagLocal, "p1", 0.0, 0.35, 1.0, 0.95, 0, 0, 0); | |
4c0d7bc7 | 1204 | p1->SetBottomMargin(0); |
1205 | p1->Draw(); | |
1206 | ||
4d4a268c | 1207 | TPad *p2 = new TPad(TString("p2")+partFlag[ipart]+chargeFlagLocal, "p2", 0.0, 0.05, 1.0, 0.35, 0, 0, 0); |
4c0d7bc7 | 1208 | p2->SetTopMargin(0); |
5f2ca1e6 | 1209 | p2->SetBottomMargin(0.3); |
4c0d7bc7 | 1210 | p2->Draw(); |
1211 | ||
1212 | Float_t scaleFonts = (0.95-0.3)/(0.3-0.05); | |
1213 | ||
1214 | // Draw spectra | |
1215 | p1->cd(); | |
1216 | p1->SetLogy(); | |
4d4a268c | 1217 | Float_t maxy = sumCharge ? 10 : 5; |
1218 | TH2F * hempty = new TH2F(TString("hempty")+long(icharge)+long(ipart),"hempty",100,0.,3, 100, 0.0015,maxy); | |
4c0d7bc7 | 1219 | hempty->SetXTitle("p_{t} (GeV/c)"); |
1220 | hempty->SetYTitle("1/N_{ev} d^{2}N / dydp_{t} (GeV/c)^{-1}"); | |
1221 | hempty->Draw(); | |
1222 | c1->SetLogy(); | |
d8e6677d | 1223 | |
4c0d7bc7 | 1224 | |
1225 | ||
d8e6677d | 1226 | TLegend * l =new TLegend(0.502874, 0.493056, 0.892241, 0.904762); |
4d4a268c | 1227 | l->SetFillColor(kWhite); |
1228 | if(sumCharge) { | |
1229 | hSpectra[iCombInStudy][ipart][kPos]->Add(hSpectra[iCombInStudy][ipart][kNeg]); // Draw pos+neg | |
1230 | } | |
ee08b77d | 1231 | hSpectra[iCombInStudy][ipart][icharge]->Draw("same"); |
d8e6677d | 1232 | l->AddEntry(hSpectra[kTOF][ipart][icharge],TString ("Data")); |
4c0d7bc7 | 1233 | for(Int_t itune = 0; itune < kNTunes; itune++){ |
4d4a268c | 1234 | if(sumCharge) hSpectraMC[itune][ipart][icharge]->Add(hSpectraMC[itune][ipart][kNeg]); // Draw pos+neg; |
4c0d7bc7 | 1235 | l->AddEntry(hSpectraMC[itune][ipart][icharge],mcTuneName[itune]); |
1236 | hSpectraMC[itune][ipart][icharge]->SetLineWidth(2); | |
d8e6677d | 1237 | AliBWTools::GetGraphFromHisto(hSpectraMC[itune][ipart][icharge])->Draw("CX"); |
4c0d7bc7 | 1238 | } |
1239 | l->Draw("same"); | |
d8e6677d | 1240 | |
4d4a268c | 1241 | TLatex * tex = new TLatex(0.6712643,2.353486,sumCharge ? Form("%s+%s",partLabel[ipart][icharge],partLabel[ipart][1]) : partLabel[ipart][icharge]); |
d8e6677d | 1242 | tex->SetTextFont(42); |
1243 | tex->SetTextSize(0.07936508); | |
1244 | tex->SetLineWidth(2); | |
1245 | tex->Draw(); | |
4c0d7bc7 | 1246 | |
1247 | // Draw ratio | |
1248 | p2->cd(); | |
4d4a268c | 1249 | Float_t maxyr = sumCharge ? 3.99 : 2.99; |
1250 | TH2F * hemptyr = new TH2F(TString("hemptyratio")+long(icharge)+long(ipart),"hempty",100,0.,3, 100, 0.01,maxyr); | |
4c0d7bc7 | 1251 | hemptyr->SetXTitle("p_{t} (GeV/c)"); |
1252 | hemptyr->SetYTitle("Data/MC"); | |
1253 | hemptyr->GetXaxis()->SetLabelSize(0.04*scaleFonts); | |
1254 | hemptyr->GetYaxis()->SetLabelSize(0.04*scaleFonts); | |
1255 | hemptyr->GetYaxis()->SetTitleSize(0.05*scaleFonts); | |
1256 | hemptyr->GetYaxis()->SetTitleOffset(1.4/scaleFonts); | |
1257 | hemptyr->GetXaxis()->SetTitleSize(0.05*scaleFonts); | |
d8e6677d | 1258 | hemptyr->GetXaxis()->SetTitleOffset(1.05); |
4c0d7bc7 | 1259 | hemptyr->SetTickLength(0.03*scaleFonts, "X"); |
1260 | hemptyr->SetTickLength(0.02*scaleFonts, "Y"); | |
1261 | // hemptyr->GetXaxis()->SetTitleOffset(1.4/scaleFonts); | |
1262 | hemptyr->GetYaxis()->SetNdivisions(5); | |
1263 | hemptyr->Draw(""); | |
1264 | ||
1265 | AliBWFunc fm; | |
1266 | for(Int_t itune = 0; itune < kNTunes; itune++){ | |
1267 | TF1* f = fm.GetHistoFunc(hSpectraMC[itune][ipart][icharge], TString("f")+mcTuneName[itune]); | |
1268 | ||
1269 | // l->AddEntry(hSpectraMC[itune][ipart][icharge],mcTuneName[itune]); | |
ee08b77d | 1270 | TH1F* hRatio = AliBWTools::DivideHistoByFunc(hSpectra[iCombInStudy][ipart][icharge],f); |
4c0d7bc7 | 1271 | hRatio->SetLineStyle(hSpectraMC[itune][ipart][icharge]->GetLineStyle()); |
1272 | hRatio->SetLineColor(hSpectraMC[itune][ipart][icharge]->GetLineColor()); | |
1273 | hRatio->SetLineWidth(hSpectraMC[itune][ipart][icharge]->GetLineWidth()); | |
4d4a268c | 1274 | hRatio->SetMarkerStyle(0); |
1275 | hRatio->Draw("same"); | |
1276 | // AliBWTools::GetGraphFromHisto(hRatio)->Draw("CX"); | |
4c0d7bc7 | 1277 | } |
1278 | ||
1279 | ||
1280 | ||
4d4a268c | 1281 | PrintCanvas(c1,printFormats); |
4c0d7bc7 | 1282 | } |
1283 | } | |
1284 | ||
1285 | ||
1286 | ||
1287 | } | |
1288 | ||
1289 | void DrawAllAndKaons() { | |
1290 | ||
22abbb27 | 1291 | |
4c0d7bc7 | 1292 | |
1293 | // gROOT->LoadMacro("ALICEWorkInProgress.C"); | |
1294 | ||
eea641de | 1295 | // gStyle->SetOptFit(0); |
1296 | gStyle->SetStatX(0.9); | |
1297 | gStyle->SetStatY(0.88); | |
1298 | gStyle->SetStatW(0.4); | |
1299 | gStyle->SetStatH(0.1); | |
4c0d7bc7 | 1300 | |
ee08b77d | 1301 | TH1F * hKaonsAllTPCTOF = (TH1F*) hSpectra[iCombInStudy][kKaon][kPos]->Clone(); |
1302 | hKaonsAllTPCTOF->Add(hSpectra[iCombInStudy][kKaon][kNeg]); | |
4c0d7bc7 | 1303 | |
1304 | TH1F * hK0Scaled = (TH1F*) hSpectra[kK0][kKaon][kPos]->Clone(); | |
1305 | hK0Scaled->Add(hSpectra[kK0][kKaon][kPos]); | |
1306 | ||
1307 | hSpectra[kKinks][kKaon][kPos]->SetMarkerStyle(25); | |
5f2ca1e6 | 1308 | hSpectra[kKinks][kKaon][kPos]->SetLineColor(4); |
4c0d7bc7 | 1309 | hSpectra[kKinks][kKaon][kPos]->SetStats(0); |
1310 | TH1F * hKinksAll = (TH1F*) hSpectra[kKinks][kKaon][kPos]->Clone(); | |
1311 | hKinksAll->Add(hSpectra[kKinks][kKaon][kNeg]); | |
1312 | ||
5f2ca1e6 | 1313 | TCanvas * c1 = new TCanvas("cKaons","cKaons",700,700); |
4c0d7bc7 | 1314 | c1->SetLogy(); |
5f2ca1e6 | 1315 | TH2F * hempty = new TH2F("hempty_allkaons","hempty",100,0.,3, 100, 1e-3,6); |
4c0d7bc7 | 1316 | hempty->SetXTitle("p_{t} (GeV/c)"); |
22abbb27 | 1317 | // hempty->SetYTitle("dN / dp_{t} (A.U.)"); |
1318 | hempty->SetYTitle("1/N_{ev} d^{2}N / dydp_{t} (GeV/c)^{-1}"); | |
4c0d7bc7 | 1319 | hempty->Draw(); |
4c0d7bc7 | 1320 | hK0Scaled->Draw("same"); |
5f2ca1e6 | 1321 | hKaonsAllTPCTOF->Draw("same"); |
4c0d7bc7 | 1322 | hKinksAll->Draw("same"); |
22abbb27 | 1323 | |
4c0d7bc7 | 1324 | |
22abbb27 | 1325 | TLegend * leg = new TLegend(0.456897, 0.71875, 0.892241, 0.936012,NULL,"brNDC"); |
5f2ca1e6 | 1326 | // leg->SetBorderSize(0); |
4c0d7bc7 | 1327 | // leg->SetLineColor(1); |
1328 | // leg->SetLineStyle(1); | |
1329 | // leg->SetLineWidth(1); | |
1330 | // leg->SetFillColor(19); | |
5f2ca1e6 | 1331 | leg->SetFillColor(0); |
1332 | TLegendEntry *entry=leg->AddEntry(hKaonsAllTPCTOF,"K^{+} + K^{-}, ITS+TPC+TOF ","lpf"); | |
1333 | entry=leg->AddEntry(hK0Scaled,"K^{0} #times 2","lpf"); | |
1334 | entry=leg->AddEntry(hKinksAll,"K^{+} + K ^{-}, Kinks","lpf"); | |
4c0d7bc7 | 1335 | leg->Draw(); |
1336 | ||
22abbb27 | 1337 | // ALICEWorkInProgress(c1,today.Data(),"#splitline{ALICE Prelimiary}{Statistical Error Only}"); |
49ad829f | 1338 | // TLatex * tex = new TLatex(0.2120805,0.01288336,"Statistical error only"); |
1339 | // tex->SetTextColor(2); | |
1340 | // tex->SetTextFont(42); | |
1341 | // tex->SetTextSize(0.03496503); | |
1342 | // tex->Draw(); | |
4c0d7bc7 | 1343 | |
1344 | c1->Update(); | |
22abbb27 | 1345 | PrintCanvas(c1, printFormats); |
4c0d7bc7 | 1346 | |
1347 | // Draw all "stable" hadrons | |
57b5143c | 1348 | Bool_t divideCanvas=kFALSE; |
1349 | TH1F * hSpectraSystError[kNHist][kNPart][kNCharge];// used to store spectra with syst error only | |
1350 | ||
4c0d7bc7 | 1351 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ |
78033eeb | 1352 | TCanvas * c1h = 0; |
1353 | if(divideCanvas) c1h = new TCanvas(TString("cAll_")+chargeFlag[icharge],TString("cAll_")+chargeFlag[icharge],1200,500); | |
57b5143c | 1354 | else { |
1355 | c1h = new TCanvas(TString("cAll_")+chargeFlag[icharge],TString("cAll_")+chargeFlag[icharge],700,700); | |
1356 | c1h->SetLogy(); | |
1357 | } | |
2210b002 | 1358 | c1h->SetLeftMargin(0.14); |
78033eeb | 1359 | if(divideCanvas) { |
1360 | c1h->Divide(3,1); | |
1361 | c1h->cd(1); | |
1362 | } | |
57b5143c | 1363 | TH2F * hemptyLoc = new TH2F(TString("hempty")+long(icharge),"hempty",100,0.,4, 100, 1e-2,5); |
2210b002 | 1364 | hemptyLoc->SetXTitle("p_{t} (GeV/c)"); |
1365 | hemptyLoc->SetYTitle("1/N_{ev} d^{2}N / dydp_{t} (GeV/c)^{-1}"); | |
1366 | hemptyLoc->GetYaxis()->SetTitleOffset(1.35); | |
1367 | hemptyLoc->GetXaxis()->SetTitleOffset(1.1); | |
57b5143c | 1368 | hemptyLoc->GetXaxis()->SetRangeUser(0,1); |
2210b002 | 1369 | hemptyLoc->Draw(); |
57b5143c | 1370 | |
1371 | for(Int_t ipart = 0; ipart < kNPart; ipart++) { | |
1372 | if(ipart == kPion) leg = new TLegend(0.558908, 0.78125, 0.889368, 0.924107, NULL,"brNDC"); | |
1373 | else if (ipart == kKaon) leg = new TLegend(0.16523, 0.177083, 0.50431, 0.31994, NULL,"brNDC"); | |
1374 | else if (ipart == kProton) leg = new TLegend(0.51, 0.177083, 0.84, 0.31994, NULL,"brNDC"); | |
1375 | leg->SetFillColor(0); | |
78033eeb | 1376 | if(divideCanvas) { |
1377 | cout << "CD!!" << ipart+1 << endl; | |
1378 | c1h->cd(ipart+1); | |
1379 | TH2F * hClone = new TH2F(*hemptyLoc); | |
1380 | hClone->GetXaxis()->SetRangeUser(0,1); | |
1381 | if(ipart == kPion) hClone->GetYaxis()->SetRangeUser(0,3.5); | |
1382 | if(ipart == kKaon) hClone->GetYaxis()->SetRangeUser(0,0.28); | |
1383 | if(ipart == kProton) hClone->GetYaxis()->SetRangeUser(0,0.1); | |
1384 | hClone->Draw(); | |
1385 | } | |
4c0d7bc7 | 1386 | for(Int_t idet = 0; idet <= kITSTPC; idet++){ |
1387 | // if (idet == kITS) continue; | |
1388 | // if (idet == kITSTPC) hSpectra[idet][ipart][icharge]->SetMarkerColor(kGreen); | |
1389 | hSpectra[idet][ipart][icharge]->SetMarkerStyle(marker[idet]); | |
57b5143c | 1390 | TH1F * hsyst = new TH1F(*hSpectra[idet][ipart][icharge]); |
1391 | hsyst->Reset(); | |
1392 | AliBWTools::GetValueAndError(hsyst, hSpectra[idet][ipart][icharge],hSystError[idet][ipart][icharge],kTRUE); | |
1393 | // hSpectra[idet][ipart][icharge]->Draw("same"); | |
1394 | hsyst->Draw("same"); | |
1395 | hSpectraSystError[idet][ipart][icharge] = hsyst; | |
4c0d7bc7 | 1396 | leg->AddEntry(hSpectra[idet][ipart][icharge],TString(partLabel[ipart][icharge])+" (" + detFlag[idet] + ")","lpf"); |
1397 | } | |
1398 | // leg->AddLine(); | |
57b5143c | 1399 | leg->Draw(); |
78033eeb | 1400 | } |
1401 | if(divideCanvas) c1h->cd(1); | |
57b5143c | 1402 | // leg->Draw(); |
2210b002 | 1403 | // ALICEWorkInProgress(c1h,today.Data(),"#splitline{ALICE Preliminary}{Statistical Error Only}"); |
1404 | PrintCanvas(c1h,printFormats); | |
4c0d7bc7 | 1405 | } |
1406 | ||
1407 | ||
eea641de | 1408 | // Draw ratios |
4c0d7bc7 | 1409 | |
eea641de | 1410 | // K-/K+ in the different detectors |
5f2ca1e6 | 1411 | TCanvas * cpm=new TCanvas("cpm","Kminus/Kplus",700,700); |
1412 | cpm->Divide(2,2); | |
1413 | cpm->cd(1); | |
1414 | TH1F* hRatioKPKM_TPC=new TH1F(*(hSpectra[kTPC][kKaon][kNeg])); | |
1415 | hRatioKPKM_TPC->SetMinimum(0.5); | |
1416 | hRatioKPKM_TPC->SetMaximum(1.5); | |
1417 | hRatioKPKM_TPC->Divide(hSpectra[kTPC][kKaon][kPos]); | |
1418 | hRatioKPKM_TPC->GetYaxis()->SetTitle("K-/K+ (TPC)"); | |
1419 | hRatioKPKM_TPC->Draw(); | |
1420 | cpm->cd(2); | |
1421 | TH1F* hRatioKPKM_ITS=new TH1F(*(hSpectra[kITS][kKaon][kNeg])); | |
1422 | hRatioKPKM_ITS->Divide(hSpectra[kITS][kKaon][kPos]); | |
1423 | hRatioKPKM_ITS->SetMinimum(0.5); | |
1424 | hRatioKPKM_ITS->SetMaximum(1.5); | |
1425 | hRatioKPKM_ITS->GetYaxis()->SetTitle("K-/K+ (ITSsa)"); | |
1426 | hRatioKPKM_ITS->Draw(""); | |
1427 | cpm->cd(3); | |
1428 | TH1F* hRatioKPKM_TOF=new TH1F(*(hSpectra[kTOF][kKaon][kNeg])); | |
1429 | hRatioKPKM_TOF->Divide(hSpectra[kTOF][kKaon][kPos]); | |
1430 | hRatioKPKM_TOF->SetMinimum(0.5); | |
1431 | hRatioKPKM_TOF->SetMaximum(1.5); | |
1432 | hRatioKPKM_TOF->GetYaxis()->SetTitle("K-/K+ (TOF)"); | |
1433 | hRatioKPKM_TOF->Draw(""); | |
1434 | cpm->cd(4); | |
1435 | TH1F* hRatioKPKM_ITSTPC=new TH1F(*(hSpectra[kITSTPC][kKaon][kNeg])); | |
1436 | hRatioKPKM_ITSTPC->Divide(hSpectra[kITSTPC][kKaon][kPos]); | |
1437 | hRatioKPKM_ITSTPC->SetMinimum(0.5); | |
1438 | hRatioKPKM_ITSTPC->SetMaximum(1.5); | |
1439 | hRatioKPKM_ITSTPC->GetYaxis()->SetTitle("K-/K+ (ITS Global)"); | |
1440 | hRatioKPKM_ITSTPC->Draw(""); | |
4c0d7bc7 | 1441 | |
eea641de | 1442 | |
1443 | // ITS/TPC | |
5f2ca1e6 | 1444 | TH1F * hRatioITSTPC[kNPart][kNCharge]; |
eea641de | 1445 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ // loop over charges |
1446 | // Create canvas | |
2210b002 | 1447 | TCanvas * c1h = new TCanvas(TString("cITSTPCRatio_")+chargeFlag[icharge],TString("cITSTPCRatio_")+chargeFlag[icharge],1200,500); |
1448 | c1h->Divide(3,1); | |
1449 | c1h->SetGridy(); | |
1450 | TH2F * hemptyLoc = new TH2F(TString("hemptyR")+long(icharge),"ITSsa/TPC ",100,0.,1., 100, 0.5,1.5); | |
1451 | hemptyLoc->SetXTitle("p_{t} (GeV/c)"); | |
1452 | hemptyLoc->SetYTitle("ITSsa / TPC"); | |
eea641de | 1453 | // Loop over particles |
5f2ca1e6 | 1454 | for(Int_t ipart = 0; ipart < kNPart; ipart++) { |
eea641de | 1455 | // Clone histo |
57b5143c | 1456 | hRatioITSTPC[ipart][icharge]=new TH1F(*hSpectraSystError[kITS][ipart][icharge]); |
1457 | Int_t nBinsITS=hSpectraSystError[kITS][ipart][icharge]->GetNbinsX(); | |
1458 | Int_t nBinsTPC=hSpectraSystError[kTPC][ipart][icharge]->GetNbinsX(); | |
eea641de | 1459 | // Loop over ITS bins, |
5f2ca1e6 | 1460 | for(Int_t iBin=1; iBin<=nBinsITS; iBin++){ |
1461 | hRatioITSTPC[ipart][icharge]->SetBinContent(iBin,0.); | |
57b5143c | 1462 | hRatioITSTPC[ipart][icharge]->SetBinError(iBin,0.); |
1463 | Float_t lowPtITS=hSpectraSystError[kITS][ipart][icharge]->GetBinLowEdge(iBin); | |
1464 | Float_t binWidITS=hSpectraSystError[kITS][ipart][icharge]->GetBinWidth(iBin); | |
eea641de | 1465 | // Loop over TPC bins and find overlapping bins to ITS |
5f2ca1e6 | 1466 | for(Int_t jBin=1; jBin<=nBinsTPC; jBin++){ |
57b5143c | 1467 | Float_t lowPtTPC=hSpectraSystError[kTPC][ipart][icharge]->GetBinLowEdge(jBin); |
1468 | Float_t binWidTPC=hSpectraSystError[kTPC][ipart][icharge]->GetBinWidth(jBin); | |
5f2ca1e6 | 1469 | if(TMath::Abs(lowPtITS-lowPtTPC)<0.001 && TMath::Abs(binWidTPC-binWidITS)<0.001){ |
57b5143c | 1470 | Float_t numer=hSpectraSystError[kITS][ipart][icharge]->GetBinContent(iBin); |
1471 | Float_t denom=hSpectraSystError[kTPC][ipart][icharge]->GetBinContent(jBin); | |
1472 | Float_t enumer=hSpectraSystError[kITS][ipart][icharge]->GetBinError(iBin); | |
1473 | Float_t edenom=hSpectraSystError[kTPC][ipart][icharge]->GetBinError(jBin); | |
5f2ca1e6 | 1474 | Double_t ratio=0.; |
1475 | Double_t eratio=0.; | |
1476 | if(numer>0. && denom>0.){ | |
1477 | ratio=numer/denom; | |
1478 | eratio=ratio*TMath::Sqrt((enumer/numer)*(enumer/numer)+(edenom/denom)*(edenom/denom)); | |
1479 | } | |
1480 | hRatioITSTPC[ipart][icharge]->SetBinContent(iBin,ratio); | |
1481 | hRatioITSTPC[ipart][icharge]->SetBinError(iBin,eratio); | |
1482 | break; | |
1483 | } | |
1484 | } | |
1485 | } | |
2210b002 | 1486 | c1h->cd(ipart+1); |
1487 | // hemptyLoc->SetStats(1); | |
1488 | // hemptyLoc->Draw(); | |
1489 | hRatioITSTPC[ipart][icharge]->SetYTitle("ITSsa/TPC"); | |
eea641de | 1490 | hRatioITSTPC[ipart][icharge]->SetStats(1); |
1491 | hRatioITSTPC[ipart][icharge]->GetYaxis()->SetRangeUser(0.5,1.5); | |
1492 | hRatioITSTPC[ipart][icharge]->Draw(""); | |
57b5143c | 1493 | TF1 * f = new TF1("fpol0","[0]"); |
1494 | f->SetParameter(0,1); | |
1495 | hRatioITSTPC[ipart][icharge]->Fit(f,"","same"); | |
eea641de | 1496 | |
5f2ca1e6 | 1497 | } |
2210b002 | 1498 | PrintCanvas(c1h,printFormats); |
5f2ca1e6 | 1499 | } |
1500 | ||
5b103cfc | 1501 | // ITS/TPC |
1502 | TH1F * hRatioITSGlobalTPC[kNPart][kNCharge]; | |
1503 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ // loop over charges | |
1504 | // Create canvas | |
2210b002 | 1505 | TCanvas * c1h = new TCanvas(TString("cITSGlobalTPCRatio_")+chargeFlag[icharge],TString("cITSGlobalTPCRatio_")+chargeFlag[icharge],1200,500); |
1506 | c1h->Divide(3,1); | |
1507 | c1h->SetGridy(); | |
1508 | TH2F * hemptyLoc = new TH2F(TString("hemptyR")+long(icharge),"ITS Global/TPC ",100,0.,1., 100, 0.5,1.5); | |
1509 | hemptyLoc->SetXTitle("p_{t} (GeV/c)"); | |
1510 | hemptyLoc->SetYTitle("ITS Global / TPC"); | |
5b103cfc | 1511 | // Loop over particles |
1512 | for(Int_t ipart = 0; ipart < kNPart; ipart++) { | |
1513 | // Clone histo | |
57b5143c | 1514 | hRatioITSGlobalTPC[ipart][icharge]=new TH1F(*hSpectraSystError[kITSTPC][ipart][icharge]); |
1515 | Int_t nBinsITS=hSpectraSystError[kITSTPC][ipart][icharge]->GetNbinsX(); | |
1516 | Int_t nBinsTPC=hSpectraSystError[kTPC][ipart][icharge]->GetNbinsX(); | |
5b103cfc | 1517 | // Loop over ITS bins, |
1518 | for(Int_t iBin=1; iBin<=nBinsITS; iBin++){ | |
1519 | hRatioITSGlobalTPC[ipart][icharge]->SetBinContent(iBin,0.); | |
57b5143c | 1520 | hRatioITSGlobalTPC[ipart][icharge]->SetBinError(iBin,0.); |
1521 | Float_t lowPtITS=hSpectraSystError[kITSTPC][ipart][icharge]->GetBinLowEdge(iBin); | |
1522 | Float_t binWidITS=hSpectraSystError[kITSTPC][ipart][icharge]->GetBinWidth(iBin); | |
5b103cfc | 1523 | // Loop over TPC bins and find overlapping bins to ITS |
1524 | for(Int_t jBin=1; jBin<=nBinsTPC; jBin++){ | |
57b5143c | 1525 | Float_t lowPtTPC=hSpectraSystError[kTPC][ipart][icharge]->GetBinLowEdge(jBin); |
1526 | Float_t binWidTPC=hSpectraSystError[kTPC][ipart][icharge]->GetBinWidth(jBin); | |
5b103cfc | 1527 | if(TMath::Abs(lowPtITS-lowPtTPC)<0.001 && TMath::Abs(binWidTPC-binWidITS)<0.001){ |
57b5143c | 1528 | Float_t numer=hSpectraSystError[kITSTPC][ipart][icharge]->GetBinContent(iBin); |
1529 | Float_t denom=hSpectraSystError[kTPC][ipart][icharge]->GetBinContent(jBin); | |
1530 | Float_t enumer=hSpectraSystError[kITSTPC][ipart][icharge]->GetBinError(iBin); | |
1531 | Float_t edenom=hSpectraSystError[kTPC][ipart][icharge]->GetBinError(jBin); | |
5b103cfc | 1532 | Double_t ratio=0.; |
1533 | Double_t eratio=0.; | |
1534 | if(numer>0. && denom>0.){ | |
1535 | ratio=numer/denom; | |
1536 | eratio=ratio*TMath::Sqrt((enumer/numer)*(enumer/numer)+(edenom/denom)*(edenom/denom)); | |
1537 | } | |
1538 | hRatioITSGlobalTPC[ipart][icharge]->SetBinContent(iBin,ratio); | |
1539 | hRatioITSGlobalTPC[ipart][icharge]->SetBinError(iBin,eratio); | |
1540 | break; | |
1541 | } | |
1542 | } | |
1543 | } | |
2210b002 | 1544 | c1h->cd(ipart+1); |
1545 | // hemptyLoc->SetStats(1); | |
1546 | // hemptyLoc->Draw(); | |
1547 | hRatioITSGlobalTPC[ipart][icharge]->SetYTitle("ITS Global/TPC"); | |
5b103cfc | 1548 | hRatioITSGlobalTPC[ipart][icharge]->SetStats(1); |
1549 | hRatioITSGlobalTPC[ipart][icharge]->GetYaxis()->SetRangeUser(0.5,1.5); | |
1550 | hRatioITSGlobalTPC[ipart][icharge]->Draw(""); | |
57b5143c | 1551 | TF1 * f = new TF1("fpol0","[0]"); |
1552 | f->SetParameter(0,1); | |
1553 | hRatioITSGlobalTPC[ipart][icharge]->Fit(f,"","same"); | |
5b103cfc | 1554 | |
1555 | } | |
2210b002 | 1556 | PrintCanvas(c1h,printFormats); |
5b103cfc | 1557 | } |
1558 | ||
eea641de | 1559 | // TOF/TPC |
5f2ca1e6 | 1560 | TH1F * hRatioTOFTPC[kNPart][kNCharge]; |
eea641de | 1561 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ // loop over charges |
1562 | // create canvas | |
1563 | TCanvas * c1t = new TCanvas(TString("cTOFTPCRatio_")+chargeFlag[icharge],TString("cTOFTPCRatio_")+chargeFlag[icharge],1200,500); | |
5f2ca1e6 | 1564 | c1t->SetGridy(); |
eea641de | 1565 | c1t->Divide(3,1); |
5f2ca1e6 | 1566 | TH2F * hemptyt = new TH2F(TString("hemptyRt")+long(icharge),"TOF/TPC ",100,0.,1., 100, 0.5,1.5); |
1567 | hemptyt->SetXTitle("p_{t} (GeV/c)"); | |
1568 | hemptyt->SetYTitle("TOF / TPC"); | |
eea641de | 1569 | // hemptyt->Draw(); |
1570 | for(Int_t ipart = 0; ipart < kNPart; ipart++) { // loop over particles | |
1571 | // Clone histo | |
57b5143c | 1572 | hRatioTOFTPC[ipart][icharge]=new TH1F(*hSpectraSystError[kTOF][ipart][icharge]); |
1573 | Int_t nBinsTOF=hSpectraSystError[kTOF][ipart][icharge]->GetNbinsX(); | |
1574 | Int_t nBinsTPC=hSpectraSystError[kTPC][ipart][icharge]->GetNbinsX(); | |
eea641de | 1575 | // Loop over TOF bins |
5f2ca1e6 | 1576 | for(Int_t iBin=1; iBin<=nBinsTOF; iBin++){ |
1577 | hRatioTOFTPC[ipart][icharge]->SetBinContent(iBin,0.); | |
57b5143c | 1578 | hRatioTOFTPC[ipart][icharge]->SetBinError(iBin,0.); |
1579 | Float_t lowPtTOF=hSpectraSystError[kTOF][ipart][icharge]->GetBinLowEdge(iBin); | |
1580 | Float_t binWidTOF=hSpectraSystError[kTOF][ipart][icharge]->GetBinWidth(iBin); | |
eea641de | 1581 | // Loop over TPC bins and find overlapping bins to ITS |
5f2ca1e6 | 1582 | for(Int_t jBin=1; jBin<=nBinsTPC; jBin++){ |
57b5143c | 1583 | Float_t lowPtTPC=hSpectraSystError[kTPC][ipart][icharge]->GetBinLowEdge(jBin); |
1584 | Float_t binWidTPC=hSpectraSystError[kTPC][ipart][icharge]->GetBinWidth(jBin); | |
5f2ca1e6 | 1585 | if(TMath::Abs(lowPtTOF-lowPtTPC)<0.001 && TMath::Abs(binWidTPC-binWidTOF)<0.001){ |
57b5143c | 1586 | Float_t numer=hSpectraSystError[kTOF][ipart][icharge]->GetBinContent(iBin); |
1587 | Float_t denom=hSpectraSystError[kTPC][ipart][icharge]->GetBinContent(jBin); | |
1588 | Float_t enumer=hSpectraSystError[kTOF][ipart][icharge]->GetBinError(iBin); | |
1589 | Float_t edenom=hSpectraSystError[kTPC][ipart][icharge]->GetBinError(jBin); | |
5f2ca1e6 | 1590 | Double_t ratio=0.; |
1591 | Double_t eratio=0.; | |
1592 | if(numer>0. && denom>0.){ | |
1593 | ratio=numer/denom; | |
1594 | eratio=ratio*TMath::Sqrt((enumer/numer)*(enumer/numer)+(edenom/denom)*(edenom/denom)); | |
1595 | } | |
1596 | hRatioTOFTPC[ipart][icharge]->SetBinContent(iBin,ratio); | |
1597 | hRatioTOFTPC[ipart][icharge]->SetBinError(iBin,eratio); | |
1598 | break; | |
1599 | } | |
1600 | } | |
1601 | } | |
eea641de | 1602 | c1t->cd(ipart+1); |
2210b002 | 1603 | hRatioTOFTPC[ipart][icharge]->SetYTitle("TOF/TPC"); |
0e6f7392 | 1604 | hRatioTOFTPC[ipart][icharge]->SetStats(1); |
1605 | hRatioTOFTPC[ipart][icharge]->GetYaxis()->SetRangeUser(0.5,1.5); | |
1606 | hRatioTOFTPC[ipart][icharge]->Draw(""); | |
57b5143c | 1607 | hRatioTOFTPC[ipart][icharge]->Print(); |
1608 | TF1 * f = new TF1("fpol0","[0]"); | |
1609 | f->SetParameter(0,1); | |
1610 | hRatioTOFTPC[ipart][icharge]->Fit(f,"","same"); | |
5f2ca1e6 | 1611 | } |
22abbb27 | 1612 | PrintCanvas(c1t,printFormats); |
5f2ca1e6 | 1613 | } |
1614 | ||
4c0d7bc7 | 1615 | } |
1616 | ||
1617 | void DrawWithJacek() { | |
1618 | ||
1619 | //1. Convert spectra to dNdeta and sum | |
1620 | TH1F * hsum = (TH1F*) htemplate->Clone(); | |
1621 | hsum->Reset(); | |
ee08b77d | 1622 | Int_t idet= iCombInStudy; |
4c0d7bc7 | 1623 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ |
1624 | for(Int_t ipart = 0; ipart < kNPart; ipart++){ | |
1625 | TH1 * h = hSpectra[idet][ipart][icharge]; | |
1626 | Int_t nbin = h->GetNbinsX(); | |
1627 | for(Int_t ibin = 1; ibin <= nbin; ibin++){ | |
1628 | Double_t pt = h->GetBinCenter(ibin); | |
1629 | Double_t mt = TMath::Sqrt(pt*pt + mass[ipart]*mass[ipart]); | |
1630 | Double_t jacobian = pt/mt; | |
1631 | h->SetBinContent(ibin,h->GetBinContent(ibin)*jacobian); | |
1632 | h->SetBinError (ibin,h->GetBinError (ibin)*jacobian); | |
1633 | Int_t ibinSum = hsum->FindBin(pt); | |
1634 | Double_t epsilon = 0.0001; | |
1635 | if ( h->GetBinContent(ibin) > 0 && | |
1636 | (TMath::Abs(h->GetBinLowEdge(ibin) - hsum->GetBinLowEdge(ibinSum)) > epsilon || | |
1637 | TMath::Abs(h->GetBinLowEdge(ibin+1) - hsum->GetBinLowEdge(ibinSum+1)) ) | |
1638 | ) { | |
1639 | cout << "DISCREPANCY IN BIN RANGES" << endl; | |
1640 | cout << pt << " " << ibinSum << " " << ibin << "; " << h->GetBinContent(ibin) << endl | |
1641 | << h->GetBinLowEdge(ibin) << "-" << h->GetBinLowEdge(ibin+1) << endl | |
1642 | << hsum->GetBinLowEdge(ibinSum) << "-" << hsum->GetBinLowEdge(ibinSum+1) << endl; | |
1643 | cout << "" << endl; | |
1644 | } | |
1645 | ||
1646 | ||
1647 | hsum->SetBinContent(ibinSum,hsum->GetBinContent(ibinSum)+h->GetBinContent(ibin)); // EROOR FIXME | |
1648 | hsum->SetBinError (ibinSum,0); | |
1649 | } | |
1650 | // hsum->Add(h); | |
1651 | } | |
1652 | } | |
1653 | ||
1654 | ||
1655 | // Load Jacek and Draw both: | |
77233d47 | 1656 | new TFile ("./Files/pt_allcharged_published.root"); |
1657 | TGraphErrors * gJacek = (TGraphErrors*) gDirectory->Get("dndpt"); | |
1658 | gJacek->SetMarkerStyle(20); | |
1659 | gJacek->SetMarkerColor(kGray); | |
1660 | gJacek->Draw("AP"); | |
1661 | gJacek->GetHistogram()->SetXTitle("p_{t} (GeV)"); | |
1662 | gJacek->GetHistogram()->SetYTitle("d^{2}N/dp_{t}d#eta (GeV^{-1})"); | |
1663 | hsum->Draw("same"); | |
1664 | ||
1665 | TGraphErrors * gRatio = AliBWTools::DivideGraphByHisto(gJacek,hsum); | |
1666 | new TCanvas(); | |
1667 | gRatio->Draw("AP"); | |
4c0d7bc7 | 1668 | |
1669 | ||
1670 | ||
1671 | } | |
1672 | ||
1673 | ||
1674 | void DrawRatioToStar() { | |
1675 | ||
1676 | // Star data | |
1677 | // gROOT->LoadMacro("StarPPSpectra.C"); | |
1678 | TGraphErrors ** gStar = StarPPSpectra(); | |
c4fe358b | 1679 | gStar[0]->SetMarkerStyle(kOpenStar); |
1680 | gStar[1]->SetMarkerStyle(kFullStar); | |
1681 | gStar[2]->SetMarkerStyle(kOpenStar); | |
1682 | gStar[3]->SetMarkerStyle(kFullStar); | |
1683 | gStar[4]->SetMarkerStyle(kOpenStar); | |
1684 | gStar[5]->SetMarkerStyle(kFullStar); | |
4c0d7bc7 | 1685 | |
1686 | // ALICE, INEL -> NSD | |
1687 | Double_t scaleYield = 3.58/3.02; // from paper 2 | |
1688 | for(Int_t ipart = 0; ipart < kNPart; ipart++){ | |
1689 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ | |
ee08b77d | 1690 | hSpectra[iCombInStudy][ipart][icharge]->Scale(scaleYield); |
4c0d7bc7 | 1691 | } |
1692 | } | |
1693 | ||
1694 | ||
1695 | TCanvas * c1 = new TCanvas("cRatioToStarNeg","cRatioToStarNeg"); | |
1696 | TH2F * hempty = new TH2F(TString("hemptyNeg"),"hemptyNeg",100,0.,1.5, 100, 0.001,1.8); | |
1697 | hempty->SetXTitle("p_{t} (GeV/c)"); | |
1698 | hempty->SetYTitle("ALICE/STAR (NSD)"); | |
c4fe358b | 1699 | hempty->Draw(); |
1700 | ||
1701 | TCanvas * c1Comp = new TCanvas("cCompToStarNeg","cCompToStarNeg"); | |
1702 | c1Comp->SetLogy(); | |
1703 | TH2F * hempty2 = new TH2F(TString("hemptyCompNeg"),"hemptyCompNeg",100,0.,1.5, 100, 0.001,10); | |
1704 | hempty2->SetXTitle("p_{t} (GeV/c)"); | |
1705 | hempty2->SetYTitle("1/N_{ev} d^{2}N / dydp_{t} (GeV/c)^{-1} [NSD]"); | |
1706 | hempty2->Draw(); | |
4c0d7bc7 | 1707 | |
1708 | TLegend *leg = new TLegend(0.6510067,0.1853147,0.8892617,0.4178322,"Negative","brNDC"); | |
1709 | leg->SetBorderSize(0); | |
1710 | leg->SetLineColor(1); | |
1711 | leg->SetLineStyle(1); | |
1712 | leg->SetLineWidth(1); | |
1713 | leg->SetFillColor(0); | |
1714 | leg->SetFillStyle(1001); | |
1715 | ||
1716 | ||
c4fe358b | 1717 | |
1718 | c1->cd(); | |
4c0d7bc7 | 1719 | TGraphErrors * g ; |
ee08b77d | 1720 | g = AliBWTools::DivideGraphByHisto(gStar[0],hSpectra[iCombInStudy][kPion][kNeg],1); |
4c0d7bc7 | 1721 | g->SetMarkerStyle(kFullCircle); |
1722 | g->SetMarkerColor(kBlack); | |
1723 | g->Draw("p"); | |
1724 | leg->AddEntry(g,"#pi^{-}","lp"); | |
ee08b77d | 1725 | g = AliBWTools::DivideGraphByHisto(gStar[2],hSpectra[iCombInStudy][kKaon][kNeg],1); |
4c0d7bc7 | 1726 | g->SetMarkerStyle(kOpenCircle); |
1727 | g->SetMarkerColor(kRed); | |
1728 | g->Draw("p"); | |
1729 | leg->AddEntry(g,"K^{-}","lp"); | |
ee08b77d | 1730 | g = AliBWTools::DivideGraphByHisto(gStar[4],hSpectra[iCombInStudy][kProton][kNeg],1); |
4c0d7bc7 | 1731 | g->SetMarkerStyle(kOpenSquare); |
1732 | g->SetMarkerColor(kBlue); | |
1733 | g->Draw("p"); | |
1734 | leg->AddEntry(g,"#bar{p}","lp"); | |
1735 | leg->Draw(); | |
1736 | ||
c4fe358b | 1737 | c1Comp->cd(); |
1738 | gStar[0]->Draw("p"); | |
1739 | hSpectra[iCombInStudy][kPion][kNeg]->Draw("same"); | |
1740 | gStar[2]->Draw("p"); | |
1741 | hSpectra[iCombInStudy][kKaon][kNeg]->Draw("same"); | |
1742 | gStar[4]->Draw("p"); | |
1743 | hSpectra[iCombInStudy][kProton][kNeg]->Draw("same"); | |
4c0d7bc7 | 1744 | |
1745 | ||
1746 | ||
1747 | TCanvas * c2 = new TCanvas("cRatioToStarPos","cRatioToStarPos"); | |
1748 | hempty->Draw(); | |
1749 | leg = new TLegend(0.6510067,0.1853147,0.8892617,0.4178322,"Positive","brNDC"); | |
1750 | leg->SetBorderSize(0); | |
1751 | leg->SetLineColor(1); | |
1752 | leg->SetLineStyle(1); | |
1753 | leg->SetLineWidth(1); | |
1754 | leg->SetFillColor(0); | |
1755 | leg->SetFillStyle(1001); | |
c4fe358b | 1756 | |
1757 | TCanvas * c2Comp = new TCanvas("cCompToStarPos","cCompToStarPos"); | |
1758 | c2Comp->SetLogy(); | |
1759 | hempty2->Draw(); | |
1760 | ||
1761 | c2->cd(); | |
4c0d7bc7 | 1762 | // TGraphErrors * g ; |
ee08b77d | 1763 | g = AliBWTools::DivideGraphByHisto(gStar[1],hSpectra[iCombInStudy][kPion][kPos],1); |
4c0d7bc7 | 1764 | g->SetMarkerStyle(kFullCircle); |
1765 | g->SetMarkerColor(kBlack); | |
1766 | g->Draw("p"); | |
1767 | leg->AddEntry(g,"#pi^{+}","lp"); | |
ee08b77d | 1768 | g = AliBWTools::DivideGraphByHisto(gStar[3],hSpectra[iCombInStudy][kKaon][kPos],1); |
4c0d7bc7 | 1769 | g->SetMarkerStyle(kOpenCircle); |
1770 | g->SetMarkerColor(kRed); | |
1771 | g->Draw("p"); | |
1772 | leg->AddEntry(g,"K^{+}","lp"); | |
ee08b77d | 1773 | g = AliBWTools::DivideGraphByHisto(gStar[5],hSpectra[iCombInStudy][kProton][kPos],1); |
4c0d7bc7 | 1774 | g->SetMarkerStyle(kOpenSquare); |
1775 | g->SetMarkerColor(kBlue); | |
1776 | g->Draw("p"); | |
1777 | leg->AddEntry(g,"p","lp"); | |
c4fe358b | 1778 | leg->Draw(); |
4c0d7bc7 | 1779 | |
1780 | ||
c4fe358b | 1781 | c2Comp->cd(); |
1782 | gStar[1]->Draw("p"); | |
1783 | hSpectra[iCombInStudy][kPion][kPos]->Draw("same"); | |
1784 | gStar[3]->Draw("p"); | |
1785 | hSpectra[iCombInStudy][kKaon][kPos]->Draw("same"); | |
1786 | gStar[5]->Draw("p"); | |
1787 | hSpectra[iCombInStudy][kProton][kPos]->Draw("same"); | |
1788 | ||
4c0d7bc7 | 1789 | |
1790 | c1->Update(); | |
1791 | c2->Update(); | |
1792 | gSystem->ProcessEvents(); | |
1793 | c1->Print(TString(c1->GetName()) + ".eps"); | |
1794 | c2->Print(TString(c2->GetName()) + ".eps"); | |
1795 | ALICEWorkInProgress(c1,today.Data(),"#splitline{ALICE Preliminary}{Statistical Error Only}"); | |
1796 | ALICEWorkInProgress(c2,today.Data(),"#splitline{ALICE Preliminary}{Statistical Error Only}"); | |
1797 | c1->Update(); | |
1798 | c2->Update(); | |
1799 | c1->Print(TString(c1->GetName()) + ".png"); | |
1800 | c2->Print(TString(c2->GetName()) + ".png"); | |
1801 | ||
1802 | ||
1803 | ||
1804 | ||
1805 | } | |
1806 | ||
1807 | ||
1808 | ||
1809 | void DrawRatios() { | |
1810 | ||
1811 | // Draws ratios of combined spectra | |
1812 | ||
4d4a268c | 1813 | // Rebin histograms (hi pt was fluctuating too much in the TOF) |
1814 | RebinHistosForRatios(); | |
1815 | ||
4c0d7bc7 | 1816 | // Compute ratios |
1817 | TH1F * hPosNegRatio[kNPart]; | |
1818 | ||
1819 | for(Int_t ipart = 0; ipart < kNPart; ipart++){ | |
ee08b77d | 1820 | hPosNegRatio[ipart] = (TH1F*) hSpectra[iCombInStudy][ipart][kPos]->Clone(); |
1821 | hPosNegRatio[ipart]->Divide(hSpectra[iCombInStudy][ipart][kNeg]); | |
4c0d7bc7 | 1822 | hPosNegRatio[ipart]->SetYTitle(TString(partLabel[ipart][kPos])+"/"+partLabel[ipart][kNeg]); |
1823 | hPosNegRatio[ipart]->SetMinimum(0.5); | |
1824 | hPosNegRatio[ipart]->SetMaximum(1.5); | |
1825 | } | |
1826 | ||
4d4a268c | 1827 | |
ee08b77d | 1828 | TH1F * hKPiRatio = (TH1F*) hSpectra[iCombInStudy][kKaon][kPos]->Clone(); |
1829 | hKPiRatio->Add(hSpectra[iCombInStudy][kKaon][kNeg]); | |
1830 | TH1F * htmp = (TH1F*) hSpectra[iCombInStudy][kPion][kPos]->Clone(); | |
1831 | htmp->Add(hSpectra[iCombInStudy][kPion][kNeg]); | |
4c0d7bc7 | 1832 | hKPiRatio->Divide(htmp); |
1833 | hKPiRatio->SetYTitle("(K^{+}+K^{-})/(#pi^{+}+#pi^{-})"); | |
4d4a268c | 1834 | // Get ratio of syst errors |
1835 | // FIXME: commented combined error | |
1836 | TH1F * hKPiRatioSyst = new TH1F(*htemplate);; | |
1837 | AliBWTools::GetValueAndError(hKPiRatioSyst,hSpectra[iCombInStudy][kKaon][kPos],hSystError[iCombInStudy][kKaon][kPos],kTRUE); | |
1838 | // AliBWTools::GetHistoCombinedErrors(hKPiRatioSyst,hSpectra[iCombInStudy][kKaon][kPos]); | |
1839 | AliBWTools::GetValueAndError(htmp,hSpectra[iCombInStudy][kKaon][kNeg],hSystError[iCombInStudy][kKaon][kNeg],kTRUE); | |
1840 | // // AliBWTools::GetHistoCombinedErrors(htmp,hSpectra[iCombInStudy][kKaon][kNeg]); | |
1841 | hKPiRatioSyst->Add(htmp); | |
1842 | AliBWTools::GetValueAndError(htmp,hSpectra[iCombInStudy][kPion][kNeg],hSystError[iCombInStudy][kPion][kNeg],kTRUE); | |
1843 | // AliBWTools::GetHistoCombinedErrors(htmp,hSpectra[iCombInStudy][kPion][kNeg]); | |
1844 | TH1F * htmp2 = new TH1F(*htemplate); | |
1845 | AliBWTools::GetValueAndError(htmp2,hSpectra[iCombInStudy][kPion][kPos],hSystError[iCombInStudy][kPion][kPos],kTRUE); | |
1846 | // AliBWTools::GetHistoCombinedErrors(htmp2,hSpectra[iCombInStudy][kPion][kPos]); | |
1847 | htmp->Add(htmp2); | |
1848 | hKPiRatioSyst->Divide(htmp); | |
1849 | hKPiRatioSyst->SetYTitle("(K^{+}+K^{-})/(#pi^{+}+#pi^{-})"); | |
1850 | // Multiply error for sqrt(2): the syst is not reduced by combining positive and negative charges | |
1851 | Int_t nbin = hKPiRatioSyst->GetNbinsX(); | |
1852 | for(Int_t ibin = 1; ibin <= nbin; ibin++){ | |
1853 | hKPiRatioSyst->SetBinError(ibin, hKPiRatioSyst->GetBinError(ibin)*TMath::Sqrt(2)); | |
1854 | } | |
1855 | ||
1856 | ||
1857 | ||
4c0d7bc7 | 1858 | |
1859 | TH1F * hKPiRatioMC[kNTunes]; | |
5f2ca1e6 | 1860 | if(showMC){ |
1861 | for(Int_t itune = 0; itune < kNTunes; itune++){ | |
1862 | hKPiRatioMC[itune] = (TH1F*) hSpectraMC[itune][kKaon][kPos]->Clone(); | |
1863 | hKPiRatioMC[itune]->Add(hSpectraMC[itune][kKaon][kNeg]); | |
2210b002 | 1864 | TH1F * htmpLoc = (TH1F*) hSpectraMC[itune][kPion][kPos]->Clone(); |
1865 | htmpLoc->Add(hSpectraMC[itune][kPion][kNeg]); | |
1866 | hKPiRatioMC[itune]->Divide(htmpLoc); | |
5f2ca1e6 | 1867 | hKPiRatioMC[itune]->SetYTitle("(K^{+}+K^{-})/(#pi^{+}+#pi^{-})"); |
1868 | } | |
4c0d7bc7 | 1869 | } |
1870 | ||
1871 | ||
1872 | ||
ee08b77d | 1873 | TH1F * hPPiRatio = (TH1F*) hSpectra[iCombInStudy][kProton][kPos]->Clone(); |
1874 | hPPiRatio->Add(hSpectra[iCombInStudy][kProton][kNeg]); | |
1875 | htmp = (TH1F*) hSpectra[iCombInStudy][kPion][kPos]->Clone(); | |
1876 | htmp->Add(hSpectra[iCombInStudy][kPion][kNeg]); | |
4c0d7bc7 | 1877 | hPPiRatio->Divide(htmp); |
1878 | hPPiRatio->SetYTitle("(p+#bar{p})/(#pi^{+}+#pi^{-})"); | |
1879 | ||
55056a39 | 1880 | TH1F * hPPiRatioMC[kNTunes]; |
5f2ca1e6 | 1881 | if(showMC){ |
5f2ca1e6 | 1882 | for(Int_t itune = 0; itune < kNTunes; itune++){ |
1883 | hPPiRatioMC[itune] = (TH1F*) hSpectraMC[itune][kProton][kPos]->Clone(); | |
1884 | hPPiRatioMC[itune]->Add(hSpectraMC[itune][kProton][kNeg]); | |
2210b002 | 1885 | TH1F * htmpLoc = (TH1F*) hSpectraMC[itune][kPion][kPos]->Clone(); |
1886 | htmpLoc->Add(hSpectraMC[itune][kPion][kNeg]); | |
1887 | hPPiRatioMC[itune]->Divide(htmpLoc); | |
5f2ca1e6 | 1888 | hPPiRatioMC[itune]->SetYTitle("(p+#bar{p})/(#pi^{+}+#pi^{-})"); |
1889 | } | |
4c0d7bc7 | 1890 | } |
1891 | ||
d8e6677d | 1892 | // Phenix ratio |
1893 | TFile * fPhenix = new TFile("./Files/phenixdata.root"); | |
1894 | TH1F * hPPiRatioPhenix = (TH1F*) fPhenix->Get("ProtonPHNX"); | |
1895 | hPPiRatioPhenix->Add((TH1F*) fPhenix->Get("PbarPHNX")); | |
1896 | htmp = (TH1F*) fPhenix->Get("PiPlusPHNX"); | |
1897 | htmp->Add((TH1F*) fPhenix->Get("PiMinusPHNX")); | |
1898 | hPPiRatioPhenix = AliBWTools::DivideHistosDifferentBins(hPPiRatioPhenix,htmp); | |
1899 | hPPiRatioPhenix->SetMarkerStyle(24); | |
f7bf9679 | 1900 | hPPiRatioPhenix->SetBinContent(14,0); |
1901 | hPPiRatioPhenix->SetBinError(14,0); | |
4c0d7bc7 | 1902 | // Draw |
1903 | // TH2F * hempty = new TH2F(TString("hempty"),"hempty",100,0.,1.5, 100, 0.001,1.8); | |
1904 | // hempty->SetXTitle("p_{t} (GeV/c)"); | |
1905 | // hempty->SetYTitle(""); | |
1906 | ||
1907 | // tmp: overlay levi fits | |
1908 | AliBWFunc * fm2 = new AliBWFunc; | |
1909 | fm2->SetVarType(AliBWFunc::kdNdpt); | |
1910 | TF1 * fLevi[kNPart] [kNCharge]; | |
1911 | fLevi[kPion] [kPos] = fm2->GetLevi (mass[0], 0.1243, 7.614785, 1.524167, "fLeviPiPlus"); | |
1912 | fLevi[kKaon] [kPos] = fm2->GetLevi (mass[1], 0.1625, 5.869318, 0.186361, "fLeviKPlus"); | |
1913 | fLevi[kProton][kPos] = fm2->GetLevi (mass[2], 0.1773, 6.918065, 0.086389, "fLeviPPlus"); | |
1914 | fLevi[kPion] [kNeg] = fm2->GetLevi (mass[0], 0.1267, 7.979582, 1.515908, "fLeviPiNeg"); | |
1915 | fLevi[kKaon] [kNeg] = fm2->GetLevi (mass[1], 0.1721, 6.927956, 0.191140, "fLeviKNeg"); | |
1916 | fLevi[kProton][kNeg] = fm2->GetLevi (mass[2], 0.1782, 8.160362, 0.082091, "fLeviPNeg"); | |
1917 | for(Int_t ipart = 0; ipart < kNPart; ipart++){ | |
1918 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ | |
1919 | fLevi[ipart][icharge]->SetRange(0,4); | |
1920 | } | |
1921 | } | |
1922 | ||
1923 | ||
1924 | for(Int_t ipart = 0; ipart < kNPart; ipart++){ | |
b7018bc7 | 1925 | TString detName = detFlag[iCombInStudy]; |
1926 | detName.ReplaceAll(" ", "_"); | |
1927 | detName.ReplaceAll("+", ""); | |
1928 | ||
5f2ca1e6 | 1929 | TCanvas * c1 = new TCanvas(TString("cRatio_")+detName+TString("_")+partFlag[ipart], TString("cRatio_")+detName+partFlag[ipart]); |
1930 | c1->SetGridy(); | |
4c0d7bc7 | 1931 | hPosNegRatio[ipart]->Draw(); |
1932 | TF1 * fRatio = new TF1 (TString("fRatio")+partFlag[ipart], TString(fLevi[ipart][kPos]->GetName())+"/"+fLevi[ipart][kNeg]->GetName()); | |
1933 | // fRatio->Draw("same"); | |
1934 | fRatio->SetRange(0,5); | |
4d4a268c | 1935 | PrintCanvas(c1,printFormats); |
4c0d7bc7 | 1936 | |
1937 | } | |
1938 | ||
1939 | ||
1940 | TCanvas * c2 = new TCanvas(TString("cRatio_KPi"),TString("cRatio_KPi")); | |
d8e6677d | 1941 | // c2->SetGridy(); |
4d4a268c | 1942 | |
1943 | hKPiRatioSyst->SetMaximum(0.62); | |
1944 | hKPiRatioSyst->SetFillColor(kGray); | |
1945 | hKPiRatioSyst->SetMarkerStyle(24); | |
1946 | // hKPiRatioSyst->Draw("e2"); | |
1947 | // TMP FIXME: ERROR FROM FP: | |
1948 | TFile * f = TFile::Open("./Files/ITSsa-systematics-20101014.root"); | |
1949 | TH1F * hsystRatio = (TH1F*) f->Get("hSystRatioKpi"); | |
1950 | TH1F * hsyst = new TH1F(*htemplate); | |
1951 | AliBWTools::GetValueAndError(hsyst,hKPiRatio,hsystRatio,1); | |
1952 | // AliBWTools::GetHistoCombinedErrors(hsyst,hKPiRatio); | |
1953 | // hsyst->Draw("samee1"); | |
1954 | // END TMP | |
1955 | hKPiRatio->SetMaximum(0.62); | |
1956 | hKPiRatio->Draw(""); | |
4c0d7bc7 | 1957 | TLegend * lMC = new TLegend(0.526846, 0.18007, 0.887584,0.407343); |
1958 | lMC->SetFillColor(kWhite); | |
4d4a268c | 1959 | TGraphErrors ** gStar = StarPPSpectra(); |
1960 | cout << gStar[6] << endl; | |
1961 | ||
1962 | gStar[6]->SetMarkerStyle(kOpenStar); | |
1963 | gStar[6]->Draw("P"); | |
1964 | ||
5f2ca1e6 | 1965 | if(showE735){ |
1966 | gROOT->LoadMacro("GetE735Ratios.C"); | |
4d4a268c | 1967 | // GetE735Ratios(0,0)->Draw("EX0,same"); |
1968 | // GetE735Ratios(0,1)->Draw("EX0,same"); | |
1969 | // GetE735Ratios(0,2)->Draw("EX0,same"); | |
5f2ca1e6 | 1970 | GetE735Ratios(0,3)->Draw("EX0,same"); |
1971 | } | |
4c0d7bc7 | 1972 | hKPiRatio->SetMarkerStyle(20); |
d8e6677d | 1973 | hKPiRatio->SetMarkerColor(kRed); |
1974 | hKPiRatio->SetLineColor(kRed); | |
4d4a268c | 1975 | //hKPiRatio->Draw("same"); |
4c0d7bc7 | 1976 | |
5f2ca1e6 | 1977 | if(showMC){ |
1978 | for(Int_t itune = 0; itune < kNTunes; itune++){ | |
1979 | lMC->AddEntry(hKPiRatioMC[itune],mcTuneName[itune]); | |
1980 | hKPiRatioMC[itune]->SetLineWidth(2); | |
1981 | hKPiRatioMC[itune]->Draw("same,chist"); | |
1982 | } | |
4c0d7bc7 | 1983 | |
5f2ca1e6 | 1984 | lMC->Draw(); |
1985 | } | |
4d4a268c | 1986 | cout << "12" << endl; |
4c0d7bc7 | 1987 | |
5f2ca1e6 | 1988 | if(showE735){ |
1989 | TLegend * l = new TLegend( 0.1879, 0.68, 0.54,0.92); | |
1990 | l->SetFillColor(kWhite); | |
1991 | l->AddEntry(hKPiRatio, "ALICE, #sqrt{s} = 900 GeV","lpf"); | |
4d4a268c | 1992 | l->AddEntry(gStar[6], "STAR, #sqrt{s} = 200 GeV", "lp"); |
1993 | // l->AddEntry(GetE735Ratios(0,0), "E735, #sqrt{s} = 300 GeV","lpf"); | |
1994 | // l->AddEntry(GetE735Ratios(0,1), "E735, #sqrt{s} = 540 GeV","lpf"); | |
1995 | // l->AddEntry(GetE735Ratios(0,2), "E735, #sqrt{s} = 1000 GeV","lpf"); | |
5f2ca1e6 | 1996 | l->AddEntry(GetE735Ratios(0,3), "E735, #sqrt{s} = 1800 GeV","lpf"); |
1997 | l->Draw(); | |
1998 | } | |
4c0d7bc7 | 1999 | |
2000 | ||
2001 | TCanvas * c3 = new TCanvas(TString("cRatio_PPi"),TString("cRatio_PPi")); | |
4d4a268c | 2002 | cout << "DRAWING" << endl; |
d8e6677d | 2003 | // c3->SetGridy(); |
2004 | hPPiRatio->SetMarkerStyle(20); | |
2005 | hPPiRatio->SetMarkerColor(kRed); | |
2006 | hPPiRatio->SetLineColor(kRed); | |
4c0d7bc7 | 2007 | hPPiRatio->Draw(); |
2008 | hPPiRatio->SetMaximum(0.39); | |
5f2ca1e6 | 2009 | if(showMC){ |
2010 | lMC = new TLegend(0.526846, 0.18007, 0.887584,0.407343); | |
2011 | lMC->SetFillColor(kWhite); | |
2012 | ||
2013 | for(Int_t itune = 0; itune < kNTunes; itune++){ | |
2014 | lMC->AddEntry(hKPiRatioMC[itune],mcTuneName[itune]); | |
55056a39 | 2015 | hPPiRatioMC[itune]->SetLineWidth(2); |
2016 | hPPiRatioMC[itune]->Draw("same,chist"); | |
5f2ca1e6 | 2017 | } |
4c0d7bc7 | 2018 | |
5f2ca1e6 | 2019 | lMC->Draw(); |
2020 | } | |
d8e6677d | 2021 | hPPiRatioPhenix->Draw("same"); |
2022 | ||
2023 | ||
2024 | TLegend * l = new TLegend(0.186242, 0.781469, 0.538591, 0.921329); | |
2025 | l->SetFillColor(kWhite); | |
2026 | l->AddEntry(hPPiRatio, "ALICE, #sqrt{s} = 900 GeV","lpf"); | |
2027 | l->AddEntry(hPPiRatioPhenix, "PHENIX, #sqrt{s} = 200 GeV","lpf"); | |
2028 | // TLegend * l = new TLegend(0.206376, 0.77972, 0.600671, 0.909091); | |
2029 | l->Draw(); | |
2030 | ||
4d4a268c | 2031 | PrintCanvas(c2,printFormats); |
2032 | PrintCanvas(c3,printFormats); | |
4c0d7bc7 | 2033 | |
2034 | } | |
2035 | ||
3725266c | 2036 | void DrawWithSyst() { |
2037 | ||
2038 | // Draws detector and combined with syst errors. | |
2039 | ||
2040 | for(Int_t idet = 0; idet < kNHist; idet++){ | |
57b5143c | 2041 | if(idet > kITSTPC && idet != iCombInStudy) continue; |
3725266c | 2042 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ |
2043 | TCanvas * c = new TCanvas(TString("cWithSyst")+detFlag[idet]+chargeFlag[icharge],TString("cWithSyst")+detFlag[idet]+chargeFlag[icharge]); | |
2044 | TH2F * hempty = new TH2F(TString("hempty")+long(icharge),"hempty",100,0.,2.9, 100, 0.0005,5); | |
2045 | hempty->SetXTitle("p_{t} (GeV/c)"); | |
2046 | hempty->SetYTitle("1/N_{ev} d^{2}N / dydp_{t} (GeV/c)^{-1}"); | |
2047 | hempty->GetYaxis()->SetTitleOffset(1.35); | |
2048 | hempty->GetXaxis()->SetTitleOffset(1.1); | |
2049 | hempty->Draw(); | |
2050 | c->SetLogy(); | |
2051 | for(Int_t ipart = 0; ipart < kNPart; ipart++){ | |
2052 | // cout << detFlag[idet] << " " << chargeFlag[icharge] << " " << partFlag[ipart] << endl; | |
2053 | ||
2054 | TString opt = ipart ? "" : "same"; | |
2055 | TH1F * hsyst = new TH1F(*htemplate); | |
2056 | AliBWTools::GetValueAndError(hsyst,hSpectra[idet][ipart][icharge],hSystError[idet][ipart][icharge],kTRUE); | |
55056a39 | 2057 | hsyst->SetFillColor(kGray); |
3725266c | 2058 | hsyst->Draw("e5,same"); |
2059 | hSpectra[idet][ipart][icharge]->Draw("same"); | |
22abbb27 | 2060 | // Do not draw syst error outside of spectra range |
2061 | Float_t lowEdge = AliBWTools::GetLowestNotEmptyBinEdge(hSpectra[idet][ipart][icharge]); | |
2062 | Float_t hiEdge = AliBWTools::GetHighestNotEmptyBinEdge(hSpectra[idet][ipart][icharge]); | |
2063 | Int_t lowBin = hSystError[idet][ipart][icharge]->FindBin(lowEdge); | |
2064 | Int_t hiBin = hSystError[idet][ipart][icharge]->FindBin(hiEdge); | |
2065 | Int_t nbin = hSystError[idet][ipart][icharge]->GetNbinsX(); | |
2210b002 | 2066 | for(Int_t ibin = 1; ibin < lowBin; ibin++){ |
22abbb27 | 2067 | hSystError[idet][ipart][icharge]->SetBinContent(ibin,0); |
2068 | hSystError[idet][ipart][icharge]->SetBinError(ibin,0); | |
2069 | } | |
2070 | for(Int_t ibin = hiBin; ibin <= nbin; ibin++){ | |
2071 | hSystError[idet][ipart][icharge]->SetBinContent(ibin,0); | |
2072 | hSystError[idet][ipart][icharge]->SetBinError(ibin,0); | |
2073 | } | |
2074 | // Draw relative syst error, convert to graph to avoid vertical bars at the edges (empty bins in histo) | |
2075 | AliBWTools::GetGraphFromHisto(hSystError[idet][ipart][icharge])->Draw("LX"); | |
2076 | ||
3725266c | 2077 | |
2078 | ||
2079 | } | |
57b5143c | 2080 | PrintCanvas(c,printFormats); |
3725266c | 2081 | } |
2082 | } | |
2083 | } | |
c4fe358b | 2084 | |
2085 | void Help() { | |
2086 | ||
2087 | cout << "Macro: void CombineSpectra(Int_t analysisType=kDoFits, Int_t fitFuncID = kFitLevi) " << endl; | |
2088 | cout << "" << endl; | |
4d4a268c | 2089 | cout << "WARNING: some function modify the histos upon execution. You are advised to quit root" << endl; |
2090 | cout << " before running a different function" << endl; | |
2091 | cout << "" << endl; | |
c4fe358b | 2092 | cout << "Possible Arguments" << endl; |
2093 | cout << "- analysisType:" << endl; | |
eea641de | 2094 | cout << " kDoFits: Fit Combined Spectra " << endl; |
2095 | cout << " kDoRatios: Particle ratios K/pi and p/pi" << endl; | |
2096 | cout << " kDoSuperposition: Compare different detectors (superimpose and ratios)" << endl; | |
2097 | cout << " kDoCompareStar: Compare combined spectra to star results" << endl; | |
2098 | cout << " kDoDrawWithModels: Compare combined spectra and models" << endl; | |
3725266c | 2099 | cout << " kDoDrawSyst: Draws spectra from individual detectors with their systematic error" << endl; |
77233d47 | 2100 | cout << " kDoCompareToAllCharged: convert to dn/deta/dpt, sum up and compare to all charged" << endl; |
eea641de | 2101 | cout << " kDoHelp: This help" << endl; |
2102 | cout << "- fitFuncID, function used to extrapolate and compute yields" << endl; | |
2103 | cout << " An analitic fit function [kFitLevi, kFitUA1, kFitPowerLaw]" << endl; | |
2104 | cout << " Or a shape from a MC moder [kFitPhojet, kFitAtlasCSC, kFitCMS6D6T, kFitPerugia0]" << endl; | |
2f7f947e | 2105 | cout << " Which is fitted to the data at low pt and used to extrapolate at low pt" << endl; |
c4fe358b | 2106 | |
2107 | ||
2108 | } | |
3725266c | 2109 | |
2110 | void PrintCanvas(TCanvas* c,const TString formats) { | |
2111 | // print a canvas in every of the given comma-separated formats | |
2112 | // ensure the canvas is updated | |
22abbb27 | 2113 | if(!doPrint) return; |
3725266c | 2114 | c->Update(); |
2115 | gSystem->ProcessEvents(); | |
2116 | TObjArray * arr = formats.Tokenize(","); | |
2117 | TIterator * iter = arr->MakeIterator(); | |
2118 | TObjString * element = 0; | |
2119 | TString name =c ->GetName(); | |
2120 | name.ReplaceAll(" ","_"); | |
2121 | name.ReplaceAll("+","Plus"); | |
2122 | name.ReplaceAll("-",""); | |
2123 | while ((element = (TObjString*) iter->Next())) { | |
2124 | c->Print(name+ "."+element->GetString()); | |
2125 | } | |
2126 | } | |
2127 | ||
22abbb27 | 2128 | |
3725266c | 2129 | void LoadLibs(){ |
2130 | ||
2131 | gSystem->Load("libTree.so"); | |
2132 | gSystem->Load("libVMC.so"); | |
2133 | gSystem->Load("libMinuit.so"); | |
2134 | gSystem->Load("libSTEERBase.so"); | |
2135 | gSystem->Load("libESD.so"); | |
2136 | gSystem->Load("libAOD.so"); | |
2137 | gSystem->Load("libANALYSIS.so"); | |
2138 | gSystem->Load("libANALYSISalice.so"); | |
2139 | gSystem->Load("libCORRFW.so"); | |
2140 | gSystem->Load("libPWG2spectra.so"); | |
2141 | ||
2142 | } | |
22abbb27 | 2143 | |
4d4a268c | 2144 | void RebinHistosForRatios() { |
2145 | // rebins all spectra histograms | |
2146 | ||
2147 | cout << "WARNING: REBINNING HISTOS BEFORE RATIOS" << endl; | |
2148 | ||
2149 | ||
2150 | const Float_t rebinBins[] = {0.05,0.1,0.12,0.14,0.16,0.18,0.20,0.25,0.30,0.35,0.4,0.45,0.5,0.55,0.6,0.65,0.7,0.75,0.8,0.85,0.9,0.95,1.0,1.1,1.2,1.3,1.4,1.6,1.8,2.4}; | |
2151 | Int_t nrebinBins=29; | |
2152 | ||
2153 | ||
2154 | ||
2155 | Int_t idet = iCombInStudy; | |
2156 | for(Int_t ipart = 0; ipart < kNPart; ipart++){ | |
2157 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ | |
2158 | ||
2159 | TH1F * h = hSpectra[idet][ipart][icharge]; | |
2160 | TString hname = h->GetName(); | |
2161 | hname += "_rebin"; | |
2162 | ||
2163 | TH1F *hrebin = new TH1F(hname, hname,nrebinBins,rebinBins); | |
2164 | hrebin->Sumw2(); | |
2165 | ||
2166 | for(Int_t i=1;i <= h->GetNbinsX();i++){ | |
2167 | Float_t x = h->GetBinCenter(i); | |
2168 | Int_t j=0; | |
2169 | while(x > rebinBins[j+1]) j++; | |
2170 | j++; | |
2171 | hrebin->AddBinContent(j,h->GetBinContent(i)); | |
2172 | ||
2173 | Float_t err = hrebin->GetBinError(j) + h->GetBinError(i)*h->GetBinError(i); | |
2174 | hrebin->SetBinError(j,err); | |
2175 | } | |
2176 | ||
2177 | for(Int_t i=1;i <= hrebin->GetNbinsX();i++){ | |
2178 | hrebin->SetBinError(i,TMath::Sqrt(hrebin->GetBinError(i))); | |
2179 | } | |
2180 | hSpectra[idet][ipart][icharge] = hrebin; | |
2181 | } | |
2182 | } | |
2183 | ||
2184 | } |