CMS 3D CMS Logo

Public Member Functions | Static Public Member Functions | Private Member Functions | Private Attributes

muonisolation::CaloExtractor Class Reference

#include <CaloExtractor.h>

Inheritance diagram for muonisolation::CaloExtractor:
reco::isodeposit::IsoDepositExtractor

List of all members.

Public Member Functions

 CaloExtractor ()
 CaloExtractor (const edm::ParameterSet &par)
virtual reco::IsoDeposit deposit (const edm::Event &ev, const edm::EventSetup &evSetup, const reco::Track &track) const
virtual void fillVetos (const edm::Event &ev, const edm::EventSetup &evSetup, const reco::TrackCollection &tracks)
virtual ~CaloExtractor ()

Static Public Member Functions

static GlobalPoint MuonAtCaloPosition (const reco::Track &muon, const double bz, const GlobalPoint &endpos, bool fixVxy=false, bool fixVz=false)
 Extrapolate muons to calorimeter-object positions.

Private Member Functions

double noiseEcal (const CaloTower &tower) const
double noiseHcal (const CaloTower &tower) const

Private Attributes

edm::InputTag theCaloTowerCollectionLabel
std::string theDepositLabel
double theDR_Max
double theDR_Veto_E
double theDR_Veto_H
double theThreshold_E
double theThreshold_H
std::vector< DetIdtheVetoCollection
double theWeight_E
double theWeight_H
bool vertexConstraintFlag_XY
bool vertexConstraintFlag_Z

Detailed Description

Definition at line 18 of file CaloExtractor.h.


Constructor & Destructor Documentation

muonisolation::CaloExtractor::CaloExtractor ( ) [inline]

Definition at line 22 of file CaloExtractor.h.

{};
CaloExtractor::CaloExtractor ( const edm::ParameterSet par)

Definition at line 28 of file CaloExtractor.cc.

                                                    :
  theCaloTowerCollectionLabel(par.getParameter<edm::InputTag>("CaloTowerCollectionLabel")),
  theDepositLabel(par.getUntrackedParameter<string>("DepositLabel")),
  theWeight_E(par.getParameter<double>("Weight_E")),
  theWeight_H(par.getParameter<double>("Weight_H")),
  theThreshold_E(par.getParameter<double>("Threshold_E")),
  theThreshold_H(par.getParameter<double>("Threshold_H")),
  theDR_Veto_E(par.getParameter<double>("DR_Veto_E")),
  theDR_Veto_H(par.getParameter<double>("DR_Veto_H")),
  theDR_Max(par.getParameter<double>("DR_Max")),
  vertexConstraintFlag_XY(par.getParameter<bool>("Vertex_Constraint_XY")),
  vertexConstraintFlag_Z(par.getParameter<bool>("Vertex_Constraint_Z"))
{
}
virtual muonisolation::CaloExtractor::~CaloExtractor ( ) [inline, virtual]

Definition at line 25 of file CaloExtractor.h.

{}

Member Function Documentation

IsoDeposit CaloExtractor::deposit ( const edm::Event ev,
const edm::EventSetup evSetup,
const reco::Track track 
) const [virtual]

make single IsoDeposit based on track as input purely virtual: have to implement in concrete implementations

make this abit faster

Implements reco::isodeposit::IsoDepositExtractor.

Definition at line 95 of file CaloExtractor.cc.

References ecalTB2006H4_GenSimDigiReco_cfg::bField, deltaR(), PV3DBase< T, PVType, FrameType >::eta(), reco::TrackBase::eta(), spr::find(), edm::EventSetup::get(), LogDebug, MuonAtCaloPosition(), noiseEcal(), noiseHcal(), PV3DBase< T, PVType, FrameType >::phi(), reco::TrackBase::phi(), reco::TrackBase::pt(), theCaloTowerCollectionLabel, theDR_Max, theDR_Veto_E, theDR_Veto_H, theThreshold_E, theThreshold_H, theVetoCollection, theWeight_E, theWeight_H, vertexConstraintFlag_XY, and vertexConstraintFlag_Z.

