]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGLF/ThermalFits/PlotRatiosForQM14.C
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Merge branch 'feature-movesplit'
[u/mrichter/AliRoot.git] / PWGLF / ThermalFits / PlotRatiosForQM14.C
CommitLineData
9d84731e 1#if !defined (__CINT__) || (defined(__MAKECINT__))
2#include <iostream>
3#include "TClonesArray.h"
4#include "AliParticleYield.h"
5#include "TH1F.h"
6#include "TCanvas.h"
046fa53e 7#include "TStyle.h"
9d84731e 8#include <fstream>
046fa53e 9#include "TLatex.h"
10#include "TLegend.h"
11#include "TList.h"
12#include "TF1.h"
13#include "AliPWGHistoTools.h"
14#include "TGraphErrors.h"
15#include "TMath.h"
79936afc 16#include "TDatabasePDG.h"
69f3eeda 17#include "TH2F.h"
18#include "TSystem.h"
15fccaa0 19#include "TPaveText.h"
b926af7d 20#include <map>
21#include <fstream>
22#include <istream>
23#include "TMarker.h"
24#include "TObjString.h"
25#include "TLegendEntry.h"
046fa53e 26
9d84731e 27#endif
28
046fa53e 29// Plots ratios for QM and saves input files for thermal models
30
2f73d64b 31#if 1 // DUMMY IFDEF USED TO HIDE PREAMBLE in EMACS
046fa53e 32
9d84731e 33enum MyParticles { kPDGPi = 211, kPDGK = 321, kPDGProton = 2212, kPDGKS0 = 310, kPDGLambda=3122, kPDGXi=3312,kPDGOmega=3334,kPDGPhi=333,kPDGKStar=313,kPDGDeuteron=1000010020,kPDGHE3 = 1000020030, kPDGHyperTriton = 1010010030, kPDGSigmaStarPlus=3224,kPDGSigmaStarMinus=3114,kPDGXiStar=3324};
34
046fa53e 35typedef enum {kStatError, kSystError, kTotalError} myerror_t;
36
37TH1F * GetHistoRatios(TClonesArray * arr, Int_t system, Float_t energy, TString centrality, const char * histotitle, Int_t icolor, Int_t imarker = kOpenSquare, Int_t errorsType = kTotalError, Float_t shift = 0) ;
b926af7d 38TH1F * GetHistoYields(TClonesArray * arr, Int_t system, Float_t energy, TString centrality, const char * histotitle, Int_t icolor, Int_t imarker = kOpenSquare, Int_t errorsType = kTotalError, Float_t shift = 0) ;
9d84731e 39void PrepareThermalModelsInputFiles(TClonesArray * arr, Int_t system, Float_t energy, TString centrality, Bool_t separateCharges=0) ;
046fa53e 40void SetStyle(Bool_t graypalette=0) ;
b926af7d 41TLegend * NewLegendQM(Double_t x1, Double_t y1, Double_t x2, Double_t y2, Bool_t isYield = 0) ;
42void DrawRatio(TString what, Bool_t isYield = kFALSE, Double_t shiftloc=0.);
15fccaa0 43
69f3eeda 44void DrawFrame(Bool_t yields = 0) ;
b926af7d 45void DrawExtrapolatedSymbolsAndLegendPbPb0010() ;
46void DrawExtrapolatedSymbolsAndLegendpPb0005() ;
2f73d64b 47void DrawMarkerKStarNoFit(Bool_t plotLegend = 0) ;
48void DrawMarkerNucleiNoFit() ;
49void DrawExtrapNotInFitpPb0005(Bool_t drawExtrap = 1) ;
50
51void DrawExtrapolatedSymbolsYieldsPbPb0010(Double_t x1=0.144578, Double_t y1=0.408249, Double_t x2=0.351406, Double_t y2=0.542403, Bool_t plotExtraploatedLegend=1);
52Float_t shiftRatioDataModel = 0;
eee43f2a 53Double_t GetGraphRatioAndStdDev(TGraphErrors * gModel, TGraphErrors * &gRatio, TGraphErrors *&gStdDev) ;
54TString particlesToExcludeFromChi2 = "";// The above method recomputes the chi2. This string is used to store values to be excluded from this calculation. Each PDG cocde has to be enclosed in [..]. This is obviously not efficient as it involves many string conversions. But efficiency is not an issue here.
4910214c 55
b926af7d 56
046fa53e 57void LoadArrays() ;
69f3eeda 58//void AddLabel(Float_t x, Float_t y, TString text);
59void myLatexDraw(TLatex *currentLatex, Float_t currentSize=0.5, Int_t currentColor=1);
60void myPaveSetup(float rRatio=0, float rRange3=0, float rRange5=0,
61 int rFillColor=0);
62void myPadSetUp(TPad *currentPad);
b926af7d 63TGraphErrors* PlotThermusYields(const char * filename, Int_t color, Int_t lineStyle,
4910214c 64 const char * tag);
b926af7d 65TGraphErrors * PlotGSIYields(const char * fileName, Int_t color=kBlack, Int_t lineStyle = kSolid,
2f73d64b 66 const char * tag ="", Bool_t isPbPb = 1);
5ba62df2 67TGraphErrors* PlotFlorenceYields(const char * filename, Int_t color, Int_t lineStyle,
68 const char * tag) ;
4910214c 69
70void AddLineToThermalLegend(TObject * obj, TString line, const char * optFirst = "L");
69f3eeda 71
2f73d64b 72void SaveCanvas(const char * name) ;
73
74
69f3eeda 75// Ratios to be draw. Remember to change the labels in DrawFrame if you change this
76const Int_t nratio = 10;
69f3eeda 77Int_t num [nratio] = {kPDGK , kPDGProton , kPDGLambda , kPDGXi , kPDGOmega , kPDGDeuteron , kPDGHE3 , kPDGHyperTriton , kPDGPhi , kPDGKStar};
78Int_t denum[nratio] = {kPDGPi , kPDGPi , kPDGKS0 , kPDGPi , kPDGPi , kPDGProton , kPDGDeuteron , kPDGPi , kPDGK , kPDGK};
15fccaa0 79Int_t isSum[nratio] = {1 , 1 , 1 , 1 , 1 , 1 , 0 , 1 , 1 , 1 };
5ba62df2 80
81// const Int_t nratio = 10;
82// Int_t num [nratio] = {kPDGK , kPDGProton , kPDGLambda , kPDGXi , kPDGOmega , kPDGDeuteron , kPDGHE3 , kPDGHyperTriton , kPDGPhi , kPDGKStar};
83// Int_t denum[nratio] = {kPDGPi , kPDGPi , kPDGKS0 , kPDGPi , kPDGPi , kPDGProton , kPDGDeuteron , kPDGPi , kPDGK , kPDGK};
84// Int_t isSum[nratio] = {1 , 1 , 1 , 1 , 1 , 1 , 0 , 1 , 1 , 1 };
85
15fccaa0 86const char * ratiosLabels[] = {"#frac{K^{+}+K^{-}}{#pi^{+}+#pi^{-}}",
87 "#frac{p+#bar{p}}{#pi^{+}+#pi^{-}}",
b926af7d 88 "#frac{2#Lambda}{K_{S}^{0}}",
15fccaa0 89 "#frac{#Xi^{-}+#Xi^{+}}{#pi^{+}+#pi^{-}}",
90 "#frac{#Omega^{-}+#Omega^{+}}{#pi^{+}+#pi^{-}}",
91 "#frac{d}{p+#bar{p}}",
92 "#frac{{}^{3}He }{d}",
93 "#frac{{}^{3}_{#Lambda}H+{}^{3}_{#Lambda}#bar{H} }{#pi^{+}+#pi^{-}}",
94 "#frac{#phi}{K^{+}+K^{-}}",
95 "#frac{K*+#bar{K}*}{K^{+}+K^{-}}",};
5ba62df2 96//static const Double_t scale[] = {1 , 3 , 0.5 , 30 , 250 , 50 , 100 , 4e5 , 2 , 1 };
97static const Double_t scale[] = {1 , 3 , 0.5 , 80 , 1000 , 50 , 100 , 4e5 , 2 , 1 };
69f3eeda 98//static const Double_t scale[] = {1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,};
b926af7d 99const Int_t npart = 12;
100Int_t particleYields [npart] = {kPDGPi ,kPDGK ,kPDGKS0, kPDGKStar, kPDGPhi, kPDGProton , kPDGLambda , kPDGXi , kPDGOmega , kPDGDeuteron, kPDGHyperTriton, kPDGHE3 };
101Int_t isSumYields[npart] = {1 ,1 ,0 , 1 , 0 , 1 ,1 ,1 ,1 ,0 , 1 , 0 };
102//Int_t isSumInputFiles[npart] = {1 ,1 ,0 , 1 , 0 , 1 ,1 ,1 ,1 ,0 , 1 , 0 };
103const char * yieldsLabels[] = {"#frac{#pi^{+}+#pi^{-}}{2}",
104 "#frac{K^{+}+K^{-}}{2}",
105 "K_{S}^{0}",
eee43f2a 106 "#frac{K*+K*}{2}",
b926af7d 107 "#phi",
108 "#frac{p+#bar{p}}{2}",
109 "#Lambda",
110 "#frac{#Xi^{-}+#Xi^{+}}{2}",
111 "#frac{#Omega^{-}+#Omega^{+}}{2}",
112 "d",
113 "#frac{{}^{3}_{#Lambda}H+{}^{3}_{#Lambda}#bar{H}}{2}",
114 "{}^{3}He",
115};
9d84731e 116
69f3eeda 117//
b926af7d 118Int_t markerNoFit = 28;
119Int_t markerExtrap = 27;
15fccaa0 120Double_t maxy = 0.5;
9d84731e 121
eee43f2a 122
046fa53e 123// Data arrays;
124TClonesArray *arrPbPb=0, *arrpp7=0, *arrpPb=0, * arrpp276=0, * arrpp900=0, * arrThermus=0;
125TClonesArray *arrSTARPbPb=0, *arrPHENIXPbPb=0, *arrBRAHMSPbPb=0;
126TClonesArray *arrSTARpp =0, *arrPHENIXpp=0;
127
69f3eeda 128//const Double_t *scaleRatios = 0;
4c7bd9d1 129Double_t *correlatedUnc = 0;
2f73d64b 130//Double_t correlatedUncLocalPbPb[14] = {0.0424 , 0.0424 , 0.041 , 0 , 0 , 0.0424 , 0.0424 , 0 , 0.05 , 0.05 };
131Double_t correlatedUncLocalPbPbOnlyKStarPhi[14] = {0. , 0. , 0. , 0 , 0 , 0. , 0. , 0 , 0.05 , 0.05 };
4c7bd9d1 132Double_t correlatedUncLocalPP [14] = {0.0424 , 0.0424 , 0.041 , 0 , 0 , 0.0424 , 0.0424 , 0 , 0.0424 , 0.0424 };
133Double_t correlatedUncZero[14] = {0,0,0,0,0,0,0,0,0,0,0,0,0,0};
134
69f3eeda 135TCanvas *myCan = 0;
b926af7d 136TPad *myPadStdDev =0;
eee43f2a 137TPad *myPadRatio =0;
b926af7d 138TPad *myPadHisto =0;
139TPad *myPadLabel =0;
140TLegend * legThermal = 0;
2f73d64b 141#endif
142
5ba62df2 143Bool_t saveCanvas = 1; // if true, the canvas is saved and copied to the analysis note folder
4c7bd9d1 144
9d84731e 145TClonesArray * PlotRatiosForQM14() {
ad431d97 146#if !(!defined (__CINT__) || (defined(__MAKECINT__)))
046fa53e 147 LoadLibs();
ad431d97 148#endif
9d84731e 149
046fa53e 150 //
151 LoadArrays();
ad431d97 152
046fa53e 153 // Uncomment stuff in this section to save the inputs for thermal models
eee43f2a 154 //#define SAVE_INPUT_THERMAL_MODEL
046fa53e 155#ifdef SAVE_INPUT_THERMAL_MODEL
8f3e9877 156 PrepareThermalModelsInputFiles(arrPbPb, AliParticleYield::kCSPbPb, 2760, "V0M0010", /*separateCharges*/1);
4910214c 157 // PrepareThermalModelsInputFiles(arrpp7, AliParticleYield::kCSpp, 7000, "", /*separateCharges*/1);
158 // PrepareThermalModelsInputFiles(arrpPb, AliParticleYield::kCSpPb, 5020, "V0A0005", /*separateCharges*/1);
159 // PrepareThermalModelsInputFiles(arrpPb, AliParticleYield::kCSpPb, 5020, "V0A2040", /*separateCharges*/1);
160 // PrepareThermalModelsInputFiles(arrpPb, AliParticleYield::kCSpPb, 5020, "V0A6080", /*separateCharges*/1);
161 // PrepareThermalModelsInputFiles(arrPbPb, AliParticleYield::kCSPbPb, 2760, "V0M6080", /* separateCharges*/1);
162 // PrepareThermalModelsInputFiles(arrPbPb, AliParticleYield::kCSPbPb, 2760, "V0M2040", /*separateCharges*/1);
69f3eeda 163
b926af7d 164 PrepareThermalModelsInputFiles(arrPbPb, AliParticleYield::kCSPbPb, 2760, "V0M0010", /*separateCharges*/0);
4910214c 165 // PrepareThermalModelsInputFiles(arrpp7, AliParticleYield::kCSpp, 7000, "", /*separateCharges*/0);
166 // PrepareThermalModelsInputFiles(arrpPb, AliParticleYield::kCSpPb, 5020, "V0A0005", /*separateCharges*/0);
167 // PrepareThermalModelsInputFiles(arrpPb, AliParticleYield::kCSpPb, 5020, "V0A2040", /*separateCharges*/0);
168 // PrepareThermalModelsInputFiles(arrpPb, AliParticleYield::kCSpPb, 5020, "V0A6080", /*separateCharges*/0);
169 // PrepareThermalModelsInputFiles(arrPbPb, AliParticleYield::kCSPbPb, 2760, "V0M6080", /*separateCharges*/0);
170 // PrepareThermalModelsInputFiles(arrPbPb, AliParticleYield::kCSPbPb, 2760, "V0M2040", /*separateCharges*/0);
9d84731e 171
046fa53e 172 return 0;
173#endif
9d84731e 174
2f73d64b 175
046fa53e 176 SetStyle();
ad431d97 177
5ba62df2 178 DrawRatio("allpp");
179 // DrawRatio("allppWithRHIC");
15fccaa0 180 // DrawRatio("PbPbWithPP7TeV");
2f73d64b 181 //DrawRatio("allsyst");
182 //DrawRatio("PbPb6080andpPb0005");
183 // DrawRatio("pp_vsRHIC");
4910214c 184 // DrawRatio("PbPb_vsRHIC");
8f3e9877 185 //DrawRatio("aliceall");
ad431d97 186
b926af7d 187
188 // Yields and FITS
5ba62df2 189 // maxy=2000;
2f73d64b 190
5ba62df2 191 // DrawRatio("fit_ReferenceFit_PbPb0010", 1);
8f3e9877 192 // DrawRatio("fit_ReferenceFit_GSIONLY_PbPb0010", 1);
193 // DrawRatio("fit_ReferenceFit_GSITHERMUS_PbPb0010",1);
5ba62df2 194 // DrawRatio("fitSHARE_NoPionsNoProtons_PbPb0010",1);
4910214c 195 // DrawRatio("fitGSI_NoPionsNoProtons_PbPb0010", 1);
5ba62df2 196 // DrawRatio("fitShare_All_PbPb0010", 1);
2f73d64b 197
4910214c 198 // DrawRatio("fitShareWithWithoutNuclei_PbPb0010", 1);
2f73d64b 199 // maxy=200;
200 // DrawRatio("fitGSI_PbPb6080", 1);
5ba62df2 201 // maxy=150;
202 // DrawRatio("fitGSI_PbPb2040", 1);
203 // maxy = 60;
8f3e9877 204 // DrawRatio("fitThermus_GammaSFree_pPb0005");
5ba62df2 205 // DrawRatio("fitShare_pPb0005");
206 // DrawRatio("fitShare_pPb0005_NoOmega", 1);
207 // maxy=20;
208 // DrawRatio("fitThermus_GammaSFree_pPb2040");
209 // maxy=9;
210 // DrawRatio("fitThermus_GammaSFree_pPb6080");
211 // maxy=9;
212 // DrawRatio("fitGSI_pp");
213 // DrawRatio("fitFlorence_pp");
2f73d64b 214
15fccaa0 215 // NewLegendQM();
046fa53e 216 return arrPbPb;
9d84731e 217}
218
046fa53e 219TH1F * GetHistoRatios(TClonesArray * arr, Int_t system, Float_t energy, TString centrality, const char * histotitle, Int_t icolor, Int_t imarker, Int_t errorType, Float_t shift) {
ad431d97 220 // FIXME: THIS SHOULD BE REVIEWED TO MAKE SURE THE PLOTS ARE LABELLED CORRECTLY
9d84731e 221
4c7bd9d1 222
69f3eeda 223 // scaleRatios = scale;
046fa53e 224 TH1F * h = new TH1F(Form("hRatio_%d_%0.0f_%s_%d",system,energy,centrality.Data(),errorType), histotitle, nratio, 1+shift, nratio+1+shift);
225
226 TClonesArray * arrRatios = new TClonesArray ("AliParticleYield");// We save to disk the precomputed ratios
227 Int_t iratioArr = 0;// Index used only to add particles to the array above
9d84731e 228
229 // Double_t isSum = -1; // if this is -1, then the sum criterion is ignored
230 for(Int_t iratio = 1; iratio <= nratio; iratio++){
79936afc 231 std::cout << "******** " << num[iratio-1] << "/" << denum[iratio-1]<< " ("<<isSum[iratio-1]<<")*******" << std::endl ;
9d84731e 232
79936afc 233 AliParticleYield * ratio = AliParticleYield::FindRatio(arr, num[iratio-1], denum[iratio-1], system, energy, centrality,isSum[iratio-1]);
234 if(ratio)
235 {
236 ratio = new AliParticleYield(*ratio); // We need to clone it to avoid a mess if we need to use this particle again later (e.g. double scalings)
237 std::cout << " " ;
238 ratio->Print("short");
239 }
240
9d84731e 241 if(!ratio) {
242 // If the ratio is not found, try to build it!
79936afc 243 std::cout << " Looking for " << num[iratio-1] << " ("<<isSum[iratio-1]<<")"<<std::endl;
9d84731e 244 AliParticleYield * part1 = AliParticleYield::FindParticle(arr, num[iratio-1], system, energy, centrality, isSum[iratio-1]);
69f3eeda 245 if(part1) {
246 part1 = new AliParticleYield(*part1); // We need to clone it to avoid a mess if we need to use this particle again later
247 if(isSum[iratio-1] && part1->IsTypeAverage()) {
248 std::cout << "Sum requested, found average, scaling x2" << std::endl;
249 part1->Scale(2.);
250 }
251 }
ad431d97 252 // Try with the !sum, if part 1 is not found
253 if(!part1) {
79936afc 254 std::cout << " Looking for " << num[iratio-1] << " ("<<!isSum[iratio-1]<<")"<<std::endl;
9d84731e 255 part1 = AliParticleYield::FindParticle(arr, num[iratio-1], system, energy, centrality,!isSum[iratio-1]);
79936afc 256 if(part1)
257 {
258 part1 = new AliParticleYield(*part1); // We need to clone it to avoid a mess if we need to use this particle again later
259 // If the sum was requested, try to recover it!
4c7bd9d1 260 if(isSum[iratio-1]) {
79936afc 261 std::cout << " Looking for " << -num[iratio-1] <<std::endl;
15fccaa0 262 // If it's lambda and ALICE, use 2L instead of L + Lbar // FIXME: do the same for deuterons?
263 if((num[iratio-1] == kPDGLambda || num[iratio-1] == kPDGDeuteron) && energy > 300) {
264 std::cout << " It's Lambda or deuteron ALICE: Scaling x2 " << std::endl;
265 part1->Print();
266 part1->Scale(2.);
267 } else {
268 AliParticleYield * part1_bar = AliParticleYield::FindParticle(arr, -num[iratio-1], system, energy, centrality,0);
269 if(part1 && part1_bar) {
270 std::cout << "Adding " << part1_bar->GetPartName() << " " << part1->GetPartName() << std::endl;
271 part1 = AliParticleYield::Add(part1, part1_bar);
272
273 } else if (TDatabasePDG::Instance()->GetParticle(-num[iratio-1])){ // Before scaling x2 check it it makes sense (antiparticle exists)
274 std::cout << "WARNING: Sum requested but not found, scaling x2 " << part1->GetName() << std::endl;
275 part1->Scale(2);
276 }
79936afc 277 }
278 } else if(part1) { // if we are here, it means the sum was *not* requested (isSum=0), but we found something with (!isSum) = 1
279 // if the not sum was requested, but the sum is found, divide by 2 so that it is comparable
280 std::cout << "WARNING: Using sum/2 for " << part1->GetName() << std::endl;
046fa53e 281
79936afc 282 part1->Scale(0.5);
283 }
ad431d97 284 }
285
9d84731e 286 }
15fccaa0 287
288
79936afc 289 std::cout << " Looking for " << denum[iratio-1] << " ("<<isSum[iratio-1]<<")"<<std::endl;
9d84731e 290 AliParticleYield * part2 = AliParticleYield::FindParticle(arr, denum[iratio-1], system, energy, centrality,isSum[iratio-1]);
69f3eeda 291 if(part2) {
292 part2 = new AliParticleYield(*part2); // We need to clone it to avoid a mess if we need to use this particle again later
293 if(isSum[iratio-1] && part2->IsTypeAverage()) {
294 std::cout << "Sum requested, found average, scaling x2" << std::endl;
295 part2->Scale(2.);
296 }
297 }
9d84731e 298 if(!part2) {// Try with the !sum, if part 2 is not found
79936afc 299 std::cout << " Looking for " << denum[iratio-1] << " ("<<!isSum[iratio-1]<<")"<<std::endl;
9d84731e 300 part2 = AliParticleYield::FindParticle(arr, denum[iratio-1], system, energy, centrality,!isSum[iratio-1]);
79936afc 301 if(part2)
302 {
303 part2 = new AliParticleYield(*part2); // We need to clone it to avoid a mess if we need to use this particle again later
4c7bd9d1 304 if(isSum[iratio-1]) {
79936afc 305 std::cout << " Looking for " << -denum[iratio-1] << std::endl;
306 AliParticleYield * part2_bar = AliParticleYield::FindParticle(arr, -denum[iratio-1], system, energy, centrality,0);
307 if(part2 && part2_bar){
308 std::cout << "Adding " << part2_bar->GetPartName() << " " << part2->GetPartName() << std::endl;
309 part2 = AliParticleYield::Add(part2, part2_bar);
4c7bd9d1 310 } else if (TDatabasePDG::Instance()->GetParticle(-denum[iratio-1])){ // Before scaling x2 check it it makes sense (antiparticle exists)
79936afc 311 std::cout << "WARNING: Sum requested but not found, scaling x2 " << part2->GetName() << std::endl;
312 part2->Scale(2);
313 }
314 } else if(part2){
315 // if the not sum was requested, but the sum is found, divide by 2 so that it is comparable
316 std::cout << "WARNING: Using sum/2 for " << part2->GetName() << std::endl;
317 part2->Scale(0.5);
318 }
319 }
ad431d97 320
9d84731e 321 }
4c7bd9d1 322 ratio = AliParticleYield::Divide(part1, part2, correlatedUnc[iratio-1], "YQ"); // Assume by that the systematics of part1 and part2 are uncorrelated.
69f3eeda 323
9d84731e 324 if(ratio) {
ad431d97 325 std::cout << "" << std::endl;
9d84731e 326 std::cout << "WARNING: building ratio " << num[iratio-1] <<"/"<<denum[iratio-1]<<": Check uncertainties!!" << std::endl;
ad431d97 327 ratio->Print("");
328 std::cout << "" << std::endl;
9d84731e 329 }
330 }
331 if(ratio){
ad431d97 332 ratio->Scale(scale[iratio-1]);
9d84731e 333 h->SetBinContent(iratio, ratio->GetYield());
046fa53e 334 if(errorType == kTotalError) {
335 h->SetBinError (iratio, ratio->GetTotalError(0/* 0 = no normalization error */));
336 } else if (errorType == kStatError) {
337 h->SetBinError (iratio, ratio->GetStatError());
338 } else if (errorType == kSystError) {
339 h->SetBinError (iratio, ratio->GetSystError());
340 } else {
341 std::cout << "ERROR: Unknown Error Type " << errorType << std::endl;
342 }
343
344 // h->GetXaxis()->SetBinLabel(iratio, Form("#splitline{%s}{%s}",Form("#times%2.2f", scale[iratio-1]), ratio->GetLatexName()));
345 h->GetXaxis()->SetBinLabel(iratio, ratio->GetLatexName());
346 new ((*arrRatios)[iratioArr++]) AliParticleYield(*ratio);
9d84731e 347 }
348 else {
349 h->GetXaxis()->SetBinLabel(iratio, Form("#frac{%d}{%d}",num[iratio-1], denum[iratio-1]));
350
351 }
79936afc 352 std::cout << "*** END OF " << num[iratio-1] << "/" << denum[iratio-1]<< " *******" << std::endl ;
353
9d84731e 354 }
355
046fa53e 356 h->GetYaxis()->SetRangeUser(0, maxy);
357 h->SetLineColor(icolor);
358 h->SetMarkerColor(icolor);
359 h->SetMarkerStyle(imarker);
69f3eeda 360 // the "if" avoids saving twice the same ratios
361 if(errorType == kSystError) AliParticleYield::SaveAsASCIIFile(arrRatios, TString("ratios_")+h->GetName());
9d84731e 362 return h;
363
364
365
366}
367
368void PrepareThermalModelsInputFiles(TClonesArray * arr, Int_t system, Float_t energy, TString centrality, Bool_t separateCharges) {
369 // If "Separate charges" is true, tries to dump both charges are dumped
370 TClonesArray * arrOut = new TClonesArray("AliParticleYield");
69f3eeda 371 TClonesArray * arrOutGSI = new TClonesArray("AliParticleYield"); // We add dummy lines to the GSI output file if needed!
f5b6253e 372;
9d84731e 373
374 Int_t ipartOut = 0; // Index for the array
69f3eeda 375 Int_t ipartOutGSI = 0; // Index for the array
376
9d84731e 377 for(Int_t ipart = 0; ipart < npart; ipart++){
378 if(!separateCharges) {
f5b6253e 379 AliParticleYield * part = AliParticleYield::FindParticle(arr, particleYields[ipart], system, energy, centrality, isSumYields[ipart]);
380 if(!part && isSumYields[ipart]) {
9d84731e 381 //Could not find the particle, but the sum was requested: build the sum!
f5b6253e 382 part = AliParticleYield::FindParticle(arr, particleYields[ipart], system, energy, centrality, 0);
383 AliParticleYield * part2 = AliParticleYield::FindParticle(arr, -particleYields[ipart], system, energy, centrality, 0);
69f3eeda 384 if(part2 && part) part = AliParticleYield::Add(part, part2);
385 else if(part) part->Scale(2.); // If we only found a particle, we can scale it by a factor 2.
9d84731e 386 else part = 0;
387 }
69f3eeda 388 // We want to save the average of particle and antiparticle in this case
389 if(part) {
f5b6253e 390 if(isSumYields[ipart] && !part->IsTypeAverage()) part->Scale(0.5); // If it's not already an average, but just a sum, divide by 2
69f3eeda 391 new((*arrOut )[ipartOut++]) AliParticleYield(*part);
392 new((*arrOutGSI)[ipartOutGSI++]) AliParticleYield(*part);
393 } else { // Add dummy particle to the GSI list
f5b6253e 394 new((*arrOutGSI)[ipartOutGSI++]) AliParticleYield(particleYields[ipart], system, energy, -10, -10, -10, -10, -10, -10, 5, 256, "DUMMY", 1, "ALICE");
69f3eeda 395 }
9d84731e 396 }
397 else {
f5b6253e 398 // ignore isSumYields and try to find both particleYields
9d84731e 399 Bool_t notFound = 0;
f5b6253e 400 AliParticleYield * part = AliParticleYield::FindParticle(arr, particleYields[ipart], system, energy, centrality, 0);
401 if(part) {
402 new((*arrOut)[ipartOut++]) AliParticleYield(*part);
403 new((*arrOutGSI)[ipartOutGSI++]) AliParticleYield(*part);
404 }
405 else {
406 // std::cout << "ADDING DUMMY part " << particleYields[ipart] << std::endl;
407 // new((*arrOutGSI)[ipartOutGSI++]) AliParticleYield(particleYields[ipart], system, energy, -10, -10, -10, -10, -10, -10, 5, 256, "DUMMY", 1, "ALICE");
408 notFound=1;
409 }
9d84731e 410 // Try to find antiparticle (-pdg code)
f5b6253e 411 part = AliParticleYield::FindParticle(arr, -particleYields[ipart], system, energy, centrality, 0);
69f3eeda 412 if(part) {
413 new((*arrOut)[ipartOut++]) AliParticleYield(*part);
414 new((*arrOutGSI)[ipartOutGSI++]) AliParticleYield(*part);
415 }
9d84731e 416 else if (notFound) {
417 // If neither charge was found, check if we at least have the sum
f5b6253e 418 part = AliParticleYield::FindParticle(arr, abs(particleYields[ipart]), system, energy, centrality, 1);
419 if (part) {
420 if(!part->IsTypeAverage()) part->Scale(0.5); // If it's a sum (not an average) divide by 2
69f3eeda 421 new((*arrOut)[ipartOut++]) AliParticleYield(*part);
422 new((*arrOutGSI)[ipartOutGSI++]) AliParticleYield(*part);
f5b6253e 423 if(TDatabasePDG::Instance()->GetParticle(-particleYields[ipart]) &&
424 (particleYields[ipart] != kPDGLambda) &&
425 (particleYields[ipart] != kPDGKStar)
426 ){// if only the sum was found, add a dummy entry to the
427 // GSI file, so that anton always has the same # of
428 // lines. However, we don't do this for the Lambda and
429 // KStar (always one of the charges or the average)
430 new((*arrOutGSI)[ipartOutGSI++]) AliParticleYield(-particleYields[ipart], system, energy, -10, -10, -10, -10, -10, -10, 5, 256, "DUMMY", 1, "ALICE");
431 }
69f3eeda 432 }
433 else {
f5b6253e 434 std::cout << "ADDING DUMMY sum " << particleYields[ipart] << std::endl;
435 new((*arrOutGSI)[ipartOutGSI++]) AliParticleYield(particleYields[ipart], system, energy, -10, -10, -10, -10, -10, -10, 5, 256, "DUMMY", 1, "ALICE");
436 if (particleYields[ipart]==kPDGHyperTriton) {
437 // If is the 3LH, add another one for the antiparticle
438 AliParticleYield(-particleYields[ipart], system, energy, -10, -10, -10, -10, -10, -10, 5, 256, "DUMMY", 1, "ALICE");
439 }
69f3eeda 440 }
9d84731e 441 }
442
443 }
444 }
445 std::cout << "Particles for thermal model fits:" << std::endl;
69f3eeda 446 arrOut->Print("short");
447 // arrOut->Print("");
9d84731e 448 std::cout << "" << std::endl;
449 // Write GSI input file
69f3eeda 450 TIter it(arrOutGSI);
9d84731e 451 AliParticleYield * part = 0;
452 ofstream fout(Form("gsi_System_%d_Energy_%0.0f_Centr_%s_BothCharges_%d", system, energy, centrality.Data(), separateCharges));
453 while ((part = (AliParticleYield*) it.Next())){
454 fout << part->GetYield() << " " << part->GetTotalError() << std::endl;
455 }
456 fout.close();
457 // Write thermus file
458 AliParticleYield::WriteThermusFile(arrOut, Form("thermus_System_%d_Energy_%0.0f_Centr_%s_BothCharges_%d", system, energy, centrality.Data(), separateCharges));
459}
ad431d97 460
461
b926af7d 462TH1F * GetHistoYields(TClonesArray * arr, Int_t system, Float_t energy, TString centrality, const char * histotitle,
463 Int_t icolor, Int_t imarker, Int_t errorsType, Float_t shift) {
ad431d97 464
b926af7d 465 TH1F * h = new TH1F(Form("hPart_%d_%0.0f_%s",system,energy,centrality.Data()), histotitle, npart, 1+shift, npart+1+shift);
ad431d97 466
ad431d97 467 for(Int_t ipart = 1; ipart <= npart; ipart++){
4910214c 468 std::cout << "----- Searching " << particleYields[ipart-1] << " -------" << std::endl;
469
b926af7d 470 AliParticleYield * part = AliParticleYield::FindParticle(arr, particleYields[ipart-1], system, energy, centrality,isSumYields[ipart-1]);
4910214c 471 if(part) {
472 std::cout << "found" << std::endl;
473 part->Print();
474 }
b926af7d 475 if(!part && isSumYields[ipart-1]) {
476 //Could not find the particle, but the sum was requested: build the sum!
477 part = AliParticleYield::FindParticle(arr, particleYields[ipart-1], system, energy, centrality, 0);
478 AliParticleYield * part2 = AliParticleYield::FindParticle(arr, -particleYields[ipart-1], system, energy, centrality, 0);
4910214c 479 if(part2 && part) {
480 std::cout << " Building sum" << std::endl;
481 part->Print();
482 part2->Print();
483 part = AliParticleYield::Add(part, part2);
484 }
485 else if(part) {
486 std::cout << "Scaling part" << std::endl;
487 part->Print();
488 part = new AliParticleYield(*part); // Always clone before scaling
489 part->Scale(2.); // If we only found a particle, we can scale it by a factor 2.
490 }
491 else part = 0;
b926af7d 492 }
493 if(!part){
494 std::cout << "Cannot find " << particleYields[ipart-1] << std::endl;
495 continue;
496 }
4910214c 497 if(isSumYields[ipart-1] && !part->IsTypeAverage()) {
498 std::cout << " scaling /2" << std::endl;
499 part = new AliParticleYield(*part); // Always clone before scaling
500 part->Scale(0.5); // take average
501 }
b926af7d 502 std::cout << " Plotting " << particleYields[ipart-1] << std::endl;
4910214c 503 part->Print();
b926af7d 504 // part->Scale(scale[ipart-1]);
ad431d97 505 h->SetBinContent(ipart, part->GetYield());
b926af7d 506 if(errorsType == kTotalError) {
507 h->SetBinError (ipart, part->GetTotalError(0/* 0 = no normalization error */));
508 } else if (errorsType == kSystError) {
509 h->SetBinError (ipart, part->GetSystError());
510 } else if (errorsType == kStatError) {
511 h->SetBinError (ipart, part->GetStatError());
512 }
513 h->GetXaxis()->SetBinLabel(ipart, part->GetLatexName());
ad431d97 514
515 }
b926af7d 516 h->SetMarkerStyle(imarker);
517 h->SetMarkerColor(icolor);
2f73d64b 518 h->SetLineColor(icolor);
4910214c 519 h->SetMarkerSize(1.4);
b926af7d 520
ad431d97 521 return h;
522}
046fa53e 523
524void LoadArrays() {
525 arrPbPb = AliParticleYield::ReadFromASCIIFile("PbPb_2760_Cascades.txt");
526 arrPbPb->AbsorbObjects( AliParticleYield::ReadFromASCIIFile("PbPb_2760_DeuHelium3.txt"));
527 arrPbPb->AbsorbObjects( AliParticleYield::ReadFromASCIIFile("PbPb_2760_Hypertriton.txt"));
528 arrPbPb->AbsorbObjects( AliParticleYield::ReadFromASCIIFile("PbPb_2760_Kstar892.txt"));
529 arrPbPb->AbsorbObjects( AliParticleYield::ReadFromASCIIFile("PbPb_2760_LambdaK0.txt"));
530 arrPbPb->AbsorbObjects( AliParticleYield::ReadFromASCIIFile("PbPb_2760_PiKaPr.txt"));
531 arrPbPb->AbsorbObjects( AliParticleYield::ReadFromASCIIFile("PbPb_2760_phi1020.txt"));
532 arrPbPb->AbsorbObjects( AliParticleYield::ReadFromASCIIFile("PbPb_2760_AveragedNumbers.txt"));
533
534 arrpp7 = AliParticleYield::ReadFromASCIIFile("pp_7000.txt");
535
536 arrpp276 = AliParticleYield::ReadFromASCIIFile("pp_2760.txt");
537 arrpp900 = AliParticleYield::ReadFromASCIIFile("pp_900.txt");
538
539 arrpPb = AliParticleYield::ReadFromASCIIFile("pPb_5020_MultiStrange.txt");
540 arrpPb->AbsorbObjects(AliParticleYield::ReadFromASCIIFile("pPb_5020_PiKaPrLamndaK0.txt"));
541 arrpPb->AbsorbObjects(AliParticleYield::ReadFromASCIIFile("pPb_5020_deuteron.txt"));
542 arrpPb->AbsorbObjects(AliParticleYield::ReadFromASCIIFile("pPb_5020_AveragedNumbers.txt"));
69f3eeda 543 arrpPb->AbsorbObjects(AliParticleYield::ReadFromASCIIFile("pPb_5020_phi.txt"));
544 arrpPb->AbsorbObjects(AliParticleYield::ReadFromASCIIFile("pPb_5020_Kstar.txt"));
046fa53e 545 arrThermus = AliParticleYield::ReadFromASCIIFile("PbPb_2760_Thermus_Boris_20140407.txt");
546
547 // RHIC data
548 arrSTARPbPb = AliParticleYield::ReadFromASCIIFile("PbPb_200_STAR-AntonQM12.txt");
549 arrPHENIXPbPb = AliParticleYield::ReadFromASCIIFile("PbPb_200_PHENIX-AntonQM12.txt");
550 arrBRAHMSPbPb = AliParticleYield::ReadFromASCIIFile("PbPb_200_BRAHMS-AntonQM12.txt");
551 arrSTARpp = AliParticleYield::ReadFromASCIIFile("pp_200_STAR.txt");
552 arrPHENIXpp = AliParticleYield::ReadFromASCIIFile("pp_200_PHENIX.txt");
553
554}
555
556void SetStyle(Bool_t graypalette) {
557 std::cout << "Setting style!" << std::endl;
558
559 gStyle->Reset("Plain");
560 gStyle->SetOptTitle(0);
561 gStyle->SetOptStat(0);
562 if(graypalette) gStyle->SetPalette(8,0);
563 else gStyle->SetPalette(1);
69f3eeda 564 gStyle->SetDrawBorder(0);
046fa53e 565 gStyle->SetCanvasColor(10);
566 gStyle->SetCanvasBorderMode(0);
69f3eeda 567 gStyle->SetPadBorderMode(0);
046fa53e 568 gStyle->SetFrameLineWidth(1);
569 gStyle->SetFrameFillColor(kWhite);
570 gStyle->SetPadColor(10);
571 gStyle->SetPadTickX(1);
572 gStyle->SetPadTickY(1);
573 gStyle->SetPadBottomMargin(0.15);
574 gStyle->SetPadLeftMargin(0.15);
575 gStyle->SetHistLineWidth(1);
576 gStyle->SetHistLineColor(kRed);
577 gStyle->SetFuncWidth(2);
578 gStyle->SetFuncColor(kGreen);
eee43f2a 579 gStyle->SetLineWidth(1);
046fa53e 580 gStyle->SetLabelSize(0.045,"yz");
581 gStyle->SetLabelSize(0.06,"x");
582 gStyle->SetLabelOffset(0.01,"y");
583 gStyle->SetLabelOffset(0.01,"x");
584 gStyle->SetLabelColor(kBlack,"xyz");
585 gStyle->SetTitleSize(0.05,"xyz");
586 gStyle->SetTitleOffset(1.25,"y");
587 gStyle->SetTitleOffset(1.2,"x");
588 gStyle->SetTitleFillColor(kWhite);
589 gStyle->SetTextSizePixels(26);
590 gStyle->SetTextFont(42);
eee43f2a 591 gStyle->SetTickLength(0.012,"X"); gStyle->SetTickLength(0.012,"Y");
046fa53e 592
593 gStyle->SetLegendBorderSize(0);
594 gStyle->SetLegendFillColor(kWhite);
595 // gStyle->SetFillColor(kWhite);
596 gStyle->SetLegendFont(42);
597
598 gStyle->SetErrorX(0);
599 gStyle->SetEndErrorSize(5);
600}
601
b926af7d 602TLegend * NewLegendQM(Double_t x1, Double_t y1, Double_t x2, Double_t y2, Bool_t isYield) {
046fa53e 603
b926af7d 604 const char * style = "p";
046fa53e 605 const char ** labels=0;
79936afc 606 // Bool_t beautify=kFALSE;
046fa53e 607 Bool_t useTitle=kTRUE;
046fa53e 608
15fccaa0 609 TLegend * l = new TLegend(x1, y1, x2, y2);
610 l->SetFillColor(kWhite);
611 l->SetTextFont(43);
612 l->SetTextSize(25);
613 l->SetBorderSize(1);
614 l->SetLineWidth(1);
b926af7d 615 l->SetMargin(0.1);
046fa53e 616 // const Int_t markers[] = {20,24,21,25,23,28,33,20,24,21,25,23,28,33};
617 // const Int_t colors[] = {1,2,3,4,6,7,8,9,10,11,1,2,3,4,6,7,8,9,10};
618
619 TList * list = gPad->GetListOfPrimitives();
620 TIterator * iter = list->MakeIterator();
621 TObject * obj = 0;
622 Int_t ilabel = -1;
623 while ((obj = (TObject*) iter->Next())){
624 if (obj->InheritsFrom("TH1") || obj->InheritsFrom("TGraph") || obj->InheritsFrom("TF1")) {
625 if( (TString(obj->GetName()) == "hframe" ) ) continue;
626 ilabel++;
627 if (labels != 0)
628 l->AddEntry(obj, labels[ilabel], style);
629 else{
630 if (useTitle) {
631 if(TString(obj->GetTitle()).Contains("NoLegend")) continue;
2f73d64b 632 TString title = obj->GetTitle();
633 if(title.Contains("p-Pb")) {
634 TObjArray * tokens = title.Tokenize(",");
635 if(tokens) {
636 TString system = ((TObjString*) tokens->At(0))->String().Strip(TString::kBoth, ' ');
637 TString centr = ((TObjString*) tokens->At(1))->String().Strip(TString::kBoth, ' ');
638 l->AddEntry(obj, system, style);
639 l->AddEntry(obj, centr, "0");
640 delete tokens;
641 }
642
643 } else {
644 l->AddEntry(obj, obj->GetTitle(), style);
645 }
646
046fa53e 647 }
648 else
649 l->AddEntry(obj, obj->GetName(), style);
650 }
79936afc 651 // if(beautify) {
046fa53e 652
79936afc 653 // if(!obj->InheritsFrom("TF1")){
654 // ((TH1*)obj)->SetLineColor(colors[ilabel]);
655 // ((TH1*)obj)->SetMarkerColor(colors[ilabel]);
656 // ((TH1*)obj)->SetMarkerStyle(markers[ilabel]);
657 // } else {
658 // ((TF1*)obj)->SetLineColor(colors[ilabel]);
659 // }
660
046fa53e 661 }
662 }
b926af7d 663 if(isYield) {
664 // Add some details on excluded stuff
eee43f2a 665 l->SetTextSize(22);
666 l->SetBorderSize(0);
667 l->GetEntry()->SetOption("0");
668 l->SetMargin(0.01);
b926af7d 669 }
670 l->Draw();
671 return l;
046fa53e 672}
673
674
15fccaa0 675
b926af7d 676void DrawRatio(TString what, Bool_t isYield, Double_t shift) {
046fa53e 677 // This is used to simplify the code above
678 // In order to draw syst error bars, we need to convert to graphs the syst errors histos
b926af7d 679 // if isYield == true plots yields rather than ratios
046fa53e 680
4c7bd9d1 681 // Sample colors
682 // const Int_t colors[] = {kBlack, kRed+1 , kBlue+1, kGreen+3, kMagenta+1, kOrange-1,kCyan+2,kYellow+2 , kWhite};
683
046fa53e 684 TClonesArray * array = 0;
685 Int_t system, color, marker;
b926af7d 686 Float_t energy = 0, shiftloc = shift;
046fa53e 687 TString centrality, label;
4c7bd9d1 688 // FIXME: move this in the different sections below
b926af7d 689 correlatedUnc = correlatedUncZero;
4c7bd9d1 690 std::cout << "Plotting " << what.Data() << std::endl;
691
046fa53e 692
693 if (what == "frame" ) {
4c7bd9d1 694 correlatedUnc = correlatedUncZero;
b926af7d 695 DrawFrame(isYield);
69f3eeda 696 // TH1 * h = GetHistoRatios(arrPbPb, AliParticleYield::kCSPbPb, 2760, "V0M0010", "NoLegend", kWhite);
697 // h->Draw();
698 // h->GetYaxis()->SetDecimals(1);
699 // h->GetYaxis()->SetNdivisions(505);
700 // h->GetXaxis()->CenterLabels(1);
701 // Int_t nratio = h->GetNbinsX();
b926af7d 702 if(!isYield) {
703 for(Int_t iratio = 0; iratio < nratio; iratio++){
704 Double_t exp = TMath::Floor(TMath::Log10(TMath::Abs(scale[iratio])));
705 Double_t man = scale[iratio] / TMath::Power(10, exp);
706 if(exp > 2) {
707 // TLatex * scaleLabel = new TLatex(iratio+1+0.2,maxy*1.01, Form("#times %0.0f 10^{%0.0f}", man, exp));
708 TLatex * scaleLabel = new TLatex(iratio+1+0.2,0.005, Form("#times %0.0f 10^{%0.0f}", man, exp));
709 scaleLabel->SetTextFont(43);
710 scaleLabel->SetTextSize(20);
711 scaleLabel->Draw();
712 } else {
713 Double_t shiftloc2 = scale[iratio] < 50 ? 0.3 : 0.2;
714 TLatex * scaleLabel = new TLatex(iratio+1+shiftloc2, 0.005, Form("#times %g", scale[iratio]));
715 scaleLabel->SetTextFont(43);
716 scaleLabel->SetTextSize(20);
717 scaleLabel->Draw();
718 }
719
720 }
046fa53e 721 }
eee43f2a 722 if(isYield) {
723 TPaveText *pt = new TPaveText(0.691767, 0.86069, 0.893574, 0.944865,"brNDC");
724 pt->SetBorderSize(0);
725 pt->SetFillColor(0);
726 pt->SetTextAlign(12);
727 pt->SetTextFont(43);
728 pt->SetTextSize(23);
729 pt->AddText("ALICE Preliminary");
730 pt->Draw();
731 }
732 else {
733 TPaveText *pt = new TPaveText( 0.176, 0.842881, 0.378514, 0.929595,"brNDC");
734 pt->SetBorderSize(0);
735 pt->SetFillColor(0);
736 pt->SetLineWidth(1);
737 pt->SetTextFont(43);
738 pt->SetTextSize(23);
739 pt->AddText("ALICE Preliminary");
740 pt->Draw();
741 }
046fa53e 742
4c7bd9d1 743 gPad->SetGridx();
046fa53e 744
745 }
746 else if (what == "PbPb_0010") {
747 array = arrPbPb;
748 system = 2; energy = 2760.; centrality = "V0M0010";
4910214c 749 label = "Pb-Pb #sqrt{s}_{NN} = 2.76 TeV, 0-10%";
2f73d64b 750 color = kRed+1;
046fa53e 751 marker = kFullCircle;
b926af7d 752 if(!shift) shiftloc = 0;
2f73d64b 753 correlatedUnc = correlatedUncLocalPbPbOnlyKStarPhi;
4c7bd9d1 754
046fa53e 755 }
756
757 else if (what == "PbPb_6080") {
79936afc 758 array = arrPbPb;
759 system = 2; energy = 2760.; centrality = "V0M6080";
4910214c 760 label = "Pb-Pb #sqrt{s}_{NN} = 2.76 TeV, 60-80%";
79936afc 761 color = kBlue+1;
762 marker = kFullCircle;
b926af7d 763 if(!shift) shiftloc = 0.0;
2f73d64b 764 correlatedUnc = correlatedUncLocalPbPbOnlyKStarPhi;
765 }
766 else if (what == "PbPb_2040") {
767 array = arrPbPb;
768 system = 2; energy = 2760.; centrality = "V0M2040";
769 label = "Pb-Pb #sqrt{s}_{NN} = 2.76 TeV, 20-40%";
770 color = kBlue+1;
771 marker = kFullCircle;
772 if(!shift) shiftloc = 0.0;
773 correlatedUnc = correlatedUncLocalPbPbOnlyKStarPhi;
046fa53e 774 }
775 else if (what == "pp7") {
79936afc 776 array = arrpp7;
777 system = 0; energy = 7000.; centrality = "";
778 label = "pp #sqrt{s} = 7 TeV";
4c7bd9d1 779 color = kMagenta+1;
79936afc 780 marker = kFullCircle;
b926af7d 781 if(!shift) shiftloc = 0.2;
fe520153 782 // correlatedUnc = correlatedUncLocalPP;
046fa53e 783 }
784 else if (what == "pp900") {
79936afc 785 array = arrpp900;
786 system = 0; energy = 900.; centrality = "";
787 label = "pp #sqrt{s} = 0.9 TeV";
4c7bd9d1 788 color = kCyan+2;
79936afc 789 marker = kFullCircle;
b926af7d 790 if(!shift) shiftloc = -0.2;
fe520153 791 // correlatedUnc = correlatedUncLocalPP;
046fa53e 792 }
793 else if (what == "pp276") {
79936afc 794 array = arrpp276;
795 system = 0; energy = 2760.; centrality = "";
796 label = "pp #sqrt{s} = 2.76 TeV";
4c7bd9d1 797 color = kYellow+2;
79936afc 798 marker = kFullCircle;
b926af7d 799 if(!shift) shiftloc = 0;
fe520153 800 // correlatedUnc = correlatedUncLocalPP;
046fa53e 801 }
802 else if (what == "pPb0005") {
79936afc 803 array = arrpPb;
4c7bd9d1 804 system = 1; energy = 5020.; centrality = "V0A0005";
2f73d64b 805 // label = "p-Pb #sqrt{s}_{NN} = 5.02 TeV, V0A 0-5%";
806 label = "p-Pb #sqrt{s}_{NN} = 5.02 TeV, V0A Multiplicity (Pb-Side) 0-5%";
79936afc 807 color = kBlack;
808 marker = kFullCircle;
b926af7d 809 if(!shift) shiftloc = -0.2;
fe520153 810 // correlatedUnc = correlatedUncLocalPP;
79936afc 811 }
2f73d64b 812 else if (what == "pPb2040") {
813 array = arrpPb;
814 system = 1; energy = 5020.; centrality = "V0A2040";
815 label = "p-Pb #sqrt{s}_{NN} = 5.02 TeV, V0A Multiplicity (Pb-Side) 20-40%";
816 color = kBlue+1;
817 marker = kFullCircle;
818 if(!shift) shiftloc = 0.;
819 // correlatedUnc = correlatedUncLocalPP;
820 }
821 else if (what == "pPb6080") {
822 array = arrpPb;
823 system = 1; energy = 5020.; centrality = "V0A6080";
824 label = "p-Pb #sqrt{s}_{NN} = 5.02 TeV, V0A Multiplicity (Pb-Side) 60-80%";
825 color = kGreen+3;
826 marker = kFullCircle;
827 if(!shift) shiftloc = 0.;
828 // correlatedUnc = correlatedUncLocalPP;
829 }
79936afc 830 else if (what == "PbPbSTAR") {
046fa53e 831 array = arrSTARPbPb;
4c7bd9d1 832 system = 2; energy = 200.; centrality = "0005";
4910214c 833 label = "STAR, Au-Au #sqrt{s}_{NN} = 0.2 TeV, 0-5%";
046fa53e 834 color = kBlack;
835 marker = kOpenStar;
b926af7d 836 if(!shift) shiftloc = +0.2;
fe520153 837 // correlatedUnc = correlatedUncZero;
79936afc 838 }
839 else if (what == "PbPbPHENIX") {
840 array = arrPHENIXPbPb;
4c7bd9d1 841 system = 2; energy = 200.; centrality = "0005";
4910214c 842 label = "PHENIX, Au-Au #sqrt{s}_{NN} = 0.2 TeV, 0-5%";
79936afc 843 color = kBlack;
844 marker = kOpenSquare;
b926af7d 845 if(!shift) shiftloc = -0.15;
fe520153 846 // correlatedUnc = correlatedUncZero;
79936afc 847 }
848 else if (what == "PbPbBRAHMS") {
849 array = arrBRAHMSPbPb;
4c7bd9d1 850 system = 2; energy = 200.; centrality = "0010";
4910214c 851 label = "BRAHMS, Au-Au #sqrt{s}_{NN} = 0.2 TeV, 0-10%";
79936afc 852 color = kBlack;
853 marker = kOpenCross;
b926af7d 854 if(!shift) shiftloc = -0.3;
fe520153 855 // correlatedUnc = correlatedUncZero;
4c7bd9d1 856 }
857 else if (what == "ppSTAR") {
858 array = arrSTARpp;
859 system = 0; energy = 200.; centrality = "";
4910214c 860 label = "STAR, pp #sqrt{s} = 0.2 TeV";
4c7bd9d1 861 color = kBlack;
862 marker = kOpenStar;
2f73d64b 863 if(!shift) shiftloc = -0.15;
fe520153 864 // correlatedUnc = correlatedUncZero;
4c7bd9d1 865 }
866 else if (what == "ppPHENIX") {
867 array = arrPHENIXpp;
868 system = 0; energy = 200.; centrality = "";
4910214c 869 label = "PHENIX, pp #sqrt{s} = 0.2 TeV";
4c7bd9d1 870 color = kBlack;
871 marker = kOpenSquare;
b926af7d 872 if(!shift) shiftloc = -0.2;
fe520153 873 // correlatedUnc = correlatedUncZero;
69f3eeda 874 }
875 // From here on, it's meta names, to draw several series of ratios
876 else if (what == "allpp"){
4910214c 877 DrawRatio("frame");
15fccaa0 878 DrawRatio("pp900");
4910214c 879 DrawRatio("pp276");
880 DrawRatio("pp7");
15fccaa0 881 array =0;
882 NewLegendQM(0.62249, 0.635734, 0.910643, 0.94673);
883
2f73d64b 884 SaveCanvas("Ratios_pponly");
15fccaa0 885 }
5ba62df2 886 else if (what == "allppWithRHIC"){
887 DrawRatio("frame");
888 DrawRatio("pp900", 0, 0.05);
889 DrawRatio("pp276", 0, 0.20);
890 DrawRatio("pp7", 0, 0.35);
891 DrawRatio("ppPHENIX", 0, -0.35);
892 DrawRatio("ppSTAR", 0, -0.15);
893
894 array =0;
895
896 NewLegendQM(0.588353, 0.636857, 0.910643, 0.948352);
897
898 SaveCanvas("Ratios_pponly_withRHIC");
899 }
15fccaa0 900 else if (what == "PbPbWithPP7TeV"){
901
902 DrawRatio("frame");
2f73d64b 903 DrawRatio("PbPb_0010", 0, 0.1);
904 DrawRatio("pp7", 0, -0.1);
15fccaa0 905 array =0;
2f73d64b 906
907
908 NewLegendQM(0.538153, 0.749397, 0.893574, 0.950362);
909 DrawExtrapolatedSymbolsAndLegendPbPb0010();
15fccaa0 910
2f73d64b 911 SaveCanvas("Ratios_withpp7tev");
15fccaa0 912
15fccaa0 913 } else if(what == "allsyst") {
914
915 DrawRatio("frame");
b926af7d 916 DrawRatio("pp7", 0, -0.2);
917 DrawRatio("pPb0005", 0, 0.00001);
4910214c 918 DrawRatio("PbPb_0010", 0, 0.2);
69f3eeda 919 array =0;
15fccaa0 920
15fccaa0 921
b926af7d 922 NewLegendQM(0.462851, 0.631722, 0.89257, 0.936697);
923 DrawExtrapolatedSymbolsAndLegendPbPb0010();
924 DrawExtrapolatedSymbolsAndLegendpPb0005();
2f73d64b 925
926 SaveCanvas("Ratios_allsystems");
15fccaa0 927
928 }else if(what =="PbPb6080andpPb0005") {
929 DrawRatio("frame");
2f73d64b 930 DrawRatio("pPb0005",0, -0.1);
931 DrawRatio("PbPb_6080", 0, 0.1);
15fccaa0 932 array=0;
933
4910214c 934 NewLegendQM( 0.46988, 0.730036, 0.910643, 0.948736);
b926af7d 935 DrawExtrapolatedSymbolsAndLegendpPb0005();
2f73d64b 936 SaveCanvas("Ratios_6080vspPb");
15fccaa0 937
938
939 }else if(what =="pp_vsRHIC") {
940 DrawRatio("frame");
941 DrawRatio("pp7");
2f73d64b 942 DrawRatio("ppPHENIX", 0, 0.00001);
15fccaa0 943 DrawRatio("ppSTAR");
15fccaa0 944 array=0;
945
b926af7d 946 NewLegendQM( 0.554217, 0.677869, 0.910643, 0.948736);
15fccaa0 947
2f73d64b 948 SaveCanvas("Ratios_vsRHIC_pp");
949
15fccaa0 950
951 } else if (what =="PbPb_vsRHIC") {
952 DrawRatio("frame");
953 DrawRatio("PbPb_0010");
954 DrawRatio("PbPbSTAR");
955 DrawRatio("PbPbPHENIX");
956 DrawRatio("PbPbBRAHMS");
957 array = 0;
958
959
b926af7d 960
961 NewLegendQM( 0.434739, 0.591593, 0.939759, 0.936697);
962 DrawExtrapolatedSymbolsAndLegendPbPb0010();
2f73d64b 963
964 SaveCanvas("Ratios_vsRHIC_PbPb");
15fccaa0 965 } else if (what == "aliceall") {
966 DrawRatio("frame");
967 DrawRatio("PbPb_0010");
968 DrawRatio("PbPb_6080");
969 DrawRatio("pPb0005");
970 DrawRatio("pp7");
971 DrawRatio("pp276");
972 DrawRatio("pp900");
046fa53e 973 }
974
b926af7d 975 // FROM HERE: IT's yields
976 else if( what == "fit_ReferenceFit_PbPb0010") {
977
978 DrawRatio("frame",1);
979 DrawRatio("PbPb_0010",1);
eee43f2a 980 particlesToExcludeFromChi2="[313]"; // do not consider K*
4910214c 981 legThermal->SetNColumns(4);
eee43f2a 982 AddLineToThermalLegend(legThermal, "Model,T (MeV), V (fm^{3}), #chi^{2}/NDF", "0");
2f73d64b 983 // FIXME: sistemare valori rapporti
984 shiftRatioDataModel=-0.2;
985 PlotThermusYields("lhc2760_final_0005_single_gc_output_gs1_wevc_nkst.txt" , kBlack , kSolid , "THERMUS 2.3 , 155 #pm 2 , 5924 #pm 543 , 23.6/9");
986 shiftRatioDataModel =0;
987 PlotGSIYields("data+therm_fit2_s2760_0-10qm14.dat" , kOrange-1 , kDashed , "GSI , 156 #pm 2 , 5330 #pm 505 , 17.4/9");
988 shiftRatioDataModel = 0.2;
989 PlotThermusYields("NEW_fit_gamma_q_s_fixed_PbPb_0010_with_nuclei.txt" , kBlue+1 , kDashDotted , "SHARE 3 , 156 #pm 3 , 4476 #pm 696 , 14.1/9");
b926af7d 990
eee43f2a 991
992 NewLegendQM(0.651606, 0.765993, 0.909639, 0.865951, 1);
4910214c 993 DrawMarkerKStarNoFit() ;
994 DrawExtrapolatedSymbolsYieldsPbPb0010();
b926af7d 995
2f73d64b 996 SaveCanvas("Fit_PbPb0010_Reference");
b926af7d 997
998 array =0;
8f3e9877 999 } else if( what == "fit_ReferenceFit_GSIONLY_PbPb0010") {
1000
1001 DrawRatio("frame",1);
1002 DrawRatio("PbPb_0010",1);
1003 particlesToExcludeFromChi2="[313]"; // do not consider K*
1004 legThermal->SetNColumns(4);
1005 AddLineToThermalLegend(legThermal, "Model,T (MeV), V (fm^{3}), #chi^{2}/NDF", "0");
1006 // FIXME: sistemare valori rapporti
1007 shiftRatioDataModel =0;
1008 PlotGSIYields("data+therm_fit2_s2760_0-10qm14.dat" , kBlack , kSolid , "GSI , 156 #pm 2 , 5330 #pm 505 , 17.4/9");
1009
1010 myPadHisto->cd();
1011 NewLegendQM(0.651606, 0.765993, 0.909639, 0.865951, 1);
1012 DrawMarkerKStarNoFit() ;
1013 DrawExtrapolatedSymbolsYieldsPbPb0010();
1014
1015 SaveCanvas("Fit_PbPb0010_Reference_GSI");
1016
1017 array =0;
1018 } else if( what == "fit_ReferenceFit_GSITHERMUS_PbPb0010") {
1019
1020 DrawRatio("frame",1);
1021 DrawRatio("PbPb_0010",1);
1022 particlesToExcludeFromChi2="[313]"; // do not consider K*
1023 legThermal->SetNColumns(4);
1024 AddLineToThermalLegend(legThermal, "Model,T (MeV), V (fm^{3}), #chi^{2}/NDF", "0");
1025 // FIXME: sistemare valori rapporti
1026 shiftRatioDataModel =-0.1;
1027 PlotThermusYields("lhc2760_final_0005_single_gc_output_gs1_wevc_nkst.txt" , kBlack , kSolid , "THERMUS 2.3 , 155 #pm 2 , 5924 #pm 543 , 23.6/9");
1028 shiftRatioDataModel =0.1;
1029 PlotGSIYields("data+therm_fit2_s2760_0-10qm14.dat" , kOrange-1 , kDashed , "GSI , 156 #pm 2 , 5330 #pm 505 , 17.4/9");
1030
1031
1032 myPadHisto->cd();
1033 NewLegendQM(0.651606, 0.765993, 0.909639, 0.865951, 1);
1034 DrawMarkerKStarNoFit() ;
1035 DrawExtrapolatedSymbolsYieldsPbPb0010();
1036
1037 SaveCanvas("Fit_PbPb0010_Reference_GSITHERMUS");
1038 // SaveCanvas("");
1039
1040 array =0;
1041 }else if( what == "fitSHARE_NoPionsNoProtons_PbPb0010") {
b926af7d 1042 array =0;
4910214c 1043
1044 DrawRatio("frame",1);
1045 DrawRatio("PbPb_0010",1);
1046 legThermal->SetNColumns(4);
1047 AddLineToThermalLegend(legThermal, "Model,T (MeV), V (fm^{3}), #chi^{2}/NDF", "0");
1048
2f73d64b 1049 particlesToExcludeFromChi2 = "[313]";
1050 shiftRatioDataModel=-0.2;
1051 PlotThermusYields("NEW_fit_gamma_q_s_fixed_PbPb_0010_with_nuclei.txt" ,kBlack , kSolid ,"SHARE 3 , 156 #pm 3 , 4476 #pm 696 , 14.1/9 ");
1052 particlesToExcludeFromChi2 = "[313][2212]";
1053 shiftRatioDataModel=0;
1054 PlotThermusYields("NEW_fit_gamma_q_s_fixed_PbPb_0010_with_nuclei_protons.txt" ,kOrange-1 , kDashed ,"SHARE 3 (no p) , 156 #pm 3 , 4520 #pm 623 , 8.2/8");
1055 particlesToExcludeFromChi2 = "[313][211]";
1056 shiftRatioDataModel=0.2;
1057 PlotThermusYields("NEW_fit_gamma_q_s_fixed_PbPb_0010_with_nuclei_pions.txt" ,kBlue+1 , kDashDotted ,"SHARE 3 (no #pi) , 157 #pm 3 , 4103 #pm 677 , 12.2/8");
4910214c 1058
2f73d64b 1059 myPadHisto->cd();
1060 NewLegendQM(0.651606, 0.765993, 0.909639, 0.865951, 1);
4910214c 1061
1062 // Add markers for additional particles not in fit
1063 DrawMarkerKStarNoFit() ;
1064 DrawExtrapolatedSymbolsYieldsPbPb0010();
b926af7d 1065
2f73d64b 1066 SaveCanvas("Fit_PbPb0010_SHARE_NoPiNoP");
4910214c 1067
1068
1069 } else if( what == "fitGSI_NoPionsNoProtons_PbPb0010") {
b926af7d 1070 array =0;
4910214c 1071
1072 DrawRatio("frame",1);
1073 DrawRatio("PbPb_0010",1);
1074 legThermal->SetNColumns(4);
1075 AddLineToThermalLegend(legThermal, "Model,T (MeV), V (fm^{3}), #chi^{2}/NDF", "0");
1076
2f73d64b 1077 particlesToExcludeFromChi2 = "[313]";
1078 shiftRatioDataModel=-0.2;
1079 PlotGSIYields("data+therm_fit2_s2760_0-10qm14.dat" , kBlack , kSolid , "GSI , 156 #pm 1.5 , 5330 #pm 505 , 18.1/9" );
1080 particlesToExcludeFromChi2 = "[313][2212]";
1081 shiftRatioDataModel=0.;
1082 PlotGSIYields("data+therm_fit2_s2760_0-10qm14NOp.dat" , kOrange-1 , kDashed , "GSI (no p) , 156 #pm 2 , 5590 #pm 330 , 7.7/8" );
1083 particlesToExcludeFromChi2 = "[313][211]";
1084 shiftRatioDataModel=0.2;
1085 PlotGSIYields("data+therm_fit2_s2760_0-10qm14NOpi.dat" , kBlue+1 , kDashDotted , "GSI (no #pi) , 157 #pm 2 , 4990 #pm 630 , 16.5/8");
4910214c 1086
2f73d64b 1087 myPadHisto->cd();
1088 NewLegendQM(0.651606, 0.765993, 0.909639, 0.865951, 1);
4910214c 1089
1090 // Add markers for additional particles not in fit
1091 DrawMarkerKStarNoFit() ;
1092 DrawExtrapolatedSymbolsYieldsPbPb0010();
2f73d64b 1093
1094 SaveCanvas("Fit_PbPb0010_GSI_NoPiNoP");
4910214c 1095
1096
1097 } else if (what == "fitShareWithWithoutNuclei_PbPb0010") {
1098 DrawRatio("frame",1);
1099 DrawRatio("PbPb_0010",1);
1100 legThermal->SetNColumns(4);
1101 AddLineToThermalLegend(legThermal, "Model,T (MeV), V (fm^{3}), #chi^{2}/NDF", "0");
1102
2f73d64b 1103 shiftRatioDataModel = -0.1;
1104 particlesToExcludeFromChi2="[313]"; // do not consider K* and nuclei
1105 PlotThermusYields("NEW_fit_gamma_q_s_fixed_PbPb_0010_with_nuclei.txt" , kBlack , kSolid , "SHARE 3 , 156 #pm 3 , 4476 #pm 696 , 14.1/9");
1106 shiftRatioDataModel = 0.1;
1107 particlesToExcludeFromChi2="[313][1000010020][1000020030][1010010030]"; // do not consider K* and nuclei
1108 PlotThermusYields("NEW_fit_gamma_q_fixed_PbPb_0010_without_nuclei.txt" , kCyan+2 , kDashDotted , "SHARE 3 (no nuclei), 156 #pm 4 , 4364 #pm 848 , 9.6/6");
4910214c 1109
2f73d64b 1110 NewLegendQM(0.651606, 0.765993, 0.909639, 0.865951, 1);
4910214c 1111
1112 DrawMarkerKStarNoFit() ;
1113 DrawExtrapolatedSymbolsYieldsPbPb0010();
2f73d64b 1114
1115 SaveCanvas("Fit_PbPb0010_SHARE_WithWoNuclei");
1116
b926af7d 1117
1118 } else if( what == "fitThermus_GammaSFree_PbPb0010") {
1119 std::cout << "MISSING DATA" << std::endl;
1120 array =0;
1121 } else if( what == "fitShare_GammaSGammaQFree_PbPb0010") {
1122 std::cout << "MISSING DATA" << std::endl;
1123 array =0;
1124 } else if( what == "fitShare_All_PbPb0010") {
1125 DrawRatio("frame",1);
1126 DrawRatio("PbPb_0010",1);
2f73d64b 1127
1128
1129 legThermal->SetNColumns(6);
1130 AddLineToThermalLegend(legThermal, "Model,T (MeV), V (fm^{3}), #gamma_{s}, #gamma_{q}, #chi^{2}/NDF", "0");
1131
1132 particlesToExcludeFromChi2="[313][1000010020][1000020030][1010010030]"; // do not consider K* and nuclei
1133 shiftRatioDataModel = -0.3;
1134 PlotThermusYields("NEW_fit_gamma_q_s_fixed_PbPb_0010_without_nuclei.txt" , kBlack , kSolid , "SHARE 3 , 156 #pm 4 , 4364 #pm 848 , 1 (fix) , 1 (fix) , 12.4/6");
1135 shiftRatioDataModel = -0.1;
1136 PlotThermusYields("NEW_fit_gamma_q_fixed_PbPb_0010_without_nuclei.txt" , kRed+1 , kDashDotted , "SHARE 3 , 155 #pm 3 , 4406 #pm 766 , 1.07 #pm 0.05 , 1 (fix) , 9.6/5");
1137 shiftRatioDataModel = 0.1;
1138 PlotThermusYields("NEW_fit_gamma_q_s_free_PbPb_0010_without_nuclei.txt" , kCyan+2 , kDashed , "SHARE 3 , 138 #pm 6 , 3064 #pm 1319 , 1.98 #pm 0.68 , 1.63 #pm 0.38 , 3.1/4");
1139 shiftRatioDataModel = 0.3;
1140 particlesToExcludeFromChi2="[313]"; // do not consider K*
1141 PlotThermusYields("NEW_fit_gamma_q_s_free_PbPb_0010_with_nuclei.txt" , kOrange - 1 , kDotted , "SHARE 3 (with nuclei) , 152 #pm 8 , 4445 #pm 743 , 1.16 #pm 0.20 , 1.06 #pm 0.12 , 9.0/7");
1142
1143
1144 NewLegendQM(0.651606, 0.765993, 0.909639, 0.865951, 1);
1145 DrawMarkerKStarNoFit() ;
1146 DrawMarkerNucleiNoFit();
1147 DrawExtrapolatedSymbolsYieldsPbPb0010();
1148
1149
1150 legThermal->SetX1(0.12249 );
1151 legThermal->SetY1(0.0454769);
1152 legThermal->SetX2(0.821285 );
1153 legThermal->SetY2(0.383481 );
1154
1155 SaveCanvas("Fit_PbPb0010_SHARE_All");
b926af7d 1156 array =0;
1157 } else if( what == "fitGSI_PbPb6080") {
4910214c 1158 DrawRatio("frame",1);
1159 DrawRatio("PbPb_6080",1);
1160 legThermal->SetNColumns(4);
1161 AddLineToThermalLegend(legThermal, "Model,T (MeV), V (fm^{3}), #chi^{2}/NDF", "0");
2f73d64b 1162 particlesToExcludeFromChi2 = "[313][1000020030][1010010030]";
1163 PlotGSIYields ("data+therm_fit2_s2760_60-80qm14.dat", kBlack, kSolid, "GSI, 157 #pm 2, 210 #pm 20, 8.2/7");
4910214c 1164
2f73d64b 1165 legThermal->SetX1(0.143574 );
1166 legThermal->SetY1(0.0731318 );
1167 legThermal->SetX2(0.659639 );
1168 legThermal->SetY2(0.245206 );
1169
1170 myPadHisto->cd();
1171 NewLegendQM(0.651606, 0.765993, 0.909639, 0.865951, 1);
4910214c 1172
2f73d64b 1173 DrawMarkerKStarNoFit(1) ;
1174
1175 // DrawExtrapolatedSymbolsYieldsPbPb0010(0.143574, 0.251352, 0.351406, 0.343535,0);
1176 //DrawExtrapolatedSymbolsYieldsPbPb0010();
1177
1178 SaveCanvas("Fit_PbPb6080_GSI");
1179 array =0;
1180 }else if( what == "fitGSI_PbPb2040") {
1181 DrawRatio("frame",1);
1182 DrawRatio("PbPb_2040",1);
1183 legThermal->SetNColumns(4);
1184 AddLineToThermalLegend(legThermal, "Model,T (MeV), V (fm^{3}), #chi^{2}/NDF", "0");
1185 particlesToExcludeFromChi2 = "[313][1000020030][1010010030]";
1186 PlotGSIYields ("data+therm_fit2_s2760_20-40qm14.dat", kBlack, kSolid, "GSI, 161 #pm 3, 1725 #pm 220, 24.6/7");
1187
1188 legThermal->SetX1(0.143574 );
1189 legThermal->SetY1(0.0731318 );
1190 legThermal->SetX2(0.659639 );
1191 legThermal->SetY2(0.245206 );
1192
1193 myPadHisto->cd();
1194 NewLegendQM(0.651606, 0.765993, 0.909639, 0.865951, 1);
1195
1196 DrawMarkerKStarNoFit(1) ;
1197 // DrawMarkerKStarNoFit() ;
1198
1199 // DrawExtrapolatedSymbolsYieldsPbPb0010(0.143574, 0.251352, 0.351406, 0.343535,0);
4910214c 1200 //DrawExtrapolatedSymbolsYieldsPbPb0010();
1201
2f73d64b 1202 SaveCanvas("Fit_PbPb2040_GSI");
b926af7d 1203 array =0;
1204 } else if( what == "fitGSI_pPb0005") {
1205 std::cout << "MISSING DATA" << std::endl;
1206 array =0;
1207 } else if( what == "fitGSI_pPb2040") {
1208 std::cout << "MISSING DATA" << std::endl;
1209 array =0;
1210 } else if( what == "fitThermus_GammaSFree_pPb0005") {
2f73d64b 1211
1212 DrawRatio("frame",1);
1213 DrawRatio("pPb0005",1, 0.00001);
1214 particlesToExcludeFromChi2="[313][1000020030][1010010030] "; // do not consider K*
1215 legThermal->SetNColumns(6);
1216 AddLineToThermalLegend(legThermal, "Model,T (MeV), #gamma_{s}, r_{C} (fm), r (fm), #chi^{2}/NDF", "0");
1217
1218 shiftRatioDataModel=-0.2;
1219 PlotThermusYields("lhc5020_final_0005_single_gc_output_gs1_wevc_nkst.txt" , kBlack , kSolid , "THERMUS 2.3 GC , 158 #pm 2 , 1 (fix) , N/A , 3.40 #pm 0.11 , 32.0/7");
1220 shiftRatioDataModel =0;
1221 PlotThermusYields("lhc5020_final_0005_single_gc_output_gsf_wevc_nkst.txt" , kOrange-1 , kDashed , "THERMUS 2.3 GC , 159 #pm 2 , 0.98 #pm 0.03 , N/A , 3.40 #pm 0.11 , 31.3/6");
1222 shiftRatioDataModel = 0.2;
1223 PlotThermusYields("lhc5020_final_0005_single_sc_output_gs1_woevc_nkst.txt", kCyan+2 , kDashDotted , "THERMUS 2.3 SC , 158 #pm 3 , 1 (fix) , 4.61 #pm 3.77 , 3.07 #pm 0.13 , 29.6/6");
1224 //lhc5020_final_0005_single_gc_output_gsf.txt
1225 //lhc5020_final_0005_single_sc_output_gs1.txt
1226
1227 NewLegendQM(0.650602, 0.694971, 0.909639, 0.865951, 1);
1228 DrawExtrapNotInFitpPb0005() ;
1229
1230 legThermal->SetX1(0.121486 );
1231 legThermal->SetY1(0.0741793);
1232 legThermal->SetX2(0.759036 );
1233 legThermal->SetY2(0.384575 );
1234
1235 SaveCanvas("Fit_pPb0005_THERMUS");
1236
1237 //DrawExtrapolatedSymbolsAndLegendPbPb0010();
1238
b926af7d 1239 array =0;
1240 } else if( what == "fitThermus_GammaSFree_pPb2040") {
2f73d64b 1241 DrawRatio("frame",1);
1242 DrawRatio("pPb2040",1, 0.00001);
1243 particlesToExcludeFromChi2="[313][1000020030][1010010030] "; // do not consider K*
1244 legThermal->SetNColumns(6);
1245 AddLineToThermalLegend(legThermal, "Model,T (MeV), #gamma_{s}, r_{C} (fm), r (fm), #chi^{2}/NDF", "0");
1246
1247 shiftRatioDataModel=-0.2;
1248 PlotThermusYields("lhc5020_final_2040_single_gc_output_gs1_wevc_nkst.txt" , kBlack , kSolid , "THERMUS 2.3 GC , 155 #pm 2 , 1 (fix) , N/A , 2.83 #pm 0.08 , 40.1/7");
1249 shiftRatioDataModel =0;
1250 PlotThermusYields("lhc5020_final_2040_single_gc_output_gsf_wevc_nkst.txt" , kOrange-1 , kDashed , "THERMUS 2.3 GC , 156 #pm 2 , 0.93 #pm 0.03 , N/A , 2.83 #pm 0.08 , 34.6/6");
1251 shiftRatioDataModel = 0.2;
1252 PlotThermusYields("lhc5020_final_2040_single_sc_output_gs1_woevc_nkst.txt", kCyan+2 , kDashDotted , "THERMUS 2.3 SC , 156 #pm 2 , 1 (fix) , 3.45 #pm 0.77 , 2.54 #pm 0.09 , 35.5/6");
1253 //lhc5020_final_0005_single_gc_output_gsf.txt
1254 //lhc5020_final_0005_single_sc_output_gs1.txt
1255
1256 NewLegendQM(0.650602, 0.694971, 0.909639, 0.865951, 1);
1257 DrawExtrapNotInFitpPb0005(0) ;
1258
1259 legThermal->SetX1(0.121486 );
1260 legThermal->SetY1(0.0741793);
1261 legThermal->SetX2(0.759036 );
1262 legThermal->SetY2(0.384575 );
1263
1264
1265
1266 SaveCanvas("Fit_pPb2040_THERMUS");
1267
b926af7d 1268 array =0;
1269 } else if( what == "fitThermus_GammaSFree_pPb6080") {
2f73d64b 1270
1271 DrawRatio("frame",1);
1272 DrawRatio("pPb6080",1, 0.00001);
1273 particlesToExcludeFromChi2="[313][1000020030][1010010030][1000010020]"; // do not consider K*
1274 legThermal->SetNColumns(6);
1275 AddLineToThermalLegend(legThermal, "Model,T (MeV), #gamma_{s}, r_{C} (fm), r (fm), #chi^{2}/NDF", "0");
1276
1277 shiftRatioDataModel=-0.2;
1278
1279 PlotThermusYields("lhc5020_final_6080_single_gc_output_gs1_wevc_nkst.txt" , kBlack , kSolid , "THERMUS 2.3 GC , 152 #pm 2 , 1 (fix) , N/A , 2.18 #pm 0.06 , 48.2/7");
1280 shiftRatioDataModel =0;
1281 PlotThermusYields("lhc5020_final_6080_single_gc_output_gsf_wevc_nkst.txt" , kOrange-1 , kDashed , "THERMUS 2.3 GC , 154 #pm 2 , 0.88 #pm 3 , N/A , 2.21 #pm 0.07 , 28.8/6");
1282 shiftRatioDataModel = 0.2;
1283 PlotThermusYields("lhc5020_final_6080_single_sc_output_gs1_woevc_nkst.txt", kCyan+2 , kDashDotted , "THERMUS 2.3 SC , 154 #pm 2 ,1 (fix) ,3.18 #pm 0.57 ,1.96 #pm 0.07 ,40.5/6");
1284
1285
1286 //lhc5020_final_0005_single_gc_output_gsf.txt
1287 //lhc5020_final_0005_single_sc_output_gs1.txt
1288
1289 NewLegendQM(0.650602, 0.694971, 0.909639, 0.865951, 1);
1290 DrawExtrapNotInFitpPb0005(0) ;
1291
1292 legThermal->SetX1(0.121486 );
1293 legThermal->SetY1(0.0268308 );
1294 legThermal->SetX2(0.758032 );
1295 legThermal->SetY2(0.337226 );
1296
1297
1298 SaveCanvas("Fit_pPb6080_THERMUS");
1299
1300
1301
b926af7d 1302 array =0;
1303 } else if( what == "fitThermus_RC_pPb0005") {
1304 std::cout << "MISSING DATA" << std::endl;
1305 array =0;
1306 } else if( what == "fitThermus_RC_pPb2040") {
1307 std::cout << "MISSING DATA" << std::endl;
1308 array =0;
1309 } else if( what == "fitThermus_RC_pPb6080") {
1310 std::cout << "MISSING DATA" << std::endl;
1311 array =0;
1312 } else if( what == "fitShare_pPb0005") {
8f3e9877 1313
1314
b926af7d 1315 DrawRatio("frame",1);
8f3e9877 1316 DrawRatio("pPb0005",1, 0.00001);
1317 legThermal->SetNColumns(6);
1318 AddLineToThermalLegend(legThermal, "Model,T (MeV), V (fm^{3}), #gamma_{s}, #gamma_{q}, #chi^{2}/NDF", "0");
1319 particlesToExcludeFromChi2="[313][1000010020][1000020030][1010010030]"; // do not consider K* and nuclei
1320
1321 shiftRatioDataModel = -0.2;
1322 PlotThermusYields("NEW_fit_gamma_q_s_fixed_pPb_0005_without_nuclei.txt" , kBlack , kSolid , "SHARE 3 , 158 #pm 3 , 121 #pm 18 , 1 (fix) , 1 (fix) , 24.3/6");
1323 shiftRatioDataModel = 0;
1324 PlotThermusYields("NEW_fit_gamma_q_fixed_pPb_0005_without_nuclei.txt" , kBlue+1 , kDashDotted , "SHARE 3 , 161 #pm 3 , 115 #pm 17 , 0.93 #pm 0.04 , 1 (fix) , 20.3/5");
1325 shiftRatioDataModel = 0.2;
1326 PlotThermusYields("NEW_fit_gamma_q_s_free_pPb_0005_without_nuclei.txt" , kOrange+2 , kDashDotted , "SHARE 3 , 144 #pm 1 , 81 #pm 25 , 1.599 #pm 0 , 1.71 #pm 0.06 , 11.4/4");
1327
1328 NewLegendQM(0.650602, 0.694971, 0.909639, 0.865951, 1);
1329 DrawExtrapNotInFitpPb0005() ;
1330 DrawMarkerNucleiNoFit();
1331 legThermal->SetX1(0.124498 );
1332 legThermal->SetY1(0.0715488);
1333 legThermal->SetX2(0.672691 );
1334 legThermal->SetY2(0.384575 );
1335
1336 SaveCanvas("Fit_pPb0005_SHARE");
1337
1338
b926af7d 1339 array =0;
1340
5ba62df2 1341 } else if( what == "fitShare_pPb0005_NoOmega") {
1342
1343
1344 DrawRatio("frame",1);
1345 DrawRatio("pPb0005",1, 0.00001);
1346 legThermal->SetNColumns(6);
1347 AddLineToThermalLegend(legThermal, "Model,T (MeV), V (fm^{3}), #gamma_{s}, #gamma_{q}, #chi^{2}/NDF", "0");
1348 particlesToExcludeFromChi2="[313][1000010020][1000020030][1010010030]"; // do not consider K* and nuclei
1349
1350 shiftRatioDataModel = -0.3;
1351 PlotThermusYields("NEW_fit_gamma_q_s_fixed_pPb_0005_without_nuclei.txt" , kBlack , kSolid , "SHARE 3 , 158 #pm 3 , 121 #pm 18 , 1 (fix) , 1 (fix), 24.3/6");
1352 shiftRatioDataModel = -.1;
1353 PlotThermusYields("NEW_fit_gamma_q_s_fixed_pPb_0005_without_nuclei_excluding_Omega.txt" , kBlue+1 , kDashDotted , "SHARE 3 (No #Omega) , 163 #pm 4 , 99 #pm 17 , 1 (fix) , 1 (fix), 15.6/5");
1354 shiftRatioDataModel = 0.1;
1355 PlotThermusYields("NEW_fit_gamma_q_fixed_pPb_0005_without_nuclei_excluding_Omega.txt" , kOrange+2 , kDashed , "SHARE 3 (No #Omega) , 163 #pm 3 , 100 #pm 16 , 0.96 #pm 0.04 , 1 (fix), 14.7/4");
1356 shiftRatioDataModel = 0.3;
1357 PlotThermusYields("NEW_fit_gamma_q_s_fixed_pPb_0005_with_nuclei_excluding_Omega.txt" , kCyan+2 , kDotted , "SHARE 3 (No #Omega + d) , 160 #pm 3 , 114 #pm 17 , 1 (fix) , 1 (fix), 18.6/6");
1358
1359 NewLegendQM(0.650602, 0.694971, 0.909639, 0.865951, 1);
1360 DrawExtrapNotInFitpPb0005() ;
1361 DrawMarkerNucleiNoFit();
1362 legThermal->SetX1(0.124498 );
1363 legThermal->SetY1(0.0715488);
1364 legThermal->SetX2(0.672691 );
1365 legThermal->SetY2(0.384575 );
1366
1367 SaveCanvas("Fit_pPb0005_SHARE_NoOmega");
1368
1369
1370 array =0;
1371
1372 } else if( what == "fitShare_pPb2040") {
b926af7d 1373 std::cout << "MISSING DATA" << std::endl;
1374 array =0;
1375 } else if( what == "fitShare_pPb6080") {
1376 std::cout << "MISSING DATA" << std::endl;
1377 array =0;
1378 } else if( what == "fitGSI_pp") {
2f73d64b 1379 DrawRatio("frame",1);
1380 DrawRatio("pp7",1, 0.00001);
1381 legThermal->SetNColumns(5);
1382 AddLineToThermalLegend(legThermal, "Model,T (MeV), V (fm^{3}), #gamma_{s}, #chi^{2}/NDF", "0");
1383 particlesToExcludeFromChi2 = "[313][1000020030][1010010030]";
1384 shiftRatioDataModel=-0.1;
1385 PlotGSIYields ("data+therm_fit2_s7000.dat", kBlack, kSolid, "GSI, 146 #pm 2, 25 #pm 2, 1 (fix), 78.2/7", 0);
1386 shiftRatioDataModel=0.1;
1387 PlotGSIYields ("data+therm_fit2_s7000gs.dat", kRed+1, kDashed, "GSI, 150 #pm 2, 23 #pm 2, 0.88, 45.6/7", 0);
1388 legThermal->SetX1(0.143574 );
1389 legThermal->SetY1(0.0731318);
1390 legThermal->SetX2(0.659639 );
1391 legThermal->SetY2(0.374263 );
1392 myPadHisto->cd();
1393
1394 NewLegendQM(0.725904, 0.761431, 0.874498, 0.841323, 1);
1395
1396 DrawMarkerKStarNoFit();
1397 TLegend * leg = new TLegend(0.143574, 0.380408, 0.351406, 0.4603, NULL, "brNDC");
1398 leg->SetBorderSize(0);
1399 leg->SetTextFont(43);
1400 leg->SetTextSize(14);
1401 leg->SetLineColor(1);
1402 leg->SetLineStyle(1);
1403 leg->SetLineWidth(1);
1404 leg->SetFillColor(0);
1405 leg->SetFillStyle(1001);
1406 TLegendEntry * entry=leg->AddEntry("TMarker","Not in fit","p");
1407 entry->SetLineColor(1);
1408 entry->SetLineStyle(1);
1409 entry->SetLineWidth(1);
1410 entry->SetMarkerColor(1);
1411 entry->SetMarkerStyle(markerNoFit);
1412 entry->SetMarkerSize(1.2);
1413 entry->SetTextFont(43);
1414 leg->Draw();
1415
1416 SaveCanvas("Fit_pp7000_GSI");
1417 //DrawExtrapNotInFitpPb0005(0);
b926af7d 1418 array =0;
5ba62df2 1419 }else if( what == "fitFlorence_pp") {
1420 DrawRatio("frame",1);
1421 DrawRatio("pp7",1, 0.00001);
1422 legThermal->SetNColumns(1);
1423 AddLineToThermalLegend(legThermal, "Model", "0");
1424 shiftRatioDataModel = -0.1;
1425 particlesToExcludeFromChi2 = "[1000020030][1010010030][1000010020]";
1426 PlotFlorenceYields("pp7000_Florence.txt", kBlack, kSolid, "Becattini et al.");
1427 shiftRatioDataModel = 0.1;
1428 particlesToExcludeFromChi2 = "[313][1000020030][1010010030]";
1429 PlotGSIYields ("data+therm_fit2_s7000gs.dat", kRed+1, kDashed, "GSI GC",0);
1430 myPadHisto->cd();
1431 NewLegendQM(0.725904, 0.761431, 0.874498, 0.841323, 1);
1432
1433 SaveCanvas("Fit_pp7000_Florence");
1434 //DrawExtrapNotInFitpPb0005(0);
1435 array =0;
b926af7d 1436 } else if( what == "fitGSI_fullCanonical") {
1437 std::cout << "MISSING DATA" << std::endl;
1438 array =0;
1439 } else if( what == "fitThermus_RC_pp") {
1440 std::cout << "MISSING DATA" << std::endl;
1441 array =0;
1442 } else if( what == "fitShare_pp") {
1443 std::cout << "MISSING DATA" << std::endl;
1444 array =0;
1445 }
1446
046fa53e 1447
15fccaa0 1448
046fa53e 1449 else {
1450 std::cout << "Unknown ratio " << what.Data() << std::endl;
15fccaa0 1451 return;
046fa53e 1452 }
1453
4c7bd9d1 1454 if(!correlatedUnc) {
1455 std::cout << "correlatedUnc not set!" << std::endl;
1456
1457 }
1458 std::cout << "CORR: " << correlatedUnc[1] << std::endl;
1459
046fa53e 1460 if(array) {
b926af7d 1461 if(isYield) {
4910214c 1462 TGraphErrors * hstat =
b926af7d 1463 AliPWGHistoTools::GetGraphFromHisto(GetHistoYields(array, system, energy, centrality, label, color, marker, kStatError, shiftloc)
4910214c 1464 ,0);
1465 hstat->SetMarkerSize(1.2);
1466 hstat->Draw("PZ");
b926af7d 1467 AliPWGHistoTools::GetGraphFromHisto(GetHistoYields(array, system, energy, centrality, label+"NoLegend", color, marker, kSystError, shiftloc)
1468 ,0)->Draw("[]");
1469 } else {
1470 AliPWGHistoTools::GetGraphFromHisto(GetHistoRatios(array, system, energy, centrality, label, color, marker, kStatError, shiftloc)
1471 ,0)->Draw("PZ");
1472 AliPWGHistoTools::GetGraphFromHisto(GetHistoRatios(array, system, energy, centrality, label+"NoLegend", color, marker, kSystError, shiftloc)
1473 ,0)->Draw("[]");
1474 }
046fa53e 1475 }
1476
1477
1478}
b926af7d 1479
1480void DrawFrame(Bool_t isYield) {
69f3eeda 1481
5ba62df2 1482 myCan = new TCanvas("myCan","ThermalFits",50,10,1000,isYield ? 950 : 650);
69f3eeda 1483 myCan->Draw();
1484 myCan->cd();
1485 // Set the Pads
eee43f2a 1486 Double_t boundaryLabels = isYield ? 0.88 : 0.85;//0.92;
b926af7d 1487 myPadLabel = new TPad("myPadLabel","myPadLabel",0.0, boundaryLabels,1.0,1.0);
69f3eeda 1488 myPadSetUp(myPadLabel);
1489 myPadLabel->Draw();
1490
eee43f2a 1491 myPadHisto = new TPad("myPadHisto","myPadHisto",0.0,isYield ? 0.4 : 0.05 ,1.0, boundaryLabels,0);
69f3eeda 1492 myPadSetUp(myPadHisto);
1493 myPadHisto->Draw();
1494
5ba62df2 1495 myPadStdDev = new TPad("myPadStdDev","myPadStdDev",0.0,0.035,1.0,0.215,0);
69f3eeda 1496 myPadSetUp(myPadStdDev);
b926af7d 1497 if(isYield) myPadStdDev->Draw();
1498 myPadStdDev->SetGridx();
1499 myPadStdDev->SetGridy();
eee43f2a 1500
1501 // This pad is for the ratios data/model
5ba62df2 1502 myPadRatio = new TPad("myPadRatio","myPadRatio",0.0,0.22,1.0,0.4,0);
eee43f2a 1503 myPadSetUp(myPadRatio);
1504 if(isYield) myPadRatio->Draw();
1505 myPadRatio->SetGridx();
1506 myPadRatio->SetGridy();
1507
1508
1509
69f3eeda 1510 myPadLabel->cd();
1511
15fccaa0 1512 double xLabelPosition[nratio] = {0.124498, 0.211847, 0.31, 0.38, 0.465, 0.575, 0.644, 0.72, 0.82, 0.905 };
b926af7d 1513 double xLabelYields[npart] = {0.115, 0.185, 0.270, 0.330, 0.422, 0.485, 0.570, 0.625, 0.695, 0.785, 0.835, 0.915};
15fccaa0 1514 double yLabelPosition = 0.40;
69f3eeda 1515
1516 // labels
b926af7d 1517 if(isYield) {
1518 for(Int_t ipart = 0; ipart < npart; ipart++){
1519 TLatex *myPart = new TLatex(xLabelYields[ipart],yLabelPosition, yieldsLabels[ipart]);
1520 myLatexDraw(myPart,20);
1521 }
1522 } else {
1523 for(Int_t iratio = 0; iratio < nratio; iratio++){
1524 TLatex *myRatio = new TLatex(xLabelPosition[iratio],yLabelPosition, ratiosLabels[iratio]);
1525 myLatexDraw(myRatio,20);
1526 }
69f3eeda 1527 }
15fccaa0 1528 // Xi's and Omega's bar (there was no way to convince root to draw it properly)
b926af7d 1529 if(isYield) {
eee43f2a 1530 TLine *line = new TLine(0.653,0.660,0.663,0.660);
b926af7d 1531 line->SetLineWidth(2);
1532 line->Draw();
eee43f2a 1533 line = new TLine(0.725,0.660,0.738,0.660);
1534 line->SetLineWidth(2);
1535 line->Draw();
1536 line = new TLine(0.362,0.660,0.374,0.660);
b926af7d 1537 line->SetLineWidth(2);
1538 line->Draw();
1539
b926af7d 1540
1541 }
1542 else {
2f73d64b 1543 TLine *line = new TLine(0.408,0.68,0.418,0.68);
b926af7d 1544 line->SetLineWidth(2);
1545 line->Draw();
2f73d64b 1546 line = new TLine(0.498,0.68,0.513,0.68);
b926af7d 1547 line->SetLineWidth(2);
1548 line->Draw();
1549
1550
1551 }
15fccaa0 1552
69f3eeda 1553
b926af7d 1554 if(isYield) {
69f3eeda 1555 myPadStdDev->cd();
1556 myPadStdDev->SetTopMargin(0.0);
1557 myPadStdDev->SetTicks(1,1);
1558
eee43f2a 1559 Float_t devMax = 4.7;
69f3eeda 1560
b926af7d 1561 TH2F *myBlankStdDev = new TH2F("myBlankStdDev","myBlankStdDev",npart,1,npart+1,10,-devMax,+devMax);
69f3eeda 1562 myBlankStdDev->GetXaxis()->SetLabelFont(43); // precision 3: size will be in pixels
1563 myBlankStdDev->GetYaxis()->SetLabelFont(43);
1564 myBlankStdDev->GetYaxis()->SetTitleFont(43);
1565 myBlankStdDev->SetLabelSize(23,"xy");
b926af7d 1566 myBlankStdDev->SetTitleSize(20,"y");
1567 myBlankStdDev->SetNdivisions(20,"x");
69f3eeda 1568 myBlankStdDev->SetNdivisions(605,"y");
1569 myBlankStdDev->SetLabelOffset(0.012,"xy");
2f73d64b 1570 myBlankStdDev->SetYTitle("(mod.-data)/#sigma_{data}");
eee43f2a 1571 myBlankStdDev->SetTitleOffset(1.65,"y");
69f3eeda 1572 myBlankStdDev->Draw();
eee43f2a 1573
2f73d64b 1574 TH2F *myBlankRatio = new TH2F("myBlankRatio","myBlankRatio",npart,1,npart+1,10,-0.8,0.8);
eee43f2a 1575 myBlankRatio->GetXaxis()->SetLabelFont(43); // precision 3: size will be in pixels
1576 myBlankRatio->GetYaxis()->SetLabelFont(43);
1577 myBlankRatio->GetYaxis()->SetTitleFont(43);
1578 myBlankRatio->SetLabelSize(23,"xy");
1579 myBlankRatio->SetTitleSize(20,"y");
1580 myBlankRatio->SetNdivisions(20,"x");
1581 myBlankRatio->SetNdivisions(605,"y");
2f73d64b 1582 myBlankRatio->SetLabelOffset(0.012,"xy");
1583 myBlankRatio->SetYTitle("(mod.-data)/mod.");
eee43f2a 1584 myBlankRatio->SetTitleOffset(1.65,"y");
1585 myPadRatio->cd();
2f73d64b 1586 // myPadRatio->SetLogy();
eee43f2a 1587 myBlankRatio->Draw();
1588
69f3eeda 1589 }
1590
1591 myPadHisto->cd();
2f73d64b 1592 myPadHisto->SetBottomMargin(isYield ? 0.01 : 0.04);
69f3eeda 1593 // myPadHisto->SetLogy();
1594 myPadHisto->SetTicks(1,1);
b926af7d 1595 TH1 *myBlankHisto = 0;
1596 if(isYield) {
1597 myBlankHisto = new TH2F("NoLegend","NoLegend",npart,1,npart+1,10, 0.00002, maxy+0.01 );
1598 gPad->SetLogy();
1599
1600 } else {
1601 myBlankHisto = new TH2F("NoLegend","NoLegend",nratio,1,nratio+1,10, 0, maxy+0.01 );
1602 }
69f3eeda 1603 myBlankHisto->GetXaxis()->SetLabelFont(43); // precision 3: size will be in pixels
1604 myBlankHisto->GetYaxis()->SetLabelFont(43);
1605 myBlankHisto->GetYaxis()->SetTitleFont(43);
2f73d64b 1606 myBlankHisto->SetLabelSize(23,"y");
1607 myBlankHisto->SetLabelSize(0,"x");
69f3eeda 1608 myBlankHisto->SetTitleSize(26,"y");
1609 myBlankHisto->SetMaximum(10);
1610 myBlankHisto->SetMinimum(0);
b926af7d 1611 myBlankHisto->SetNdivisions(isYield? 20 :10,"x");
69f3eeda 1612 myBlankHisto->SetNdivisions(505,"y");
b926af7d 1613 if(isYield) myBlankHisto->SetYTitle("d#it{N}/d#it{y}");
69f3eeda 1614 myBlankHisto->SetLabelOffset(0.012,"xy");
eee43f2a 1615 myBlankHisto->SetTitleOffset(isYield ? 1.3 : 1,"y");
69f3eeda 1616 myBlankHisto->Draw();
b926af7d 1617
1618 if(isYield) {
eee43f2a 1619 legThermal = new TLegend(0.144578, 0.0702247, 0.659639, 0.383226);
b926af7d 1620 legThermal->SetBorderSize(1);
1621 legThermal->SetTextFont(43);
2f73d64b 1622 legThermal->SetTextSize(18);
4910214c 1623 // legThermal->SetNColumns(6);
b926af7d 1624 legThermal->SetFillColor(0);
1625 legThermal->SetLineWidth(1);
1626 legThermal->Draw();
b926af7d 1627 }
69f3eeda 1628}
1629
1630void myLatexDraw(TLatex *currentLatex, Float_t currentSize, Int_t currentColor){
15fccaa0 1631 currentLatex->SetTextFont(43);
1632 currentLatex->SetTextSizePixels(Int_t(currentSize));
1633 // currentLatex->SetTextAngle(0);
69f3eeda 1634 currentLatex->Draw();
1635 return;
1636}
1637
1638void myPaveSetup(float rRatio, float rRange3, float rRange5,
1639 int rFillColor){
1640
1641 float cHiRange = 0, cLoRange = 0;
1642 if(rRange3<rRange5) {cHiRange=rRange5;cLoRange=rRange3;}
1643 else {cHiRange=rRange3;cLoRange=rRange5;}
1644
1645 TPave *cPave= new TPave(rRatio-0.25,cLoRange,rRatio+0.25,cHiRange,0,"br");
1646 cPave->SetFillColor(rFillColor);
1647 cPave->SetLineColor(1);
1648 cPave->Draw();
1649}
1650
1651void myPadSetUp(TPad *currentPad){
1652 currentPad->SetLeftMargin(0.10);
1653 currentPad->SetRightMargin(0.02);
1654 currentPad->SetTopMargin(0.02);
1655 currentPad->SetBottomMargin(0.02);
1656 return;
1657}
b926af7d 1658TGraphErrors* PlotThermusYields(const char * filename, Int_t color, Int_t lineStyle,
4910214c 1659 const char * tag) {
b926af7d 1660
1661 Int_t lw = lineStyle == kSolid ? 2 : 3; // Set line width
1662
1663
1664 std::map<Int_t,Double_t> mapYields;
1665 std::map<Int_t,Double_t> mapStdDev;
1666
1667 Int_t pdg;
1668 Double_t yield, stddev;
1669 ifstream thermusFile(filename);
1670 TString line;
eee43f2a 1671 Double_t chi2 = 0;
b926af7d 1672 std::cout << "---"<<tag<<"---" << std::endl;
1673
1674 // std::istream is(thermusFile);
1675 // Read the std dev and the ratio in 2 maps, then plot them in a graph.
1676 while(line.ReadLine(thermusFile, kTRUE)) {
1677 if(line.BeginsWith("#")) continue;
1678 TObjArray * tokens = line.Tokenize(" \t");
1679 if(tokens->GetEntries() != 3) continue;// not a line with data
1680 // thermusFile >> pdg >> yield >> stddev;
1681 pdg = ((TObjString*)tokens->At(0))->String().Atof();
1682 yield = ((TObjString*)tokens->At(1))->String().Atof();
1683 stddev = ((TObjString*)tokens->At(2))->String().Atof();
1684
1685 if( thermusFile.eof() ) break;
1686 std::cout << "PDG " << pdg << " " << yield << " " << stddev << std::endl;
eee43f2a 1687
b926af7d 1688 mapYields[TMath::Abs(pdg)] += yield;
1689 mapStdDev[TMath::Abs(pdg)] += stddev;
1690 if(pdg < 0) {
1691 // found the antiparticle: now compute the mean
1692 mapYields[TMath::Abs(pdg)] /=2;
1693 mapStdDev[TMath::Abs(pdg)] /=2;
1694 }
1695 delete tokens;
1696 }
eee43f2a 1697
b926af7d 1698 // Now plot
1699 TGraphErrors * gThermus = new TGraphErrors;
1700 TGraphErrors * gThermusStdDev = new TGraphErrors;
1701 for(Int_t ipart = 0; ipart < npart; ipart++){
1702 gThermus->SetPoint(ipart, ipart+1.5, mapYields[particleYields[ipart]]);
1703 gThermus->SetPointError(ipart, 0.3, 0);
1704 gThermusStdDev->SetPoint(ipart, ipart+1.5, mapStdDev[particleYields[ipart]]);
1705 gThermusStdDev->SetPointError(ipart, 0.3, 0);
1706 }
1707
1708 myPadHisto->cd();
1709 gThermus->Draw("PZ");
1710 gThermus->SetLineWidth(lw);
1711
1712 gThermus->SetLineColor(color);
1713 gThermus->SetLineStyle(lineStyle);
1714 gThermus->SetTitle("NoLegend");
1715
1716
2f73d64b 1717 // myPadStdDev->cd();
1718 // gThermusStdDev->Draw("PZ");
1719 // gThermusStdDev->SetLineWidth(lw);
1720 // gThermusStdDev->SetLineColor(color);
1721 // gThermusStdDev->SetLineStyle(lineStyle);
eee43f2a 1722 TGraphErrors* gStdDev2 = 0;
1723 TGraphErrors* gRatio = 0;
1724 std::cout << "CHI2: "
1725 << GetGraphRatioAndStdDev(gThermus, gRatio, gStdDev2)
1726 << std::endl;
1727 myPadRatio->cd();
1728 gRatio->Draw("PZ");
1729 myPadStdDev->cd();
1730 gStdDev2->Draw("PZ");
b926af7d 1731 myPadHisto->cd();
4910214c 1732 AddLineToThermalLegend(gThermus, tag, "l");
b926af7d 1733 return gThermus;
1734}
1735
5ba62df2 1736TGraphErrors* PlotFlorenceYields(const char * filename, Int_t color, Int_t lineStyle,
1737 const char * tag) {
1738
1739 Int_t lw = lineStyle == kSolid ? 2 : 3; // Set line width
1740
1741
1742 std::map<Int_t,Double_t> mapYields;
1743 std::map<Int_t,Double_t> mapStdDev;
1744
1745 Int_t pdg;
1746 Double_t yield, yieldbar, stddev;
1747 ifstream thermusFile(filename);
1748 TString line;
1749 std::cout << "---"<<tag<<"---" << std::endl;
1750
1751 // std::istream is(thermusFile);
1752 // Read the std dev and the ratio in 2 maps, then plot them in a graph.
1753 while(line.ReadLine(thermusFile, kTRUE)) {
1754 if(line.BeginsWith("#")) continue;
1755 TObjArray * tokens = line.Tokenize(" \t");
1756 if(tokens->GetEntries() != 4) continue;// not a line with data
1757 // thermusFile >> pdg >> yield >> stddev;
1758 pdg = ((TObjString*)tokens->At(0))->String().Atof();
1759 yield = ((TObjString*)tokens->At(2))->String().Atof();
1760 yieldbar = ((TObjString*)tokens->At(3))->String().Atof(); // Antiparticle yield
1761 if(pdg == 0) {
1762 // not a line with data
1763 delete tokens;
1764 continue ;
1765 }
1766 if( thermusFile.eof() ) break;
1767 std::cout << "PDG " << pdg << " " << yield << " " << yieldbar << std::endl;
1768
1769 mapYields[TMath::Abs(pdg)] += yieldbar ? ((yield+yieldbar)/2) : yield; // If the antiparticle exists, use the average
1770
1771 delete tokens;
1772 }
1773
1774 // Now plot
1775 TGraphErrors * gFlorence = new TGraphErrors;
1776 for(Int_t ipart = 0; ipart < npart; ipart++){
1777 gFlorence->SetPoint(ipart, ipart+1.5, mapYields[particleYields[ipart]]);
1778 gFlorence->SetPointError(ipart, 0.3, 0);
1779 }
1780
1781 myPadHisto->cd();
1782 gFlorence->Draw("PZ");
1783 gFlorence->SetLineWidth(lw);
1784
1785 gFlorence->SetLineColor(color);
1786 gFlorence->SetLineStyle(lineStyle);
1787 gFlorence->SetTitle("NoLegend");
1788
1789
1790
1791 TGraphErrors* gStdDev2 = 0;
1792 TGraphErrors* gRatio = 0;
1793 std::cout << "CHI2: "
1794 << GetGraphRatioAndStdDev(gFlorence, gRatio, gStdDev2)
1795 << std::endl;
1796 myPadRatio->cd();
1797 gRatio->Draw("PZ");
1798 myPadStdDev->cd();
1799 gStdDev2->Draw("PZ");
1800 myPadHisto->cd();
1801 AddLineToThermalLegend(gFlorence, tag, "l");
1802 return gFlorence;
1803}
b926af7d 1804
1805
1806TGraphErrors* PlotGSIYields(const char * filename, Int_t color, Int_t lineStyle,
2f73d64b 1807 const char * tag, Bool_t isPbPb) {
4910214c 1808
1809 // tag is a comma separated list of elements to be added to the legend as diferent columns
b926af7d 1810
1811 Int_t lw = lineStyle == kSolid ? 2 : 3; // Set line width
1812
2f73d64b 1813
b926af7d 1814 const Int_t pdgPbPb0010[] = {211, -211, 321, -321, 310, 313, 333, 2212, -2212, 3122, 3312, -3312, 3334, -3334, 1000010020, 1000020030, 1010010030, -1010010030};
2f73d64b 1815 const Int_t pdgpp[] = {211, 321, 310, 313, 333, 2212, 3122, 3312, 3334, 1000010020};
1816
1817 const Int_t *pdgOrder = isPbPb ? pdgPbPb0010 : pdgpp;
b926af7d 1818
1819 std::map<Int_t,Double_t> mapYields;
1820 std::map<Int_t,Double_t> mapStdDev;
1821 std::map<Int_t,Double_t> mapUncert;
1822 std::map<Int_t,Double_t> mapData;
1823
1824 Double_t data, uncert, model;
1825
1826 ifstream gsiFile(filename);
1827 // std::istream is(thermusFile);
2f73d64b 1828 std::cout << "----- "<<tag<<" -----" << std::endl;
b926af7d 1829
1830 // Read the std dev and the ratio in 2 maps, then plot them in a graph.
1831 Int_t ipart = 0;
1832 while(gsiFile) {
1833 gsiFile >> data >> uncert >> model;
1834 if( gsiFile.eof() ) break;
2f73d64b 1835 Int_t pdg = pdgOrder[ipart];
eee43f2a 1836 std::cout << "PDG " << pdg << " " << data << std::endl;
b926af7d 1837 mapYields[TMath::Abs(pdg)] += model;
1838 mapUncert[TMath::Abs(pdg)] += uncert;
1839 mapData[TMath::Abs(pdg)] += data;
1840 if(pdg < 0) {
1841 // found the antiparticle: now compute the mean
1842 mapYields[TMath::Abs(pdg)] /=2;
1843 mapData[TMath::Abs(pdg)] /=2;
1844 mapUncert[TMath::Abs(pdg)] /=2;
1845
1846 }
1847 ipart++;
1848 }
1849
1850 // Now plot
1851 TGraphErrors * gGsi = new TGraphErrors;
1852 TGraphErrors * gGsiStdDev = new TGraphErrors;
1853 std::cout << "PDG \tmodel\tdata\tuncert\tstddev" << std::endl; // header
1854 for(Int_t ipart = 0; ipart < npart; ipart++){ // Here we use npart, as this is what we wnat to plot!
1855 Int_t pdg = particleYields[ipart];
1856 mapStdDev[TMath::Abs(pdg)] = ( mapYields[TMath::Abs(pdg)] - mapData[TMath::Abs(pdg)]) / mapUncert[TMath::Abs(pdg)] ;
b926af7d 1857 std::cout << "PDG " << pdg <<"\t"
1858 << mapYields[TMath::Abs(pdg)] << "\t" << mapData[TMath::Abs(pdg)] <<"\t"
1859 << mapUncert[TMath::Abs(pdg)] << "\t" << mapStdDev[TMath::Abs(pdg)]
1860 << std::endl;
1861
2f73d64b 1862 if(!mapYields[particleYields[ipart]]) mapYields[particleYields[ipart]] = -10;
1863 if(!mapStdDev[particleYields[ipart]]) mapStdDev[particleYields[ipart]] = -10;
1864
1865
b926af7d 1866 gGsi->SetPoint(ipart, ipart+1.5, mapYields[particleYields[ipart]]);
1867 gGsi->SetPointError(ipart, 0.3, 0);
1868 gGsiStdDev->SetPoint(ipart, ipart+1.5, mapStdDev[particleYields[ipart]]);
1869 gGsiStdDev->SetPointError(ipart, 0.3, 0);
1870 }
1871
eee43f2a 1872
1873
b926af7d 1874 myPadHisto->cd();
1875 gGsi->Draw("PZ");
1876 gGsi->SetLineWidth(lw);
15fccaa0 1877
b926af7d 1878 gGsi->SetLineColor(color);
1879 gGsi->SetLineStyle(lineStyle);
1880 gGsi->SetTitle("NoLegend");
1881
1882
2f73d64b 1883 // myPadStdDev->cd();
1884 // gGsiStdDev->Draw("PZ");
1885 // gGsiStdDev->SetLineWidth(lw);
1886 // gGsiStdDev->SetLineColor(color);
1887 // gGsiStdDev->SetLineStyle(lineStyle);
1888 // myPadHisto->cd();
b926af7d 1889
eee43f2a 1890 TGraphErrors* gStdDev2 = 0;
1891 TGraphErrors* gRatio = 0;
1892 std::cout << "CHI2: "
1893 << GetGraphRatioAndStdDev(gGsi, gRatio, gStdDev2)
1894 << std::endl;
1895 myPadRatio->cd();
1896 gRatio->Draw("PZ");
1897 myPadStdDev->cd();
1898 gStdDev2->Draw("PZ");
1899
4910214c 1900 AddLineToThermalLegend(gGsi, tag, "L");
b926af7d 1901
1902 return gGsi;
1903}
1904
1905
1906void DrawExtrapolatedSymbolsAndLegendPbPb0010() {
1907 myPadHisto->cd();
1908
1909
1910 TLegend *leg = new TLegend( 0.149598, 0.782203, 0.415663, 0.858447,NULL,"pNDC");
1911 leg->SetBorderSize(0);
1912 leg->SetTextFont(43);
1913 leg->SetTextSize(18);
1914 leg->SetLineColor(1);
1915 leg->SetLineStyle(1);
1916 leg->SetLineWidth(2);
1917 leg->SetFillColor(0);
1918 leg->SetFillStyle(1001);
1919 leg->SetMargin(0.1);
1920
1921 TLegendEntry *entry=leg->AddEntry("TMarker","Extrapolated (Pb-Pb 0-10%)","p");
1922 entry->SetLineColor(1);
1923 entry->SetLineStyle(1);
1924 entry->SetLineWidth(1);
1925 entry->SetMarkerColor(1);
1926 entry->SetMarkerStyle(27);
1927 entry->SetMarkerSize(1.2);
1928 entry->SetTextFont(43);
1929 leg->Draw();
1930 myPadLabel->cd();
1931 // Markers for extrapolated points
2f73d64b 1932
b926af7d 1933 TMarker *marker = new TMarker(0.666667,0.111825,27);
1934 marker->SetMarkerStyle(27);
1935 marker->SetMarkerSize(1.2);
1936 marker->Draw();
2f73d64b 1937 marker = new TMarker(0.920683,0.111825,27);
b926af7d 1938 marker->SetMarkerStyle(27);
1939 marker->SetMarkerSize(1.2);
1940 marker->Draw();
1941
2f73d64b 1942 // BR for 3He
1943
1944 myPadLabel->cd();
1945 TLatex * tex = new TLatex(0.73,0.05,"BR = 25%");
1946 tex->SetNDC();
1947 tex->SetTextFont(43);
1948 tex->SetTextSize(14);
1949 tex->SetLineWidth(2);
1950 tex->Draw();
1951
1952
1953
b926af7d 1954}
1955void DrawExtrapolatedSymbolsAndLegendpPb0005() {
1956 myPadHisto->cd();
1957
1958
1959 TLegend *leg = new TLegend( 0.149598, 0.709972, 0.415663, 0.786216,NULL,"pNDC");
1960 leg->SetBorderSize(0);
1961 leg->SetTextFont(43);
1962 leg->SetTextSize(18);
1963 leg->SetLineColor(1);
1964 leg->SetLineStyle(1);
1965 leg->SetLineWidth(2);
1966 leg->SetFillColor(0);
1967 leg->SetFillStyle(1001);
1968 leg->SetMargin(0.1);
1969 TLegendEntry *entry=leg->AddEntry("TMarker","Extrapolated (p-Pb 0-5%)","p");
1970 entry->SetLineColor(1);
1971 entry->SetLineStyle(1);
1972 entry->SetLineWidth(1);
1973 entry->SetMarkerColor(1);
1974 entry->SetMarkerStyle(28);
1975 entry->SetMarkerSize(1.2);
1976 entry->SetTextFont(43);
1977 leg->Draw();
1978 myPadLabel->cd();
1979
2f73d64b 1980 TMarker *marker = new TMarker(0.590361,0.111825,28);
b926af7d 1981 marker->SetMarkerStyle(28);
1982 marker->SetMarkerSize(1.2);
1983 marker->Draw();
2f73d64b 1984 marker = new TMarker(0.938755,0.111825,28);
b926af7d 1985 marker->SetMarkerStyle(28);
1986 marker->SetMarkerSize(1.2);
1987 marker->Draw();
1988}
4910214c 1989
1990void AddLineToThermalLegend(TObject * obj, TString line, const char * optFirst) {
1991
1992 // This should be a comma-separated list of text to be added to the
1993 // columns. If the number of entries does not match the numer of
1994 // columns, it returns an error
1995 TObjArray * tokens = line.Tokenize(",");
1996 if(tokens->GetEntries() != legThermal->GetNColumns()) {
1997 std::cout << "Wrong number of columns (" << tokens->GetEntries() << ","<<legThermal->GetNColumns()<<") for the thermal legend, not adding " << line.Data() << std::endl;
1998 return;
1999 }
2000
2001 TIter iter(tokens);
2002 TObjString * col;
2003 Bool_t first = 1;
2004 while((col = (TObjString*)iter.Next())) {
2005 // Add entry removing whitespaces
2006 legThermal->AddEntry(obj, col->String().Strip(TString::kBoth, ' ').Data(), first ? optFirst : "0");
2007 if (first) first = 0;
2008 }
2009}
2010
2f73d64b 2011void DrawExtrapolatedSymbolsYieldsPbPb0010(Double_t x1, Double_t y1, Double_t x2, Double_t y2, Bool_t plotExtraploatedLegend){
4910214c 2012 // Markers for extrapolated points
2013 myPadLabel->cd();
2f73d64b 2014 TMarker * marker = new TMarker(0.36245,0.111825,markerExtrap);
4910214c 2015 marker->SetMarkerStyle(markerExtrap);
2016 marker->SetMarkerSize(1.2);
2017 marker->Draw();
2f73d64b 2018 marker = new TMarker(0.945783,0.111825,markerExtrap);
4910214c 2019 marker->SetMarkerStyle(markerExtrap);
2020 marker->SetMarkerSize(1.2);
2021 marker->Draw();
2f73d64b 2022 // BR for 3He
2023 myPadHisto->cd();
2024 TLatex * tex = new TLatex(11.15, 2.5e-5,"BR = 25%");
2025 tex->SetTextFont(43);
2026 tex->SetTextSize(14);
2027 tex->SetLineWidth(2);
2028 tex->Draw();
2029
4910214c 2030 myPadHisto->cd();
eee43f2a 2031
2032
2f73d64b 2033
2034 TLegend * leg = new TLegend(x1,y1,x2,y2,NULL,"brNDC");
eee43f2a 2035 leg->SetBorderSize(0);
2036 leg->SetTextFont(43);
2037 leg->SetTextSize(14);
2038 leg->SetLineColor(1);
2039 leg->SetLineStyle(1);
2040 leg->SetLineWidth(1);
2041 leg->SetFillColor(0);
2042 leg->SetFillStyle(1001);
2043 TLegendEntry * entry=leg->AddEntry("TMarker","Not in fit","p");
2044 entry->SetLineColor(1);
2045 entry->SetLineStyle(1);
2046 entry->SetLineWidth(1);
2047 entry->SetMarkerColor(1);
2048 entry->SetMarkerStyle(markerNoFit);
2049 entry->SetMarkerSize(1.2);
2050 entry->SetTextFont(43);
2f73d64b 2051 if(plotExtraploatedLegend) {
2052 entry=leg->AddEntry("TMarker","Extrapolated","p");
2053 entry->SetLineColor(1);
2054 entry->SetLineStyle(1);
2055 entry->SetLineWidth(1);
2056 entry->SetMarkerColor(1);
2057 entry->SetMarkerStyle(markerExtrap);
2058 entry->SetMarkerSize(1.2);
2059 entry->SetTextFont(43);
2060 }
eee43f2a 2061 leg->Draw();
2062 leg->SetMargin(0.1);
2063
2064
4910214c 2065}
2066
2f73d64b 2067void DrawMarkerKStarNoFit(Bool_t plotLegend) {
4910214c 2068 myPadLabel->cd();
2f73d64b 2069 TMarker *marker = new TMarker(0.344378,0.111825,markerNoFit);
2070 marker->SetMarkerStyle(markerNoFit);
2071 marker->SetMarkerSize(1.2);
2072 marker->Draw();
2073 myPadHisto->cd();
2074
2075 if(plotLegend) {
2076
2077
2078 TLegend * leg = new TLegend(0.126506, 0.253051, 0.335341, 0.345118,NULL,"brNDC");
2079 leg->SetBorderSize(0);
2080 leg->SetTextFont(43);
2081 leg->SetTextSize(14);
2082 leg->SetLineColor(1);
2083 leg->SetLineStyle(1);
2084 leg->SetLineWidth(1);
2085 leg->SetFillColor(0);
2086 leg->SetFillStyle(1001);
2087 TLegendEntry * entry=leg->AddEntry("TMarker","Not in fit","p");
2088 entry->SetLineColor(1);
2089 entry->SetLineStyle(1);
2090 entry->SetLineWidth(1);
2091 entry->SetMarkerColor(1);
2092 entry->SetMarkerStyle(markerNoFit);
2093 entry->SetMarkerSize(1.2);
2094 entry->SetTextFont(43);
2095 leg->Draw();
2096 }
2097}
2098
2099void DrawMarkerNucleiNoFit() {
2100 myPadLabel->cd();
2101 TMarker *marker = new TMarker(0.928715,0.111825,markerNoFit);
2102 marker->SetMarkerStyle(markerNoFit);
2103 marker->SetMarkerSize(1.2);
2104 marker->Draw();
2105
2106 marker = new TMarker(0.791751,0.111825,markerNoFit);
2107 marker->SetMarkerStyle(markerNoFit);
2108 marker->SetMarkerSize(1.2);
2109 marker->Draw();
2110
2111 marker = new TMarker(0.866466,0.111825,markerNoFit);
2112 marker->SetMarkerStyle(markerNoFit);
4910214c 2113 marker->SetMarkerSize(1.2);
2114 marker->Draw();
2115 myPadHisto->cd();
2f73d64b 2116
4910214c 2117}
eee43f2a 2118
2119Double_t GetGraphRatioAndStdDev(TGraphErrors * gModel, TGraphErrors * &gRatio, TGraphErrors *&gStdDev) {
2120
2121 // Plots the ratios data/model
2122
2123 // I recomputed the stddev and the chi2 here, because some values
2124 // changed slightly due to typos in some of the input files with
2125 // respect to the ones used for those fits.
2126
eee43f2a 2127
2128 if(!gModel) {
2129 std::cout << "EMPTY MODEL" << std::endl;
2130 return 0;
2131 }
2132 // 0. set up the graphs which we will need
2133 TGraphErrors * gTemp = 0;
2134 TGraphErrors * gStat = 0;
2135 TGraphErrors * gSyst = 0;
2136 // the cloning below takes care of the style
2137 gRatio = (TGraphErrors*)gModel->Clone();
2138 gRatio->Clear();
2139 gStdDev = (TGraphErrors*)gModel->Clone();
2140 gStdDev->Clear();
2141
2142 // 1. Find the data graphs. We need both the stat and the syst, since we want the ratio with the total uncertainty (or not?)
2143 TVirtualPad * currentPad = gPad;
2144 myPadHisto->cd();
2145 TList * list = myPadHisto->GetListOfPrimitives();
2146 TIterator * iter = list->MakeIterator();
2147
2148 while ((gTemp = (TGraphErrors*) iter->Next())){
2149 if(gTemp->InheritsFrom("TGraphErrors")) {
2150 // Found a graph, it is the correct one?
2151 TString name = gTemp->GetName();
2152 TString title = gTemp->GetTitle();
2153 std::cout << "name " << name.Data() << std::endl;
2154
2155 if (name.Contains("hPart")) {
2156 // ok, it's the data
2157 if (title.Contains("NoLegend")) gSyst = gTemp; // it's the syst error FIXME: it would be better to add the error type to the data
2158 else gStat = gTemp;
2159 }
2160 if(gStat && gSyst) {
2161 std::cout << "Cannot find gStat or gSyst (" << gStat << "," << gSyst <<")" << std::endl;
2162 break; // found both stat and syst uncertainties
2163 }
2164 }
2165 }
2166
2167
2168
2169 // Compute the ratio, the stddev and the chi2
2170 Int_t npoint = gModel->GetN(); // We are sure that data and model have the same number of points in the same order, because they were created using the particleYields array (FIXME: is this also true for GSI?)
2171 Double_t chi2 = 0;
2172 for(Int_t ipoint = 0; ipoint < npoint; ipoint++){
2173
2174 Double_t yield = gStat->GetY()[ipoint];
2175 Double_t stat = gStat->GetEY()[ipoint];
2176 Double_t syst = gSyst->GetEY()[ipoint];
2177 Double_t error = TMath::Sqrt(stat*stat+syst*syst);
2178 Double_t model = gModel->GetY()[ipoint];
2179 Double_t width = gModel->GetEX()[ipoint];
2f73d64b 2180 Double_t ratio = (model-yield)/model;
eee43f2a 2181 Double_t stddev = (model-yield)/error;
2182 if(!particlesToExcludeFromChi2.Contains(Form("[%d]", particleYields[ipoint]))) chi2 += TMath::Power(stddev,2);
2183 else std::cout << "Ecluding PDG "<< particleYields[ipoint] <<" from chi2 calculation" << std::endl;
2f73d64b 2184 Double_t errorRatio = error/model;
2185 gRatio->SetPoint(ipoint, gModel->GetX()[ipoint]+shiftRatioDataModel, ratio);
eee43f2a 2186 gRatio->SetPointError(ipoint, 0, errorRatio);
2f73d64b 2187 if(model) {
2188 gStdDev->SetPoint(ipoint, gModel->GetX()[ipoint], stddev);
2189 gStdDev->SetPointError(ipoint, width, 0);
2190 }
2191 // The commented block down here was used to compare to the estimate provided directly by Boris and Benjamin
2192 // gStdDev->SetPoint(ipoint, gModel->GetX()[ipoint]-0.4, stddev);
2193 // gStdDev->SetPointError(ipoint, 0.2, 0);
eee43f2a 2194 // std::cout << "PDG " << particleYields[ipoint] <<"\t"
2195 // << model << "\t" << yield <<"\t"
2196 // << error << "\t" << stddev
2197 // << std::endl;
2198 }
2199
2200 gRatio->SetLineStyle(kSolid);
2f73d64b 2201 gRatio->SetMarkerStyle(kOpenSquare);
eee43f2a 2202 gRatio->SetMarkerColor(gRatio->GetLineColor());
2203 gRatio->SetLineWidth(2);
2204
eee43f2a 2205 currentPad->cd();
2206 return chi2;
2207}
2f73d64b 2208
2209
2210void SaveCanvas(const char * name) {
2211
2212 if(!saveCanvas) return;
2213
2214 std::cout << "Saving " << name << ".{eps,pdf,root,C}" << std::endl;
2215
2216 myCan->Update();
2217 gSystem->ProcessEvents();
2218 myCan->Print(Form("%s.eps",name));
2219 myCan->Print(Form("%s.root",name));
2220 myCan->Print(Form("%s.C",name));
2221 gSystem->Exec(Form("epstopdf %s.eps", name));
2222 gSystem->Exec(Form("if [ \"$USER\" = \"mfloris\" ]; then cp %s.{eps,pdf,root,C} /Users/mfloris/Documents/PapersNTalks/ALICE/ThermalFits/img/; fi ",name));
2223
2224}
2225
2226void DrawExtrapNotInFitpPb0005(Bool_t drawExtrap) {
2227 myPadLabel->cd();
2228 TMarker *marker = new TMarker(0.344378,0.111825,28);
2229 marker->SetMarkerStyle(28);
2230 marker->SetMarkerSize(1.2);
2231 marker->Draw();
2232 // marker = new TMarker(0.7851406,0.111825,28);
2233 // marker->SetMarkerStyle(28);
2234 // marker->SetMarkerSize(1.2);
2235 // marker->Draw();
2236 if(drawExtrap) {
2237 marker = new TMarker(0.364458,0.111825,27);
2238 marker->SetMarkerStyle(27);
2239 marker->SetMarkerSize(1.2);
2240 marker->Draw();
2241 marker = new TMarker(0.792,0.111825,27);
2242 marker->SetMarkerStyle(27);
2243 marker->SetMarkerSize(1.2);
2244 marker->Draw();
2245 }
2246 myPadHisto->cd();
2247
2248
2249
2250 TLegend * leg = 0;
2251
2252 if (drawExtrap) leg = new TLegend(0.123494, 0.400358, 0.331325, 0.534512,NULL,"brNDC");
2253 else leg = new TLegend(0.123494, 0.395097, 0.332329, 0.474011, NULL, "brNDC");
2254 leg->SetBorderSize(0);
2255 leg->SetTextFont(43);
2256 leg->SetTextSize(14);
2257 leg->SetLineColor(1);
2258 leg->SetLineStyle(1);
2259 leg->SetLineWidth(1);
2260 leg->SetFillColor(0);
2261 leg->SetFillStyle(1001);
2262 TLegendEntry * entry=leg->AddEntry("TMarker","Not in fit","p");
2263 entry->SetLineColor(1);
2264 entry->SetLineStyle(1);
2265 entry->SetLineWidth(1);
2266 entry->SetMarkerColor(1);
2267 entry->SetMarkerStyle(markerNoFit);
2268 entry->SetMarkerSize(1.2);
2269 entry->SetTextFont(43);
2270
2271 if(drawExtrap) {
2272 entry=leg->AddEntry("TMarker","Extrapolated","p");
2273 entry->SetLineColor(1);
2274 entry->SetLineStyle(1);
2275 entry->SetLineWidth(1);
2276 entry->SetMarkerColor(1);
2277 entry->SetMarkerStyle(markerExtrap);
2278 entry->SetMarkerSize(1.2);
2279 entry->SetTextFont(43);
2280 }
2281 leg->Draw();
2282 leg->SetMargin(0.1);
2283}