Back to home page

EIC code displayed by LXR
 
 

    

File indexing completed on 2025-01-18 09:27:06 

0001 
0002 // Jet Benchmarks
0003 // Author: B. Page (bpage@bnl.gov)
0004 
0005 // To run: root -b -l -q jetBenchmarks.C'("path/to/files",true)' Second argument prints histograms
0006 
0007 void jetBenchmarks(TString infile="root://dtn-eic.jlab.org//work/eic2/EPIC/RECO/24.05.0/epic_craterlake/DIS/NC/18x275/minQ2=10/pythia8NCDIS_18x275_minQ2=10_beamEffects_xAngle=-0.025_hiDiv_5.2501.eicrecon.tree.edm4eic.root", bool PRINT=false, TString outputDir="")
0008 {
0009   // Input
0010   TChain *mychain = new TChain("events");
0011   mychain->Add(infile);
0012 
0013   const int seabornRed = TColor::GetColor(213, 94, 0);
0014 
0015 // Seaborn Green: #009E73 -> (0, 158, 115)
0016 const int seabornGreen = TColor::GetColor(0, 158, 115);
0017 
0018 // Seaborn Blue: #56B4E9 -> (86, 180, 233)
0019 const int seabornBlue = TColor::GetColor(100, 149, 237);
0020 
0021   // Output
0022   //TFile *ofile = TFile::Open("test_24-05-0.hist.root","RECREATE");
0023 
0024   // TTreeReader
0025   TTreeReader tree_reader(mychain);
0026 
0027   // Set Delta R Cut
0028   float DELTARCUT = 0.05;
0029 
0030   // Reco Jets
0031   TTreeReaderArray<int> recoType = {tree_reader, "ReconstructedChargedJets.type"};
0032   TTreeReaderArray<float> recoNRG = {tree_reader, "ReconstructedChargedJets.energy"};
0033   TTreeReaderArray<int> recoPDG = {tree_reader, "ReconstructedChargedJets.PDG"};
0034   TTreeReaderArray<float> recoMomX = {tree_reader, "ReconstructedChargedJets.momentum.x"};
0035   TTreeReaderArray<float> recoMomY = {tree_reader, "ReconstructedChargedJets.momentum.y"};
0036   TTreeReaderArray<float> recoMomZ = {tree_reader, "ReconstructedChargedJets.momentum.z"};
0037   TTreeReaderArray<float> recoM = {tree_reader, "ReconstructedChargedJets.mass"};
0038   TTreeReaderArray<unsigned int> partsBegin = {tree_reader, "ReconstructedChargedJets.particles_begin"};
0039   TTreeReaderArray<unsigned int> partsEnd = {tree_reader, "ReconstructedChargedJets.particles_end"};
0040 
0041   TTreeReaderArray<int> recoClustIndex = {tree_reader, "_ReconstructedChargedJets_clusters.index"};
0042   TTreeReaderArray<int> recoTrackIndex = {tree_reader, "_ReconstructedChargedJets_tracks.index"};
0043   TTreeReaderArray<int> recoPartIndex = {tree_reader, "_ReconstructedChargedJets_particles.index"};
0044 
0045   // Reconstructed Particles
0046   TTreeReaderArray<float> recoPartMomX = {tree_reader, "ReconstructedChargedParticles.momentum.x"};
0047   TTreeReaderArray<float> recoPartMomY = {tree_reader, "ReconstructedChargedParticles.momentum.y"};
0048   TTreeReaderArray<float> recoPartMomZ = {tree_reader, "ReconstructedChargedParticles.momentum.z"};
0049   TTreeReaderArray<float> recoPartM = {tree_reader, "ReconstructedChargedParticles.mass"};
0050   TTreeReaderArray<int> recoPartPDG = {tree_reader, "ReconstructedChargedParticles.PDG"};
0051   TTreeReaderArray<float> recoPartNRG = {tree_reader, "ReconstructedChargedParticles.energy"};
0052 
0053   // Generated Jets
0054   TTreeReaderArray<int> genType = {tree_reader, "GeneratedChargedJets.type"};
0055   TTreeReaderArray<float> genNRG = {tree_reader, "GeneratedChargedJets.energy"};
0056   TTreeReaderArray<int> genPDG = {tree_reader, "GeneratedChargedJets.PDG"};
0057   TTreeReaderArray<float> genMomX = {tree_reader, "GeneratedChargedJets.momentum.x"};
0058   TTreeReaderArray<float> genMomY = {tree_reader, "GeneratedChargedJets.momentum.y"};
0059   TTreeReaderArray<float> genMomZ = {tree_reader, "GeneratedChargedJets.momentum.z"};
0060   TTreeReaderArray<float> genM = {tree_reader, "GeneratedChargedJets.mass"};
0061   TTreeReaderArray<unsigned int> genPartsBegin = {tree_reader, "GeneratedChargedJets.particles_begin"};
0062   TTreeReaderArray<unsigned int> genPartsEnd = {tree_reader, "GeneratedChargedJets.particles_end"};
0063   
0064   TTreeReaderArray<int> genPartIndex = {tree_reader, "_GeneratedChargedJets_particles.index"};
0065   //TTreeReaderArray<int> genChargedIndex = {tree_reader, "GeneratedChargedParticles_objIdx.index"};
0066   
0067   // MC
0068   //TTreeReaderArray<int> mcGenStat = {tree_reader, "MCParticles.generatorStatus"};
0069   TTreeReaderArray<float> mcMomX = {tree_reader, "GeneratedParticles.momentum.x"};
0070   TTreeReaderArray<float> mcMomY = {tree_reader, "GeneratedParticles.momentum.y"};
0071   TTreeReaderArray<float> mcMomZ = {tree_reader, "GeneratedParticles.momentum.z"};
0072   TTreeReaderArray<float> mcM = {tree_reader, "GeneratedParticles.mass"};
0073   TTreeReaderArray<int> pdg = {tree_reader, "GeneratedParticles.PDG"};
0074 
0075   // Define Histograms
0076   TH1D *counter = new TH1D("counter","",10,0.,10.);
0077 
0078   
0079   // Reco
0080   TH1D *numRecoChargedJetsECutHist = new TH1D("numRecoChargedJetsECut","",20,0.,20.);
0081   TH1D *recoChargedJetEHist = new TH1D("recoChargedJetE","",300,0.,300.);
0082   TH1D *recoChargedJetEtaECutHist = new TH1D("recoChargedJetEtaECut","",60,-3.,3.);
0083   TH2D *recoChargedJetEvsEtaHist = new TH2D("recoChargedJetEvsEta","",60,-3.,3.,300,0.,300.);
0084   TH2D *recoChargedJetPhiVsEtaECutHist = new TH2D("recoChargedJetPhiVsEtaECut","",60,-3.,3.,100,-TMath::Pi(),TMath::Pi());
0085 
0086   TH1D *numRecoChargedJetsECutNoElecHist = new TH1D("numRecoChargedJetsECutNoElec","",20,0.,20.);
0087   TH1D *recoChargedJetENoElecHist = new TH1D("recoChargedJetENoElec","",300,0.,300.);
0088   TH1D *recoChargedJetEtaECutNoElecHist = new TH1D("recoChargedJetEtaECutNoElec","",60,-3.,3.);
0089   TH2D *recoChargedJetEvsEtaNoElecHist = new TH2D("recoChargedJetEvsEtaNoElec","",60,-3.,3.,300,0.,300.);
0090   TH2D *recoChargedJetPhiVsEtaECutNoElecHist = new TH2D("recoChargedJetPhiVsEtaECutNoElec","",60,-3.,3.,100,-TMath::Pi(),TMath::Pi());
0091 
0092   TH1D *numRecoChargedJetPartsHist = new TH1D("numRecoChargedJetParts","",20,0.,20.);
0093   TH1D *recoChargedJetPartEHist = new TH1D("recoChargedJetPartE","",500,0.,100.);
0094   TH1D *recoChargedJetPartEtaHist = new TH1D("recoChargedJetPartEta","",80,-4.,4.);
0095   TH2D *recoChargedJetPartEvsEtaHist = new TH2D("recoChargedJetPartEvsEta","",80,-4.,4.,500,0.,100.);
0096   TH2D *recoChargedJetPartPhiVsEtaHist = new TH2D("recoChargedJetPartPhiVsEta","",80,-4.,4.,100,-TMath::Pi(),TMath::Pi());
0097 
0098   TH1D *numRecoChargedJetPartsNoElecHist = new TH1D("numRecoChargedJetPartsNoElec","",20,0.,20.);
0099   TH1D *recoChargedJetPartENoElecHist = new TH1D("recoChargedJetPartENoElec","",500,0.,100.);
0100   TH1D *recoChargedJetPartEtaNoElecHist = new TH1D("recoChargedJetPartEtaNoElec","",80,-4.,4.);
0101   TH2D *recoChargedJetPartEvsEtaNoElecHist = new TH2D("recoChargedJetPartEvsEtaNoElec","",80,-4.,4.,500,0.,100.);
0102   TH2D *recoChargedJetPartPhiVsEtaNoElecHist = new TH2D("recoChargedJetPartPhiVsEtaNoElec","",80,-4.,4.,100,-TMath::Pi(),TMath::Pi());
0103 
0104   TH1D *recoChargedJetPartPairwiseDeltaRHist = new TH1D("recoChargedJetPartPairwiseDeltaRHist","",5000,0.,5.);
0105 
0106   // Gen
0107   TH1D *numGenChargedJetsECutHist = new TH1D("numGenChargedJetsECut","",20,0.,20.);
0108   TH1D *genChargedJetEHist = new TH1D("genChargedJetE","",300,0.,300.);
0109   TH1D *genChargedJetEtaECutHist = new TH1D("genChargedJetEtaECut","",60,-3.,3.);
0110   TH2D *genChargedJetEvsEtaHist = new TH2D("genChargedJetEvsEta","",60,-3.,3.,300,0.,300.);
0111   TH2D *genChargedJetPhiVsEtaECutHist = new TH2D("genChargedJetPhiVsEtaECut","",60,-3.,3.,100,-TMath::Pi(),TMath::Pi());
0112 
0113   TH1D *numGenChargedJetsECutNoElecHist = new TH1D("numGenChargedJetsECutNoElec","",20,0.,20.);
0114   TH1D *genChargedJetENoElecHist = new TH1D("genChargedJetENoElec","",300,0.,300.);
0115   TH1D *genChargedJetEtaECutNoElecHist = new TH1D("genChargedJetEtaECutNoElec","",60,-3.,3.);
0116   TH2D *genChargedJetEvsEtaNoElecHist = new TH2D("genChargedJetEvsEtaNoElec","",60,-3.,3.,300,0.,300.);
0117   TH2D *genChargedJetPhiVsEtaECutNoElecHist = new TH2D("genChargedJetPhiVsEtaECutNoElec","",60,-3.,3.,100,-TMath::Pi(),TMath::Pi());
0118 
0119   TH1D *numGenChargedJetPartsHist = new TH1D("numGenChargedJetParts","",20,0.,20.);
0120   TH1D *genChargedJetPartEHist = new TH1D("genChargedJetPartE","",500,0.,100.);
0121   TH1D *genChargedJetPartEtaHist = new TH1D("genChargedJetPartEta","",80,-4.,4.);
0122   TH2D *genChargedJetPartEvsEtaHist = new TH2D("genChargedJetPartEvsEta","",80,-4.,4.,500,0.,100.);
0123   TH2D *genChargedJetPartPhiVsEtaHist = new TH2D("genChargedJetPartPhiVsEta","",80,-4.,4.,100,-TMath::Pi(),TMath::Pi());
0124 
0125   TH1D *numGenChargedJetPartsNoElecHist = new TH1D("numGenChargedJetPartsNoElec","",20,0.,20.);
0126   TH1D *genChargedJetPartENoElecHist = new TH1D("genChargedJetPartENoElec","",500,0.,100.);
0127   TH1D *genChargedJetPartEtaNoElecHist = new TH1D("genChargedJetPartEtaNoElec","",80,-4.,4.);
0128   TH2D *genChargedJetPartEvsEtaNoElecHist = new TH2D("genChargedJetPartEvsEtaNoElec","",80,-4.,4.,500,0.,100.);
0129   TH2D *genChargedJetPartPhiVsEtaNoElecHist = new TH2D("genChargedJetPartPhiVsEtaNoElec","",80,-4.,4.,100,-TMath::Pi(),TMath::Pi());
0130 
0131   TH1D *genChargedJetPartPairwiseDeltaRHist = new TH1D("genChargedJetPartPairwiseDeltaRHist","",5000,0.,5.);
0132 
0133   // Matched
0134   TH1D *matchJetDeltaRHist = new TH1D("matchJetDeltaR","",5000,0.,5.);
0135   TH1D *matchJetDeltaRBackHist = new TH1D("matchJetDeltaRBack","",5000,0.,5.);
0136   TH2D *recoVsGenChargedJetEtaHist = new TH2D("recoVsGenChargedJetEta","",80,-4.,4.,80,-4.,4.);
0137   TH2D *recoVsGenChargedJetPhiHist = new TH2D("recoVsGenChargedJetPhi","",100,-TMath::Pi(),TMath::Pi(),100,-TMath::Pi(),TMath::Pi());
0138   TH2D *recoVsGenChargedJetEHist = new TH2D("recoVsGenChargedJetE","",100,0.,100.,100,0.,100.);
0139   TH2D *recoVsGenChargedJetENoDRHist = new TH2D("recoVsGenChargedJetENoDRHist","",100,0.,100.,100,0.,100.);
0140   TH2D *recoVsGenChargedJetENoDupHist = new TH2D("recoVsGenChargedJetENoDup","",100,0.,100.,100,0.,100.);
0141 
0142   TH2D *jetResVsEtaHist = new TH2D("jetResVsEta","",80,-4.,4.,10000,-10.,10.);
0143   TH2D *jetResVsEHist = new TH2D("jetResVsE","",100,0.,100.,10000,-10.,10.);
0144   TH2D *jetResVsENegEtaHist = new TH2D("jetResVsENegEta","",20,0.,100.,10000,-10.,10.);
0145   TH2D *jetResVsEMidEtaHist = new TH2D("jetResVsEMidEta","",20,0.,100.,10000,-10.,10.);
0146   TH2D *jetResVsEPosEtaHist = new TH2D("jetResVsEPosEta","",20,0.,100.,10000,-10.,10.);
0147 
0148   TH2D *jetResVsENegEtaNoDupHist = new TH2D("jetResVsENegEtaNoDup","",20,0.,100.,10000,-10.,10.);
0149   TH2D *jetResVsEMidEtaNoDupHist = new TH2D("jetResVsEMidEtaNoDup","",20,0.,100.,10000,-10.,10.);
0150   TH2D *jetResVsEPosEtaNoDupHist = new TH2D("jetResVsEPosEtaNoDup","",20,0.,100.,10000,-10.,10.);
0151 
0152 
0153   // Loop Through Events
0154   int NEVENTS = 0;
0155   while(tree_reader.Next()) {
0156 
0157     if(NEVENTS%10000 == 0) cout << "Events Processed: " << NEVENTS << endl;
0158 
0159     counter->Fill(0);
0160 
0161     //////////////////////////////////////////////////////////////////////////
0162     //////////////////////  Analyze Reconstructed Jets  //////////////////////
0163     //////////////////////////////////////////////////////////////////////////
0164     int numRecoChargedJets = 0;
0165     int numRecoChargedJetsNoElec = 0;
0166     for(unsigned int i=0; i<recoType.GetSize(); i++)
0167       {
0168     TVector3 jetMom(recoMomX[i],recoMomY[i],recoMomZ[i]);
0169 
0170     counter->Fill(3);
0171 
0172     // Place eta cut to avoid edges of tracking acceptance
0173     if(TMath::Abs(jetMom.PseudoRapidity()) > 2.5) continue;
0174 
0175     // Place a minimum energy condition for several plots
0176     bool ECut = recoNRG[i] > 5.0;
0177 
0178     if(ECut) numRecoChargedJets++; 
0179 
0180     recoChargedJetEHist->Fill(recoNRG[i]);
0181     if(ECut) recoChargedJetEtaECutHist->Fill(jetMom.PseudoRapidity());
0182     recoChargedJetEvsEtaHist->Fill(jetMom.PseudoRapidity(),recoNRG[i]);
0183     if(ECut) recoChargedJetPhiVsEtaECutHist->Fill(jetMom.PseudoRapidity(),jetMom.Phi());
0184 
0185     bool noElectron = true;
0186     if(ECut)
0187       {
0188         // Find Jets with Electrons
0189         for(unsigned int m=partsBegin[i]; m<partsEnd[i]; m++)
0190           {
0191         if(recoPartM[recoPartIndex[m]] > 0.00050 && recoPartM[recoPartIndex[m]] < 0.00052)
0192           noElectron = false;
0193           }
0194 
0195         for(unsigned int j=partsBegin[i]; j<partsEnd[i]; j++)
0196           {
0197         // partsbegin and partsEnd specify the entries from _ReconstructedChargedJets_particles.index that make up the jet
0198         // _ReconstructedChargedJets_particles.index stores the ReconstructedChargedParticles index of the jet constituent
0199         double mX = recoPartMomX[recoPartIndex[j]];
0200         double mY = recoPartMomY[recoPartIndex[j]];
0201         double mZ = recoPartMomZ[recoPartIndex[j]];
0202         double mM = recoPartM[recoPartIndex[j]];
0203         double tmpE = TMath::Sqrt(mX*mX + mY*mY + mZ*mZ + mM*mM);
0204         
0205         TVector3 partMom(mX,mY,mZ);
0206         
0207         recoChargedJetPartEHist->Fill(tmpE);
0208         recoChargedJetPartEtaHist->Fill(partMom.PseudoRapidity());
0209         recoChargedJetPartEvsEtaHist->Fill(partMom.PseudoRapidity(),tmpE);
0210         recoChargedJetPartPhiVsEtaHist->Fill(partMom.PseudoRapidity(),partMom.Phi());
0211 
0212         if(noElectron)
0213           {
0214             recoChargedJetPartENoElecHist->Fill(tmpE);
0215             recoChargedJetPartEtaNoElecHist->Fill(partMom.PseudoRapidity());
0216             recoChargedJetPartEvsEtaNoElecHist->Fill(partMom.PseudoRapidity(),tmpE);
0217             recoChargedJetPartPhiVsEtaNoElecHist->Fill(partMom.PseudoRapidity(),partMom.Phi());
0218           }
0219 
0220         // Pairwise Distance Between Constituents
0221         if(j<(partsEnd[i]-1))
0222           {
0223             for(unsigned int k=j+1; k<partsEnd[i]; k++)
0224               {
0225             double mXB = recoPartMomX[recoPartIndex[k]];
0226             double mYB = recoPartMomY[recoPartIndex[k]];
0227             double mZB = recoPartMomZ[recoPartIndex[k]];
0228 
0229             TVector3 partMomB(mXB,mYB,mZB);
0230 
0231             double dEta = partMom.PseudoRapidity() - partMomB.PseudoRapidity();
0232             double dPhi = TVector2::Phi_mpi_pi(partMom.Phi() - partMomB.Phi());
0233             double dR = TMath::Sqrt(dEta*dEta + dPhi*dPhi);
0234 
0235             recoChargedJetPartPairwiseDeltaRHist->Fill(dR);
0236               }
0237           }
0238           }
0239         numRecoChargedJetPartsHist->Fill(partsEnd[i] - partsBegin[i]);
0240         if(noElectron) numRecoChargedJetPartsNoElecHist->Fill(partsEnd[i] - partsBegin[i]);
0241       }
0242 
0243     // No Electrons
0244     if(noElectron)
0245       {
0246         recoChargedJetENoElecHist->Fill(recoNRG[i]);
0247         if(ECut) recoChargedJetEtaECutNoElecHist->Fill(jetMom.PseudoRapidity());
0248         recoChargedJetEvsEtaNoElecHist->Fill(jetMom.PseudoRapidity(),recoNRG[i]);
0249         if(ECut) recoChargedJetPhiVsEtaECutNoElecHist->Fill(jetMom.PseudoRapidity(),jetMom.Phi());
0250 
0251         if(ECut) numRecoChargedJetsNoElec++; 
0252       }
0253       }
0254     numRecoChargedJetsECutHist->Fill(numRecoChargedJets);
0255     numRecoChargedJetsECutNoElecHist->Fill(numRecoChargedJetsNoElec);
0256 
0257     //////////////////////////////////////////////////////////////////////////
0258     ////////////////////////  Analyze Generator Jets  ////////////////////////
0259     //////////////////////////////////////////////////////////////////////////
0260     int numGenChargedJets = 0;
0261     int numGenChargedJetsNoElec = 0;
0262     for(unsigned int i=0; i<genType.GetSize(); i++)
0263       {
0264     TVector3 jetMom(genMomX[i],genMomY[i],genMomZ[i]);
0265 
0266     counter->Fill(4);
0267 
0268     // Place eta cut to avoid edges of tracking acceptance
0269     if(TMath::Abs(jetMom.PseudoRapidity()) > 2.5) continue;
0270 
0271     // Place a minimum energy condition for several plots
0272     bool ECut = genNRG[i] > 5.0;
0273 
0274     if(ECut) numGenChargedJets++; 
0275 
0276     genChargedJetEHist->Fill(genNRG[i]);
0277     if(ECut) genChargedJetEtaECutHist->Fill(jetMom.PseudoRapidity());
0278     genChargedJetEvsEtaHist->Fill(jetMom.PseudoRapidity(),genNRG[i]);
0279     if(ECut) genChargedJetPhiVsEtaECutHist->Fill(jetMom.PseudoRapidity(),jetMom.Phi());
0280 
0281     bool noElectron = true;
0282     if(ECut)
0283       {
0284         // Find Jets with Electrons
0285         for(unsigned int m=genPartsBegin[i]; m<genPartsEnd[i]; m++)
0286           {
0287         if(mcM[genPartIndex[m]] > 0.00050 && mcM[genPartIndex[m]] < 0.00052)
0288           noElectron = false;
0289           }
0290 
0291         for(unsigned int j=genPartsBegin[i]; j<genPartsEnd[i]; j++)
0292           {
0293         // partsbegin and partsEnd specify the entries from _ReconstructedChargedJets_particles.index that make up the jet
0294         // _ReconstructedChargedJets_particles.index stores the ReconstructedChargedParticles index of the jet constituent
0295         double mX = mcMomX[genPartIndex[j]];
0296         double mY = mcMomY[genPartIndex[j]];
0297         double mZ = mcMomZ[genPartIndex[j]];
0298         double mM = mcM[genPartIndex[j]];
0299         double tmpE = TMath::Sqrt(mX*mX + mY*mY + mZ*mZ + mM*mM);
0300         
0301         TVector3 partMom(mX,mY,mZ);
0302         
0303         genChargedJetPartEHist->Fill(tmpE);
0304         genChargedJetPartEtaHist->Fill(partMom.PseudoRapidity());
0305         genChargedJetPartEvsEtaHist->Fill(partMom.PseudoRapidity(),tmpE);
0306         genChargedJetPartPhiVsEtaHist->Fill(partMom.PseudoRapidity(),partMom.Phi());
0307 
0308         if(noElectron)
0309           {
0310             genChargedJetPartENoElecHist->Fill(tmpE);
0311             genChargedJetPartEtaNoElecHist->Fill(partMom.PseudoRapidity());
0312             genChargedJetPartEvsEtaNoElecHist->Fill(partMom.PseudoRapidity(),tmpE);
0313             genChargedJetPartPhiVsEtaNoElecHist->Fill(partMom.PseudoRapidity(),partMom.Phi());
0314           }
0315 
0316         // Pairwise Distance Between Constituents
0317         if(j<(genPartsEnd[i]-1))
0318           {
0319             for(unsigned int k=j+1; k<genPartsEnd[i]; k++)
0320               {
0321             double mXB = mcMomX[genPartIndex[k]];
0322             double mYB = mcMomY[genPartIndex[k]];
0323             double mZB = mcMomZ[genPartIndex[k]];
0324 
0325             TVector3 partMomB(mXB,mYB,mZB);
0326 
0327             double dEta = partMom.PseudoRapidity() - partMomB.PseudoRapidity();
0328             double dPhi = TVector2::Phi_mpi_pi(partMom.Phi() - partMomB.Phi());
0329             double dR = TMath::Sqrt(dEta*dEta + dPhi*dPhi);
0330 
0331             genChargedJetPartPairwiseDeltaRHist->Fill(dR);
0332               }
0333           }
0334           }
0335         numGenChargedJetPartsHist->Fill(genPartsEnd[i] - genPartsBegin[i]);
0336         if(noElectron) numGenChargedJetPartsNoElecHist->Fill(genPartsEnd[i] - genPartsBegin[i]);
0337       }
0338 
0339     // No Electrons
0340     if(noElectron)
0341       {
0342         genChargedJetENoElecHist->Fill(genNRG[i]);
0343         if(ECut) genChargedJetEtaECutNoElecHist->Fill(jetMom.PseudoRapidity());
0344         genChargedJetEvsEtaNoElecHist->Fill(jetMom.PseudoRapidity(),genNRG[i]);
0345         if(ECut) genChargedJetPhiVsEtaECutNoElecHist->Fill(jetMom.PseudoRapidity(),jetMom.Phi());
0346 
0347         if(ECut) numGenChargedJetsNoElec++; 
0348       }
0349       }
0350     numGenChargedJetsECutHist->Fill(numGenChargedJets);
0351     numGenChargedJetsECutNoElecHist->Fill(numGenChargedJetsNoElec);
0352 
0353     
0354     //////////////////////////////////////////////////////////////////////////
0355     /////////////////////////////  Matched Jets  /////////////////////////////
0356     //////////////////////////////////////////////////////////////////////////
0357     for(unsigned int i=0; i<genType.GetSize(); i++)
0358       {
0359     TVector3 jetMom(genMomX[i],genMomY[i],genMomZ[i]);
0360 
0361     // Place eta cut to avoid edges of tracking acceptance
0362     //if(TMath::Abs(jetMom.PseudoRapidity()) > 2.5) continue;
0363 
0364     // Place a minimum energy condition
0365     //if(genNRG[i] < 5.0) continue;
0366     
0367     // Don't Look at Electron Jets
0368     bool hasElectron = false;
0369     // Find Jets with Electrons
0370     for(unsigned int m=genPartsBegin[i]; m<genPartsEnd[i]; m++)
0371       {
0372         if(mcM[genPartIndex[m]] > 0.00050 && mcM[genPartIndex[m]] < 0.00052)
0373           hasElectron = true;
0374       }
0375     //if(hasElectron) continue;
0376 
0377     // Find Matching Reconstructed Jet
0378     double minDeltaR = 999.;
0379     int minIndex = -1;
0380     for(unsigned int j=0; j<recoType.GetSize(); j++)
0381       {
0382         TVector3 recoMom(recoMomX[j],recoMomY[j],recoMomZ[j]);
0383 
0384         double dEta = jetMom.PseudoRapidity() - recoMom.PseudoRapidity();
0385         double dPhi = TVector2::Phi_mpi_pi(jetMom.Phi() - recoMom.Phi());
0386         double dR = TMath::Sqrt(dEta*dEta + dPhi*dPhi);
0387 
0388         if(dR < minDeltaR)
0389           {
0390         minDeltaR = dR;
0391         minIndex = j;
0392           }
0393       }
0394 
0395     // Do Backwards Match
0396     double minDeltaRBack = 999.;
0397     double minIndexBack = -1;
0398     if(minIndex > -1)
0399       {
0400         TVector3 recoMatchMom(recoMomX[minIndex],recoMomY[minIndex],recoMomZ[minIndex]);
0401         for(unsigned int j=0; j<genType.GetSize(); j++)
0402           {
0403         TVector3 genMom(genMomX[j],genMomY[j],genMomZ[j]);
0404         
0405         double dEta = recoMatchMom.PseudoRapidity() - genMom.PseudoRapidity();
0406         double dPhi = TVector2::Phi_mpi_pi(recoMatchMom.Phi() - genMom.Phi());
0407         double dR = TMath::Sqrt(dEta*dEta + dPhi*dPhi);
0408         
0409         if(dR < minDeltaRBack)
0410           {
0411             minDeltaRBack = dR;
0412             minIndexBack = j;
0413           }
0414           }
0415       }
0416 
0417     // Look at Best Match
0418     if(genNRG[i] > 5.0 && TMath::Abs(jetMom.PseudoRapidity()) < 2.5 && minIndex > -1 && !hasElectron) matchJetDeltaRHist->Fill(minDeltaR);
0419     if(genNRG[i] > 5.0 && TMath::Abs(jetMom.PseudoRapidity()) < 2.5 && minIndex > -1) matchJetDeltaRBackHist->Fill(minDeltaR);
0420     if(minIndex > -1 && genNRG[i] > 5.0 && TMath::Abs(jetMom.PseudoRapidity()) < 2.5 && !hasElectron)
0421       {
0422         TVector3 recoMatchMom(recoMomX[minIndex],recoMomY[minIndex],recoMomZ[minIndex]);
0423 
0424         recoVsGenChargedJetENoDRHist->Fill(genNRG[i],recoNRG[minIndex]);
0425 
0426         if(minDeltaR < DELTARCUT)
0427           {
0428         recoVsGenChargedJetEtaHist->Fill(jetMom.PseudoRapidity(),recoMatchMom.PseudoRapidity());
0429         recoVsGenChargedJetPhiHist->Fill(jetMom.Phi(),recoMatchMom.Phi());
0430         recoVsGenChargedJetEHist->Fill(genNRG[i],recoNRG[minIndex]);
0431         
0432         double jetERes = (recoNRG[minIndex] - genNRG[i])/genNRG[i];
0433         
0434         jetResVsEtaHist->Fill(jetMom.PseudoRapidity(),jetERes);
0435         jetResVsEHist->Fill(genNRG[i],jetERes);
0436         if(jetMom.PseudoRapidity() > -2.5 && jetMom.PseudoRapidity() < -1.0)
0437           jetResVsENegEtaHist->Fill(genNRG[i],jetERes);
0438         if(jetMom.PseudoRapidity() > -1.0 && jetMom.PseudoRapidity() < 1.0)
0439           jetResVsEMidEtaHist->Fill(genNRG[i],jetERes);
0440         if(jetMom.PseudoRapidity() > 1.0 && jetMom.PseudoRapidity() < 2.5)
0441           jetResVsEPosEtaHist->Fill(genNRG[i],jetERes);
0442         
0443         // Check for Duplicate Tracks
0444         bool noDuplicate = true;
0445         for(unsigned int j=partsBegin[minIndex]; j<partsEnd[minIndex]; j++)
0446           {
0447             double mX = recoPartMomX[recoPartIndex[j]];
0448             double mY = recoPartMomY[recoPartIndex[j]];
0449             double mZ = recoPartMomZ[recoPartIndex[j]];
0450             double mM = recoPartM[recoPartIndex[j]];
0451             double tmpE = TMath::Sqrt(mX*mX + mY*mY + mZ*mZ + mM*mM);
0452             
0453             TVector3 partMom(mX,mY,mZ);
0454             
0455             // Pairwise Distance Between Constituents
0456             if(j<(partsEnd[minIndex]-1))
0457               {
0458             for(unsigned int k=j+1; k<partsEnd[minIndex]; k++)
0459               {
0460                 double mXB = recoPartMomX[recoPartIndex[k]];
0461                 double mYB = recoPartMomY[recoPartIndex[k]];
0462                 double mZB = recoPartMomZ[recoPartIndex[k]];
0463                 
0464                 TVector3 partMomB(mXB,mYB,mZB);
0465                 
0466                 double dEta = partMom.PseudoRapidity() - partMomB.PseudoRapidity();
0467                 double dPhi = TVector2::Phi_mpi_pi(partMom.Phi() - partMomB.Phi());
0468                 double dR = TMath::Sqrt(dEta*dEta + dPhi*dPhi);
0469                 
0470                 if(dR < 0.02) noDuplicate = false;
0471               }
0472               }
0473           }
0474 
0475         if(noDuplicate)
0476           {
0477             recoVsGenChargedJetENoDupHist->Fill(genNRG[i],recoNRG[minIndex]);
0478 
0479             if(jetMom.PseudoRapidity() > -2.5 && jetMom.PseudoRapidity() < -1.0)
0480               jetResVsENegEtaNoDupHist->Fill(genNRG[i],jetERes);
0481             if(jetMom.PseudoRapidity() > -1.0 && jetMom.PseudoRapidity() < 1.0)
0482               jetResVsEMidEtaNoDupHist->Fill(genNRG[i],jetERes);
0483             if(jetMom.PseudoRapidity() > 1.0 && jetMom.PseudoRapidity() < 2.5)
0484               jetResVsEPosEtaNoDupHist->Fill(genNRG[i],jetERes);
0485           }
0486           }
0487       }
0488       }
0489 
0490     NEVENTS++;
0491   }
0492 
0493   
0494   gStyle->SetOptStat(0);
0495   ////////////////////////  Reconstructed Jets Plots  ////////////////////////
0496   // Reco Number
0497   TCanvas *c1 = new TCanvas("c1","Number Reco Jets",800,600);
0498   c1->Clear();
0499   c1->Divide(1,1);
0500 
0501   c1->cd(1);
0502   numRecoChargedJetsECutHist->Draw("HIST");
0503   numRecoChargedJetsECutNoElecHist->SetLineColor(seabornRed);
0504   numRecoChargedJetsECutNoElecHist->Draw("HISTSAME");
0505   numRecoChargedJetsECutHist->SetLineWidth(2); // Set line width to 2 (adjust as needed)
0506 numRecoChargedJetsECutNoElecHist->SetLineWidth(2);
0507   numRecoChargedJetsECutHist->SetTitle("Reconstructed Jets per Event (|eta| < 2.5 && E > 5);Number");
0508 
0509 TLegend *legend1 = new TLegend(0.7, 0.7, 0.9, 0.9); // Adjust the coordinates as needed
0510 legend1->AddEntry(numRecoChargedJetsECutHist, "With Electrons", "l");
0511 legend1->AddEntry(numRecoChargedJetsECutNoElecHist, "No Electrons", "l");
0512 legend1->Draw();
0513 
0514 
0515 
0516   gPad->SetLogy();
0517   if(PRINT) c1->Print(outputDir+"numberRecoJets.png"); // Number of reconstructed jets per event with energy > 5 GeV and Abs(eta) < 2.5
0518    delete c1;
0519 
0520   // Reco Energy
0521   TCanvas *c2 = new TCanvas("c2","Reco Jet Energy",800,600);
0522   c2->Clear();
0523   c2->Divide(1,1);
0524 
0525   c2->cd(1);
0526   recoChargedJetEHist->Draw("HIST");
0527   recoChargedJetENoElecHist->SetLineColor(seabornRed);
0528   recoChargedJetENoElecHist->Draw("HISTSAME");
0529 
0530   recoChargedJetEHist->SetLineWidth(2);
0531   recoChargedJetENoElecHist->SetLineWidth(2);
0532   recoChargedJetEHist->SetTitle("Reconstructed Jet Energy (|eta| < 2.5);Energy [GeV]");
0533 
0534 TLegend *legend2 = new TLegend(0.7, 0.7, 0.9, 0.9); // Adjust the coordinates as needed
0535 legend2->AddEntry(recoChargedJetEHist, "With Electrons", "l");
0536 legend2->AddEntry(recoChargedJetENoElecHist, "No Electrons", "l");
0537 legend2->Draw();
0538 
0539   gPad->SetLogy();
0540   if(PRINT) c2->Print(outputDir+"recoJetEnergy.png"); // Energy spectrum of reconstructed jets with Abs(eta) < 2.5
0541 
0542     delete c2;
0543   // Reco Eta
0544   TCanvas *c3 = new TCanvas("c3","Reco Jet Eta",800,600);
0545   c3->Clear();
0546   c3->Divide(1,1);
0547 
0548   c3->cd(1);
0549   recoChargedJetEtaECutHist->Draw("HIST");
0550   recoChargedJetEtaECutNoElecHist->SetLineColor(seabornRed);
0551   recoChargedJetEtaECutNoElecHist->Draw("HISTSAME");
0552 
0553   recoChargedJetEtaECutHist->SetLineWidth(2);
0554   recoChargedJetEtaECutNoElecHist->SetLineWidth(2);
0555   recoChargedJetEtaECutHist->SetTitle("Reconstructed Jet Eta (E > 5);Eta");
0556 
0557 //add legend 
0558 TLegend *legend3 = new TLegend(0.7, 0.7, 0.9, 0.9); // Adjust the coordinates as needed
0559 legend3->AddEntry(recoChargedJetEtaECutHist, "With Electrons", "l");
0560 legend3->AddEntry(recoChargedJetEtaECutNoElecHist, "No Electrons", "l");
0561 legend3->Draw();
0562 
0563   gPad->SetLogy();
0564   if(PRINT) c3->Print(outputDir+"recoJetEta.png"); // Eta spectrum of reconstructed jets with energy > 5 GeV
0565     delete c3;
0566   // Reco E Vs Eta
0567   TCanvas *c4 = new TCanvas("c4","Reco Jet E Vs Eta",800,600);
0568   c4->Clear();
0569   c4->Divide(1,1);
0570 
0571   c4->cd(1);
0572   recoChargedJetEvsEtaHist->Draw("COLZ");
0573   recoChargedJetEvsEtaHist->SetTitle("Reconstructed Jet Energy Vs Eta;Eta;Energy [GeV]");
0574   gPad->SetLogz();
0575   if(PRINT) c4->Print(outputDir+"recoJetEnergyvsEta.png"); // Energy vs eta of reconstructed jets
0576 
0577   // Reco Phi Vs Eta
0578   TCanvas *c5 = new TCanvas("c5","Reco Jet Phi Vs Eta",800,600);
0579   c5->Clear();
0580   c5->Divide(1,1);
0581 
0582   c5->cd(1);
0583   recoChargedJetPhiVsEtaECutHist->Draw("COLZ");
0584   recoChargedJetPhiVsEtaECutHist->SetTitle("Reconstructed Jet Phi Vs Eta (E > 5);Eta;Phi");
0585   gPad->SetLogz();
0586   if(PRINT) c5->Print(outputDir+"recoJetPhivsEta.png"); // Phi vs eta of reconstructed jets
0587 
0588   // Num Particles Per Reco Jet
0589   TCanvas *c6 = new TCanvas("c6","Number Constituents Per Reco Jet",800,600);
0590   c6->Clear();
0591   c6->Divide(1,1);
0592 
0593   c6->cd(1);
0594   numRecoChargedJetPartsHist->Draw("HIST");
0595   numRecoChargedJetPartsNoElecHist->SetLineColor(seabornRed);
0596   numRecoChargedJetPartsNoElecHist->Draw("HISTSAME");
0597 
0598   numRecoChargedJetPartsHist->SetLineWidth(2);
0599   numRecoChargedJetPartsNoElecHist->SetLineWidth(2);
0600   numRecoChargedJetPartsHist->SetTitle("Number of Constituents Per Reco Jet;Number of Constituents");
0601 
0602 TLegend *legend6 = new TLegend(0.7, 0.7, 0.9, 0.9); // Adjust the coordinates as needed
0603 legend6->AddEntry(numRecoChargedJetPartsHist, "With Electrons", "l");
0604 legend6->AddEntry(numRecoChargedJetPartsNoElecHist, "No Electrons", "l");
0605 legend6->Draw();
0606 
0607   gPad->SetLogy();
0608   if(PRINT) c6->Print(outputDir+"numConstituentsPerRecoJet.png"); // Number of constituents in reconstructed jets
0609 
0610   // Reco Part Energy
0611   TCanvas *c7 = new TCanvas("c7","Reco Jet Constituent Energy",800,600);
0612   c7->Clear();
0613   c7->Divide(1,1);
0614 
0615   c7->cd(1);
0616   recoChargedJetPartEHist->Draw("HIST");
0617   recoChargedJetPartENoElecHist->SetLineColor(seabornRed);
0618   recoChargedJetPartENoElecHist->Draw("HISTSAME");
0619 
0620   recoChargedJetPartEHist->SetLineWidth(2);
0621   recoChargedJetPartENoElecHist->SetLineWidth(2);
0622   recoChargedJetPartEHist->SetTitle("Reconstructed Jet Constituent Energy;Energy [GeV]");
0623 
0624   TLegend *legend7 = new TLegend(0.7, 0.7, 0.9, 0.9); // Adjust the coordinates as needed
0625   legend7->AddEntry(recoChargedJetPartEHist, "With Electrons", "l");
0626   legend7->AddEntry(recoChargedJetPartENoElecHist, "No Electrons", "l");
0627   legend7->Draw();
0628 
0629   gPad->SetLogy();
0630   if(PRINT) c7->Print(outputDir+"recoJetConstituentEnergy.png"); // Energy of reconstructed jet constituents
0631 
0632   // Reco Part Eta
0633   TCanvas *c8 = new TCanvas("c8","Reco Jet Constituent Eta",800,600);
0634   c8->Clear();
0635   c8->Divide(1,1);
0636 
0637   c8->cd(1);
0638   recoChargedJetPartEtaHist->Draw("HIST");
0639   recoChargedJetPartEtaNoElecHist->SetLineColor(seabornRed);
0640   recoChargedJetPartEtaNoElecHist->Draw("HISTSAME");
0641 
0642   recoChargedJetPartEtaHist->SetLineWidth(2);
0643   recoChargedJetPartEtaNoElecHist->SetLineWidth(2);
0644 
0645   recoChargedJetPartEtaHist->SetTitle("Reconstructed Jet Constituent Eta;Eta");
0646 
0647   TLegend *legend8 = new TLegend(0.7, 0.7, 0.9, 0.9); // Adjust the coordinates as needed
0648   legend8->AddEntry(recoChargedJetPartEtaHist, "With Electrons", "l");
0649   legend8->AddEntry(recoChargedJetPartEtaNoElecHist, "No Electrons", "l");
0650   legend8->Draw();
0651 
0652   gPad->SetLogy();
0653   if(PRINT) c8->Print(outputDir+"recoJetConstituentEta.png"); // Eta of reconstructed jet constituents
0654 
0655   // Reco Part E Vs Eta
0656   TCanvas *c9 = new TCanvas("c9","Reco Jet Constituent E Vs Eta",800,600);
0657   c9->Clear();
0658   c9->Divide(1,1);
0659 
0660   c9->cd(1);
0661   recoChargedJetPartEvsEtaHist->Draw("COLZ");
0662   recoChargedJetPartEvsEtaHist->SetTitle("Reconstructed Jet Constituent Energy Vs Eta;Eta;Energy [GeV]");
0663   gPad->SetLogz();
0664   if(PRINT) c9->Print(outputDir+"recoJetConstituentEnergyVsEta.png"); // Energy vs eta of reconstructed jet constituents
0665 
0666   // Reco Part Phi Vs Eta
0667   TCanvas *c10 = new TCanvas("c10","Reco Jet Constituent Phi Vs Eta",800,600);
0668   c10->Clear();
0669   c10->Divide(1,1);
0670 
0671   c10->cd(1);
0672   recoChargedJetPartPhiVsEtaHist->Draw("COLZ");
0673   recoChargedJetPartPhiVsEtaHist->SetTitle("Reconstructed Jet Constituent Phi Vs Eta;Eta;Phi");
0674   gPad->SetLogz();
0675   if(PRINT) c10->Print(outputDir+"recoJetConstituentPhiVsEta.png"); // Phi vs eta of reconstructed jet constituents
0676 
0677   // Reco Constituent Pairwise delta R
0678   TCanvas *c11 = new TCanvas("c11","Reco Jet Constituent Pairwise Delta R",800,600);
0679   c11->Clear();
0680   c11->Divide(1,1);
0681 
0682   c11->cd(1);
0683   recoChargedJetPartPairwiseDeltaRHist->Draw("COLZ");
0684   recoChargedJetPartPairwiseDeltaRHist->SetTitle("Pairwise Constituent Delta R;Delta R");
0685   recoChargedJetPartPairwiseDeltaRHist->GetXaxis()->SetRangeUser(0,0.5);
0686   gPad->SetLogy();
0687   if(PRINT) c11->Print(outputDir+"recoJetConstituentPairwiseDR.png"); // Distance between each pair of constituents in reconstructed jets
0688 
0689   // Reco E Vs Eta No Electron Jets
0690   TCanvas *c12 = new TCanvas("c12","Reco Jet E Vs Eta (No Electrons)",800,600);
0691   c12->Clear();
0692   c12->Divide(1,1);
0693 
0694   c12->cd(1);
0695   recoChargedJetEvsEtaNoElecHist->Draw("COLZ");
0696   recoChargedJetEvsEtaNoElecHist->SetTitle("Reconstructed Jet Energy Vs Eta (No Electrons);Eta;Energy [GeV]");
0697   gPad->SetLogz();
0698   if(PRINT) c12->Print(outputDir+"recoJetEnergyVsEtaNoElectron.png"); // Reconstructed jet energy - no jets containing electrons included
0699 
0700   // Reco Phi Vs Eta No Electron Jets
0701   TCanvas *c13 = new TCanvas("c13","Reco Jet Phi Vs Eta (No Electrons)",800,600);
0702   c13->Clear();
0703   c13->Divide(1,1);
0704 
0705   c13->cd(1);
0706   recoChargedJetPhiVsEtaECutNoElecHist->Draw("COLZ");
0707   recoChargedJetPhiVsEtaECutNoElecHist->SetTitle("Reconstructed Jet Phi Vs Eta (E > 5) (No Electrons);Eta;Phi");
0708   gPad->SetLogz();
0709   if(PRINT) c13->Print(outputDir+"recoJetPhiVsEtaNoElectron.png"); // Reconstructed Jet phi vs eta - no jets containing electrons included
0710 
0711   // Reco Part E Vs Eta No Electron Jets
0712   TCanvas *c14 = new TCanvas("c14","Reco Jet Constituent E Vs Eta (No Electrons)",800,600);
0713   c14->Clear();
0714   c14->Divide(1,1);
0715 
0716   c14->cd(1);
0717   recoChargedJetPartEvsEtaNoElecHist->Draw("COLZ");
0718   recoChargedJetPartEvsEtaNoElecHist->SetTitle("Reconstructed Jet Constituent Energy Vs Eta (No Electrons);Eta;Energy [GeV]");
0719   gPad->SetLogz();
0720   if(PRINT) c14->Print(outputDir+"recoJetConstituentEnergyVsEtaNoElectron.png"); // Reconstructed jet constituent energy vs eta - no jets containing electrons included
0721 
0722   // Reco Part Phi Vs Eta No Electron Jets
0723   TCanvas *c15 = new TCanvas("c15","Reco Jet Constituent Phi Vs Eta (No Electrons)",800,600);
0724   c15->Clear();
0725   c15->Divide(1,1);
0726 
0727   c15->cd(1);
0728   recoChargedJetPartPhiVsEtaNoElecHist->Draw("COLZ");
0729   recoChargedJetPartPhiVsEtaNoElecHist->SetTitle("Reconstructed Jet Constituent Phi Vs Eta (No Electrons);Eta;Phi");
0730   gPad->SetLogz();
0731   if(PRINT) c15->Print(outputDir+"recoJetConstituentPhiVsEtaNoElectron.png"); // Reconstructed jet constituent phi vs eta - no jets containing electrons included
0732 
0733   
0734   ////////////////////////  Generated Jets Plots  ////////////////////////
0735   // Gen Number
0736   TCanvas *c16 = new TCanvas("c16","Number Gen Jets",800,600);
0737   c16->Clear();
0738   c16->Divide(1,1);
0739 
0740   c16->cd(1);
0741   numGenChargedJetsECutHist->Draw("HIST");
0742   numGenChargedJetsECutNoElecHist->SetLineColor(seabornRed);
0743   numGenChargedJetsECutNoElecHist->Draw("HISTSAME");
0744 
0745   numGenChargedJetsECutHist->SetLineWidth(2);
0746   numGenChargedJetsECutNoElecHist->SetLineWidth(2);
0747 
0748   numGenChargedJetsECutHist->SetTitle("Generator Jets per Event (|eta| < 2.5 && E > 5);Number");
0749 
0750   TLegend *legend16 = new TLegend(0.7, 0.7, 0.9, 0.9); // Adjust the coordinates as needed
0751   legend16->AddEntry(numGenChargedJetsECutHist, "With Electrons", "l");
0752   legend16->AddEntry(numGenChargedJetsECutNoElecHist, "No Electrons", "l");
0753   legend16->Draw();
0754 
0755   gPad->SetLogy();
0756   if(PRINT) c16->Print(outputDir+"numberGenJets.png"); // Number of generator jets per event with energy > 5 GeV and Abs(eta) < 2.5
0757 
0758   // Gen Energy
0759   TCanvas *c17 = new TCanvas("c17","Gen Jet Energy",800,600);
0760   c17->Clear();
0761   c17->Divide(1,1);
0762 
0763   c17->cd(1);
0764   genChargedJetEHist->Draw("HIST");
0765   genChargedJetENoElecHist->SetLineColor(seabornRed);
0766   genChargedJetENoElecHist->Draw("HISTSAME");
0767 
0768   genChargedJetEHist->SetLineWidth(2);
0769   genChargedJetENoElecHist->SetLineWidth(2);
0770 
0771   genChargedJetEHist->SetTitle("Generator Jet Energy (|eta| < 2.5);Energy [GeV]");
0772   TLegend *legend17 = new TLegend(0.7, 0.7, 0.9, 0.9); // Adjust the coordinates as needed
0773   legend17->AddEntry(genChargedJetEHist, "With Electrons", "l");
0774   legend17->AddEntry(genChargedJetENoElecHist, "No Electrons", "l");
0775   legend17->Draw();
0776 
0777   gPad->SetLogy();
0778   if(PRINT) c17->Print(outputDir+"genJetEnergy.png"); // Energy spectrum of generated jets with Abs(eta) < 2.5
0779 
0780   // Gen Eta
0781   TCanvas *c18 = new TCanvas("c18","Gen Jet Eta",800,600);
0782   c18->Clear();
0783   c18->Divide(1,1);
0784 
0785   c18->cd(1);
0786   genChargedJetEtaECutHist->Draw("HIST");
0787   genChargedJetEtaECutNoElecHist->SetLineColor(seabornRed);
0788   genChargedJetEtaECutNoElecHist->Draw("HISTSAME");
0789 
0790   genChargedJetEtaECutHist->SetLineWidth(2);
0791   genChargedJetEtaECutNoElecHist->SetLineWidth(2);
0792 
0793   genChargedJetEtaECutHist->SetTitle("Generator Jet Eta (E > 5);Eta");
0794 
0795   TLegend *legend18 = new TLegend(0.7, 0.7, 0.9, 0.9); // Adjust the coordinates as needed
0796   legend18->AddEntry(genChargedJetEtaECutHist, "With Electrons", "l");
0797   legend18->AddEntry(genChargedJetEtaECutNoElecHist, "No Electrons", "l");
0798   legend18->Draw();
0799 
0800   gPad->SetLogy();
0801   if(PRINT) c18->Print(outputDir+"genJetEta.png"); // Eta spectrum of generator jets with energy > 5 GeV
0802 
0803   // Gen E Vs Eta
0804   TCanvas *c19 = new TCanvas("c19","Gen Jet E Vs Eta",800,600);
0805   c19->Clear();
0806   c19->Divide(1,1);
0807 
0808   c19->cd(1);
0809   genChargedJetEvsEtaHist->Draw("COLZ");
0810   genChargedJetEvsEtaHist->SetTitle("Generator Jet Energy Vs Eta;Eta;Energy [GeV]");
0811   gPad->SetLogz();
0812   if(PRINT) c19->Print(outputDir+"genJetEnergyvsEta.png"); // Energy vs eta of generator jets
0813 
0814   // Gen Phi Vs Eta
0815   TCanvas *c20 = new TCanvas("c20","Gen Jet Phi Vs Eta",800,600);
0816   c20->Clear();
0817   c20->Divide(1,1);
0818 
0819   c20->cd(1);
0820   genChargedJetPhiVsEtaECutHist->Draw("COLZ");
0821   genChargedJetPhiVsEtaECutHist->SetTitle("Generator Jet Phi Vs Eta (E > 5);Eta;Phi");
0822   gPad->SetLogz();
0823   if(PRINT) c20->Print(outputDir+"genJetPhiVsEta.png"); // Phi vs eta of generator jets
0824 
0825   // Num Particles Per Gen Jet
0826   TCanvas *c21 = new TCanvas("c21","Number Constituents Per Gen Jet",800,600);
0827   c21->Clear();
0828   c21->Divide(1,1);
0829 
0830   c21->cd(1);
0831   numGenChargedJetPartsHist->Draw("HIST");
0832   numGenChargedJetPartsNoElecHist->SetLineColor(seabornRed);
0833   numGenChargedJetPartsNoElecHist->Draw("HISTSAME");
0834 
0835   numGenChargedJetPartsHist->SetLineWidth(2);
0836   numGenChargedJetPartsNoElecHist->SetLineWidth(2);
0837 
0838   numGenChargedJetPartsHist->SetTitle("Number of Constituents Per Gen Jet;Number of Constituents");
0839 
0840   TLegend *legend21 = new TLegend(0.7, 0.7, 0.9, 0.9); // Adjust the coordinates as needed
0841   legend21->AddEntry(numGenChargedJetPartsHist, "With Electrons", "l");
0842   legend21->AddEntry(numGenChargedJetPartsNoElecHist, "No Electrons", "l");
0843   legend21->Draw();
0844   gPad->SetLogy();
0845   if(PRINT) c21->Print(outputDir+"numConstituentsPerGenJet.png"); // Number of constituents in generator jets
0846 
0847   // Gen Part Energy
0848   TCanvas *c22 = new TCanvas("c22","Gen Jet Constituent Energy",800,600);
0849   c22->Clear();
0850   c22->Divide(1,1);
0851 
0852   c22->cd(1);
0853   genChargedJetPartEHist->Draw("HIST");
0854   genChargedJetPartENoElecHist->SetLineColor(seabornRed);
0855   genChargedJetPartENoElecHist->Draw("HISTSAME");
0856 
0857   genChargedJetPartEHist->SetLineWidth(2);
0858   genChargedJetPartENoElecHist->SetLineWidth(2);
0859 
0860   genChargedJetPartEHist->SetTitle("Generator Jet Constituent Energy;Energy [GeV]");
0861 
0862   TLegend *legend22 = new TLegend(0.7, 0.7, 0.9, 0.9); // Adjust the coordinates as needed
0863   legend22->AddEntry(genChargedJetPartEHist, "With Electrons", "l");
0864   legend22->AddEntry(genChargedJetPartENoElecHist, "No Electrons", "l");
0865   legend22->Draw();
0866 
0867   gPad->SetLogy();
0868   if(PRINT) c22->Print(outputDir+"genJetConstituentEnergy.png"); // Energy of generator jet constituents
0869 
0870   // Gen Part Eta
0871   TCanvas *c23 = new TCanvas("c23","Gen Jet Constituent Eta",800,600);
0872   c23->Clear();
0873   c23->Divide(1,1);
0874 
0875   c23->cd(1);
0876   genChargedJetPartEtaHist->Draw("HIST");
0877   genChargedJetPartEtaNoElecHist->SetLineColor(seabornRed);
0878   genChargedJetPartEtaNoElecHist->Draw("HISTSAME");
0879 
0880   genChargedJetPartEtaHist->SetLineWidth(2);
0881   genChargedJetPartEtaNoElecHist->SetLineWidth(2);
0882 
0883   genChargedJetPartEtaHist->SetTitle("Generator Jet Constituent Eta;Eta");
0884 
0885   TLegend *legend23 = new TLegend(0.7, 0.7, 0.9, 0.9); // Adjust the coordinates as needed
0886   legend23->AddEntry(genChargedJetPartEtaHist, "With Electrons", "l");
0887   legend23->AddEntry(genChargedJetPartEtaNoElecHist, "No Electrons", "l");
0888   legend23->Draw();
0889 
0890   gPad->SetLogy();
0891   if(PRINT) c23->Print(outputDir+"genJetConstituentEta.png"); // Eta of generator jet constituents
0892 
0893   // Gen Part E Vs Eta
0894   TCanvas *c24 = new TCanvas("c24","Gen Jet Constituent E Vs Eta",800,600);
0895   c24->Clear();
0896   c24->Divide(1,1);
0897 
0898   c24->cd(1);
0899   genChargedJetPartEvsEtaHist->Draw("COLZ");
0900   genChargedJetPartEvsEtaHist->SetTitle("Generator Jet Constituent Energy Vs Eta;Eta;Energy [GeV]");
0901   gPad->SetLogz();
0902   if(PRINT) c24->Print(outputDir+"genJetConstituentEnergyVsEta.png"); // Energy vs eta of generator jet constituents
0903 
0904   // Gen Part Phi Vs Eta
0905   TCanvas *c25 = new TCanvas("c25","Gen Jet Constituent Phi Vs Eta",800,600);
0906   c25->Clear();
0907   c25->Divide(1,1);
0908 
0909   c25->cd(1);
0910   genChargedJetPartPhiVsEtaHist->Draw("COLZ");
0911   genChargedJetPartPhiVsEtaHist->SetTitle("Generator Jet Constituent Phi Vs Eta;Eta;Phi");
0912   gPad->SetLogz();
0913   if(PRINT) c25->Print(outputDir+"genJetConstituentPhiVsEta.png"); // Phi vs eta of generator jet constituents
0914 
0915   // Gen Constituent Pairwise delta R
0916   TCanvas *c26 = new TCanvas("c26","Gen Jet Constituent Pairwise Delta R",800,600);
0917   c26->Clear();
0918   c26->Divide(1,1);
0919 
0920   c26->cd(1);
0921   genChargedJetPartPairwiseDeltaRHist->Draw("COLZ");
0922   genChargedJetPartPairwiseDeltaRHist->SetTitle("Generator Jet Pairwise Constituent Delta R;Delta R");
0923   genChargedJetPartPairwiseDeltaRHist->GetXaxis()->SetRangeUser(0,0.5);
0924   gPad->SetLogy();
0925   if(PRINT) c26->Print(outputDir+"genJetConstituentPairwiseDR.png"); // Distance between each pair of constituents in generator jets
0926 
0927   // Gen E Vs Eta No Electron Jets
0928   TCanvas *c27 = new TCanvas("c27","Gen Jet E Vs Eta (No Electrons)",800,600);
0929   c27->Clear();
0930   c27->Divide(1,1);
0931 
0932   c27->cd(1);
0933   genChargedJetEvsEtaNoElecHist->Draw("COLZ");
0934   genChargedJetEvsEtaNoElecHist->SetTitle("Generator Jet Energy Vs Eta (No Electrons);Eta;Energy [GeV]");
0935   gPad->SetLogz();
0936   if(PRINT) c27->Print(outputDir+"genJetEnergyVsEtaNoElectron.png"); // Generator jet energy vs eta - no jets containing electrons included
0937 
0938   // Gen Phi Vs Eta No Electron Jets
0939   TCanvas *c28 = new TCanvas("c28","Gen Jet Phi Vs Eta (No Electrons)",800,600);
0940   c28->Clear();
0941   c28->Divide(1,1);
0942 
0943   c28->cd(1);
0944   genChargedJetPhiVsEtaECutNoElecHist->Draw("COLZ");
0945   genChargedJetPhiVsEtaECutNoElecHist->SetTitle("Generator Jet Phi Vs Eta (E > 5) (No Electrons);Eta;Phi");
0946   gPad->SetLogz();
0947   if(PRINT) c28->Print(outputDir+"genJetPhiVsEtaNoElectron.png"); // Generator Jet phi vs eta - no jets containing electrons included
0948 
0949   // Gen Part E Vs Eta No Electron Jets
0950   TCanvas *c29 = new TCanvas("c29","Gen Jet Constituent E Vs Eta (No Electrons)",800,600);
0951   c29->Clear();
0952   c29->Divide(1,1);
0953 
0954   c29->cd(1);
0955   genChargedJetPartEvsEtaNoElecHist->Draw("COLZ");
0956   genChargedJetPartEvsEtaNoElecHist->SetTitle("Generator Jet Constituent Energy Vs Eta (No Electrons);Eta;Energy [GeV]");
0957   gPad->SetLogz();
0958   if(PRINT) c29->Print(outputDir+"genJetConstituentEnergyVsEtaNoElectron.png"); // Generator jet constituent energy vs eta - no jets containing electrons included
0959 
0960   // Gen Part Phi Vs Eta No Electron Jets
0961   TCanvas *c30 = new TCanvas("c30","Gen Jet Constituent Phi Vs Eta (No Electrons)",800,600);
0962   c30->Clear();
0963   c30->Divide(1,1);
0964 
0965   c30->cd(1);
0966   genChargedJetPartPhiVsEtaNoElecHist->Draw("COLZ");
0967   genChargedJetPartPhiVsEtaNoElecHist->SetTitle("Generator Jet Constituent Phi Vs Eta (No Electrons);Eta;Phi");
0968   gPad->SetLogz();
0969   //c30->Print(outputDir+"recoJetEvsEta.png");
0970   if(PRINT) c30->Print(outputDir+"genJetConstituentPhiVsEtaNoElectron.png"); // Generator jet constituent phi vs eta - no jets containing electrons included
0971 
0972   
0973   ////////////////////////  Matched Jets Plots  ////////////////////////
0974   // Matched Delta R
0975   TCanvas *c31 = new TCanvas("c31","Gen - Reco Delta R",800,600);
0976   c31->Clear();
0977   c31->Divide(1,1);
0978 
0979   c31->cd(1);
0980   matchJetDeltaRHist->Draw("HIST");
0981   matchJetDeltaRBackHist->SetLineColor(seabornRed);
0982   //matchJetDeltaRBackHist->Draw("HISTSAME");
0983   matchJetDeltaRHist->SetTitle("Matched Gen - Reco Jet Delta R;Delta R");
0984   gPad->SetLogy();
0985   if(PRINT) c31->Print(outputDir+"genRecoJetDeltaR.png"); // Distance between closest generated and reconstructed jet pair
0986 
0987   // Matched Reco Vs Gen Eta
0988   TCanvas *c32 = new TCanvas("c32","Reco Vs Gen Eta",800,600);
0989   c32->Clear();
0990   c32->Divide(1,1);
0991 
0992   c32->cd(1);
0993   recoVsGenChargedJetEtaHist->Draw("COLZ");
0994   recoVsGenChargedJetEtaHist->SetTitle("Reconstructed Vs Generator Jet Eta;Gen Eta;Reco Eta");
0995   gPad->SetLogz();
0996   if(PRINT) c32->Print(outputDir+"matchedRecoVsGenJetEta.png"); // Matched Reconstructed Vs Generator Jet eta
0997 
0998   // Matched Reco Vs Gen Phi
0999   TCanvas *c33 = new TCanvas("c33","Reco Vs Gen Phi",800,600);
1000   c33->Clear();
1001   c33->Divide(1,1);
1002 
1003   c33->cd(1);
1004   recoVsGenChargedJetPhiHist->Draw("COLZ");
1005   recoVsGenChargedJetPhiHist->SetTitle("Reconstructed Vs Generator Jet Phi;Gen Phi;Reco Phi");
1006   gPad->SetLogz();
1007   if(PRINT) c33->Print(outputDir+"matchedRecoVsGenJetPhi.png"); // Matched reconstructed vs generator jet phi
1008 
1009   // Matched Reco Vs Gen Energy
1010   TCanvas *c34 = new TCanvas("c34","Reco Vs Gen Energy",800,600);
1011   c34->Clear();
1012   c34->Divide(1,1);
1013 
1014   TF1 *f1_34 = new TF1("f1_34","1.0*x + 0.0",1,100);
1015   TF1 *f2_34 = new TF1("f2_34","2.0*x + 0.0",1,100);
1016   TF1 *f3_34 = new TF1("f3_34","3.0*x + 0.0",1,100);
1017 
1018   c34->cd(1);
1019   recoVsGenChargedJetEHist->Draw("COLZ");
1020   recoVsGenChargedJetEHist->SetTitle("Reconstructed Vs Generator Jet Energy;Gen E;Reco E");
1021   f1_34->Draw("SAME");
1022   f2_34->Draw("SAME");
1023   f3_34->Draw("SAME");
1024   gPad->SetLogz();
1025   if(PRINT) c34->Print(outputDir+"matchedRecoVsGenJetEnergy.png"); // Matched reconstructed vs generator jet energy
1026 
1027   // Jet Res Vs Gen Eta
1028   TCanvas *c35 = new TCanvas("c35","Jet Res Vs Gen Eta",800,600);
1029   c35->Clear();
1030   c35->Divide(1,1);
1031 
1032   c35->cd(1);
1033   jetResVsEtaHist->Draw("COLZ");
1034   jetResVsEtaHist->SetTitle("(Reco - Gen)/Gen Jet Energy Vs Gen Eta;Gen Eta;Res");
1035   gPad->SetLogz();
1036   if(PRINT) c35->Print(outputDir+"matchedJetResolutionVsEta.png"); // Matched jet resolution vs generator jet eta
1037 
1038   // Jet Res Vs Gen E
1039   TCanvas *c36 = new TCanvas("c36","Jet Res Vs Gen E",800,600);
1040   c36->Clear();
1041   c36->Divide(1,1);
1042 
1043   c36->cd(1);
1044   jetResVsEHist->Draw("COLZ");
1045   jetResVsEHist->SetTitle("(Reco - Gen)/Gen Jet Energy Vs Gen Energy;Gen E;Res");
1046   gPad->SetLogz();
1047   if(PRINT) c36->Print(outputDir+"matchedJetResolutionVsEnergy.png"); // Matched jet resolution vs generator jet energy
1048 
1049   // Jet Res Vs Gen E Neg Eta
1050   TCanvas *c37 = new TCanvas("c37","Jet Res Vs Gen E (-2.5 < eta < -1.0)",800,600);
1051   c37->Clear();
1052   c37->Divide(1,1);
1053 
1054   c37->cd(1);
1055   jetResVsENegEtaNoDupHist->Draw("COLZ");
1056   jetResVsENegEtaNoDupHist->SetTitle("(Reco - Gen)/Gen Jet Energy Vs Gen Energy (-2.5 < eta < -1.0) No Duplicate;Gen E;Res");
1057   gPad->SetLogz();
1058   if(PRINT) c37->Print(outputDir+"matchedJetResolutionVsEnergyNegEta.png"); // Matched jet resolution vs generator jet energy -2.5 < eta < -1.0
1059 
1060   // Jet Res Vs Gen E Mid Eta
1061   TCanvas *c38 = new TCanvas("c38","Jet Res Vs Gen E (-1.0 < eta < 1.0)",800,600);
1062   c38->Clear();
1063   c38->Divide(1,1);
1064 
1065   c38->cd(1);
1066   jetResVsEMidEtaNoDupHist->Draw("COLZ");
1067   jetResVsEMidEtaNoDupHist->SetTitle("(Reco - Gen)/Gen Jet Energy Vs Gen Energy (-1.0 < eta < 1.0) No Duplicate;Gen E;Res");
1068   gPad->SetLogz();
1069   if(PRINT) c38->Print(outputDir+"matchedJetResolutionVsEnergyMidEta.png"); // Matched jet resolution vs generator jet energy -1.0 < eta < 1.0
1070     delete c38;
1071   // Jet Res Vs Gen E Pos Eta
1072   TCanvas *c39 = new TCanvas("c39","Jet Res Vs Gen E (1.0 < eta < 2.5)",800,600);
1073   c39->Clear();
1074   c39->Divide(1,1);
1075 
1076   c39->cd(1);
1077   jetResVsEPosEtaNoDupHist->Draw("COLZ");
1078   jetResVsEPosEtaNoDupHist->SetTitle("(Reco - Gen)/Gen Jet Energy Vs Gen Energy (1.0 < eta < 2.5) No Duplicate;Gen E;Res");
1079   gPad->SetLogz();
1080   if(PRINT) c39->Print(outputDir+"matchedJetResolutionVsEnergyPosEta.png"); // Matched jet resolution vs generator jet energy 1.0 < eta < 2.5
1081   delete c39;
1082 
1083   
1084   // Generate Resolution Plots
1085   const int BINS = 20;
1086   double binCent[BINS];
1087   double jesVsENeg[BINS];
1088   double jesVsEMid[BINS];
1089   double jesVsEPos[BINS];
1090   double jerVsENeg[BINS];
1091   double jerVsEMid[BINS];
1092   double jerVsEPos[BINS];
1093 
1094   std::fill(std::begin(binCent), std::end(binCent), -999.);
1095   std::fill(std::begin(jesVsENeg), std::end(jesVsENeg), -999.);
1096   std::fill(std::begin(jesVsEMid), std::end(jesVsEMid), -999.);
1097   std::fill(std::begin(jesVsEPos), std::end(jesVsEPos), -999.);
1098   std::fill(std::begin(jerVsENeg), std::end(jerVsENeg), -999.);
1099   std::fill(std::begin(jerVsEMid), std::end(jerVsEMid), -999.);
1100   std::fill(std::begin(jerVsEPos), std::end(jerVsEPos), -999.);
1101 
1102   TH1D *pxA = jetResVsENegEtaNoDupHist->ProjectionX("pxA",1,10000);
1103   for(int i=0; i<BINS; i++)
1104     {
1105       binCent[i] = pxA->GetBinCenter(i+1);
1106     }
1107 
1108   TCanvas *c40 = new TCanvas("c40","Negative Rapidity Fit Results",800,600);
1109   c40->Clear();
1110   c40->Divide(5,4);
1111 
1112   TH1D *hA[20];
1113   for(int i=1; i<21; i++)
1114     {
1115       hA[i-1] = (TH1D *)jetResVsENegEtaNoDupHist->ProjectionY(Form("projYA_%d",i),i,i);
1116 
1117       TF1 *myGausA = new TF1("myGausA","gaus",-0.5,0.5);
1118       myGausA->SetParameters(hA[i-1]->GetMaximum(),0.0,0.01);
1119 
1120       c40->cd(i);
1121       //hA[i-1]->Draw("HIST");
1122       hA[i-1]->Fit("myGausA","B","",-0.5,0.5);
1123       hA[i-1]->GetXaxis()->SetRangeUser(-1,1);
1124       gPad->SetLogy();
1125 
1126       if(hA[i-1]->GetEntries() > 2)
1127     {
1128       auto fA = hA[i-1]->GetFunction("myGausA");
1129 
1130       if(fA->GetParError(2)/fA->GetParameter(2) < 0.5)
1131         {
1132           jesVsENeg[i-1] = fA->GetParameter(1);
1133           jerVsENeg[i-1] = fA->GetParameter(2);
1134         }
1135       //cout << fA->GetParameter(0) << " " << fA->GetParameter(1) << " " << fA->GetParameter(2) << endl;
1136       //cout << fA->GetParError(0) << " " << fA->GetParError(1) << " " << fA->GetParError(2) << endl;
1137     }
1138     }
1139   if(PRINT) c40->Print(outputDir+"matchedJetResolutionVsEnergyNegEtaFitSummary.png"); // Matched jet resolution vs generator jet energy -2.5 < eta < -1.0 fits
1140     delete c40;
1141 
1142   TCanvas *c41 = new TCanvas("c41","Mid Rapidity Fit Results",800,600);
1143   c41->Clear();
1144   c41->Divide(5,4);
1145 
1146   TH1D *hB[20];
1147   for(int i=1; i<21; i++)
1148     {
1149       hB[i-1] = (TH1D *)jetResVsEMidEtaNoDupHist->ProjectionY(Form("projYB_%d",i),i,i);
1150 
1151       TF1 *myGausB = new TF1("myGausB","gaus",-0.5,0.5);
1152       myGausB->SetParameters(hB[i-1]->GetMaximum(),0.0,0.01);
1153 
1154       c41->cd(i);
1155       //hA[i-1]->Draw("HIST");
1156       hB[i-1]->Fit("myGausB","B","",-0.5,0.5);
1157       hB[i-1]->GetXaxis()->SetRangeUser(-1,1);
1158       gPad->SetLogy();
1159 
1160       if(hB[i-1]->GetEntries() > 2)
1161     {
1162       auto fB = hB[i-1]->GetFunction("myGausB");
1163 
1164       if(fB->GetParError(2)/fB->GetParameter(2) < 0.5)
1165         {
1166           jesVsEMid[i-1] = fB->GetParameter(1);
1167           jerVsEMid[i-1] = fB->GetParameter(2);
1168         }
1169       //cout << fB->GetParameter(0) << " " << fB->GetParameter(1) << " " << fB->GetParameter(2) << endl;
1170       //cout << fB->GetParError(0) << " " << fB->GetParError(1) << " " << fB->GetParError(2) << endl;
1171     }
1172     }
1173   if(PRINT) c41->Print(outputDir+"matchedJetResolutionVsEnergyMidEtaFitSummary.png"); // Matched jet resolution vs generator jet energy -1.0 < eta < 1.0 fits
1174     delete c41;
1175 
1176   TCanvas *c42 = new TCanvas("c42","Positive Rapidity Fit Results",800,600);
1177   c42->Clear();
1178   c42->Divide(5,4);
1179 
1180   TH1D *hC[20];
1181   for(int i=1; i<21; i++)
1182     {
1183       hC[i-1] = (TH1D *)jetResVsEPosEtaNoDupHist->ProjectionY(Form("projYC_%d",i),i,i);
1184 
1185       TF1 *myGausC = new TF1("myGausC","gaus",-0.5,0.5);
1186       myGausC->SetParameters(hC[i-1]->GetMaximum(),0.0,0.01);
1187 
1188       c42->cd(i);
1189       //hA[i-1]->Draw("HIST");
1190       hC[i-1]->Fit("myGausC","B","",-0.5,0.5);
1191       hC[i-1]->GetXaxis()->SetRangeUser(-1,1);
1192       gPad->SetLogy();
1193 
1194       if(hC[i-1]->GetEntries() > 2)
1195     {
1196       auto fC = hC[i-1]->GetFunction("myGausC");
1197 
1198       if(fC->GetParError(2)/fC->GetParameter(2) < 0.5)
1199         {
1200           jesVsEPos[i-1] = fC->GetParameter(1);
1201           jerVsEPos[i-1] = fC->GetParameter(2);
1202         }
1203       //cout << fC->GetParameter(0) << " " << fC->GetParameter(1) << " " << fC->GetParameter(2) << endl;
1204       //cout << fC->GetParError(0) << " " << fC->GetParError(1) << " " << fC->GetParError(2) << endl;
1205     }
1206     }
1207   if(PRINT) c42->Print(outputDir+"matchedJetResolutionVsEnergyPosEtaFitSummary.png"); // Matched jet resolution vs generator jet energy 1.0 < eta < 2.5 fits
1208     delete c42;
1209   TCanvas *c43 = new TCanvas("c43","Positive JES/JER",800,600);
1210   c43->Clear();
1211   c43->Divide(1,1);
1212 
1213   TGraph *gJESvsENeg = new TGraph(BINS,binCent,jesVsENeg);
1214   TGraph *gJERvsENeg = new TGraph(BINS,binCent,jerVsENeg);
1215 
1216   TGraph *gJESvsEMid = new TGraph(BINS,binCent,jesVsEMid);
1217   TGraph *gJERvsEMid = new TGraph(BINS,binCent,jerVsEMid);
1218 
1219   TGraph *gJESvsEPos = new TGraph(BINS,binCent,jesVsEPos);
1220   TGraph *gJERvsEPos = new TGraph(BINS,binCent,jerVsEPos);
1221 
1222   TH2D *test43 = new TH2D("test43","Jet Energy Scale / Resolution Vs Eta;True Eta;JES/JER",1,0.,100.,1,-0.2,0.2);
1223   test43->Draw();
1224 
1225   c43->cd(1);
1226   gJERvsENeg->Draw("*");
1227   gJERvsENeg->SetMarkerStyle(21);
1228   gJERvsENeg->SetMarkerSize(1);
1229   gJERvsENeg->SetMarkerColor(seabornBlue);
1230 
1231   gJESvsENeg->Draw("*");
1232   gJESvsENeg->SetMarkerStyle(26);
1233   gJESvsENeg->SetMarkerSize(1);
1234   gJESvsENeg->SetMarkerColor(seabornBlue);
1235 
1236   gJERvsEMid->Draw("*");
1237   gJERvsEMid->SetMarkerStyle(21);
1238   gJERvsEMid->SetMarkerSize(1);
1239   gJERvsEMid->SetMarkerColor(seabornRed);
1240 
1241   gJESvsEMid->Draw("*");
1242   gJESvsEMid->SetMarkerStyle(26);
1243   gJESvsEMid->SetMarkerSize(1);
1244   gJESvsEMid->SetMarkerColor(seabornRed);
1245 
1246   gJERvsEPos->Draw("*");
1247   gJERvsEPos->SetMarkerStyle(21);
1248   gJERvsEPos->SetMarkerSize(1);
1249   gJERvsEPos->SetMarkerColor(seabornGreen);
1250 
1251   gJESvsEPos->Draw("*");
1252   gJESvsEPos->SetMarkerStyle(26);
1253   gJESvsEPos->SetMarkerSize(1);
1254   gJESvsEPos->SetMarkerColor(seabornGreen);
1255 
1256   TLegend *legend = new TLegend(0.7,0.7,0.9,0.9); 
1257 legend->AddEntry(gJERvsENeg, "JER, (-2.5 < #eta < -1)", "p");
1258 legend->AddEntry(gJESvsENeg, "JES, (-2.5 < #eta < -1)","p");
1259 legend->AddEntry(gJERvsEMid, "JER, (-1 < #eta < 1)", "p");
1260 legend->AddEntry(gJESvsEMid, "JES, (-1 < #eta < 1)", "p");
1261 legend->AddEntry(gJERvsEPos, "JER, (1 < #eta < 2.5) ", "p");
1262 legend->AddEntry(gJESvsEPos, "JES, (1 < #eta < 2.5)", "p");
1263 legend->Draw();
1264 
1265   if(PRINT) c43->Print(outputDir+"matchedJetScaleResolutionSummary.png"); // Matched jet JER/JES summary
1266     delete c43;
1267 
1268 delete mychain;
1269 }
[ 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!