{
  IsoDeposit dep(muon.eta(), muon.phi() );
  LogDebug("Muon|RecoMuon|L2MuonIsolationProducer")
          << " >>> Muon: pt " << muon.pt()
          << " eta " << muon.eta()
          << " phi " << muon.phi();

  Handle<CaloTowerCollection> towers;
  event.getByLabel(theCaloTowerCollectionLabel,towers);

  edm::ESHandle<CaloGeometry> caloGeom;
  eventSetup.get<CaloGeometryRecord>().get(caloGeom);

  edm::ESHandle<MagneticField> bField;
  eventSetup.get<IdealMagneticFieldRecord>().get(bField);
  double bz = bField->inInverseGeV(GlobalPoint(0.,0.,0.)).z();

  CaloTowerCollection::const_iterator cal;
  CaloTowerCollection::const_iterator calEnd(towers->end());
  for ( cal = towers->begin(); cal != calEnd; ++cal ) {
      double dEta = fabs(muon.eta()-cal->eta());
      if (fabs(dEta) > theDR_Max) continue;
    
      double deltar0 = reco::deltaR(muon,*cal);
      if (deltar0>theDR_Max) continue;

      double etecal = cal->emEt();
      double eecal = cal->emEnergy();
      bool doEcal = theWeight_E>0 && etecal>theThreshold_E && eecal>3*noiseEcal(*cal);
      double ethcal = cal->hadEt();
      double ehcal = cal->hadEnergy();
      bool doHcal = theWeight_H>0 && ethcal>theThreshold_H && ehcal>3*noiseHcal(*cal);
      if ((!doEcal) && (!doHcal)) continue;

      DetId calId = cal->id();
      GlobalPoint endpos = caloGeom->getPosition(calId);
      GlobalPoint muatcal = MuonAtCaloPosition(muon,bz,endpos,vertexConstraintFlag_XY, vertexConstraintFlag_Z);
      double deltar = reco::deltaR(muatcal,endpos);

      if (deltar<theDR_Veto_H) {
              dep.setVeto(IsoDeposit::Veto(reco::isodeposit::Direction(muatcal.eta(), muatcal.phi()), theDR_Veto_H));
      }

      if (doEcal) {
            if (deltar<theDR_Veto_E) { 
                  double calodep = theWeight_E*etecal;
                  if (doHcal) calodep += theWeight_H*ethcal;
                  dep.addCandEnergy(calodep);
                    LogDebug("Muon|RecoMuon|L2MuonIsolationProducer")
                    << " >>> Calo deposit inside veto (with ECAL): deltar " << deltar
                    << " calodep " << calodep
                    << " ecaldep " << etecal
                    << " hcaldep " << ethcal
                    << " eta " << cal->eta()
                    << " phi " << cal->phi();
                  continue;
            }
      } else {
            if (deltar<theDR_Veto_H) { 
                  dep.addCandEnergy(theWeight_H*ethcal);
                    LogDebug("Muon|RecoMuon|L2MuonIsolationProducer")
                    << " >>> Calo deposit inside veto (no ECAL): deltar " << deltar
                    << " calodep " << theWeight_H*ethcal
                    << " eta " << cal->eta()
                    << " phi " << cal->phi();
                  continue;
            }
      }

      if (std::find(theVetoCollection.begin(), theVetoCollection.end()
                  , calId)!=theVetoCollection.end()) {
            LogDebug("Muon|RecoMuon|L2MuonIsolationProducer")
            << " >>> Deposits belongs to other track: deltar, etecal, ethcal= " 
            << deltar << ", " << etecal << ", " << ethcal;
            continue;
      }

      if (doEcal) {
            if (deltar>theDR_Veto_E) { 
                  double calodep = theWeight_E*etecal;
                  if (doHcal) calodep += theWeight_H*ethcal;
                  dep.addDeposit(reco::isodeposit::Direction(endpos.eta(), endpos.phi()),calodep);
                    LogDebug("Muon|RecoMuon|L2MuonIsolationProducer")
                    << " >>> Calo deposit (with ECAL): deltar " << deltar
                    << " calodep " << calodep
                    << " ecaldep " << etecal
                    << " hcaldep " << ethcal
                    << " eta " << cal->eta()
                    << " phi " << cal->phi();
            }
      } else {
            if (deltar>theDR_Veto_H) { 
                    dep.addDeposit(reco::isodeposit::Direction(endpos.eta(), endpos.phi()),theWeight_H*ethcal);
                    LogDebug("Muon|RecoMuon|L2MuonIsolationProducer")
                    << " >>> Calo deposit (no ECAL): deltar " << deltar
                    << " calodep " << theWeight_H*ethcal
                    << " eta " << cal->eta()
                    << " phi " << cal->phi();
            }
      }
  }

  return dep;

}
void CaloExtractor::fillVetos ( const edm::Event ev,
const edm::EventSetup evSetup,
const reco::TrackCollection tracks 
) [virtual]

