MuonEnergyDepositAnalyzer.cc

Go to the documentation of this file.

/*
 *  See header file for a description of this class.
 *
 *  \author G. Mila - INFN Torino
 */

#include "DQMOffline/Muon/interface/MuonEnergyDepositAnalyzer.h"

#include "DataFormats/Common/interface/Handle.h"

#include "FWCore/MessageLogger/interface/MessageLogger.h"

#include "DataFormats/TrackReco/interface/Track.h"
#include "DataFormats/TrackReco/interface/TrackFwd.h"
#include "TrackingTools/TransientTrack/interface/TransientTrackBuilder.h"
#include "TrackingTools/Records/interface/TransientTrackRecord.h"
#include "TrackingTools/TransientTrack/interface/TransientTrack.h"

#include <cmath>
#include <string>
using namespace std;
using namespace edm;

MuonEnergyDepositAnalyzer::MuonEnergyDepositAnalyzer(const edm::ParameterSet& pSet){
  parameters = pSet;
  
  // the services
  theService = new MuonServiceProxy(parameters.getParameter<ParameterSet>("ServiceParameters"));

  theMuonCollectionLabel_ = consumes<reco::MuonCollection>(parameters.getParameter<InputTag>("MuonCollection"));

  AlgoName = parameters.getParameter<std::string>("AlgoName");
  emNoBin = parameters.getParameter<int>("emSizeBin");
  emNoMin = parameters.getParameter<double>("emSizeMin");
  emNoMax = parameters.getParameter<double>("emSizeMax");

  hadNoBin = parameters.getParameter<int>("hadSizeBin");
  hadNoMin = parameters.getParameter<double>("hadSizeMin");
  hadNoMax = parameters.getParameter<double>("hadSizeMax");
  
  emS9NoBin = parameters.getParameter<int>("emS9SizeBin");
  emS9NoMin = parameters.getParameter<double>("emS9SizeMin");
  emS9NoMax = parameters.getParameter<double>("emS9SizeMax");
  
  hadS9NoBin = parameters.getParameter<int>("hadS9SizeBin");
  hadS9NoMin = parameters.getParameter<double>("hadS9SizeMin");
  hadS9NoMax = parameters.getParameter<double>("hadS9SizeMax");

  hoNoBin = parameters.getParameter<int>("hoSizeBin");
  hoNoMin = parameters.getParameter<double>("hoSizeMin");
  hoNoMax = parameters.getParameter<double>("hoSizeMax");
  
  hoS9NoBin = parameters.getParameter<int>("hoS9SizeBin");
  hoS9NoMin = parameters.getParameter<double>("hoS9SizeMin");
  hoS9NoMax = parameters.getParameter<double>("hoS9SizeMax");
}
MuonEnergyDepositAnalyzer::~MuonEnergyDepositAnalyzer() {   
  delete theService;  
}
void MuonEnergyDepositAnalyzer::bookHistograms(DQMStore::IBooker & ibooker,
                           edm::Run const & /*iRun*/,
                           edm::EventSetup const & /* iSetup */){
  ibooker.cd();
  ibooker.setCurrentFolder("Muons/MuonEnergyDepositAnalyzer");

  std::string histname = "ecalDepositedEnergyBarrel_";
  ecalDepEnergyBarrel = ibooker.book1D(histname+AlgoName, "Energy deposited in the ECAL barrel cells", emNoBin, emNoMin, emNoMax);
  ecalDepEnergyBarrel->setAxisTitle("GeV");
  
  histname = "ecalDepositedEnergyEndcap_";
  ecalDepEnergyEndcap = ibooker.book1D(histname+AlgoName, "Energy deposited in the ECAL endcap cells", emNoBin, emNoMin, emNoMax);
  ecalDepEnergyEndcap->setAxisTitle("GeV");

  histname = "ecalS9DepositedEnergyBarrel_";
  ecalS9DepEnergyBarrel = ibooker.book1D(histname+AlgoName, "Energy deposited in the ECAL barrel 3*3 towers", emS9NoBin, emS9NoMin, emS9NoMax);
  ecalS9DepEnergyBarrel->setAxisTitle("GeV");
  histname = "ecalS9DepositedEnergyEndcap_";
  ecalS9DepEnergyEndcap = ibooker.book1D(histname+AlgoName, "Energy deposited in the ECAL endcap 3*3 towers", emS9NoBin, emS9NoMin, emS9NoMax);
  ecalS9DepEnergyEndcap->setAxisTitle("GeV");
  histname = "ecalS9PointingMuDepositedEnergy_Glb_";
  ecalS9PointingMuDepEnergy_Glb = ibooker.book1D(histname+AlgoName, "Pointing glb muons energy deposited in the ECAL 3*3 towers", emS9NoBin, emS9NoMin, emS9NoMax);
  ecalS9PointingMuDepEnergy_Glb->setAxisTitle("GeV"); 
  histname = "ecalS9PointingMuDepositedEnergy_Tk_";
  ecalS9PointingMuDepEnergy_Tk = ibooker.book1D(histname+AlgoName, "Pointing tk muons energy deposited in the ECAL 3*3 towers", emS9NoBin, emS9NoMin, emS9NoMax);
  ecalS9PointingMuDepEnergy_Tk->setAxisTitle("GeV"); 
  histname = "ecalS9PointingMuDepositedEnergy_Sta_";
  ecalS9PointingMuDepEnergy_Sta = ibooker.book1D(histname+AlgoName, "Pointing sta muons energy deposited in the ECAL 3*3 towers", emS9NoBin, emS9NoMin, emS9NoMax);
  ecalS9PointingMuDepEnergy_Sta->setAxisTitle("GeV"); 
  
  histname = "hadDepositedEnergyBarrel_";
  hcalDepEnergyBarrel = ibooker.book1D(histname+AlgoName, "Energy deposited in the HCAL barrel cells", hadNoBin, hadNoMin, hadNoMax);
  hcalDepEnergyBarrel->setAxisTitle("GeV");
  histname = "hadDepositedEnergyEndcap_";
  hcalDepEnergyEndcap = ibooker.book1D(histname+AlgoName, "Energy deposited in the HCAL endcap cells", hadNoBin, hadNoMin, hadNoMax);
  hcalDepEnergyEndcap->setAxisTitle("GeV");
  
  histname = "hadS9DepositedEnergyBarrel_";
  hcalS9DepEnergyBarrel = ibooker.book1D(histname+AlgoName, "Energy deposited in the HCAL barrel 3*3 towers", hadS9NoBin, hadS9NoMin, hadS9NoMax);
  hcalS9DepEnergyBarrel->setAxisTitle("GeV");
  histname = "hadS9DepositedEnergyEndcap_";
  hcalS9DepEnergyEndcap = ibooker.book1D(histname+AlgoName, "Energy deposited in the HCAL endcap 3*3 towers", hadS9NoBin, hadS9NoMin, hadS9NoMax);
  hcalS9DepEnergyEndcap->setAxisTitle("GeV");
  histname = "hadS9PointingMuDepositedEnergy_Glb_";
  hcalS9PointingMuDepEnergy_Glb = ibooker.book1D(histname+AlgoName, "Pointing glb muons energy deposited in the HCAL endcap 3*3 towers", hadS9NoBin, hadS9NoMin, hadS9NoMax);
  hcalS9PointingMuDepEnergy_Glb->setAxisTitle("GeV");
  histname = "hadS9PointingMuDepositedEnergy_Tk_";
  hcalS9PointingMuDepEnergy_Tk = ibooker.book1D(histname+AlgoName, "Pointing tk muons energy deposited in the HCAL endcap 3*3 towers", hadS9NoBin, hadS9NoMin, hadS9NoMax);
  hcalS9PointingMuDepEnergy_Tk->setAxisTitle("GeV");
  histname = "hadS9PointingMuDepositedEnergy_Sta_";
  hcalS9PointingMuDepEnergy_Sta = ibooker.book1D(histname+AlgoName, "Pointing sta muons energy deposited in the HCAL endcap 3*3 towers", hadS9NoBin, hadS9NoMin, hadS9NoMax);
  hcalS9PointingMuDepEnergy_Sta->setAxisTitle("GeV");

  histname = "hoDepositedEnergy_";
  hoDepEnergy = ibooker.book1D(histname+AlgoName, "Energy deposited in the HO cells", hoNoBin, hoNoMin, hoNoMax);
  hoDepEnergy->setAxisTitle("GeV");

  histname = "hoS9DepositedEnergy_";
  hoS9DepEnergy = ibooker.book1D(histname+AlgoName, "Energy deposited in the HO 3*3 towers", hoS9NoBin, hoS9NoMin, hoS9NoMax);
  hoS9DepEnergy->setAxisTitle("GeV");
  histname = "hoS9PointingMuDepositedEnergy_Glb_";
  hoS9PointingMuDepEnergy_Glb = ibooker.book1D(histname+AlgoName, "Pointing glb muons energy deposited in the HO 3*3 towers", hoS9NoBin, hoS9NoMin, hoS9NoMax);
  hoS9PointingMuDepEnergy_Glb->setAxisTitle("GeV");
  histname = "hoS9PointingMuDepositedEnergy_Tk_";
  hoS9PointingMuDepEnergy_Tk = ibooker.book1D(histname+AlgoName, "Pointing tk muons energy deposited in the HO 3*3 towers", hoS9NoBin, hoS9NoMin, hoS9NoMax);
  hoS9PointingMuDepEnergy_Tk->setAxisTitle("GeV");
  histname = "hoS9PointingMuDepositedEnergy_Sta_";
  hoS9PointingMuDepEnergy_Sta = ibooker.book1D(histname+AlgoName, "Pointing sta muons energy deposited in the HO 3*3 towers", hoS9NoBin, hoS9NoMin, hoS9NoMax);
  hoS9PointingMuDepEnergy_Sta->setAxisTitle("GeV");
}
void MuonEnergyDepositAnalyzer::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup){
  LogTrace(metname)<<"[MuonEnergyDepositAnalyzer] Filling the histos";
  theService->update(iSetup);
  
  // Take the muon container
  edm::Handle<reco::MuonCollection> muons;
  iEvent.getByToken(theMuonCollectionLabel_,muons);
  
  if (!muons.isValid()) return;
  
  for (reco::MuonCollection::const_iterator recoMu = muons->begin(); recoMu!=muons->end(); ++recoMu){
    // get all the mu energy deposits
    reco::MuonEnergy muEnergy = recoMu->calEnergy();
  
    // energy deposited in ECAL
    LogTrace(metname) << "Energy deposited in ECAL: "<<muEnergy.em;
    if (fabs(recoMu->eta()) > 1.479) 
      ecalDepEnergyEndcap->Fill(muEnergy.em);
    else
      ecalDepEnergyBarrel->Fill(muEnergy.em);

    // energy deposited in HCAL
    LogTrace(metname) << "Energy deposited in HCAL: "<<muEnergy.had;
    if (fabs(recoMu->eta()) > 1.4)
      hcalDepEnergyEndcap->Fill(muEnergy.had);
    else
      hcalDepEnergyBarrel->Fill(muEnergy.had);
  
    // energy deposited in HO
    LogTrace(metname) << "Energy deposited in HO: "<<muEnergy.ho;
    if (fabs(recoMu->eta()) < 1.26)
      hoDepEnergy->Fill(muEnergy.ho);
  
    // energy deposited in ECAL in 3*3 towers
    LogTrace(metname) << "Energy deposited in ECAL: "<<muEnergy.emS9;
    if (fabs(recoMu->eta()) > 1.479) 
      ecalS9DepEnergyEndcap->Fill(muEnergy.emS9);
    else
      ecalS9DepEnergyBarrel->Fill(muEnergy.emS9);
     
    // energy deposited in HCAL in 3*3 crystals
    LogTrace(metname) << "Energy deposited in HCAL: "<<muEnergy.hadS9;
    if (fabs(recoMu->eta()) > 1.4)
      hcalS9DepEnergyEndcap->Fill(muEnergy.hadS9);
    else
      hcalS9DepEnergyBarrel->Fill(muEnergy.hadS9);
  
    // energy deposited in HO in 3*3 crystals
    LogTrace(metname) << "Energy deposited in HO: "<<muEnergy.hoS9;
    if (fabs(recoMu->eta()) < 1.26)
      hoS9DepEnergy->Fill(muEnergy.hoS9);
  
    // plot for energy tests
    edm::ESHandle<TransientTrackBuilder> theB;
    iSetup.get<TransientTrackRecord>().get("TransientTrackBuilder",theB);
    reco::TransientTrack TransTrack;
  
    if(recoMu->isGlobalMuon())
      TransTrack = theB->build(recoMu->globalTrack());
    else if((recoMu->isTrackerMuon() || recoMu->isRPCMuon()))
      TransTrack = theB->build(recoMu->innerTrack());
    else if(recoMu->isStandAloneMuon())
      TransTrack = theB->build(recoMu->outerTrack());
    else
      continue;

    TrajectoryStateOnSurface TSOS;
    TSOS = TransTrack.impactPointState();
    // section for vertex pointing muon
    if((abs(TSOS.globalPosition().z())<30) && (abs(TSOS.globalPosition().perp())<20)){
      // GLB muon
      if(recoMu->isGlobalMuon()){
    ecalS9PointingMuDepEnergy_Glb->Fill(muEnergy.emS9);
    hcalS9PointingMuDepEnergy_Glb->Fill(muEnergy.hadS9);
    hoS9PointingMuDepEnergy_Glb->Fill(muEnergy.hoS9);
      }
      // TK muon
      else if(recoMu->isTrackerMuon()){
    ecalS9PointingMuDepEnergy_Tk->Fill(muEnergy.emS9);
    hcalS9PointingMuDepEnergy_Tk->Fill(muEnergy.hadS9);
    hoS9PointingMuDepEnergy_Tk->Fill(muEnergy.hoS9);
      }
      // STA muon
      else if(recoMu->isStandAloneMuon()){
    ecalS9PointingMuDepEnergy_Sta->Fill(muEnergy.emS9);
    hcalS9PointingMuDepEnergy_Sta->Fill(muEnergy.hadS9);
    hoS9PointingMuDepEnergy_Sta->Fill(muEnergy.hoS9);
      } 
    }
  }

}