Skip to content

Commit

Permalink
Merge pull request #20130 from bsunanda/Run2-hcx142
Browse files Browse the repository at this point in the history 
Run2-hcx142 Make codes compatible with current HCAL geometry
  • Loading branch information
@cmsbuild
cmsbuild committed Aug 18, 2017
2 parents d4519c4 + 3843fdf commit 28b86df
Show file tree
Hide file tree
Showing 8 changed files with 141 additions and 151 deletions.
47 changes: 21 additions & 26 deletions Calibration/HcalIsolatedTrackReco/interface/EcalIsolatedParticleCandidateProducer.h
View file
Expand Up @@ -2,46 +2,41 @@
#define Calibration_EcalIsolatedParticleCandidateProducer_h

#include "FWCore/Framework/interface/Frameworkfwd.h"
#include "FWCore/Framework/interface/EDProducer.h"
#include "FWCore/Framework/interface/global/EDProducer.h"

#include "FWCore/Framework/interface/Event.h"
#include "FWCore/Framework/interface/MakerMacros.h"

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

#include "Geometry/CaloGeometry/interface/CaloGeometry.h"

#include "DataFormats/EcalRecHit/interface/EcalRecHitCollections.h"
#include "DataFormats/HLTReco/interface/TriggerFilterObjectWithRefs.h"
#include "DataFormats/HcalIsolatedTrack/interface/IsolatedPixelTrackCandidate.h"
//
// class decleration
//

