SCorrelatorQAMaker/src/SReadLambdaJetTree.cc

0001 // ----------------------------------------------------------------------------
0002 // 'ReadLambdaJetTree.cc'
0003 // Derek Anderson
0004 // 02.19.2024
0005 //
0006 // SCorrelatorQAMaker plugin to read lambda-tagged jet
0007 // trees and draw plots and calculate the total no.
0008 // of lambdas.
0009 // ----------------------------------------------------------------------------
0010
0011 #define SCORRELATORQAMAKER_SREADLAMBDAJETTREE_CC
0012
0013 // plugin definition
0014 #include "SReadLambdaJetTree.h"
0015
0016 // make common namespaces implicit
0017 using namespace std;
0018
0019
0020
0021 namespace SColdQcdCorrelatorAnalysis {
0022
0023   // SReadLambdaJetTree public methods ----------------------------------------
0024
0025   void SReadLambdaJetTree::Init() {
0026
0027     // announce start
0028     cout << "\n  Starting lambda jet tree reader!" << endl;
0029
0030     // initialize output
0031     InitOutput();
0032
0033     // run internal routines
0034     InitInput();
0035     InitTree();
0036     InitHists();
0037     return;
0038
0039   }  // end 'Init()'
0040
0041
0042
0043   void SReadLambdaJetTree::Analyze() {
0044
0045     DoAnalysis();
0046     return;
0047
0048   }  // end 'Analyze()'
0049
0050
0051
0052   void SReadLambdaJetTree::End() {
0053
0054     // run internal routines
0055     SetHistogramStyles();
0056     SaveOutput();
0057     CloseInput();
0058
0059     // close output
0060     CloseOutput();
0061
0062     // announce end and exit
0063     cout << "  Finished lambda jet tree reader!\n" << endl;
0064     return;
0065
0066   }  // end 'End()'
0067
0068
0069
0070   // SReadLambdaJetTree internal methods --------------------------------------
0071
0072   void SReadLambdaJetTree::InitInput() {
0073
0074     // open input files
0075     m_fInput = new TFile(m_config.inFileName.data(), "read");
0076     if (!m_fInput) {
0077       cerr << "PANIC: couldn't open input file:\n"
0078            << "       \"" << m_config.inFileName << "\"\n"
0079            << endl;
0080       assert(m_fInput);
0081     }
0082     cout << "    Opened input file:\n"
0083          << "      \"" << m_config.inFileName << "\""
0084          << endl;
0085
0086     // grab input tree
0087     m_tInput = (TTree*) m_fInput -> Get(m_config.inTreeName.data());
0088     if (!m_tInput) {
0089       cerr << "PANIC: couldn't grab input tree \"" << m_config.inTreeName << "\"!\n" << endl;
0090       assert(m_tInput);
0091     }
0092     cout << "    Grabbed input tree \"" << m_config.inTreeName << "\"." << endl;
0093
0094     // exit routine
0095     return;
0096
0097   }  // end 'InitInput()'
0098
0099
0100
0101   void SReadLambdaJetTree::InitTree() {
0102
0103     // set branch addresses
0104     m_tInput -> SetMakeClass(1);
0105     m_tInput -> SetBranchAddress("EvtNJets",       &m_evtNJets,       &m_brEvtNJets);
0106     m_tInput -> SetBranchAddress("EvtNLambdas",    &m_evtNLambdas,    &m_brEvtNLambdas);
0107     m_tInput -> SetBranchAddress("EvtNTaggedJets", &m_evtNTaggedJets, &m_brEvtNTaggedJets);
0108     m_tInput -> SetBranchAddress("EvtNChrgPars",   &m_evtNChrgPars,   &m_brEvtNChrgPars);
0109     m_tInput -> SetBranchAddress("EvtNNeuPars",    &m_evtNNeuPars,    &m_brEvtNNeuPars);
0110     m_tInput -> SetBranchAddress("EvtSumEPar",     &m_evtSumEPar,     &m_brEvtSumEPar);
0111     m_tInput -> SetBranchAddress("EvtVtxX",        &m_evtVtxX,        &m_brEvtVtxX);
0112     m_tInput -> SetBranchAddress("EvtVtxY",        &m_evtVtxY,        &m_brEvtVtxY);
0113     m_tInput -> SetBranchAddress("EvtVtxZ",        &m_evtVtxZ,        &m_brEvtVtxZ);
0114     m_tInput -> SetBranchAddress("PartonA_ID",     &m_partonA_ID,     &m_brPartonA_ID);
0115     m_tInput -> SetBranchAddress("PartonB_ID",     &m_partonB_ID,     &m_brPartonB_ID);
0116     m_tInput -> SetBranchAddress("PartonA_Px",     &m_partonA_Px,     &m_brPartonA_Px);
0117     m_tInput -> SetBranchAddress("PartonA_Py",     &m_partonA_Py,     &m_brPartonA_Py);
0118     m_tInput -> SetBranchAddress("PartonA_Pz",     &m_partonA_Pz,     &m_brPartonA_Pz);
0119     m_tInput -> SetBranchAddress("PartonA_E",      &m_partonA_E,      &m_brPartonA_E);
0120     m_tInput -> SetBranchAddress("PartonB_Px",     &m_partonB_Px,     &m_brPartonB_Px);
0121     m_tInput -> SetBranchAddress("PartonB_Py",     &m_partonB_Py,     &m_brPartonB_Py);
0122     m_tInput -> SetBranchAddress("PartonB_Pz",     &m_partonB_Pz,     &m_brPartonB_Pz);
0123     m_tInput -> SetBranchAddress("PartonB_E",      &m_partonB_E,      &m_brPartonB_E);
0124     m_tInput -> SetBranchAddress("LambdaID",       &m_lambdaID,       &m_brLambdaID);
0125     m_tInput -> SetBranchAddress("LambdaPID",      &m_lambdaPID,      &m_brLambdaPID);
0126     m_tInput -> SetBranchAddress("LambdaJetID",    &m_lambdaJetID,    &m_brLambdaJetID);
0127     m_tInput -> SetBranchAddress("LambdaEmbedID",  &m_lambdaEmbedID,  &m_brLambdaEmbedID);
0128     m_tInput -> SetBranchAddress("LambdaZ",        &m_lambdaZ,        &m_brLambdaZ);
0129     m_tInput -> SetBranchAddress("LambdaDr",       &m_lambdaDr,       &m_brLambdaDr);
0130     m_tInput -> SetBranchAddress("LambdaEnergy",   &m_lambdaEnergy,   &m_brLambdaEnergy);
0131     m_tInput -> SetBranchAddress("LambdaPt",       &m_lambdaPt,       &m_brLambdaPt);
0132     m_tInput -> SetBranchAddress("LambdaEta",      &m_lambdaEta,      &m_brLambdaEta);
0133     m_tInput -> SetBranchAddress("LambdaPhi",      &m_lambdaPhi,      &m_brLambdaPhi);
0134     m_tInput -> SetBranchAddress("JetHasLambda",   &m_jetHasLambda,   &m_brJetHasLambda);
0135     m_tInput -> SetBranchAddress("JetNCst",        &m_jetNCst,        &m_brJetNCst);
0136     m_tInput -> SetBranchAddress("JetID",          &m_jetID,          &m_brJetID);
0137     m_tInput -> SetBranchAddress("JetE",           &m_jetE,           &m_brJetE);
0138     m_tInput -> SetBranchAddress("JetPt",          &m_jetPt,          &m_brJetPt);
0139     m_tInput -> SetBranchAddress("JetEta",         &m_jetEta,         &m_brJetEta);
0140     m_tInput -> SetBranchAddress("JetPhi",         &m_jetPhi,         &m_brJetPhi);
0141     m_tInput -> SetBranchAddress("CstID",          &m_cstID,          &m_brCstID);
0142     m_tInput -> SetBranchAddress("CstPID",         &m_cstPID,         &m_brCstPID);
0143     m_tInput -> SetBranchAddress("CstJetID",       &m_cstJetID,       &m_brCstJetID);
0144     m_tInput -> SetBranchAddress("CstEmbedID",     &m_cstEmbedID,     &m_brCstEmbedID);
0145     m_tInput -> SetBranchAddress("CstZ",           &m_cstZ,           &m_brCstZ);
0146     m_tInput -> SetBranchAddress("CstDr",          &m_cstDr,          &m_brCstDr);
0147     m_tInput -> SetBranchAddress("CstEnergy",      &m_cstEnergy,      &m_brCstEnergy);
0148     m_tInput -> SetBranchAddress("CstPt",          &m_cstPt,          &m_brCstPt);
0149     m_tInput -> SetBranchAddress("CstEta",         &m_cstEta,         &m_brCstEta);
0150     m_tInput -> SetBranchAddress("CstPhi",         &m_cstPhi,         &m_brCstPhi);
0151     cout << "    Initialized input tree." << endl;
0152
0153     // exit routine
0154     return;
0155
0156   }  // end 'InitTree()'
0157
0158
0159
0160   void SReadLambdaJetTree::InitHists() {
0161
0162     // make sure sumw2 is on
0163     TH1::SetDefaultSumw2(true);
0164     TH2::SetDefaultSumw2(true);
0165
0166     // create event histograms
0167     for (auto evtNameAndTitle : m_hist.vecEvtNameAndTitles) {
0168       vecHistEvt.push_back(
0169         new TH1D(
0170           evtNameAndTitle.first.data(),
0171           evtNameAndTitle.second.data(),
0172           m_hist.nNumBins,
0173           m_hist.rNumBins.first,
0174           m_hist.rNumBins.second
0175         )
0176       );
0177     }
0178
0179     // create jet/lambda histograms
0180     vecHist1D.resize( m_hist.vecTypeNames.size() );
0181     vecHist2D.resize( m_hist.vecTypeNames.size() );
0182     for (size_t iType = 0; iType < m_hist.vecTypeNames.size(); iType++) {
0183
0184       // loop over base variables
0185       vecHist1D[iType].resize( m_hist.vecBaseNames.size() );
0186       vecHist2D[iType].resize( m_hist.vecBaseNames.size() );
0187       for (size_t iVar = 0; iVar < m_hist.vecBaseNames.size(); iVar++) {
0188
0189         // make 1d name and title
0190         const string sName1D  = m_hist.vecBaseNames[iVar] + "_" + m_hist.vecTypeNames[iType];
0191         const string sTitle1D = ";" + get<0>(m_hist.vecAxisDef[iVar]) + ";counts";
0192
0193         // create 1d histogram
0194         vecHist1D[iType][iVar] = new TH1D(
0195           sName1D.data(),
0196           sTitle1D.data(),
0197           get<1>(m_hist.vecAxisDef[iVar]),
0198           get<2>(m_hist.vecAxisDef[iVar]).first,
0199           get<2>(m_hist.vecAxisDef[iVar]).second
0200         );
0201
0202         // loop over vs variables
0203         vecHist2D[iType][iVar].resize( m_hist.vecVsMods.size() );
0204         for (size_t iVs = 0; iVs < m_hist.vecVsMods.size(); iVs++) {
0205
0206           // make 2d name and title
0207           const string sName2D  = m_hist.vecBaseNames[iVar] + m_hist.vecVsMods[iVs] + "_" + m_hist.vecTypeNames[iType];
0208           const string sTitle2D = ";" + get<0>(m_hist.vecVsDef[iVs]) + ";" + get<0>(m_hist.vecAxisDef[iVar]) + ";counts";
0209
0210           // create 2d histogram
0211           vecHist2D[iType][iVar][iVs] = new TH2D(
0212             sName2D.data(),
0213             sTitle2D.data(),
0214             get<1>(m_hist.vecVsDef[iVs]),
0215             get<2>(m_hist.vecVsDef[iVs]).first,
0216             get<2>(m_hist.vecVsDef[iVs]).second,
0217             get<1>(m_hist.vecAxisDef[iVar]),
0218             get<2>(m_hist.vecAxisDef[iVar]).first,
0219             get<2>(m_hist.vecAxisDef[iVar]).second
0220           );
0221         }  // end vs loop
0222       }  // end variable loop
0223     }  // end type loop
0224     cout << "    Initialized histograms." << endl;
0225
0226     // exit routine
0227     return;
0228
0229   }  // end 'InitHists()'
0230
0231
0232
0233   void SReadLambdaJetTree::DoAnalysis() {
0234
0235     const int64_t nEvents = m_tInput -> GetEntries();
0236     cout << "    Beginning event loop: " << nEvents << " to process" << endl;
0237
0238     // for counting types
0239     vector<uint64_t> nTot(m_hist.vecTypeNames.size(), 0);
0240     vector<uint64_t> nEvt(m_hist.vecTypeNames.size(), 0);
0241
0242     int64_t nBytes = 0;
0243     for (int64_t iEvt = 0; iEvt < nEvents; iEvt++) {
0244
0245       // grab event
0246       const int64_t bytes = m_tInput -> GetEntry(iEvt);
0247       if (bytes < 0) {
0248         cerr << "WARNING: issue with event " << iEvt << "! Aborting event loop!" << endl;
0249         break;
0250       } else {
0251         nBytes += bytes;
0252       }
0253
0254       // announce progress
0255       const int64_t iProg = iEvt + 1;
0256       if (iProg == nEvents) {
0257         cout << "      Processing entry " << iEvt << "/" << nEvents << "..." << endl;
0258       } else {
0259         cout << "      Processing entry " << iEvt << "/" << nEvents << "...\r" << flush;
0260       }
0261
0262       // reset per-event counters
0263       for (size_t iType = 0; iType < nEvt.size(); iType++) {
0264         nEvt[iType] = 0;
0265       }
0266
0267       // get total no. of jets, lambdas in vectors
0268       const size_t nVecJets = m_jetPt    -> size();
0269       const size_t nVecLams = m_lambdaPt -> size();
0270
0271       // identify highest pt jet
0272       bool     foundTopPt = false;
0273       double   ptTop      = 0.;
0274       uint64_t iTopPt     = 0;
0275       uint64_t nLamTop    = 0;
0276       for (size_t iJet = 0; iJet < nVecJets; iJet++) {
0277
0278         // make sure jet satisfies cuts
0279         const bool isGoodJet = IsGoodJet(m_jetPt -> at(iJet), m_jetEta -> at(iJet));
0280         if (!isGoodJet) continue;
0281
0282         // find associated lambda(s)
0283         uint64_t nLamJet    = 0;
0284         for (size_t iLam = 0; iLam < nVecLams; iLam++) {
0285
0286           // check if lambda is associated
0287           const bool isAssocLam = IsAssociatedLambda(m_lambdaJetID -> at(iLam), m_jetID -> at(iJet));
0288           if (!isAssocLam) {
0289             continue;
0290           }
0291           ++nLamJet;
0292         }  // end lambda loop
0293
0294         if (m_jetPt -> at(iJet) > ptTop) {
0295           ptTop      = m_jetPt -> at(iJet);
0296           iTopPt     = iJet;
0297           nLamTop    = nLamJet;
0298           foundTopPt = true;
0299         }
0300       }  // end 1st jet loop
0301
0302       // fill highest pt histograms
0303       if (foundTopPt) {
0304         Hist hTopPtJet = {
0305           .eta  = m_jetEta -> at(iTopPt),
0306           .ene  = m_jetE   -> at(iTopPt),
0307           .pt   = m_jetPt  -> at(iTopPt),
0308           .df   = 0.,
0309           .dh   = 0.,
0310           .dr   = 0.,
0311           .z    = 1.,
0312           .nlam = (double) nLamTop,
0313           .ncst = (double) m_jetNCst -> at(iTopPt),
0314           .plam = (double) nLamTop / (double) m_jetNCst -> at(iTopPt)
0315         };
0316         VsVar vsTopPtJet = {
0317           .eta = m_jetEta -> at(iTopPt),
0318           .ene = m_jetE   -> at(iTopPt),
0319           .pt  = m_jetPt  -> at(iTopPt),
0320           .df  = 0.,
0321           .dh  = 0.
0322         };
0323         FillHist1D(Type::HJet, hTopPtJet);
0324         FillHist2D(Type::HJet, hTopPtJet, vsTopPtJet);
0325       }
0326
0327       // loop over lambdas
0328       for (size_t iLam = 0; iLam < nVecLams; iLam++) {
0329
0330         // make sure lambda satisfies cuts
0331         const bool isGoodLam = IsGoodLambda(m_lambdaPt -> at(iLam), m_lambdaEta -> at(iLam));
0332         if (!isGoodLam) continue;
0333
0334         // do calculations
0335         const double dfLam = GetDeltaPhi(m_lambdaPhi -> at(iLam), m_jetPhi -> at(iTopPt));
0336         const double dhLam = GetDeltaEta(m_lambdaEta -> at(iLam), m_jetEta -> at(iTopPt));
0337
0338         // fill general lambda histograms
0339         Hist hLambda = {
0340           .eta  = m_lambdaEta    -> at(iLam),
0341           .ene  = m_lambdaEnergy -> at(iLam),
0342           .pt   = m_lambdaPt     -> at(iLam),
0343           .df   = dfLam,
0344           .dh   = dhLam,
0345           .dr   = m_lambdaDr -> at(iLam),
0346           .z    = m_lambdaZ  -> at(iLam),
0347           .nlam = 1,
0348           .ncst = 0,
0349           .plam = 1
0350         };
0351         VsVar vsLambda = {
0352           .eta = m_lambdaEta    -> at(iLam),
0353           .ene = m_lambdaEnergy -> at(iLam),
0354           .pt  = m_lambdaPt     -> at(iLam),
0355           .df  = dfLam,
0356           .dh  = dhLam
0357         };
0358         FillHist1D(Type::Lam, hLambda);
0359         FillHist2D(Type::Lam, hLambda, vsLambda);
0360         ++nEvt[Type::Lam];
0361         ++nTot[Type::Lam];
0362
0363         const bool isLeadLam = IsLeadingLambda(m_lambdaZ -> at(iLam));
0364         if (isLeadLam) {
0365           FillHist1D(Type::LLam, hLambda);
0366           FillHist2D(Type::LLam, hLambda, vsLambda);
0367           ++nEvt[Type::LLam];
0368           ++nTot[Type::LLam];
0369         }
0370       }  // end lambda loop
0371
0372       // loop over jets
0373       for (size_t iJet = 0; iJet < nVecJets; iJet++) {
0374
0375         // make sure jet satisfies cuts
0376         const bool isGoodJet = IsGoodJet(m_jetPt -> at(iJet), m_jetEta -> at(iJet));
0377         if (!isGoodJet) continue;
0378
0379         // do calculations
0380         const double dfJet = GetDeltaPhi(m_jetPhi -> at(iJet), m_jetPhi -> at(iTopPt));
0381         const double dhJet = GetDeltaEta(m_jetEta -> at(iJet), m_jetEta -> at(iTopPt));
0382
0383         // find associated lambda(s)
0384         bool     hasLambda  = false;
0385         bool     hasLeadLam = false;
0386         uint64_t nLamJet    = 0;
0387         for (size_t iLam = 0; iLam < nVecLams; iLam++) {
0388
0389           // check if lambda is associated
0390           const bool isAssocLam = IsAssociatedLambda(m_lambdaJetID -> at(iLam), m_jetID -> at(iJet));
0391           if (!isAssocLam) {
0392             continue;
0393           } else {
0394             hasLambda = true;
0395           }
0396           ++nLamJet;
0397
0398           // check if lambda is leading
0399           const bool isLeadLam = IsLeadingLambda(m_lambdaZ -> at(iLam));
0400           if (isLeadLam) hasLeadLam = true;
0401
0402         }  // end lambda loop
0403         vecHistEvt.at(Evt::NLamJet) -> Fill(nLamJet);
0404
0405         // fill general jet histograms
0406         Hist hJet = {
0407           .eta  = m_jetEta -> at(iJet),
0408           .ene  = m_jetE   -> at(iJet),
0409           .pt   = m_jetPt  -> at(iJet),
0410           .df   = dfJet,
0411           .dh   = dhJet,
0412           .dr   = 0.,
0413           .z    = 1.,
0414           .nlam = (double) nLamJet,
0415           .ncst = (double) m_jetNCst -> at(iJet),
0416           .plam = (double) nLamJet / (double) m_jetNCst -> at(iJet)
0417         };
0418         VsVar vsJet = {
0419           .eta = m_jetEta -> at(iJet),
0420           .ene = m_jetE   -> at(iJet),
0421           .pt  = m_jetPt  -> at(iJet),
0422           .df  = dfJet,
0423           .dh  = dhJet
0424         };
0425         FillHist1D(Type::Jet, hJet);
0426         FillHist2D(Type::Jet, hJet, vsJet);
0427         ++nEvt[Type::Jet];
0428         ++nTot[Type::Jet];
0429
0430         // if no lambdas, continue
0431         //   otherwise fill hists
0432         if (!hasLambda) {
0433           continue;
0434         } else {
0435           FillHist1D(Type::LJet, hJet);
0436           FillHist2D(Type::LJet, hJet, vsJet);
0437           ++nEvt[Type::LJet];
0438           ++nTot[Type::LJet];
0439         }
0440
0441         // fill multi-lambda jet histogmras
0442         if (nLamJet >= 2) {
0443
0444           // fill jet histograms
0445           FillHist1D(Type::MLJet, hJet);
0446           FillHist2D(Type::MLJet, hJet, vsJet);
0447           ++nEvt[Type::MLJet];
0448           ++nTot[Type::MLJet];
0449
0450           // loop over lambdas in jet
0451           for (size_t iLam = 0; iLam < nVecLams; iLam++) {
0452
0453             // do calculations
0454             const double dfLam = GetDeltaPhi(m_lambdaPhi -> at(iLam), m_jetPhi -> at(iTopPt));
0455             const double dhLam = GetDeltaEta(m_lambdaEta -> at(iLam), m_jetEta -> at(iTopPt));
0456
0457             // check if lambda is associated
0458             const bool isAssocLam = IsAssociatedLambda(m_lambdaJetID -> at(iLam), m_jetID -> at(iJet));
0459             if (!isAssocLam) {
0460               continue;
0461             }
0462
0463             // fill lambda in multi-lambda jet histograms
0464             Hist hLamInMLJet = {
0465               .eta  = m_lambdaEta    -> at(iLam),
0466               .ene  = m_lambdaEnergy -> at(iLam),
0467               .pt   = m_lambdaPt     -> at(iLam),
0468               .df   = dfLam,
0469               .dh   = dhLam,
0470               .dr   = m_lambdaDr -> at(iLam),
0471               .z    = m_lambdaZ  -> at(iLam),
0472               .nlam = 1,
0473               .ncst = 0,
0474               .plam = 1
0475             };
0476             VsVar vsLamInMLJet = {
0477               .eta = m_lambdaEta    -> at(iLam),
0478               .ene = m_lambdaEnergy -> at(iLam),
0479               .pt  = m_lambdaPt     -> at(iLam),
0480               .df  = dfLam,
0481               .dh  = dhLam
0482             };
0483             FillHist1D(Type::MLJetLam, hLamInMLJet);
0484             FillHist2D(Type::MLJetLam, hLamInMLJet, vsLamInMLJet);
0485           }  // end lambda loop
0486         }  // end if (nLamJet >= 2)
0487
0488         // fill jet w/ leading lambda histograms
0489         if (hasLeadLam) {
0490           FillHist1D(Type::LLJet, hJet);
0491           FillHist2D(Type::LLJet, hJet, vsJet);
0492           ++nEvt[Type::LLJet];
0493           ++nTot[Type::LLJet];
0494         }
0495       }  // end 2nd jet loop
0496
0497       // fill event histograms
0498       vecHistEvt.at(Evt::NJet)      -> Fill(nEvt[Type::Jet]);
0499       vecHistEvt.at(Evt::NTagJet)   -> Fill(nEvt[Type::LJet]);
0500       vecHistEvt.at(Evt::NLeadJet)  -> Fill(nEvt[Type::LLJet]);
0501       vecHistEvt.at(Evt::NMultiJet) -> Fill(nEvt[Type::MLJet]);
0502       vecHistEvt.at(Evt::NLam)      -> Fill(nEvt[Type::Lam]);
0503       vecHistEvt.at(Evt::NLeadLam)  -> Fill(nEvt[Type::LLam]);
0504
0505     }  // end event loop
0506     cout << "    Event loop finished!\n"
0507          << "      nLambda          = " << nTot[Type::Lam]   << "\n"
0508          << "      nLeadLambda      = " << nTot[Type::LLam]  << "\n"
0509          << "      nJet             = " << nTot[Type::Jet]   << "\n"
0510          << "      nTaggedJets      = " << nTot[Type::LJet]  << "\n"
0511          << "      nLeadTagJets     = " << nTot[Type::LLJet] << "\n"
0512          << "      nMultiLambdaJets = " << nTot[Type::MLJet]
0513          << endl;
0514
0515     // exit routine
0516     return;
0517
0518   }  // end 'DoAnalysis()'
0519
0520
0521
0522   void SReadLambdaJetTree::SetHistogramStyles() {
0523
0524     for (auto hEvt : vecHistEvt) {
0525       hEvt -> SetLineStyle(m_config.defLineStyle);
0526       hEvt -> SetLineColor(m_config.defHistColor);
0527       hEvt -> SetFillStyle(m_config.defFillStyle);
0528       hEvt -> SetFillColor(m_config.defHistColor);
0529       hEvt -> SetMarkerStyle(m_config.defMarkStyle);
0530       hEvt -> SetMarkerColor(m_config.defHistColor);
0531       hEvt -> GetXaxis() -> CenterTitle(m_config.centerTitle);
0532       hEvt -> GetXaxis() -> SetTitleFont(m_config.defHistFont);
0533       hEvt -> GetXaxis() -> SetTitleSize(m_config.defTitleX);
0534       hEvt -> GetXaxis() -> SetTitleOffset(m_config.defOffX);
0535       hEvt -> GetXaxis() -> SetLabelSize(m_config.defLabelX);
0536       hEvt -> GetYaxis() -> CenterTitle(m_config.centerTitle);
0537       hEvt -> GetYaxis() -> SetTitleFont(m_config.defHistFont);
0538       hEvt -> GetYaxis() -> SetTitleSize(m_config.defTitleY);
0539       hEvt -> GetYaxis() -> SetTitleOffset(m_config.defOffY);
0540       hEvt -> GetYaxis() -> SetLabelSize(m_config.defLabelY);
0541     }
0542     for (auto type : vecHist1D) {
0543       for (auto h1D : type) {
0544         h1D -> SetLineStyle(m_config.defLineStyle);
0545         h1D -> SetLineColor(m_config.defHistColor);
0546         h1D -> SetFillStyle(m_config.defFillStyle);
0547         h1D -> SetFillColor(m_config.defHistColor);
0548         h1D -> SetMarkerStyle(m_config.defMarkStyle);
0549         h1D -> SetMarkerColor(m_config.defHistColor);
0550         h1D -> GetXaxis() -> CenterTitle(m_config.centerTitle);
0551         h1D -> GetXaxis() -> SetTitleFont(m_config.defHistFont);
0552         h1D -> GetXaxis() -> SetTitleSize(m_config.defTitleX);
0553         h1D -> GetXaxis() -> SetTitleOffset(m_config.defOffX);
0554         h1D -> GetXaxis() -> SetLabelSize(m_config.defLabelX);
0555         h1D -> GetYaxis() -> CenterTitle(m_config.centerTitle);
0556         h1D -> GetYaxis() -> SetTitleFont(m_config.defHistFont);
0557         h1D -> GetYaxis() -> SetTitleSize(m_config.defTitleY);
0558         h1D -> GetYaxis() -> SetTitleOffset(m_config.defOffY);
0559         h1D -> GetYaxis() -> SetLabelSize(m_config.defLabelY);
0560       }
0561     }
0562     for (auto type : vecHist2D) {
0563       for (auto var : type) {
0564         for (auto h2D : var) {
0565           h2D -> SetLineStyle(m_config.defLineStyle);
0566           h2D -> SetLineColor(m_config.defHistColor);
0567           h2D -> SetFillStyle(m_config.defFillStyle);
0568           h2D -> SetFillColor(m_config.defHistColor);
0569           h2D -> SetMarkerStyle(m_config.defMarkStyle);
0570           h2D -> SetMarkerColor(m_config.defHistColor);
0571           h2D -> GetXaxis() -> CenterTitle(m_config.centerTitle);
0572           h2D -> GetXaxis() -> SetTitleFont(m_config.defHistFont);
0573           h2D -> GetXaxis() -> SetTitleSize(m_config.defTitleX);
0574           h2D -> GetXaxis() -> SetTitleOffset(m_config.defOffX);
0575           h2D -> GetXaxis() -> SetLabelSize(m_config.defLabelX);
0576           h2D -> GetYaxis() -> CenterTitle(m_config.centerTitle);
0577           h2D -> GetYaxis() -> SetTitleFont(m_config.defHistFont);
0578           h2D -> GetYaxis() -> SetTitleSize(m_config.defTitleY);
0579           h2D -> GetYaxis() -> SetTitleOffset(m_config.defOffY);
0580           h2D -> GetYaxis() -> SetLabelSize(m_config.defLabelY);
0581           h2D -> GetZaxis() -> CenterTitle(m_config.centerTitle);
0582           h2D -> GetZaxis() -> SetTitleFont(m_config.defHistFont);
0583           h2D -> GetZaxis() -> SetTitleSize(m_config.defTitleZ);
0584           h2D -> GetZaxis() -> SetTitleOffset(m_config.defOffZ);
0585           h2D -> GetZaxis() -> SetLabelSize(m_config.defLabelZ);
0586         }
0587       }
0588     }
0589     cout << "    Set histogram styles." << endl;
0590
0591     // exit internal routine
0592     return;
0593
0594   }  // end 'SetHistogramStyles()'
0595
0596
0597
0598   void SReadLambdaJetTree::SaveOutput() {
0599
0600     m_outDir -> cd();
0601     for (auto hEvt : vecHistEvt) {
0602       hEvt -> Write();
0603     }
0604     for (auto type : vecHist1D) {
0605       for (auto h1D : type) {
0606         h1D -> Write();
0607       }
0608     }
0609     for (auto type : vecHist2D) {
0610       for (auto var : type) {
0611         for (auto h2D : var) {
0612           h2D -> Write();
0613         }
0614       }
0615     }
0616     cout << "    Saved histograms." << endl;
0617
0618     // exit internal routine
0619     return;
0620
0621   }  // end 'SaveOutput()'
0622
0623
0624
0625   void SReadLambdaJetTree::CloseInput() {
0626
0627     m_fInput -> cd();
0628     m_fInput -> Close();
0629     return;
0630
0631   }  // end 'CloseInput()'
0632
0633
0634
0635   void SReadLambdaJetTree::FillHist1D(const int type, Hist hist) {
0636
0637     vecHist1D.at(type)[Var::Eta]  -> Fill(hist.eta);
0638     vecHist1D.at(type)[Var::Ene]  -> Fill(hist.ene);
0639     vecHist1D.at(type)[Var::Pt]   -> Fill(hist.pt);
0640     vecHist1D.at(type)[Var::DPhi] -> Fill(hist.df);
0641     vecHist1D.at(type)[Var::DEta] -> Fill(hist.dh);
0642     vecHist1D.at(type)[Var::Dr]   -> Fill(hist.dr);
0643     vecHist1D.at(type)[Var::Z]    -> Fill(hist.z);
0644     vecHist1D.at(type)[Var::NL]   -> Fill(hist.nlam);
0645     vecHist1D.at(type)[Var::NC]   -> Fill(hist.ncst);
0646     vecHist1D.at(type)[Var::RLC]  -> Fill(hist.plam);
0647     return;
0648
0649   }  // end 'FillHist1D(int, Hist)'
0650
0651
0652
0653   void SReadLambdaJetTree::FillHist2D(const int type, Hist hist, VsVar vs) {
0654
0655     // fill vs. eta histograms
0656     vecHist2D.at(type)[Var::Eta][Mod::VsEta]  -> Fill(vs.eta, hist.eta);
0657     vecHist2D.at(type)[Var::Ene][Mod::VsEta]  -> Fill(vs.eta, hist.ene);
0658     vecHist2D.at(type)[Var::Pt][Mod::VsEta]   -> Fill(vs.eta, hist.pt);
0659     vecHist2D.at(type)[Var::DPhi][Mod::VsEta] -> Fill(vs.eta, hist.df);
0660     vecHist2D.at(type)[Var::DEta][Mod::VsEta] -> Fill(vs.eta, hist.dh);
0661     vecHist2D.at(type)[Var::Dr][Mod::VsEta]   -> Fill(vs.eta, hist.dr);
0662     vecHist2D.at(type)[Var::Z][Mod::VsEta]    -> Fill(vs.eta, hist.z);
0663     vecHist2D.at(type)[Var::NL][Mod::VsEta]   -> Fill(vs.eta, hist.nlam);
0664     vecHist2D.at(type)[Var::NC][Mod::VsEta]    -> Fill(vs.eta, hist.ncst);
0665     vecHist2D.at(type)[Var::RLC][Mod::VsEta]  -> Fill(vs.eta, hist.plam);
0666
0667     // fill vs. energy histograms
0668     vecHist2D.at(type)[Var::Eta][Mod::VsEne]  -> Fill(vs.ene, hist.eta);
0669     vecHist2D.at(type)[Var::Ene][Mod::VsEne]  -> Fill(vs.ene, hist.ene);
0670     vecHist2D.at(type)[Var::Pt][Mod::VsEne]   -> Fill(vs.ene, hist.pt);
0671     vecHist2D.at(type)[Var::DPhi][Mod::VsEne] -> Fill(vs.ene, hist.df);
0672     vecHist2D.at(type)[Var::DEta][Mod::VsEne] -> Fill(vs.ene, hist.dh);
0673     vecHist2D.at(type)[Var::Dr][Mod::VsEne]   -> Fill(vs.ene, hist.dr);
0674     vecHist2D.at(type)[Var::Z][Mod::VsEne]    -> Fill(vs.ene, hist.z);
0675     vecHist2D.at(type)[Var::NL][Mod::VsEne]   -> Fill(vs.ene, hist.nlam);
0676     vecHist2D.at(type)[Var::NC][Mod::VsEne]    -> Fill(vs.ene, hist.ncst);
0677     vecHist2D.at(type)[Var::RLC][Mod::VsEne]  -> Fill(vs.ene, hist.plam);
0678
0679     // fill vs. pt histograms
0680     vecHist2D.at(type)[Var::Eta][Mod::VsPt]  -> Fill(vs.pt, hist.eta);
0681     vecHist2D.at(type)[Var::Ene][Mod::VsPt]  -> Fill(vs.pt, hist.ene);
0682     vecHist2D.at(type)[Var::Pt][Mod::VsPt]   -> Fill(vs.pt, hist.pt);
0683     vecHist2D.at(type)[Var::DPhi][Mod::VsPt] -> Fill(vs.pt, hist.df);
0684     vecHist2D.at(type)[Var::DEta][Mod::VsPt] -> Fill(vs.pt, hist.dh);
0685     vecHist2D.at(type)[Var::Dr][Mod::VsPt]   -> Fill(vs.pt, hist.dr);
0686     vecHist2D.at(type)[Var::Z][Mod::VsPt]    -> Fill(vs.pt, hist.z);
0687     vecHist2D.at(type)[Var::NL][Mod::VsPt]   -> Fill(vs.pt, hist.nlam);
0688     vecHist2D.at(type)[Var::NC][Mod::VsPt]    -> Fill(vs.pt, hist.ncst);
0689     vecHist2D.at(type)[Var::RLC][Mod::VsPt]  -> Fill(vs.pt, hist.plam);
0690
0691     // fill vs. delta-phi histograms
0692     vecHist2D.at(type)[Var::Eta][Mod::VsDPhi]  -> Fill(vs.df, hist.eta);
0693     vecHist2D.at(type)[Var::Ene][Mod::VsDPhi]  -> Fill(vs.df, hist.ene);
0694     vecHist2D.at(type)[Var::Pt][Mod::VsDPhi]   -> Fill(vs.df, hist.pt);
0695     vecHist2D.at(type)[Var::DPhi][Mod::VsDPhi] -> Fill(vs.df, hist.df);
0696     vecHist2D.at(type)[Var::DEta][Mod::VsDPhi] -> Fill(vs.df, hist.dh);
0697     vecHist2D.at(type)[Var::Dr][Mod::VsDPhi]   -> Fill(vs.df, hist.dr);
0698     vecHist2D.at(type)[Var::Z][Mod::VsDPhi]    -> Fill(vs.df, hist.z);
0699     vecHist2D.at(type)[Var::NL][Mod::VsDPhi]   -> Fill(vs.df, hist.nlam);
0700     vecHist2D.at(type)[Var::NC][Mod::VsDPhi]    -> Fill(vs.df, hist.ncst);
0701     vecHist2D.at(type)[Var::RLC][Mod::VsDPhi]  -> Fill(vs.df, hist.plam);
0702
0703     // fill vs. delta-eta histograms
0704     vecHist2D.at(type)[Var::Eta][Mod::VsDEta]  -> Fill(vs.dh, hist.eta);
0705     vecHist2D.at(type)[Var::Ene][Mod::VsDEta]  -> Fill(vs.dh, hist.ene);
0706     vecHist2D.at(type)[Var::Pt][Mod::VsDEta]   -> Fill(vs.dh, hist.pt);
0707     vecHist2D.at(type)[Var::DPhi][Mod::VsDEta] -> Fill(vs.dh, hist.df);
0708     vecHist2D.at(type)[Var::DEta][Mod::VsDEta] -> Fill(vs.dh, hist.dh);
0709     vecHist2D.at(type)[Var::Dr][Mod::VsDEta]   -> Fill(vs.dh, hist.dr);
0710     vecHist2D.at(type)[Var::Z][Mod::VsDEta]    -> Fill(vs.dh, hist.z);
0711     vecHist2D.at(type)[Var::NL][Mod::VsDEta]   -> Fill(vs.dh, hist.nlam);
0712     vecHist2D.at(type)[Var::NC][Mod::VsDEta]    -> Fill(vs.dh, hist.ncst);
0713     vecHist2D.at(type)[Var::RLC][Mod::VsDEta]  -> Fill(vs.dh, hist.plam);
0714     return;
0715
0716   }  // end 'FillHist2D(int, Hist, VsVar)'
0717
0718
0719
0720   bool SReadLambdaJetTree::IsGoodJet(const double pt, const double eta) {
0721
0722     const bool isGoodPt  = (pt > m_config.ptJetMin);
0723     const bool isGoodEta = (abs(eta) < m_config.etaJetMax);
0724     const bool isGoodJet = (isGoodPt && isGoodEta);
0725     return isGoodJet;
0726
0727   }  // end 'IsGoodJet(double, double)'
0728
0729
0730
0731   bool SReadLambdaJetTree::IsGoodLambda(const double pt, const double eta) {
0732
0733     const bool isGoodPt  = (pt > m_config.ptLamMin);
0734     const bool isGoodEta = (abs(eta) < m_config.etaLamMax);
0735     const bool isGoodLam = (isGoodPt && isGoodEta);
0736     return isGoodLam;
0737
0738   }  // end 'IsGoodLambda(double, double)'
0739
0740
0741
0742   bool SReadLambdaJetTree::IsLeadingLambda(const double z) {
0743
0744     const bool isLeadLam = (z > m_config.zLeadMin);
0745     return isLeadLam;
0746
0747   }  // end 'IsLeadingLambda(double)'
0748
0749
0750
0751   bool SReadLambdaJetTree::IsAssociatedLambda(const int idLam, const int idJet) {
0752
0753     const bool isAssoc = (idLam == idJet);
0754     return isAssoc;
0755
0756   }  // end 'IsAssociatedLambda(int, int)'
0757
0758
0759
0760   double SReadLambdaJetTree::GetDeltaPhi(const double phiA, const double phiB) {
0761
0762     double dPhi = phiA - phiB;
0763     if (dPhi < m_const.minDPhi) dPhi += TMath::TwoPi();
0764     if (dPhi > m_const.maxDPhi) dPhi -= TMath::TwoPi();
0765     return dPhi;
0766
0767   }  // end 'GetDeltaPhi(double, double)'
0768
0769
0770
0771   double SReadLambdaJetTree::GetDeltaEta(const double etaA, const double etaB) {
0772
0773     const double dEta = etaA - etaB;
0774     return dEta;
0775
0776   }  // end 'GetDeltaEta(double, double)'
0777
0778 }  // end SColdQcdCorrelatorAnalysis namespace
0779
0780 // end ------------------------------------------------------------------------