]>
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" |
5f2ca1e6 | 29 | #endif |
4c0d7bc7 | 30 | |
31 | using namespace std; | |
32 | ||
33 | ||
34 | // A bunch of useful enums and constants | |
35 | enum {kPion=0,kKaon,kProton,kNPart}; | |
ee08b77d | 36 | enum {kTPC=0,kTOF,kITS,kITSTPC,kK0,kKinks,kCombTOFTPC,kCombAll,kNHist};// "k0" listed here as a kind of PID method... |
4c0d7bc7 | 37 | const Int_t kNDet = kITS+2; |
38 | enum {kPos=0,kNeg,kNCharge}; | |
39 | enum {kPhojet=0,kPyTuneAtlasCSC, kPyTuneCMS6D6T, kPyTunePerugia0, kNTunes} ; | |
40 | enum {kFitLevi=0, kFitUA1, kFitPowerLaw, | |
41 | kFitPhojet, kFitAtlasCSC, kFitCMS6D6T, kFitPerugia0, | |
42 | kNFit}; | |
c4fe358b | 43 | enum {kDoFits=0, kDoRatios, kDoSuperposition, kDoDrawWithModels, kDoCompareToStar, kDoHelp}; |
4c0d7bc7 | 44 | |
45 | // flags, labels and names | |
46 | const char * partFlag[] = {"Pion", "Kaon", "Proton"}; | |
b7018bc7 | 47 | const char * detFlag[] = {"TPC", "TOF", "ITS", "ITS Global", "K0", "Kinks", "Combined TOF + TPC", "Combined TOF + TPC + ITS"}; |
4c0d7bc7 | 48 | const char * chargeFlag[] = {"Pos", "Neg"}; |
49 | const char * chargeLabel[] = {"Positive", "Negative"}; | |
50 | const char * partLabel[kNPart][kNCharge] = {{"#pi^{+}", "#pi^{-}"}, | |
51 | // {"K^{+} (#times 2)", "K^{-} (#times 2)"}, | |
52 | {"K^{+}", "K^{-}"}, | |
53 | {"p" , "#bar{p}"}}; | |
5f2ca1e6 | 54 | const char * partLatex[kNPart][kNCharge] = {{"$\\pi^{+}$", "$\\pi^{-}$"}, |
55 | // {"K^{+} (#times 2)", "K^{-} (#times 2)"}, | |
56 | {"$K^{+}$", "$K^{-}$"}, | |
57 | {"$p$" , "$\\bar{p}$"}}; | |
4c0d7bc7 | 58 | const char * mcTuneName[] = {"Phojet", |
59 | "Pythia - CSC 306", | |
60 | "Pythia - D6T 109", | |
61 | "Pythia - Perugia0 - 320", }; | |
62 | const char * funcName[] = { "Levi", "UA1", "PowerLaw", "Phojet", "AtlasCSC", "CMS6D6T", "Perugia0"}; | |
63 | ||
64 | // Style | |
65 | //const Int_t marker[] = {25,24,28,20,21}; // set for highlithining marek | |
66 | const Int_t marker[] = {20,24,25,28,21}; // standard set | |
67 | const Int_t color [] = {1,2,4}; | |
68 | ||
69 | const Int_t mcLineColor[] = {kGreen,kRed,kBlue,kBlack}; | |
70 | const Int_t mcLineStyle[] = {1,2,3,4}; | |
71 | ||
72 | ||
73 | // Template needed to combine different detectors | |
ee08b77d | 74 | 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}; |
75 | Int_t nbinsTempl=34; | |
4c0d7bc7 | 76 | |
77 | TH1F * htemplate = new TH1F("htemplate", "htemplate",nbinsTempl, templBins ); | |
78 | ||
79 | // Globals | |
80 | TH1F * hSpectra[kNHist][kNPart][kNCharge]; | |
81 | TH1F * hSpectraMC[kNTunes][kNPart][kNCharge]; | |
82 | Double_t mass[kNPart]; | |
83 | ||
84 | // Functions: | |
85 | // Loading | |
86 | void LoadSpectra() ; | |
87 | void LoadMC() ; | |
88 | ||
89 | // Additional tasks (may be uncommented below) | |
90 | void DrawStar(Int_t icharge); | |
91 | void GetITSResiduals(); | |
92 | void DrawWithModels() ; | |
93 | void DrawAllAndKaons(); | |
94 | void DrawWithJacek(); | |
95 | void DrawRatioToStar(); | |
96 | void DrawRatios(); | |
5f2ca1e6 | 97 | void FitCombined(); |
c4fe358b | 98 | void Help(); |
4c0d7bc7 | 99 | |
100 | // External stuff | |
5f2ca1e6 | 101 | void ALICEWorkInProgress(TCanvas *c,TString today="11/05/2010", TString label = "ALICE performance"){ |
102 | ||
103 | TPad *myPadLogo = new TPad("myPadLogo", "Pad for ALICE Logo",0.72,0.79,0.82,0.89); | |
104 | myPadLogo->SetFillColor(0); | |
105 | myPadLogo->SetBorderMode(0); | |
106 | myPadLogo->SetBorderSize(2); | |
107 | myPadLogo->SetFrameBorderMode(0); | |
108 | myPadLogo->SetLeftMargin(0.0); | |
109 | myPadLogo->SetTopMargin(0.0); | |
110 | myPadLogo->SetBottomMargin(0.0); | |
111 | myPadLogo->SetRightMargin(0.0); | |
112 | myPadLogo->SetFillStyle(0); | |
113 | myPadLogo->Draw(); | |
114 | myPadLogo->cd(); | |
115 | TASImage *myAliceLogo = new TASImage("alice_logo.png"); | |
116 | myAliceLogo->Draw(); | |
117 | c->cd(); | |
118 | TPaveText* t1=new TPaveText(0.65,0.73,0.89,0.78,"NDC"); | |
119 | t1->SetFillStyle(0); | |
120 | t1->SetBorderSize(0); | |
121 | t1->AddText(0.,0.,label); | |
122 | t1->SetTextColor(kRed); | |
123 | t1->SetTextFont(42); | |
124 | t1->Draw(); | |
125 | TPaveText* t2=new TPaveText(0.65,0.65,0.89,0.7,"NDC"); | |
126 | t2->SetFillStyle(0); | |
127 | t2->SetBorderSize(0); | |
128 | t2->SetTextColor(kRed); | |
129 | t2->SetTextFont(52); | |
130 | t2->AddText(0.,0.,today.Data()); | |
131 | t2->Draw(); | |
132 | } | |
4c0d7bc7 | 133 | |
134 | // Used to tag plots | |
135 | TDatime dt; | |
136 | TString today = ""; | |
137 | ||
138 | ||
139 | ||
140 | // Switches | |
141 | Bool_t convertToMT = 0; | |
b7018bc7 | 142 | Bool_t doPrint = 1; |
4c0d7bc7 | 143 | Bool_t scaleKaons = kFALSE; |
c4fe358b | 144 | Bool_t drawStar = kTRUE; // Overlay star when doing fits |
4c0d7bc7 | 145 | Bool_t correctSecondaries = 1; |
146 | Bool_t correctGeantFlukaXS = 1; | |
5f2ca1e6 | 147 | Int_t iCombInStudy = kCombAll; //kCombTOFTPC |
c4fe358b | 148 | Int_t fitFuncID = kFitLevi; |
5f2ca1e6 | 149 | Bool_t showMC=kTRUE; |
150 | Bool_t showE735=kTRUE; | |
eea641de | 151 | Bool_t useSecCorrFromDCA=kFALSE; |
4c0d7bc7 | 152 | |
eea641de | 153 | void CombineSpectra(Int_t analysisType=kDoHelp, Int_t locfitFuncID = kFitLevi) { //kDoSuperposition;//kDoDrawWithModels;// kDoFits; //kDoRatios; |
4c0d7bc7 | 154 | |
155 | // This macro is used to combine the 900 GeV spectra from 2009 | |
156 | ||
c4fe358b | 157 | fitFuncID=locfitFuncID; |
158 | ||
4c0d7bc7 | 159 | // Load stuff |
160 | gSystem->Load("libTree.so"); | |
161 | gSystem->Load("libVMC.so"); | |
162 | gSystem->Load("libMinuit.so"); | |
163 | gSystem->Load("libSTEERBase.so"); | |
164 | gSystem->Load("libESD.so"); | |
165 | gSystem->Load("libAOD.so"); | |
166 | gSystem->Load("libANALYSIS.so"); | |
167 | gSystem->Load("libANALYSISalice.so"); | |
168 | gSystem->Load("libCORRFW.so"); | |
169 | gSystem->Load("libPWG2spectra.so"); | |
170 | ||
eea641de | 171 | // Print Help and quit |
172 | if (analysisType == kDoHelp) { | |
173 | Help(); | |
174 | return; | |
175 | } | |
176 | ||
177 | ||
4c0d7bc7 | 178 | // Set date |
179 | today = today + long(dt.GetDay()) +"/" + long(dt.GetMonth()) +"/"+ long(dt.GetYear()); | |
180 | ||
181 | ||
182 | // Set Masses | |
183 | mass[kPion] = TDatabasePDG::Instance()->GetParticle("pi+")->Mass(); | |
184 | mass[kKaon] = TDatabasePDG::Instance()->GetParticle("K+")->Mass(); | |
185 | mass[kProton] = TDatabasePDG::Instance()->GetParticle("proton")->Mass(); | |
186 | ||
187 | // Load histos | |
188 | LoadSpectra(); | |
189 | LoadMC(); | |
190 | ||
191 | // Additional tasks | |
192 | //DrawStar(icharge); | |
193 | // GetITSResiduals(); | |
5f2ca1e6 | 194 | if(analysisType==kDoSuperposition) DrawAllAndKaons(); |
195 | else if(analysisType==kDoDrawWithModels) DrawWithModels() ; | |
4c0d7bc7 | 196 | //DrawWithJacek(); |
c4fe358b | 197 | else if(analysisType==kDoCompareToStar) DrawRatioToStar(); |
5f2ca1e6 | 198 | else if(analysisType==kDoRatios) DrawRatios(); |
199 | else if(analysisType==kDoFits) FitCombined(); | |
b7018bc7 | 200 | return; |
5f2ca1e6 | 201 | } |
4c0d7bc7 | 202 | |
203 | ||
5f2ca1e6 | 204 | void FitCombined() { |
eea641de | 205 | gStyle->SetOptTitle(0); |
206 | gStyle->SetOptStat(0); | |
207 | ||
4c0d7bc7 | 208 | // Draw combined & Fit |
209 | AliBWFunc * fm = new AliBWFunc; | |
210 | fm->SetVarType(AliBWFunc::kdNdpt); | |
211 | if (convertToMT) fm->SetVarType(AliBWFunc::kOneOverMtdNdmt); | |
212 | ||
213 | // table to print results | |
214 | AliLatexTable table(10,"c|cccccccc"); | |
215 | if (fitFuncID == kFitLevi) { | |
5f2ca1e6 | 216 | table.InsertCustomRow("Part & Yield & Yield (FIT) & T Slope & n & $\\chi^2$/NDF & Min X & Frac Above & $\\langle p_{t} \\rangle$ & $\\langle p_{t}^{2} \\rangle$ \\\\"); |
4c0d7bc7 | 217 | } else if (fitFuncID == kFitPowerLaw) { |
5f2ca1e6 | 218 | table.InsertCustomRow("Part & Yield & Norm & n & pt0 & $\\chi^2$/NDF & Min X & Frac Above & $\\langle p_{t} \\rangle$ & $\\langle p_{t}^{2} \\rangle$ \\\\"); |
4c0d7bc7 | 219 | } else { |
5f2ca1e6 | 220 | table.InsertCustomRow("Part & Yield & Par0 & Par2 & Par1 & $\\chi^2$/NDF & Min X & Frac Above & $\\langle p_{t} \\rangle$ & $\\langle p_{t}^{2} \\rangle$ \\\\"); |
4c0d7bc7 | 221 | |
222 | } | |
223 | ||
224 | AliLatexTable tempTable(4,"c|ccc"); | |
225 | tempTable.InsertCustomRow("Part & Yield & Yield Below & Frac Above\\\\"); | |
226 | ||
227 | ||
5f2ca1e6 | 228 | TH1F* hRatiosToFit[kNPart][kNCharge]; |
4c0d7bc7 | 229 | // Fit all |
230 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ | |
5f2ca1e6 | 231 | |
232 | TCanvas * c2 = new TCanvas(TString("cCombined")+chargeFlag[icharge]+"_"+funcName[fitFuncID], TString("cCombined")+chargeFlag[icharge],700,700); | |
233 | c2->SetTickx(); | |
234 | c2->SetTicky(); | |
235 | c2->SetLeftMargin(0.14); | |
236 | TCanvas * c2r = new TCanvas(TString("cCombinedRatio")+chargeFlag[icharge]+"_"+funcName[fitFuncID], TString("cCombinedRatio")+chargeFlag[icharge],700,700); | |
5f2ca1e6 | 237 | c2->cd(); |
238 | gPad->SetLogy(); | |
239 | TH2F * hempty = new TH2F(TString("hempty")+long(icharge),"hempty",100,0.,2.9, 100, 0.0005,5); | |
4c0d7bc7 | 240 | hempty->SetXTitle("p_{t} (GeV/c)"); |
241 | hempty->SetYTitle("1/N_{ev} d^{2}N / dydp_{t} (GeV/c)^{-1}"); | |
5f2ca1e6 | 242 | hempty->GetYaxis()->SetTitleOffset(1.35); |
243 | hempty->GetXaxis()->SetTitleOffset(1.1); | |
4c0d7bc7 | 244 | hempty->Draw(); |
5f2ca1e6 | 245 | c2r->cd(); |
246 | gPad->SetGridy(); | |
247 | TH2F * hemptyR = new TH2F(TString("hemptyR")+long(icharge),"hemptyR",100,0.,2.9, 100, 0.5,1.5); | |
248 | hemptyR->SetXTitle("p_{t} (GeV/c)"); | |
249 | hemptyR->SetYTitle("Data/Fit"); | |
250 | hemptyR->Draw(); | |
251 | ||
252 | TLegend * l = new TLegend(0.516779, 0.7, 0.89094 ,0.916084, chargeLabel[icharge]); | |
4c0d7bc7 | 253 | l->SetFillColor(kWhite); |
5f2ca1e6 | 254 | l->SetTextSize(0.035); |
255 | TPaveText* tf=new TPaveText(0.18,0.14,0.56,0.29,"NDC"); | |
256 | if(fitFuncID == kFitLevi){ | |
257 | tf->AddText("#frac{dN}{dp_{t}} #propto p_{t} #left(1+#frac{#sqrt{m^{2}+p_{t}^{2}} -m}{nT} #right)^{-n}"); | |
258 | // tf->SetNDC(); | |
259 | tf->SetTextFont(12); | |
260 | tf->SetTextSize(0.035); | |
261 | } | |
4c0d7bc7 | 262 | for(Int_t ipart = 0; ipart < kNPart; ipart++){ |
5f2ca1e6 | 263 | printf(" ----- Fit %s %s ------\n",partFlag[ipart],chargeFlag[icharge]); |
4c0d7bc7 | 264 | Float_t fitmin = 0; |
265 | Float_t fitmax = 3; | |
266 | ||
267 | // Get functions | |
268 | TF1 * func = 0; | |
3472a634 | 269 | Int_t normPar = -1; |
4c0d7bc7 | 270 | if(fitFuncID == kFitLevi) { |
3472a634 | 271 | normPar = 0; // The levi is normalized by parameter 0 |
4c0d7bc7 | 272 | if (ipart == kPion) |
273 | func = fm->GetLevi(mass[ipart], 0.12, 7, 1.5); | |
274 | if (ipart == kKaon) | |
275 | func = fm->GetLevi(mass[ipart], 0.17, 7, 0.17); | |
276 | if (ipart == kProton) | |
5f2ca1e6 | 277 | func = fm->GetLevi(mass[ipart], 0.15, 8.5, 0.09); |
278 | } | |
4c0d7bc7 | 279 | else if(fitFuncID == kFitUA1) func = fm->GetUA1(mass[ipart],0.2,1.25,1.25,0.2,1.5); |
280 | else if(fitFuncID == kFitPowerLaw) { | |
281 | if (ipart == kPion) | |
282 | func = fm->GetPowerLaw(1.0, 8.6, 7); | |
283 | if (ipart == kKaon) | |
284 | func = fm->GetPowerLaw(3.0, 12, 2.6); | |
285 | if (ipart == kProton) | |
286 | func = fm->GetPowerLaw(24, 72, 0.8); | |
287 | } | |
288 | else if(fitFuncID == kFitPhojet) func = fm->GetHistoFunc(hSpectraMC[kPhojet] [ipart][icharge]); | |
289 | else if(fitFuncID == kFitAtlasCSC) func = fm->GetHistoFunc(hSpectraMC[kPyTuneAtlasCSC][ipart][icharge]); | |
290 | else if(fitFuncID == kFitPerugia0) func = fm->GetHistoFunc(hSpectraMC[kPyTunePerugia0][ipart][icharge]); | |
291 | else if(fitFuncID == kFitCMS6D6T) func = fm->GetHistoFunc(hSpectraMC[kPyTuneCMS6D6T] [ipart][icharge]); | |
292 | else { | |
293 | cout << "Unknown fit ID " << fitFuncID << endl; | |
294 | return; | |
295 | } | |
296 | ||
297 | if(fitFuncID >= kFitPhojet){ | |
298 | fitmin = 0.0; | |
299 | fitmax = 1.0; | |
300 | } | |
301 | ||
ee08b77d | 302 | if(!AliBWTools::Fit(hSpectra[iCombInStudy][ipart][icharge],func,fitmin,fitmax)) { |
4c0d7bc7 | 303 | cout << " FIT ERROR " << endl; |
304 | return; | |
305 | } | |
5f2ca1e6 | 306 | c2->cd(); |
ee08b77d | 307 | hSpectra[iCombInStudy][ipart][icharge]->Draw("same"); |
5f2ca1e6 | 308 | TF1* fitfunc=(TF1*)hSpectra[iCombInStudy][ipart][icharge]->GetListOfFunctions()->At(0); |
309 | fitfunc->Draw("same"); | |
310 | fitfunc->SetRange(0,4); | |
311 | fitfunc->SetLineColor(hSpectra[iCombInStudy][ipart][icharge]->GetLineColor()); | |
c4fe358b | 312 | if(drawStar) DrawStar(icharge); |
313 | hRatiosToFit[ipart][icharge]=(TH1F*)hSpectra[iCombInStudy][ipart][icharge]->Clone(Form("hRatio%s%s",chargeFlag[icharge],partFlag[icharge])); // Ratio data/fit | |
5f2ca1e6 | 314 | for(Int_t iBin=1; iBin<hSpectra[iCombInStudy][ipart][icharge]->GetNbinsX(); iBin++){ |
315 | Double_t lowLim=hSpectra[iCombInStudy][ipart][icharge]->GetBinLowEdge(iBin); | |
316 | Double_t highLim=hSpectra[iCombInStudy][ipart][icharge]->GetBinLowEdge(iBin+1); | |
317 | Double_t contFunc=fitfunc->Integral(lowLim,highLim)/(highLim-lowLim); | |
318 | Double_t ratio=hSpectra[iCombInStudy][ipart][icharge]->GetBinContent(iBin)/contFunc; | |
319 | Double_t eratio=hSpectra[iCombInStudy][ipart][icharge]->GetBinError(iBin)/contFunc; | |
320 | hRatiosToFit[ipart][icharge]->SetBinContent(iBin,ratio); | |
321 | hRatiosToFit[ipart][icharge]->SetBinError(iBin,eratio); | |
322 | } | |
323 | // hSpectra[iCombInStudy][ipart][icharge]->GetListOfFunctions()->At(0)->Draw("same"); | |
324 | // ((TF1*)hSpectra[iCombInStudy][ipart][icharge]->GetListOfFunctions()->At(0))->SetRange(0,4); | |
325 | // ((TF1*)hSpectra[iCombInStudy][ipart][icharge]->GetListOfFunctions()->At(0))->SetLineColor(hSpectra[iCombInStudy][ipart][icharge]->GetLineColor()); | |
326 | c2->Update(); | |
ee08b77d | 327 | l->AddEntry(hSpectra[iCombInStudy][ipart][icharge], |
4c0d7bc7 | 328 | scaleKaons && ipart == kKaon ? |
329 | (TString(partLabel[ipart][icharge])+" #times 2").Data() | |
330 | : partLabel[ipart][icharge]); | |
ee08b77d | 331 | // TF1 * fClone = (TF1*) hSpectra[iCombInStudy][ipart][icharge]->GetListOfFunctions()->At(0)->Clone(); |
332 | // hSpectra[iCombInStudy][ipart][icharge]->GetListOfFunctions()->Add(fClone); | |
4c0d7bc7 | 333 | // fClone->SetLineStyle(kDashed); |
334 | // fClone->SetRange(0,100); | |
335 | // fClone->Draw("same"); | |
336 | ||
337 | // populate table | |
338 | // Float_t yield = func->Integral(0.45,1.05); | |
339 | // Float_t yieldE = func->IntegralError(0.45,1.05); | |
340 | ||
341 | Float_t yield = func->Integral(0.,100); | |
342 | //Float_t yieldE = func->IntegralError(0.,100); | |
343 | ||
344 | Double_t yieldTools, yieldETools; | |
345 | Double_t partialYields[3],partialYieldsErrors[3]; | |
ee08b77d | 346 | AliBWTools::GetYield(hSpectra[iCombInStudy][ipart][icharge], func, yieldTools, yieldETools, |
4c0d7bc7 | 347 | 0, 100, partialYields,partialYieldsErrors); |
348 | Double_t tslope = func->GetParameter(2); | |
349 | Double_t tslopeE = func->GetParError(2); | |
350 | ||
5f2ca1e6 | 351 | table.SetNextCol(partLatex[ipart][icharge]); |
4c0d7bc7 | 352 | //table.SetNextCol(yield,yieldE,-4); |
353 | table.SetNextCol(yieldTools, yieldETools,-4); | |
354 | table.SetNextCol(func->GetParameter(0)); | |
355 | table.SetNextCol(tslope,tslopeE,-4); | |
5f2ca1e6 | 356 | table.SetNextCol(func->GetParameter(1),func->GetParError(1)); |
4c0d7bc7 | 357 | table.SetNextCol(Form("%2.2f/%d",func->GetChisquare(),func->GetNDF())); |
ee08b77d | 358 | Float_t lowestPoint = AliBWTools::GetLowestNotEmptyBinEdge(hSpectra[iCombInStudy][ipart][icharge]); |
4c0d7bc7 | 359 | //Float_t lowestPoint = AliBWTools::GetLowestNotEmptyBinEdge(hSpectra[kITS][ipart][icharge]); |
360 | Float_t yieldAbove = func->Integral(lowestPoint,100); | |
361 | table.SetNextCol(lowestPoint,-2); | |
362 | table.SetNextCol(yieldAbove/yield,-2); | |
363 | Float_t mean, meane; | |
364 | Float_t mean2, mean2e; | |
3472a634 | 365 | AliBWTools::GetMean (func, mean, meane , 0.,100., normPar); |
366 | AliBWTools::GetMeanSquare(func, mean2, mean2e, 0.,100., normPar); | |
4c0d7bc7 | 367 | table.SetNextCol(mean, meane ,-4); |
368 | table.SetNextCol(mean2, mean2e,-4); | |
369 | ||
370 | // fMean2->IntegralError(0,100)/func->Integral(0,100),-7); | |
371 | table.InsertRow(); | |
372 | ||
373 | ||
374 | /// TEMP TABLE | |
375 | tempTable.SetNextCol(partLabel[ipart][icharge]); | |
376 | tempTable.SetNextCol(yieldTools, yieldETools, -4); | |
377 | tempTable.SetNextCol(partialYields[1], partialYieldsErrors[1], -4); | |
378 | tempTable.SetNextCol(yieldAbove/yield,-2); | |
379 | tempTable.InsertRow(); | |
5f2ca1e6 | 380 | c2r->cd(); |
381 | hRatiosToFit[ipart][icharge]->Draw("esame"); | |
382 | ||
4c0d7bc7 | 383 | } |
5f2ca1e6 | 384 | c2->cd(); |
385 | l->Draw(); | |
386 | c2r->cd(); | |
4c0d7bc7 | 387 | l->Draw(); |
388 | if (doPrint) { | |
5f2ca1e6 | 389 | c2->cd(); |
4c0d7bc7 | 390 | c2->Update(); |
391 | gSystem->ProcessEvents(); | |
5f2ca1e6 | 392 | tf->Draw(); |
4c0d7bc7 | 393 | c2->Print(TString(c2->GetName()) + ".eps"); |
394 | ALICEWorkInProgress(c2,"","#splitline{ALICE Preliminary}{Statistical Error Only}"); | |
5f2ca1e6 | 395 | c2->Update(); |
4c0d7bc7 | 396 | c2->Print(TString(c2->GetName()) + ".png"); |
5f2ca1e6 | 397 | c2r->Update(); |
398 | gSystem->ProcessEvents(); | |
399 | c2r->Print(TString(c2r->GetName()) + ".eps"); | |
400 | c2r->Print(TString(c2r->GetName()) + ".png"); | |
4c0d7bc7 | 401 | } |
402 | ||
5f2ca1e6 | 403 | |
4c0d7bc7 | 404 | } |
405 | ||
406 | ||
3472a634 | 407 | table.PrintTable("ASCII"); |
4c0d7bc7 | 408 | |
409 | cout << "" << endl; | |
410 | tempTable.PrintTable("ASCII"); | |
411 | ||
5f2ca1e6 | 412 | |
413 | ||
4c0d7bc7 | 414 | } |
415 | ||
416 | void LoadSpectra() { | |
417 | ||
418 | TFile * f=0; | |
419 | ||
420 | // Load | |
421 | ||
422 | ||
423 | // TOF | |
5f2ca1e6 | 424 | // Load Efficiencies |
425 | f = new TFile("./Files/effhistos.root"); | |
426 | TH1D * hEffTrackTOF[kNPart][kNCharge]; | |
427 | TH1D * hEffMatchTOF[kNPart][kNCharge]; | |
428 | hEffTrackTOF[kPion] [kPos] = (TH1D*) f->Get("hpitrk_pos"); | |
429 | hEffTrackTOF[kKaon] [kPos] = (TH1D*) f->Get("hkatrk_pos"); | |
430 | hEffTrackTOF[kProton][kPos] = (TH1D*) f->Get("hprtrk_pos"); | |
431 | hEffMatchTOF[kPion] [kPos] = (TH1D*) f->Get("hpieff_pos"); | |
432 | hEffMatchTOF[kKaon] [kPos] = (TH1D*) f->Get("hkaeff_pos"); | |
433 | hEffMatchTOF[kProton][kPos] = (TH1D*) f->Get("hpreff_pos"); | |
434 | hEffTrackTOF[kPion] [kNeg] = (TH1D*) f->Get("hpitrk_neg"); | |
435 | hEffTrackTOF[kKaon] [kNeg] = (TH1D*) f->Get("hkatrk_neg"); | |
436 | hEffTrackTOF[kProton][kNeg] = (TH1D*) f->Get("hprtrk_neg"); | |
437 | hEffMatchTOF[kPion] [kNeg] = (TH1D*) f->Get("hpieff_neg"); | |
438 | hEffMatchTOF[kKaon] [kNeg] = (TH1D*) f->Get("hkaeff_neg"); | |
439 | hEffMatchTOF[kProton][kNeg] = (TH1D*) f->Get("hpreff_neg"); | |
440 | ||
4c0d7bc7 | 441 | // f = new TFile("./Files/spectra-pos-y.root"); |
5f2ca1e6 | 442 | f = new TFile("./Files/spectraRaw-pos-y.root"); |
4c0d7bc7 | 443 | hSpectra[kTOF][kPion] [kPos]= (TH1F*) f->Get("hpi"); |
444 | hSpectra[kTOF][kProton][kPos]= (TH1F*) f->Get("hpr"); | |
445 | hSpectra[kTOF][kKaon] [kPos]= (TH1F*) f->Get("hka"); | |
446 | hSpectra[kTOF][kPion] [kPos]->SetName("hpiPos"); | |
447 | hSpectra[kTOF][kProton][kPos]->SetName("hprPos"); | |
448 | hSpectra[kTOF][kKaon] [kPos]->SetName("hkaPos"); | |
449 | //f = new TFile("./Files/spectra-neg-y.root"); | |
5f2ca1e6 | 450 | f = new TFile("./Files/spectraRaw-neg-y.root"); |
4c0d7bc7 | 451 | hSpectra[kTOF][kPion] [kNeg]= (TH1F*) f->Get("hpi"); |
452 | hSpectra[kTOF][kProton][kNeg]= (TH1F*) f->Get("hpr"); | |
453 | hSpectra[kTOF][kKaon] [kNeg]= (TH1F*) f->Get("hka"); | |
454 | hSpectra[kTOF][kPion] [kNeg]->SetName("hpiNeg"); | |
455 | hSpectra[kTOF][kProton][kNeg]->SetName("hprNeg"); | |
456 | hSpectra[kTOF][kKaon] [kNeg]->SetName("hkaNeg"); | |
457 | ||
5f2ca1e6 | 458 | // Divide for efficiency |
459 | hSpectra[kTOF][kPion] [kPos]->Divide(hEffTrackTOF[kPion] [kPos]); | |
460 | hSpectra[kTOF][kKaon] [kPos]->Divide(hEffTrackTOF[kKaon] [kPos]); | |
461 | hSpectra[kTOF][kProton][kPos]->Divide(hEffTrackTOF[kProton][kPos]); | |
462 | hSpectra[kTOF][kPion] [kPos]->Divide(hEffMatchTOF[kPion] [kPos]); | |
463 | hSpectra[kTOF][kKaon] [kPos]->Divide(hEffMatchTOF[kKaon] [kPos]); | |
464 | hSpectra[kTOF][kProton][kPos]->Divide(hEffMatchTOF[kProton][kPos]); | |
465 | hSpectra[kTOF][kPion] [kNeg]->Divide(hEffTrackTOF[kPion] [kNeg]); | |
466 | hSpectra[kTOF][kKaon] [kNeg]->Divide(hEffTrackTOF[kKaon] [kNeg]); | |
467 | hSpectra[kTOF][kProton][kNeg]->Divide(hEffTrackTOF[kProton][kNeg]); | |
468 | hSpectra[kTOF][kPion] [kNeg]->Divide(hEffMatchTOF[kPion] [kNeg]); | |
469 | hSpectra[kTOF][kKaon] [kNeg]->Divide(hEffMatchTOF[kKaon] [kNeg]); | |
470 | hSpectra[kTOF][kProton][kNeg]->Divide(hEffMatchTOF[kProton][kNeg]); | |
471 | ||
472 | ||
473 | // Clean UP TOF spectra, removing unwanted points | |
4c0d7bc7 | 474 | cout << "Cleaning Up TOF spectra" << endl; |
475 | Int_t nbin = hSpectra[kTOF][kKaon][kPos]->GetNbinsX(); | |
476 | for(Int_t ibin = 1; ibin <= nbin; ibin++){ | |
477 | Float_t pt = hSpectra[kTOF][kKaon][kPos]->GetBinCenter(ibin); | |
478 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ | |
479 | if (pt > 2.35) { | |
480 | hSpectra[kTOF][kKaon][icharge]->SetBinContent(ibin,0); | |
481 | hSpectra[kTOF][kKaon][icharge]->SetBinError (ibin,0); | |
482 | hSpectra[kTOF][kProton][icharge]->SetBinContent(ibin,0); | |
483 | hSpectra[kTOF][kProton][icharge]->SetBinError (ibin,0); | |
484 | } | |
485 | } | |
486 | } | |
487 | // cout << "Scaling TOF to TPC" << endl; | |
488 | // // Scale protons to TPC | |
489 | // hSpectra[kTOF][kProton][kPos]->Scale(1./1.05); | |
490 | // // Scale pbar so that pbar/p is compatible with Panos | |
491 | // hSpectra[kTOF][kProton][kNeg]->Scale(1./1.1); | |
492 | ||
493 | ||
494 | ||
5f2ca1e6 | 495 | // ITS SA |
496 | f = new TFile("./Files/ITSsaSpectraCorr_20100727.root"); | |
4c0d7bc7 | 497 | hSpectra[kITS][kPion] [kPos]= (TH1F*) f->Get("hSpectraPos0"); |
498 | hSpectra[kITS][kKaon] [kPos]= (TH1F*) f->Get("hSpectraPos1"); | |
499 | hSpectra[kITS][kProton][kPos]= (TH1F*) f->Get("hSpectraPos2"); | |
500 | hSpectra[kITS][kPion] [kNeg]= (TH1F*) f->Get("hSpectraNeg0"); | |
501 | hSpectra[kITS][kKaon] [kNeg]= (TH1F*) f->Get("hSpectraNeg1"); | |
502 | hSpectra[kITS][kProton][kNeg]= (TH1F*) f->Get("hSpectraNeg2"); | |
503 | ||
504 | for(Int_t ipart = 0; ipart < kNPart; ipart++){ | |
505 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ | |
506 | Int_t nbin = hSpectra[kITS][ipart][icharge]->GetNbinsX(); | |
4c0d7bc7 | 507 | for(Int_t ibin = 1; ibin <= nbin; ibin++){ |
508 | if(hSpectra[kITS][ipart][icharge]->GetBinContent(ibin) < 0 ) { | |
509 | hSpectra[kITS][ipart][icharge]->SetBinContent(ibin,0); | |
510 | hSpectra[kITS][ipart][icharge]->SetBinError (ibin,0); | |
511 | } | |
5f2ca1e6 | 512 | if(ipart == kProton && ibin==9){ |
513 | 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)); | |
514 | hSpectra[kITS][ipart][icharge]->SetBinContent(ibin,0); | |
515 | hSpectra[kITS][ipart][icharge]->SetBinError (ibin,0); | |
516 | } | |
4c0d7bc7 | 517 | // if ((ipart == kKaon && ibin >= 12) || (ipart == kProton && ibin >= 20)) { |
518 | // hSpectra[kITS][ipart][icharge]->SetBinContent(ibin,0); | |
519 | // hSpectra[kITS][ipart][icharge]->SetBinError (ibin,0); | |
520 | // } | |
521 | } | |
522 | ||
523 | } | |
524 | } | |
525 | ||
eea641de | 526 | if(useSecCorrFromDCA){ |
527 | TFile* fseccorr = new TFile("./Files/CorrFac-SecProtonsFromDCA-ITSsa.root"); | |
528 | TH1F* hcorrp=(TH1F*)fseccorr->Get("hSecPCorrFromDCA"); | |
529 | TH1F* hcorrpbar=(TH1F*)fseccorr->Get("hSecPbarCorrFromDCA"); | |
530 | hSpectra[kITS][kProton][kPos]->Multiply(hcorrp); | |
531 | hSpectra[kITS][kProton][kNeg]->Multiply(hcorrpbar); | |
532 | fseccorr->Close(); | |
533 | } | |
4c0d7bc7 | 534 | |
535 | // ITS + TPC (Marek) | |
5f2ca1e6 | 536 | f = TFile::Open("./Files/SpectraCorrectedITSBeforeProtons20100720.root"); |
4c0d7bc7 | 537 | TList * list = (TList*) gDirectory->Get("output"); |
538 | hSpectra[kITSTPC][kPion] [kPos]= (TH1F*) list->FindObject("Pions"); | |
539 | hSpectra[kITSTPC][kKaon] [kPos]= (TH1F*) list->FindObject("Kaons"); | |
540 | hSpectra[kITSTPC][kProton][kPos]= (TH1F*) list->FindObject("Protons"); | |
541 | hSpectra[kITSTPC][kPion] [kNeg]= (TH1F*) list->FindObject("AntiPions"); | |
542 | hSpectra[kITSTPC][kKaon] [kNeg]= (TH1F*) list->FindObject("AntiKaons"); | |
543 | hSpectra[kITSTPC][kProton][kNeg]= (TH1F*) list->FindObject("AntiProtons"); | |
544 | ||
545 | // TPC | |
4c0d7bc7 | 546 | f = new TFile("./Files/protonSpectra_20100615.root"); |
547 | hSpectra[kTPC][kProton][kPos]= AliBWTools::GetHistoFromGraph((TGraphErrors*)f->Get("protonPosClassic"),htemplate); | |
548 | hSpectra[kTPC][kProton][kNeg]= AliBWTools::GetHistoFromGraph((TGraphErrors*)f->Get("protonNegClassic"),htemplate); | |
549 | f = new TFile("./Files/pionSpectra_20100615.root"); | |
550 | hSpectra[kTPC][kPion][kPos]= AliBWTools::GetHistoFromGraph((TGraphErrors*)f->Get("pionPosClassic"),htemplate); | |
551 | hSpectra[kTPC][kPion][kNeg]= AliBWTools::GetHistoFromGraph((TGraphErrors*)f->Get("pionNegClassic"),htemplate); | |
552 | f = new TFile("./Files/kaonSpectra_20100615.root"); | |
553 | hSpectra[kTPC][kKaon][kPos]= AliBWTools::GetHistoFromGraph((TGraphErrors*)f->Get("kaonPosClassic"),htemplate); | |
554 | hSpectra[kTPC][kKaon][kNeg]= AliBWTools::GetHistoFromGraph((TGraphErrors*)f->Get("kaonNegClassic"),htemplate); | |
555 | ||
556 | // Clean UP TPC spectra, removing unwanted points | |
557 | cout << "Cleaning Up TPC spectra" << endl; | |
558 | nbin = hSpectra[kTPC][kKaon][kPos]->GetNbinsX(); | |
559 | for(Int_t ibin = 0; ibin < nbin; ibin++){ | |
560 | Float_t pt = hSpectra[kTPC][kKaon][kPos]->GetBinCenter(ibin); | |
5f2ca1e6 | 561 | if (pt > 0.45){ // || pt<0.25) { |
4c0d7bc7 | 562 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ |
563 | hSpectra[kTPC][kKaon][icharge]->SetBinContent(ibin,0); | |
564 | hSpectra[kTPC][kKaon][icharge]->SetBinError (ibin,0); | |
565 | } | |
566 | } | |
567 | if (pt < 0.45) { | |
568 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ | |
569 | hSpectra[kTPC][kProton][icharge]->SetBinContent(ibin,0); | |
570 | hSpectra[kTPC][kProton][icharge]->SetBinError (ibin,0); | |
571 | } | |
572 | } | |
573 | } | |
574 | ||
575 | ||
576 | ||
577 | ||
578 | ||
579 | // K0s | |
5f2ca1e6 | 580 | f = new TFile ("./Files/PtSpectraCorrectedK0sOff.root"); |
4c0d7bc7 | 581 | // hSpectra[kK0][kKaon][kPos] = (TH1F*) AliBWTools::GetdNdPtFromOneOverPt((TH1*) gDirectory->Get("hSpectraOff")); |
582 | hSpectra[kK0][kKaon][kPos] = (TH1F*) gDirectory->Get("hSpectraOff"); | |
583 | // hSpectra[kK0][kKaon][kPos]->Scale(2*TMath::Pi()); | |
584 | // hSpectra[kK0][kKaon][kPos]->Scale(1./272463); | |
585 | hSpectra[kK0][kKaon][kNeg] = hSpectra[kK0][kKaon][kPos]; | |
586 | ||
587 | // Kinks: TO BE FIXED WITH POSITIVES AND NEGATIVES | |
588 | // f = new TFile ("./Files/PtAllKaonKinkRap6Apr24.root"); | |
589 | f = new TFile ("./Files/PtKaonKinkJune13AllPN_20100615.root"); | |
590 | hSpectra[kKinks][kKaon][kPos] = (TH1F*)gDirectory->Get("fptallKPA"); | |
591 | hSpectra[kKinks][kKaon][kNeg] = (TH1F*)gDirectory->Get("fptallKNA"); | |
5f2ca1e6 | 592 | hSpectra[kKinks][kKaon][kPos]->Scale(0.5/0.7); // different rapidity range for kinks |
593 | hSpectra[kKinks][kKaon][kNeg]->Scale(0.5/0.7); // different rapidity range for kinks | |
594 | hSpectra[kKinks][kKaon][kPos]->Scale(276004./263345.); // different N of events | |
595 | hSpectra[kKinks][kKaon][kNeg]->Scale(276004./263345.); // different N of events | |
4c0d7bc7 | 596 | |
5f2ca1e6 | 597 | // Apply correction factors |
4c0d7bc7 | 598 | // Secondaries for protons |
5f2ca1e6 | 599 | |
4c0d7bc7 | 600 | f = new TFile ("./Files/corrFactorProtons_20100615.root"); |
601 | TH1F * hCorrSecondaries = (TH1F*) gDirectory->Get("corrFactorProtons"); | |
602 | if(correctSecondaries) { | |
603 | cout << "CORRECTING SECONDARIES" << endl; | |
604 | ||
605 | for(Int_t idet = 0; idet <= kTOF; idet++){ // TPC and TOF only | |
606 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ | |
607 | Int_t ipart = kProton; | |
608 | TH1* h = hSpectra[idet][ipart][icharge]; // lighter notation below | |
609 | if (h){ | |
610 | Int_t nbins = h->GetNbinsX(); | |
611 | for(Int_t ibin = 0; ibin < nbins; ibin++){ | |
612 | // Float_t pt = convertToMT ? TMath::Sqrt(h->GetBinCenter(ibin)*h->GetBinCenter(ibin)-mass[kProton]*mass[kProton]) : h->GetBinCenter(ibin); | |
613 | Float_t pt = h->GetBinCenter(ibin); | |
614 | if (icharge == kNeg) pt = -pt; | |
615 | Int_t binCorrection = hCorrSecondaries->FindBin(pt); | |
616 | Float_t correction = hCorrSecondaries->GetBinContent(binCorrection); | |
617 | // cout << pt << " " << correction << endl; | |
618 | ||
c4fe358b | 619 | if (correction) {// If the bin is empty this is a 0 |
4c0d7bc7 | 620 | h->SetBinContent(ibin,h->GetBinContent(ibin)/correction); |
621 | h->SetBinError (ibin,h->GetBinError (ibin)/correction); | |
622 | } else if (h->GetBinContent(ibin) > 0) { // If we are skipping a non-empty bin, we notify the user | |
623 | cout << "Not correcting bin "<<ibin << " for protons secondaries, " << detFlag[idet] << " " << chargeFlag[icharge] << endl; | |
624 | cout << " Bin content: " << h->GetBinContent(ibin) << endl; | |
625 | ||
626 | } | |
627 | } | |
628 | } | |
629 | } | |
630 | } | |
631 | } | |
632 | // geant/fluka absorption | |
633 | if(correctGeantFlukaXS) { | |
634 | cout << "CORRECTING GEANT3/FLUKA" << endl; | |
eea641de | 635 | // common file for Kaons |
636 | TFile *fFlukakaons = TFile::Open("./Files/correctionForCrossSection.321.root"); | |
637 | TH1D * hCorrFlukakaon[kNCharge]; | |
638 | hCorrFlukakaon[kPos] = (TH1D*)fFlukakaons->Get("gHistCorrectionForCrossSectionParticles"); | |
639 | hCorrFlukakaon[kNeg] = (TH1D*)fFlukakaons->Get("gHistCorrectionForCrossSectionAntiParticles"); | |
640 | ||
641 | for(Int_t idet = 0; idet < kNDet; idet++){ | |
642 | if( idet != kITS) continue; // comment to use fluka for kaons for all dets | |
643 | if (idet == kTOF) continue; // TOF already corrected | |
644 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ | |
645 | Int_t ipart = kKaon; | |
646 | TH1 * h = hSpectra[kITS][ipart][icharge]; // only ITS sa | |
647 | if (h){ | |
648 | Int_t nbins = h->GetNbinsX(); | |
649 | Int_t nbinsy=hCorrFlukakaon[icharge]->GetNbinsY(); | |
650 | for(Int_t ibin = 0; ibin < nbins; ibin++){ | |
651 | Float_t pt = h->GetBinCenter(ibin); | |
652 | Float_t minPtCorrection = hCorrFlukakaon[icharge]->GetXaxis()->GetBinLowEdge(1); | |
653 | Float_t maxPtCorrection = hCorrFlukakaon[icharge]->GetXaxis()->GetBinLowEdge(nbinsy+1); | |
654 | if (pt < minPtCorrection) pt = minPtCorrection+0.0001; | |
655 | if (pt > maxPtCorrection) pt = maxPtCorrection; | |
656 | Float_t correction = hCorrFlukakaon[icharge]->GetBinContent(hCorrFlukakaon[icharge]->GetXaxis()->FindBin(pt)); | |
657 | if (correction != 0) {// If the bin is empty this is a 0 | |
658 | h->SetBinContent(ibin,h->GetBinContent(ibin)*correction); | |
659 | h->SetBinError (ibin,h->GetBinError (ibin)*correction); | |
660 | } else if (h->GetBinContent(ibin) > 0) { // If we are skipping a non-empty bin, we notify the user | |
661 | cout << "Fluka/GEANT: Not correcting bin "<<ibin << " for protons secondaries, ITS, " << chargeFlag[icharge] << endl; | |
662 | cout << " Bin content: " << h->GetBinContent(ibin) << endl; | |
663 | } | |
664 | } | |
665 | } | |
666 | } | |
667 | } | |
668 | ||
669 | // ITS specific file for protons/antiprotons | |
5f2ca1e6 | 670 | TFile* fITS = new TFile ("./Files/correctionForCrossSectionITS_20100719.root"); |
671 | TH2D * hCorrFlukaITS[kNCharge]; | |
672 | hCorrFlukaITS[kPos] = (TH2D*)fITS->Get("gHistCorrectionForCrossSectionProtons"); | |
673 | hCorrFlukaITS[kNeg] = (TH2D*)fITS->Get("gHistCorrectionForCrossSectionAntiProtons"); | |
674 | ||
675 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ | |
676 | Int_t ipart = kProton; | |
677 | TH1 * h = hSpectra[kITS][ipart][icharge]; // only ITS sa | |
678 | if (h){ | |
679 | Int_t nbins = h->GetNbinsX(); | |
680 | Int_t nbinsy=hCorrFlukaITS[icharge]->GetNbinsY(); | |
681 | for(Int_t ibin = 0; ibin < nbins; ibin++){ | |
682 | Float_t pt = h->GetBinCenter(ibin); | |
683 | Float_t minPtCorrection = hCorrFlukaITS[icharge]->GetYaxis()->GetBinLowEdge(1); | |
684 | Float_t maxPtCorrection = hCorrFlukaITS[icharge]->GetYaxis()->GetBinLowEdge(nbinsy+1); | |
685 | if (pt < minPtCorrection) pt = minPtCorrection+0.0001; | |
686 | if (pt > maxPtCorrection) pt = maxPtCorrection; | |
687 | Float_t correction = hCorrFlukaITS[icharge]->GetBinContent(1,hCorrFlukaITS[icharge]->GetYaxis()->FindBin(pt)); | |
688 | if (correction != 0) {// If the bin is empty this is a 0 | |
689 | h->SetBinContent(ibin,h->GetBinContent(ibin)*correction); | |
690 | h->SetBinError (ibin,h->GetBinError (ibin)*correction); | |
691 | } else if (h->GetBinContent(ibin) > 0) { // If we are skipping a non-empty bin, we notify the user | |
692 | cout << "Fluka/GEANT: Not correcting bin "<<ibin << " for protons secondaries, ITS, " << chargeFlag[icharge] << endl; | |
693 | cout << " Bin content: " << h->GetBinContent(ibin) << endl; | |
694 | } | |
695 | } | |
696 | } | |
697 | } | |
698 | ||
4c0d7bc7 | 699 | |
700 | f = new TFile ("./Files/correctionForCrossSection_20100615.root"); | |
701 | TH2D * hCorrFluka[kNCharge]; | |
702 | hCorrFluka[kPos] = (TH2D*)gDirectory->Get("gHistCorrectionForCrossSectionProtons"); | |
703 | hCorrFluka[kNeg] = (TH2D*)gDirectory->Get("gHistCorrectionForCrossSectionAntiProtons"); | |
704 | for(Int_t idet = 0; idet < kNDet; idet++){ | |
705 | if (idet == kITS) continue; // skip ITS sa | |
706 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ | |
707 | Int_t ipart = kProton; | |
708 | TH1 * h = hSpectra[idet][ipart][icharge]; // lighter notation below | |
709 | if (h){ | |
710 | Int_t nbins = h->GetNbinsX(); | |
711 | for(Int_t ibin = 0; ibin < nbins; ibin++){ | |
712 | // Float_t pt = convertToMT ? | |
713 | // TMath::Sqrt(h->GetBinCenter(ibin)*h->GetBinCenter(ibin)-mass[kProton]*mass[kProton]) : | |
714 | // h->GetBinCenter(ibin); | |
715 | Float_t pt = h->GetBinCenter(ibin); | |
716 | Float_t minPtCorrection = hCorrFluka[icharge]->GetYaxis()->GetBinLowEdge(1); | |
717 | Float_t maxPtCorrection = hCorrFluka[icharge]->GetYaxis()->GetBinLowEdge(h->GetNbinsY()+1); | |
718 | if (pt < minPtCorrection) pt = minPtCorrection+0.0001; | |
719 | if (pt > maxPtCorrection) pt = maxPtCorrection; | |
720 | Float_t correction = hCorrFluka[icharge]->GetBinContent(1,hCorrFluka[icharge]->GetYaxis()->FindBin(pt)); | |
4c0d7bc7 | 721 | |
5f2ca1e6 | 722 | // already in the efficiency correction (F. Noferini) |
723 | if (idet == kTOF && icharge == kNeg) { | |
724 | correction = 1; | |
725 | // // Apply parametrized correction computed by francesco | |
726 | // // Fitted panos correction and using momentum at the outer radius of the TPC | |
4c0d7bc7 | 727 | |
5f2ca1e6 | 728 | // Float_t ptav = pt; // Just to use the same name francesco uses... |
4c0d7bc7 | 729 | |
5f2ca1e6 | 730 | // // from pT constrained at P.V. (ptav) to pT TPC outer (ptTPCout) |
731 | // Float_t ptTPCout=ptav*(1-6.81059e-01*TMath::Exp(-ptav*4.20094)); | |
4c0d7bc7 | 732 | |
5f2ca1e6 | 733 | // // traking correction (fit to Panos) |
734 | // Float_t antiprotonEC = 1 - 0.129758 *TMath::Exp(-ptav*0.679612); | |
4c0d7bc7 | 735 | |
5f2ca1e6 | 736 | // // TOF matching efficiency correction (derived from Panos one scaled for M.B.(TOF)/M.B.(TPC)). |
737 | // Float_t antiprotonEC2 = TMath::Power(1 - 0.129758*TMath::Exp(-ptTPCout*0.679612),0.07162/0.03471); | |
738 | // correction = antiprotonEC * antiprotonEC2; | |
739 | } | |
740 | ||
4c0d7bc7 | 741 | // cout << icharge<< " " << h->GetBinCenter(ibin) << " " << pt << " " << correction << endl; |
742 | if (correction != 0) {// If the bin is empty this is a 0 | |
743 | h->SetBinContent(ibin,h->GetBinContent(ibin)*correction); | |
744 | h->SetBinError (ibin,h->GetBinError (ibin)*correction); | |
745 | // if (idet == kTOF) { | |
746 | // cout << "CORRECTING TOF TWICE" << endl; | |
747 | // h->SetBinContent(ibin,h->GetBinContent(ibin)*correction); | |
748 | // h->SetBinError (ibin,h->GetBinError (ibin)*correction); | |
749 | // } | |
750 | } else if (h->GetBinContent(ibin) > 0) { // If we are skipping a non-empty bin, we notify the user | |
751 | cout << "Fluka/GEANT: Not correcting bin "<<ibin << " for protons secondaries, " << detFlag[idet] << " " << chargeFlag[icharge] << endl; | |
752 | cout << " Bin content: " << h->GetBinContent(ibin) << endl; | |
753 | } | |
754 | ||
755 | } | |
756 | ||
757 | } | |
758 | } | |
759 | } | |
760 | } | |
761 | ||
762 | ||
763 | // Set style and scale | |
764 | for(Int_t idet = 0; idet < kNDet; idet++){ | |
765 | for(Int_t ipart = 0; ipart < kNPart; ipart++){ | |
766 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ | |
767 | if (hSpectra[idet][ipart][icharge]){ | |
768 | hSpectra[idet][ipart][icharge]->SetStats(0); // disable stats | |
769 | hSpectra[idet][ipart][icharge]->SetMarkerColor (color[ipart] ); | |
770 | hSpectra[idet][ipart][icharge]->SetLineColor (color[ipart] ); | |
771 | hSpectra[idet][ipart][icharge]->SetMarkerStyle (marker[ipart]); | |
772 | hSpectra[idet][ipart][icharge]->SetXTitle("p_{t} (GeV/c)"); | |
773 | hSpectra[idet][ipart][icharge]->SetYTitle("1/N_{ev} dN/dp_{t} (GeV/c)^{-1}"); | |
774 | if (convertToMT) { | |
775 | TH1F * htmp = (TH1F*) AliBWTools::GetOneOverPtdNdPt(hSpectra[idet][ipart][icharge]); | |
776 | hSpectra[idet][ipart][icharge] = (TH1F*)AliBWTools::GetdNdmtFromdNdpt(htmp,mass[ipart]); | |
777 | hSpectra[idet][ipart][icharge]->SetXTitle("m_{t} (GeV)"); | |
778 | hSpectra[idet][ipart][icharge]->SetYTitle("1/N_{ev} 1/m_{t} dN/dp_{t} (GeV)^{-1}"); | |
779 | } | |
780 | // if (idet == kTOF || idet == kTPC) { | |
781 | // hSpectra[idet][ipart][icharge]->Scale(1./236763); | |
782 | // } | |
783 | // if (idet == kITS){ | |
784 | // hSpectra[idet][ipart][icharge]->Scale(202945./236763); | |
785 | // } | |
786 | if (scaleKaons && ipart == kKaon) { | |
787 | hSpectra[idet][ipart][icharge]->Scale(2.); | |
788 | } | |
789 | } else { | |
790 | printf("Cannot load %s,%s,%s\n",detFlag[idet],partFlag[ipart],chargeFlag[icharge]); | |
791 | } | |
792 | } | |
793 | } | |
794 | } | |
795 | ||
796 | ||
797 | ||
798 | ||
ee08b77d | 799 | |
800 | // Create fake weights for the mean; To be update once we have syst errors | |
801 | TH1F * hWeights[3]; | |
5f2ca1e6 | 802 | const Double_t kWeights[3] = {0.1,0.1,0.1}; // TPC, TOF, ITS |
ee08b77d | 803 | for(Int_t idet = 0; idet <= kITS ; idet++){ |
804 | hWeights[idet] = (TH1F*) hSpectra[idet][kPion][kPos]->Clone(); | |
805 | Int_t nbin = hWeights[idet]->GetNbinsX(); | |
806 | for(Int_t ibin = 0; ibin < nbin; ibin++){ | |
807 | hWeights[idet]->SetBinError(ibin, kWeights[idet]); | |
808 | } | |
809 | } | |
810 | ||
811 | ||
812 | ||
4c0d7bc7 | 813 | // Combine detectors |
814 | for(Int_t ipart = 0; ipart < kNPart; ipart++){ | |
815 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ | |
816 | TH1F * htemplLocal = htemplate; // If we are converting to 1/mt dNdmt we need to convert the template as well... | |
817 | ||
818 | if (convertToMT) { | |
819 | TH1F * htmp = (TH1F*) AliBWTools::GetOneOverPtdNdPt(htemplate); | |
820 | htemplLocal = (TH1F*)AliBWTools::GetdNdmtFromdNdpt(htmp,mass[ipart]); | |
821 | ||
822 | } | |
823 | hSpectra[kCombTOFTPC][ipart][icharge] = AliBWTools::CombineHistos(hSpectra[kTOF][ipart][icharge], | |
824 | hSpectra[kTPC][ipart][icharge], | |
825 | htemplLocal,1.);; | |
ee08b77d | 826 | hSpectra[kCombAll][ipart][icharge] = |
827 | AliBWTools::Combine3HistosWithErrors(hSpectra[kITS][ipart][icharge], | |
828 | hSpectra[kTPC][ipart][icharge], | |
829 | hSpectra[kTOF][ipart][icharge], | |
830 | hWeights[kITS], | |
831 | hWeights[kTPC], | |
832 | hWeights[kTOF], | |
833 | htemplLocal,1);; | |
4c0d7bc7 | 834 | // if (convertToMT) { |
835 | // TH1F * htmp = (TH1F*) AliBWTools::GetOneOverPtdNdPt(hSpectra[kCombTOFTPC][ipart][icharge]); | |
836 | // hSpectra[kCombTOFTPC][ipart][icharge] = (TH1F*)AliBWTools::GetdNdmtFromdNdpt(htmp,mass[ipart]); | |
837 | // hSpectra[kCombTOFTPC][ipart][icharge]->SetXTitle("m_{t} (GeV)"); | |
838 | // hSpectra[kCombTOFTPC][ipart][icharge]->SetYTitle("1/N_{ev} 1/m_{t} dN/dp_{t} (GeV)^{-1}"); | |
839 | // } | |
840 | } | |
841 | } | |
842 | ||
843 | ||
844 | // Scale for the number of inelastic collisions and correct for | |
845 | // efficiency losses due to physics selection: | |
846 | ||
847 | Double_t effPhysSel[kNPart]; | |
848 | effPhysSel[kPion] = 1.012; | |
849 | effPhysSel[kKaon] = 1.013; | |
850 | effPhysSel[kProton] = 1.014; | |
851 | ||
852 | ||
853 | for(Int_t idet = 0; idet < kNHist; idet++){ | |
854 | for(Int_t ipart = 0; ipart < kNPart; ipart++){ | |
855 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ | |
856 | if(hSpectra[idet][ipart][icharge]) { | |
857 | // cout << "Scaling!" << endl; | |
5f2ca1e6 | 858 | if(idet!=kK0){ |
859 | hSpectra[idet][ipart][icharge]->Scale(1.*effPhysSel[ipart]/278366.15); // Scale PhysSel tutti? // FIXME | |
860 | } | |
4c0d7bc7 | 861 | } |
862 | } | |
863 | } | |
864 | } | |
865 | ||
866 | ||
867 | } | |
868 | ||
869 | void LoadMC() { | |
870 | ||
871 | TFile * f = 0; | |
872 | const char * evClass= "INEL"; | |
873 | const char * files[] = {"./Files/dndeta_Phojet_10M_900GeV.root", | |
874 | "./Files/dndeta_AtlasCSC306_10M_900GeV.root", | |
875 | "./Files/dndeta_CMSD6T109_10M_900GeV.root", | |
876 | "./Files/dndeta_Perugia0320_10M_900GeV.root", }; | |
877 | ||
878 | // Phojet | |
879 | for(Int_t itune = 0; itune < kNTunes; itune++){ | |
880 | f = new TFile(files[itune]); | |
881 | for(Int_t ipart = 0; ipart < kNPart; ipart++){ | |
882 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ | |
883 | hSpectraMC[itune][ipart][icharge] = (TH1F*) f->Get(Form("fHistPtID_%s_%s%s",evClass,partFlag[ipart],icharge==0 ? "Pos" : "Neg")); | |
884 | } | |
885 | } | |
886 | } | |
887 | ||
888 | ||
889 | // Set style | |
890 | for(Int_t itune = 0; itune < kNTunes; itune++){ | |
891 | for(Int_t ipart = 0; ipart < kNPart; ipart++){ | |
892 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ | |
893 | if (hSpectraMC[itune][ipart][icharge]){ | |
894 | hSpectraMC[itune][ipart][icharge]->SetName(TString(hSpectraMC[itune][ipart][icharge]->GetName())+mcTuneName[itune]); | |
895 | hSpectraMC[itune][ipart][icharge]->SetMarkerColor (mcLineColor[itune] ); | |
896 | hSpectraMC[itune][ipart][icharge]->SetLineColor (mcLineColor[itune] ); | |
897 | hSpectraMC[itune][ipart][icharge]->SetLineStyle (mcLineStyle[itune] ); | |
898 | hSpectraMC[itune][ipart][icharge]->SetMarkerStyle (1); | |
899 | if (convertToMT) { | |
900 | TH1F * htmp = (TH1F*)AliBWTools::GetOneOverPtdNdPt(hSpectraMC[itune][ipart][icharge]); | |
901 | hSpectraMC[itune][ipart][icharge] = (TH1F*)AliBWTools::GetdNdmtFromdNdpt(htmp,mass[ipart]); | |
902 | hSpectraMC[itune][ipart][icharge]->SetXTitle("m_{t} (GeV)"); | |
903 | hSpectraMC[itune][ipart][icharge]->SetYTitle("1/N_{ev} 1/m_{t} dN/dp_{t} (GeV)^{-1}"); | |
904 | } | |
905 | ||
906 | } else { | |
907 | printf("Cannot load MC # %d,%s,%s\n",itune,partFlag[ipart],chargeFlag[icharge]); | |
908 | } | |
909 | } | |
910 | } | |
911 | } | |
912 | ||
913 | } | |
914 | ||
915 | void DrawStar(Int_t icharge) { | |
916 | // cout << icharge << endl; | |
917 | ||
918 | // gROOT->LoadMacro("StarPPSpectra.C"); | |
919 | TGraphErrors ** gStar = StarPPSpectra(); | |
920 | ||
921 | for(Int_t istar = 0; istar < 6; istar++){ | |
922 | gStar[istar]->SetMarkerStyle(kOpenStar); | |
923 | if (icharge==kPos && (istar%2) ) gStar[istar]->Draw("P"); | |
924 | else if (icharge==kNeg && !(istar%2) ) gStar[istar]->Draw("P"); | |
925 | else cout << "Skipping Star " << istar << endl; | |
926 | } | |
927 | } | |
928 | ||
929 | void GetITSResiduals() { | |
930 | ||
931 | ||
932 | for(Int_t ipart = 0; ipart < kNPart; ipart++){ | |
933 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ | |
934 | // cout << "1 " << ipart << " " << icharge << endl; | |
935 | ||
936 | // gSystem->ProcessEvents(); | |
937 | // gSystem->Sleep(1000); | |
938 | TF1* f = (TF1*) hSpectra[kCombTOFTPC][ipart][icharge]->GetListOfFunctions()->At(0); | |
939 | TH1F * hres1, *hres2; | |
940 | AliBWTools::GetResiduals(hSpectra[kITS][ipart][icharge], f, &hres1, &hres2); | |
941 | ||
942 | TCanvas * c1 = new TCanvas(TString(partFlag[ipart])+"_"+chargeFlag[icharge],TString(partFlag[ipart])+"_"+chargeFlag[icharge]); | |
943 | c1->SetLogy(); | |
944 | hSpectra[kCombTOFTPC][ipart][icharge]->Draw(); | |
945 | hSpectra[kITS][ipart][icharge]->SetMarkerStyle(24); | |
946 | hSpectra[kCombTOFTPC][ipart][icharge]->SetMarkerStyle(20); | |
947 | hSpectra[kITS][ipart][icharge]->Draw("same"); | |
948 | hSpectra[kCombTOFTPC][ipart][icharge]->GetListOfFunctions()->At(0)->Draw("same"); | |
949 | TLegend* l = new TLegend( 0.182886, 0.192308, 0.505034,0.384615, TString(partLabel[ipart][icharge])+" "+chargeFlag[icharge]); | |
950 | l->AddEntry(hSpectra[kCombTOFTPC][ipart][icharge], "TOF+TPC"); | |
951 | l->AddEntry(hSpectra[kITS][ipart][icharge], "ITS"); | |
952 | l->AddEntry(f, "Fit to TOF+TPC"); | |
953 | l->Draw(); | |
954 | ||
955 | TCanvas * c2 = new TCanvas(TString(partFlag[ipart])+"_"+chargeFlag[icharge]+"_res", | |
956 | TString(partFlag[ipart])+"_"+chargeFlag[icharge]+"_res"); | |
957 | c2->SetGridy(); | |
958 | hres2->SetMinimum(-1); | |
959 | hres2->SetMaximum(1); | |
960 | hres2->Draw(); | |
961 | hres2->GetYaxis()->SetTitleOffset(1.2); | |
962 | Float_t x = AliBWTools::GetLowestNotEmptyBinEdge(hSpectra[kCombTOFTPC][ipart][icharge]); | |
963 | TLine * line = new TLine(x,-1,x,1); | |
964 | line->SetLineStyle(kDashed); | |
965 | line->Draw("same"); | |
966 | ||
967 | if (doPrint) { | |
968 | c1->Update(); | |
969 | c2->Update(); | |
970 | gSystem->ProcessEvents(); | |
971 | c1->Print(TString(c1->GetName()) + ".png"); | |
972 | c2->Print(TString(c2->GetName()) + ".png"); | |
973 | } | |
974 | } | |
975 | } | |
976 | } | |
977 | ||
978 | void DrawWithModels() { | |
979 | ||
980 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ | |
981 | ||
982 | // Draw with models | |
983 | for(Int_t ipart = 0; ipart < kNPart; ipart++){ | |
984 | // Multipad canvas | |
5f2ca1e6 | 985 | TCanvas * c1 = new TCanvas(TString("cSpectra")+partFlag[ipart]+chargeFlag[icharge],TString("cSpectra")+partFlag[ipart]+chargeFlag[icharge],700,700); |
986 | TPad *p1 = new TPad(TString("p1")+partFlag[ipart]+chargeFlag[icharge], "p1", 0.0, 0.35, 1.0, 0.95, 0, 0, 0); | |
4c0d7bc7 | 987 | p1->SetBottomMargin(0); |
988 | p1->Draw(); | |
989 | ||
5f2ca1e6 | 990 | TPad *p2 = new TPad(TString("p2")+partFlag[ipart]+chargeFlag[icharge], "p2", 0.0, 0.05, 1.0, 0.35, 0, 0, 0); |
4c0d7bc7 | 991 | p2->SetTopMargin(0); |
5f2ca1e6 | 992 | p2->SetBottomMargin(0.3); |
4c0d7bc7 | 993 | p2->Draw(); |
994 | ||
995 | Float_t scaleFonts = (0.95-0.3)/(0.3-0.05); | |
996 | ||
997 | // Draw spectra | |
998 | p1->cd(); | |
999 | p1->SetLogy(); | |
5f2ca1e6 | 1000 | TH2F * hempty = new TH2F(TString("hempty")+long(icharge)+long(ipart),"hempty",100,0.,4, 100, 0.0015,5); |
4c0d7bc7 | 1001 | hempty->SetXTitle("p_{t} (GeV/c)"); |
1002 | hempty->SetYTitle("1/N_{ev} d^{2}N / dydp_{t} (GeV/c)^{-1}"); | |
1003 | hempty->Draw(); | |
1004 | c1->SetLogy(); | |
1005 | ||
1006 | ||
1007 | TLegend * l =new TLegend( 0.543624, 0.431818, 0.892617,0.629371); | |
1008 | l->SetFillColor(kWhite); | |
ee08b77d | 1009 | hSpectra[iCombInStudy][ipart][icharge]->Draw("same"); |
4c0d7bc7 | 1010 | l->AddEntry(hSpectra[kTOF][ipart][icharge],TString ("Data: ")+partLabel[ipart][icharge]); |
1011 | for(Int_t itune = 0; itune < kNTunes; itune++){ | |
1012 | l->AddEntry(hSpectraMC[itune][ipart][icharge],mcTuneName[itune]); | |
1013 | hSpectraMC[itune][ipart][icharge]->SetLineWidth(2); | |
1014 | hSpectraMC[itune][ipart][icharge]->Draw("same,chist"); | |
1015 | } | |
1016 | l->Draw("same"); | |
1017 | ||
1018 | // Draw ratio | |
1019 | p2->cd(); | |
5f2ca1e6 | 1020 | TH2F * hemptyr = new TH2F(TString("hemptyratio")+long(icharge)+long(ipart),"hempty",100,0.,4, 100, 0.01,2.99); |
4c0d7bc7 | 1021 | hemptyr->SetXTitle("p_{t} (GeV/c)"); |
1022 | hemptyr->SetYTitle("Data/MC"); | |
1023 | hemptyr->GetXaxis()->SetLabelSize(0.04*scaleFonts); | |
1024 | hemptyr->GetYaxis()->SetLabelSize(0.04*scaleFonts); | |
1025 | hemptyr->GetYaxis()->SetTitleSize(0.05*scaleFonts); | |
1026 | hemptyr->GetYaxis()->SetTitleOffset(1.4/scaleFonts); | |
1027 | hemptyr->GetXaxis()->SetTitleSize(0.05*scaleFonts); | |
1028 | hemptyr->SetTickLength(0.03*scaleFonts, "X"); | |
1029 | hemptyr->SetTickLength(0.02*scaleFonts, "Y"); | |
1030 | // hemptyr->GetXaxis()->SetTitleOffset(1.4/scaleFonts); | |
1031 | hemptyr->GetYaxis()->SetNdivisions(5); | |
1032 | hemptyr->Draw(""); | |
1033 | ||
1034 | AliBWFunc fm; | |
1035 | for(Int_t itune = 0; itune < kNTunes; itune++){ | |
1036 | TF1* f = fm.GetHistoFunc(hSpectraMC[itune][ipart][icharge], TString("f")+mcTuneName[itune]); | |
1037 | ||
1038 | // l->AddEntry(hSpectraMC[itune][ipart][icharge],mcTuneName[itune]); | |
ee08b77d | 1039 | TH1F* hRatio = AliBWTools::DivideHistoByFunc(hSpectra[iCombInStudy][ipart][icharge],f); |
4c0d7bc7 | 1040 | hRatio->SetLineStyle(hSpectraMC[itune][ipart][icharge]->GetLineStyle()); |
1041 | hRatio->SetLineColor(hSpectraMC[itune][ipart][icharge]->GetLineColor()); | |
1042 | hRatio->SetLineWidth(hSpectraMC[itune][ipart][icharge]->GetLineWidth()); | |
1043 | hRatio->Draw("lhist,same"); | |
1044 | } | |
1045 | ||
1046 | ||
1047 | ||
1048 | if(doPrint) { | |
1049 | c1->Update(); | |
1050 | gSystem->ProcessEvents(); | |
1051 | c1->Print(TString(c1->GetName())+".eps"); | |
1052 | ALICEWorkInProgress(c1,"","#splitline{ALICE Preliminary}{Statistical Error Only}"); | |
1053 | c1->Print(TString(c1->GetName())+".png"); | |
1054 | ||
1055 | } | |
1056 | } | |
1057 | } | |
1058 | ||
1059 | ||
1060 | ||
1061 | } | |
1062 | ||
1063 | void DrawAllAndKaons() { | |
1064 | ||
1065 | ||
1066 | // gROOT->LoadMacro("ALICEWorkInProgress.C"); | |
1067 | ||
eea641de | 1068 | // gStyle->SetOptFit(0); |
1069 | gStyle->SetStatX(0.9); | |
1070 | gStyle->SetStatY(0.88); | |
1071 | gStyle->SetStatW(0.4); | |
1072 | gStyle->SetStatH(0.1); | |
4c0d7bc7 | 1073 | |
ee08b77d | 1074 | TH1F * hKaonsAllTPCTOF = (TH1F*) hSpectra[iCombInStudy][kKaon][kPos]->Clone(); |
1075 | hKaonsAllTPCTOF->Add(hSpectra[iCombInStudy][kKaon][kNeg]); | |
4c0d7bc7 | 1076 | |
1077 | TH1F * hK0Scaled = (TH1F*) hSpectra[kK0][kKaon][kPos]->Clone(); | |
1078 | hK0Scaled->Add(hSpectra[kK0][kKaon][kPos]); | |
1079 | ||
1080 | hSpectra[kKinks][kKaon][kPos]->SetMarkerStyle(25); | |
5f2ca1e6 | 1081 | hSpectra[kKinks][kKaon][kPos]->SetLineColor(4); |
4c0d7bc7 | 1082 | hSpectra[kKinks][kKaon][kPos]->SetStats(0); |
1083 | TH1F * hKinksAll = (TH1F*) hSpectra[kKinks][kKaon][kPos]->Clone(); | |
1084 | hKinksAll->Add(hSpectra[kKinks][kKaon][kNeg]); | |
1085 | ||
5f2ca1e6 | 1086 | TCanvas * c1 = new TCanvas("cKaons","cKaons",700,700); |
4c0d7bc7 | 1087 | c1->SetLogy(); |
5f2ca1e6 | 1088 | TH2F * hempty = new TH2F("hempty_allkaons","hempty",100,0.,3, 100, 1e-3,6); |
4c0d7bc7 | 1089 | hempty->SetXTitle("p_{t} (GeV/c)"); |
1090 | hempty->SetYTitle("dN / dp_{t} (A.U.)"); | |
1091 | hempty->Draw(); | |
4c0d7bc7 | 1092 | hK0Scaled->Draw("same"); |
5f2ca1e6 | 1093 | hKaonsAllTPCTOF->Draw("same"); |
4c0d7bc7 | 1094 | hKinksAll->Draw("same"); |
1095 | ||
1096 | TLegend * leg = new TLegend(0.2013423,0.2255245,0.5503356,0.4335664,NULL,"brNDC"); | |
5f2ca1e6 | 1097 | // leg->SetBorderSize(0); |
4c0d7bc7 | 1098 | // leg->SetLineColor(1); |
1099 | // leg->SetLineStyle(1); | |
1100 | // leg->SetLineWidth(1); | |
1101 | // leg->SetFillColor(19); | |
5f2ca1e6 | 1102 | leg->SetFillColor(0); |
1103 | TLegendEntry *entry=leg->AddEntry(hKaonsAllTPCTOF,"K^{+} + K^{-}, ITS+TPC+TOF ","lpf"); | |
1104 | entry=leg->AddEntry(hK0Scaled,"K^{0} #times 2","lpf"); | |
1105 | entry=leg->AddEntry(hKinksAll,"K^{+} + K ^{-}, Kinks","lpf"); | |
4c0d7bc7 | 1106 | leg->Draw(); |
1107 | ||
5f2ca1e6 | 1108 | ALICEWorkInProgress(c1,today.Data(),"#splitline{ALICE Prelimiary}{Statistical Error Only}"); |
4c0d7bc7 | 1109 | TLatex * tex = new TLatex(0.2120805,0.01288336,"Statistical error only"); |
1110 | tex->SetTextColor(2); | |
1111 | tex->SetTextFont(42); | |
1112 | tex->SetTextSize(0.03496503); | |
1113 | tex->Draw(); | |
1114 | ||
1115 | c1->Update(); | |
1116 | if(doPrint) c1->Print(TString(c1->GetName())+".png"); | |
1117 | ||
1118 | // Draw all "stable" hadrons | |
1119 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ | |
5f2ca1e6 | 1120 | TCanvas * c1 = new TCanvas(TString("cAll_")+chargeFlag[icharge],TString("cAll_")+chargeFlag[icharge],700,700); |
4c0d7bc7 | 1121 | c1->SetLogy(); |
5f2ca1e6 | 1122 | c1->SetLeftMargin(0.14); |
4c0d7bc7 | 1123 | TH2F * hempty = new TH2F(TString("hempty")+long(icharge),"hempty",100,0.,4, 100, 1e-4,10); |
1124 | hempty->SetXTitle("p_{t} (GeV/c)"); | |
5f2ca1e6 | 1125 | hempty->SetYTitle("1/N_{ev} d^{2}N / dydp_{t} (GeV/c)^{-1}"); |
1126 | hempty->GetYaxis()->SetTitleOffset(1.35); | |
1127 | hempty->GetXaxis()->SetTitleOffset(1.1); | |
4c0d7bc7 | 1128 | hempty->Draw(); |
5f2ca1e6 | 1129 | leg = new TLegend( 0.645973, 0.2, 0.892617,0.636364, NULL,"brNDC"); |
1130 | leg->SetFillColor(0); | |
4c0d7bc7 | 1131 | for(Int_t ipart = 0; ipart < kNPart; ipart++) { |
1132 | for(Int_t idet = 0; idet <= kITSTPC; idet++){ | |
1133 | // if (idet == kITS) continue; | |
1134 | // if (idet == kITSTPC) hSpectra[idet][ipart][icharge]->SetMarkerColor(kGreen); | |
1135 | hSpectra[idet][ipart][icharge]->SetMarkerStyle(marker[idet]); | |
1136 | hSpectra[idet][ipart][icharge]->Draw("same"); | |
1137 | leg->AddEntry(hSpectra[idet][ipart][icharge],TString(partLabel[ipart][icharge])+" (" + detFlag[idet] + ")","lpf"); | |
1138 | } | |
1139 | // leg->AddLine(); | |
1140 | } | |
1141 | leg->Draw(); | |
5f2ca1e6 | 1142 | ALICEWorkInProgress(c1,today.Data(),"#splitline{ALICE Preliminary}{Statistical Error Only}"); |
4c0d7bc7 | 1143 | c1->Update(); |
1144 | if(doPrint) c1->Print(TString(c1->GetName())+".png"); | |
1145 | } | |
1146 | ||
1147 | ||
eea641de | 1148 | // Draw ratios |
4c0d7bc7 | 1149 | |
eea641de | 1150 | // K-/K+ in the different detectors |
5f2ca1e6 | 1151 | TCanvas * cpm=new TCanvas("cpm","Kminus/Kplus",700,700); |
1152 | cpm->Divide(2,2); | |
1153 | cpm->cd(1); | |
1154 | TH1F* hRatioKPKM_TPC=new TH1F(*(hSpectra[kTPC][kKaon][kNeg])); | |
1155 | hRatioKPKM_TPC->SetMinimum(0.5); | |
1156 | hRatioKPKM_TPC->SetMaximum(1.5); | |
1157 | hRatioKPKM_TPC->Divide(hSpectra[kTPC][kKaon][kPos]); | |
1158 | hRatioKPKM_TPC->GetYaxis()->SetTitle("K-/K+ (TPC)"); | |
1159 | hRatioKPKM_TPC->Draw(); | |
1160 | cpm->cd(2); | |
1161 | TH1F* hRatioKPKM_ITS=new TH1F(*(hSpectra[kITS][kKaon][kNeg])); | |
1162 | hRatioKPKM_ITS->Divide(hSpectra[kITS][kKaon][kPos]); | |
1163 | hRatioKPKM_ITS->SetMinimum(0.5); | |
1164 | hRatioKPKM_ITS->SetMaximum(1.5); | |
1165 | hRatioKPKM_ITS->GetYaxis()->SetTitle("K-/K+ (ITSsa)"); | |
1166 | hRatioKPKM_ITS->Draw(""); | |
1167 | cpm->cd(3); | |
1168 | TH1F* hRatioKPKM_TOF=new TH1F(*(hSpectra[kTOF][kKaon][kNeg])); | |
1169 | hRatioKPKM_TOF->Divide(hSpectra[kTOF][kKaon][kPos]); | |
1170 | hRatioKPKM_TOF->SetMinimum(0.5); | |
1171 | hRatioKPKM_TOF->SetMaximum(1.5); | |
1172 | hRatioKPKM_TOF->GetYaxis()->SetTitle("K-/K+ (TOF)"); | |
1173 | hRatioKPKM_TOF->Draw(""); | |
1174 | cpm->cd(4); | |
1175 | TH1F* hRatioKPKM_ITSTPC=new TH1F(*(hSpectra[kITSTPC][kKaon][kNeg])); | |
1176 | hRatioKPKM_ITSTPC->Divide(hSpectra[kITSTPC][kKaon][kPos]); | |
1177 | hRatioKPKM_ITSTPC->SetMinimum(0.5); | |
1178 | hRatioKPKM_ITSTPC->SetMaximum(1.5); | |
1179 | hRatioKPKM_ITSTPC->GetYaxis()->SetTitle("K-/K+ (ITS Global)"); | |
1180 | hRatioKPKM_ITSTPC->Draw(""); | |
4c0d7bc7 | 1181 | |
eea641de | 1182 | |
1183 | // ITS/TPC | |
5f2ca1e6 | 1184 | TH1F * hRatioITSTPC[kNPart][kNCharge]; |
eea641de | 1185 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ // loop over charges |
1186 | // Create canvas | |
1187 | TCanvas * c1 = new TCanvas(TString("cITSTPCRatio_")+chargeFlag[icharge],TString("cITSTPCRatio_")+chargeFlag[icharge],1200,500); | |
1188 | c1->Divide(3,1); | |
5f2ca1e6 | 1189 | c1->SetGridy(); |
1190 | TH2F * hempty = new TH2F(TString("hemptyR")+long(icharge),"ITSsa/TPC ",100,0.,1., 100, 0.5,1.5); | |
1191 | hempty->SetXTitle("p_{t} (GeV/c)"); | |
1192 | hempty->SetYTitle("ITSsa / TPC"); | |
eea641de | 1193 | // Loop over particles |
5f2ca1e6 | 1194 | for(Int_t ipart = 0; ipart < kNPart; ipart++) { |
eea641de | 1195 | // Clone histo |
5f2ca1e6 | 1196 | hRatioITSTPC[ipart][icharge]=new TH1F(*hSpectra[kITS][ipart][icharge]); |
1197 | Int_t nBinsITS=hSpectra[kITS][ipart][icharge]->GetNbinsX(); | |
1198 | Int_t nBinsTPC=hSpectra[kTPC][ipart][icharge]->GetNbinsX(); | |
eea641de | 1199 | // Loop over ITS bins, |
5f2ca1e6 | 1200 | for(Int_t iBin=1; iBin<=nBinsITS; iBin++){ |
1201 | hRatioITSTPC[ipart][icharge]->SetBinContent(iBin,0.); | |
1202 | hRatioITSTPC[ipart][icharge]->SetBinContent(iBin,0.); | |
1203 | Float_t lowPtITS=hSpectra[kITS][ipart][icharge]->GetBinLowEdge(iBin); | |
1204 | Float_t binWidITS=hSpectra[kITS][ipart][icharge]->GetBinWidth(iBin); | |
eea641de | 1205 | // Loop over TPC bins and find overlapping bins to ITS |
5f2ca1e6 | 1206 | for(Int_t jBin=1; jBin<=nBinsTPC; jBin++){ |
1207 | Float_t lowPtTPC=hSpectra[kTPC][ipart][icharge]->GetBinLowEdge(jBin); | |
1208 | Float_t binWidTPC=hSpectra[kTPC][ipart][icharge]->GetBinWidth(jBin); | |
1209 | if(TMath::Abs(lowPtITS-lowPtTPC)<0.001 && TMath::Abs(binWidTPC-binWidITS)<0.001){ | |
1210 | Float_t numer=hSpectra[kITS][ipart][icharge]->GetBinContent(iBin); | |
1211 | Float_t denom=hSpectra[kTPC][ipart][icharge]->GetBinContent(jBin); | |
1212 | Float_t enumer=hSpectra[kITS][ipart][icharge]->GetBinError(iBin); | |
1213 | Float_t edenom=hSpectra[kTPC][ipart][icharge]->GetBinError(jBin); | |
1214 | Double_t ratio=0.; | |
1215 | Double_t eratio=0.; | |
1216 | if(numer>0. && denom>0.){ | |
1217 | ratio=numer/denom; | |
1218 | eratio=ratio*TMath::Sqrt((enumer/numer)*(enumer/numer)+(edenom/denom)*(edenom/denom)); | |
1219 | } | |
1220 | hRatioITSTPC[ipart][icharge]->SetBinContent(iBin,ratio); | |
1221 | hRatioITSTPC[ipart][icharge]->SetBinError(iBin,eratio); | |
1222 | break; | |
1223 | } | |
1224 | } | |
1225 | } | |
eea641de | 1226 | c1->cd(ipart+1); |
1227 | // hempty->SetStats(1); | |
1228 | // hempty->Draw(); | |
1229 | hRatioITSTPC[ipart][icharge]->SetStats(1); | |
1230 | hRatioITSTPC[ipart][icharge]->GetYaxis()->SetRangeUser(0.5,1.5); | |
1231 | hRatioITSTPC[ipart][icharge]->Draw(""); | |
1232 | hRatioITSTPC[ipart][icharge]->Fit("pol0","","same"); | |
1233 | ||
5f2ca1e6 | 1234 | } |
1235 | if(doPrint) c1->Print(TString(c1->GetName())+".png"); | |
1236 | } | |
1237 | ||
eea641de | 1238 | // TOF/TPC |
5f2ca1e6 | 1239 | TH1F * hRatioTOFTPC[kNPart][kNCharge]; |
eea641de | 1240 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ // loop over charges |
1241 | // create canvas | |
1242 | TCanvas * c1t = new TCanvas(TString("cTOFTPCRatio_")+chargeFlag[icharge],TString("cTOFTPCRatio_")+chargeFlag[icharge],1200,500); | |
5f2ca1e6 | 1243 | c1t->SetGridy(); |
eea641de | 1244 | c1t->Divide(3,1); |
5f2ca1e6 | 1245 | TH2F * hemptyt = new TH2F(TString("hemptyRt")+long(icharge),"TOF/TPC ",100,0.,1., 100, 0.5,1.5); |
1246 | hemptyt->SetXTitle("p_{t} (GeV/c)"); | |
1247 | hemptyt->SetYTitle("TOF / TPC"); | |
eea641de | 1248 | // hemptyt->Draw(); |
1249 | for(Int_t ipart = 0; ipart < kNPart; ipart++) { // loop over particles | |
1250 | // Clone histo | |
5f2ca1e6 | 1251 | hRatioTOFTPC[ipart][icharge]=new TH1F(*hSpectra[kTOF][ipart][icharge]); |
1252 | Int_t nBinsTOF=hSpectra[kTOF][ipart][icharge]->GetNbinsX(); | |
1253 | Int_t nBinsTPC=hSpectra[kTPC][ipart][icharge]->GetNbinsX(); | |
eea641de | 1254 | // Loop over TOF bins |
5f2ca1e6 | 1255 | for(Int_t iBin=1; iBin<=nBinsTOF; iBin++){ |
1256 | hRatioTOFTPC[ipart][icharge]->SetBinContent(iBin,0.); | |
1257 | hRatioTOFTPC[ipart][icharge]->SetBinContent(iBin,0.); | |
1258 | Float_t lowPtTOF=hSpectra[kTOF][ipart][icharge]->GetBinLowEdge(iBin); | |
1259 | Float_t binWidTOF=hSpectra[kTOF][ipart][icharge]->GetBinWidth(iBin); | |
eea641de | 1260 | // Loop over TPC bins and find overlapping bins to ITS |
5f2ca1e6 | 1261 | for(Int_t jBin=1; jBin<=nBinsTPC; jBin++){ |
1262 | Float_t lowPtTPC=hSpectra[kTPC][ipart][icharge]->GetBinLowEdge(jBin); | |
1263 | Float_t binWidTPC=hSpectra[kTPC][ipart][icharge]->GetBinWidth(jBin); | |
1264 | if(TMath::Abs(lowPtTOF-lowPtTPC)<0.001 && TMath::Abs(binWidTPC-binWidTOF)<0.001){ | |
1265 | Float_t numer=hSpectra[kTOF][ipart][icharge]->GetBinContent(iBin); | |
1266 | Float_t denom=hSpectra[kTPC][ipart][icharge]->GetBinContent(jBin); | |
1267 | Float_t enumer=hSpectra[kTOF][ipart][icharge]->GetBinError(iBin); | |
1268 | Float_t edenom=hSpectra[kTPC][ipart][icharge]->GetBinError(jBin); | |
1269 | Double_t ratio=0.; | |
1270 | Double_t eratio=0.; | |
1271 | if(numer>0. && denom>0.){ | |
1272 | ratio=numer/denom; | |
1273 | eratio=ratio*TMath::Sqrt((enumer/numer)*(enumer/numer)+(edenom/denom)*(edenom/denom)); | |
1274 | } | |
1275 | hRatioTOFTPC[ipart][icharge]->SetBinContent(iBin,ratio); | |
1276 | hRatioTOFTPC[ipart][icharge]->SetBinError(iBin,eratio); | |
1277 | break; | |
1278 | } | |
1279 | } | |
1280 | } | |
eea641de | 1281 | c1t->cd(ipart+1); |
1282 | hRatioITSTPC[ipart][icharge]->SetStats(1); | |
1283 | hRatioITSTPC[ipart][icharge]->GetYaxis()->SetRangeUser(0.5,1.5); | |
1284 | hRatioITSTPC[ipart][icharge]->Draw(""); | |
1285 | hRatioITSTPC[ipart][icharge]->Fit("pol0","","same"); | |
5f2ca1e6 | 1286 | } |
1287 | if(doPrint) c1t->Print(TString(c1t->GetName())+".png"); | |
1288 | } | |
1289 | ||
4c0d7bc7 | 1290 | } |
1291 | ||
1292 | void DrawWithJacek() { | |
1293 | ||
1294 | //1. Convert spectra to dNdeta and sum | |
1295 | TH1F * hsum = (TH1F*) htemplate->Clone(); | |
1296 | hsum->Reset(); | |
ee08b77d | 1297 | Int_t idet= iCombInStudy; |
4c0d7bc7 | 1298 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ |
1299 | for(Int_t ipart = 0; ipart < kNPart; ipart++){ | |
1300 | TH1 * h = hSpectra[idet][ipart][icharge]; | |
1301 | Int_t nbin = h->GetNbinsX(); | |
1302 | for(Int_t ibin = 1; ibin <= nbin; ibin++){ | |
1303 | Double_t pt = h->GetBinCenter(ibin); | |
1304 | Double_t mt = TMath::Sqrt(pt*pt + mass[ipart]*mass[ipart]); | |
1305 | Double_t jacobian = pt/mt; | |
1306 | h->SetBinContent(ibin,h->GetBinContent(ibin)*jacobian); | |
1307 | h->SetBinError (ibin,h->GetBinError (ibin)*jacobian); | |
1308 | Int_t ibinSum = hsum->FindBin(pt); | |
1309 | Double_t epsilon = 0.0001; | |
1310 | if ( h->GetBinContent(ibin) > 0 && | |
1311 | (TMath::Abs(h->GetBinLowEdge(ibin) - hsum->GetBinLowEdge(ibinSum)) > epsilon || | |
1312 | TMath::Abs(h->GetBinLowEdge(ibin+1) - hsum->GetBinLowEdge(ibinSum+1)) ) | |
1313 | ) { | |
1314 | cout << "DISCREPANCY IN BIN RANGES" << endl; | |
1315 | cout << pt << " " << ibinSum << " " << ibin << "; " << h->GetBinContent(ibin) << endl | |
1316 | << h->GetBinLowEdge(ibin) << "-" << h->GetBinLowEdge(ibin+1) << endl | |
1317 | << hsum->GetBinLowEdge(ibinSum) << "-" << hsum->GetBinLowEdge(ibinSum+1) << endl; | |
1318 | cout << "" << endl; | |
1319 | } | |
1320 | ||
1321 | ||
1322 | hsum->SetBinContent(ibinSum,hsum->GetBinContent(ibinSum)+h->GetBinContent(ibin)); // EROOR FIXME | |
1323 | hsum->SetBinError (ibinSum,0); | |
1324 | } | |
1325 | // hsum->Add(h); | |
1326 | } | |
1327 | } | |
1328 | ||
1329 | ||
1330 | // Load Jacek and Draw both: | |
5f2ca1e6 | 1331 | // new TFile ("./Files/dNdPt_Data_Points_ALICE_900GeV.root"); |
1332 | // TGraphErrors * gJacek = (TGraphErrors*) gDirectory->Get("inel"); | |
1333 | // gJacek->Draw("AP"); | |
1334 | // hsum->Draw("same"); | |
4c0d7bc7 | 1335 | |
5f2ca1e6 | 1336 | // TGraphErrors * gRatio = AliBWTools::DivideGraphByHisto(gJacek,hsum); |
4c0d7bc7 | 1337 | |
1338 | ||
5f2ca1e6 | 1339 | // new TCanvas(); |
1340 | // gRatio->Draw("AP"); | |
4c0d7bc7 | 1341 | |
1342 | ||
1343 | ||
1344 | } | |
1345 | ||
1346 | ||
1347 | void DrawRatioToStar() { | |
1348 | ||
1349 | // Star data | |
1350 | // gROOT->LoadMacro("StarPPSpectra.C"); | |
1351 | TGraphErrors ** gStar = StarPPSpectra(); | |
c4fe358b | 1352 | gStar[0]->SetMarkerStyle(kOpenStar); |
1353 | gStar[1]->SetMarkerStyle(kFullStar); | |
1354 | gStar[2]->SetMarkerStyle(kOpenStar); | |
1355 | gStar[3]->SetMarkerStyle(kFullStar); | |
1356 | gStar[4]->SetMarkerStyle(kOpenStar); | |
1357 | gStar[5]->SetMarkerStyle(kFullStar); | |
4c0d7bc7 | 1358 | |
1359 | // ALICE, INEL -> NSD | |
1360 | Double_t scaleYield = 3.58/3.02; // from paper 2 | |
1361 | for(Int_t ipart = 0; ipart < kNPart; ipart++){ | |
1362 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ | |
ee08b77d | 1363 | hSpectra[iCombInStudy][ipart][icharge]->Scale(scaleYield); |
4c0d7bc7 | 1364 | } |
1365 | } | |
1366 | ||
1367 | ||
1368 | TCanvas * c1 = new TCanvas("cRatioToStarNeg","cRatioToStarNeg"); | |
1369 | TH2F * hempty = new TH2F(TString("hemptyNeg"),"hemptyNeg",100,0.,1.5, 100, 0.001,1.8); | |
1370 | hempty->SetXTitle("p_{t} (GeV/c)"); | |
1371 | hempty->SetYTitle("ALICE/STAR (NSD)"); | |
c4fe358b | 1372 | hempty->Draw(); |
1373 | ||
1374 | TCanvas * c1Comp = new TCanvas("cCompToStarNeg","cCompToStarNeg"); | |
1375 | c1Comp->SetLogy(); | |
1376 | TH2F * hempty2 = new TH2F(TString("hemptyCompNeg"),"hemptyCompNeg",100,0.,1.5, 100, 0.001,10); | |
1377 | hempty2->SetXTitle("p_{t} (GeV/c)"); | |
1378 | hempty2->SetYTitle("1/N_{ev} d^{2}N / dydp_{t} (GeV/c)^{-1} [NSD]"); | |
1379 | hempty2->Draw(); | |
4c0d7bc7 | 1380 | |
1381 | TLegend *leg = new TLegend(0.6510067,0.1853147,0.8892617,0.4178322,"Negative","brNDC"); | |
1382 | leg->SetBorderSize(0); | |
1383 | leg->SetLineColor(1); | |
1384 | leg->SetLineStyle(1); | |
1385 | leg->SetLineWidth(1); | |
1386 | leg->SetFillColor(0); | |
1387 | leg->SetFillStyle(1001); | |
1388 | ||
1389 | ||
c4fe358b | 1390 | |
1391 | c1->cd(); | |
4c0d7bc7 | 1392 | TGraphErrors * g ; |
ee08b77d | 1393 | g = AliBWTools::DivideGraphByHisto(gStar[0],hSpectra[iCombInStudy][kPion][kNeg],1); |
4c0d7bc7 | 1394 | g->SetMarkerStyle(kFullCircle); |
1395 | g->SetMarkerColor(kBlack); | |
1396 | g->Draw("p"); | |
1397 | leg->AddEntry(g,"#pi^{-}","lp"); | |
ee08b77d | 1398 | g = AliBWTools::DivideGraphByHisto(gStar[2],hSpectra[iCombInStudy][kKaon][kNeg],1); |
4c0d7bc7 | 1399 | g->SetMarkerStyle(kOpenCircle); |
1400 | g->SetMarkerColor(kRed); | |
1401 | g->Draw("p"); | |
1402 | leg->AddEntry(g,"K^{-}","lp"); | |
ee08b77d | 1403 | g = AliBWTools::DivideGraphByHisto(gStar[4],hSpectra[iCombInStudy][kProton][kNeg],1); |
4c0d7bc7 | 1404 | g->SetMarkerStyle(kOpenSquare); |
1405 | g->SetMarkerColor(kBlue); | |
1406 | g->Draw("p"); | |
1407 | leg->AddEntry(g,"#bar{p}","lp"); | |
1408 | leg->Draw(); | |
1409 | ||
c4fe358b | 1410 | c1Comp->cd(); |
1411 | gStar[0]->Draw("p"); | |
1412 | hSpectra[iCombInStudy][kPion][kNeg]->Draw("same"); | |
1413 | gStar[2]->Draw("p"); | |
1414 | hSpectra[iCombInStudy][kKaon][kNeg]->Draw("same"); | |
1415 | gStar[4]->Draw("p"); | |
1416 | hSpectra[iCombInStudy][kProton][kNeg]->Draw("same"); | |
4c0d7bc7 | 1417 | |
1418 | ||
1419 | ||
1420 | TCanvas * c2 = new TCanvas("cRatioToStarPos","cRatioToStarPos"); | |
1421 | hempty->Draw(); | |
1422 | leg = new TLegend(0.6510067,0.1853147,0.8892617,0.4178322,"Positive","brNDC"); | |
1423 | leg->SetBorderSize(0); | |
1424 | leg->SetLineColor(1); | |
1425 | leg->SetLineStyle(1); | |
1426 | leg->SetLineWidth(1); | |
1427 | leg->SetFillColor(0); | |
1428 | leg->SetFillStyle(1001); | |
c4fe358b | 1429 | |
1430 | TCanvas * c2Comp = new TCanvas("cCompToStarPos","cCompToStarPos"); | |
1431 | c2Comp->SetLogy(); | |
1432 | hempty2->Draw(); | |
1433 | ||
1434 | c2->cd(); | |
4c0d7bc7 | 1435 | // TGraphErrors * g ; |
ee08b77d | 1436 | g = AliBWTools::DivideGraphByHisto(gStar[1],hSpectra[iCombInStudy][kPion][kPos],1); |
4c0d7bc7 | 1437 | g->SetMarkerStyle(kFullCircle); |
1438 | g->SetMarkerColor(kBlack); | |
1439 | g->Draw("p"); | |
1440 | leg->AddEntry(g,"#pi^{+}","lp"); | |
ee08b77d | 1441 | g = AliBWTools::DivideGraphByHisto(gStar[3],hSpectra[iCombInStudy][kKaon][kPos],1); |
4c0d7bc7 | 1442 | g->SetMarkerStyle(kOpenCircle); |
1443 | g->SetMarkerColor(kRed); | |
1444 | g->Draw("p"); | |
1445 | leg->AddEntry(g,"K^{+}","lp"); | |
ee08b77d | 1446 | g = AliBWTools::DivideGraphByHisto(gStar[5],hSpectra[iCombInStudy][kProton][kPos],1); |
4c0d7bc7 | 1447 | g->SetMarkerStyle(kOpenSquare); |
1448 | g->SetMarkerColor(kBlue); | |
1449 | g->Draw("p"); | |
1450 | leg->AddEntry(g,"p","lp"); | |
c4fe358b | 1451 | leg->Draw(); |
4c0d7bc7 | 1452 | |
1453 | ||
c4fe358b | 1454 | c2Comp->cd(); |
1455 | gStar[1]->Draw("p"); | |
1456 | hSpectra[iCombInStudy][kPion][kPos]->Draw("same"); | |
1457 | gStar[3]->Draw("p"); | |
1458 | hSpectra[iCombInStudy][kKaon][kPos]->Draw("same"); | |
1459 | gStar[5]->Draw("p"); | |
1460 | hSpectra[iCombInStudy][kProton][kPos]->Draw("same"); | |
1461 | ||
4c0d7bc7 | 1462 | |
1463 | c1->Update(); | |
1464 | c2->Update(); | |
1465 | gSystem->ProcessEvents(); | |
1466 | c1->Print(TString(c1->GetName()) + ".eps"); | |
1467 | c2->Print(TString(c2->GetName()) + ".eps"); | |
1468 | ALICEWorkInProgress(c1,today.Data(),"#splitline{ALICE Preliminary}{Statistical Error Only}"); | |
1469 | ALICEWorkInProgress(c2,today.Data(),"#splitline{ALICE Preliminary}{Statistical Error Only}"); | |
1470 | c1->Update(); | |
1471 | c2->Update(); | |
1472 | c1->Print(TString(c1->GetName()) + ".png"); | |
1473 | c2->Print(TString(c2->GetName()) + ".png"); | |
1474 | ||
1475 | ||
1476 | ||
1477 | ||
1478 | } | |
1479 | ||
1480 | ||
1481 | ||
1482 | void DrawRatios() { | |
1483 | ||
1484 | // Draws ratios of combined spectra | |
1485 | ||
1486 | // Compute ratios | |
1487 | TH1F * hPosNegRatio[kNPart]; | |
1488 | ||
1489 | for(Int_t ipart = 0; ipart < kNPart; ipart++){ | |
ee08b77d | 1490 | hPosNegRatio[ipart] = (TH1F*) hSpectra[iCombInStudy][ipart][kPos]->Clone(); |
1491 | hPosNegRatio[ipart]->Divide(hSpectra[iCombInStudy][ipart][kNeg]); | |
4c0d7bc7 | 1492 | hPosNegRatio[ipart]->SetYTitle(TString(partLabel[ipart][kPos])+"/"+partLabel[ipart][kNeg]); |
1493 | hPosNegRatio[ipart]->SetMinimum(0.5); | |
1494 | hPosNegRatio[ipart]->SetMaximum(1.5); | |
1495 | } | |
1496 | ||
ee08b77d | 1497 | TH1F * hKPiRatio = (TH1F*) hSpectra[iCombInStudy][kKaon][kPos]->Clone(); |
1498 | hKPiRatio->Add(hSpectra[iCombInStudy][kKaon][kNeg]); | |
1499 | TH1F * htmp = (TH1F*) hSpectra[iCombInStudy][kPion][kPos]->Clone(); | |
1500 | htmp->Add(hSpectra[iCombInStudy][kPion][kNeg]); | |
4c0d7bc7 | 1501 | hKPiRatio->Divide(htmp); |
1502 | hKPiRatio->SetYTitle("(K^{+}+K^{-})/(#pi^{+}+#pi^{-})"); | |
1503 | ||
1504 | TH1F * hKPiRatioMC[kNTunes]; | |
5f2ca1e6 | 1505 | if(showMC){ |
1506 | for(Int_t itune = 0; itune < kNTunes; itune++){ | |
1507 | hKPiRatioMC[itune] = (TH1F*) hSpectraMC[itune][kKaon][kPos]->Clone(); | |
1508 | hKPiRatioMC[itune]->Add(hSpectraMC[itune][kKaon][kNeg]); | |
1509 | TH1F * htmp = (TH1F*) hSpectraMC[itune][kPion][kPos]->Clone(); | |
1510 | htmp->Add(hSpectraMC[itune][kPion][kNeg]); | |
1511 | hKPiRatioMC[itune]->Divide(htmp); | |
1512 | hKPiRatioMC[itune]->SetYTitle("(K^{+}+K^{-})/(#pi^{+}+#pi^{-})"); | |
1513 | } | |
4c0d7bc7 | 1514 | } |
1515 | ||
1516 | ||
1517 | ||
ee08b77d | 1518 | TH1F * hPPiRatio = (TH1F*) hSpectra[iCombInStudy][kProton][kPos]->Clone(); |
1519 | hPPiRatio->Add(hSpectra[iCombInStudy][kProton][kNeg]); | |
1520 | htmp = (TH1F*) hSpectra[iCombInStudy][kPion][kPos]->Clone(); | |
1521 | htmp->Add(hSpectra[iCombInStudy][kPion][kNeg]); | |
4c0d7bc7 | 1522 | hPPiRatio->Divide(htmp); |
1523 | hPPiRatio->SetYTitle("(p+#bar{p})/(#pi^{+}+#pi^{-})"); | |
1524 | ||
5f2ca1e6 | 1525 | if(showMC){ |
1526 | TH1F * hPPiRatioMC[kNTunes]; | |
1527 | for(Int_t itune = 0; itune < kNTunes; itune++){ | |
1528 | hPPiRatioMC[itune] = (TH1F*) hSpectraMC[itune][kProton][kPos]->Clone(); | |
1529 | hPPiRatioMC[itune]->Add(hSpectraMC[itune][kProton][kNeg]); | |
1530 | TH1F * htmp = (TH1F*) hSpectraMC[itune][kPion][kPos]->Clone(); | |
1531 | htmp->Add(hSpectraMC[itune][kPion][kNeg]); | |
1532 | hPPiRatioMC[itune]->Divide(htmp); | |
1533 | hPPiRatioMC[itune]->SetYTitle("(p+#bar{p})/(#pi^{+}+#pi^{-})"); | |
1534 | } | |
4c0d7bc7 | 1535 | } |
1536 | ||
1537 | // Draw | |
1538 | // TH2F * hempty = new TH2F(TString("hempty"),"hempty",100,0.,1.5, 100, 0.001,1.8); | |
1539 | // hempty->SetXTitle("p_{t} (GeV/c)"); | |
1540 | // hempty->SetYTitle(""); | |
1541 | ||
1542 | // tmp: overlay levi fits | |
1543 | AliBWFunc * fm2 = new AliBWFunc; | |
1544 | fm2->SetVarType(AliBWFunc::kdNdpt); | |
1545 | TF1 * fLevi[kNPart] [kNCharge]; | |
1546 | fLevi[kPion] [kPos] = fm2->GetLevi (mass[0], 0.1243, 7.614785, 1.524167, "fLeviPiPlus"); | |
1547 | fLevi[kKaon] [kPos] = fm2->GetLevi (mass[1], 0.1625, 5.869318, 0.186361, "fLeviKPlus"); | |
1548 | fLevi[kProton][kPos] = fm2->GetLevi (mass[2], 0.1773, 6.918065, 0.086389, "fLeviPPlus"); | |
1549 | fLevi[kPion] [kNeg] = fm2->GetLevi (mass[0], 0.1267, 7.979582, 1.515908, "fLeviPiNeg"); | |
1550 | fLevi[kKaon] [kNeg] = fm2->GetLevi (mass[1], 0.1721, 6.927956, 0.191140, "fLeviKNeg"); | |
1551 | fLevi[kProton][kNeg] = fm2->GetLevi (mass[2], 0.1782, 8.160362, 0.082091, "fLeviPNeg"); | |
1552 | for(Int_t ipart = 0; ipart < kNPart; ipart++){ | |
1553 | for(Int_t icharge = 0; icharge < kNCharge; icharge++){ | |
1554 | fLevi[ipart][icharge]->SetRange(0,4); | |
1555 | } | |
1556 | } | |
1557 | ||
1558 | ||
1559 | for(Int_t ipart = 0; ipart < kNPart; ipart++){ | |
b7018bc7 | 1560 | TString detName = detFlag[iCombInStudy]; |
1561 | detName.ReplaceAll(" ", "_"); | |
1562 | detName.ReplaceAll("+", ""); | |
1563 | ||
5f2ca1e6 | 1564 | TCanvas * c1 = new TCanvas(TString("cRatio_")+detName+TString("_")+partFlag[ipart], TString("cRatio_")+detName+partFlag[ipart]); |
1565 | c1->SetGridy(); | |
4c0d7bc7 | 1566 | hPosNegRatio[ipart]->Draw(); |
1567 | TF1 * fRatio = new TF1 (TString("fRatio")+partFlag[ipart], TString(fLevi[ipart][kPos]->GetName())+"/"+fLevi[ipart][kNeg]->GetName()); | |
1568 | // fRatio->Draw("same"); | |
1569 | fRatio->SetRange(0,5); | |
1570 | if (doPrint) { | |
1571 | c1->Update(); | |
1572 | gSystem->ProcessEvents(); | |
b7018bc7 | 1573 | c1->Print(TString(c1->GetName()) + ".png"); |
4c0d7bc7 | 1574 | } |
1575 | ||
1576 | } | |
1577 | ||
1578 | ||
1579 | TCanvas * c2 = new TCanvas(TString("cRatio_KPi"),TString("cRatio_KPi")); | |
5f2ca1e6 | 1580 | c2->SetGridy(); |
4c0d7bc7 | 1581 | hKPiRatio->Draw(); |
1582 | TLegend * lMC = new TLegend(0.526846, 0.18007, 0.887584,0.407343); | |
1583 | lMC->SetFillColor(kWhite); | |
1584 | ||
5f2ca1e6 | 1585 | if(showE735){ |
1586 | gROOT->LoadMacro("GetE735Ratios.C"); | |
1587 | GetE735Ratios(0,0)->Draw("EX0,same"); | |
1588 | GetE735Ratios(0,1)->Draw("EX0,same"); | |
1589 | GetE735Ratios(0,2)->Draw("EX0,same"); | |
1590 | GetE735Ratios(0,3)->Draw("EX0,same"); | |
1591 | } | |
4c0d7bc7 | 1592 | hKPiRatio->SetMarkerStyle(20); |
1593 | hKPiRatio->Draw("same"); | |
1594 | ||
5f2ca1e6 | 1595 | if(showMC){ |
1596 | for(Int_t itune = 0; itune < kNTunes; itune++){ | |
1597 | lMC->AddEntry(hKPiRatioMC[itune],mcTuneName[itune]); | |
1598 | hKPiRatioMC[itune]->SetLineWidth(2); | |
1599 | hKPiRatioMC[itune]->Draw("same,chist"); | |
1600 | } | |
4c0d7bc7 | 1601 | |
5f2ca1e6 | 1602 | lMC->Draw(); |
1603 | } | |
4c0d7bc7 | 1604 | |
5f2ca1e6 | 1605 | if(showE735){ |
1606 | TLegend * l = new TLegend( 0.1879, 0.68, 0.54,0.92); | |
1607 | l->SetFillColor(kWhite); | |
1608 | l->AddEntry(hKPiRatio, "ALICE, #sqrt{s} = 900 GeV","lpf"); | |
1609 | l->AddEntry(GetE735Ratios(0,0), "E735, #sqrt{s} = 300 GeV","lpf"); | |
1610 | l->AddEntry(GetE735Ratios(0,1), "E735, #sqrt{s} = 540 GeV","lpf"); | |
1611 | l->AddEntry(GetE735Ratios(0,2), "E735, #sqrt{s} = 1000 GeV","lpf"); | |
1612 | l->AddEntry(GetE735Ratios(0,3), "E735, #sqrt{s} = 1800 GeV","lpf"); | |
1613 | l->Draw(); | |
1614 | } | |
4c0d7bc7 | 1615 | |
1616 | ||
1617 | TCanvas * c3 = new TCanvas(TString("cRatio_PPi"),TString("cRatio_PPi")); | |
5f2ca1e6 | 1618 | c3->SetGridy(); |
4c0d7bc7 | 1619 | hPPiRatio->Draw(); |
1620 | hPPiRatio->SetMaximum(0.39); | |
5f2ca1e6 | 1621 | if(showMC){ |
1622 | lMC = new TLegend(0.526846, 0.18007, 0.887584,0.407343); | |
1623 | lMC->SetFillColor(kWhite); | |
1624 | ||
1625 | for(Int_t itune = 0; itune < kNTunes; itune++){ | |
1626 | lMC->AddEntry(hKPiRatioMC[itune],mcTuneName[itune]); | |
1627 | hKPiRatioMC[itune]->SetLineWidth(2); | |
1628 | hKPiRatioMC[itune]->Draw("same,chist"); | |
1629 | } | |
4c0d7bc7 | 1630 | |
5f2ca1e6 | 1631 | lMC->Draw(); |
1632 | } | |
4c0d7bc7 | 1633 | |
1634 | if (doPrint) { | |
1635 | c2->Update(); | |
1636 | gSystem->ProcessEvents(); | |
5f2ca1e6 | 1637 | c2->Print(TString(c2->GetName()) + ".png"); |
4c0d7bc7 | 1638 | c2->Print(TString(c2->GetName()) + ".eps"); |
1639 | c3->Update(); | |
1640 | gSystem->ProcessEvents(); | |
5f2ca1e6 | 1641 | c3->Print(TString(c3->GetName()) + ".png"); |
4c0d7bc7 | 1642 | c3->Print(TString(c3->GetName()) + ".eps"); |
1643 | } | |
1644 | ||
1645 | ||
1646 | } | |
1647 | ||
c4fe358b | 1648 | |
1649 | void Help() { | |
1650 | ||
1651 | cout << "Macro: void CombineSpectra(Int_t analysisType=kDoFits, Int_t fitFuncID = kFitLevi) " << endl; | |
1652 | cout << "" << endl; | |
1653 | ||
1654 | cout << "Possible Arguments" << endl; | |
1655 | cout << "- analysisType:" << endl; | |
eea641de | 1656 | cout << " kDoFits: Fit Combined Spectra " << endl; |
1657 | cout << " kDoRatios: Particle ratios K/pi and p/pi" << endl; | |
1658 | cout << " kDoSuperposition: Compare different detectors (superimpose and ratios)" << endl; | |
1659 | cout << " kDoCompareStar: Compare combined spectra to star results" << endl; | |
1660 | cout << " kDoDrawWithModels: Compare combined spectra and models" << endl; | |
1661 | cout << " kDoHelp: This help" << endl; | |
1662 | cout << "- fitFuncID, function used to extrapolate and compute yields" << endl; | |
1663 | cout << " An analitic fit function [kFitLevi, kFitUA1, kFitPowerLaw]" << endl; | |
1664 | cout << " Or a shape from a MC moder [kFitPhojet, kFitAtlasCSC, kFitCMS6D6T, kFitPerugia0]" << endl; | |
1665 | cout << " Which is fitted to the data at low pt and used to extrapolate at low pt" | |
c4fe358b | 1666 | |
1667 | ||
1668 | } |