EICMCAnalysis/analysis_eicsmear/eicsmear_dvmp_plot.C

0001 #include <vector>
0002 int eicsmear_dvmp_plot()
0003 {
0004
0005   /* PRELIMINARY ROOT STUFF */
0006
0007
0008   /* Loading libraries and setting sphenix style */
0009   gSystem->Load("libeicsmear");
0010   gROOT->LoadMacro("/sphenix/user/gregtom3/SBU/research/macros/macros/sPHENIXStyle/sPhenixStyle.C");
0011   SetsPhenixStyle();
0012
0013   /* Open TFile */
0014   TFile *input = new TFile("/sphenix/user/gregtom3/data/Fall2018/LOI_smear_studies/DVMP.root","READ");
0015
0016   /* Open histograms in TFile */
0017   TH1F *h1_invariant_mass = (TH1F*)input->Get("h1_invariant_mass");
0018   TH1F *h1_invariant_mass_comb = (TH1F*)input->Get("h1_invariant_mass_comb");
0019   TH1F *h1_invariant_mass_alone = (TH1F*)input->Get("h1_invariant_mass_alone");
0020
0021
0022   TCanvas *c_inv = new TCanvas("c_inv","",700,500);
0023   h1_invariant_mass_comb->SetLineColor(kGray);
0024   h1_invariant_mass_comb->SetFillColor(kGray);
0025   h1_invariant_mass->GetXaxis()->SetTitle("e^{-}e^{+} Invariant Mass [GeV/c^{2}]");
0026   h1_invariant_mass->GetYaxis()->SetTitle("Entries");
0027
0028   //Scale for Luminosity
0029   h1_invariant_mass->Scale(10/1.55);
0030   h1_invariant_mass_comb->Scale(10/1.55);
0031   h1_invariant_mass_alone->Scale(10/1.55);
0032
0033
0034
0035
0036   // Legend text
0037   auto legend = new TLegend(0.6,0.65,0.8,0.75);
0038   legend->SetTextSize(0.03);
0039   legend->AddEntry(h1_invariant_mass,"All e^{+}e^{-} pairs","l");
0040   legend->AddEntry(h1_invariant_mass_comb,"Combinatorial background","l");
0041
0042   //h1_invariant_mass->Draw();
0043   //h1_invariant_mass_comb->Draw("same");
0044   //h1_invariant_mass->Draw("same");
0045   legend->Draw("same");
0046   gPad->SetLogy();
0047   gPad->RedrawAxis();
0048   c_inv->SaveAs("inv_mass_poly.eps");
0049
0050   TCanvas *c_sing = new TCanvas("c_sing","",700,500);
0051   double x[100]; double y[100]; double ex[100]; double ey[100];
0052   for(int entry = 0; entry<100; entry++)
0053     {
0054       x[entry] = h1_invariant_mass->GetXaxis()->GetBinCenter(entry);
0055       if(entry==0)
0056     y[entry]=h1_invariant_mass->GetBinContent(1);
0057       else
0058     y[entry] = h1_invariant_mass->GetBinContent(entry);
0059       ey[entry] = sqrt(h1_invariant_mass->GetBinContent(entry));
0060       ex[entry] = 0;
0061
0062     }
0063
0064   auto gr = new TGraph(100,x,y);
0065   auto gr_error = new TGraphErrors(100,x,y,ex,ey);
0066   gr_error->SetMarkerSize(1);
0067   gr_error->SetMarkerColor(kBlack);
0068   gr->SetMarkerSize(.5);
0069   gr->SetMarkerColor(kBlue);
0070
0071   gr->Draw("AP");
0072
0073   gr->GetXaxis()->SetTitle("e^{+}e^{-} Invariant Mass [GeV/c^{2}]");
0074   gr->GetYaxis()->SetTitle("Entries");
0075   //gr->GetYaxis()->SetTitleOffset(1.6);
0076   auto th2 = new TLatex(3.16,38000,"Luminosity 10 fb^{-1}");
0077   th2->SetTextSize(0.04);
0078   th2->Draw("same");
0079   h1_invariant_mass_comb->Draw("same");
0080   c_sing->SaveAs("inv_mass_sing.eps");
0081
0082   TCanvas *c_nils = new TCanvas("c_nils","",700,500);
0083   //gr_error->Draw("ap");
0084   h1_invariant_mass_alone->SetMarkerColor(kWhite);
0085   h1_invariant_mass_alone->SetLineColor(kWhite);
0086   gr_error->GetXaxis()->SetRangeUser(2.95,3.25);
0087   gr_error->GetXaxis()->SetTitle("e^{+}e^{-} Invariant Mass [GeV/c^{2}]");
0088   gr_error->GetYaxis()->SetTitle("Entries");
0089   gr_error->GetYaxis()->SetNdivisions(505);
0090   gr_error->Draw("aP");
0091   // h1_invariant_mass_alone->GetYaxis()->SetRangeUser(0,50000);
0092   /*h1_invariant_mass_alone->Draw();
0093   h1_invariant_mass_comb->SetLineWidth(0);
0094   h1_invariant_mass_comb->Draw("same");
0095   gr->Draw("p");
0096   gPad->RedrawAxis();
0097   */
0098   //h1_invariant_mass_alone->Draw();
0099   using namespace RooFit;
0100   //double hmin = h1_invariant_mass_alone->GetXaxis()->GetXmin();
0101   //double hmax = h1_invariant_mass_alone->GetXaxis()->GetXmax();
0102
0103   double hmin = 3.075;
0104   double hmax = 3.125;
0105   // Declare observable x
0106   RooRealVar xx("xx","e^{+}e^{-} Invariant Mass [GeV/c^{2}]",hmin,hmax);
0107   RooDataHist dh("dh","Entries",xx,h1_invariant_mass);
0108   RooPlot* frame = xx.frame("Z mass");
0109   //h1_invariant_mass->Draw();
0110   h1_invariant_mass->GetYaxis()->SetNdivisions(505);
0111   //frame->Draw("same");
0112   //dh.plotOn(frame,MarkerColor(4),MarkerSize(2),MarkerStyle(20));  //this will show histogram data points on canvas
0113   //dh.statOn(frame);  //this will display hist stat on canvas
0114
0115   RooRealVar mean("mean","mean",3.096, 3, 1200.0);
0116   RooRealVar width("width","width",0.0926, 0.0900, 0.1);
0117   RooRealVar sigma("sigma","sigma",0.02, 0.0, 100.0);
0118   RooRealVar alpha("alpha","alpha",0.02, -1000.0, 1000.0);
0119   RooRealVar n("n","n",3, -5, 5);
0120   RooGaussian gauss("gauss","gauss",xx,mean,sigma);
0121   //RooCBShape gauss("gauss","gauss",xx,mean,sigma,alpha,n);
0122   //RooBreitWigner gauss("gauss","gauss",xx,mean,sigma);
0123   //RooVoigtian gauss("gauss","gauss",xx,mean,width,sigma);
0124   //h1_invariant_mass->SetLineColor(kWhite);
0125   //h1_invariant_mass->Draw("same");
0126   RooFitResult* filters = gauss.fitTo(dh,"qr");
0127   //h1_invariant_mass_alone->SetLineColor(kWhite);
0128   //h1_invariant_mass_alone->Draw("same");
0129   frame->GetYaxis()->SetTitle("Entries");
0130   frame->GetYaxis()->SetNdivisions(505);
0131
0132
0133   TF1 *gfit = new TF1("Gaussian","gaus",hmin,hmax); // Create the fit function
0134   h1_invariant_mass->Fit("Gaussian","RQ0"); //
0135   gfit->SetLineColor(kBlue);
0136   gfit->Draw("same");
0137
0138   //gauss.plotOn(frame,LineColor(4));//this will show fit overlay on canvas
0139   //gauss.paramOn(frame); //this will display the fit parameters on canvas
0140   dh.plotOn(frame,MarkerColor(4),MarkerSize(1),MarkerStyle(20));  //this will show histogram data points on canvas
0141   //frame->Draw();
0142   //h1_invariant_mass_comb->Draw("same");
0143   frame->Draw("same");
0144   auto legend = new TLegend(0.2,0.79,0.4,0.89);
0145   legend->SetTextSize(0.04);
0146   legend->AddEntry(h1_invariant_mass,"All e^{+}e^{-} pairs","l");
0147   legend->AddEntry(gr,"Bin points","p");
0148   legend->AddEntry(gfit,"Fitted Gaussian","l");
0149   //legend->AddEntry(h1_invariant_mass_comb,"Combinatorial Background","f");
0150   //legend->Draw("same");
0151   // Luminosity Text
0152   auto th1 = new TLatex(3.14,130000,"Luminosity 10 fb^{-1}");
0153   th1->SetTextSize(0.04);
0154   th1->Draw("same");
0155
0156   auto th3 = new TLatex(3.14,150000,"18x275 e+p");
0157   th3->SetTextSize(0.04);
0158   th3->Draw("same");
0159
0160
0161   auto th4 = new TLatex(3.14,140000,"DVMP Sartre");
0162   th4->SetTextSize(0.04);
0163   th4->Draw("same");
0164
0165   auto th5 = new TLatex(3.14,160000,"#bf{#it{sPHENIX}} Simulation");
0166   th5->SetTextSize(0.04);
0167   th5->Draw("same");
0168
0169
0170
0171   c_nils->SaveAs("dvmp_gaussian.eps");
0172
0173
0174   return 0;
0175 }
0176
0177 Double_t mybw(Double_t* x, Double_t* par)
0178 {
0179   Double_t arg1 = 14.0/22.0; // 2 over pi
0180   Double_t arg2 = par[1]*par[1]*par[2]*par[2]; //Gamma=par[1]  M=par[2]
0181   Double_t arg3 = ((x[0]*x[0]) - (par[2]*par[2]))*((x[0]*x[0]) - (par[2]*par[2]));
0182   Double_t arg4 = x[0]*x[0]*x[0]*x[0]*((par[1]*par[1])/(par[2]*par[2]));
0183   return par[0]*arg1*arg2/(arg3 + arg4);
0184 }