#include <RealisticSimClusterMapper.h>

Inheritance diagram for RealisticSimClusterMapper:

Public Member Functions
void	buildClusters (const edm::Handle< reco::PFRecHitCollection > &, const std::vector< bool > &, const std::vector< bool > &, reco::PFClusterCollection &) override

RealisticSimClusterMapper &	operator= (const RealisticSimClusterMapper &)=delete

	RealisticSimClusterMapper (const edm::ParameterSet &conf, edm::ConsumesCollector &sumes)

	RealisticSimClusterMapper (const RealisticSimClusterMapper &)=delete

void	update (const edm::EventSetup &) final

void	updateEvent (const edm::Event &) final

	~RealisticSimClusterMapper () override

Public Member Functions inherited from InitialClusteringStepBase
	InitialClusteringStepBase (const edm::ParameterSet &conf, edm::ConsumesCollector &sumes)

	InitialClusteringStepBase (const ICSB &)=delete

std::ostream &	operator<< (std::ostream &o) const

ICSB &	operator= (const ICSB &)=delete

void	reset ()

virtual	~InitialClusteringStepBase ()=default

Private Attributes
const float	calibMaxEta_ = 3.0

const float	calibMinEta_ = 1.4

std::vector< double >	egammaCalib_

const float	exclusiveFraction_ = 0.7f

std::vector< double >	hadronCalib_

const float	invisibleFraction_ = 0.3f

const float	maxDforTimingSquared_ = 4.0f

const float	maxDistance_ = 10.f

const bool	maxDistanceFilter_ = false

const unsigned int	minNHitsforTiming_ = 3

hgcal::RecHitTools	rhtools_

edm::Handle< SimClusterCollection >	simClusterH_

edm::EDGetTokenT< SimClusterCollection >	simClusterToken_

const float	timeOffset_

const bool	useMCFractionsForExclEnergy_ = false

Additional Inherited Members
Protected Types inherited from InitialClusteringStepBase
typedef std::tuple< std::vector< int >, std::vector< double >, std::vector< double > >	I3tuple

Protected Member Functions inherited from InitialClusteringStepBase
reco::PFRecHitRef	makeRefhit (const edm::Handle< reco::PFRecHitCollection > &h, const unsigned i) const

Protected Attributes inherited from InitialClusteringStepBase
const std::unordered_map< std::string, int >	_layerMap

unsigned	_nClustersFound

unsigned	_nSeeds

std::unordered_map< int, I3tuple >	_thresholds

Detailed Description

Definition at line 14 of file RealisticSimClusterMapper.h.

Constructor & Destructor Documentation

RealisticSimClusterMapper::RealisticSimClusterMapper	(	const edm::ParameterSet &	conf,
		edm::ConsumesCollector &	sumes
	)

inline

Definition at line 16 of file RealisticSimClusterMapper.h.

References edm::ConsumesCollector::consumes(), egammaCalib_, edm::ParameterSet::getParameter(), hadronCalib_, and simClusterToken_.

Referenced by ~RealisticSimClusterMapper().

       : InitialClusteringStepBase(conf, sumes),
         invisibleFraction_(conf.getParameter<double>("invisibleFraction")),
         exclusiveFraction_(conf.getParameter<double>("exclusiveFraction")),
         maxDistanceFilter_(conf.getParameter<bool>("maxDistanceFilter")),
         maxDistance_(conf.getParameter<double>("maxDistance")),
         maxDforTimingSquared_(conf.getParameter<double>("maxDforTimingSquared")),
         timeOffset_(conf.getParameter<double>("timeOffset")),
         minNHitsforTiming_(conf.getParameter<unsigned int>("minNHitsforTiming")),
         useMCFractionsForExclEnergy_(conf.getParameter<bool>("useMCFractionsForExclEnergy")),
         calibMinEta_(conf.getParameter<double>("calibMinEta")),
         calibMaxEta_(conf.getParameter<double>("calibMaxEta")) {
     simClusterToken_ = sumes.consumes<SimClusterCollection>(conf.getParameter<edm::InputTag>("simClusterSrc"));
     hadronCalib_ = conf.getParameter<std::vector<double> >("hadronCalib");
     egammaCalib_ = conf.getParameter<std::vector<double> >("egammaCalib");
   }

RealisticSimClusterMapper::~RealisticSimClusterMapper ( )

inlineoverride

Definition at line 33 of file RealisticSimClusterMapper.h.

References buildClusters(), operator=(), RealisticSimClusterMapper(), update(), and updateEvent().

33 {}

RealisticSimClusterMapper::RealisticSimClusterMapper ( const RealisticSimClusterMapper & )

delete

Member Function Documentation

void RealisticSimClusterMapper::buildClusters	(	const edm::Handle< reco::PFRecHitCollection > &	input,
		const std::vector< bool > &	rechitMask,
		const std::vector< bool > &	seedable,
		reco::PFClusterCollection &	output
	)

overridevirtual

Implements InitialClusteringStepBase.

Definition at line 48 of file RealisticSimClusterMapper.cc.