class EcalIsolatedParticleCandidateProducer : public edm::EDProducer {
public:
explicit EcalIsolatedParticleCandidateProducer(const edm::ParameterSet&);
~EcalIsolatedParticleCandidateProducer();

private:

const CaloGeometry* geo;

double InConeSize_;
double OutConeSize_;
double hitCountEthr_;
double hitEthr_;

edm::EDGetTokenT<l1extra::L1JetParticleCollection> tok_l1tau_;
edm::EDGetTokenT<trigger::TriggerFilterObjectWithRefs> tok_hlt_;
edm::EDGetTokenT<EcalRecHitCollection> tok_EB_;
edm::EDGetTokenT<EcalRecHitCollection> tok_EE_;

virtual void beginJob() ;
virtual void produce(edm::Event&, const edm::EventSetup&);
virtual void endJob() ;
class EcalIsolatedParticleCandidateProducer : public edm::global::EDProducer<> {
public:
explicit EcalIsolatedParticleCandidateProducer(const edm::ParameterSet&);
~EcalIsolatedParticleCandidateProducer();

private:
void produce(edm::StreamID, edm::Event&, const edm::EventSetup&) const override;
void beginJob() override;
void endJob() override;

double InConeSize_;
double OutConeSize_;
double hitCountEthr_;
double hitEthr_;

edm::EDGetTokenT<l1extra::L1JetParticleCollection> tok_l1tau_;
edm::EDGetTokenT<trigger::TriggerFilterObjectWithRefs> tok_hlt_;
edm::EDGetTokenT<EcalRecHitCollection> tok_EB_;
edm::EDGetTokenT<EcalRecHitCollection> tok_EE_;

// ----------member data ---------------------------
// ----------member data ---------------------------
};

#endif
196 changes: 88 additions & 108 deletions Calibration/HcalIsolatedTrackReco/src/EcalIsolatedParticleCandidateProducer.cc
View file
Expand Up @@ -18,6 +18,7 @@


// system include files
#include <cmath>
#include <memory>

// user include files
Expand All @@ -27,15 +28,12 @@

#include "DataFormats/Common/interface/Ref.h"
#include "DataFormats/DetId/interface/DetId.h"
#include "Geometry/Records/interface/CaloGeometryRecord.h"
#include "Geometry/CaloGeometry/interface/CaloGeometry.h"
#include "Geometry/CaloTopology/interface/EcalBarrelTopology.h"

#include "Geometry/Records/interface/CaloGeometryRecord.h"

#include "Calibration/HcalIsolatedTrackReco/interface/EcalIsolatedParticleCandidateProducer.h"




EcalIsolatedParticleCandidateProducer::EcalIsolatedParticleCandidateProducer(const edm::ParameterSet& conf)
{
InConeSize_ = conf.getParameter<double>("EcalInnerConeSize");
Expand Down Expand Up @@ -68,142 +66,124 @@ EcalIsolatedParticleCandidateProducer::~EcalIsolatedParticleCandidateProducer()

// ------------ method called to produce the data ------------
void
EcalIsolatedParticleCandidateProducer::produce(edm::Event& iEvent, const edm::EventSetup& iSetup)
{

using namespace edm;
EcalIsolatedParticleCandidateProducer::produce(edm::StreamID, edm::Event& iEvent, const edm::EventSetup& iSetup) const {

// std::cout<<"get tau"<<std::endl;

Handle<l1extra::L1JetParticleCollection> l1Taus;
edm::Handle<l1extra::L1JetParticleCollection> l1Taus;
iEvent.getByToken(tok_l1tau_,l1Taus);

// std::cout<<"get geom"<<std::endl;

ESHandle<CaloGeometry> pG;
edm::ESHandle<CaloGeometry> pG;
iSetup.get<CaloGeometryRecord>().get(pG);
geo = pG.product();
const CaloGeometry* geo = pG.product();

// std::cout<<" get ec rechit"<<std::endl;

Handle<EcalRecHitCollection> ecalEB;
edm::Handle<EcalRecHitCollection> ecalEB;
iEvent.getByToken(tok_EB_,ecalEB);

Handle<EcalRecHitCollection> ecalEE;
edm::Handle<EcalRecHitCollection> ecalEE;
iEvent.getByToken(tok_EE_,ecalEE);

// std::cout<<"get l1 trig obj"<<std::endl;

Handle<trigger::TriggerFilterObjectWithRefs> l1trigobj;
edm::Handle<trigger::TriggerFilterObjectWithRefs> l1trigobj;
iEvent.getByToken(tok_hlt_, l1trigobj);

std::vector< edm::Ref<l1extra::L1JetParticleCollection> > l1tauobjref;
std::vector< edm::Ref<l1extra::L1JetParticleCollection> > l1jetobjref;

l1trigobj->getObjects(trigger::TriggerL1TauJet, l1tauobjref);
l1trigobj->getObjects(trigger::TriggerL1CenJet, l1jetobjref);
std::vector< edm::Ref<l1t::TauBxCollection> > l1tauobjref;
std::vector< edm::Ref<l1t::JetBxCollection> > l1jetobjref;
l1trigobj->getObjects(trigger::TriggerL1Tau, l1tauobjref);
l1trigobj->getObjects(trigger::TriggerL1Jet, l1jetobjref);

double ptTriggered=-10;
double etaTriggered=-100;
double phiTriggered=-100;

// std::cout<<"find highest pT triggered obj"<<std::endl;

for (unsigned int p=0; p<l1tauobjref.size(); p++)
{
if (l1tauobjref[p]->pt()>ptTriggered)
{
ptTriggered=l1tauobjref[p]->pt();
phiTriggered=l1tauobjref[p]->phi();
etaTriggered=l1tauobjref[p]->eta();
}
}
for (unsigned int p=0; p<l1jetobjref.size(); p++)
{
if (l1jetobjref[p]->pt()>ptTriggered)
{
ptTriggered=l1jetobjref[p]->pt();
phiTriggered=l1jetobjref[p]->phi();
etaTriggered=l1jetobjref[p]->eta();
}
}
for (unsigned int p=0; p<l1tauobjref.size(); p++) {
if (l1tauobjref[p]->pt()>ptTriggered) {
ptTriggered=l1tauobjref[p]->pt();
phiTriggered=l1tauobjref[p]->phi();
etaTriggered=l1tauobjref[p]->eta();
}
}
for (unsigned int p=0; p<l1jetobjref.size(); p++) {
if (l1jetobjref[p]->pt()>ptTriggered) {
ptTriggered=l1jetobjref[p]->pt();
phiTriggered=l1jetobjref[p]->phi();
etaTriggered=l1jetobjref[p]->eta();
}
}

auto iptcCollection = std::make_unique<reco::IsolatedPixelTrackCandidateCollection>();

// std::cout<<"loop over l1taus"<<std::endl;

for (l1extra::L1JetParticleCollection::const_iterator tit=l1Taus->begin(); tit!=l1Taus->end(); tit++)
{
double dphi=fabs(tit->phi()-phiTriggered);
if (dphi>3.1415926535) dphi=2*3.1415926535-dphi;
double Rseed=sqrt(pow(etaTriggered-tit->eta(),2)+dphi*dphi);
if (Rseed<1.2) continue;
int nhitOut=0;
int nhitIn=0;
double OutEnergy=0;
double InEnergy=0;
for (l1extra::L1JetParticleCollection::const_iterator tit=l1Taus->begin(); tit!=l1Taus->end(); tit++) {
double dphi=fabs(tit->phi()-phiTriggered);
if (dphi>M_PI) dphi=2*M_PI-dphi;
double Rseed=sqrt(pow(etaTriggered-tit->eta(),2)+dphi*dphi);
if (Rseed<1.2) continue;
int nhitOut=0;
int nhitIn=0;
double OutEnergy=0;
double InEnergy=0;
// std::cout<<" loops over rechits"<<std::endl;
for (EcalRecHitCollection::const_iterator eItr=ecalEB->begin(); eItr!=ecalEB->end(); eItr++)
{
double phiD, R;
GlobalPoint pos = geo->getPosition(eItr->detid());
double phihit = pos.phi();
double etahit = pos.eta();
phiD=fabs(phihit-tit->phi());
if (phiD>3.1415926535) phiD=2*3.1415926535-phiD;
R=sqrt(pow(etahit-tit->eta(),2)+phiD*phiD);
for (EcalRecHitCollection::const_iterator eItr=ecalEB->begin(); eItr!=ecalEB->end(); eItr++) {
double phiD, R;
const GlobalPoint& pos = geo->getPosition(eItr->detid());
double phihit = pos.phi();
double etahit = pos.eta();
phiD=fabs(phihit-tit->phi());
if (phiD>M_PI) phiD=2*M_PI-phiD;
R=sqrt(pow(etahit-tit->eta(),2)+phiD*phiD);

if (R<OutConeSize_&&R>InConeSize_&&eItr->energy()>hitCountEthr_)
{
nhitOut++;
}
if (R<InConeSize_&&eItr->energy()>hitCountEthr_)
{
nhitIn++;
}

if (R<OutConeSize_&&R>InConeSize_&&eItr->energy()>hitEthr_)
{
OutEnergy+=eItr->energy();
}
if (R<InConeSize_&&eItr->energy()>hitEthr_)
{
InEnergy+=eItr->energy();
}

}

for (EcalRecHitCollection::const_iterator eItr=ecalEE->begin(); eItr!=ecalEE->end(); eItr++)
{
double phiD, R;
GlobalPoint pos = geo->getPosition(eItr->detid());
double phihit = pos.phi();
double etahit = pos.eta();
phiD=fabs(phihit-tit->phi());
if (phiD>3.1415926535) phiD=2*3.1415926535-phiD;
R=sqrt(pow(etahit-tit->eta(),2)+phiD*phiD);
if (R<OutConeSize_&&R>InConeSize_&&eItr->energy()>hitCountEthr_)
{
nhitOut++;
}
if (R<InConeSize_&&eItr->energy()>hitCountEthr_)
{
nhitIn++;
}
if (R<OutConeSize_&&R>InConeSize_&&eItr->energy()>hitEthr_)
{
OutEnergy+=eItr->energy();
}
if (R<InConeSize_&&eItr->energy()>hitEthr_)
{
InEnergy+=eItr->energy();
}

}
if (R<OutConeSize_&&R>InConeSize_&&eItr->energy()>hitCountEthr_) {
nhitOut++;
}
if (R<InConeSize_&&eItr->energy()>hitCountEthr_) {
nhitIn++;
}

if (R<OutConeSize_&&R>InConeSize_&&eItr->energy()>hitEthr_) {
OutEnergy+=eItr->energy();
}
if (R<InConeSize_&&eItr->energy()>hitEthr_) {
InEnergy+=eItr->energy();
}

}

for (EcalRecHitCollection::const_iterator eItr=ecalEE->begin(); eItr!=ecalEE->end(); eItr++) {
double phiD, R;
const GlobalPoint& pos = geo->getPosition(eItr->detid());
double phihit = pos.phi();
double etahit = pos.eta();
phiD=fabs(phihit-tit->phi());
if (phiD>M_PI) phiD=2*M_PI-phiD;
R=sqrt(pow(etahit-tit->eta(),2)+phiD*phiD);
if (R<OutConeSize_&&R>InConeSize_&&eItr->energy()>hitCountEthr_) {
nhitOut++;
}
if (R<InConeSize_&&eItr->energy()>hitCountEthr_) {
nhitIn++;
}
if (R<OutConeSize_&&R>InConeSize_&&eItr->energy()>hitEthr_) {
OutEnergy+=eItr->energy();
}
if (R<InConeSize_&&eItr->energy()>hitEthr_) {
InEnergy+=eItr->energy();
}

}
// std::cout<<"create and push_back candidate"<<std::endl;
reco::IsolatedPixelTrackCandidate newca(l1extra::L1JetParticleRef(l1Taus,tit-l1Taus->begin()), InEnergy, OutEnergy, nhitIn, nhitOut);
iptcCollection->push_back(newca);
}
reco::IsolatedPixelTrackCandidate newca(l1extra::L1JetParticleRef(l1Taus,tit-l1Taus->begin()), InEnergy, OutEnergy, nhitIn, nhitOut);
iptcCollection->push_back(newca);
}



Expand Down
5 changes: 3 additions & 2 deletions Calibration/IsolatedParticles/BuildFile.xml
View file
Expand Up @@ -27,11 +27,12 @@
<use name="SimTracker/TrackerHitAssociation"/>
<use name="SimGeneral/HepPDTRecord"/>
<use name="Geometry/CaloGeometry"/>
<use name="Geometry/CaloTopology"/>
<use name="Geometry/CommonDetUnit"/>
<use name="Geometry/TrackerGeometryBuilder"/>
<use name="Geometry/EcalAlgo"/>
<use name="Geometry/CaloTopology"/>
<use name="Geometry/HcalTowerAlgo"/>
<use name="Geometry/Records"/>
<use name="Geometry/TrackerGeometryBuilder"/>
<use name="RecoLocalCalo/EcalRecAlgos"/>
<use name="TrackingTools/TrajectoryState"/>
<use name="TrackingTools/GeomPropagators"/>
Expand Down
4 changes: 2 additions & 2 deletions Calibration/IsolatedParticles/interface/CaloSimInfo.icc
View file
Expand Up @@ -307,7 +307,7 @@ namespace spr{
std::cout << " " << matchedId[it];
std::cout << std::endl;
}
if (matchedId.size() > 0) {
if (!matchedId.empty()) {
typename T::const_iterator ihit;
for (ihit=hits->begin(); ihit!=hits->end(); ihit++) {

Expand Down Expand Up @@ -348,7 +348,7 @@ namespace spr{
std::cout << " " << matchedId[it];
std::cout << std::endl;
}
if (matchedId.size() > 0) {
if (!matchedId.empty()) {
typename T::const_iterator ihit;
for (ihit=hitsEB->begin(); ihit!=hitsEB->end(); ihit++) {

Expand Down
9 changes: 7 additions & 2 deletions Calibration/IsolatedParticles/src/CaloSimInfo.cc
View file
@@ -1,7 +1,9 @@
#include "DataFormats/DetId/interface/DetId.h"
#include "DataFormats/EcalDetId/interface/EcalSubdetector.h"
#include "DataFormats/HcalDetId/interface/HcalSubdetector.h"
#include "Calibration/IsolatedParticles/interface/CaloConstants.h"
#include "Calibration/IsolatedParticles/interface/CaloSimInfo.h"
#include "Geometry/HcalTowerAlgo/interface/HcalGeometry.h"

#include "CLHEP/Units/PhysicalConstants.h"
#include "CLHEP/Units/SystemOfUnits.h"
Expand All @@ -14,13 +16,16 @@ namespace spr{

double timeOfFlight(DetId id, const CaloGeometry* geo, bool debug) {

double R = geo->getPosition(id).mag();
GlobalPoint point = (id.det() == DetId::Hcal) ?
((HcalGeometry*)(geo->getSubdetectorGeometry(id)))->getPosition(id) :
geo->getPosition(id);
double R = point.mag();
double tmp = R/CLHEP::c_light/CLHEP::ns;
#ifdef EDM_ML_DEBUG
if (debug) {
DetId::Detector det = id.det();
int subdet = id.subdetId();
double eta = geo->getPosition(id).eta();
double eta = point.eta();
double theta = 2.0*atan(exp(-std::abs(eta)));
double dist = 0;
if (det == DetId::Ecal) {
Expand Down

0 comments on commit 28b86df

Please sign in to comment.