fill vetoes: to exclude deposits at IsoDeposit creation stage check concrete extractors if it's no-op !

make this abit faster

Implements reco::isodeposit::IsoDepositExtractor.

Definition at line 43 of file CaloExtractor.cc.

References ecalTB2006H4_GenSimDigiReco_cfg::bField, deltaR(), edm::EventSetup::get(), RPCpg::mu, MuonAtCaloPosition(), noiseEcal(), noiseHcal(), theCaloTowerCollectionLabel, theDR_Max, theDR_Veto_E, theDR_Veto_H, theThreshold_E, theThreshold_H, theVetoCollection, theWeight_E, theWeight_H, vertexConstraintFlag_XY, and vertexConstraintFlag_Z.

{
  theVetoCollection.clear();

  Handle<CaloTowerCollection> towers;
  event.getByLabel(theCaloTowerCollectionLabel,towers);

  edm::ESHandle<CaloGeometry> caloGeom;
  eventSetup.get<CaloGeometryRecord>().get(caloGeom);

  edm::ESHandle<MagneticField> bField;
  eventSetup.get<IdealMagneticFieldRecord>().get(bField);
  double bz = bField->inInverseGeV(GlobalPoint(0.,0.,0.)).z();

  TrackCollection::const_iterator mu;
  TrackCollection::const_iterator muEnd(muons.end());
  
  CaloTowerCollection::const_iterator cal;
  CaloTowerCollection::const_iterator calEnd(towers->end());
  
  for ( mu = muons.begin(); mu != muEnd; ++mu ) {
    for ( cal = towers->begin(); cal != calEnd; ++cal ) {
      double dEta = fabs(mu->eta()-cal->eta());
      if (fabs(dEta) > theDR_Max) continue;

            double deltar0 = reco::deltaR(*mu,*cal);
            if (deltar0>theDR_Max) continue;

            double etecal = cal->emEt();
            double eecal = cal->emEnergy();
            bool doEcal = theWeight_E>0 && etecal>theThreshold_E && eecal>3*noiseEcal(*cal);
            double ethcal = cal->hadEt();
            double ehcal = cal->hadEnergy();
            bool doHcal = theWeight_H>0 && ethcal>theThreshold_H && ehcal>3*noiseHcal(*cal);
            if ((!doEcal) && (!doHcal)) continue;

            DetId calId = cal->id();
            GlobalPoint endpos = caloGeom->getPosition(calId);
            GlobalPoint muatcal = MuonAtCaloPosition(*mu,bz,endpos, vertexConstraintFlag_XY, vertexConstraintFlag_Z);
            double deltar = reco::deltaR(muatcal,endpos);

            if (doEcal) {
                  if (deltar<theDR_Veto_E) theVetoCollection.push_back(calId);
            } else {
                  if (deltar<theDR_Veto_H) theVetoCollection.push_back(calId);
            }
      }
  }
     
}
GlobalPoint CaloExtractor::MuonAtCaloPosition ( const reco::Track muon,
const double  bz,
const GlobalPoint endpos,
bool  fixVxy = false,
bool  fixVz = false 
) [static]

