Reference Manual

Go to the documentation of this file.
00001 /*
00002  *  See header file for a description of this class.
00003  *
00004  *  $Date: 2011/08/12 15:37:12 $
00005  *  $Revision: 1.26 $
00006  *  \author F. Chlebana - Fermilab
00007  */
00008 
00009 #include "DQMOffline/JetMET/interface/JetAnalyzer.h"
00010 #include "DataFormats/Common/interface/Handle.h"
00011 
00012 #include "DataFormats/JetReco/interface/CaloJetCollection.h"
00013 #include "DataFormats/JetReco/interface/CaloJet.h"
00014 
00015 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00016 
00017 
00018 #include <string>
00019 using namespace edm;
00020 
00021 // ***********************************************************
00022 JetAnalyzer::JetAnalyzer(const edm::ParameterSet& pSet) {
00023   
00024   parameters   = pSet;
00025   _leadJetFlag = 0;
00026   _JetLoPass   = 0;
00027   _JetHiPass   = 0;
00028   _ptThreshold = 20.;
00029   _asymmetryThirdJetCut = 5.;
00030   _balanceThirdJetCut   = 0.2; 
00031   _n90HitsMin =0;
00032   _fHPDMax=1.;
00033   _resEMFMin=0.;
00034   _n90HitsMinLoose =0;
00035   _fHPDMaxLoose=1.;
00036   _resEMFMinLoose=0.;
00037   _n90HitsMinTight =0;
00038   _fHPDMaxTight=1.;
00039   _resEMFMinTight=0.;
00040   _sigmaEtaMinTight=-999.;
00041   _sigmaPhiMinTight=-999.;
00042 
00043 } 
00044   
00045 // ***********************************************************
00046 JetAnalyzer::~JetAnalyzer() { }
00047 
00048 
00049 // ***********************************************************
00050 void JetAnalyzer::beginJob(DQMStore * dbe) {
00051     
00052   jetname = "jetAnalyzer";
00053   
00054   LogTrace(jetname)<<"[JetAnalyzer] Parameters initialization";
00055   dbe->setCurrentFolder("JetMET/Jet/"+_source);
00056 
00057   jetME = dbe->book1D("jetReco", "jetReco", 3, 1, 4);
00058   jetME->setBinLabel(1,"CaloJets",1);
00059 
00060   //
00061   jetID = new reco::helper::JetIDHelper(parameters.getParameter<ParameterSet>("JetIDParams"));
00062   //
00063 
00064   fillJIDPassFrac = parameters.getParameter<int>("fillJIDPassFrac");
00065   makedijetselection = parameters.getParameter<int>("makedijetselection");
00066 
00067   // monitoring of eta parameter
00068   etaBin = parameters.getParameter<int>("etaBin");
00069   etaMin = parameters.getParameter<double>("etaMin");
00070   etaMax = parameters.getParameter<double>("etaMax");
00071 
00072   // monitoring of phi paramater
00073   phiBin = parameters.getParameter<int>("phiBin");
00074   phiMin = parameters.getParameter<double>("phiMin");
00075   phiMax = parameters.getParameter<double>("phiMax");
00076 
00077   // monitoring of the transverse momentum
00078   ptBin = parameters.getParameter<int>("ptBin");
00079   ptMin = parameters.getParameter<double>("ptMin");
00080   ptMax = parameters.getParameter<double>("ptMax");
00081 
00082   // 
00083   eBin = parameters.getParameter<int>("eBin");
00084   eMin = parameters.getParameter<double>("eMin");
00085   eMax = parameters.getParameter<double>("eMax");
00086 
00087   // 
00088   pBin = parameters.getParameter<int>("pBin");
00089   pMin = parameters.getParameter<double>("pMin");
00090   pMax = parameters.getParameter<double>("pMax");
00091 
00092   //
00093   _ptThreshold = parameters.getParameter<double>("ptThreshold");
00094   _asymmetryThirdJetCut = parameters.getParameter<double>("asymmetryThirdJetCut");
00095   _balanceThirdJetCut   = parameters.getParameter<double>("balanceThirdJetCut");
00096   _n90HitsMin = parameters.getParameter<int>("n90HitsMin");
00097   _fHPDMax = parameters.getParameter<double>("fHPDMax");
00098   _resEMFMin = parameters.getParameter<double>("resEMFMin");
00099   _sigmaEtaMinTight = parameters.getParameter<double>("sigmaEtaMinTight");
00100   _sigmaPhiMinTight = parameters.getParameter<double>("sigmaPhiMinTight");
00101 
00102   _n90HitsMinLoose = parameters.getParameter<int>("n90HitsMinLoose");
00103   _fHPDMaxLoose = parameters.getParameter<double>("fHPDMaxLoose");
00104   _resEMFMinLoose = parameters.getParameter<double>("resEMFMinLoose");
00105   _n90HitsMinTight = parameters.getParameter<int>("n90HitsMinTight");
00106   _fHPDMaxTight = parameters.getParameter<double>("fHPDMaxTight");
00107   _resEMFMinTight = parameters.getParameter<double>("resEMFMinTight");
00108 
00109   
00110   // Generic Jet Parameters
00111   mPt                      = dbe->book1D("Pt",  "Pt", ptBin, ptMin, ptMax);
00112   mEta                     = dbe->book1D("Eta", "Eta", etaBin, etaMin, etaMax);
00113   mPhi                     = dbe->book1D("Phi", "Phi", phiBin, phiMin, phiMax);
00114   mConstituents            = dbe->book1D("Constituents", "# of Constituents", 50, 0, 100);
00115   mHFrac                   = dbe->book1D("HFrac", "HFrac", 120, -0.1, 1.1);
00116   mEFrac                   = dbe->book1D("EFrac", "EFrac", 120, -0.1, 1.1);
00117   //
00118   //mE                       = dbe->book1D("E", "E", eBin, eMin, eMax);
00119   //mP                       = dbe->book1D("P", "P", pBin, pMin, pMax);
00120   //  mMass                    = dbe->book1D("Mass", "Mass", 100, 0, 25);
00121   //
00122   mPhiVSEta                     = dbe->book2D("PhiVSEta", "PhiVSEta", 50, etaMin, etaMax, 24, phiMin, phiMax);
00123   if(makedijetselection!=1){
00124     mPt_1                    = dbe->book1D("Pt1", "Pt1", 20, 0, 100);   
00125     mPt_2                    = dbe->book1D("Pt2", "Pt2", 60, 0, 300);   
00126     mPt_3                    = dbe->book1D("Pt3", "Pt3", 100, 0, 5000);
00127     // Low and high pt trigger paths
00128     mPt_Lo                  = dbe->book1D("Pt_Lo", "Pt (Pass Low Pt Jet Trigger)", 20, 0, 100);   
00129     //mEta_Lo                 = dbe->book1D("Eta_Lo", "Eta (Pass Low Pt Jet Trigger)", etaBin, etaMin, etaMax);
00130     mPhi_Lo                 = dbe->book1D("Phi_Lo", "Phi (Pass Low Pt Jet Trigger)", phiBin, phiMin, phiMax);
00131     
00132     mPt_Hi                  = dbe->book1D("Pt_Hi", "Pt (Pass Hi Pt Jet Trigger)", 60, 0, 300);   
00133     mEta_Hi                 = dbe->book1D("Eta_Hi", "Eta (Pass Hi Pt Jet Trigger)", etaBin, etaMin, etaMax);
00134     mPhi_Hi                 = dbe->book1D("Phi_Hi", "Phi (Pass Hi Pt Jet Trigger)", phiBin, phiMin, phiMax);
00135     mNJets                   = dbe->book1D("NJets", "Number of Jets", 100, 0, 100);
00136 
00137     //mPt_Barrel_Lo            = dbe->book1D("Pt_Barrel_Lo", "Pt Barrel (Pass Low Pt Jet Trigger)", 20, 0, 100);   
00138     //mPhi_Barrel_Lo           = dbe->book1D("Phi_Barrel_Lo", "Phi Barrel (Pass Low Pt Jet Trigger)", phiBin, phiMin, phiMax);
00139     mConstituents_Barrel     = dbe->book1D("Constituents_Barrel", "Constituents Barrel above", 50, 0, 100);
00140     mHFrac_Barrel            = dbe->book1D("HFrac_Barrel", "HFrac Barrel", 100, 0, 1);
00141     mEFrac_Barrel            = dbe->book1D("EFrac_Barrel", "EFrac Barrel", 110, -0.05, 1.05);
00142     
00143     //mPt_EndCap_Lo            = dbe->book1D("Pt_EndCap_Lo", "Pt EndCap (Pass Low Pt Jet Trigger)", 20, 0, 100);   
00144     //mPhi_EndCap_Lo           = dbe->book1D("Phi_EndCap_Lo", "Phi EndCap (Pass Low Pt Jet Trigger)", phiBin, phiMin, phiMax);
00145     mConstituents_EndCap     = dbe->book1D("Constituents_EndCap", "Constituents EndCap", 50, 0, 100);
00146     mHFrac_EndCap            = dbe->book1D("HFrac_Endcap", "HFrac EndCap", 100, 0, 1);
00147     mEFrac_EndCap            = dbe->book1D("EFrac_Endcap", "EFrac EndCap", 110, -0.05, 1.05);
00148     
00149     //mPt_Forward_Lo           = dbe->book1D("Pt_Forward_Lo", "Pt Forward (Pass Low Pt Jet Trigger)", 20, 0, 100);  
00150     //mPhi_Forward_Lo          = dbe->book1D("Phi_Forward_Lo", "Phi Forward (Pass Low Pt Jet Trigger)", phiBin, phiMin, phiMax);
00151     mConstituents_Forward    = dbe->book1D("Constituents_Forward", "Constituents Forward", 50, 0, 100);
00152     mHFrac_Forward           = dbe->book1D("HFrac_Forward", "HFrac Forward", 100, 0, 1);
00153     mEFrac_Forward           = dbe->book1D("EFrac_Forward", "EFrac Forward", 110, -0.05, 1.05);
00154     
00155     mPt_Barrel_Hi            = dbe->book1D("Pt_Barrel_Hi", "Pt Barrel (Pass Hi Pt Jet Trigger)", 60, 0, 300);   
00156     mPhi_Barrel_Hi           = dbe->book1D("Phi_Barrel_Hi", "Phi Barrel (Pass Hi Pt Jet Trigger)", phiBin, phiMin, phiMax);
00157     //mConstituents_Barrel_Hi  = dbe->book1D("Constituents_Barrel_Hi", "Constituents Barrel (Pass Hi Pt Jet Trigger)", 50, 0, 100);
00158     //mHFrac_Barrel_Hi         = dbe->book1D("HFrac_Barrel_Hi", "HFrac Barrel (Pass Hi Pt Jet Trigger)", 100, 0, 1);
00159     
00160     mPt_EndCap_Hi            = dbe->book1D("Pt_EndCap_Hi", "Pt EndCap (Pass Hi Pt Jet Trigger)", 60, 0, 300);  
00161     mPhi_EndCap_Hi           = dbe->book1D("Phi_EndCap_Hi", "Phi EndCap (Pass Hi Pt Jet Trigger)", phiBin, phiMin, phiMax);
00162     //mConstituents_EndCap_Hi  = dbe->book1D("Constituents_EndCap_Hi", "Constituents EndCap (Pass Hi Pt Jet Trigger)", 50, 0, 100);
00163     //mHFrac_EndCap_Hi         = dbe->book1D("HFrac_EndCap_Hi", "HFrac EndCap (Pass Hi Pt Jet Trigger)", 100, 0, 1);
00164     
00165     mPt_Forward_Hi           = dbe->book1D("Pt_Forward_Hi", "Pt Forward (Pass Hi Pt Jet Trigger)", 60, 0, 300);  
00166     mPhi_Forward_Hi          = dbe->book1D("Phi_Forward_Hi", "Phi Forward (Pass Hi Pt Jet Trigger)", phiBin, phiMin, phiMax);
00167     //mConstituents_Forward_Hi = dbe->book1D("Constituents_Forward_Hi", "Constituents Forward (Pass Hi Pt Jet Trigger)", 50, 0, 100);
00168     //mHFrac_Forward_Hi        = dbe->book1D("HFrac_Forward_Hi", "HFrac Forward (Pass Hi Pt Jet Trigger)", 100, 0, 1);
00169     
00170     mPhi_Barrel              = dbe->book1D("Phi_Barrel", "Phi_Barrel", phiBin, phiMin, phiMax);
00171     //mE_Barrel                = dbe->book1D("E_Barrel", "E_Barrel", eBin, eMin, eMax);
00172     mPt_Barrel               = dbe->book1D("Pt_Barrel", "Pt_Barrel", ptBin, ptMin, ptMax);
00173     
00174     mPhi_EndCap              = dbe->book1D("Phi_EndCap", "Phi_EndCap", phiBin, phiMin, phiMax);
00175     //mE_EndCap                = dbe->book1D("E_EndCap", "E_EndCap", eBin, eMin, 2*eMax);
00176     mPt_EndCap               = dbe->book1D("Pt_EndCap", "Pt_EndCap", ptBin, ptMin, ptMax);
00177     
00178     mPhi_Forward             = dbe->book1D("Phi_Forward", "Phi_Forward", phiBin, phiMin, phiMax);
00179     //mE_Forward               = dbe->book1D("E_Forward", "E_Forward", eBin, eMin, 4*eMax);
00180     mPt_Forward              = dbe->book1D("Pt_Forward", "Pt_Forward", ptBin, ptMin, ptMax);
00181     
00182     // Leading Jet Parameters
00183     mEtaFirst                = dbe->book1D("EtaFirst", "EtaFirst", 100, -5, 5);
00184     mPhiFirst                = dbe->book1D("PhiFirst", "PhiFirst", 70, -3.5, 3.5);
00185     //mEFirst                  = dbe->book1D("EFirst", "EFirst", 100, 0, 1000);
00186     mPtFirst                 = dbe->book1D("PtFirst", "PtFirst", 100, 0, 500);
00187     if(fillJIDPassFrac==1){//fillJIDPassFrac defines a collection of cleaned jets, for which we will want to fill the cleaning passing fraction
00188       mLooseJIDPassFractionVSeta      = dbe->bookProfile("LooseJIDPassFractionVSeta","LooseJIDPassFractionVSeta",etaBin, etaMin, etaMax,0.,1.2);
00189       mLooseJIDPassFractionVSpt       = dbe->bookProfile("LooseJIDPassFractionVSpt","LooseJIDPassFractionVSpt",ptBin, ptMin, ptMax,0.,1.2);
00190       mTightJIDPassFractionVSeta      = dbe->bookProfile("TightJIDPassFractionVSeta","TightJIDPassFractionVSeta",etaBin, etaMin, etaMax,0.,1.2);
00191       mTightJIDPassFractionVSpt       = dbe->bookProfile("TightJIDPassFractionVSpt","TightJIDPassFractionVSpt",ptBin, ptMin, ptMax,0.,1.2);
00192 
00193 
00194     }
00195   }
00196   // CaloJet specific
00197   mMaxEInEmTowers         = dbe->book1D("MaxEInEmTowers", "MaxEInEmTowers", 100, 0, 100);
00198   mMaxEInHadTowers        = dbe->book1D("MaxEInHadTowers", "MaxEInHadTowers", 100, 0, 100);
00199   if(makedijetselection!=1) {
00200     mHadEnergyInHO          = dbe->book1D("HadEnergyInHO", "HadEnergyInHO", 100, 0, 10);
00201     mHadEnergyInHB          = dbe->book1D("HadEnergyInHB", "HadEnergyInHB", 100, 0, 50);
00202     mHadEnergyInHF          = dbe->book1D("HadEnergyInHF", "HadEnergyInHF", 100, 0, 50);
00203     mHadEnergyInHE          = dbe->book1D("HadEnergyInHE", "HadEnergyInHE", 100, 0, 100);
00204     mEmEnergyInEB           = dbe->book1D("EmEnergyInEB", "EmEnergyInEB", 100, 0, 50);
00205     mEmEnergyInEE           = dbe->book1D("EmEnergyInEE", "EmEnergyInEE", 100, 0, 50);
00206     mEmEnergyInHF           = dbe->book1D("EmEnergyInHF", "EmEnergyInHF", 120, -20, 100);
00207   }
00208   mDPhi                   = dbe->book1D("DPhi", "dPhi btw the two leading jets", 100, 0., acos(-1.));
00209   
00210   //JetID variables
00211   
00212   mresEMF                 = dbe->book1D("resEMF", "resEMF", 50, 0., 1.);
00213   mN90Hits                = dbe->book1D("N90Hits", "N90Hits", 50, 0., 50);
00214   mfHPD                   = dbe->book1D("fHPD", "fHPD", 50, 0., 1.);
00215   mfRBX                   = dbe->book1D("fRBX", "fRBX", 50, 0., 1.);
00216 
00217   //  msigmaEta                   = dbe->book1D("sigmaEta", "sigmaEta", 50, 0., 1.);
00218   //  msigmaPhi                   = dbe->book1D("sigmaPhi", "sigmaPhi", 50, 0., 0.5);
00219   
00220   if(makedijetselection==1) {
00221     mDijetAsymmetry                   = dbe->book1D("DijetAsymmetry", "DijetAsymmetry", 100, -1., 1.);
00222     mDijetBalance                     = dbe->book1D("DijetBalance",   "DijetBalance",   100, -2., 2.);
00223     if (fillJIDPassFrac==1) {
00224       mLooseJIDPassFractionVSeta  = dbe->bookProfile("LooseJIDPassFractionVSeta","LooseJIDPassFractionVSeta",50, -3., 3.,0.,1.2);
00225       mLooseJIDPassFractionVSpt   = dbe->bookProfile("LooseJIDPassFractionVSpt","LooseJIDPassFractionVSpt",ptBin, ptMin, ptMax,0.,1.2);
00226       mTightJIDPassFractionVSeta  = dbe->bookProfile("TightJIDPassFractionVSeta","TightJIDPassFractionVSeta",50, -3., 3.,0.,1.2);
00227       mTightJIDPassFractionVSpt   = dbe->bookProfile("TightJIDPassFractionVSpt","TightJIDPassFractionVSpt",ptBin, ptMin, ptMax,0.,1.2);
00228     }
00229   }
00230 }
00231 
00232 
00233 //void JetAnalyzer::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup, 
00234 //                        const edm::TriggerResults& triggerResults,
00235 //                        const reco::CaloJet& jet) {
00236 
00237 
00238 // ***********************************************************
00239 void JetAnalyzer::endJob() {
00240   delete jetID;
00241 }
00242 
00243 
00244 // ***********************************************************
00245 void JetAnalyzer::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup, 
00246                           const reco::CaloJetCollection& caloJets) {
00247   int numofjets=0;
00248   double  fstPhi=0.;
00249   double  sndPhi=0.;
00250   double  diff = 0.;
00251   double  corr = 0.;
00252   double  dphi = -999. ;
00253   bool thiscleaned=false;
00254   bool Loosecleaned=false;
00255   bool Tightcleaned=false;
00256   bool thisemfclean=true;
00257   bool emfcleanLoose=true;
00258   bool emfcleanTight=true;
00259 
00260   srand( iEvent.id().event() % 10000);
00261   
00262   if(makedijetselection==1){
00263     //Dijet selection - careful: the pT is uncorrected!
00264     //if(makedijetselection==1 && caloJets.size()>=2){
00265     if(caloJets.size()>=2){
00266       double  dphiDJ = -999. ;
00267       bool emfcleanLooseFirstJet=true;
00268       bool emfcleanLooseSecondJet=true;
00269       bool emfcleanTightFirstJet=true;
00270       bool emfcleanTightSecondJet=true;
00271       bool LoosecleanedFirstJet = false;
00272       bool LoosecleanedSecondJet = false;
00273       bool TightcleanedFirstJet = false;
00274       bool TightcleanedSecondJet = false;
00275       //both jets pass pt threshold
00276       if ((caloJets.at(0)).pt() > _ptThreshold && (caloJets.at(1)).pt() > _ptThreshold ) {
00277         if(fabs((caloJets.at(0)).eta())<3. && fabs((caloJets.at(1)).eta())<3. ){
00278           //calculate dphi
00279           dphiDJ = fabs((caloJets.at(0)).phi()-(caloJets.at(1)).phi());
00280           if (dphiDJ > 3.14) dphiDJ=fabs(dphiDJ -6.28 );
00281           //fill DPhi histo (before cutting)
00282           if (mDPhi) mDPhi->Fill (dphiDJ);
00283           //dphi cut
00284           if(fabs(dphiDJ)>2.1){
00285             //JetID 
00286             emfcleanLooseFirstJet=true;
00287             emfcleanTightFirstJet=true;
00288             emfcleanLooseSecondJet=true;
00289             emfcleanTightSecondJet=true;
00290             //jetID for first jet
00291             jetID->calculate(iEvent, (caloJets.at(0)));
00292             if(jetID->restrictedEMF()<_resEMFMinLoose && fabs((caloJets.at(0)).eta())<2.6) emfcleanLooseFirstJet=false;
00293             if(jetID->restrictedEMF()<_resEMFMinTight && fabs((caloJets.at(0)).eta())<2.6) emfcleanTightFirstJet=false;
00294             if(jetID->n90Hits()>=_n90HitsMinLoose && jetID->fHPD()<_fHPDMaxLoose && emfcleanLooseFirstJet) LoosecleanedFirstJet=true;
00295             if(jetID->n90Hits()>=_n90HitsMinTight && jetID->fHPD()<_fHPDMaxTight && sqrt((caloJets.at(0)).etaetaMoment())>_sigmaEtaMinTight && sqrt((caloJets.at(0)).phiphiMoment())>_sigmaPhiMinTight && emfcleanTightFirstJet) TightcleanedFirstJet=true;
00296             //fill the JID variables histograms BEFORE you cut on them
00297             if (mN90Hits)         mN90Hits->Fill (jetID->n90Hits());
00298             if (mfHPD)            mfHPD->Fill (jetID->fHPD());
00299             if (mresEMF)         mresEMF->Fill (jetID->restrictedEMF());
00300             if (mfRBX)            mfRBX->Fill (jetID->fRBX());
00301 
00302             //jetID for second jet
00303             jetID->calculate(iEvent, (caloJets.at(1)));
00304             if(jetID->restrictedEMF()<_resEMFMinLoose && fabs((caloJets.at(1)).eta())<2.6) emfcleanLooseSecondJet=false;
00305             if(jetID->restrictedEMF()<_resEMFMinTight && fabs((caloJets.at(1)).eta())<2.6) emfcleanTightSecondJet=false;
00306             if(jetID->n90Hits()>=_n90HitsMinLoose && jetID->fHPD()<_fHPDMaxLoose && emfcleanLooseSecondJet) LoosecleanedSecondJet=true;
00307             if(jetID->n90Hits()>=_n90HitsMinTight && jetID->fHPD()<_fHPDMaxTight && sqrt((caloJets.at(1)).etaetaMoment())>_sigmaEtaMinTight && sqrt((caloJets.at(1)).phiphiMoment())>_sigmaPhiMinTight && emfcleanTightSecondJet) TightcleanedSecondJet=true;
00308             //fill the JID variables histograms BEFORE you cut on them
00309             if (mN90Hits)         mN90Hits->Fill (jetID->n90Hits());
00310             if (mfHPD)            mfHPD->Fill (jetID->fHPD());
00311             if (mresEMF)         mresEMF->Fill (jetID->restrictedEMF());
00312             if (mfRBX)            mfRBX->Fill (jetID->fRBX());
00313 
00314             if(LoosecleanedFirstJet && LoosecleanedSecondJet) { //only if both jets are (loose) cleaned
00315               //fill histos for first jet
00316               if (mPt)   mPt->Fill ((caloJets.at(0)).pt());
00317               if (mEta)  mEta->Fill ((caloJets.at(0)).eta());
00318               if (mPhi)  mPhi->Fill ((caloJets.at(0)).phi());
00319               if (mPhiVSEta) mPhiVSEta->Fill((caloJets.at(0)).eta(),(caloJets.at(0)).phi());
00320               if (mConstituents) mConstituents->Fill ((caloJets.at(0)).nConstituents());
00321               if (mHFrac)        mHFrac->Fill ((caloJets.at(0)).energyFractionHadronic());
00322               if (mEFrac)        mEFrac->Fill ((caloJets.at(0)).emEnergyFraction());
00323               //if (mE)    mE->Fill ((caloJets.at(0)).energy());
00324               //if (mP)    mP->Fill ((caloJets.at(0)).p());
00325               //if (mMass) mMass->Fill ((caloJets.at(0)).mass());
00326               if (mMaxEInEmTowers)  mMaxEInEmTowers->Fill ((caloJets.at(0)).maxEInEmTowers());
00327               if (mMaxEInHadTowers) mMaxEInHadTowers->Fill ((caloJets.at(0)).maxEInHadTowers());
00328               if (mN90Hits)         mN90Hits->Fill (jetID->n90Hits());
00329               if (mfHPD)            mfHPD->Fill (jetID->fHPD());
00330               if (mresEMF)         mresEMF->Fill (jetID->restrictedEMF());
00331               if (mfRBX)            mfRBX->Fill (jetID->fRBX());
00332               //sigmaeta and sigmaphi only used in the tight selection.
00333               //fill the histos for them AFTER the loose selection 
00334               //  if (msigmaEta)  msigmaEta->Fill(sqrt((caloJets.at(0)).etaetaMoment()));
00335               //  if (msigmaPhi)  msigmaPhi->Fill(sqrt((caloJets.at(0)).phiphiMoment()));
00336               //fill histos for second jet
00337               if (mPt)   mPt->Fill ((caloJets.at(1)).pt());
00338               if (mEta)  mEta->Fill ((caloJets.at(1)).eta());
00339               if (mPhi)  mPhi->Fill ((caloJets.at(1)).phi());
00340               if (mPhiVSEta) mPhiVSEta->Fill((caloJets.at(1)).eta(),(caloJets.at(1)).phi());
00341               if (mConstituents) mConstituents->Fill ((caloJets.at(1)).nConstituents());
00342               if (mHFrac)        mHFrac->Fill ((caloJets.at(1)).energyFractionHadronic());
00343               if (mEFrac)        mEFrac->Fill ((caloJets.at(1)).emEnergyFraction());
00344               //if (mE)    mE->Fill ((caloJets.at(1)).energy());
00345               //if (mP)    mP->Fill ((caloJets.at(1)).p());
00346               //if (mMass) mMass->Fill ((caloJets.at(1)).mass());
00347               if (mMaxEInEmTowers)  mMaxEInEmTowers->Fill ((caloJets.at(1)).maxEInEmTowers());
00348               if (mMaxEInHadTowers) mMaxEInHadTowers->Fill ((caloJets.at(1)).maxEInHadTowers());
00349               //sigmaeta and sigmaphi only used in the tight selection.
00350               //fill the histos for them AFTER the loose selection 
00351               //  if (msigmaEta)  msigmaEta->Fill(sqrt((caloJets.at(1)).etaetaMoment()));
00352               //  if (msigmaPhi)  msigmaPhi->Fill(sqrt((caloJets.at(1)).phiphiMoment()));
00353 
00354             }
00355             //let's see how many of these jets passed the JetID cleaning
00356             if(fillJIDPassFrac==1) {
00357               if(LoosecleanedFirstJet) {
00358                 mLooseJIDPassFractionVSeta->Fill((caloJets.at(0)).eta(),1.);
00359                 mLooseJIDPassFractionVSpt->Fill((caloJets.at(0)).pt(),1.);
00360               } else  {
00361                 mLooseJIDPassFractionVSeta->Fill((caloJets.at(0)).eta(),0.);
00362                 mLooseJIDPassFractionVSpt->Fill((caloJets.at(0)).pt(),0.);
00363               }
00364               if(LoosecleanedSecondJet) {
00365                 mLooseJIDPassFractionVSeta->Fill((caloJets.at(1)).eta(),1.);
00366                 mLooseJIDPassFractionVSpt->Fill((caloJets.at(1)).pt(),1.);
00367               } else  {
00368                 mLooseJIDPassFractionVSeta->Fill((caloJets.at(1)).eta(),0.);
00369                 mLooseJIDPassFractionVSpt->Fill((caloJets.at(1)).pt(),0.);
00370               }
00371               //TIGHT JID
00372               if(TightcleanedFirstJet) {
00373                 mTightJIDPassFractionVSeta->Fill((caloJets.at(0)).eta(),1.);
00374                 mTightJIDPassFractionVSpt->Fill((caloJets.at(0)).pt(),1.);
00375               } else  {
00376                 mTightJIDPassFractionVSeta->Fill((caloJets.at(0)).eta(),0.);
00377                 mTightJIDPassFractionVSpt->Fill((caloJets.at(0)).pt(),0.);
00378               }
00379               if(TightcleanedSecondJet) {
00380                 mTightJIDPassFractionVSeta->Fill((caloJets.at(1)).eta(),1.);
00381                 mTightJIDPassFractionVSpt->Fill((caloJets.at(1)).pt(),1.);
00382               } else  {
00383                 mTightJIDPassFractionVSeta->Fill((caloJets.at(1)).eta(),0.);
00384                 mTightJIDPassFractionVSpt->Fill((caloJets.at(1)).pt(),0.);
00385               }
00386 
00387             }//if fillJIDPassFrac
00388           }// FABS DPHI < 2.1
00389         }// fabs eta < 3
00390       }// pt jets > threshold
00391       //now do the dijet balance and asymmetry calculations
00392       if (fabs(caloJets.at(0).eta() < 1.4)) {
00393         double pt_dijet = (caloJets.at(0).pt() + caloJets.at(1).pt())/2;
00394         
00395         double dPhi = fabs((caloJets.at(0)).phi()-(caloJets.at(1)).phi());
00396         if (dPhi > 3.14) dPhi=fabs(dPhi -6.28 );
00397         
00398         if (dPhi > 2.7) {
00399           double pt_probe;
00400           double pt_barrel;
00401           int jet1, jet2;
00402 
00403           int randJet = rand() % 2;
00404 
00405           if (fabs(caloJets.at(1).eta() < 1.4)) {
00406             if (randJet) {
00407               jet1 = 0;
00408               jet2 = 1;
00409             }
00410             else {
00411               jet1 = 1;
00412               jet2 = 0;
00413             }
00414           
00415             /***Di-Jet Asymmetry****
00416              * leading jets eta < 1.4
00417              * leading jets dphi > 2.7
00418              * pt_third jet < threshold
00419              * A = (pt_1 - pt_2)/(pt_1 + pt_2)
00420              * jets 1 and two are randomly ordered
00421              */
00422             bool thirdJetCut = true;
00423             for (unsigned int third = 2; third < caloJets.size(); ++third) 
00424               if (caloJets.at(third).pt() > _asymmetryThirdJetCut) 
00425                 thirdJetCut = false;
00426             if (thirdJetCut) {
00427               double dijetAsymmetry = (caloJets.at(jet1).pt() - caloJets.at(jet2).pt()) / (caloJets.at(jet1).pt() + caloJets.at(jet2).pt());
00428               mDijetAsymmetry->Fill(dijetAsymmetry);
00429             }// end restriction on third jet pt in asymmetry calculation
00430               
00431           }
00432           else {
00433             jet1 = 0;
00434             jet2 = 1;
00435           }
00436           
00437           pt_barrel = caloJets.at(jet1).pt();
00438           pt_probe  = caloJets.at(jet2).pt();
00439           
00440           //dijet balance cuts
00441           /***Di-Jet Balance****
00442            * pt_dijet = (pt_probe+pt_barrel)/2
00443            * leading jets dphi > 2.7
00444            * reject evnets where pt_third/pt_dijet > 0.2
00445            * pv selection
00446            * B = (pt_probe - pt_barrel)/pt_dijet
00447            * select probe randomly from 2 jets if both leading jets are in the barrel
00448            */
00449           bool thirdJetCut = true;
00450           for (unsigned int third = 2; third < caloJets.size(); ++third) 
00451             if (caloJets.at(third).pt()/pt_dijet > _balanceThirdJetCut) 
00452               thirdJetCut = false;
00453           if (thirdJetCut) {
00454             double dijetBalance = (pt_probe - pt_barrel) / pt_dijet;
00455             mDijetBalance->Fill(dijetBalance);
00456           }// end restriction on third jet pt ratio in balance calculation
00457         }// dPhi > 2.7
00458       }// leading jet eta cut for asymmetry and balance calculations
00459     }//jet size >= 2
00460   }// do dijet selection
00461   else {
00462     for (reco::CaloJetCollection::const_iterator jet = caloJets.begin(); jet!=caloJets.end(); ++jet) {
00463       LogTrace(jetname)<<"[JetAnalyzer] Analyze Calo Jet";
00464       Loosecleaned=false;
00465       Tightcleaned=false;
00466       if (jet == caloJets.begin()) {
00467         fstPhi = jet->phi();
00468         _leadJetFlag = 1;
00469       } else {
00470         _leadJetFlag = 0;
00471       }
00472       if (jet == (caloJets.begin()+1)) sndPhi = jet->phi();
00473       //jetID
00474       jetID->calculate(iEvent, *jet);
00475       //minimal (uncorrected!) pT cut
00476       if (jet->pt() > _ptThreshold) {
00477         //  if (msigmaEta)  msigmaEta->Fill(sqrt(jet->etaetaMoment()));
00478         //  if (msigmaPhi)  msigmaPhi->Fill(sqrt(jet->phiphiMoment()));
00479         //cleaning to use for filling histograms
00480         thisemfclean=true;
00481         if(jetID->restrictedEMF()<_resEMFMin && fabs(jet->eta())<2.6) thisemfclean=false;
00482         if(jetID->n90Hits()>=_n90HitsMin && jetID->fHPD()<_fHPDMax && thisemfclean) thiscleaned=true;
00483         //loose and tight cleaning, used to fill the JetIDPAssFraction histos
00484         if(jetID->n90Hits()>=_n90HitsMinLoose && jetID->fHPD()<_fHPDMaxLoose && emfcleanLoose) Loosecleaned=true;
00485         if(jetID->n90Hits()>=_n90HitsMinTight && jetID->fHPD()<_fHPDMaxTight && sqrt(jet->etaetaMoment())>_sigmaEtaMinTight && sqrt(jet->phiphiMoment())>_sigmaPhiMinTight && emfcleanTight) Tightcleaned=true;
00486 
00487         if(fillJIDPassFrac==1) {
00488           if(Loosecleaned) {
00489             mLooseJIDPassFractionVSeta->Fill(jet->eta(),1.);
00490             mLooseJIDPassFractionVSpt->Fill(jet->pt(),1.);
00491           } else {
00492             mLooseJIDPassFractionVSeta->Fill(jet->eta(),0.);
00493             mLooseJIDPassFractionVSpt->Fill(jet->pt(),0.);
00494           }
00495           //TIGHT
00496           if(Tightcleaned) {
00497             mTightJIDPassFractionVSeta->Fill(jet->eta(),1.);
00498             mTightJIDPassFractionVSpt->Fill(jet->pt(),1.);
00499           } else {
00500             mTightJIDPassFractionVSeta->Fill(jet->eta(),0.);
00501             mTightJIDPassFractionVSpt->Fill(jet->pt(),0.);
00502           }
00503         }
00504         //eventually we could define the "cleaned" flag differently for e.g. HF
00505         if(thiscleaned) {
00506           numofjets++ ;
00507           jetME->Fill(1);      
00508           
00509           // Leading jet
00510           // Histograms are filled once per event
00511           if (_leadJetFlag == 1) { 
00512             if (mEtaFirst) mEtaFirst->Fill (jet->eta());
00513             if (mPhiFirst) mPhiFirst->Fill (jet->phi());
00514             //if (mEFirst)   mEFirst->Fill (jet->energy());
00515             if (mPtFirst)  mPtFirst->Fill (jet->pt());
00516           }
00517           // --- Passed the low pt jet trigger
00518           if (_JetLoPass == 1) {
00519           /*  if (fabs(jet->eta()) <= 1.3) {
00520               if (mPt_Barrel_Lo)           mPt_Barrel_Lo->Fill(jet->pt());
00521               if (mEta_Lo)          mEta_Lo->Fill(jet->eta());
00522               if (mPhi_Barrel_Lo)          mPhi_Barrel_Lo->Fill(jet->phi());
00523             }
00524             if ( (fabs(jet->eta()) > 1.3) && (fabs(jet->eta()) <= 3) ) {
00525               if (mPt_EndCap_Lo)           mPt_EndCap_Lo->Fill(jet->pt());
00526               if (mEta_Lo)          mEta_Lo->Fill(jet->eta());
00527               if (mPhi_EndCap_Lo)          mPhi_EndCap_Lo->Fill(jet->phi());
00528             }
00529             if (fabs(jet->eta()) > 3.0) {
00530               if (mPt_Forward_Lo)           mPt_Forward_Lo->Fill(jet->pt());
00531               if (mEta_Lo)          mEta_Lo->Fill(jet->eta());
00532               if (mPhi_Forward_Lo)          mPhi_Forward_Lo->Fill(jet->phi());
00533             } */
00534             //if (mEta_Lo) mEta_Lo->Fill (jet->eta());
00535             if (mPhi_Lo) mPhi_Lo->Fill (jet->phi());
00536             if (mPt_Lo)  mPt_Lo->Fill (jet->pt());
00537           }
00538           
00539           // --- Passed the high pt jet trigger
00540           if (_JetHiPass == 1) {
00541             if (fabs(jet->eta()) <= 1.3) {
00542               if (mPt_Barrel_Hi && jet->pt()>100.)           mPt_Barrel_Hi->Fill(jet->pt());
00543               if (mEta_Hi && jet->pt()>100.)          mEta_Hi->Fill(jet->eta());
00544               if (mPhi_Barrel_Hi)          mPhi_Barrel_Hi->Fill(jet->phi());
00545               //if (mConstituents_Barrel_Hi) mConstituents_Barrel_Hi->Fill(jet->nConstituents());       
00546               //if (mHFrac_Barrel_Hi)        mHFrac_Barrel_Hi->Fill(jet->energyFractionHadronic());     
00547             }
00548             if ( (fabs(jet->eta()) > 1.3) && (fabs(jet->eta()) <= 3) ) {
00549               if (mPt_EndCap_Hi && jet->pt()>100.)           mPt_EndCap_Hi->Fill(jet->pt());
00550               if (mEta_Hi && jet->pt()>100.)          mEta_Hi->Fill(jet->eta());
00551               if (mPhi_EndCap_Hi)          mPhi_EndCap_Hi->Fill(jet->phi());
00552               //if (mConstituents_EndCap_Hi) mConstituents_EndCap_Hi->Fill(jet->nConstituents());       
00553               //if (mHFrac_EndCap_Hi)        mHFrac_EndCap_Hi->Fill(jet->energyFractionHadronic());     
00554             }
00555             if (fabs(jet->eta()) > 3.0) {
00556               if (mPt_Forward_Hi && jet->pt()>100.)           mPt_Forward_Hi->Fill(jet->pt());
00557               if (mEta_Hi && jet->pt()>100.)          mEta_Hi->Fill(jet->eta());
00558               if (mPhi_Forward_Hi)          mPhi_Forward_Hi->Fill(jet->phi());
00559               //if (mConstituents_Forward_Hi) mConstituents_Forward_Hi->Fill(jet->nConstituents());     
00560               //if (mHFrac_Forward_Hi)        mHFrac_Forward_Hi->Fill(jet->energyFractionHadronic());   
00561             }
00562             
00563             if (mEta_Hi && jet->pt()>100.) mEta_Hi->Fill (jet->eta());
00564             if (mPhi_Hi) mPhi_Hi->Fill (jet->phi());
00565             if (mPt_Hi)  mPt_Hi->Fill (jet->pt());
00566           }
00567           
00568           if (mPt)   mPt->Fill (jet->pt());
00569           if (mPt_1) mPt_1->Fill (jet->pt());
00570           if (mPt_2) mPt_2->Fill (jet->pt());
00571           if (mPt_3) mPt_3->Fill (jet->pt());
00572           if (mEta)  mEta->Fill (jet->eta());
00573           if (mPhi)  mPhi->Fill (jet->phi());
00574           
00575           if (mPhiVSEta) mPhiVSEta->Fill(jet->eta(),jet->phi());
00576           
00577           if (mConstituents) mConstituents->Fill (jet->nConstituents());
00578           if (mHFrac)        mHFrac->Fill (jet->energyFractionHadronic());
00579           if (mEFrac)        mEFrac->Fill (jet->emEnergyFraction());
00580           
00581           if (fabs(jet->eta()) <= 1.3) {
00582             if (mPt_Barrel)   mPt_Barrel->Fill (jet->pt());
00583             if (mPhi_Barrel)  mPhi_Barrel->Fill (jet->phi());
00584             //if (mE_Barrel)    mE_Barrel->Fill (jet->energy());
00585       if (mConstituents_Barrel)    mConstituents_Barrel->Fill(jet->nConstituents());    
00586       if (mHFrac_Barrel)           mHFrac_Barrel->Fill(jet->energyFractionHadronic());  
00587       if (mEFrac_Barrel)           mEFrac_Barrel->Fill(jet->emEnergyFraction());        
00588           }
00589           if ( (fabs(jet->eta()) > 1.3) && (fabs(jet->eta()) <= 3) ) {
00590             if (mPt_EndCap)   mPt_EndCap->Fill (jet->pt());
00591             if (mPhi_EndCap)  mPhi_EndCap->Fill (jet->phi());
00592             //if (mE_EndCap)    mE_EndCap->Fill (jet->energy());
00593       if (mConstituents_EndCap)    mConstituents_EndCap->Fill(jet->nConstituents());    
00594       if (mHFrac_EndCap)           mHFrac_EndCap->Fill(jet->energyFractionHadronic());
00595       if (mEFrac_EndCap)           mEFrac_EndCap->Fill(jet->emEnergyFraction());        
00596           }
00597           if (fabs(jet->eta()) > 3.0) {
00598             if (mPt_Forward)   mPt_Forward->Fill (jet->pt());
00599             if (mPhi_Forward)  mPhi_Forward->Fill (jet->phi());
00600             //if (mE_Forward)    mE_Forward->Fill (jet->energy());
00601       if (mConstituents_Forward)    mConstituents_Forward->Fill(jet->nConstituents());  
00602       if (mHFrac_Forward)           mHFrac_Forward->Fill(jet->energyFractionHadronic());
00603       if (mEFrac_Forward)           mEFrac_Forward->Fill(jet->emEnergyFraction());      
00604           }
00605           
00606           //if (mE)    mE->Fill (jet->energy());
00607           //if (mP)    mP->Fill (jet->p());
00608           // if (mMass) mMass->Fill (jet->mass());
00609           
00610           if (mMaxEInEmTowers)  mMaxEInEmTowers->Fill (jet->maxEInEmTowers());
00611           if (mMaxEInHadTowers) mMaxEInHadTowers->Fill (jet->maxEInHadTowers());
00612           
00613           if (mHadEnergyInHO)   mHadEnergyInHO->Fill (jet->hadEnergyInHO());
00614           if (mHadEnergyInHB)   mHadEnergyInHB->Fill (jet->hadEnergyInHB());
00615           if (mHadEnergyInHF)   mHadEnergyInHF->Fill (jet->hadEnergyInHF());
00616           if (mHadEnergyInHE)   mHadEnergyInHE->Fill (jet->hadEnergyInHE());
00617           if (mEmEnergyInEB)    mEmEnergyInEB->Fill (jet->emEnergyInEB());
00618           if (mEmEnergyInEE)    mEmEnergyInEE->Fill (jet->emEnergyInEE());
00619           if (mEmEnergyInHF)    mEmEnergyInHF->Fill (jet->emEnergyInHF());
00620           
00621           if (mN90Hits)         mN90Hits->Fill (jetID->n90Hits());
00622           if (mfHPD)            mfHPD->Fill (jetID->fHPD());
00623           if (mresEMF)         mresEMF->Fill (jetID->restrictedEMF());
00624           if (mfRBX)            mfRBX->Fill (jetID->fRBX());
00625           
00626           //calculate correctly the dphi
00627           if(numofjets>1) {
00628             diff = fabs(fstPhi - sndPhi);
00629             corr = 2*acos(-1.) - diff;
00630             if(diff < acos(-1.)) { 
00631               dphi = diff; 
00632             } else { 
00633               dphi = corr;
00634             }
00635           }
00636         }
00637       }//pt cut
00638     }
00639     if (mNJets) mNJets->Fill (numofjets);
00640     if (mDPhi && dphi>-998.) mDPhi->Fill (dphi);
00641   }//not dijet
00642 }
00643
CMSSW_4_4_3_patch1/src/DQMOffline/JetMET/src/JetAnalyzer.cc