Referenced by ~RealisticSimClusterMapper().

                                                                                {
   const SimClusterCollection& simClusters = *simClusterH_;
   auto const& hits = *input;
   RealisticHitToClusterAssociator realisticAssociator;
   const int numberOfLayers = rhtools_.getGeometryType() == 0 ? rhtools_.getLayer(ForwardSubdetector::ForwardEmpty)
                                                              : rhtools_.getLayer(DetId::Forward);
   realisticAssociator.init(hits.size(), simClusters.size(), numberOfLayers + 1);
   // for quick indexing back to hit energy
   std::unordered_map<uint32_t, size_t> detIdToIndex(hits.size());
   for (uint32_t i = 0; i < hits.size(); ++i) {
     detIdToIndex[hits[i].detId()] = i;
     auto ref = makeRefhit(input, i);
     const auto& hitPos = rhtools_.getPosition(ref->detId());
 
     realisticAssociator.insertHitPosition(hitPos.x(), hitPos.y(), hitPos.z(), i);
     realisticAssociator.insertHitEnergy(ref->energy(), i);
     realisticAssociator.insertLayerId(rhtools_.getLayerWithOffset(ref->detId()), i);
   }
 
   for (unsigned int ic = 0; ic < simClusters.size(); ++ic) {
     const auto& sc = simClusters[ic];
     const auto& hitsAndFractions = sc.hits_and_fractions();
     for (const auto& hAndF : hitsAndFractions) {
       auto itr = detIdToIndex.find(hAndF.first);
       if (itr == detIdToIndex.end()) {
         continue;  // hit wasn't saved in reco or did not pass the SNR threshold
       }
       auto hitId = itr->second;
       auto ref = makeRefhit(input, hitId);
       float fraction = hAndF.second;
       float associatedEnergy = fraction * ref->energy();
       realisticAssociator.insertSimClusterIdAndFraction(ic, fraction, hitId, associatedEnergy);
     }
   }
   realisticAssociator.computeAssociation(
       exclusiveFraction_, useMCFractionsForExclEnergy_, rhtools_.lastLayerEE(), rhtools_.lastLayerFH());
   realisticAssociator.findAndMergeInvisibleClusters(invisibleFraction_, exclusiveFraction_);
   realisticAssociator.findCentersOfGravity();
   if (maxDistanceFilter_)
     realisticAssociator.filterHitsByDistance(maxDistance_);
 
   const auto& realisticClusters = realisticAssociator.realisticClusters();
   unsigned int nClusters = realisticClusters.size();
   float bin_norm = 1. / (calibMaxEta_ - calibMinEta_);
   float egamma_bin_norm = egammaCalib_.size() * bin_norm;
   float hadron_bin_norm = hadronCalib_.size() * bin_norm;
   for (unsigned ic = 0; ic < nClusters; ++ic) {
     float highest_energy = 0.0f;
     output.emplace_back();
     reco::PFCluster& back = output.back();
     edm::Ref<std::vector<reco::PFRecHit> > seed;
     float energyCorrection = 1.f;
     float timeRealisticSC = -99.;
     if (realisticClusters[ic].isVisible()) {
       int pdgId = simClusters[ic].pdgId();
       auto abseta = std::abs(simClusters[ic].eta());
       if ((abseta >= calibMinEta_) and (abseta <= calibMaxEta_))  //protecting range
       {
         if ((isEGamma(pdgId) or isPi0(pdgId)) and !egammaCalib_.empty()) {
           unsigned int etabin = std::floor(((abseta - calibMinEta_) * egamma_bin_norm));
           energyCorrection = egammaCalib_[etabin];
         } else if (isHadron(pdgId) and !(isPi0(pdgId)) and
                    !hadronCalib_
                         .empty())  // this function is expensive.. should we treat as hadron everything which is not egamma?
         {
           unsigned int etabin = std::floor(((abseta - calibMinEta_) * hadron_bin_norm));
           energyCorrection = hadronCalib_[etabin];
         }
       }
       std::vector<float> timeHits;
       const auto& hitsIdsAndFractions = realisticClusters[ic].hitsIdsAndFractions();
       for (const auto& idAndF : hitsIdsAndFractions) {
         auto fraction = idAndF.second;
         if (fraction > 0.f) {
           auto ref = makeRefhit(input, idAndF.first);
           back.addRecHitFraction(reco::PFRecHitFraction(ref, fraction));
           const float hit_energy = fraction * ref->energy();
           if (hit_energy > highest_energy || highest_energy == 0.0) {
             highest_energy = hit_energy;
             seed = ref;
           }
           //select hits good for timing
           if (ref->time() > -1.) {
             int rhLayer = rhtools_.getLayerWithOffset(ref->detId());
             std::array<float, 3> scPosition = realisticClusters[ic].getCenterOfGravity(rhLayer);
             float distanceSquared =
                 std::pow((ref->position().x() - scPosition[0]), 2) + std::pow((ref->position().y() - scPosition[1]), 2);
             if (distanceSquared < maxDforTimingSquared_) {
               timeHits.push_back(ref->time() - timeOffset_);
             }
           }
         }
       }
       //assign time if minimum number of hits
       if (timeHits.size() >= minNHitsforTiming_)
         timeRealisticSC = hgcalsimclustertime::fixSizeHighestDensity(timeHits);
     }
     if (!back.hitsAndFractions().empty()) {
       back.setSeed(seed->detId());
       back.setEnergy(realisticClusters[ic].getEnergy());
       back.setCorrectedEnergy(energyCorrection * realisticClusters[ic].getEnergy());  //applying energy correction
       back.setTime(timeRealisticSC);
     } else {
       back.setSeed(-1);
       back.setEnergy(0.f);
     }
   }
 }