/** A macro for applying a pt cut on jets and finding invariant mass. This examples creates two histograms. The first histogram shows the invariant mass of two jets, and the second histogram shows the pt of jets, after applying a cut. */ #include "DataFormats/FWLite/interface/Handle.h" #include "DataFormats/FWLite/interface/Event.h" #include "TH1.h" #include "TFile.h" #include "TCanvas.h" #include "TLegend.h" #include #if !defined(__CINT__) && !defined(__MAKECINT__) #include "DataFormats/PatCandidates/interface/Jet.h" #endif #include using namespace std; void example() { TFile * file = TFile::Open("file:SM.root"); // TFile * file = TFile::Open("file:BSM.root"); TH1D * hist_invarMass_SC5 = new TH1D("hist_invarMass_SC5", "Invariant Mass of the 4-vector sum of Two Jets", 100, 0, 1000); TH1D * hist_pt1_SC5 = new TH1D("hist_pt1_SC5", "Pt of highest Pt jet", 100, 0, 500); TH1D * hist_pt2_SC5 = new TH1D("hist_pt2_SC5", "Pt of 2nd highest Pt jet", 100, 0, 500); TH1D * hist_eta1_SC5 = new TH1D("hist_eta1_SC5", "Eta of highest Pt jet", 100, -5, 5); TH1D * hist_eta2_SC5 = new TH1D("hist_eta2_SC5", "Eta of 2nd highest Pt jet", 100, -5, 5); TH1D * hist_invarMass_SC7 = new TH1D("hist_invarMass_SC7", "Invariant Mass of the 4-vector sum of Two Jets", 100, 0, 1000); TH1D * hist_pt1_SC7 = new TH1D("hist_pt1_SC7", "Pt of highest Pt jet", 100, 0, 500); TH1D * hist_pt2_SC7 = new TH1D("hist_pt2_SC7", "Pt of 2nd highest Pt jet", 100, 0, 500); TH1D * hist_eta1_SC7 = new TH1D("hist_eta1_SC7", "Eta of highest Pt jet", 100, -5, 5); TH1D * hist_eta2_SC7 = new TH1D("hist_eta2_SC7", "Eta of 2nd highest Pt jet", 100, -5, 5); TH1D * hist_invarMass_KT4 = new TH1D("hist_invarMass_KT4", "Invariant Mass of the 4-vector sum of Two Jets", 100, 0, 1000); TH1D * hist_pt1_KT4 = new TH1D("hist_pt1_KT4", "Pt of highest Pt jet", 100, 0, 500); TH1D * hist_pt2_KT4 = new TH1D("hist_pt2_KT4", "Pt of 2nd highest Pt jet", 100, 0, 500); TH1D * hist_eta1_KT4 = new TH1D("hist_eta1_KT4", "Eta of highest Pt jet", 100, -5, 5); TH1D * hist_eta2_KT4= new TH1D("hist_eta2_KT4", "Eta of 2nd highest Pt jet", 100, -5, 5); TH1D * hist_invarMass_KT6 = new TH1D("hist_invarMass_KT6", "Invariant Mass of the 4-vector sum of Two Jets", 100, 0, 1000); TH1D * hist_pt1_KT6 = new TH1D("hist_pt1_KT6", "Pt of highest Pt jet", 100, 0, 500); TH1D * hist_pt2_KT6 = new TH1D("hist_pt2_KT6", "Pt of 2nd highest Pt jet", 100, 0, 500); TH1D * hist_eta1_KT6 = new TH1D("hist_eta1_KT6", "Eta of highest Pt jet", 100, -5, 5); TH1D * hist_eta2_KT6= new TH1D("hist_eta2_KT6", "Eta of 2nd highest Pt jet", 100, -5, 5); fwlite::Event ev(file); for( ev.toBegin(); ! ev.atEnd(); ++ev) { fwlite::Handle > h_jet_SC5; h_jet_SC5.getByLabel(ev,"selectedLayer1JetsSC5"); if ( !h_jet_SC5.isValid() ) continue; vector const & jets_SC5= *h_jet_SC5; fwlite::Handle > h_jet_SC7; h_jet_SC7.getByLabel(ev,"selectedLayer1JetsSC7"); if ( !h_jet_SC7.isValid() ) continue; vector const & jets_SC7= *h_jet_SC7; fwlite::Handle > h_jet_KT4; h_jet_KT4.getByLabel(ev,"selectedLayer1JetsKT4"); if ( !h_jet_KT4.isValid() ) continue; vector const & jets_KT4= *h_jet_KT4; fwlite::Handle > h_jet_KT6; h_jet_KT6.getByLabel(ev,"selectedLayer1JetsKT6"); if ( !h_jet_KT6.isValid() ) continue; vector const & jets_KT6= *h_jet_KT6; // Applies a cut on the collection of jets and prints the rest to the histogram // vector::const_iterator iter; // for ( iter = jets.begin(); iter != jets.end(); ++iter) { // if (iter->pt() > 10) // hist_jets->Fill( iter->correctedJet("raw").pt() ); // } // end jet for loop // if at least more than one jet is reconstructed: if ((h_jet_SC5->size() > 1) && (h_jet_SC7->size() > 1) && (h_jet_KT4->size() > 1) && (h_jet_KT6->size() > 1)){ // in the next lines you can replace "raw" by the other stages of corrections: // "off" - L1Offset offset correction // "rel" - L2Relative relative inter eta correction // "abs" - L3Absolute absolute pt correction // Note: if you don't use e.g. correctedJet("xxx").pt but jet.pt, // you will get by default the "abs" correction // In addition, you can correct for: // "emf" - L4Emf correction as a function of the jet emf // "had" - L5Flavour hadron level correction for gluons, light quarks, charm, beauty // "ue" - underlying event correction for gluons, light quarks, charm, beauty // "part" - L7Parton parton level correction for gluons, light quarks, charm, beauty // One more sideremark: by choosing one correction, the corrections mentioned above // your chosen correction are applied automatically. pat::Jet jets0_SC5=jets_SC5[0].correctedJet("raw"); pat::Jet jets1_SC5=jets_SC5[1].correctedJet("raw"); pat::Jet jets0_SC7=jets_SC7[0].correctedJet("raw"); pat::Jet jets1_SC7=jets_SC7[1].correctedJet("raw"); pat::Jet jets0_KT4=jets_KT4[0].correctedJet("raw"); pat::Jet jets1_KT4=jets_KT4[1].correctedJet("raw"); pat::Jet jets0_KT6=jets_KT6[0].correctedJet("raw"); pat::Jet jets1_KT6=jets_KT6[1].correctedJet("raw"); TLorentzVector jet1_SC5(jets0_SC5.px(), jets0_SC5.py(), jets0_SC5.pz(), jets0_SC5.energy() ); TLorentzVector jet2_SC5(jets1_SC5.px(), jets1_SC5.py(), jets1_SC5.pz(), jets1_SC5.energy() ); TLorentzVector jet1_SC7(jets0_SC7.px(), jets0_SC7.py(), jets0_SC7.pz(), jets0_SC7.energy() ); TLorentzVector jet2_SC7(jets1_SC7.px(), jets1_SC7.py(), jets1_SC7.pz(), jets1_SC7.energy() ); TLorentzVector jet1_KT4(jets0_KT4.px(), jets0_KT4.py(), jets0_KT4.pz(), jets0_KT4.energy() ); TLorentzVector jet2_KT4(jets1_KT4.px(), jets1_KT4.py(), jets1_KT4.pz(), jets1_KT4.energy() ); TLorentzVector jet1_KT6(jets0_KT6.px(), jets0_KT6.py(), jets0_KT6.pz(), jets0_KT6.energy() ); TLorentzVector jet2_KT6(jets1_KT6.px(), jets1_KT6.py(), jets1_KT6.pz(), jets1_KT6.energy() ); double invm_SC5 = (jet1_SC5 + jet2_SC5).Mag(); double invm_SC7 = (jet1_SC7 + jet2_SC7).Mag(); double invm_KT4 = (jet1_KT4 + jet2_KT4).Mag(); double invm_KT6 = (jet1_KT6 + jet2_KT6).Mag(); //PT hist_pt1_SC5->Fill(jets0_SC5.pt()); hist_pt2_SC5->Fill(jets1_SC5.pt()); hist_pt1_SC7->Fill(jets0_SC7.pt()); hist_pt2_SC7->Fill(jets1_SC7.pt()); hist_pt1_KT4->Fill(jets0_KT4.pt()); hist_pt2_KT4->Fill(jets1_KT4.pt()); hist_pt1_KT6->Fill(jets0_KT6.pt()); hist_pt2_KT6->Fill(jets1_KT6.pt()); //ETA hist_eta1_SC5->Fill(jets0_SC5.eta()); hist_eta2_SC5->Fill(jets1_SC5.eta()); hist_eta1_SC7->Fill(jets0_SC7.eta()); hist_eta2_SC7->Fill(jets1_SC7.eta()); hist_eta1_KT4->Fill(jets0_KT4.eta()); hist_eta2_KT4->Fill(jets1_KT4.eta()); hist_eta1_KT6->Fill(jets0_KT6.eta()); hist_eta2_KT6->Fill(jets1_KT6.eta()); //INV MASS hist_invarMass_SC5->Fill(invm_SC5); hist_invarMass_SC7->Fill(invm_SC7); hist_invarMass_KT4->Fill(invm_KT4); hist_invarMass_KT6->Fill(invm_KT6); } //end if loop } //end of event loop // TCanvas *invariantMass = new TCanvas("invariantMass", "invariantMass"); // invariantMass->Divide(3,5); // invariantMass->cd(1); // hist_pt1->Draw(); // invariantMass->cd(2); // hist_pt1_SC5->Draw(); // invariantMass->cd(3); // hist_pt1_KT4->Draw(); // invariantMass->cd(4); // hist_pt2->Draw(); // invariantMass->cd(5); // hist_pt2_SC5->Draw(); // invariantMass->cd(6); // hist_pt2_KT4->Draw(); // invariantMass->cd(7); // hist_eta1->Draw(); // invariantMass->cd(8); // hist_eta1_SC5->Draw(); // invariantMass->cd(9); // hist_eta1_KT4->Draw(); // invariantMass->cd(10); // hist_eta2->Draw(); // invariantMass->cd(11); // hist_eta2_SC5->Draw(); // invariantMass->cd(12); // hist_eta2_KT4->Draw(); // invariantMass->cd(13); // hist_invarMass->Draw(); // invariantMass->cd(14); // hist_invarMass_SC5->Draw(); // invariantMass->cd(15); // hist_invarMass_KT4->Draw(); // invariantMass->cd(6); // hist_delta_phi->Draw(); TCanvas *invariantMass = new TCanvas("invariantMass", "invariantMass"); invariantMass->Divide(2,2); invariantMass->cd(1); hist_invarMass_SC5->Draw(); invariantMass->cd(2); hist_invarMass_SC7->Draw(); invariantMass->cd(3); hist_invarMass_KT4->Draw(); invariantMass->cd(4); hist_invarMass_KT6->Draw(); TCanvas *pT = new TCanvas("pT", "pT"); pT->Divide(2,4); pT->cd(1); hist_pt1_SC5->Draw(); pT->cd(2); hist_pt2_SC5->Draw(); pT->cd(3); hist_pt1_SC7->Draw(); pT->cd(4); hist_pt2_SC7->Draw(); pT->cd(5); hist_pt1_KT4->Draw(); pT->cd(6); hist_pt2_KT4->Draw(); pT->cd(7); hist_pt1_KT6->Draw(); pT->cd(8); hist_pt2_KT6->Draw(); TCanvas *eta = new TCanvas("eta", "eta"); eta->Divide(2,4); eta->cd(1); hist_eta1_SC5->Draw(); eta->cd(2); hist_eta2_SC5->Draw(); eta->cd(3); hist_eta1_SC7->Draw(); eta->cd(4); hist_eta2_SC7->Draw(); eta->cd(5); hist_eta1_KT4->Draw(); eta->cd(6); hist_eta2_KT4->Draw(); eta->cd(7); hist_eta1_KT6->Draw(); eta->cd(8); hist_eta2_KT6->Draw(); }