Extrapolate muons to calorimeter-object positions.

Definition at line 203 of file CaloExtractor.cc.

References reco::TrackBase::charge(), funct::cos(), reco::TrackBase::covariance(), reco::TrackBase::dsz(), reco::TrackBase::dxy(), reco::TrackBase::dz(), reco::TrackBase::i_dsz, reco::TrackBase::i_dxy, reco::TrackBase::i_lambda, reco::TrackBase::i_phi, reco::TrackBase::i_qoverp, normalizedPhi(), reco::TrackBase::p(), reco::TrackBase::phi(), funct::pow(), reco::TrackBase::pt(), reco::TrackBase::qoverp(), funct::sin(), funct::tan(), theta(), reco::TrackBase::theta(), reco::TrackBase::vx(), reco::TrackBase::vy(), PV3DBase< T, PVType, FrameType >::x(), x, PV3DBase< T, PVType, FrameType >::y(), detailsBasic3DVector::y, PV3DBase< T, PVType, FrameType >::z(), and z.

Referenced by deposit(), and fillVetos().

                                                                                                                                    {
      double qoverp= muon.qoverp();
      double cur = bz*muon.charge()/muon.pt();
      double phi0 = muon.phi();
      double dca = muon.dxy();
      double theta = muon.theta();
      double dz = muon.dz();

      //LogDebug("Muon|RecoMuon|L2MuonIsolationProducer")
          //<< " Pt(GeV): " <<  muon.pt()
          //<< ", phi0 " <<  muon.phi0()
             //<< ", eta " <<  muon.eta();
      //LogDebug("Muon|RecoMuon|L2MuonIsolationProducer")
          //<< " d0 " <<  muon.d0()
          //<< ", dz " <<  muon.dz();
      //LogDebug("Muon|RecoMuon|L2MuonIsolationProducer")
          //<< " rhocal " <<  endpos.perp()
          //<< ", zcal " <<  endpos.z();

      if (fixVxy && fixVz) {
            // Note that here we assume no correlation between XY and Z projections
            // This should be a reasonable approximation for our purposes
            double errd02 = muon.covariance(muon.i_dxy,muon.i_dxy);
            if (pow(muon.dxy(),2)<4*errd02) {
                  phi0 -= muon.dxy()*muon.covariance(muon.i_dxy,muon.i_phi)
                                     /errd02;
                  cur -= muon.dxy()*muon.covariance(muon.i_dxy,muon.i_qoverp)
                                     /errd02 * (cur/qoverp);
                  dca = 0;
            } 
            double errdsz2 = muon.covariance(muon.i_dsz,muon.i_dsz);
            if (pow(muon.dsz(),2)<4*errdsz2) {
                  theta += muon.dsz()*muon.covariance(muon.i_dsz,muon.i_lambda)
                                     /errdsz2;
                  dz = 0;
            } 
      } else if (fixVxy) {
            double errd02 = muon.covariance(muon.i_dxy,muon.i_dxy);
            if (pow(muon.dxy(),2)<4*errd02) {
                  phi0  -= muon.dxy()*muon.covariance(muon.i_dxy,muon.i_phi)
                                     /errd02;
                  cur -= muon.dxy()*muon.covariance(muon.i_dxy,muon.i_qoverp)
                                     /errd02 * (cur/qoverp);
                  theta += muon.dxy()*muon.covariance(muon.i_dxy,muon.i_lambda)
                                     /errd02;
                  dz    -= muon.dxy()*muon.covariance(muon.i_dxy,muon.i_dsz)
                                     /errd02 * muon.p()/muon.pt();
                  dca = 0;
            } 
      } else if (fixVz) {
            double errdsz2 = muon.covariance(muon.i_dsz,muon.i_dsz);
            if (pow(muon.dsz(),2)<4*errdsz2) {
                  theta += muon.dsz()*muon.covariance(muon.i_dsz,muon.i_lambda)
                                     /errdsz2;
                  phi0  -= muon.dsz()*muon.covariance(muon.i_dsz,muon.i_phi)
                                     /errdsz2;
                  cur -= muon.dsz()*muon.covariance(muon.i_dsz,muon.i_qoverp)
                                     /errdsz2 * (cur/qoverp);
                  dca   -= muon.dsz()*muon.covariance(muon.i_dsz,muon.i_dxy)
                                     /errdsz2;
                  dz = 0;
            } 
      }

      double sphi0 = sin(phi0);
      double cphi0 = cos(phi0);

      double xsin =  endpos.x()*sphi0 - endpos.y()*cphi0;
      double xcos =  endpos.x()*cphi0 + endpos.y()*sphi0;
      double fcdca = fabs(1-cur*dca);
      double phif = atan2( fcdca*sphi0-cur*endpos.x()
                         , fcdca*cphi0+cur*endpos.y());
      double tphif2 = tan(0.5*(phif-phi0));
      double dcaf = dca + xsin + xcos*tphif2;

      double x = endpos.x() - dcaf*sin(phif);
      double y = endpos.y() + dcaf*cos(phif);

      double deltas =  (x-muon.vx())*cphi0 + (y-muon.vy())*sphi0;
      double deltaphi = normalizedPhi(phif-phi0);
      if (deltaphi!=0) deltas = deltas*deltaphi/sin(deltaphi);

      double z =dz;
      double tantheta = tan(theta);
      if (tantheta!=0) {
            z += deltas/tan(theta);
      } else {
            z = endpos.z();
      }

      return GlobalPoint(x,y,z);
}
double CaloExtractor::noiseEcal ( const CaloTower tower) const [private]

