File indexing completed on 2023-03-17 11:14:51
0001 #include "TH1F.h"
0002 #include "TProfile.h"
0003 #include "TProfile.h"
0004 #include "TCanvas.h"
0005 #include "TFile.h"
0006 #include "TStyle.h"
0007 #include "TROOT.h"
0008
0009 #include <iostream>
0010
0011
0012
0013 class histos
0014 {
0015 public:
0016 histos(TFile * inputFile)
0017 {
0018 mass = dynamic_cast<TH1F*>(inputFile->Get("hRecBestResAllEvents_Mass"));
0019 TDirectory * MassPdir = dynamic_cast<TDirectory*>(inputFile->Get("Mass_P"));
0020 massProb = dynamic_cast<TProfile*>(MassPdir->Get("Mass_PProf"));
0021
0022
0023 TDirectory * massFinePdir = dynamic_cast<TDirectory*>(inputFile->Get("Mass_fine_P"));
0024 massFineProb = dynamic_cast<TProfile*>(massFinePdir->Get("Mass_fine_PProf"));
0025
0026
0027 likePt = dynamic_cast<TProfile*>(inputFile->Get("hLikeVSMu_LikelihoodVSPt_prof"));
0028 likePhi = dynamic_cast<TProfile*>(inputFile->Get("hLikeVSMu_LikelihoodVSPhi_prof"));
0029 likeEta = dynamic_cast<TProfile*>(inputFile->Get("hLikeVSMu_LikelihoodVSEta_prof"));
0030 }
0031 TH1F * mass;
0032 TProfile * massProb;
0033 TProfile * massFineProb;
0034 TProfile * likePt;
0035 TProfile * likePhi;
0036 TProfile * likeEta;
0037 };
0038
0039
0040 int getXbins(const TH1 * h, const double & xMin, const double & xMax)
0041 {
0042
0043 const TAxis * xAxis = h->GetXaxis();
0044 double xAxisMin = xAxis->GetXmin();
0045
0046 double binWidth = xAxis->GetBinWidth(1);
0047 int xMinBin = int((xMin - xAxisMin)/binWidth) + 1;
0048 int xMaxBin = int((xMax - xAxisMin)/binWidth) + 1;
0049
0050
0051
0052 return( xMaxBin - xMinBin );
0053 }
0054
0055
0056 void drawMasses(const double ResMass, const double ResHalfWidth, histos & h, const int ires, const int rebin = 1)
0057 {
0058 TH1F * mass = (TH1F*)h.mass->Clone();
0059 TProfile * massProb = 0;
0060 mass->Rebin(rebin);
0061
0062 if( ires == 0 ) massProb = (TProfile*)h.massProb->Clone();
0063 else massProb = (TProfile*)h.massFineProb->Clone();
0064 mass->SetAxisRange(ResMass - ResHalfWidth, ResMass + ResHalfWidth);
0065 mass->SetLineColor(kRed);
0066
0067
0068 int massXbins = getXbins( mass, (ResMass - ResHalfWidth), (ResMass + ResHalfWidth) );
0069 int massProbXbins = getXbins( massProb, (ResMass - ResHalfWidth), (ResMass + ResHalfWidth) );
0070
0071 if( massProbXbins > massXbins && massXbins != 0 ) {
0072 std::cout << "massProbXbins("<<massProbXbins<<") > " << "massXbins("<<massXbins<<")" << std::endl;
0073 std::cout << "massProb = " << massProb << std::endl;
0074 std::cout << "mass = " << mass << std::endl;
0075 massProb->Rebin(massProbXbins/massXbins);
0076 }
0077 else if( massXbins > massProbXbins && massProbXbins != 0 ) {
0078 std::cout << "massXbins("<<massXbins<<") > " << "massProbXbins("<<massProbXbins<<")" << std::endl;
0079 std::cout << "massProb = " << massProb << std::endl;
0080 std::cout << "mass = " << mass << std::endl;
0081 mass->Rebin(massXbins/massProbXbins);
0082 }
0083
0084 mass->SetAxisRange(ResMass - ResHalfWidth, ResMass + ResHalfWidth);
0085 massProb->SetAxisRange(ResMass - ResHalfWidth, ResMass + ResHalfWidth);
0086
0087 double massProbIntegral = massProb->Integral("width");
0088
0089 double normFactor = 1.;
0090 if( massProbIntegral != 0 ) normFactor = (mass->Integral("width"))/massProbIntegral;
0091
0092
0093 massProb->SetLineColor(kBlue);
0094 massProb->Scale(normFactor);
0095 if( ires == 3 ) {
0096 std::cout << "massProbIntegralWidth = " << massProb->Integral("width") << std::endl;
0097 std::cout << "massIntegralWidth = " << mass->Integral("width") << std::endl;
0098 std::cout << "massProbIntegral = " << massProb->Integral() << std::endl;
0099 std::cout << "massIntegral = " << mass->Integral() << std::endl;
0100 }
0101
0102 mass->SetMarkerColor(kRed);
0103
0104 if( (massProb->GetMaximum()) > mass->GetMaximum() ) {
0105 massProb->Draw("PE");
0106 mass->Draw("PESAMEHISTO");
0107 }
0108 else {
0109 mass->Draw("PE");
0110 massProb->Draw("PESAMEHISTO");
0111 }
0112 }
0113
0114
0115 void drawLike(TProfile * likeHisto1, TProfile * likeHisto2)
0116 {
0117 likeHisto1->SetMinimum(-8);
0118 likeHisto1->SetMaximum(-4);
0119 likeHisto1->Draw("");
0120 likeHisto2->SetLineColor(kRed);
0121 likeHisto2->Draw("SAME");
0122 }
0123
0124
0125 void Plot_mass(const TString & fileNameBefore = "0", const TString & fileNameAfter = "1", const int rebin = 1) {
0126
0127 gROOT->SetBatch(true);
0128 gStyle->SetOptStat ("111111");
0129 gStyle->SetOptFit (1);
0130
0131 double ResHalfWidth[6] = {20., 0.5, 0.5, 0.5, 0.2, 0.2};
0132
0133
0134 double ResMass[6] = {90.986, 10.3552, 10.0233, 9.4603, 3.68609, 3.0969};
0135
0136 TFile * inputFile1 = new TFile(fileNameBefore+"_MuScleFit.root", "READ");
0137 TFile * inputFile2 = new TFile(fileNameAfter+"_MuScleFit.root", "READ");
0138 histos histos1(inputFile1);
0139 histos histos2(inputFile2);
0140
0141 TFile * outputFile = new TFile("plotMassOutput.root", "RECREATE");
0142
0143 TCanvas * Allres = new TCanvas ("Allres", "All resonances", 600, 600);
0144 Allres->Divide (2,3);
0145 for (int ires=0; ires<6; ires++) {
0146 Allres->cd(ires+1);
0147 drawMasses(ResMass[ires], ResHalfWidth[ires], histos1, ires, rebin);
0148 }
0149 Allres->Write();
0150 TCanvas * AllresTogether = new TCanvas ("AllresTogether", "All resonances together", 600, 600);
0151 AllresTogether->Divide(2,2);
0152
0153 AllresTogether->cd(1);
0154 drawMasses(ResMass[5], 2, histos1, 3, rebin);
0155
0156 AllresTogether->cd(2);
0157 drawMasses(ResMass[2], 2, histos1, 3, rebin);
0158
0159 AllresTogether->cd(3);
0160 drawMasses(6., 6, histos1, 3, rebin);
0161
0162 AllresTogether->cd(4);
0163 drawMasses(50., 50, histos1, 0, rebin);
0164 AllresTogether->Write();
0165
0166 TCanvas * Allres2 = new TCanvas ("Allres2", "All resonances", 600, 600);
0167 Allres2->Divide (2,3);
0168 for (int ires=0; ires<6; ires++) {
0169 Allres2->cd(ires+1);
0170 drawMasses(ResMass[ires], ResHalfWidth[ires], histos2, ires, rebin);
0171 }
0172 Allres2->Write();
0173 TCanvas * AllresTogether2 = new TCanvas ("AllresTogether2", "All resonances together", 600, 600);
0174 AllresTogether2->Divide(2,2);
0175
0176 AllresTogether2->cd(1);
0177 drawMasses(ResMass[5], 2, histos2, 3, rebin);
0178
0179 AllresTogether2->cd(2);
0180 drawMasses(ResMass[2], 2, histos2, 3, rebin);
0181
0182 AllresTogether2->cd(3);
0183 drawMasses(6., 6, histos2, 3, rebin);
0184
0185 AllresTogether2->cd(4);
0186 drawMasses(50., 50, histos2, 0, rebin);
0187 AllresTogether2->Write();
0188
0189 TCanvas * LR = new TCanvas ("LR", "Likelihood and resolution before and after corrections", 600, 600 );
0190 LR->Divide (2,3);
0191
0192 LR->cd(1);
0193 drawLike(histos1.likePt, histos2.likePt);
0194 LR->cd(3);
0195 drawLike(histos1.likeEta, histos2.likeEta);
0196 LR->cd(5);
0197 drawLike(histos1.likePhi, histos2.likePhi);
0198
0199 LR->Write();
0200
0201 TCanvas * G = new TCanvas ("G", "Mass before and after corrections", 600, 600 );
0202 G->Divide (3,3);
0203
0204 TH1F * MuZ = histos1.mass;
0205 TH1F * McZ = histos2.mass;
0206
0207 G->cd(1);
0208 MuZ->SetAxisRange (2.8, 3.4);
0209 MuZ->SetLineColor (kBlue);
0210 McZ->SetLineColor (kRed);
0211 MuZ->DrawCopy ("HISTO");
0212 McZ->DrawCopy ("SAME");
0213
0214 G->cd(2);
0215 MuZ->SetAxisRange (9., 10.);
0216 McZ->SetAxisRange (9., 10.);
0217 MuZ->DrawCopy ("HISTO");
0218 McZ->DrawCopy ("SAME");
0219
0220 G->cd(3);
0221 MuZ->SetAxisRange (70., 110.);
0222 McZ->SetAxisRange (70., 110.);
0223 MuZ->DrawCopy ("HISTO");
0224 McZ->DrawCopy ("SAME");
0225
0226 G->cd(4);
0227 MuZ->SetAxisRange (2.8, 3.4);
0228 MuZ->Fit ("gaus", "", "", 3., 3.2);
0229 MuZ->DrawCopy();
0230 G->cd(7);
0231 McZ->SetAxisRange (2.8, 3.4);
0232 McZ->Fit ("gaus", "", "", 3., 3.2);
0233 McZ->DrawCopy();
0234 G->cd(5);
0235 MuZ->SetAxisRange (9., 10.);
0236 MuZ->Fit ("gaus", "", "", 9., 10.);
0237 MuZ->DrawCopy();
0238 G->cd(8);
0239 McZ->SetAxisRange (9., 10.);
0240 McZ->Fit ("gaus", "", "", 9., 10.);
0241 McZ->DrawCopy();
0242 G->cd(6);
0243 MuZ->SetAxisRange (70., 110.);
0244 MuZ->Fit ("gaus", "", "", 70., 110.);
0245 MuZ->DrawCopy();
0246 G->cd(9);
0247 McZ->SetAxisRange (70., 110.);
0248 McZ->Fit ("gaus", "", "", 70., 110.);
0249 McZ->DrawCopy();
0250
0251 G->Write();
0252
0253 TCanvas * Spectrum = new TCanvas ("Spectrum", "Mass before and after corrections", 600, 600 );
0254 Spectrum->cd(1);
0255 MuZ->SetAxisRange(1.,15.);
0256 McZ->SetAxisRange(1.,15.);
0257 McZ->SetLineColor(kRed);
0258 McZ->DrawCopy("HISTO");
0259 Spectrum->cd(2);
0260 MuZ->SetLineColor(kBlack);
0261 MuZ->DrawCopy("SAMEHISTO");
0262
0263 Spectrum->Write();
0264
0265 outputFile->Close();
0266 }
0267