Back to home page

sPhenix code displayed by LXR
 
 

    

File indexing completed on 2025-08-06 08:11:52 

0001 // ----------------------------------------------------------------------------
0002 // 'MakeNewMatcherPlots.cxx'
0003 // Derek Anderson
0004 // 08.24.2023
0005 //
0006 // Short macro to make plots from the
0007 // output of the ClusMatchTree module.
0008 // ----------------------------------------------------------------------------
0009 
0010 #include <string>
0011 #include <vector>
0012 #include <utility>
0013 #include <iostream>
0014 #include "TH1.h"
0015 #include "TH2.h"
0016 #include "TFile.h"
0017 #include "TTree.h"
0018 #include "TError.h"
0019 
0020 using namespace std;
0021 
0022 
0023 
0024 void MakeNewMatcherPlots() {
0025 
0026   // io parameters
0027   const string sOutput("newMatcherPlots_oneMatchPerParticle_oddFrac02120.pt10n1evt500pim.d18m1y2024.root");
0028   const string sInTrue("input/merged/sPhenixG4_oneMatchPerParticle_newMatcher.pt10num1evt500pim.d4m1y2024.root");
0029   const string sInReco("input/merged/sPhenixG4_oneMatchPerParticle_newMatcher.pt10num1evt500pim.d4m1y2024.root");
0030   const string sTreeTrue("T");
0031   const string sTreeReco("T");
0032 
0033   // weird track parameters
0034   const pair<float, float> oddPtFrac = {0.2, 1.2};
0035 
0036   // lower verbosity
0037   gErrorIgnoreLevel = kError;
0038   cout << "\n  Staring new matcher plot script..." << endl;
0039 
0040   // open files
0041   TFile* fOutput = new TFile(sOutput.data(), "recreate");
0042   TFile* fInTrue = new TFile(sInTrue.data(), "read");
0043   TFile* fInReco = new TFile(sInReco.data(), "read");
0044   if (!fOutput || !fInTrue || !fInReco) {
0045      cerr << "PANIC: couldn't open a file!\n"
0046           << "       fOutput = " << fOutput << "fInTrue = " << fInTrue << ", fInReco = " << fInReco << "\n"
0047           << endl;
0048     return;
0049   }
0050   cout << "    Opened files." << endl;
0051 
0052   // grab input tree
0053   TTree* tInTrue = (TTree*) fInTrue -> Get(sTreeTrue.data());
0054   TTree* tInReco = (TTree*) fInReco -> Get(sTreeReco.data());
0055   if (!tInTrue || !tInReco) {
0056     cerr << "PANIC: couldn't grab input tree!\n"
0057          << "       sTreeTrue = " << sTreeTrue << ", tInTrue = " << tInTrue << "\n"
0058          << "       sTreeReco = " << sTreeReco << ", tInReco = " << tInReco << "\n"
0059          << endl;
0060     return;
0061   }
0062   cout << "    Grabbed input trees." << endl;
0063 
0064   // declare input leaves (ignore cluster branches for now)
0065   Int_t   tru_event;
0066   Int_t   tru_nphg4_part;
0067   Float_t tru_centrality;
0068   Int_t   tru_ntrackmatches;
0069   Int_t   tru_nphg4;
0070   Int_t   tru_nsvtx;
0071   Int_t   tru_trackid;
0072   Bool_t  tru_is_G4track;
0073   Bool_t  tru_is_Svtrack;
0074   Bool_t  tru_is_matched;
0075   Float_t tru_trkpt;
0076   Float_t tru_trketa;
0077   Float_t tru_trkphi;
0078   Int_t   tru_nclus;
0079   Int_t   tru_nclustpc;
0080   Int_t   tru_nclusmvtx;
0081   Int_t   tru_nclusintt;
0082   Float_t tru_matchrat;
0083   Float_t tru_matchrat_intt;
0084   Float_t tru_matchrat_mvtx;
0085   Float_t tru_matchrat_tpc;
0086 
0087   Int_t   rec_event;
0088   Int_t   rec_nphg4_part;
0089   Float_t rec_centrality;
0090   Int_t   rec_ntrackmatches;
0091   Int_t   rec_nphg4;
0092   Int_t   rec_nsvtx;
0093   Int_t   rec_trackid;
0094   Bool_t  rec_is_G4track;
0095   Bool_t  rec_is_Svtrack;
0096   Bool_t  rec_is_matched;
0097   Float_t rec_trkpt;
0098   Float_t rec_trketa;
0099   Float_t rec_trkphi;
0100   Int_t   rec_nclus;
0101   Int_t   rec_nclustpc;
0102   Int_t   rec_nclusmvtx;
0103   Int_t   rec_nclusintt;
0104   Float_t rec_matchrat;
0105   Float_t rec_matchrat_intt;
0106   Float_t rec_matchrat_mvtx;
0107   Float_t rec_matchrat_tpc;
0108 
0109   // set branch addresses (ignore cluster branches for now)
0110   tInTrue -> SetBranchAddress("event",         &tru_event);
0111   tInTrue -> SetBranchAddress("nphg4_part",    &tru_nphg4_part);
0112   tInTrue -> SetBranchAddress("centrality",    &tru_centrality);
0113   tInTrue -> SetBranchAddress("ntrackmatches", &tru_ntrackmatches);
0114   tInTrue -> SetBranchAddress("nphg4",         &tru_nphg4);
0115   tInTrue -> SetBranchAddress("nsvtx",         &tru_nsvtx);
0116   tInTrue -> SetBranchAddress("trackid",       &tru_trackid);
0117   tInTrue -> SetBranchAddress("is_G4track",    &tru_is_G4track);
0118   tInTrue -> SetBranchAddress("is_Svtrack",    &tru_is_Svtrack);
0119   tInTrue -> SetBranchAddress("is_matched",    &tru_is_matched);
0120   tInTrue -> SetBranchAddress("trkpt",         &tru_trkpt);
0121   tInTrue -> SetBranchAddress("trketa",        &tru_trketa);
0122   tInTrue -> SetBranchAddress("trkphi",        &tru_trkphi);
0123   tInTrue -> SetBranchAddress("nclus",         &tru_nclus);
0124   tInTrue -> SetBranchAddress("nclustpc",      &tru_nclustpc);
0125   tInTrue -> SetBranchAddress("nclusmvtx",     &tru_nclusmvtx);
0126   tInTrue -> SetBranchAddress("nclusintt",     &tru_nclusintt);
0127   tInTrue -> SetBranchAddress("matchrat",      &tru_matchrat);
0128   tInTrue -> SetBranchAddress("matchrat_intt", &tru_matchrat_intt);
0129   tInTrue -> SetBranchAddress("matchrat_mvtx", &tru_matchrat_mvtx);
0130   tInTrue -> SetBranchAddress("matchrat_tpc",  &tru_matchrat_tpc);
0131 
0132   tInReco -> SetBranchAddress("event",         &rec_event);
0133   tInReco -> SetBranchAddress("nphg4_part",    &rec_nphg4_part);
0134   tInReco -> SetBranchAddress("centrality",    &rec_centrality);
0135   tInReco -> SetBranchAddress("ntrackmatches", &rec_ntrackmatches);
0136   tInReco -> SetBranchAddress("nphg4",         &rec_nphg4);
0137   tInReco -> SetBranchAddress("nsvtx",         &rec_nsvtx);
0138   tInReco -> SetBranchAddress("trackid",       &rec_trackid);
0139   tInReco -> SetBranchAddress("is_G4track",    &rec_is_G4track);
0140   tInReco -> SetBranchAddress("is_Svtrack",    &rec_is_Svtrack);
0141   tInReco -> SetBranchAddress("is_matched",    &rec_is_matched);
0142   tInReco -> SetBranchAddress("trkpt",         &rec_trkpt);
0143   tInReco -> SetBranchAddress("trketa",        &rec_trketa);
0144   tInReco -> SetBranchAddress("trkphi",        &rec_trkphi);
0145   tInReco -> SetBranchAddress("nclus",         &rec_nclus);
0146   tInReco -> SetBranchAddress("nclustpc",      &rec_nclustpc);
0147   tInReco -> SetBranchAddress("nclusmvtx",     &rec_nclusmvtx);
0148   tInReco -> SetBranchAddress("nclusintt",     &rec_nclusintt);
0149   tInReco -> SetBranchAddress("matchrat",      &rec_matchrat);
0150   tInReco -> SetBranchAddress("matchrat_intt", &rec_matchrat_intt);
0151   tInReco -> SetBranchAddress("matchrat_mvtx", &rec_matchrat_mvtx);
0152   tInReco -> SetBranchAddress("matchrat_tpc",  &rec_matchrat_tpc);
0153   cout << "    Set input branches." << endl;
0154 
0155   // histogram indices
0156   enum Var {
0157     NTot,
0158     NIntt,
0159     NMvtx,
0160     NTpc,
0161     RTot,
0162     RIntt,
0163     RMvtx,
0164     RTpc,
0165     Phi,
0166     Eta,
0167     Pt,
0168     Frac,
0169     Qual,
0170     PtErr,
0171     EtaErr,
0172     PhiErr,
0173     PtRes,
0174     EtaRes,
0175     PhiRes
0176   };
0177   enum Type {
0178     Truth,
0179     Track,
0180     Weird,
0181     Normal
0182   };
0183   enum Comp {
0184     VsTruthPt,
0185     VsNumTpc
0186   };
0187 
0188   // output histogram base names
0189   const vector<vector<string>> vecHistBase = {
0190     {"hTruthNumTot",  "hTrackNumTot",  "hWeirdNumTot",  "hNormNumTot"},
0191     {"hTruthNumIntt", "hTrackNumIntt", "hWeirdNumIntt", "hNormNumIntt"},
0192     {"hTruthNumMvtx", "hTrackNumMvtx", "hWeirdNumMvtx", "hNormNumMvtx"},
0193     {"hTruthNumTpc",  "hTrackNumTpc",  "hWeirdNumTpc",  "hNormNumTpc"},
0194     {"hTruthRatTot",  "hTrackRatTot",  "hWeirdRatTot",  "hNormRatTot"},
0195     {"hTruthRatIntt", "hTrackRatIntt", "hWeirdRatIntt", "hNormRatIntt"},
0196     {"hTruthRatMvtx", "hTrackRatMvtx", "hWeirdRatMvtx", "hNormRatMvtx"},
0197     {"hTruthRatTpc",  "hTrackRatTpc",  "hWeirdRatTpc",  "hNormRatTpc"},
0198     {"hTruthPhi",     "hTrackPhi",     "hWeirdPhi",     "hNormPhi"},
0199     {"hTruthEta",     "hTrackEta",     "hWeirdEta",     "hNormEta"},
0200     {"hTruthPt",      "hTrackPt",      "hWeirdPt",      "hNormPt"},
0201     {"hTruthFrac",    "hTrackFrac",    "hWeirdFrac",    "hNormFrac"},
0202     {"hTruthQual",    "hTrackQual",    "hWeirdQual",    "hNormQual"},
0203     {"hTruthPtErr",   "hTrackPtErr",   "hWeirdPtErr",   "hNormPtErr"},
0204     {"hTruthEtaErr",  "hTrackEtaErr",  "hWeirdEtaErr",  "hNormEtaErr"},
0205     {"hTruthPhiErr",  "hTrackPhiErr",  "hWeirdPhiErr",  "hNormPhiErr"},
0206     {"hTruthPtRes",   "hTrackPtRes",   "hWeirdPtRes",   "hNormPtRes"},
0207     {"hTruthEtaRes",  "hTrackEtaRes",  "hWeirdEtaRes",  "hNormEtaRes"},
0208     {"hTruthPhiRes",  "hTrackPhiRes",  "hWeirdPhiRes",  "hNormPhiRes"}
0209   };
0210   const size_t nHistRows  = vecHistBase.size();
0211   const size_t nHistTypes = vecHistBase[0].size();
0212 
0213   // 2D histogram name modifiers
0214   const vector<string> vecVsHistModifiers = {
0215     "VsTruthPt",
0216     "VsNumTpc"
0217   };
0218   const size_t nVsHistMods = vecVsHistModifiers.size();
0219 
0220   // axis titles
0221   const string         sCount("counts");
0222   const vector<string> vecBaseAxisVars = {
0223     "N^{tot} = N_{hit}^{mvtx} + N_{hit}^{intt} + N_{clust}^{tpc}",
0224     "N_{hit}^{intt}",
0225     "N_{hit}^{mvtx}",
0226     "N_{clust}^{tpc}",
0227     "N_{reco}^{tot} / N_{true}^{tot}",
0228     "N_{reco}^{intt} / N_{true}^{intt}",
0229     "N_{reco}^{mvtx} / N_{true}^{mvtx}",
0230     "N_{reco}^{tpc} / N_{true}^{tpc}",
0231     "#varphi",
0232     "#eta",
0233     "p_{T} [GeV/c]",
0234     "p_{T}^{reco} / p_{T}^{true}"
0235   };
0236   const vector<string> vecVsAxisVars = {
0237     "p_{T}^{true} [GeV/c]",
0238     "N_{clust}^{tpc}"
0239   };
0240 
0241   // output histogram no. of bins
0242   const uint32_t nNumBins  = 101;
0243   const uint32_t nRatBins  = 120;
0244   const uint32_t nEtaBins  = 80;
0245   const uint32_t nPhiBins  = 360;
0246   const uint32_t nPtBins   = 202;
0247   const uint32_t nFracBins = 220;
0248   const uint32_t nQualBins = ;
0249   const uint32_t nResBins  = ;
0250 
0251   // output histogram bin ranges
0252   const pair<float, float> xNumBins  = {-0.5,  100.5};
0253   const pair<float, float> xRatBins  = {-0.5,  5.5};
0254   const pair<float, float> xEtaBins  = {-2.,   2.};
0255   const pair<float, float> xPhiBins  = {-3.15, 3.15};
0256   const pair<float, float> xPtBins   = {-1.,   100.};
0257   const pair<float, float> xFracBins = {-0.5,  10.5};
0258   const pair<float, float> xQualBins = {};
0259   const pair<float, float> xResBins  = {};
0260 
0261   // output histogram base binning definitions
0262   vector<tuple<uint32_t, pair<float, float>>> vecBaseHistBins = {
0263     make_tuple(nNumBins,  xNumBins),
0264     make_tuple(nNumBins,  xNumBins),
0265     make_tuple(nNumBins,  xNumBins),
0266     make_tuple(nNumBins,  xNumBins),
0267     make_tuple(nRatBins,  xRatBins),
0268     make_tuple(nRatBins,  xRatBins),
0269     make_tuple(nRatBins,  xRatBins),
0270     make_tuple(nRatBins,  xRatBins),
0271     make_tuple(nPhiBins,  xPhiBins),
0272     make_tuple(nEtaBins,  xEtaBins),
0273     make_tuple(nPtBins,   xPtBins),
0274     make_tuple(nFracBins, xFracBins),
0275     make_tuple(nResBins,  xResBins)
0276   };
0277 
0278   // output 2D histogram x-axis binning
0279   vector<tuple<uint32_t, pair<float, float>>> vecVsHistBins = {
0280     make_tuple(nPtBins,  xPtBins),
0281     make_tuple(nNumBins, xNumBins)
0282   };
0283 
0284   // make 1D histograms
0285   vector<vector<TH1D*>> vecHist1D(nHistRows);
0286   for (size_t iHistRow = 0; iHistRow < nHistRows; iHistRow++) {
0287     for (const string sHistBase : vecHistBase[iHistRow]) {
0288       vecHist1D[iHistRow].push_back(
0289         new TH1D(
0290           sHistBase.data(),
0291           "",
0292           get<0>(vecBaseHistBins[iHistRow]),
0293           get<1>(vecBaseHistBins[iHistRow]).first,
0294           get<1>(vecBaseHistBins[iHistRow]).second
0295         )
0296       );
0297     }  // end type loop
0298   }  // end row loop
0299 
0300   // make 2D histograms
0301   vector<vector<vector<TH2D*>>> vecHist2D(nHistRows, vector<vector<TH2D*>>(nHistTypes));
0302   for (size_t iHistRow = 0; iHistRow < nHistRows; iHistRow++) {
0303     for (size_t iHistType = 0; iHistType < nHistTypes; iHistType++) {
0304       for (size_t iVsHistMod = 0; iVsHistMod < nVsHistMods; iVsHistMod++) {
0305         const string sHistName2D = vecHistBase[iHistRow][iHistType] + vecVsHistModifiers[iVsHistMod];
0306         vecHist2D[iHistRow][iHistType].push_back(
0307           new TH2D(
0308             sHistName2D.data(),
0309             "",
0310             get<0>(vecVsHistBins[iVsHistMod]),
0311             get<1>(vecVsHistBins[iVsHistMod]).first,
0312             get<1>(vecVsHistBins[iVsHistMod]).second,
0313             get<0>(vecBaseHistBins[iHistRow]),
0314             get<1>(vecBaseHistBins[iHistRow]).first,
0315             get<1>(vecBaseHistBins[iHistRow]).second
0316           )
0317         );
0318       }  // end modifier loop
0319     }  // end type loop
0320   }  // end row loop
0321 
0322   // set errors
0323   for (auto histRow1D : vecHist1D) {
0324     for (auto hist1D : histRow1D) {
0325       hist1D -> Sumw2();
0326     }
0327   }
0328 
0329   for (auto histRow2D : vecHist2D) {
0330     for (auto histType2D : histRow2D) {
0331       for (auto hist2D: histType2D) {
0332         hist2D -> Sumw2();
0333       }
0334     }
0335   }
0336 
0337   // grab no. of entries
0338   const int64_t nTrueEntries = tInTrue -> GetEntries();
0339   const int64_t nRecoEntries = tInReco -> GetEntries(); 
0340   cout << "    Beginning truth particle loop: " << nTrueEntries << " to process" << endl;
0341 
0342   // loop over truth particles
0343   int64_t nTrueBytes = 0;
0344   for (int64_t iTrueEntry = 0; iTrueEntry < nTrueEntries; iTrueEntry++) {
0345 
0346     // grab truth particle entry
0347     const int64_t trueBytes = tInTrue -> GetEntry(iTrueEntry);
0348     if (trueBytes < 0) {
0349       cerr << "PANIC: issue with entry " << iTrueEntry << "! Aborting loop!\n" << endl;
0350       break;
0351     } else {
0352       nTrueBytes += trueBytes;
0353     }
0354 
0355     const int64_t iTrueProg = iTrueEntry + 1;
0356     if (iTrueProg == nTrueEntries) {
0357       cout << "      Processing entry " << iTrueProg << "/" << nTrueEntries << "..." << endl;
0358     } else {
0359       cout << "      Processing entry " << iTrueProg << "/" << nTrueEntries << "...\r" << flush;
0360     }
0361 
0362     // select truth particles
0363     if (!tru_is_G4track) continue;
0364 
0365     // fill truth 1D histograms
0366     vecHist1D[Var::NTot][Type::Truth]  -> Fill(tru_nclus);
0367     vecHist1D[Var::NIntt][Type::Truth] -> Fill(tru_nclusintt);
0368     vecHist1D[Var::NMvtx][Type::Truth] -> Fill(tru_nclusmvtx);
0369     vecHist1D[Var::NTpc][Type::Truth]  -> Fill(tru_nclustpc);
0370     vecHist1D[Var::RTot][Type::Truth]  -> Fill(1.);
0371     vecHist1D[Var::RIntt][Type::Truth] -> Fill(1.);
0372     vecHist1D[Var::RMvtx][Type::Truth] -> Fill(1.);
0373     vecHist1D[Var::RTpc][Type::Truth]  -> Fill(1.);
0374     vecHist1D[Var::Phi][Type::Truth]   -> Fill(tru_trkphi);
0375     vecHist1D[Var::Eta][Type::Truth]   -> Fill(tru_trketa);
0376     vecHist1D[Var::Pt][Type::Truth]    -> Fill(tru_trkpt);
0377     vecHist1D[Var::Frac][Type::Truth]  -> Fill(1.);
0378 
0379     // fill truth 2D histograms
0380     vecHist2D[Var::NTot][Type::Truth][Comp::VsTruthPt]  -> Fill(tru_trkpt, tru_nclus);
0381     vecHist2D[Var::NIntt][Type::Truth][Comp::VsTruthPt] -> Fill(tru_trkpt, tru_nclusintt);
0382     vecHist2D[Var::NMvtx][Type::Truth][Comp::VsTruthPt] -> Fill(tru_trkpt, tru_nclusmvtx);
0383     vecHist2D[Var::NTpc][Type::Truth][Comp::VsTruthPt]  -> Fill(tru_trkpt, tru_nclustpc);
0384     vecHist2D[Var::RTot][Type::Truth][Comp::VsTruthPt]  -> Fill(tru_trkpt, 1.);
0385     vecHist2D[Var::RIntt][Type::Truth][Comp::VsTruthPt] -> Fill(tru_trkpt, 1.);
0386     vecHist2D[Var::RMvtx][Type::Truth][Comp::VsTruthPt] -> Fill(tru_trkpt, 1.);
0387     vecHist2D[Var::RTpc][Type::Truth][Comp::VsTruthPt]  -> Fill(tru_trkpt, 1.);
0388     vecHist2D[Var::Phi][Type::Truth][Comp::VsTruthPt]   -> Fill(tru_trkpt, tru_trkphi);
0389     vecHist2D[Var::Eta][Type::Truth][Comp::VsTruthPt]   -> Fill(tru_trkpt, tru_trketa);
0390     vecHist2D[Var::Pt][Type::Truth][Comp::VsTruthPt]    -> Fill(tru_trkpt, tru_trkpt);
0391     vecHist2D[Var::Frac][Type::Truth][Comp::VsTruthPt]  -> Fill(tru_trkpt, 1.);
0392 
0393     vecHist2D[Var::NTot][Type::Truth][Comp::VsNumTpc]  -> Fill(tru_nclustpc, tru_nclus);
0394     vecHist2D[Var::NIntt][Type::Truth][Comp::VsNumTpc] -> Fill(tru_nclustpc, tru_nclusintt);
0395     vecHist2D[Var::NMvtx][Type::Truth][Comp::VsNumTpc] -> Fill(tru_nclustpc, tru_nclusmvtx);
0396     vecHist2D[Var::NTpc][Type::Truth][Comp::VsNumTpc]  -> Fill(tru_nclustpc, tru_nclustpc);
0397     vecHist2D[Var::RTot][Type::Truth][Comp::VsNumTpc]  -> Fill(tru_nclustpc, 1.);
0398     vecHist2D[Var::RIntt][Type::Truth][Comp::VsNumTpc] -> Fill(tru_nclustpc, 1.);
0399     vecHist2D[Var::RMvtx][Type::Truth][Comp::VsNumTpc] -> Fill(tru_nclustpc, 1.);
0400     vecHist2D[Var::RTpc][Type::Truth][Comp::VsNumTpc]  -> Fill(tru_nclustpc, 1.);
0401     vecHist2D[Var::Phi][Type::Truth][Comp::VsNumTpc]   -> Fill(tru_nclustpc, tru_trkphi);
0402     vecHist2D[Var::Eta][Type::Truth][Comp::VsNumTpc]   -> Fill(tru_nclustpc, tru_trketa);
0403     vecHist2D[Var::Pt][Type::Truth][Comp::VsNumTpc]    -> Fill(tru_nclustpc, tru_trkpt);
0404     vecHist2D[Var::Frac][Type::Truth][Comp::VsNumTpc]  -> Fill(tru_nclustpc, 1.);
0405   }  // end truth particle loop
0406 
0407   // announce next entry loop
0408   cout << "    Finished truth particle loop.\n"
0409        << "    Beginning reconstructed track loop: " << nRecoEntries << " to process"
0410        << endl;
0411 
0412   // loop over reco tracks
0413   // TODO identify matched truth particles
0414   int64_t nRecoBytes = 0;
0415   for (int64_t iRecoEntry = 0; iRecoEntry < nRecoEntries; iRecoEntry++) {
0416 
0417     // grab reco track entry
0418     const int64_t recoBytes = tInReco -> GetEntry(iRecoEntry);
0419     if (recoBytes < 0) {
0420       cerr << "PANIC: issue with entry " << iRecoEntry << "! Aborting loop!\n" << endl;
0421       break;
0422     } else {
0423       nRecoBytes += recoBytes;
0424     }
0425 
0426     const int64_t iRecoProg = iRecoEntry + 1;
0427     if (iRecoProg == nRecoEntries) {
0428       cout << "      Processing entry " << iRecoProg << "/" << nRecoEntries << "..." << endl;
0429     } else {
0430       cout << "      Processing entry " << iRecoProg << "/" << nRecoEntries << "...\r" << flush;
0431     }
0432 
0433     // select only tracks matched to truth particle
0434     if (!rec_is_matched || !rec_is_Svtrack) continue;
0435 
0436     // fill all matched reco 1D histograms
0437     // FIXME actually calculate pt fraction
0438     vecHist1D[Var::NTot][Type::Track]  -> Fill(rec_nclus);
0439     vecHist1D[Var::NIntt][Type::Track] -> Fill(rec_nclusintt);
0440     vecHist1D[Var::NMvtx][Type::Track] -> Fill(rec_nclusmvtx);
0441     vecHist1D[Var::NTpc][Type::Track]  -> Fill(rec_nclustpc);
0442     vecHist1D[Var::RTot][Type::Track]  -> Fill(rec_matchrat);
0443     vecHist1D[Var::RIntt][Type::Track] -> Fill(rec_matchrat_intt);
0444     vecHist1D[Var::RMvtx][Type::Track] -> Fill(rec_matchrat_mvtx);
0445     vecHist1D[Var::RTpc][Type::Track]  -> Fill(rec_matchrat_tpc);
0446     vecHist1D[Var::Phi][Type::Track]   -> Fill(rec_trkphi);
0447     vecHist1D[Var::Eta][Type::Track]   -> Fill(rec_trketa);
0448     vecHist1D[Var::Pt][Type::Track]    -> Fill(rec_trkpt);
0449     vecHist1D[Var::Frac][Type::Track]  -> Fill(1.);
0450 
0451     // fill all matched reco 2D histograms
0452     // FIXME use actual truth pt & ntpc of matched particle
0453     vecHist2D[Var::NTot][Type::Track][Comp::VsTruthPt]  -> Fill(rec_trkpt, rec_nclus);
0454     vecHist2D[Var::NIntt][Type::Track][Comp::VsTruthPt] -> Fill(rec_trkpt, rec_nclusintt);
0455     vecHist2D[Var::NMvtx][Type::Track][Comp::VsTruthPt] -> Fill(rec_trkpt, rec_nclusmvtx);
0456     vecHist2D[Var::NTpc][Type::Track][Comp::VsTruthPt]  -> Fill(rec_trkpt, rec_nclustpc);
0457     vecHist2D[Var::RTot][Type::Track][Comp::VsTruthPt]  -> Fill(rec_trkpt, rec_matchrat);
0458     vecHist2D[Var::RIntt][Type::Track][Comp::VsTruthPt] -> Fill(rec_trkpt, rec_matchrat_intt);
0459     vecHist2D[Var::RMvtx][Type::Track][Comp::VsTruthPt] -> Fill(rec_trkpt, rec_matchrat_mvtx);
0460     vecHist2D[Var::RTpc][Type::Track][Comp::VsTruthPt]  -> Fill(rec_trkpt, rec_matchrat_tpc);
0461     vecHist2D[Var::Phi][Type::Track][Comp::VsTruthPt]   -> Fill(rec_trkpt, rec_trkphi);
0462     vecHist2D[Var::Eta][Type::Track][Comp::VsTruthPt]   -> Fill(rec_trkpt, rec_trketa);
0463     vecHist2D[Var::Pt][Type::Track][Comp::VsTruthPt]    -> Fill(rec_trkpt, rec_trkpt);
0464     vecHist2D[Var::Frac][Type::Track][Comp::VsTruthPt]  -> Fill(rec_trkpt, 1.);
0465 
0466     vecHist2D[Var::NTot][Type::Track][Comp::VsNumTpc]  -> Fill(rec_nclustpc, rec_nclus);
0467     vecHist2D[Var::NIntt][Type::Track][Comp::VsNumTpc] -> Fill(rec_nclustpc, rec_nclusintt);
0468     vecHist2D[Var::NMvtx][Type::Track][Comp::VsNumTpc] -> Fill(rec_nclustpc, rec_nclusmvtx);
0469     vecHist2D[Var::NTpc][Type::Track][Comp::VsNumTpc]  -> Fill(rec_nclustpc, rec_nclustpc);
0470     vecHist2D[Var::RTot][Type::Track][Comp::VsNumTpc]  -> Fill(rec_nclustpc, rec_matchrat);
0471     vecHist2D[Var::RIntt][Type::Track][Comp::VsNumTpc] -> Fill(rec_nclustpc, rec_matchrat_intt);
0472     vecHist2D[Var::RMvtx][Type::Track][Comp::VsNumTpc] -> Fill(rec_nclustpc, rec_matchrat_mvtx);
0473     vecHist2D[Var::RTpc][Type::Track][Comp::VsNumTpc]  -> Fill(rec_nclustpc, rec_matchrat_tpc);
0474     vecHist2D[Var::Phi][Type::Track][Comp::VsNumTpc]   -> Fill(rec_nclustpc, rec_trkphi);
0475     vecHist2D[Var::Eta][Type::Track][Comp::VsNumTpc]   -> Fill(rec_nclustpc, rec_trketa);
0476     vecHist2D[Var::Pt][Type::Track][Comp::VsNumTpc]    -> Fill(rec_nclustpc, rec_trkpt);
0477     vecHist2D[Var::Frac][Type::Track][Comp::VsNumTpc]  -> Fill(rec_nclustpc, 1.);
0478 
0479     // fill weird and normal matched reco 1D histograms
0480     // FIXME actually check if is a weird track
0481     const bool isWeirdTrack = true;
0482     if (isWeirdTrack) {
0483       vecHist1D[Var::NTot][Type::Weird]  -> Fill(rec_nclus);
0484       vecHist1D[Var::NIntt][Type::Weird] -> Fill(rec_nclusintt);
0485       vecHist1D[Var::NMvtx][Type::Weird] -> Fill(rec_nclusmvtx);
0486       vecHist1D[Var::NTpc][Type::Weird]  -> Fill(rec_nclustpc);
0487       vecHist1D[Var::RTot][Type::Weird]  -> Fill(rec_matchrat);
0488       vecHist1D[Var::RIntt][Type::Weird] -> Fill(rec_matchrat_intt);
0489       vecHist1D[Var::RMvtx][Type::Weird] -> Fill(rec_matchrat_mvtx);
0490       vecHist1D[Var::RTpc][Type::Weird]  -> Fill(rec_matchrat_tpc);
0491       vecHist1D[Var::Phi][Type::Weird]   -> Fill(rec_trkphi);
0492       vecHist1D[Var::Eta][Type::Weird]   -> Fill(rec_trketa);
0493       vecHist1D[Var::Pt][Type::Weird]    -> Fill(rec_trkpt);
0494       vecHist1D[Var::Frac][Type::Weird]  -> Fill(1.);
0495 
0496       vecHist2D[Var::NTot][Type::Weird][Comp::VsTruthPt]  -> Fill(rec_trkpt, rec_nclus);
0497       vecHist2D[Var::NIntt][Type::Weird][Comp::VsTruthPt] -> Fill(rec_trkpt, rec_nclusintt);
0498       vecHist2D[Var::NMvtx][Type::Weird][Comp::VsTruthPt] -> Fill(rec_trkpt, rec_nclusmvtx);
0499       vecHist2D[Var::NTpc][Type::Weird][Comp::VsTruthPt]  -> Fill(rec_trkpt, rec_nclustpc);
0500       vecHist2D[Var::RTot][Type::Weird][Comp::VsTruthPt]  -> Fill(rec_trkpt, rec_matchrat);
0501       vecHist2D[Var::RIntt][Type::Weird][Comp::VsTruthPt] -> Fill(rec_trkpt, rec_matchrat_intt);
0502       vecHist2D[Var::RMvtx][Type::Weird][Comp::VsTruthPt] -> Fill(rec_trkpt, rec_matchrat_mvtx);
0503       vecHist2D[Var::RTpc][Type::Weird][Comp::VsTruthPt]  -> Fill(rec_trkpt, rec_matchrat_tpc);
0504       vecHist2D[Var::Phi][Type::Weird][Comp::VsTruthPt]   -> Fill(rec_trkpt, rec_trkphi);
0505       vecHist2D[Var::Eta][Type::Weird][Comp::VsTruthPt]   -> Fill(rec_trkpt, rec_trketa);
0506       vecHist2D[Var::Pt][Type::Weird][Comp::VsTruthPt]    -> Fill(rec_trkpt, rec_trkpt);
0507       vecHist2D[Var::Frac][Type::Weird][Comp::VsTruthPt]  -> Fill(rec_trkpt, 1.);
0508 
0509       vecHist2D[Var::NTot][Type::Weird][Comp::VsNumTpc]  -> Fill(rec_nclustpc, rec_nclus);
0510       vecHist2D[Var::NIntt][Type::Weird][Comp::VsNumTpc] -> Fill(rec_nclustpc, rec_nclusintt);
0511       vecHist2D[Var::NMvtx][Type::Weird][Comp::VsNumTpc] -> Fill(rec_nclustpc, rec_nclusmvtx);
0512       vecHist2D[Var::NTpc][Type::Weird][Comp::VsNumTpc]  -> Fill(rec_nclustpc, rec_nclustpc);
0513       vecHist2D[Var::RTot][Type::Weird][Comp::VsNumTpc]  -> Fill(rec_nclustpc, rec_matchrat);
0514       vecHist2D[Var::RIntt][Type::Weird][Comp::VsNumTpc] -> Fill(rec_nclustpc, rec_matchrat_intt);
0515       vecHist2D[Var::RMvtx][Type::Weird][Comp::VsNumTpc] -> Fill(rec_nclustpc, rec_matchrat_mvtx);
0516       vecHist2D[Var::RTpc][Type::Weird][Comp::VsNumTpc]  -> Fill(rec_nclustpc, rec_matchrat_tpc);
0517       vecHist2D[Var::Phi][Type::Weird][Comp::VsNumTpc]   -> Fill(rec_nclustpc, rec_trkphi);
0518       vecHist2D[Var::Eta][Type::Weird][Comp::VsNumTpc]   -> Fill(rec_nclustpc, rec_trketa);
0519       vecHist2D[Var::Pt][Type::Weird][Comp::VsNumTpc]    -> Fill(rec_nclustpc, rec_trkpt);
0520       vecHist2D[Var::Frac][Type::Weird][Comp::VsNumTpc]  -> Fill(rec_nclustpc, 1.);
0521     } else {
0522       vecHist1D[Var::NTot][Type::Normal]  -> Fill(rec_nclus);
0523       vecHist1D[Var::NIntt][Type::Normal] -> Fill(rec_nclusintt);
0524       vecHist1D[Var::NMvtx][Type::Normal] -> Fill(rec_nclusmvtx);
0525       vecHist1D[Var::NTpc][Type::Normal]  -> Fill(rec_nclustpc);
0526       vecHist1D[Var::RTot][Type::Normal]  -> Fill(rec_matchrat);
0527       vecHist1D[Var::RIntt][Type::Normal] -> Fill(rec_matchrat_intt);
0528       vecHist1D[Var::RMvtx][Type::Normal] -> Fill(rec_matchrat_mvtx);
0529       vecHist1D[Var::RTpc][Type::Normal]  -> Fill(rec_matchrat_tpc);
0530       vecHist1D[Var::Phi][Type::Normal]   -> Fill(rec_trkphi);
0531       vecHist1D[Var::Eta][Type::Normal]   -> Fill(rec_trketa);
0532       vecHist1D[Var::Pt][Type::Normal]    -> Fill(rec_trkpt);
0533       vecHist1D[Var::Frac][Type::Normal]  -> Fill(1.);
0534 
0535       vecHist2D[Var::NTot][Type::Normal][Comp::VsTruthPt]  -> Fill(rec_trkpt, rec_nclus);
0536       vecHist2D[Var::NIntt][Type::Normal][Comp::VsTruthPt] -> Fill(rec_trkpt, rec_nclusintt);
0537       vecHist2D[Var::NMvtx][Type::Normal][Comp::VsTruthPt] -> Fill(rec_trkpt, rec_nclusmvtx);
0538       vecHist2D[Var::NTpc][Type::Normal][Comp::VsTruthPt]  -> Fill(rec_trkpt, rec_nclustpc);
0539       vecHist2D[Var::RTot][Type::Normal][Comp::VsTruthPt]  -> Fill(rec_trkpt, rec_matchrat);
0540       vecHist2D[Var::RIntt][Type::Normal][Comp::VsTruthPt] -> Fill(rec_trkpt, rec_matchrat_intt);
0541       vecHist2D[Var::RMvtx][Type::Normal][Comp::VsTruthPt] -> Fill(rec_trkpt, rec_matchrat_mvtx);
0542       vecHist2D[Var::RTpc][Type::Normal][Comp::VsTruthPt]  -> Fill(rec_trkpt, rec_matchrat_tpc);
0543       vecHist2D[Var::Phi][Type::Normal][Comp::VsTruthPt]   -> Fill(rec_trkpt, rec_trkphi);
0544       vecHist2D[Var::Eta][Type::Normal][Comp::VsTruthPt]   -> Fill(rec_trkpt, rec_trketa);
0545       vecHist2D[Var::Pt][Type::Normal][Comp::VsTruthPt]    -> Fill(rec_trkpt, rec_trkpt);
0546       vecHist2D[Var::Frac][Type::Normal][Comp::VsTruthPt]  -> Fill(rec_trkpt, 1.);
0547 
0548       vecHist2D[Var::NTot][Type::Normal][Comp::VsNumTpc]  -> Fill(rec_nclustpc, rec_nclus);
0549       vecHist2D[Var::NIntt][Type::Normal][Comp::VsNumTpc] -> Fill(rec_nclustpc, rec_nclusintt);
0550       vecHist2D[Var::NMvtx][Type::Normal][Comp::VsNumTpc] -> Fill(rec_nclustpc, rec_nclusmvtx);
0551       vecHist2D[Var::NTpc][Type::Normal][Comp::VsNumTpc]  -> Fill(rec_nclustpc, rec_nclustpc);
0552       vecHist2D[Var::RTot][Type::Normal][Comp::VsNumTpc]  -> Fill(rec_nclustpc, rec_matchrat);
0553       vecHist2D[Var::RIntt][Type::Normal][Comp::VsNumTpc] -> Fill(rec_nclustpc, rec_matchrat_intt);
0554       vecHist2D[Var::RMvtx][Type::Normal][Comp::VsNumTpc] -> Fill(rec_nclustpc, rec_matchrat_mvtx);
0555       vecHist2D[Var::RTpc][Type::Normal][Comp::VsNumTpc]  -> Fill(rec_nclustpc, rec_matchrat_tpc);
0556       vecHist2D[Var::Phi][Type::Normal][Comp::VsNumTpc]   -> Fill(rec_nclustpc, rec_trkphi);
0557       vecHist2D[Var::Eta][Type::Normal][Comp::VsNumTpc]   -> Fill(rec_nclustpc, rec_trketa);
0558       vecHist2D[Var::Pt][Type::Normal][Comp::VsNumTpc]    -> Fill(rec_nclustpc, rec_trkpt);
0559       vecHist2D[Var::Frac][Type::Normal][Comp::VsNumTpc]  -> Fill(rec_nclustpc, 1.);
0560     }
0561   }  // end reco track loop
0562   cout << "    Finished reconstructed track loop."  << endl;
0563 
0564   // set axis titles
0565   size_t iVar = 0;
0566   for (auto histRow1D : vecHist1D) {
0567     for (auto hist1D : histRow1D) {
0568       hist1D -> GetXaxis() -> SetTitle(vecBaseAxisVars.at(iVar).data());
0569       hist1D -> GetYaxis() -> SetTitle(sCount.data());
0570     }
0571     ++iVar;
0572   }
0573   iVar = 0;
0574 
0575   size_t iComp = 0;
0576   for (auto histRow2D : vecHist2D) {
0577     for (auto histType2D : histRow2D) {
0578       iComp = 0;
0579       for (auto hist2D : histType2D) {
0580         hist2D -> GetXaxis() -> SetTitle(vecVsAxisVars.at(iComp).data());
0581         hist2D -> GetYaxis() -> SetTitle(vecBaseAxisVars.at(iVar).data());
0582         hist2D -> GetZaxis() -> SetTitle(sCount.data());
0583         ++iComp;
0584       }
0585     }
0586     ++iVar;
0587   }
0588   cout << "    Set axis titles." << endl;
0589 
0590   // save histograms
0591   fOutput -> cd();
0592   for (auto histRow1D : vecHist1D) {
0593     for (auto hist1D : histRow1D) {
0594       hist1D -> Write();
0595     }
0596   }
0597   for (auto histRow2D : vecHist2D) {
0598     for (auto histType2D : histRow2D) {
0599       for (auto hist2D: histType2D) {
0600         hist2D -> Write();
0601       }
0602     }
0603   }
0604   cout << "    Saved histograms." << endl;
0605 
0606   // close files
0607   fOutput -> cd();
0608   fOutput -> Close();
0609   fInTrue -> cd();
0610   fInTrue -> Close();
0611   fInReco -> cd();
0612   fInReco -> Close();
0613 
0614   // exit
0615   cout << "  Finished new matcher plot script!\n" << endl;
0616   return;
0617 
0618 }
0619 
0620 // end ------------------------------------------------------------------------
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.3.7 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!