1 void CompareCorrectedSpectra(const char *datafilea, const char *datafileb){
3 gROOT->SetStyle("Plain");
5 gStyle->SetOptStat(1111);
6 gStyle->SetPadBorderMode(0);
7 gStyle->SetCanvasColor(10);
8 gStyle->SetPadLeftMargin(0.13);
9 gStyle->SetPadRightMargin(0.13);
13 TFile *filea = TFile::Open(datafilea);
14 TGraphErrors *correctedSpectrum_a = (TGraphErrors *) filea->Get("AlltogetherSpectrum");
16 TH1D *histo = (TH1D*) filea->Get("RatioUnfoldingAlltogetherSpectrum");
17 histo->SetName("historatio");
18 TH1D *histoa = (TH1D*) histo->Clone();
21 TH1D *histob = (TH1D*) histo->Clone();
26 TFile *fileb = TFile::Open(datafileb);
27 TGraphErrors *correctedSpectrum_b = (TGraphErrors *) fileb->Get("AlltogetherSpectrum");
31 correctedSpectrum_a->SetMarkerStyle(24);
32 correctedSpectrum_a->SetMarkerColor(1);
33 correctedSpectrum_a->SetLineColor(1);
35 correctedSpectrum_b->SetMarkerStyle(27);
36 correctedSpectrum_b->SetMarkerColor(4);
37 correctedSpectrum_b->SetLineColor(4);
41 TCanvas *c1 = new TCanvas("CorrectedSpectrum","CorrectedSpectrum",800,800);
43 TH1D *total = new TH1D("total","",1,0.38,4.3);
44 total->SetMaximum(1.0);
45 total->SetMinimum(1.0e-09);
46 total->SetXTitle("p_{T} [GeV/c]");
47 total->SetYTitle("1/2#pip_{T} d^{2}N/dp_{T}dy [GeV/c]^{-2}, |#eta| < 0.8");
48 total->SetTitleOffset(1.5,"Y");
50 gPad->SetLeftMargin(0.13);
53 gPad->SetFillColor(10);
54 gPad->SetFrameFillColor(10);
55 gStyle->SetOptStat(0);
56 gStyle->SetOptFit(1111);
57 gStyle->SetOptTitle(0);
58 correctedSpectrum_a->Draw("Psame");
59 correctedSpectrum_b->Draw("Psame");
60 TLegend *leg = new TLegend(0.25,0.825,0.48,0.875);
61 leg->AddEntry(correctedSpectrum_a,"Corrected spectrum a","lep");
62 leg->AddEntry(correctedSpectrum_b,"Corrected spectrum b","lep");
74 double xa,ya,xb,yb,eya,exa,eyb,exb;
75 Int_t npointsa = correctedSpectrum_a->GetN();
76 Int_t npointsb = correctedSpectrum_b->GetN();
77 if(npointsa != npointsb) {
78 printf("Problem the two spectra have not the same number of points");
81 for(Int_t k = 0; k < npointsa; k++){
82 correctedSpectrum_a->GetPoint(k,xa,ya);
83 correctedSpectrum_b->GetPoint(k,xb,yb);
85 Double_t centerhisto = histoa->GetBinCenter(k+1);
86 //printf("bin center %f and center %f\n",centerhisto,xa);
87 histoa->SetBinContent(k+1,ya);
88 histob->SetBinContent(k+1,yb);
90 if(TMath::Abs(xa-xb) > 0.0001) {
91 printf("Problem the two spectra have not the same number of points");
94 eya = correctedSpectrum_a->GetErrorY(k);
95 exa = correctedSpectrum_a->GetErrorX(k);
96 eyb = correctedSpectrum_b->GetErrorY(k);
97 exb = correctedSpectrum_b->GetErrorX(k);
100 if(yb > 0.0) y[k] = ya/yb;
102 Double_t error = 0.0;
103 if((yb > 0.0) && (ya > 0.0)) {
104 error = TMath::Sqrt(TMath::Abs(eya*eya-eyb*eyb));
107 histoa->SetBinError(k+1,eya);
108 histob->SetBinError(k+1,eyb);
112 histo->Divide(histoa,histob,1.0,1.0,"B");
114 TGraph *gratio = new TGraph(npointsa,&x[0],&y[0]);
115 gratio->SetName("RatioOfSpectra");
117 TGraphErrors *gdiff = new TGraphErrors(npointsa,&x[0],&ry[0],&rex[0],&rey[0]);
118 gdiff->SetName("Difference");
120 TCanvas *c2 = new TCanvas("ratio","ratio",800,800);
122 TH1D *ratio = new TH1D("ratio","",1,0.38,10.0);
123 ratio->SetMaximum(1.5);
124 ratio->SetMinimum(0.5);
125 ratio->SetXTitle("p_{T} [GeV/c]");
126 ratio->SetYTitle("Ratio of corrected spectra a/b");
127 ratio->SetTitleOffset(1.5,"Y");
129 histo->SetLineColor(2);
131 gPad->SetLeftMargin(0.13);
135 gPad->SetFillColor(10);
136 gPad->SetFrameFillColor(10);
137 gStyle->SetOptStat(0);
138 gStyle->SetOptFit(1111);
139 gStyle->SetOptTitle(0);
140 gratio->SetMarkerStyle(24);
141 gratio->SetMarkerColor(4);
142 gratio->SetLineColor(4);
145 TCanvas *c3 = new TCanvas("Difference","Difference",800,800);
147 TH1D *diff = new TH1D("difference","",1,0.38,4.3);
148 diff->SetMaximum(1.0e-04);
149 diff->SetMinimum(-1.0e-04);
150 diff->SetXTitle("p_{T} [GeV/c]");
151 diff->SetYTitle("Difference of corrected spectra a-b");
152 diff->SetTitleOffset(1.5,"Y");
154 gPad->SetLeftMargin(0.13);
158 gPad->SetFillColor(10);
159 gPad->SetFrameFillColor(10);
160 gStyle->SetOptStat(0);
161 gStyle->SetOptFit(1111);
162 gStyle->SetOptTitle(0);
163 gdiff->SetMarkerStyle(24);
164 gdiff->SetMarkerColor(4);
165 gdiff->SetLineColor(4);