Definition at line 296 of file CaloExtractor.cc.

References eta(), and reco::LeafCandidate::eta().

Referenced by deposit(), and fillVetos().

                                                            {
      double noise = 0.04;
      double eta = tower.eta();
      if (fabs(eta)>1.479) noise = 0.15;
      return noise;
}
double CaloExtractor::noiseHcal ( const CaloTower tower) const [private]

Definition at line 303 of file CaloExtractor.cc.

Referenced by deposit(), and fillVetos().

                                                            {
      double noise = 0.2;
      return noise;
}

Member Data Documentation

Definition at line 35 of file CaloExtractor.h.

Referenced by deposit(), and fillVetos().

Definition at line 38 of file CaloExtractor.h.

Definition at line 47 of file CaloExtractor.h.

Referenced by deposit(), and fillVetos().

Definition at line 45 of file CaloExtractor.h.

Referenced by deposit(), and fillVetos().

Definition at line 46 of file CaloExtractor.h.

Referenced by deposit(), and fillVetos().

Definition at line 43 of file CaloExtractor.h.

Referenced by deposit(), and fillVetos().

Definition at line 44 of file CaloExtractor.h.

Referenced by deposit(), and fillVetos().

Definition at line 52 of file CaloExtractor.h.

Referenced by deposit(), and fillVetos().

Definition at line 41 of file CaloExtractor.h.

Referenced by deposit(), and fillVetos().

Definition at line 42 of file CaloExtractor.h.

Referenced by deposit(), and fillVetos().

Definition at line 48 of file CaloExtractor.h.

Referenced by deposit(), and fillVetos().

Definition at line 49 of file CaloExtractor.h.

Referenced by deposit(), and fillVetos().