macros/wEEC_3D_unfolding_kinematics/draw_misses_fakes.C

0001 #include "analysisHelper.h"
0002
0003 void draw_misses_fakes(float towCut = 0.25, Mode mode = Mode::kFull)
0004 {
0005
0006     std::string towCutStr = std::format("{}",towCut);
0007     std::replace(towCutStr.begin(), towCutStr.end(), '.', 'p');
0008     std::string outDir = std::format("/sphenix/tg/tg01/jets/bkimelman/VandyDSTs_wEEC_3D_unfolding_kinematics_Apr30_2026_{}/",towCutStr);
0009
0010     TFile *fResp = new TFile(std::format("{}/response-all-fullClosure.root",outDir).c_str(),"READ");
0011     if(!fResp)
0012     {
0013         //std::cerr << "Response file " << respFile << " does not exit" << std::endl;
0014         return;
0015     }
0016
0017     //std::vector<TH1D*> hTruth(nDphi, nullptr);
0018     //std::vector<TH1D*> hReco(nDphi, nullptr);
0019     //std::vector<TH1D*> hMiss(nDphi, nullptr);
0020     //std::vector<TH1D*> hFake(nDphi, nullptr);
0021
0022     TH1D *hMisses = new TH1D("misses","Misses;#Delta#phi;Miss Rate",nDphi, &dPhiBins[0]);
0023     hMisses->SetMarkerStyle(20);
0024     hMisses->SetMarkerColor(kRed+1);
0025
0026     TH1D *hFakes = new TH1D("fakes","Fakes;#Delta#phi;Fake Rate",nDphi, &dPhiBins[0]);
0027     hFakes->SetMarkerStyle(21);
0028     hFakes->SetMarkerColor(kBlue+1);
0029
0030     vector<TGraphAsymmErrors*> hTruth_pw(nDphi, nullptr);
0031     vector<TGraphAsymmErrors*> hReco_pw(nDphi, nullptr);
0032     vector<TGraphAsymmErrors*> hMisses_pw(nDphi, nullptr);
0033     vector<TGraphAsymmErrors*> hFakes_pw(nDphi, nullptr);
0034     for(int k=0; k<nDphi; ++k) {
0035         hTruth_pw[k] = new TGraphAsymmErrors();
0036         hTruth_pw[k]->SetMarkerStyle(20);
0037         hTruth_pw[k]->SetMarkerColor(kBlue+1);
0038         hTruth_pw[k]->SetLineColor(kBlue+1);
0039
0040         hReco_pw[k] = new TGraphAsymmErrors();
0041         hReco_pw[k]->SetMarkerStyle(20);
0042         hReco_pw[k]->SetMarkerColor(kBlue+1);
0043         hReco_pw[k]->SetLineColor(kBlue+1);
0044
0045         hMisses_pw[k] = new TGraphAsymmErrors();
0046         hMisses_pw[k]->SetMarkerStyle(21);
0047         hMisses_pw[k]->SetMarkerColor(kRed+1);
0048         hMisses_pw[k]->SetLineColor(kRed+1);
0049
0050         hFakes_pw[k] = new TGraphAsymmErrors();
0051         hFakes_pw[k]->SetMarkerStyle(21);
0052         hFakes_pw[k]->SetMarkerColor(kRed+1);
0053         hFakes_pw[k]->SetLineColor(kRed+1);
0054     }
0055
0056     vector<double> tMin(nDphi, 1e20);
0057     vector<double> tMax(nDphi, 0.0);
0058     vector<double> rMin(nDphi, 1e20);
0059     vector<double> rMax(nDphi, 0.0);
0060
0061     for(int k=0; k<nDphi; ++k) {
0062         TH1D *hTruth = (TH1D*)fResp->Get(std::format("hWEEC3D_truth_{}",k).c_str());
0063         hTruth->SetDirectory(0);
0064         double truthIntegral = hTruth->Integral();
0065
0066         TH1D *hReco = (TH1D*)fResp->Get(std::format("hWEEC3D_meas_{}",k).c_str());
0067         hReco->SetDirectory(0);
0068         double recoIntegral = hReco->Integral();
0069
0070         TH1D *hMiss = (TH1D*)fResp->Get(std::format("hWEEC3D_misses_{}",k).c_str());
0071         hMiss->SetDirectory(0);
0072         double missIntegral = hMiss->Integral();
0073
0074         TH1D *hFake = (TH1D*)fResp->Get(std::format("hWEEC3D_fakes_{}",k).c_str());
0075         hFake->SetDirectory(0);
0076         double fakeIntegral = hFake->Integral();
0077
0078         double missRate = (truthIntegral > 0) ? missIntegral / truthIntegral : 0.0;
0079         hMisses->SetBinContent(k+1, missRate);
0080
0081         double fakeRate = (recoIntegral > 0) ? fakeIntegral / recoIntegral : 0.0;
0082         hFakes->SetBinContent(k+1, fakeRate);
0083
0084         int iLreco = FindBin(35, recoLeadPtBins);
0085         int iSreco = FindBin(35, recoSublPtBins);
0086         int iLtrue = FindBin(35, trueLeadPtBins);
0087         int iStrue = FindBin(35, trueSublPtBins);
0088
0089         for(int pw=0; pw<nPairWeight; pw++) {
0090             int rf = RecoFlat3DIndex(iLreco, iSreco, pw);
0091             int tf = TrueFlat3DIndex(iLtrue, iStrue, pw);
0092
0093             double truthPW = hTruth->GetBinContent(tf);
0094             double recoPW = hReco->GetBinContent(rf);
0095             double missPW = hMiss->GetBinContent(tf);
0096             double fakePW = hFake->GetBinContent(rf);
0097
0098             double pwBinCenter = 0.5*(pairWeightBins[pw] + pairWeightBins[pw+1]);
0099             double pWELow = pwBinCenter - pairWeightBins[pw];
0100             double pWEHigh = pairWeightBins[pw+1] - pwBinCenter;
0101
0102             if(truthPW > 0)
0103             {
0104                 if(truthPW < tMin[k]) tMin[k] = truthPW;
0105                 if(truthPW > tMax[k]) tMax[k] = truthPW;
0106                 hTruth_pw[k]->AddPoint(pwBinCenter, truthPW);
0107                 hTruth_pw[k]->SetPointError(hTruth_pw[k]->GetN()-1, pWELow, pWEHigh, 0.0, 0.0);
0108                 hMisses_pw[k]->AddPoint(0.5*(pairWeightBins[pw] + pairWeightBins[pw+1]), missPW / truthPW);
0109                 hMisses_pw[k]->SetPointError(hMisses_pw[k]->GetN()-1, pWELow, pWEHigh, 0.0, 0.0);
0110             }
0111             if(recoPW > 0)
0112             {
0113                 if(recoPW < rMin[k]) rMin[k] = recoPW;
0114                 if(recoPW > rMax[k]) rMax[k] = recoPW;
0115                 hReco_pw[k]->AddPoint(pwBinCenter, recoPW);
0116                 hReco_pw[k]->SetPointError(hReco_pw[k]->GetN()-1, pWELow, pWEHigh, 0.0, 0.0);
0117                 hFakes_pw[k]->AddPoint(pwBinCenter, fakePW / recoPW);
0118                 hFakes_pw[k]->SetPointError(hFakes_pw[k]->GetN()-1, pWELow, pWEHigh, 0.0, 0.0);
0119             }
0120         }
0121
0122         delete hTruth; delete hReco;
0123         delete hMiss; delete hFake;
0124     }
0125
0126     TCanvas *c1 = new TCanvas();
0127     gStyle->SetOptStat(0);
0128
0129     TH2D *tmp = new TH2D("tmp","Fake and Miss Rates;#Delta#phi;Rate",1,0,TMath::Pi(),1,0,1);
0130     tmp->Draw();
0131
0132     TLine *l25 = new TLine(0,0.25,TMath::Pi(),0.25);
0133     l25->SetLineColor(kBlack);
0134     l25->SetLineStyle(2);
0135     l25->Draw("same");
0136
0137     TLine *l5 = new TLine(0,0.5,TMath::Pi(),0.5);
0138     l5->SetLineColor(kBlack);
0139     l5->SetLineStyle(2);
0140     l5->Draw("same");
0141
0142     TLine *l75 = new TLine(0,0.75,TMath::Pi(),0.75);
0143     l75->SetLineColor(kBlack);
0144     l75->SetLineStyle(2);
0145     l75->Draw("same");
0146
0147     hMisses->Draw("PSAME");
0148     hFakes->Draw("PSAME");
0149
0150     TLegend *leg = new TLegend(0.7,0.7,0.85,0.85);
0151     leg->AddEntry(hMisses,"Misses","P");
0152     leg->AddEntry(hFakes,"Fakes","P");
0153     leg->Draw();
0154
0155     c1->SaveAs(std::format("{}/Plots/misses_and_fakes-{}-{}.png",outDir,ModeLabel(mode),towCutStr).c_str());
0156
0157     TCanvas *c2 = new TCanvas("c2","",2*957,957);
0158     //TH2D *tmp2 = new TH2D("tmp2","Fake and Miss Rates;#Delta#phi;Rate",1,1e-6,2.0,1,0,1);
0159     c2->SaveAs(std::format("{}/Plots/misses_and_fakes_dPhi-{}-{}.pdf[",outDir,ModeLabel(mode),towCutStr).c_str());
0160
0161     for(int k=0; k<nDphi; ++k) {
0162         c2->Clear();
0163         c2->Divide(2,1);
0164         //c2->cd(1)->SetLogx();
0165         //c2->cd(1)->SetLogy();
0166
0167         c2->cd(1);
0168         auto *p1L = new TPad("p1L","",0,0,1,1);
0169         p1L->SetLogy();
0170         p1L->SetLogx();
0171         p1L->Draw();
0172         p1L->cd();
0173
0174         TH1F *hlT = p1L->DrawFrame(1e-6,0.1*tMin[k],2.0,10*tMax[k]);
0175         hlT->SetXTitle("#Delta#phi");
0176         hlT->SetYTitle("x-sec counts");
0177         hlT->SetTitle(std::format("Truth Distribution and Miss Rates p_{{T,tower}}>{} GeV. #Delta#phi bin {} {:.2f}-{:.2f} 30<p_{{T,lead}}<40, 30<p_{{T,sub}}<40",towCut,k,dPhiBins[k],dPhiBins[k+1]).c_str());
0178         hlT->GetYaxis()->SetAxisColor(kBlue+1);
0179         hlT->GetYaxis()->SetLabelColor(kBlue+1);
0180         hlT->GetYaxis()->SetTitleColor(kBlue+1);
0181         hlT->Draw();
0182
0183         hTruth_pw[k]->Draw("PSAME");
0184
0185         c2->cd(1);
0186         TPad *p1R = new TPad("p1R","",0,0,1,1);
0187         p1R->SetFillStyle(4000);
0188         p1R->SetFillColor(0);
0189         p1R->SetFrameFillStyle(4000);
0190         p1R->SetFrameLineColor(0);
0191         p1R->SetFrameBorderMode(0);
0192         p1R->SetLogx();
0193         p1R->Draw();
0194         p1R->cd();
0195
0196         TH1F *hrT = p1R->DrawFrame(1e-6,0,2.0,1.0);
0197         hrT->GetXaxis()->SetLabelOffset(99);
0198         hrT->GetYaxis()->SetLabelOffset(99);
0199         hrT->GetXaxis()->SetAxisColor(0);
0200         hrT->GetXaxis()->SetTickLength(0);
0201         hrT->GetYaxis()->SetAxisColor(0);
0202         hrT->GetYaxis()->SetTickLength(0);
0203
0204         TGaxis *abT = new TGaxis(1e-6,0,2.0,0,1e-6,2,510,"G");
0205         abT->SetTickLength(0);
0206         abT->SetLabelOffset(99);
0207         abT->Draw();
0208
0209         TGaxis *atT = new TGaxis(1e-6,1,2.0,1,1e-6,2,0,"G");
0210         atT->SetTickLength(0);
0211         atT->SetLabelOffset(99);
0212         atT->Draw();
0213
0214         TGaxis *alT = new TGaxis(1e-6,0,1e-6,1,0.1*tMin[k],10*tMax[k],510,"G");
0215         alT->SetTickLength(0);
0216         alT->SetLabelOffset(99);
0217         alT->SetLineColor(kBlue+1);
0218         alT->Draw();
0219
0220         TGaxis *arT = new TGaxis(2.0,0,2.0,1.0,0,1,510,"+L");
0221         arT->SetLineColor(kRed+1);
0222         arT->SetLabelColor(kRed+1);
0223         arT->SetTitleColor(kRed+1);
0224         arT->SetTitle("Miss Rate");
0225         arT->Draw();
0226
0227         hMisses_pw[k]->Draw("PSAME");
0228
0229
0230         //
0231         c2->cd(2);
0232         auto *p2L = new TPad("p2L","",0,0,1,1);
0233         p2L->SetLogy();
0234         p2L->SetLogx();
0235         p2L->Draw();
0236         p2L->cd();
0237
0238         TH1F *hlR = p2L->DrawFrame(1e-6,0.1*rMin[k],2.0,10*rMax[k]);
0239         hlR->SetXTitle("#Delta#phi");
0240         hlR->SetYTitle("x-sec counts");
0241         hlR->SetTitle(std::format("Reco Distribution and Fake Rates p_{{T,tower}}>{} GeV. #Delta#phi bin {} {:.2f}-{:.2f} 30<p_{{T,lead}}<40, 30<p_{{T,sub}}<40",towCut,k,dPhiBins[k],dPhiBins[k+1]).c_str());
0242         hlR->GetYaxis()->SetAxisColor(kBlue+1);
0243         hlR->GetYaxis()->SetLabelColor(kBlue+1);
0244         hlR->GetYaxis()->SetTitleColor(kBlue+1);
0245         hlR->Draw();
0246
0247         hReco_pw[k]->Draw("PSAME");
0248
0249         c2->cd(2);
0250         TPad *p2R = new TPad("p2R","",0,0,1,1);
0251         p2R->SetFillStyle(4000);
0252         p2R->SetFillColor(0);
0253         p2R->SetFrameFillStyle(4000);
0254         p2R->SetFrameLineColor(0);
0255         p2R->SetFrameBorderMode(0);
0256         p2R->SetLogx();
0257         p2R->Draw();
0258         p2R->cd();
0259
0260         TH1F *hrR = p2R->DrawFrame(1e-6,0,2.0,1.0);
0261         hrR->GetXaxis()->SetLabelOffset(99);
0262         hrR->GetYaxis()->SetLabelOffset(99);
0263         hrR->GetXaxis()->SetAxisColor(0);
0264         hrR->GetXaxis()->SetTickLength(0);
0265         hrR->GetYaxis()->SetAxisColor(0);
0266         hrR->GetYaxis()->SetTickLength(0);
0267
0268         TGaxis *abR = new TGaxis(1e-6,0,2.0,0,1e-6,2,510,"G");
0269         abR->SetTickLength(0);
0270         abR->SetLabelOffset(99);
0271         abR->Draw();
0272
0273         TGaxis *atR = new TGaxis(1e-6,1,2.0,1,1e-6,2,0,"G");
0274         atR->SetTickLength(0);
0275         atR->SetLabelOffset(99);
0276         atR->Draw();
0277
0278         TGaxis *alR = new TGaxis(1e-6,0,1e-6,1,0.1*rMin[k],10*rMax[k],510,"G");
0279         alR->SetTickLength(0);
0280         alR->SetLabelOffset(99);
0281         alR->SetLineColor(kBlue+1);
0282         alR->Draw();
0283
0284         TGaxis *arR = new TGaxis(2.0,0,2.0,1.0,0,1,510,"+L");
0285         arR->SetLineColor(kRed+1);
0286         arR->SetLabelColor(kRed+1);
0287         arR->SetTitleColor(kRed+1);
0288         arR->SetTitle("Fake Rate");
0289         arR->Draw();
0290
0291
0292         hFakes_pw[k]->Draw("PSAME");
0293
0294         c2->SaveAs(std::format("{}/Plots/misses_and_fakes_dPhi-{}-{}.pdf",outDir,ModeLabel(mode),towCutStr).c_str());
0295
0296     }
0297     c2->SaveAs(std::format("{}/Plots/misses_and_fakes_dPhi-{}-{}.pdf]",outDir,ModeLabel(mode),towCutStr).c_str());
0298
0299
0300
0301 }