CMS 3D CMS Logo
List of all members | Public Member Functions | Private Member Functions | Private Attributes | Friends
PFAlgo Class Reference

#include <PFAlgo.h>

Public Member Functions

void checkCleaning (const reco::PFRecHitCollection &cleanedHF)
 Check HF Cleaning. More...
 
reco::PFCandidateCollection & getCleanedCandidates ()
 
PFMuonAlgo * getPFMuonAlgo ()
 
reco::PFCandidateCollection makeConnectedCandidates ()
 
 PFAlgo (double nSigmaECAL, double nSigmaHCAL, double nSigmaHFEM, double nSigmaHFHAD, std::vector< double > resolHF_square, PFEnergyCalibration &calibration, PFEnergyCalibrationHF &thepfEnergyCalibrationHF, const edm::ParameterSet &pset)
 constructor More...
 
void reconstructParticles (const reco::PFBlockHandle &blockHandle, PFEGammaFilters const *pfegamma)
 reconstruct particles More...
 
void setCandConnectorParameters (bool bCorrect, bool bCalibPrimary, double dptRel_PrimaryTrack, double dptRel_MergedTrack, double ptErrorSecondary, const std::vector< double > &nuclCalibFactors)
 
void setCandConnectorParameters (const edm::ParameterSet &iCfgCandConnector)
 
void setDisplacedVerticesParameters (bool rejectTracks_Bad, bool rejectTracks_Step45, bool usePFNuclearInteractions, bool usePFConversions, bool usePFDecays, double dptRel_DispVtx)
 
void setEGammaCollections (const edm::View< reco::PFCandidate > &pfEgammaCandidates, const edm::ValueMap< reco::GsfElectronRef > &valueMapGedElectrons, const edm::ValueMap< reco::PhotonRef > &valueMapGedPhotons)
 
void setEGammaParameters (bool use_EGammaFilters, bool useProtectionsForJetMET)
 
void setEGElectronCollection (const reco::GsfElectronCollection &egelectrons)
 
void setHOTag (bool ho)
 
void setMuonHandle (const edm::Handle< reco::MuonCollection > &muons)
 
void setPFVertexParameters (bool useVertex, reco::VertexCollection const &primaryVertices)
 
void setPostHFCleaningParameters (bool postHFCleaning, const edm::ParameterSet &pfHFCleaningParams)
 

Private Member Functions

void associatePSClusters (unsigned iEcal, reco::PFBlockElement::Type psElementType, const reco::PFBlock &block, const edm::OwnVector< reco::PFBlockElement > &elements, const reco::PFBlock::LinkData &linkData, std::vector< bool > &active, std::vector< double > &psEne)
 Associate PS clusters to a given ECAL cluster, and return their energy. More...
 
bool checkAndReconstructSecondaryInteraction (const reco::PFBlockRef &blockref, const edm::OwnVector< reco::PFBlockElement > &elements, bool isActive, int iElement)
 
bool checkGoodTrackDeadHcal (const reco::TrackRef &trackRef, bool hasDeadHcal)
 
bool checkHasDeadHcal (const std::multimap< double, unsigned > &hcalElems, const std::vector< bool > &deadArea)
 
void conversionAlgo (const edm::OwnVector< reco::PFBlockElement > &elements, std::vector< bool > &active)
 
void createCandidatesECAL (const reco::PFBlock &block, reco::PFBlock::LinkData &linkData, const edm::OwnVector< reco::PFBlockElement > &elements, std::vector< bool > &active, const reco::PFBlockRef &blockref, ElementIndices &inds, std::vector< bool > &deadArea)
 
void createCandidatesHCAL (const reco::PFBlock &block, reco::PFBlock::LinkData &linkData, const edm::OwnVector< reco::PFBlockElement > &elements, std::vector< bool > &active, const reco::PFBlockRef &blockref, ElementIndices &inds, std::vector< bool > &deadArea)
 
void createCandidatesHCALUnlinked (const reco::PFBlock &block, reco::PFBlock::LinkData &linkData, const edm::OwnVector< reco::PFBlockElement > &elements, std::vector< bool > &active, const reco::PFBlockRef &blockref, ElementIndices &inds, std::vector< bool > &deadArea)
 
void createCandidatesHF (const reco::PFBlock &block, reco::PFBlock::LinkData &linkData, const edm::OwnVector< reco::PFBlockElement > &elements, std::vector< bool > &active, const reco::PFBlockRef &blockref, ElementIndices &inds)
 
int decideType (const edm::OwnVector< reco::PFBlockElement > &elements, const reco::PFBlockElement::Type type, std::vector< bool > &active, ElementIndices &inds, std::vector< bool > &deadArea, unsigned int iEle)
 
void egammaFilters (const reco::PFBlockRef &blockref, std::vector< bool > &active, PFEGammaFilters const *pfegamma)
 
void elementLoop (const reco::PFBlock &block, reco::PFBlock::LinkData &linkData, const edm::OwnVector< reco::PFBlockElement > &elements, std::vector< bool > &active, const reco::PFBlockRef &blockref, ElementIndices &inds, std::vector< bool > &deadArea)
 
double hfEnergyResolution (double clusterEnergy) const
 
bool isFromSecInt (const reco::PFBlockElement &eTrack, std::string order) const
 
double neutralHadronEnergyResolution (double clusterEnergy, double clusterEta) const
 todo: use PFClusterTools for this More...
 
double nSigmaHCAL (double clusterEnergy, double clusterEta) const
 
double nSigmaHFEM (double clusterEnergy) const
 
double nSigmaHFHAD (double clusterEnergy) const
 
void postCleaning ()
 
void processBlock (const reco::PFBlockRef &blockref, std::list< reco::PFBlockRef > &hcalBlockRefs, std::list< reco::PFBlockRef > &ecalBlockRefs, PFEGammaFilters const *pfegamma)
 
unsigned reconstructCluster (const reco::PFCluster &cluster, double particleEnergy, bool useDirection=false, double particleX=0., double particleY=0., double particleZ=0.)
 
unsigned reconstructTrack (const reco::PFBlockElement &elt, bool allowLoose=false)
 
bool recoTracksNotHCAL (const reco::PFBlock &block, reco::PFBlock::LinkData &linkData, const edm::OwnVector< reco::PFBlockElement > &elements, const reco::PFBlockRef &blockref, std::vector< bool > &active, bool goodTrackDeadHcal, bool hasDeadHcal, unsigned int iTrack, std::multimap< double, unsigned > &ecalElems, reco::TrackRef &trackRef)
 
void relinkTrackToHcal (const reco::PFBlock &block, std::multimap< double, unsigned > &ecalElems, std::multimap< double, unsigned > &hcalElems, const std::vector< bool > &active, reco::PFBlock::LinkData &linkData, unsigned int iTrack)
 
void setHcalDepthInfo (reco::PFCandidate &cand, const reco::PFCluster &cluster) const
 

Private Attributes

PFEnergyCalibration & calibration_
 
PFCandConnector connector_
 
double dptRel_DispVtx_
 
std::vector< double > factors45_
 
float goodPixelTrackDeadHcal_chi2n_
 
float goodPixelTrackDeadHcal_dxy_
 
float goodPixelTrackDeadHcal_dz_
 
int goodPixelTrackDeadHcal_maxLost3Hit_
 
int goodPixelTrackDeadHcal_maxLost4Hit_
 
float goodPixelTrackDeadHcal_maxPt_
 
float goodPixelTrackDeadHcal_minEta_
 
float goodPixelTrackDeadHcal_ptErrRel_
 
float goodTrackDeadHcal_chi2n_
 
float goodTrackDeadHcal_dxy_
 
int goodTrackDeadHcal_layers_
 
float goodTrackDeadHcal_ptErrRel_
 Variables for track cleaning in bad HCal areas. More...
 
float goodTrackDeadHcal_validFr_
 
double maxDeltaPhiPt_
 
double maxSignificance_
 
double minDeltaMet_
 
double minHFCleaningPt_
 
double minSignificance_
 
double minSignificanceReduction_
 
std::vector< double > muonECAL_
 
edm::Handle< reco::MuonCollection > muonHandle_
 
std::vector< double > muonHCAL_
 Variables for muons and fakes. More...
 
std::vector< double > muonHO_
 
const double nSigmaECAL_
 number of sigma to judge energy excess in ECAL More...
 
const double nSigmaEConstHCAL = 100.
 
const double nSigmaEConstHFEM = 100.
 
const double nSigmaEConstHFHAD = 100.
 
const double nSigmaHCAL_
 number of sigma to judge energy excess in HCAL More...
 
const double nSigmaHFEM_
 number of sigma to judge energy excess in HF More...
 
const double nSigmaHFHAD_
 
double nSigmaTRACK_
 
int nVtx_
 
std::unique_ptr< reco::PFCandidateCollection > pfCandidates_
 
reco::PFCandidateCollection pfCleanedCandidates_
 
const edm::View< reco::PFCandidate > * pfEgammaCandidates_
 
std::unique_ptr< PFMuonAlgo > pfmu_
 
bool postHFCleaning_
 
bool postMuonCleaning_
 
reco::Vertex primaryVertex_
 
double ptError_
 
bool rejectTracks_Bad_
 
bool rejectTracks_Step45_
 
const std::vector< double > resolHF_square_
 
PFEnergyCalibrationHF & thepfEnergyCalibrationHF_
 
double useBestMuonTrack_
 
bool useEGammaFilters_
 Variables for NEW EGAMMA selection. More...
 
bool useHO_
 
bool usePFConversions_
 
bool usePFDecays_
 
bool usePFMuonMomAssign_
 
bool usePFNuclearInteractions_
 
bool useProtectionsForJetMET_
 
bool useVertices_ = false
 
const edm::ValueMap< reco::GsfElectronRef > * valueMapGedElectrons_
 
const edm::ValueMap< reco::PhotonRef > * valueMapGedPhotons_
 

Friends

std::ostream & operator<< (std::ostream &out, const PFAlgo &algo)
 

Detailed Description

Definition at line 53 of file PFAlgo.h.

Constructor & Destructor Documentation

◆ PFAlgo()

PFAlgo::PFAlgo ( double  nSigmaECAL,
double  nSigmaHCAL,
double  nSigmaHFEM,
double  nSigmaHFHAD,
std::vector< double >  resolHF_square,
PFEnergyCalibration &  calibration,
PFEnergyCalibrationHF &  thepfEnergyCalibrationHF,
const edm::ParameterSet &  pset 
)

constructor

Definition at line 15 of file PFAlgo.cc.

23  : pfCandidates_(new PFCandidateCollection),
24  nSigmaECAL_(nSigmaECAL),
25  nSigmaHCAL_(nSigmaHCAL),
26  nSigmaHFEM_(nSigmaHFEM),
27  nSigmaHFHAD_(nSigmaHFHAD),
28  resolHF_square_(resolHF_square),
29  calibration_(calibration),
30  thepfEnergyCalibrationHF_(thepfEnergyCalibrationHF),
31  connector_() {
32  const edm::ParameterSet pfMuonAlgoParams = pset.getParameter<edm::ParameterSet>("PFMuonAlgoParameters");
33  bool postMuonCleaning = pset.getParameter<bool>("postMuonCleaning");
34  pfmu_ = std::make_unique<PFMuonAlgo>(pfMuonAlgoParams, postMuonCleaning);
35 
36  // HF resolution parameters
37  assert(resolHF_square_.size() == 3); // make sure that stochastic, constant, noise (i.e. three) terms are specified.
38 
39  // Muon parameters
40  muonHCAL_ = pset.getParameter<std::vector<double>>("muon_HCAL");
41  muonECAL_ = pset.getParameter<std::vector<double>>("muon_ECAL");
42  muonHO_ = pset.getParameter<std::vector<double>>("muon_HO");
43  assert(muonHCAL_.size() == 2 && muonECAL_.size() == 2 && muonHO_.size() == 2);
44  nSigmaTRACK_ = pset.getParameter<double>("nsigma_TRACK");
45  ptError_ = pset.getParameter<double>("pt_Error");
46  factors45_ = pset.getParameter<std::vector<double>>("factors_45");
47  assert(factors45_.size() == 2);
48 
49  // Bad Hcal Track Parameters
50  goodTrackDeadHcal_ptErrRel_ = pset.getParameter<double>("goodTrackDeadHcal_ptErrRel");
51  goodTrackDeadHcal_chi2n_ = pset.getParameter<double>("goodTrackDeadHcal_chi2n");
52  goodTrackDeadHcal_layers_ = pset.getParameter<uint32_t>("goodTrackDeadHcal_layers");
53  goodTrackDeadHcal_validFr_ = pset.getParameter<double>("goodTrackDeadHcal_validFr");
54  goodTrackDeadHcal_dxy_ = pset.getParameter<double>("goodTrackDeadHcal_dxy");
55 
56  goodPixelTrackDeadHcal_minEta_ = pset.getParameter<double>("goodPixelTrackDeadHcal_minEta");
57  goodPixelTrackDeadHcal_maxPt_ = pset.getParameter<double>("goodPixelTrackDeadHcal_maxPt");
58  goodPixelTrackDeadHcal_ptErrRel_ = pset.getParameter<double>("goodPixelTrackDeadHcal_ptErrRel");
59  goodPixelTrackDeadHcal_chi2n_ = pset.getParameter<double>("goodPixelTrackDeadHcal_chi2n");
60  goodPixelTrackDeadHcal_maxLost3Hit_ = pset.getParameter<int32_t>("goodPixelTrackDeadHcal_maxLost3Hit");
61  goodPixelTrackDeadHcal_maxLost4Hit_ = pset.getParameter<int32_t>("goodPixelTrackDeadHcal_maxLost4Hit");
62  goodPixelTrackDeadHcal_dxy_ = pset.getParameter<double>("goodPixelTrackDeadHcal_dxy");
63  goodPixelTrackDeadHcal_dz_ = pset.getParameter<double>("goodPixelTrackDeadHcal_dz");
64 }

References cms::cuda::assert(), factors45_, goodPixelTrackDeadHcal_chi2n_, goodPixelTrackDeadHcal_dxy_, goodPixelTrackDeadHcal_dz_, goodPixelTrackDeadHcal_maxLost3Hit_, goodPixelTrackDeadHcal_maxLost4Hit_, goodPixelTrackDeadHcal_maxPt_, goodPixelTrackDeadHcal_minEta_, goodPixelTrackDeadHcal_ptErrRel_, goodTrackDeadHcal_chi2n_, goodTrackDeadHcal_dxy_, goodTrackDeadHcal_layers_, goodTrackDeadHcal_ptErrRel_, goodTrackDeadHcal_validFr_, muonECAL_, muonHCAL_, muonHO_, nSigmaTRACK_, pfmu_, HLT_FULL_cff::postMuonCleaning, muonDTDigis_cfi::pset, ptError_, and resolHF_square_.

Member Function Documentation

◆ associatePSClusters()

void PFAlgo::associatePSClusters ( unsigned  iEcal,
reco::PFBlockElement::Type  psElementType,
const reco::PFBlock &  block,
const edm::OwnVector< reco::PFBlockElement > &  elements,
const reco::PFBlock::LinkData &  linkData,
std::vector< bool > &  active,
std::vector< double > &  psEne 
)
private

Associate PS clusters to a given ECAL cluster, and return their energy.

Definition at line 3417 of file PFAlgo.cc.

3423  {
3424  // Find all PS clusters with type psElement associated to ECAL cluster iEcal,
3425  // within all PFBlockElement "elements" of a given PFBlock "block"
3426  // psElement can be reco::PFBlockElement::PS1 or reco::PFBlockElement::PS2
3427  // Returns a vector of PS cluster energies, and updates the "active" vector.
3428 
3429  // Find all PS clusters linked to the iEcal cluster
3430  std::multimap<double, unsigned> sortedPS;
3431  block.associatedElements(iEcal, linkData, sortedPS, psElementType, reco::PFBlock::LINKTEST_ALL);
3432 
3433  // Loop over these PS clusters
3434  double totalPS = 0.;
3435  for (auto const& ps : sortedPS) {
3436  // CLuster index and distance to iEcal
3437  unsigned iPS = ps.second;
3438  // double distPS = ps.first;
3439 
3440  // Ignore clusters already in use
3441  if (!active[iPS])
3442  continue;
3443 
3444  // Check that this cluster is not closer to another ECAL cluster
3445  std::multimap<double, unsigned> sortedECAL;
3446  block.associatedElements(iPS, linkData, sortedECAL, reco::PFBlockElement::ECAL, reco::PFBlock::LINKTEST_ALL);
3447  unsigned jEcal = sortedECAL.begin()->second;
3448  if (jEcal != iEcal)
3449  continue;
3450 
3451  // Update PS energy
3452  PFBlockElement::Type pstype = elements[iPS].type();
3453  assert(pstype == psElementType);
3454  PFClusterRef psclusterref = elements[iPS].clusterRef();
3455  assert(!psclusterref.isNull());
3456  totalPS += psclusterref->energy();
3457  psEne[0] += psclusterref->energy();
3458  active[iPS] = false;
3459  }
3460 }

References cms::cuda::assert(), groupFilesInBlocks::block, reco::PFBlockElement::ECAL, bookConverter::elements, edm::Ref< C, T, F >::isNull(), and reco::PFBlock::LINKTEST_ALL.

Referenced by recoTracksNotHCAL().

◆ checkAndReconstructSecondaryInteraction()

bool PFAlgo::checkAndReconstructSecondaryInteraction ( const reco::PFBlockRef &  blockref,
const edm::OwnVector< reco::PFBlockElement > &  elements,
bool  isActive,
int  iElement 
)
private

Definition at line 847 of file PFAlgo.cc.

850  {
851  bool ret = isActive;
852  if (isActive && isFromSecInt(elements[iElement], "primary")) {
853  bool isPrimaryTrack =
854  elements[iElement].displacedVertexRef(PFBlockElement::T_TO_DISP)->displacedVertexRef()->isTherePrimaryTracks();
855  if (isPrimaryTrack) {
856  LogTrace("PFAlgo|elementLoop") << "Primary Track reconstructed alone";
857 
858  unsigned tmpi = reconstructTrack(elements[iElement]);
859  (*pfCandidates_)[tmpi].addElementInBlock(blockref, iElement);
860  ret = false;
861  }
862  }
863 
864  return ret;
865 }

References bookConverter::elements, isFromSecInt(), LogTrace, reconstructTrack(), and runTheMatrix::ret.

Referenced by elementLoop().

◆ checkCleaning()

void PFAlgo::checkCleaning ( const reco::PFRecHitCollection &  cleanedHF)

Check HF Cleaning.

Definition at line 3587 of file PFAlgo.cc.

3587  {
3588  // No hits to recover, leave.
3589  if (cleanedHits.empty())
3590  return;
3591 
3592  //Compute met and met significance (met/sqrt(SumEt))
3593  double metX = 0.;
3594  double metY = 0.;
3595  double sumet = 0;
3596  std::vector<unsigned int> hitsToBeAdded;
3597  for (auto const& pfc : *pfCandidates_) {
3598  metX += pfc.px();
3599  metY += pfc.py();
3600  sumet += pfc.pt();
3601  }
3602  double met2 = metX * metX + metY * metY;
3603  double met2_Original = met2;
3604  // Select events with large MET significance.
3605  // double significance = std::sqrt(met2/sumet);
3606  // double significanceCor = significance;
3607  double metXCor = metX;
3608  double metYCor = metY;
3609  double sumetCor = sumet;
3610  double met2Cor = met2;
3611  bool next = true;
3612  unsigned iCor = 1E9;
3613 
3614  // Find the cleaned hit with the largest effect on the MET
3615  while (next) {
3616  double metReduc = -1.;
3617  // Loop on the candidates
3618  for (unsigned i = 0; i < cleanedHits.size(); ++i) {
3619  const PFRecHit& hit = cleanedHits[i];
3620  double length = std::sqrt(hit.position().mag2());
3621  double px = hit.energy() * hit.position().x() / length;
3622  double py = hit.energy() * hit.position().y() / length;
3623  double pt = std::sqrt(px * px + py * py);
3624 
3625  // Check that it is not already scheculed to be cleaned
3626  bool skip = false;
3627  for (unsigned int hitIdx : hitsToBeAdded) {
3628  if (i == hitIdx)
3629  skip = true;
3630  if (skip)
3631  break;
3632  }
3633  if (skip)
3634  continue;
3635 
3636  // Now add the candidate to the MET
3637  double metXInt = metX + px;
3638  double metYInt = metY + py;
3639  double sumetInt = sumet + pt;
3640  double met2Int = metXInt * metXInt + metYInt * metYInt;
3641 
3642  // And check if it could contribute to a MET reduction
3643  if (met2Int < met2Cor) {
3644  metXCor = metXInt;
3645  metYCor = metYInt;
3646  metReduc = (met2 - met2Int) / met2Int;
3647  met2Cor = met2Int;
3648  sumetCor = sumetInt;
3649  // significanceCor = std::sqrt(met2Cor/sumetCor);
3650  iCor = i;
3651  }
3652  }
3653  //
3654  // If the MET must be significanly reduced, schedule the candidate to be added
3655  //
3656  if (metReduc > minDeltaMet_) {
3657  hitsToBeAdded.push_back(iCor);
3658  metX = metXCor;
3659  metY = metYCor;
3660  sumet = sumetCor;
3661  met2 = met2Cor;
3662  } else {
3663  // Otherwise just stop the loop
3664  next = false;
3665  }
3666  }
3667  //
3668  // At least 10 GeV MET reduction
3669  if (std::sqrt(met2_Original) - std::sqrt(met2) > 5.) {
3670  LogTrace("PFAlgo|checkCleaning") << hitsToBeAdded.size() << " hits were re-added ";
3671  LogTrace("PFAlgo|checkCleaning") << "MET reduction = " << std::sqrt(met2_Original) << " -> " << std::sqrt(met2Cor)
3672  << " = " << std::sqrt(met2Cor) - std::sqrt(met2_Original);
3673  LogTrace("PFAlgo|checkCleaning") << "Added after cleaning check : ";
3674  for (unsigned int hitIdx : hitsToBeAdded) {
3675  const PFRecHit& hit = cleanedHits[hitIdx];
3676  PFCluster cluster(hit.layer(), hit.energy(), hit.position().x(), hit.position().y(), hit.position().z());
3677  reconstructCluster(cluster, hit.energy());
3678  LogTrace("PFAlgo|checkCleaning") << pfCandidates_->back() << ". time = " << hit.time();
3679  }
3680  }
3681 } //PFAlgo::checkCleaning

References mps_fire::i, LogTrace, minDeltaMet_, GetRecoTauVFromDQM_MC_cff::next, pfCandidates_, DiDispStaMuonMonitor_cfi::pt, multPhiCorr_741_25nsDY_cfi::px, multPhiCorr_741_25nsDY_cfi::py, reconstructCluster(), optionsL1T::skip, mathSSE::sqrt(), hit::x, hit::y, and hit::z.

Referenced by PFProducer::produce().

◆ checkGoodTrackDeadHcal()

bool PFAlgo::checkGoodTrackDeadHcal ( const reco::TrackRef &  trackRef,
bool  hasDeadHcal 
)
private

Definition at line 894 of file PFAlgo.cc.

894  {
895  bool goodTrackDeadHcal = false;
896  if (hasDeadHcal) {
897  goodTrackDeadHcal = (trackRef->ptError() < goodTrackDeadHcal_ptErrRel_ * trackRef->pt() &&
898  trackRef->normalizedChi2() < goodTrackDeadHcal_chi2n_ &&
899  trackRef->hitPattern().trackerLayersWithMeasurement() >= goodTrackDeadHcal_layers_ &&
900  trackRef->validFraction() > goodTrackDeadHcal_validFr_ &&
901  std::abs(trackRef->dxy(primaryVertex_.position())) < goodTrackDeadHcal_dxy_);
902  // now we add an extra block for tracks at high |eta|
903  if (!goodTrackDeadHcal && std::abs(trackRef->eta()) > goodPixelTrackDeadHcal_minEta_ && // high eta
904  trackRef->hitPattern().pixelLayersWithMeasurement() >= 3 && // pixel track
905  trackRef->hitPattern().trackerLayersWithoutMeasurement(HitPattern::TRACK_HITS) == 0 &&
906  trackRef->hitPattern().trackerLayersWithoutMeasurement(HitPattern::MISSING_INNER_HITS) == 0 &&
907  trackRef->hitPattern().trackerLayersWithoutMeasurement(HitPattern::MISSING_OUTER_HITS) <=
908  (trackRef->hitPattern().pixelLayersWithMeasurement() > 3 ? goodPixelTrackDeadHcal_maxLost4Hit_
909  : goodPixelTrackDeadHcal_maxLost3Hit_) &&
910  trackRef->normalizedChi2() < goodPixelTrackDeadHcal_chi2n_ && // tighter cut
911  std::abs(trackRef->dxy(primaryVertex_.position())) < goodPixelTrackDeadHcal_dxy_ &&
912  std::abs(trackRef->dz(primaryVertex_.position())) < goodPixelTrackDeadHcal_dz_ &&
913  trackRef->ptError() < goodPixelTrackDeadHcal_ptErrRel_ * trackRef->pt() && // sanity
914  trackRef->pt() < goodPixelTrackDeadHcal_maxPt_) { // sanity
915  goodTrackDeadHcal = true;
916  // FIXME: may decide to do something to the track pT
917  }
918  //if (!goodTrackDeadHcal && trackRef->hitPattern().trackerLayersWithMeasurement() == 4 && trackRef->validFraction() == 1
919  LogTrace("PFAlgo|elementLoop")
920  << " track pt " << trackRef->pt() << " +- " << trackRef->ptError() << " layers valid "
921  << trackRef->hitPattern().trackerLayersWithMeasurement() << ", lost "
922  << trackRef->hitPattern().trackerLayersWithoutMeasurement(HitPattern::TRACK_HITS) << ", lost outer "
923  << trackRef->hitPattern().trackerLayersWithoutMeasurement(HitPattern::MISSING_OUTER_HITS) << ", lost inner "
924  << trackRef->hitPattern().trackerLayersWithoutMeasurement(HitPattern::MISSING_INNER_HITS) << "(valid fraction "
925  << trackRef->validFraction() << ")"
926  << " chi2/ndf " << trackRef->normalizedChi2() << " |dxy| " << std::abs(trackRef->dxy(primaryVertex_.position()))
927  << " +- " << trackRef->dxyError() << " |dz| " << std::abs(trackRef->dz(primaryVertex_.position())) << " +- "
928  << trackRef->dzError() << (goodTrackDeadHcal ? " passes " : " fails ") << "quality cuts";
929  }
930  return goodTrackDeadHcal;
931 }

References funct::abs(), goodPixelTrackDeadHcal_chi2n_, goodPixelTrackDeadHcal_dxy_, goodPixelTrackDeadHcal_dz_, goodPixelTrackDeadHcal_maxLost3Hit_, goodPixelTrackDeadHcal_maxLost4Hit_, goodPixelTrackDeadHcal_maxPt_, goodPixelTrackDeadHcal_minEta_, goodPixelTrackDeadHcal_ptErrRel_, goodTrackDeadHcal_chi2n_, goodTrackDeadHcal_dxy_, goodTrackDeadHcal_layers_, goodTrackDeadHcal_ptErrRel_, goodTrackDeadHcal_validFr_, LogTrace, reco::Vertex::position(), and primaryVertex_.

Referenced by elementLoop().

◆ checkHasDeadHcal()

bool PFAlgo::checkHasDeadHcal ( const std::multimap< double, unsigned > &  hcalElems,
const std::vector< bool > &  deadArea 
)
private

Definition at line 867 of file PFAlgo.cc.

867  {
868  // there's 3 possible options possible here, in principle:
869  // 1) flag everything that may be associated to a dead hcal marker
870  // 2) flag everything whose closest hcal link is a dead hcal marker
871  // 3) flag only things that are linked only to dead hcal marker
872  // in our first test we go for (2)
873  //--- option (1) --
874  //bool hasDeadHcal = false;
875  //for (auto it = hcalElems.begin(), ed = hcalElems.end(); it != ed; /*NOTE NO ++it HERE */ ) {
876  // if (deadArea[it->second]) { hasDeadHcal = true; it = hcalElems.erase(it); } // std::multimap::erase returns iterator to next
877  // else ++it;
878  //}
879  //--- option (2) --
880  bool hasDeadHcal = false;
881  if (!hcalElems.empty() && deadArea[hcalElems.begin()->second]) {
882  hasDeadHcal = true;
883  }
884  //--- option (3) --
885  //bool hasDeadHcal = true;
886  //for (auto it = hcalElems.begin(), ed = hcalElems.end(); it != ed; /*NOTE NO ++it HERE */ ) {
887  // if (deadArea[it->second]) { it = hcalElems.erase(it); } // std::multimap::erase returns iterator to next
888  // else { hasDeadHcal = false; }
889  //}
890  return hasDeadHcal;
891 }

Referenced by elementLoop().

◆ conversionAlgo()

void PFAlgo::conversionAlgo ( const edm::OwnVector< reco::PFBlockElement > &  elements,
std::vector< bool > &  active 
)
private

Definition at line 362 of file PFAlgo.cc.

362  {
363  LogTrace("PFAlgo|conversionAlgo") << "start of function PFAlgo::conversionAlgo(), elements.size()="
364  << elements.size();
365  for (unsigned iEle = 0; iEle < elements.size(); iEle++) {
366  PFBlockElement::Type type = elements[iEle].type();
367  if (type == PFBlockElement::TRACK) {
368  LogTrace("PFAlgo|conversionAlgo") << "elements[" << iEle << "].type() == TRACK, active[" << iEle
369  << "]=" << active[iEle];
370  if (elements[iEle].trackRef()->algo() == reco::TrackBase::conversionStep) {
371  active[iEle] = false;
372  }
373  if (elements[iEle].trackRef()->quality(reco::TrackBase::highPurity)) {
374  LogTrace("PFAlgo|conversionAlgo") << "Track is high purity";
375  continue;
376  }
377  const auto* trackRef = dynamic_cast<const reco::PFBlockElementTrack*>((&elements[iEle]));
378  if (!(trackRef->trackType(reco::PFBlockElement::T_FROM_GAMMACONV))) {
379  LogTrace("PFAlgo|conversionAlgo") << "!trackType(T_FROM_GAMMACONV)";
380  continue;
381  }
382  if (!elements[iEle].convRefs().empty()) {
383  active[iEle] = false;
384  }
385  LogTrace("PFAlgo|conversionAlgo") << "active[iEle=" << iEle << "]=" << active[iEle];
386  }
387  }
388  LogTrace("PFAlgo|conversionAlgo") << "end of function PFAlgo::conversionAlgo";
389 }

References reco::TrackBase::conversionStep, bookConverter::elements, relativeConstraints::empty, reco::TrackBase::highPurity, LogTrace, quality, and reco::PFBlockElement::T_FROM_GAMMACONV.

Referenced by processBlock().

◆ createCandidatesECAL()

void PFAlgo::createCandidatesECAL ( const reco::PFBlock &  block,
reco::PFBlock::LinkData &  linkData,
const edm::OwnVector< reco::PFBlockElement > &  elements,
std::vector< bool > &  active,
const reco::PFBlockRef &  blockref,
ElementIndices &  inds,
std::vector< bool > &  deadArea 
)
private

Definition at line 3004 of file PFAlgo.cc.

3010  {
3011  LogTrace("PFAlgo|createCandidatesECAL")
3012  << "start of function PFAlgo::createCandidatesECAL(), ecalIs.size()=" << inds.ecalIs.size();
3013 
3014  // --------------- loop ecal ------------------
3015 
3016  // for each ecal element iEcal = ecalIs[i] in turn:
3017 
3018  for (unsigned i = 0; i < inds.ecalIs.size(); i++) {
3019  unsigned iEcal = inds.ecalIs[i];
3020 
3021  LogTrace("PFAlgo|createCandidatesECAL") << "elements[" << iEcal << "]=" << elements[iEcal];
3022 
3023  if (!active[iEcal]) {
3024  LogTrace("PFAlgo|createCandidatesECAL") << "iEcal=" << iEcal << " not active";
3025  continue;
3026  }
3027 
3028  PFBlockElement::Type type = elements[iEcal].type();
3029  assert(type == PFBlockElement::ECAL);
3030 
3031  PFClusterRef clusterref = elements[iEcal].clusterRef();
3032  assert(!clusterref.isNull());
3033 
3034  active[iEcal] = false;
3035 
3036  float ecalEnergy = clusterref->correctedEnergy();
3037  // float ecalEnergy = calibration_.energyEm( clusterref->energy() );
3038  double particleEnergy = ecalEnergy;
3039 
3040  auto& cand = (*pfCandidates_)[reconstructCluster(*clusterref, particleEnergy)];
3041 
3042  cand.setEcalEnergy(clusterref->energy(), ecalEnergy);
3043  cand.setHcalEnergy(0., 0.);
3044  cand.setHoEnergy(0., 0.);
3045  cand.setPs1Energy(0.);
3046  cand.setPs2Energy(0.);
3047  cand.addElementInBlock(blockref, iEcal);
3048 
3049  } // end loop on ecal elements iEcal = ecalIs[i]
3050  LogTrace("PFAlgo|createCandidatesECAL") << "end of function PFALgo::createCandidatesECAL";
3051 }

References cms::cuda::assert(), ECAL, ElementIndices::ecalIs, bookConverter::elements, mps_fire::i, edm::Ref< C, T, F >::isNull(), LogTrace, and reconstructCluster().

Referenced by processBlock().

◆ createCandidatesHCAL()

void PFAlgo::createCandidatesHCAL ( const reco::PFBlock &  block,
reco::PFBlock::LinkData &  linkData,
const edm::OwnVector< reco::PFBlockElement > &  elements,
std::vector< bool > &  active,
const reco::PFBlockRef &  blockref,
ElementIndices &  inds,
std::vector< bool > &  deadArea 
)
private

Definition at line 1694 of file PFAlgo.cc.

1700  {
1701  LogTrace("PFAlgo|createCandidatesHCAL")
1702  << "start of function PFAlgo::createCandidatesHCAL, inds.hcalIs.size()=" << inds.hcalIs.size();
1703 
1704  // --------------- loop hcal ------------------
1705 
1706  for (unsigned iHcal : inds.hcalIs) {
1707  PFBlockElement::Type type = elements[iHcal].type();
1708 
1709  assert(type == PFBlockElement::HCAL);
1710 
1711  LogTrace("PFAlgo|createCandidatesHCAL") << "elements[" << iHcal << "]=" << elements[iHcal];
1712 
1713  // associated tracks
1714  std::multimap<double, unsigned> sortedTracks;
1715  block.associatedElements(iHcal, linkData, sortedTracks, reco::PFBlockElement::TRACK, reco::PFBlock::LINKTEST_ALL);
1716 
1717  std::multimap<unsigned, std::pair<double, unsigned>> associatedEcals;
1718 
1719  std::map<unsigned, std::pair<double, double>> associatedPSs;
1720 
1721  std::multimap<double, std::pair<unsigned, bool>> associatedTracks;
1722 
1723  // A temporary maps for ECAL satellite clusters
1724  std::multimap<double, std::tuple<unsigned, ::math::XYZVector, double>>
1725  ecalSatellites; // last element (double) : correction under the egamma hypothesis
1726  std::tuple<unsigned, ::math::XYZVector, double> fakeSatellite(iHcal, ::math::XYZVector(0., 0., 0.), 1.);
1727  ecalSatellites.emplace(-1., fakeSatellite);
1728 
1729  std::multimap<unsigned, std::pair<double, unsigned>> associatedHOs;
1730 
1731  PFClusterRef hclusterref = elements[iHcal].clusterRef();
1732  assert(!hclusterref.isNull());
1733 
1734  //if there is no track attached to that HCAL, then do not
1735  //reconstruct an HCAL alone, check if it can be recovered
1736  //first
1737 
1738  // if no associated tracks, create a neutral hadron
1739  //if(sortedTracks.empty() ) {
1740 
1741  if (sortedTracks.empty()) {
1742  LogTrace("PFAlgo|createCandidatesHCAL") << "\tno associated tracks, keep for later";
1743  continue;
1744  }
1745 
1746  // Lock the HCAL cluster
1747  active[iHcal] = false;
1748 
1749  // in the following, tracks are associated to this hcal cluster.
1750  // will look for an excess of energy in the calorimeters w/r to
1751  // the charged energy, and turn this excess into a neutral hadron or
1752  // a photon.
1753 
1754  LogTrace("PFAlgo|createCandidatesHCAL") << sortedTracks.size() << " associated tracks:";
1755 
1756  double totalChargedMomentum = 0;
1757  double sumpError2 = 0.;
1758  double totalHO = 0.;
1759  double totalEcal = 0.;
1760  double totalEcalEGMCalib = 0.;
1761  double totalHcal = hclusterref->energy();
1762  vector<double> hcalP;
1763  vector<double> hcalDP;
1764  vector<unsigned> tkIs;
1765  double maxDPovP = -9999.;
1766 
1767  //Keep track of how much energy is assigned to calorimeter-vs-track energy/momentum excess
1768  vector<unsigned> chargedHadronsIndices;
1769  vector<unsigned> chargedHadronsInBlock;
1770  double mergedNeutralHadronEnergy = 0;
1771  double mergedPhotonEnergy = 0;
1772  double muonHCALEnergy = 0.;
1773  double muonECALEnergy = 0.;
1774  double muonHCALError = 0.;
1775  double muonECALError = 0.;
1776  unsigned nMuons = 0;
1777 
1778  ::math::XYZVector photonAtECAL(0., 0., 0.);
1779  std::vector<std::tuple<unsigned, ::math::XYZVector, double>>
1780  ecalClusters; // last element (double) : correction under the egamma hypothesis
1781  double sumEcalClusters = 0;
1782  ::math::XYZVector hadronDirection(
1783  hclusterref->position().X(), hclusterref->position().Y(), hclusterref->position().Z());
1784  hadronDirection = hadronDirection.Unit();
1785  ::math::XYZVector hadronAtECAL = totalHcal * hadronDirection;
1786 
1787  // Loop over all tracks associated to this HCAL cluster
1788  for (auto const& trk : sortedTracks) {
1789  unsigned iTrack = trk.second;
1790 
1791  // Check the track has not already been used (e.g., in electrons, conversions...)
1792  if (!active[iTrack])
1793  continue;
1794  // Sanity check 1
1795  PFBlockElement::Type type = elements[iTrack].type();
1796  assert(type == reco::PFBlockElement::TRACK);
1797  // Sanity check 2
1798  reco::TrackRef trackRef = elements[iTrack].trackRef();
1799  assert(!trackRef.isNull());
1800 
1801  // The direction at ECAL entrance
1802  const ::math::XYZPointF& chargedPosition =
1803  dynamic_cast<const reco::PFBlockElementTrack*>(&elements[iTrack])->positionAtECALEntrance();
1804  ::math::XYZVector chargedDirection(chargedPosition.X(), chargedPosition.Y(), chargedPosition.Z());
1805  chargedDirection = chargedDirection.Unit();
1806 
1807  // look for ECAL elements associated to iTrack (associated to iHcal)
1808  std::multimap<double, unsigned> sortedEcals;
1809  block.associatedElements(iTrack, linkData, sortedEcals, reco::PFBlockElement::ECAL, reco::PFBlock::LINKTEST_ALL);
1810 
1811  LogTrace("PFAlgo|createCandidatesHCAL") << "number of Ecal elements linked to this track: " << sortedEcals.size();
1812 
1813  // look for HO elements associated to iTrack (associated to iHcal)
1814  std::multimap<double, unsigned> sortedHOs;
1815  if (useHO_) {
1816  block.associatedElements(iTrack, linkData, sortedHOs, reco::PFBlockElement::HO, reco::PFBlock::LINKTEST_ALL);
1817  }
1818  LogTrace("PFAlgo|createCandidatesHCAL")
1819  << "PFAlgo : number of HO elements linked to this track: " << sortedHOs.size();
1820 
1821  // Create a PF Candidate right away if the track is a tight muon
1822  reco::MuonRef muonRef = elements[iTrack].muonRef();
1823 
1824  bool thisIsAMuon = PFMuonAlgo::isMuon(elements[iTrack]);
1825  bool thisIsAnIsolatedMuon = PFMuonAlgo::isIsolatedMuon(elements[iTrack]);
1826  bool thisIsALooseMuon = false;
1827 
1828  if (!thisIsAMuon) {
1829  thisIsALooseMuon = PFMuonAlgo::isLooseMuon(elements[iTrack]);
1830  }
1831 
1832  if (thisIsAMuon) {
1833  LogTrace("PFAlgo|createCandidatesHCAL") << "This track is identified as a muon - remove it from the stack";
1834  LogTrace("PFAlgo|createCandidatesHCAL") << elements[iTrack];
1835 
1836  // Create a muon.
1837 
1838  unsigned tmpi = reconstructTrack(elements[iTrack]);
1839 
1840  (*pfCandidates_)[tmpi].addElementInBlock(blockref, iTrack);
1841  (*pfCandidates_)[tmpi].addElementInBlock(blockref, iHcal);
1842  double muonHcal = std::min(muonHCAL_[0] + muonHCAL_[1], totalHcal);
1843 
1844  // if muon is isolated and muon momentum exceeds the calo energy, absorb the calo energy
1845  // rationale : there has been a EM showering by the muon in the calorimeter - or the coil -
1846  // and we don't want to double count the energy
1847  bool letMuonEatCaloEnergy = false;
1848 
1849  if (thisIsAnIsolatedMuon) {
1850  // The factor 1.3 is the e/pi factor in HCAL...
1851  double totalCaloEnergy = totalHcal / 1.30;
1852  unsigned iEcal = 0;
1853  if (!sortedEcals.empty()) {
1854  iEcal = sortedEcals.begin()->second;
1855  PFClusterRef eclusterref = elements[iEcal].clusterRef();
1856  totalCaloEnergy += eclusterref->energy();
1857  }
1858 
1859  if (useHO_) {
1860  // The factor 1.3 is assumed to be the e/pi factor in HO, too.
1861  unsigned iHO = 0;
1862  if (!sortedHOs.empty()) {
1863  iHO = sortedHOs.begin()->second;
1864  PFClusterRef eclusterref = elements[iHO].clusterRef();
1865  totalCaloEnergy += eclusterref->energy() / 1.30;
1866  }
1867  }
1868 
1869  if ((pfCandidates_->back()).p() > totalCaloEnergy)
1870  letMuonEatCaloEnergy = true;
1871  }
1872 
1873  if (letMuonEatCaloEnergy)
1874  muonHcal = totalHcal;
1875  double muonEcal = 0.;
1876  unsigned iEcal = 0;
1877  if (!sortedEcals.empty()) {
1878  iEcal = sortedEcals.begin()->second;
1879  PFClusterRef eclusterref = elements[iEcal].clusterRef();
1880  (*pfCandidates_)[tmpi].addElementInBlock(blockref, iEcal);
1881  muonEcal = std::min(muonECAL_[0] + muonECAL_[1], eclusterref->energy());
1882  if (letMuonEatCaloEnergy)
1883  muonEcal = eclusterref->energy();
1884  // If the muon expected energy accounts for the whole ecal cluster energy, lock the ecal cluster
1885  if (eclusterref->energy() - muonEcal < 0.2)
1886  active[iEcal] = false;
1887  (*pfCandidates_)[tmpi].setEcalEnergy(eclusterref->energy(), muonEcal);
1888  }
1889  unsigned iHO = 0;
1890  double muonHO = 0.;
1891  if (useHO_) {
1892  if (!sortedHOs.empty()) {
1893  iHO = sortedHOs.begin()->second;
1894  PFClusterRef hoclusterref = elements[iHO].clusterRef();
1895  (*pfCandidates_)[tmpi].addElementInBlock(blockref, iHO);
1896  muonHO = std::min(muonHO_[0] + muonHO_[1], hoclusterref->energy());
1897  if (letMuonEatCaloEnergy)
1898  muonHO = hoclusterref->energy();
1899  // If the muon expected energy accounts for the whole HO cluster energy, lock the HO cluster
1900  if (hoclusterref->energy() - muonHO < 0.2)
1901  active[iHO] = false;
1902  (*pfCandidates_)[tmpi].setHcalEnergy(totalHcal, muonHcal);
1903  (*pfCandidates_)[tmpi].setHoEnergy(hoclusterref->energy(), muonHO);
1904  }
1905  } else {
1906  (*pfCandidates_)[tmpi].setHcalEnergy(totalHcal, muonHcal);
1907  }
1908  setHcalDepthInfo((*pfCandidates_)[tmpi], *hclusterref);
1909 
1910  if (letMuonEatCaloEnergy) {
1911  muonHCALEnergy += totalHcal;
1912  if (useHO_)
1913  muonHCALEnergy += muonHO;
1914  muonHCALError += 0.;
1915  muonECALEnergy += muonEcal;
1916  muonECALError += 0.;
1917  photonAtECAL -= muonEcal * chargedDirection;
1918  hadronAtECAL -= totalHcal * chargedDirection;
1919  if (!sortedEcals.empty())
1920  active[iEcal] = false;
1921  active[iHcal] = false;
1922  if (useHO_ && !sortedHOs.empty())
1923  active[iHO] = false;
1924  } else {
1925  // Estimate of the energy deposit & resolution in the calorimeters
1926  muonHCALEnergy += muonHCAL_[0];
1927  muonHCALError += muonHCAL_[1] * muonHCAL_[1];
1928  if (muonHO > 0.) {
1929  muonHCALEnergy += muonHO_[0];
1930  muonHCALError += muonHO_[1] * muonHO_[1];
1931  }
1932  muonECALEnergy += muonECAL_[0];
1933  muonECALError += muonECAL_[1] * muonECAL_[1];
1934  // ... as well as the equivalent "momentum" at ECAL entrance
1935  photonAtECAL -= muonECAL_[0] * chargedDirection;
1936  hadronAtECAL -= muonHCAL_[0] * chargedDirection;
1937  }
1938 
1939  // Remove it from the stack
1940  active[iTrack] = false;
1941  // Go to next track
1942  continue;
1943  }
1944 
1945  //
1946 
1947  LogTrace("PFAlgo|createCandidatesHCAL") << "elements[iTrack=" << iTrack << "]=" << elements[iTrack];
1948 
1949  // introduce tracking errors
1950  double trackMomentum = trackRef->p();
1951  totalChargedMomentum += trackMomentum;
1952 
1953  // If the track is not a tight muon, but still resembles a muon
1954  // keep it for later for blocks with too large a charged energy
1955  if (thisIsALooseMuon && !thisIsAMuon)
1956  nMuons += 1;
1957 
1958  // ... and keep anyway the pt error for possible fake rejection
1959  // ... blow up errors of 5th and 4th iteration, to reject those
1960  // ... tracks first (in case it's needed)
1961  double dpt = trackRef->ptError();
1962  double blowError = PFTrackAlgoTools::errorScale(trackRef->algo(), factors45_);
1963  // except if it is from an interaction
1964  bool isPrimaryOrSecondary = isFromSecInt(elements[iTrack], "all");
1965 
1966  if (isPrimaryOrSecondary)
1967  blowError = 1.;
1968 
1969  std::pair<unsigned, bool> tkmuon(iTrack, thisIsALooseMuon);
1970  associatedTracks.emplace(-dpt * blowError, tkmuon);
1971 
1972  // Also keep the total track momentum error for comparison with the calo energy
1973  double dp = trackRef->qoverpError() * trackMomentum * trackMomentum;
1974  sumpError2 += dp * dp;
1975 
1976  bool connectedToEcal = false; // Will become true if there is at least one ECAL cluster connected
1977  if (!sortedEcals.empty()) { // start case: at least one ecal element associated to iTrack
1978 
1979  // Loop over all connected ECAL clusters
1980  for (auto const& ecal : sortedEcals) {
1981  unsigned iEcal = ecal.second;
1982  double dist = ecal.first;
1983 
1984  // Ignore ECAL cluters already used
1985  if (!active[iEcal]) {
1986  LogTrace("PFAlgo|createCandidatesHCAL") << "cluster locked";
1987  continue;
1988  }
1989 
1990  // Sanity checks
1991  PFBlockElement::Type type = elements[iEcal].type();
1992  assert(type == PFBlockElement::ECAL);
1993  PFClusterRef eclusterref = elements[iEcal].clusterRef();
1994  assert(!eclusterref.isNull());
1995 
1996  // Check if this ECAL is not closer to another track - ignore it in that case
1997  std::multimap<double, unsigned> sortedTracksEcal;
1998  block.associatedElements(
1999  iEcal, linkData, sortedTracksEcal, reco::PFBlockElement::TRACK, reco::PFBlock::LINKTEST_ALL);
2000  unsigned jTrack = sortedTracksEcal.begin()->second;
2001  if (jTrack != iTrack)
2002  continue;
2003 
2004  double distEcal = block.dist(jTrack, iEcal, linkData, reco::PFBlock::LINKTEST_ALL);
2005 
2006  float ecalEnergyCalibrated = eclusterref->correctedEnergy(); // calibrated based on the egamma hypothesis
2007  float ecalEnergy = eclusterref->energy();
2008  ::math::XYZVector photonDirection(
2009  eclusterref->position().X(), eclusterref->position().Y(), eclusterref->position().Z());
2010  photonDirection = photonDirection.Unit();
2011 
2012  if (!connectedToEcal) { // This is the closest ECAL cluster - will add its energy later
2013 
2014  LogTrace("PFAlgo|createCandidatesHCAL") << "closest: " << elements[iEcal];
2015 
2016  connectedToEcal = true;
2017  // PJ 1st-April-09 : To be done somewhere !!! (Had to comment it, but it is needed)
2018  // currentChargedHadron.addElementInBlock( blockref, iEcal );
2019 
2020  // KH: we don't know if this satellite is due to egamma or hadron shower. use raw energy for PF hadron calibration._ store also calibration constant.
2021  double ecalCalibFactor = (ecalEnergy > 1E-9) ? ecalEnergyCalibrated / ecalEnergy : 1.;
2022  std::tuple<unsigned, ::math::XYZVector, double> satellite(
2023  iEcal, ecalEnergy * photonDirection, ecalCalibFactor);
2024  ecalSatellites.emplace(-1., satellite);
2025 
2026  } else { // Keep satellite clusters for later
2027 
2028  // KH: same as above.
2029  double ecalCalibFactor = (ecalEnergy > 1E-9) ? ecalEnergyCalibrated / ecalEnergy : 1.;
2030  std::tuple<unsigned, ::math::XYZVector, double> satellite(
2031  iEcal, ecalEnergy * photonDirection, ecalCalibFactor);
2032  ecalSatellites.emplace(dist, satellite);
2033  }
2034 
2035  std::pair<double, unsigned> associatedEcal(distEcal, iEcal);
2036  associatedEcals.emplace(iTrack, associatedEcal);
2037 
2038  } // End loop ecal associated to iTrack
2039  } // end case: at least one ecal element associated to iTrack
2040 
2041  if (useHO_ && !sortedHOs.empty()) { // start case: at least one ho element associated to iTrack
2042 
2043  // Loop over all connected HO clusters
2044  for (auto const& ho : sortedHOs) {
2045  unsigned iHO = ho.second;
2046  double distHO = ho.first;
2047 
2048  // Ignore HO cluters already used
2049  if (!active[iHO]) {
2050  LogTrace("PFAlgo|createCandidatesHCAL") << "cluster locked";
2051  continue;
2052  }
2053 
2054  // Sanity checks
2055  PFBlockElement::Type type = elements[iHO].type();
2056  assert(type == PFBlockElement::HO);
2057  PFClusterRef hoclusterref = elements[iHO].clusterRef();
2058  assert(!hoclusterref.isNull());
2059 
2060  // Check if this HO is not closer to another track - ignore it in that case
2061  std::multimap<double, unsigned> sortedTracksHO;
2062  block.associatedElements(
2063  iHO, linkData, sortedTracksHO, reco::PFBlockElement::TRACK, reco::PFBlock::LINKTEST_ALL);
2064  unsigned jTrack = sortedTracksHO.begin()->second;
2065  if (jTrack != iTrack)
2066  continue;
2067 
2068  // double chi2HO = block.chi2(jTrack,iHO,linkData,
2069  // reco::PFBlock::LINKTEST_ALL);
2070  //double distHO = block.dist(jTrack,iHO,linkData,
2071  // reco::PFBlock::LINKTEST_ALL);
2072 
2073  // Increment the total energy by the energy of this HO cluster
2074  totalHO += hoclusterref->energy();
2075  active[iHO] = false;
2076  // Keep track for later reference in the PFCandidate.
2077  std::pair<double, unsigned> associatedHO(distHO, iHO);
2078  associatedHOs.emplace(iTrack, associatedHO);
2079 
2080  } // End loop ho associated to iTrack
2081  } // end case: at least one ho element associated to iTrack
2082 
2083  } // end loop on tracks associated to hcal element iHcal
2084 
2085  // Include totalHO in totalHCAL for the time being (it will be calibrated as HCAL energy)
2086  totalHcal += totalHO;
2087 
2088  // test compatibility between calo and tracker. //////////////
2089 
2090  double caloEnergy = 0.;
2091  double slopeEcal = 1.0;
2092  double calibEcal = 0.;
2093  double calibHcal = 0.;
2094  hadronDirection = hadronAtECAL.Unit();
2095 
2096  // Determine the expected calo resolution from the total charged momentum
2097  double caloResolution = neutralHadronEnergyResolution(totalChargedMomentum, hclusterref->positionREP().Eta());
2098  caloResolution *= totalChargedMomentum;
2099  // Account for muons
2100  caloResolution = std::sqrt(caloResolution * caloResolution + muonHCALError + muonECALError);
2101  totalEcal -= std::min(totalEcal, muonECALEnergy);
2102  totalEcalEGMCalib -= std::min(totalEcalEGMCalib, muonECALEnergy);
2103  totalHcal -= std::min(totalHcal, muonHCALEnergy);
2104  if (totalEcal < 1E-9)
2105  photonAtECAL = ::math::XYZVector(0., 0., 0.);
2106  if (totalHcal < 1E-9)
2107  hadronAtECAL = ::math::XYZVector(0., 0., 0.);
2108 
2109  // Loop over all ECAL satellites, starting for the closest to the various tracks
2110  // and adding other satellites until saturation of the total track momentum
2111  // Note : for code simplicity, the first element of the loop is the HCAL cluster
2112  // with 0 energy in the ECAL
2113  for (auto const& ecalSatellite : ecalSatellites) {
2114  // Add the energy of this ECAL cluster
2115  double previousCalibEcal = calibEcal;
2116  double previousCalibHcal = calibHcal;
2117  double previousCaloEnergy = caloEnergy;
2118  double previousSlopeEcal = slopeEcal;
2119  ::math::XYZVector previousHadronAtECAL = hadronAtECAL;
2120  //
2121  totalEcal +=
2122  sqrt(std::get<1>(ecalSatellite.second).Mag2()); // KH: raw ECAL energy for input to PF hadron calibration
2123  totalEcalEGMCalib += sqrt(std::get<1>(ecalSatellite.second).Mag2()) *
2124  std::get<2>(ecalSatellite.second); // KH: calibrated ECAL energy under the egamma hypothesis
2125  photonAtECAL += std::get<1>(ecalSatellite.second) *
2126  std::get<2>(ecalSatellite.second); // KH: calibrated ECAL energy under the egamma hypothesis
2127  calibEcal = std::max(0., totalEcal); // KH: preparing for hadron calibration
2128  calibHcal = std::max(0., totalHcal);
2129  hadronAtECAL = calibHcal * hadronDirection;
2130  // Calibrate ECAL and HCAL energy under the hadron hypothesis.
2131  calibration_.energyEmHad(totalChargedMomentum,
2132  calibEcal,
2133  calibHcal,
2134  hclusterref->positionREP().Eta(),
2135  hclusterref->positionREP().Phi());
2136  caloEnergy = calibEcal + calibHcal;
2137  if (totalEcal > 0.)
2138  slopeEcal = calibEcal / totalEcal;
2139 
2140  hadronAtECAL = calibHcal * hadronDirection;
2141 
2142  // Continue looping until all closest clusters are exhausted and as long as
2143  // the calorimetric energy does not saturate the total momentum.
2144  if (ecalSatellite.first < 0. || caloEnergy - totalChargedMomentum <= 0.) {
2145  LogTrace("PFAlgo|createCandidatesHCAL")
2146  << "\t\t\tactive, adding " << std::get<1>(ecalSatellite.second) << " to ECAL energy, and locking";
2147  active[std::get<0>(ecalSatellite.second)] = false;
2148  double clusterEnergy =
2149  sqrt(std::get<1>(ecalSatellite.second).Mag2()) *
2150  std::get<2>(ecalSatellite.second); // KH: ECAL energy calibrated under the egamma hypothesis
2151  if (clusterEnergy > 50) { // KH: used to split energetic ecal clusters (E>50 GeV)
2152  ecalClusters.push_back(ecalSatellite.second);
2153  sumEcalClusters += clusterEnergy;
2154  }
2155  continue;
2156  }
2157 
2158  // Otherwise, do not consider the last cluster examined and exit.
2159  // active[is->second.first] = true;
2160  totalEcal -= sqrt(std::get<1>(ecalSatellite.second).Mag2());
2161  totalEcalEGMCalib -= sqrt(std::get<1>(ecalSatellite.second).Mag2()) * std::get<2>(ecalSatellite.second);
2162  photonAtECAL -= std::get<1>(ecalSatellite.second) * std::get<2>(ecalSatellite.second);
2163  calibEcal = previousCalibEcal;
2164  calibHcal = previousCalibHcal;
2165  hadronAtECAL = previousHadronAtECAL;
2166  slopeEcal = previousSlopeEcal;
2167  caloEnergy = previousCaloEnergy;
2168 
2169  break;
2170  }
2171 
2172  // Sanity check !
2173  assert(caloEnergy >= 0);
2174 
2175  // And now check for hadronic energy excess...
2176 
2177  //colin: resolution should be measured on the ecal+hcal case.
2178  // however, the result will be close.
2179  // double caloResolution = neutralHadronEnergyResolution( caloEnergy );
2180  // caloResolution *= caloEnergy;
2181  // PJ The resolution is on the expected charged calo energy !
2182  //double caloResolution = neutralHadronEnergyResolution( totalChargedMomentum, hclusterref->positionREP().Eta());
2183  //caloResolution *= totalChargedMomentum;
2184  // that of the charged particles linked to the cluster!
2185 
2187  if (totalChargedMomentum - caloEnergy > nSigmaTRACK_ * caloResolution) {
2188  // First consider loose muons
2189  if (nMuons > 0) {
2190  for (auto const& trk : associatedTracks) {
2191  // Only muons
2192  if (!trk.second.second)
2193  continue;
2194 
2195  const unsigned int iTrack = trk.second.first;
2196  // Only active tracks
2197  if (!active[iTrack])
2198  continue;
2199 
2200  const double trackMomentum = elements[trk.second.first].trackRef()->p();
2201 
2202  // look for ECAL elements associated to iTrack (associated to iHcal)
2203  std::multimap<double, unsigned> sortedEcals;
2204  block.associatedElements(
2205  iTrack, linkData, sortedEcals, reco::PFBlockElement::ECAL, reco::PFBlock::LINKTEST_ALL);
2206  std::multimap<double, unsigned> sortedHOs;
2207  block.associatedElements(iTrack, linkData, sortedHOs, reco::PFBlockElement::HO, reco::PFBlock::LINKTEST_ALL);
2208 
2209  //Here allow for loose muons!
2210  auto& muon = (*pfCandidates_)[reconstructTrack(elements[iTrack], true)];
2211 
2212  muon.addElementInBlock(blockref, iTrack);
2213  muon.addElementInBlock(blockref, iHcal);
2214  const double muonHcal = std::min(muonHCAL_[0] + muonHCAL_[1], totalHcal - totalHO);
2215  double muonHO = 0.;
2216  muon.setHcalEnergy(totalHcal, muonHcal);
2217  if (!sortedEcals.empty()) {
2218  const unsigned int iEcal = sortedEcals.begin()->second;
2219  PFClusterRef eclusterref = elements[iEcal].clusterRef();
2220  muon.addElementInBlock(blockref, iEcal);
2221  const double muonEcal = std::min(muonECAL_[0] + muonECAL_[1], eclusterref->energy());
2222  muon.setEcalEnergy(eclusterref->energy(), muonEcal);
2223  }
2224  if (useHO_ && !sortedHOs.empty()) {
2225  const unsigned int iHO = sortedHOs.begin()->second;
2226  PFClusterRef hoclusterref = elements[iHO].clusterRef();
2227  muon.addElementInBlock(blockref, iHO);
2228  muonHO = std::min(muonHO_[0] + muonHO_[1], hoclusterref->energy());
2229  muon.setHcalEnergy(max(totalHcal - totalHO, 0.0), muonHcal);
2230  muon.setHoEnergy(hoclusterref->energy(), muonHO);
2231  }
2232  setHcalDepthInfo(muon, *hclusterref);
2233  // Remove it from the block
2234  const ::math::XYZPointF& chargedPosition =
2235  dynamic_cast<const reco::PFBlockElementTrack*>(&elements[trk.second.first])->positionAtECALEntrance();
2236  ::math::XYZVector chargedDirection(chargedPosition.X(), chargedPosition.Y(), chargedPosition.Z());
2237  chargedDirection = chargedDirection.Unit();
2238  totalChargedMomentum -= trackMomentum;
2239  // Update the calo energies
2240  if (totalEcal > 0.)
2241  calibEcal -= std::min(calibEcal, muonECAL_[0] * calibEcal / totalEcal);
2242  if (totalHcal > 0.)
2243  calibHcal -= std::min(calibHcal, muonHCAL_[0] * calibHcal / totalHcal);
2244  totalEcal -= std::min(totalEcal, muonECAL_[0]);
2245  totalHcal -= std::min(totalHcal, muonHCAL_[0]);
2246  if (totalEcal > muonECAL_[0])
2247  photonAtECAL -= muonECAL_[0] * chargedDirection;
2248  if (totalHcal > muonHCAL_[0])
2249  hadronAtECAL -= muonHCAL_[0] * calibHcal / totalHcal * chargedDirection;
2250  caloEnergy = calibEcal + calibHcal;
2251  muonHCALEnergy += muonHCAL_[0];
2252  muonHCALError += muonHCAL_[1] * muonHCAL_[1];
2253  if (muonHO > 0.) {
2254  muonHCALEnergy += muonHO_[0];
2255  muonHCALError += muonHO_[1] * muonHO_[1];
2256  if (totalHcal > 0.) {
2257  calibHcal -= std::min(calibHcal, muonHO_[0] * calibHcal / totalHcal);
2258  totalHcal -= std::min(totalHcal, muonHO_[0]);
2259  }
2260  }
2261  muonECALEnergy += muonECAL_[0];
2262  muonECALError += muonECAL_[1] * muonECAL_[1];
2263  active[iTrack] = false;
2264  // Stop the loop whenever enough muons are removed
2265  //Commented out: Keep looking for muons since they often come in pairs -Matt
2266  //if ( totalChargedMomentum < caloEnergy ) break;
2267  }
2268  // New calo resolution.
2269  caloResolution = neutralHadronEnergyResolution(totalChargedMomentum, hclusterref->positionREP().Eta());
2270  caloResolution *= totalChargedMomentum;
2271  caloResolution = std::sqrt(caloResolution * caloResolution + muonHCALError + muonECALError);
2272  }
2273  }
2274 
2275 #ifdef EDM_ML_DEBUG
2276  LogTrace("PFAlgo|createCandidatesHCAL") << "\tBefore Cleaning ";
2277  LogTrace("PFAlgo|createCandidatesHCAL") << "\tCompare Calo Energy to total charged momentum ";
2278  LogTrace("PFAlgo|createCandidatesHCAL") << "\t\tsum p = " << totalChargedMomentum << " +- " << sqrt(sumpError2);
2279  LogTrace("PFAlgo|createCandidatesHCAL") << "\t\tsum ecal = " << totalEcal;
2280  LogTrace("PFAlgo|createCandidatesHCAL") << "\t\tsum hcal = " << totalHcal;
2281  LogTrace("PFAlgo|createCandidatesHCAL") << "\t\t => Calo Energy = " << caloEnergy << " +- " << caloResolution;
2282  LogTrace("PFAlgo|createCandidatesHCAL")
2283  << "\t\t => Calo Energy- total charged momentum = " << caloEnergy - totalChargedMomentum << " +- "
2284  << sqrt(sumpError2 + caloResolution * caloResolution);
2285 #endif
2286 
2287  // Second consider bad tracks (if still needed after muon removal)
2288  unsigned corrTrack = 10000000;
2289  double corrFact = 1.;
2290 
2291  if (rejectTracks_Bad_ && totalChargedMomentum - caloEnergy > nSigmaTRACK_ * caloResolution) {
2292  for (auto const& trk : associatedTracks) {
2293  const unsigned iTrack = trk.second.first;
2294  // Only active tracks
2295  if (!active[iTrack])
2296  continue;
2297  const reco::TrackRef& trackref = elements[trk.second.first].trackRef();
2298 
2299  const double dptRel = fabs(trk.first) / trackref->pt() * 100;
2300  const bool isSecondary = isFromSecInt(elements[iTrack], "secondary");
2301  const bool isPrimary = isFromSecInt(elements[iTrack], "primary");
2302 
2303  if (isSecondary && dptRel < dptRel_DispVtx_)
2304  continue;
2305  // Consider only bad tracks
2306  if (fabs(trk.first) < ptError_)
2307  break;
2308  // What would become the block charged momentum if this track were removed
2309  const double wouldBeTotalChargedMomentum = totalChargedMomentum - trackref->p();
2310  // Reject worst tracks, as long as the total charged momentum
2311  // is larger than the calo energy
2312 
2313  if (wouldBeTotalChargedMomentum > caloEnergy) {
2314  if (isSecondary) {
2315  LogTrace("PFAlgo|createCandidatesHCAL")
2316  << "In bad track rejection step dptRel = " << dptRel << " dptRel_DispVtx_ = " << dptRel_DispVtx_;
2317  LogTrace("PFAlgo|createCandidatesHCAL")
2318  << "The calo energy would be still smaller even without this track but it is attached to a NI";
2319  }
2320 
2321  if (isPrimary || (isSecondary && dptRel < dptRel_DispVtx_))
2322  continue;
2323  active[iTrack] = false;
2324  totalChargedMomentum = wouldBeTotalChargedMomentum;
2325  LogTrace("PFAlgo|createCandidatesHCAL")
2326  << "\tElement " << elements[iTrack] << " rejected (dpt = " << -trk.first
2327  << " GeV/c, algo = " << trackref->algo() << ")";
2328  // Just rescale the nth worst track momentum to equalize the calo energy
2329  } else {
2330  if (isPrimary)
2331  break;
2332  corrTrack = iTrack;
2333  corrFact = (caloEnergy - wouldBeTotalChargedMomentum) / elements[trk.second.first].trackRef()->p();
2334  if (trackref->p() * corrFact < 0.05) {
2335  corrFact = 0.;
2336  active[iTrack] = false;
2337  }
2338  totalChargedMomentum -= trackref->p() * (1. - corrFact);
2339  LogTrace("PFAlgo|createCandidatesHCAL")
2340  << "\tElement " << elements[iTrack] << " (dpt = " << -trk.first << " GeV/c, algo = " << trackref->algo()
2341  << ") rescaled by " << corrFact << " Now the total charged momentum is " << totalChargedMomentum;
2342  break;
2343  }
2344  }
2345  }
2346 
2347  // New determination of the calo and track resolution avec track deletion/rescaling.
2348  caloResolution = neutralHadronEnergyResolution(totalChargedMomentum, hclusterref->positionREP().Eta());
2349  caloResolution *= totalChargedMomentum;
2350  caloResolution = std::sqrt(caloResolution * caloResolution + muonHCALError + muonECALError);
2351 
2352  // Check if the charged momentum is still very inconsistent with the calo measurement.
2353  // In this case, just drop all tracks from 4th and 5th iteration linked to this block
2354 
2355  if (rejectTracks_Step45_ && sortedTracks.size() > 1 &&
2356  totalChargedMomentum - caloEnergy > nSigmaTRACK_ * caloResolution) {
2357  for (auto const& trk : associatedTracks) {
2358  unsigned iTrack = trk.second.first;
2359  reco::TrackRef trackref = elements[iTrack].trackRef();
2360  if (!active[iTrack])
2361  continue;
2362 
2363  double dptRel = fabs(trk.first) / trackref->pt() * 100;
2364  bool isPrimaryOrSecondary = isFromSecInt(elements[iTrack], "all");
2365 
2366  if (isPrimaryOrSecondary && dptRel < dptRel_DispVtx_)
2367  continue;
2368 
2369  if (PFTrackAlgoTools::step5(trackref->algo())) {
2370  active[iTrack] = false;
2371  totalChargedMomentum -= trackref->p();
2372 
2373  LogTrace("PFAlgo|createCandidatesHCAL")
2374  << "\tElement " << elements[iTrack] << " rejected (dpt = " << -trk.first
2375  << " GeV/c, algo = " << trackref->algo() << ")";
2376  }
2377  }
2378  }
2379 
2380  // New determination of the calo and track resolution avec track deletion/rescaling.
2381  caloResolution = neutralHadronEnergyResolution(totalChargedMomentum, hclusterref->positionREP().Eta());
2382  caloResolution *= totalChargedMomentum;
2383  caloResolution = std::sqrt(caloResolution * caloResolution + muonHCALError + muonECALError);
2384 
2385  // Make PF candidates with the remaining tracks in the block
2386  sumpError2 = 0.;
2387  for (auto const& trk : associatedTracks) {
2388  unsigned iTrack = trk.second.first;
2389  if (!active[iTrack])
2390  continue;
2391  reco::TrackRef trackRef = elements[iTrack].trackRef();
2392  double trackMomentum = trackRef->p();
2393  double dp = trackRef->qoverpError() * trackMomentum * trackMomentum;
2394  unsigned tmpi = reconstructTrack(elements[iTrack]);
2395 
2396  (*pfCandidates_)[tmpi].addElementInBlock(blockref, iTrack);
2397  (*pfCandidates_)[tmpi].addElementInBlock(blockref, iHcal);
2398  setHcalDepthInfo((*pfCandidates_)[tmpi], *hclusterref);
2399  auto myEcals = associatedEcals.equal_range(iTrack);
2400  for (auto ii = myEcals.first; ii != myEcals.second; ++ii) {
2401  unsigned iEcal = ii->second.second;
2402  if (active[iEcal])
2403  continue;
2404  (*pfCandidates_)[tmpi].addElementInBlock(blockref, iEcal);
2405  }
2406 
2407  if (useHO_) {
2408  auto myHOs = associatedHOs.equal_range(iTrack);
2409  for (auto ii = myHOs.first; ii != myHOs.second; ++ii) {
2410  unsigned iHO = ii->second.second;
2411  if (active[iHO])
2412  continue;
2413  (*pfCandidates_)[tmpi].addElementInBlock(blockref, iHO);
2414  }
2415  }
2416 
2417  if (iTrack == corrTrack) {
2418  if (corrFact < 0.)
2419  corrFact = 0.; // protect against negative scaling
2420  (*pfCandidates_)[tmpi].rescaleMomentum(corrFact);
2421  trackMomentum *= corrFact;
2422  }
2423  chargedHadronsIndices.push_back(tmpi);
2424  chargedHadronsInBlock.push_back(iTrack);
2425  active[iTrack] = false;
2426  hcalP.push_back(trackMomentum);
2427  hcalDP.push_back(dp);
2428  if (dp / trackMomentum > maxDPovP)
2429  maxDPovP = dp / trackMomentum;
2430  sumpError2 += dp * dp;
2431  }
2432 
2433  // The total uncertainty of the difference Calo-Track
2434  double totalError = sqrt(sumpError2 + caloResolution * caloResolution);
2435 
2436 #ifdef EDM_ML_DEBUG
2437  LogTrace("PFAlgo|createCandidatesHCAL")
2438  << "\tCompare Calo Energy to total charged momentum " << endl
2439  << "\t\tsum p = " << totalChargedMomentum << " +- " << sqrt(sumpError2) << endl
2440  << "\t\tsum ecal = " << totalEcal << endl
2441  << "\t\tsum hcal = " << totalHcal << endl
2442  << "\t\t => Calo Energy = " << caloEnergy << " +- " << caloResolution << endl
2443  << "\t\t => Calo Energy- total charged momentum = " << caloEnergy - totalChargedMomentum << " +- "
2444  << totalError;
2445 #endif
2446 
2447  /* */
2448 
2450  double nsigma = nSigmaHCAL(totalChargedMomentum, hclusterref->positionREP().Eta());
2451  //double nsigma = nSigmaHCAL(caloEnergy,hclusterref->positionREP().Eta());
2452  if (abs(totalChargedMomentum - caloEnergy) < nsigma * totalError) {
2453  // deposited caloEnergy compatible with total charged momentum
2454  // if tracking errors are large take weighted average
2455 
2456 #ifdef EDM_ML_DEBUG
2457  LogTrace("PFAlgo|createCandidatesHCAL")
2458  << "\t\tcase 1: COMPATIBLE "
2459  << "|Calo Energy- total charged momentum| = " << abs(caloEnergy - totalChargedMomentum) << " < " << nsigma
2460  << " x " << totalError;
2461  if (maxDPovP < 0.1)
2462  LogTrace("PFAlgo|createCandidatesHCAL") << "\t\t\tmax DP/P = " << maxDPovP << " less than 0.1: do nothing ";
2463  else
2464  LogTrace("PFAlgo|createCandidatesHCAL")
2465  << "\t\t\tmax DP/P = " << maxDPovP << " > 0.1: take weighted averages ";
2466 #endif
2467 
2468  // if max DP/P < 10% do nothing
2469  if (maxDPovP > 0.1) {
2470  // for each track associated to hcal
2471  // int nrows = tkIs.size();
2472  int nrows = chargedHadronsIndices.size();
2473  TMatrixTSym<double> a(nrows);
2474  TVectorD b(nrows);
2475  TVectorD check(nrows);
2476  double sigma2E = caloResolution * caloResolution;
2477  for (int i = 0; i < nrows; i++) {
2478  double sigma2i = hcalDP[i] * hcalDP[i];
2479  LogTrace("PFAlgo|createCandidatesHCAL")
2480  << "\t\t\ttrack associated to hcal " << i << " P = " << hcalP[i] << " +- " << hcalDP[i];
2481  a(i, i) = 1. / sigma2i + 1. / sigma2E;
2482  b(i) = hcalP[i] / sigma2i + caloEnergy / sigma2E;
2483  for (int j = 0; j < nrows; j++) {
2484  if (i == j)
2485  continue;
2486  a(i, j) = 1. / sigma2E;
2487  } // end loop on j
2488  } // end loop on i
2489 
2490  // solve ax = b
2491  TDecompChol decomp(a);
2492  bool ok = false;
2493  TVectorD x = decomp.Solve(b, ok);
2494  // for each track create a PFCandidate track
2495  // with a momentum rescaled to weighted average
2496  if (ok) {
2497  for (int i = 0; i < nrows; i++) {
2498  // unsigned iTrack = trackInfos[i].index;
2499  unsigned ich = chargedHadronsIndices[i];
2500  double rescaleFactor = x(i) / hcalP[i];
2501  if (rescaleFactor < 0.)
2502  rescaleFactor = 0.; // protect against negative scaling
2503  (*pfCandidates_)[ich].rescaleMomentum(rescaleFactor);
2504 
2505  LogTrace("PFAlgo|createCandidatesHCAL")
2506  << "\t\t\told p " << hcalP[i] << " new p " << x(i) << " rescale " << rescaleFactor;
2507  }
2508  } else {
2509  edm::LogError("PFAlgo::createCandidatesHCAL") << "TDecompChol.Solve returned ok=false";
2510  assert(0);
2511  }
2512  }
2513  }
2514 
2516  else if (caloEnergy > totalChargedMomentum) {
2517  //case 2: caloEnergy > totalChargedMomentum + nsigma*totalError
2518  //there is an excess of energy in the calos
2519  //create a neutral hadron or a photon
2520 
2521  double eNeutralHadron = caloEnergy - totalChargedMomentum;
2522  double ePhoton = (caloEnergy - totalChargedMomentum) /
2523  slopeEcal; // KH: this slopeEcal is computed based on ECAL energy under the hadron hypothesis,
2524  // thought we are creating photons.
2525  // This is a fuzzy case, but it should be better than corrected twice under both egamma and hadron hypotheses.
2526 
2527 #ifdef EDM_ML_DEBUG
2528  if (!sortedTracks.empty()) {
2529  LogTrace("PFAlgo|createCandidatesHCAL") << "\t\tcase 2: NEUTRAL DETECTION " << caloEnergy << " > " << nsigma
2530  << "x" << totalError << " + " << totalChargedMomentum;
2531  LogTrace("PFAlgo|createCandidatesHCAL") << "\t\tneutral activity detected: " << endl
2532  << "\t\t\t photon = " << ePhoton << endl
2533  << "\t\t\tor neutral hadron = " << eNeutralHadron;
2534 
2535  LogTrace("PFAlgo|createCandidatesHCAL") << "\t\tphoton or hadron ?";
2536  }
2537 
2538  if (sortedTracks.empty())
2539  LogTrace("PFAlgo|createCandidatesHCAL") << "\t\tno track -> hadron ";
2540  else
2541  LogTrace("PFAlgo|createCandidatesHCAL")
2542  << "\t\t" << sortedTracks.size() << " tracks -> check if the excess is photonic or hadronic";
2543 #endif
2544 
2545  double ratioMax = 0.;
2546  reco::PFClusterRef maxEcalRef;
2547  unsigned maxiEcal = 9999;
2548 
2549  // for each track associated to hcal: iterator IE ie :
2550 
2551  LogTrace("PFAlgo|createCandidatesHCAL") << "loop over sortedTracks.size()=" << sortedTracks.size();
2552  for (auto const& trk : sortedTracks) {
2553  unsigned iTrack = trk.second;
2554 
2555  PFBlockElement::Type type = elements[iTrack].type();
2556  assert(type == reco::PFBlockElement::TRACK);
2557 
2558  reco::TrackRef trackRef = elements[iTrack].trackRef();
2559  assert(!trackRef.isNull());
2560 
2561  auto iae = associatedEcals.find(iTrack);
2562 
2563  if (iae == associatedEcals.end())
2564  continue;
2565 
2566  // double distECAL = iae->second.first;
2567  unsigned iEcal = iae->second.second;
2568 
2569  PFBlockElement::Type typeEcal = elements[iEcal].type();
2570  assert(typeEcal == reco::PFBlockElement::ECAL);
2571 
2572  reco::PFClusterRef clusterRef = elements[iEcal].clusterRef();
2573  assert(!clusterRef.isNull());
2574 
2575  double pTrack = trackRef->p();
2576  double eECAL = clusterRef->energy();
2577  double eECALOverpTrack = eECAL / pTrack;
2578 
2579  if (eECALOverpTrack > ratioMax) {
2580  ratioMax = eECALOverpTrack;
2581  maxEcalRef = clusterRef;
2582  maxiEcal = iEcal;
2583  }
2584 
2585  } // end loop on tracks associated to hcal: iterator IE ie
2586 
2587  std::vector<reco::PFClusterRef> pivotalClusterRef;
2588  std::vector<unsigned> iPivotal;
2589  std::vector<double> particleEnergy, ecalEnergy, hcalEnergy, rawecalEnergy, rawhcalEnergy;
2590  std::vector<::math::XYZVector> particleDirection;
2591 
2592  // If the excess is smaller than the ecal energy, assign the whole
2593  // excess to photons
2594  if (ePhoton < totalEcal || eNeutralHadron - calibEcal < 1E-10) {
2595  if (!maxEcalRef.isNull()) {
2596  // So the merged photon energy is,
2597  mergedPhotonEnergy = ePhoton;
2598  }
2599  } else {
2600  // Otherwise assign the whole ECAL energy to the photons
2601  if (!maxEcalRef.isNull()) {
2602  // So the merged photon energy is,
2603  mergedPhotonEnergy = totalEcalEGMCalib; // KH: use calibrated ECAL energy under the egamma hypothesis
2604  }
2605  // ... and assign the remaining excess to neutral hadrons using the direction of ecal clusters
2606  mergedNeutralHadronEnergy = eNeutralHadron - calibEcal;
2607  }
2608 
2609  if (mergedPhotonEnergy > 0) {
2610  // Split merged photon into photons for each energetic ecal cluster (necessary for jet substructure reconstruction)
2611  // make only one merged photon if less than 2 ecal clusters
2612  // KH: this part still needs review, after using non-corrected ECAL energy for PF hadron calibrations
2613  if (ecalClusters.size() <= 1) {
2614  ecalClusters.clear();
2615  ecalClusters.emplace_back(
2616  maxiEcal,
2617  photonAtECAL,
2618  1.); // KH: calibration factor of 1, which should be ok as long as sumEcalClusters is consistent with photonAtECAL in this case
2619  sumEcalClusters = sqrt(photonAtECAL.Mag2());
2620  }
2621  for (auto const& pae : ecalClusters) {
2622  const double clusterEnergyCalibrated =
2623  sqrt(std::get<1>(pae).Mag2()) *
2624  std::get<2>(
2625  pae); // KH: calibrated under the egamma hypothesis. Note: sumEcalClusters is normally calibrated under egamma hypothesis
2626  particleEnergy.push_back(mergedPhotonEnergy * clusterEnergyCalibrated / sumEcalClusters);
2627  particleDirection.push_back(std::get<1>(pae));
2628  ecalEnergy.push_back(mergedPhotonEnergy * clusterEnergyCalibrated / sumEcalClusters);
2629  hcalEnergy.push_back(0.);
2630  rawecalEnergy.push_back(totalEcal);
2631  rawhcalEnergy.push_back(0.);
2632  pivotalClusterRef.push_back(elements[std::get<0>(pae)].clusterRef());
2633  iPivotal.push_back(std::get<0>(pae));
2634  }
2635  } // mergedPhotonEnergy > 0
2636 
2637  if (mergedNeutralHadronEnergy > 1.0) {
2638  // Split merged neutral hadrons according to directions of energetic ecal clusters (necessary for jet substructure reconstruction)
2639  // make only one merged neutral hadron if less than 2 ecal clusters
2640  if (ecalClusters.size() <= 1) {
2641  ecalClusters.clear();
2642  ecalClusters.emplace_back(
2643  iHcal,
2644  hadronAtECAL,
2645  1.); // KH: calibration factor of 1, which should be ok as long as sumEcalClusters is consistent with photonAtECAL
2646  sumEcalClusters = sqrt(hadronAtECAL.Mag2());
2647  }
2648  for (auto const& pae : ecalClusters) {
2649  const double clusterEnergyCalibrated =
2650  sqrt(std::get<1>(pae).Mag2()) *
2651  std::get<2>(
2652  pae); // KH: calibrated under the egamma hypothesis. Note: sumEcalClusters is normally calibrated under egamma hypothesis
2653  particleEnergy.push_back(mergedNeutralHadronEnergy * clusterEnergyCalibrated / sumEcalClusters);
2654  particleDirection.push_back(std::get<1>(pae));
2655  ecalEnergy.push_back(0.);
2656  hcalEnergy.push_back(mergedNeutralHadronEnergy * clusterEnergyCalibrated / sumEcalClusters);
2657  rawecalEnergy.push_back(0.);
2658  rawhcalEnergy.push_back(totalHcal);
2659  pivotalClusterRef.push_back(hclusterref);
2660  iPivotal.push_back(iHcal);
2661  }
2662  } //mergedNeutralHadronEnergy > 1.0
2663 
2664  // reconstructing a merged neutral
2665  // the type of PFCandidate is known from the
2666  // reference to the pivotal cluster.
2667 
2668  for (unsigned iPivot = 0; iPivot < iPivotal.size(); ++iPivot) {
2669  if (particleEnergy[iPivot] < 0.)
2670  edm::LogWarning("PFAlgo|createCandidatesHCAL")
2671  << "ALARM = Negative energy for iPivot=" << iPivot << ", " << particleEnergy[iPivot];
2672 
2673  const bool useDirection = true;
2674  auto& neutral = (*pfCandidates_)[reconstructCluster(*pivotalClusterRef[iPivot],
2675  particleEnergy[iPivot],
2676  useDirection,
2677  particleDirection[iPivot].X(),
2678  particleDirection[iPivot].Y(),
2679  particleDirection[iPivot].Z())];
2680 
2681  neutral.setEcalEnergy(rawecalEnergy[iPivot], ecalEnergy[iPivot]);
2682  if (!useHO_) {
2683  neutral.setHcalEnergy(rawhcalEnergy[iPivot], hcalEnergy[iPivot]);
2684  neutral.setHoEnergy(0., 0.);
2685  } else { // useHO_
2686  if (rawhcalEnergy[iPivot] == 0.) { // photons should be here
2687  neutral.setHcalEnergy(0., 0.);
2688  neutral.setHoEnergy(0., 0.);
2689  } else {
2690  neutral.setHcalEnergy(max(rawhcalEnergy[iPivot] - totalHO, 0.0),
2691  hcalEnergy[iPivot] * max(1. - totalHO / rawhcalEnergy[iPivot], 0.));
2692  neutral.setHoEnergy(totalHO, totalHO * hcalEnergy[iPivot] / rawhcalEnergy[iPivot]);
2693  }
2694  }
2695  neutral.setPs1Energy(0.);
2696  neutral.setPs2Energy(0.);
2697  neutral.set_mva_nothing_gamma(-1.);
2698  // neutral.addElement(&elements[iPivotal]);
2699  // neutral.addElementInBlock(blockref, iPivotal[iPivot]);
2700  neutral.addElementInBlock(blockref, iHcal);
2701  for (unsigned iTrack : chargedHadronsInBlock) {
2702  neutral.addElementInBlock(blockref, iTrack);
2703  // Assign the position of the track at the ECAL entrance
2704  const ::math::XYZPointF& chargedPosition =
2705  dynamic_cast<const reco::PFBlockElementTrack*>(&elements[iTrack])->positionAtECALEntrance();
2706  neutral.setPositionAtECALEntrance(chargedPosition);
2707 
2708  auto myEcals = associatedEcals.equal_range(iTrack);
2709  for (auto ii = myEcals.first; ii != myEcals.second; ++ii) {
2710  unsigned iEcal = ii->second.second;
2711  if (active[iEcal])
2712  continue;
2713  neutral.addElementInBlock(blockref, iEcal);
2714  }
2715  }
2716  }
2717 
2718  } // excess of energy
2719 
2720  // will now share the hcal energy between the various charged hadron
2721  // candidates, taking into account the potential neutral hadrons
2722 
2723  //JB: The question is: we've resolved the merged photons cleanly, but how should
2724  //the remaining hadrons be assigned the remaining ecal energy?
2725  //*Temporary solution*: follow HCAL example with fractions...
2726 
2727  // remove the energy of the potential neutral hadron
2728  double totalHcalEnergyCalibrated = std::max(calibHcal - mergedNeutralHadronEnergy, 0.);
2729  // similarly for the merged photons
2730  double totalEcalEnergyCalibrated = std::max(calibEcal - mergedPhotonEnergy, 0.);
2731  // share between the charged hadrons
2732 
2733  //COLIN can compute this before
2734  // not exactly equal to sum p, this is sum E
2735  double chargedHadronsTotalEnergy = 0;
2736  for (unsigned index : chargedHadronsIndices) {
2737  reco::PFCandidate& chargedHadron = (*pfCandidates_)[index];
2738  chargedHadronsTotalEnergy += chargedHadron.energy();
2739  }
2740 
2741  for (unsigned index : chargedHadronsIndices) {
2742  reco::PFCandidate& chargedHadron = (*pfCandidates_)[index];
2743  float fraction = chargedHadron.energy() / chargedHadronsTotalEnergy;
2744 
2745  if (!useHO_) {
2746  chargedHadron.setHcalEnergy(fraction * totalHcal, fraction * totalHcalEnergyCalibrated);
2747  chargedHadron.setHoEnergy(0., 0.);
2748  } else {
2749  chargedHadron.setHcalEnergy(fraction * max(totalHcal - totalHO, 0.0),
2750  fraction * totalHcalEnergyCalibrated * (1. - totalHO / totalHcal));
2751  chargedHadron.setHoEnergy(fraction * totalHO, fraction * totalHO * totalHcalEnergyCalibrated / totalHcal);
2752  }
2753  //JB: fixing up (previously omitted) setting of ECAL energy gouzevit
2754  chargedHadron.setEcalEnergy(fraction * totalEcal, fraction * totalEcalEnergyCalibrated);
2755  }
2756 
2757  // Finally treat unused ecal satellites as photons.
2758  for (auto const& ecalSatellite : ecalSatellites) {
2759  // Ignore satellites already taken
2760  unsigned iEcal = std::get<0>(ecalSatellite.second);
2761  if (!active[iEcal])
2762  continue;
2763 
2764  // Sanity checks again (well not useful, this time!)
2765  PFBlockElement::Type type = elements[iEcal].type();
2766  assert(type == PFBlockElement::ECAL);
2767  PFClusterRef eclusterref = elements[iEcal].clusterRef();
2768  assert(!eclusterref.isNull());
2769 
2770  // Lock the cluster
2771  active[iEcal] = false;
2772 
2773  // Find the associated tracks
2774  std::multimap<double, unsigned> assTracks;
2775  block.associatedElements(iEcal, linkData, assTracks, reco::PFBlockElement::TRACK, reco::PFBlock::LINKTEST_ALL);
2776 
2777  // Create a photon
2778  double ecalClusterEnergyCalibrated =
2779  sqrt(std::get<1>(ecalSatellite.second).Mag2()) *
2780  std::get<2>(
2781  ecalSatellite.second); // KH: calibrated under the egamma hypothesis (rawEcalClusterEnergy * calibration)
2782  auto& cand = (*pfCandidates_)[reconstructCluster(*eclusterref, ecalClusterEnergyCalibrated)];
2783  cand.setEcalEnergy(eclusterref->energy(), ecalClusterEnergyCalibrated);
2784  cand.setHcalEnergy(0., 0.);
2785  cand.setHoEnergy(0., 0.);
2786  cand.setPs1Energy(associatedPSs[iEcal].first);
2787  cand.setPs2Energy(associatedPSs[iEcal].second);
2788  cand.addElementInBlock(blockref, iEcal);
2789  cand.addElementInBlock(blockref, sortedTracks.begin()->second);
2790 
2791  if (fabs(eclusterref->energy() - sqrt(std::get<1>(ecalSatellite.second).Mag2())) > 1e-3 ||
2792  fabs(eclusterref->correctedEnergy() - ecalClusterEnergyCalibrated) > 1e-3)
2793  edm::LogWarning("PFAlgo:processBlock")
2794  << "ecalCluster vs ecalSatellites look inconsistent (eCluster E, calibE, ecalSatellite E, calib E): "
2795  << eclusterref->energy() << " " << eclusterref->correctedEnergy() << " "
2796  << sqrt(std::get<1>(ecalSatellite.second).Mag2()) << " " << ecalClusterEnergyCalibrated;
2797 
2798  } // ecalSatellites
2799 
2800  } // hcalIs
2801  // end loop on hcal element iHcal= hcalIs[i]
2802  LogTrace("PFAlgo|createCandidatesHCAL") << "end of function PFAlgo::createCandidatesHCAL";
2803 }

References a, funct::abs(), cms::cuda::assert(), b, groupFilesInBlocks::block, calibration_, muons2muons_cfi::chargedHadron, RPCNoise_example::check, Calorimetry_cff::dp, dptRel_DispVtx_, MillePedeFileConverter_cfg::e, ECAL, reco::PFBlockElement::ECAL, RecoEcal_cff::ecalClusters, bookConverter::elements, PFEnergyCalibration::energyEmHad(), PFTrackAlgoTools::errorScale(), factors45_, first, HLT_FULL_cff::fraction, HCAL, ElementIndices::hcalIs, photonIsolationHIProducer_cfi::ho, DigiToRawDM_cff::HO, reco::PFBlockElement::HO, mps_fire::i, cuy::ii, isFromSecInt(), PFMuonAlgo::isIsolatedMuon(), PFMuonAlgo::isLooseMuon(), PFMuonAlgo::isMuon(), edm::Ref< C, T, F >::isNull(), phase2tkutil::isPrimary(), dqmiolumiharvest::j, reco::PFBlock::LINKTEST_ALL, LogTrace, SiStripPI::max, min(), muonECAL_, muonHCAL_, muonHO_, neutralHadronEnergyResolution(), custom_jme_cff::nMuons, nSigmaHCAL(), nSigmaTRACK_, convertSQLiteXML::ok, AlCaHLTBitMon_ParallelJobs::p, pfCandidates_, ptError_, bTagCommon_cff::ratioMax, reconstructCluster(), reconstructTrack(), rejectTracks_Bad_, rejectTracks_Step45_, edm::second(), setHcalDepthInfo(), mathSSE::sqrt(), PFTrackAlgoTools::step5(), reco::PFBlockElement::TRACK, reco::btau::trackMomentum, useHO_, x, X, BeamSpotPI::Y, and BeamSpotPI::Z.

Referenced by processBlock().

◆ createCandidatesHCALUnlinked()

void PFAlgo::createCandidatesHCALUnlinked ( const reco::PFBlock &  block,
reco::PFBlock::LinkData &  linkData,
const edm::OwnVector< reco::PFBlockElement > &  elements,
std::vector< bool > &  active,
const reco::PFBlockRef &  blockref,
ElementIndices &  inds,
std::vector< bool > &  deadArea 
)
private

Definition at line 2805 of file PFAlgo.cc.

2811  {
2812  // Processing the remaining HCAL clusters
2813  LogTrace("PFAlgo|createCandidatesHCALUnlinked")
2814  << "start of function PFAlgo::createCandidatesHCALUnlinked, hcalIs.size()=" << inds.hcalIs.size();
2815 
2816  // --------------- loop remaining hcal ------------------
2817 
2818  for (unsigned iHcal : inds.hcalIs) {
2819  // Keep ECAL and HO elements for reference in the PFCandidate
2820  std::vector<unsigned> ecalRefs;
2821  std::vector<unsigned> hoRefs;
2822 
2823  LogTrace("PFAlgo|createCandidatesHCALUnlinked") << elements[iHcal] << " ";
2824 
2825  if (!active[iHcal]) {
2826  LogTrace("PFAlgo|createCandidatesHCALUnlinked") << "not active " << iHcal;
2827  continue;
2828  }
2829 
2830  // Find the ECAL elements linked to it
2831  std::multimap<double, unsigned> ecalElems;
2832  block.associatedElements(iHcal, linkData, ecalElems, reco::PFBlockElement::ECAL, reco::PFBlock::LINKTEST_ALL);
2833 
2834  // Loop on these ECAL elements
2835  float totalEcal = 0.;
2836  float ecalMax = 0.;
2837  reco::PFClusterRef eClusterRef;
2838  for (auto const& ecal : ecalElems) {
2839  unsigned iEcal = ecal.second;
2840  double dist = ecal.first;
2841  PFBlockElement::Type type = elements[iEcal].type();
2842  assert(type == PFBlockElement::ECAL);
2843 
2844  // Check if already used
2845  if (!active[iEcal])
2846  continue;
2847 
2848  // Check the distance (one HCALPlusECAL tower, roughly)
2849  // if ( dist > 0.15 ) continue;
2850 
2851  //COLINFEB16
2852  // what could be done is to
2853  // - link by rechit.
2854  // - take in the neutral hadron all the ECAL clusters
2855  // which are within the same CaloTower, according to the distance,
2856  // except the ones which are closer to another HCAL cluster.
2857  // - all the other ECAL linked to this HCAL are photons.
2858  //
2859  // about the closest HCAL cluster.
2860  // it could maybe be easier to loop on the ECAL clusters first
2861  // to cut the links to all HCAL clusters except the closest, as is
2862  // done in the first track loop. But maybe not!
2863  // or add an helper function to the PFAlgo class to ask
2864  // if a given element is the closest of a given type to another one?
2865 
2866  // Check if not closer from another free HCAL
2867  std::multimap<double, unsigned> hcalElems;
2868  block.associatedElements(iEcal, linkData, hcalElems, reco::PFBlockElement::HCAL, reco::PFBlock::LINKTEST_ALL);
2869 
2870  const bool isClosest = std::none_of(hcalElems.begin(), hcalElems.end(), [&](auto const& hcal) {
2871  return active[hcal.second] && hcal.first < dist;
2872  });
2873 
2874  if (!isClosest)
2875  continue;
2876 
2877 #ifdef EDM_ML_DEBUG
2878  LogTrace("PFAlgo|createCandidatesHCALUnlinked")
2879  << "\telement " << elements[iEcal] << " linked with dist " << dist;
2880  LogTrace("PFAlgo|createCandidatesHCALUnlinked") << "Added to HCAL cluster to form a neutral hadron";
2881 #endif
2882 
2883  reco::PFClusterRef eclusterRef = elements[iEcal].clusterRef();
2884  assert(!eclusterRef.isNull());
2885 
2886  // KH: use raw ECAL energy for PF hadron calibration_.
2887  double ecalEnergy = eclusterRef->energy(); // ecalEnergy = eclusterRef->correctedEnergy();
2888 
2889  totalEcal += ecalEnergy;
2890  if (ecalEnergy > ecalMax) {
2891  ecalMax = ecalEnergy;
2892  eClusterRef = eclusterRef;
2893  }
2894 
2895  ecalRefs.push_back(iEcal);
2896  active[iEcal] = false;
2897 
2898  } // End loop ECAL
2899 
2900  // Now find the HO clusters linked to the HCAL cluster
2901  double totalHO = 0.;
2902  double hoMax = 0.;
2903  //unsigned jHO = 0;
2904  if (useHO_) {
2905  std::multimap<double, unsigned> hoElems;
2906  block.associatedElements(iHcal, linkData, hoElems, reco::PFBlockElement::HO, reco::PFBlock::LINKTEST_ALL);
2907 
2908  // Loop on these HO elements
2909  // double totalHO = 0.;
2910  // double hoMax = 0.;
2911  // unsigned jHO = 0;
2912  reco::PFClusterRef hoClusterRef;
2913  for (auto const& ho : hoElems) {
2914  unsigned iHO = ho.second;
2915  double dist = ho.first;
2916  PFBlockElement::Type type = elements[iHO].type();
2917  assert(type == PFBlockElement::HO);
2918 
2919  // Check if already used
2920  if (!active[iHO])
2921  continue;
2922 
2923  // Check the distance (one HCALPlusHO tower, roughly)
2924  // if ( dist > 0.15 ) continue;
2925 
2926  // Check if not closer from another free HCAL
2927  std::multimap<double, unsigned> hcalElems;
2928  block.associatedElements(iHO, linkData, hcalElems, reco::PFBlockElement::HCAL, reco::PFBlock::LINKTEST_ALL);
2929 
2930  const bool isClosest = std::none_of(hcalElems.begin(), hcalElems.end(), [&](auto const& hcal) {
2931  return active[hcal.second] && hcal.first < dist;
2932  });
2933 
2934  if (!isClosest)
2935  continue;
2936 
2937 #ifdef EDM_ML_DEBUG
2938  if (useHO_) {
2939  LogTrace("PFAlgo|createCandidatesHCALUnlinked")
2940  << "\telement " << elements[iHO] << " linked with dist " << dist;
2941  LogTrace("PFAlgo|createCandidatesHCALUnlinked") << "Added to HCAL cluster to form a neutral hadron";
2942  }
2943 #endif
2944 
2945  reco::PFClusterRef hoclusterRef = elements[iHO].clusterRef();
2946  assert(!hoclusterRef.isNull());
2947 
2948  double hoEnergy =
2949  hoclusterRef->energy(); // calibration_.energyEm(*hoclusterRef,ps1Ene,ps2Ene,crackCorrection);
2950 
2951  totalHO += hoEnergy;
2952  if (hoEnergy > hoMax) {
2953  hoMax = hoEnergy;
2954  hoClusterRef = hoclusterRef;
2955  //jHO = iHO;
2956  }
2957 
2958  hoRefs.push_back(iHO);
2959  active[iHO] = false;
2960 
2961  } // End loop HO
2962  }
2963 
2964  PFClusterRef hclusterRef = elements[iHcal].clusterRef();
2965  assert(!hclusterRef.isNull());
2966 
2967  // HCAL energy
2968  double totalHcal = hclusterRef->energy();
2969  // Include the HO energy
2970  if (useHO_)
2971  totalHcal += totalHO;
2972 
2973  // Calibration
2974  double calibEcal = totalEcal > 0. ? totalEcal : 0.;
2975  double calibHcal = std::max(0., totalHcal);
2976  if (hclusterRef->layer() == PFLayer::HF_HAD || hclusterRef->layer() == PFLayer::HF_EM) {
2977  calibEcal = totalEcal;
2978  } else {
2979  calibration_.energyEmHad(
2980  -1., calibEcal, calibHcal, hclusterRef->positionREP().Eta(), hclusterRef->positionREP().Phi());
2981  }
2982 
2983  auto& cand = (*pfCandidates_)[reconstructCluster(*hclusterRef, calibEcal + calibHcal)];
2984 
2985  cand.setEcalEnergy(totalEcal, calibEcal);
2986  if (!useHO_) {
2987  cand.setHcalEnergy(totalHcal, calibHcal);
2988  cand.setHoEnergy(0., 0.);
2989  } else {
2990  cand.setHcalEnergy(max(totalHcal - totalHO, 0.0), calibHcal * (1. - totalHO / totalHcal));
2991  cand.setHoEnergy(totalHO, totalHO * calibHcal / totalHcal);
2992  }
2993  cand.setPs1Energy(0.);
2994  cand.setPs2Energy(0.);
2995  cand.addElementInBlock(blockref, iHcal);
2996  for (auto const& ec : ecalRefs)
2997  cand.addElementInBlock(blockref, ec);
2998  for (auto const& ho : hoRefs)
2999  cand.addElementInBlock(blockref, ho);
3000 
3001  } //loop hcal elements
3002 }

References cms::cuda::assert(), groupFilesInBlocks::block, calibration_, ECAL, reco::PFBlockElement::ECAL, bookConverter::elements, PFEnergyCalibration::energyEmHad(), reco::PFBlockElement::HCAL, hltEgammaHLTExtra_cfi::hcal, ElementIndices::hcalIs, PFLayer::HF_EM, PFLayer::HF_HAD, photonIsolationHIProducer_cfi::ho, DigiToRawDM_cff::HO, reco::PFBlockElement::HO, edm::Ref< C, T, F >::isNull(), reco::PFBlock::LINKTEST_ALL, LogTrace, SiStripPI::max, reconstructCluster(), and useHO_.

Referenced by processBlock().

◆ createCandidatesHF()

void PFAlgo::createCandidatesHF ( const reco::PFBlock &  block,
reco::PFBlock::LinkData &  linkData,
const edm::OwnVector< reco::PFBlockElement > &  elements,
std::vector< bool > &  active,
const reco::PFBlockRef &  blockref,
ElementIndices &  inds 
)
private

Definition at line 1251 of file PFAlgo.cc.

1256  {
1257  LogTrace("PFAlgo|createCandidatesHF") << "starting function PFAlgo::createCandidatesHF";
1258 
1259  bool trackInBlock = !inds.trackIs.empty();
1260  // inds.trackIs can be empty, even if there are tracks in this block,
1261  // but what we want to check is if this block has any track including inactive ones
1262  if (!trackInBlock)
1263  for (unsigned iEle = 0; iEle < elements.size(); iEle++) {
1264  PFBlockElement::Type type = elements[iEle].type();
1265  if (type == PFBlockElement::TRACK) {
1266  trackInBlock = true;
1267  break;
1268  }
1269  }
1270  // there is at least one HF element in this block.
1271  // in case of no track, all elements must be HF
1272  if (!trackInBlock)
1273  assert(inds.hfEmIs.size() + inds.hfHadIs.size() == elements.size());
1274 
1275  //
1276  // Dealing with a block with at least one track
1277  // Occasionally, there are only inactive tracks and multiple HF clusters. Consider such blocks too.
1278  //
1279  if (trackInBlock) { // count any tracks (not only active tracks)
1280  // sorted tracks associated with a HfHad cluster
1281  std::multimap<double, unsigned> sortedTracks;
1282  std::multimap<double, unsigned> sortedTracksActive; // only active ones
1283  // HfEms associated with tracks linked to a HfHad cluster
1284  std::multimap<unsigned, std::pair<double, unsigned>> associatedHfEms;
1285  // Temporary map for HfEm satellite clusters
1286  std::multimap<double, std::pair<unsigned, double>> hfemSatellites;
1287 
1288  //
1289  // Loop over active HfHad clusters
1290  //
1291  for (unsigned iHfHad : inds.hfHadIs) {
1292  PFBlockElement::Type type = elements[iHfHad].type();
1293  assert(type == PFBlockElement::HFHAD);
1294 
1295  PFClusterRef hclusterRef = elements[iHfHad].clusterRef();
1296  assert(!hclusterRef.isNull());
1297 
1298  sortedTracks.clear();
1299  sortedTracksActive.clear();
1300  associatedHfEms.clear();
1301  hfemSatellites.clear();
1302 
1303  // Look for associated tracks
1304  block.associatedElements(
1305  iHfHad, linkData, sortedTracks, reco::PFBlockElement::TRACK, reco::PFBlock::LINKTEST_ALL);
1306 
1307  LogTrace("PFAlgo|createCandidatesHF") << "elements[" << iHfHad << "]=" << elements[iHfHad];
1308 
1309  if (sortedTracks.empty()) {
1310  LogTrace("PFAlgo|createCandidatesHCF") << "\tno associated tracks, keep for later";
1311  continue;
1312  }
1313 
1314  // Lock the HFHAD cluster
1315  active[iHfHad] = false;
1316 
1317  LogTrace("PFAlgo|createCandidatesHF") << sortedTracks.size() << " associated tracks:";
1318 
1319  double totalChargedMomentum = 0.;
1320  double sumpError2 = 0.;
1321 
1322  //
1323  // Loop over all tracks associated to this HFHAD cluster
1324  //
1325  for (auto const& trk : sortedTracks) {
1326  unsigned iTrack = trk.second;
1327 
1328  // Check the track has not already been used
1329  if (!active[iTrack])
1330  continue;
1331  // Sanity check 1
1332  PFBlockElement::Type type = elements[iTrack].type();
1333  assert(type == reco::PFBlockElement::TRACK);
1334  // Sanity check 2
1335  reco::TrackRef trackRef = elements[iTrack].trackRef();
1336  assert(!trackRef.isNull());
1337 
1338  // Introduce tracking errors
1339  double trackMomentum = trackRef->p();
1340  totalChargedMomentum += trackMomentum;
1341 
1342  // Also keep the total track momentum error for comparison with the calo energy
1343  double dp = trackRef->qoverpError() * trackMomentum * trackMomentum;
1344  sumpError2 += dp * dp;
1345 
1346  // Store active tracks for 2nd loop to create charged hadrons
1347  sortedTracksActive.emplace(trk);
1348 
1349  // look for HFEM elements associated to iTrack (associated to iHfHad)
1350  std::multimap<double, unsigned> sortedHfEms;
1351  block.associatedElements(
1352  iTrack, linkData, sortedHfEms, reco::PFBlockElement::HFEM, reco::PFBlock::LINKTEST_ALL);
1353 
1354  LogTrace("PFAlgo|createCandidatesHF") << "number of HfEm elements linked to this track: " << sortedHfEms.size();
1355 
1356  bool connectedToHfEm = false; // Will become true if there is at least one HFEM cluster connected
1357 
1358  //
1359  // Loop over all HFEM clusters connected to iTrack
1360  //
1361  for (auto const& hfem : sortedHfEms) {
1362  unsigned iHfEm = hfem.second;
1363  double dist = hfem.first;
1364 
1365  // Ignore HFEM cluters already used
1366  if (!active[iHfEm]) {
1367  LogTrace("PFAlgo|createCandidatesHF") << "cluster locked";
1368  continue;
1369  }
1370 
1371  // Sanity checks
1372  PFBlockElement::Type type = elements[iHfEm].type();
1373  assert(type == PFBlockElement::HFEM);
1374  PFClusterRef eclusterRef = elements[iHfEm].clusterRef();
1375  assert(!eclusterRef.isNull());
1376 
1377  // Check if this HFEM is not closer to another track - ignore it in that case
1378  std::multimap<double, unsigned> sortedTracksHfEm;
1379  block.associatedElements(
1380  iHfEm, linkData, sortedTracksHfEm, reco::PFBlockElement::TRACK, reco::PFBlock::LINKTEST_ALL);
1381  unsigned jTrack = sortedTracksHfEm.begin()->second;
1382  if (jTrack != iTrack)
1383  continue;
1384 
1385  double distHfEm = block.dist(jTrack, iHfEm, linkData, reco::PFBlock::LINKTEST_ALL);
1386  double hfemEnergy = eclusterRef->energy();
1387 
1388  if (!connectedToHfEm) { // This is the closest HFEM cluster - will add its energy later
1389 
1390  LogTrace("PFAlgo|createCandidatesHF") << "closest: " << elements[iHfEm];
1391  connectedToHfEm = true;
1392  std::pair<unsigned, double> satellite(iHfEm, hfemEnergy);
1393  hfemSatellites.emplace(-1., satellite);
1394 
1395  } else { // Keep satellite clusters for later
1396 
1397  // KH: same as above.
1398  std::pair<unsigned, double> satellite(iHfEm, hfemEnergy);
1399  hfemSatellites.emplace(dist, satellite);
1400  }
1401 
1402  std::pair<double, unsigned> associatedHfEm(distHfEm, iHfEm);
1403  associatedHfEms.emplace(iTrack, associatedHfEm);
1404 
1405  } // End loop hfem associated to iTrack
1406  } // sortedTracks
1407 
1408  // HfHad energy
1409  double uncalibratedenergyHfHad = hclusterRef->energy();
1410  double energyHfHad = uncalibratedenergyHfHad;
1411  if (thepfEnergyCalibrationHF_.getcalibHF_use()) {
1412  energyHfHad = thepfEnergyCalibrationHF_.energyHad( // HAD only calibration
1413  uncalibratedenergyHfHad,
1414  hclusterRef->positionREP().Eta(),
1415  hclusterRef->positionREP().Phi());
1416  }
1417  double calibFactorHfHad = (uncalibratedenergyHfHad > 0.) ? energyHfHad / uncalibratedenergyHfHad : 1.;
1418 
1419  // HfEm energy
1420  double energyHfEmTmp = 0.;
1421  double uncalibratedenergyHfEmTmp = 0.;
1422  double energyHfEm = 0.;
1423  double uncalibratedenergyHfEm = 0.;
1424 
1425  // estimated HF resolution and track p error
1426  double caloResolution = hfEnergyResolution(totalChargedMomentum);
1427  caloResolution *= totalChargedMomentum;
1428  double totalError = sqrt(caloResolution * caloResolution + sumpError2);
1429  double nsigmaHFEM = nSigmaHFEM(totalChargedMomentum);
1430  double nsigmaHFHAD = nSigmaHFHAD(totalChargedMomentum);
1431 
1432  // Handle case that no active track gets associated to HfHad cluster
1433  if (sortedTracksActive.empty()) {
1434  // look for HFEM elements associated to iHfHad
1435  std::multimap<double, unsigned> sortedHfEms;
1436  std::multimap<double, unsigned> sortedHfEmsActive;
1437  block.associatedElements(
1438  iHfHad, linkData, sortedHfEms, reco::PFBlockElement::HFEM, reco::PFBlock::LINKTEST_ALL);
1439  //
1440  // If iHfHad is connected to HFEM cluster, Loop over all of them
1441  //
1442  if (!sortedHfEms.empty()) {
1443  for (auto const& hfem : sortedHfEms) {
1444  unsigned iHfEm = hfem.second;
1445  // Ignore HFEM cluters already used
1446  if (!active[iHfEm])
1447  continue;
1448  sortedHfEmsActive.emplace(hfem);
1449  PFClusterRef eclusterRef = elements[iHfEm].clusterRef();
1450  assert(!eclusterRef.isNull());
1451  double hfemEnergy = eclusterRef->energy();
1452  uncalibratedenergyHfEm += hfemEnergy;
1453  energyHfEm = uncalibratedenergyHfEm;
1454  if (thepfEnergyCalibrationHF_.getcalibHF_use()) {
1455  energyHfEm = thepfEnergyCalibrationHF_.energyEmHad(
1456  uncalibratedenergyHfEm, 0.0, eclusterRef->positionREP().Eta(), eclusterRef->positionREP().Phi());
1457  energyHfHad = thepfEnergyCalibrationHF_.energyEmHad(
1458  0.0, uncalibratedenergyHfHad, hclusterRef->positionREP().Eta(), hclusterRef->positionREP().Phi());
1459  } // calib true
1460  } // loop over sortedHfEm
1461  } // if !sortedHfEms.empty()
1462  //
1463  // Create HF candidates
1464  unsigned tmpi = reconstructCluster(*hclusterRef, energyHfEm + energyHfHad);
1465  (*pfCandidates_)[tmpi].setHcalEnergy(uncalibratedenergyHfHad, energyHfHad);
1466  (*pfCandidates_)[tmpi].setEcalEnergy(uncalibratedenergyHfEm, energyHfEm);
1467  (*pfCandidates_)[tmpi].addElementInBlock(blockref, iHfHad);
1468  for (auto const& hfem : sortedHfEmsActive) {
1469  unsigned iHfEm = hfem.second;
1470  (*pfCandidates_)[tmpi].addElementInBlock(blockref, iHfEm);
1471  active[iHfEm] = false;
1472  }
1473 
1474  } // if sortedTracksActive.empty() ends
1475  //
1476  // Active tracks are associated.
1477  // Create HFHAD candidates from excess energy w.r.t. tracks
1478  else if ((energyHfHad - totalChargedMomentum) > nsigmaHFHAD * totalError) { // HfHad is excessive
1479  assert(energyHfEm == 0.);
1480  // HfHad candidate from excess
1481  double energyHfHadExcess = max(energyHfHad - totalChargedMomentum, 0.);
1482  double uncalibratedenergyHfHadExcess = energyHfHadExcess / calibFactorHfHad;
1483  unsigned tmpi = reconstructCluster(*hclusterRef, energyHfHadExcess);
1484  (*pfCandidates_)[tmpi].setHcalEnergy(uncalibratedenergyHfHadExcess, energyHfHadExcess);
1485  (*pfCandidates_)[tmpi].setEcalEnergy(0., 0.);
1486  (*pfCandidates_)[tmpi].addElementInBlock(blockref, iHfHad);
1487  energyHfHad = max(energyHfHad - energyHfHadExcess, 0.);
1488  uncalibratedenergyHfHad = max(uncalibratedenergyHfHad - uncalibratedenergyHfHadExcess, 0.);
1489  }
1490  //
1491  // If there is a room for HFEM satellites to get associated,
1492  // loop over all HFEM satellites, starting for the closest to the various tracks
1493  // and adding other satellites until saturation of the total track momentum
1494  //
1495  else {
1496  for (auto const& hfemSatellite : hfemSatellites) {
1497  //
1498  uncalibratedenergyHfEmTmp += std::get<1>(hfemSatellite.second); // KH: raw HFEM energy
1499  energyHfEmTmp = uncalibratedenergyHfEmTmp;
1500  double energyHfHadTmp = uncalibratedenergyHfHad; // now to test hfhad calibration with EM+HAD cases
1501  unsigned iHfEm = std::get<0>(hfemSatellite.second);
1502  PFClusterRef eclusterRef = elements[iHfEm].clusterRef();
1503  assert(!eclusterRef.isNull());
1504  if (thepfEnergyCalibrationHF_.getcalibHF_use()) {
1505  energyHfEmTmp = thepfEnergyCalibrationHF_.energyEmHad(
1506  uncalibratedenergyHfEmTmp, 0.0, eclusterRef->positionREP().Eta(), eclusterRef->positionREP().Phi());
1507  energyHfHadTmp = thepfEnergyCalibrationHF_.energyEmHad(
1508  0.0, uncalibratedenergyHfHad, hclusterRef->positionREP().Eta(), hclusterRef->positionREP().Phi());
1509  }
1510 
1511  double caloEnergyTmp = energyHfEmTmp + energyHfHadTmp;
1512  double calibFactorHfEm = (uncalibratedenergyHfEmTmp > 0.) ? energyHfEmTmp / uncalibratedenergyHfEmTmp : 1.;
1513 
1514  // Continue looping until all closest clusters are exhausted and as long as
1515  // the calorimetric energy does not saturate the total momentum.
1516  if (hfemSatellite.first < 0. || caloEnergyTmp < totalChargedMomentum) {
1517  LogTrace("PFAlgo|createCandidatesHF")
1518  << "\t\t\tactive, adding " << std::get<1>(hfemSatellite.second) << " to HFEM energy, and locking";
1519  active[std::get<0>(hfemSatellite.second)] = false;
1520  // HfEm is excessive (possible for the first hfemSatellite)
1521  if (hfemSatellite.first < 0. && (caloEnergyTmp - totalChargedMomentum) > nsigmaHFEM * totalError) {
1522  // HfEm candidate from excess
1523  double energyHfEmExcess = max(caloEnergyTmp - totalChargedMomentum, 0.);
1524  double uncalibratedenergyHfEmExcess = energyHfEmExcess / calibFactorHfEm;
1525  unsigned tmpi = reconstructCluster(*eclusterRef, energyHfEmExcess);
1526  (*pfCandidates_)[tmpi].setEcalEnergy(uncalibratedenergyHfEmExcess, energyHfEmExcess);
1527  (*pfCandidates_)[tmpi].setHcalEnergy(0, 0.);
1528  (*pfCandidates_)[tmpi].addElementInBlock(blockref, iHfEm);
1529  energyHfEmTmp = max(energyHfEmTmp - energyHfEmExcess, 0.);
1530  uncalibratedenergyHfEmTmp = max(uncalibratedenergyHfEmTmp - uncalibratedenergyHfEmExcess, 0.);
1531  }
1532  energyHfEm = energyHfEmTmp;
1533  uncalibratedenergyHfEm = uncalibratedenergyHfEmTmp;
1534  energyHfHad = energyHfHadTmp;
1535  continue;
1536  }
1537  break;
1538  } // loop over hfemsattelites ends
1539  } // if HFHAD is excessive or not
1540 
1541  //
1542  // Loop over all tracks associated to this HFHAD cluster *again* in order to produce charged hadrons
1543  //
1544  for (auto const& trk : sortedTracksActive) {
1545  unsigned iTrack = trk.second;
1546 
1547  // Sanity check
1548  reco::TrackRef trackRef = elements[iTrack].trackRef();
1549  assert(!trackRef.isNull());
1550 
1551  //
1552  // Reconstructing charged hadrons
1553  //
1554  unsigned tmpi = reconstructTrack(elements[iTrack]);
1555  active[iTrack] = false;
1556  (*pfCandidates_)[tmpi].addElementInBlock(blockref, iHfHad);
1557  auto myHfEms = associatedHfEms.equal_range(iTrack);
1558  for (auto ii = myHfEms.first; ii != myHfEms.second; ++ii) {
1559  unsigned iHfEm = ii->second.second;
1560  if (active[iHfEm])
1561  continue;
1562  (*pfCandidates_)[tmpi].addElementInBlock(blockref, iHfEm);
1563  }
1564  double frac = 0.;
1565  if (totalChargedMomentum)
1566  frac = trackRef->p() / totalChargedMomentum;
1567  (*pfCandidates_)[tmpi].setEcalEnergy(uncalibratedenergyHfEm * frac, energyHfEm * frac);
1568  (*pfCandidates_)[tmpi].setHcalEnergy(uncalibratedenergyHfHad * frac, energyHfHad * frac);
1569 
1570  } // sortedTracks loop ends
1571 
1572  } // iHfHad element loop ends
1573 
1574  //
1575  // Loop over remaining active HfEm clusters
1576  //
1577  for (unsigned iHfEm = 0; iHfEm < elements.size(); iHfEm++) {
1578  PFBlockElement::Type type = elements[iHfEm].type();
1579  if (type == PFBlockElement::HFEM && active[iHfEm]) {
1580  reco::PFClusterRef eclusterRef = elements[iHfEm].clusterRef();
1581  double energyHF = 0.;
1582  double uncalibratedenergyHF = 0.;
1583  unsigned tmpi = 0;
1584  // do EM-only calibration here
1585  energyHF = eclusterRef->energy();
1586  uncalibratedenergyHF = energyHF;
1587  if (thepfEnergyCalibrationHF_.getcalibHF_use()) {
1588  energyHF = thepfEnergyCalibrationHF_.energyEm(
1589  uncalibratedenergyHF, eclusterRef->positionREP().Eta(), eclusterRef->positionREP().Phi());
1590  }
1591  tmpi = reconstructCluster(*eclusterRef, energyHF);
1592  (*pfCandidates_)[tmpi].setEcalEnergy(uncalibratedenergyHF, energyHF);
1593  (*pfCandidates_)[tmpi].setHcalEnergy(0., 0.);
1594  (*pfCandidates_)[tmpi].addElementInBlock(blockref, iHfEm);
1595  active[iHfEm] = false;
1596  LogTrace("PFAlgo|createCandidatesHF") << "HF EM alone from blocks with tracks! " << energyHF;
1597  }
1598  } // remaining active HfEm cluster loop ends
1599 
1600  } // if-statement for blocks including tracks ends here
1601  //
1602  // -----------------------------------------------
1603  // From here, traditional PF HF candidate creation
1604  // -----------------------------------------------
1605  //
1606  else if (elements.size() == 1) {
1607  //Auguste: HAD-only calibration here
1608  reco::PFClusterRef clusterRef = elements[0].clusterRef();
1609  double energyHF = 0.;
1610  double uncalibratedenergyHF = 0.;
1611  unsigned tmpi = 0;
1612  switch (clusterRef->layer()) {
1613  case PFLayer::HF_EM:
1614  // do EM-only calibration here
1615  energyHF = clusterRef->energy();
1616  uncalibratedenergyHF = energyHF;
1617  if (thepfEnergyCalibrationHF_.getcalibHF_use()) {
1618  energyHF = thepfEnergyCalibrationHF_.energyEm(
1619  uncalibratedenergyHF, clusterRef->positionREP().Eta(), clusterRef->positionREP().Phi());
1620  }
1621  tmpi = reconstructCluster(*clusterRef, energyHF);
1622  (*pfCandidates_)[tmpi].setEcalEnergy(uncalibratedenergyHF, energyHF);
1623  (*pfCandidates_)[tmpi].setHcalEnergy(0., 0.);
1624  (*pfCandidates_)[tmpi].setHoEnergy(0., 0.);
1625  (*pfCandidates_)[tmpi].setPs1Energy(0.);
1626  (*pfCandidates_)[tmpi].setPs2Energy(0.);
1627  (*pfCandidates_)[tmpi].addElementInBlock(blockref, inds.hfEmIs[0]);
1628  LogTrace("PFAlgo|createCandidatesHF") << "HF EM alone ! " << energyHF;
1629  break;
1630  case PFLayer::HF_HAD:
1631  // do HAD-only calibration here
1632  energyHF = clusterRef->energy();
1633  uncalibratedenergyHF = energyHF;
1634  if (thepfEnergyCalibrationHF_.getcalibHF_use()) {
1635  energyHF = thepfEnergyCalibrationHF_.energyHad(
1636  uncalibratedenergyHF, clusterRef->positionREP().Eta(), clusterRef->positionREP().Phi());
1637  }
1638  tmpi = reconstructCluster(*clusterRef, energyHF);
1639  (*pfCandidates_)[tmpi].setHcalEnergy(uncalibratedenergyHF, energyHF);
1640  (*pfCandidates_)[tmpi].setEcalEnergy(0., 0.);
1641  (*pfCandidates_)[tmpi].setHoEnergy(0., 0.);
1642  (*pfCandidates_)[tmpi].setPs1Energy(0.);
1643  (*pfCandidates_)[tmpi].setPs2Energy(0.);
1644  (*pfCandidates_)[tmpi].addElementInBlock(blockref, inds.hfHadIs[0]);
1645  LogTrace("PFAlgo|createCandidatesHF") << "HF Had alone ! " << energyHF;
1646  break;
1647  default:
1648  assert(0);
1649  }
1650  } else if (elements.size() == 2) {
1651  //Auguste: EM + HAD calibration here
1652  reco::PFClusterRef c0 = elements[0].clusterRef();
1653  reco::PFClusterRef c1 = elements[1].clusterRef();
1654  // 2 HF elements. Must be in each layer.
1655  reco::PFClusterRef cem = (c0->layer() == PFLayer::HF_EM ? c0 : c1);
1656  reco::PFClusterRef chad = (c1->layer() == PFLayer::HF_HAD ? c1 : c0);
1657 
1658  if (cem->layer() != PFLayer::HF_EM || chad->layer() != PFLayer::HF_HAD) {
1659  edm::LogError("PFAlgo::createCandidatesHF") << "Error: 2 elements, but not 1 HFEM and 1 HFHAD";
1660  edm::LogError("PFAlgo::createCandidatesHF") << block;
1661  assert(0);
1662  // assert( c1->layer()== PFLayer::HF_EM &&
1663  // c0->layer()== PFLayer::HF_HAD );
1664  }
1665  // do EM+HAD calibration here
1666  double energyHfEm = cem->energy();
1667  double energyHfHad = chad->energy();
1668  double uncalibratedenergyHfEm = energyHfEm;
1669  double uncalibratedenergyHfHad = energyHfHad;
1670  if (thepfEnergyCalibrationHF_.getcalibHF_use()) {
1671  energyHfEm = thepfEnergyCalibrationHF_.energyEmHad(
1672  uncalibratedenergyHfEm, 0.0, c0->positionREP().Eta(), c0->positionREP().Phi());
1673  energyHfHad = thepfEnergyCalibrationHF_.energyEmHad(
1674  0.0, uncalibratedenergyHfHad, c1->positionREP().Eta(), c1->positionREP().Phi());
1675  }
1676  auto& cand = (*pfCandidates_)[reconstructCluster(*chad, energyHfEm + energyHfHad)];
1677  cand.setEcalEnergy(uncalibratedenergyHfEm, energyHfEm);
1678  cand.setHcalEnergy(uncalibratedenergyHfHad, energyHfHad);
1679  cand.setHoEnergy(0., 0.);
1680  cand.setPs1Energy(0.);
1681  cand.setPs2Energy(0.);
1682  cand.addElementInBlock(blockref, inds.hfEmIs[0]);
1683  cand.addElementInBlock(blockref, inds.hfHadIs[0]);
1684  LogTrace("PFAlgo|createCandidatesHF") << "HF EM+HAD found ! " << energyHfEm << " " << energyHfHad;
1685  } else {
1686  // Unusual blocks including HF elements, but do not fit any of the above categories
1687  edm::LogWarning("PFAlgo::createCandidatesHF")
1688  << "Warning: HF, but n elem different from 1 or 2 or >=3 or !trackIs.empty()";
1689  edm::LogWarning("PFAlgo::createCandidatesHF") << block;
1690  }
1691  LogTrace("PFAlgo|createCandidateHF") << "end of function PFAlgo::createCandidateHF";
1692 }

References cms::cuda::assert(), groupFilesInBlocks::block, fftjetpileupestimator_calo_uncalib_cfi::c0, alignmentValidation::c1, Calorimetry_cff::dp, bookConverter::elements, PFEnergyCalibrationHF::energyEm(), PFEnergyCalibrationHF::energyEmHad(), PFEnergyCalibrationHF::energyHad(), HcalResponse_cfi::energyHF, DivergingColor::frac, PFEnergyCalibrationHF::getcalibHF_use(), PFLayer::HF_EM, PFLayer::HF_HAD, reco::PFBlockElement::HFEM, ElementIndices::hfEmIs, hfEnergyResolution(), ElementIndices::hfHadIs, cuy::ii, edm::Ref< C, T, F >::isNull(), reco::PFBlock::LINKTEST_ALL, LogTrace, SiStripPI::max, nSigmaHFEM(), nSigmaHFHAD(), reconstructCluster(), reconstructTrack(), mathSSE::sqrt(), thepfEnergyCalibrationHF_, reco::PFBlockElement::TRACK, ElementIndices::trackIs, and reco::btau::trackMomentum.

Referenced by processBlock().

◆ decideType()

int PFAlgo::decideType ( const edm::OwnVector< reco::PFBlockElement > &  elements,
const reco::PFBlockElement::Type  type,
std::vector< bool > &  active,
ElementIndices &  inds,
std::vector< bool > &  deadArea,
unsigned int  iEle 
)
private

Definition at line 1193 of file PFAlgo.cc.

1198  {
1199  switch (type) {
1200  case PFBlockElement::TRACK:
1201  if (active[iEle]) {
1202  inds.trackIs.push_back(iEle);
1203  LogTrace("PFAlgo|decideType") << "TRACK, stored index, continue";
1204  }
1205  break;
1206  case PFBlockElement::ECAL:
1207  if (active[iEle]) {
1208  inds.ecalIs.push_back(iEle);
1209  LogTrace("PFAlgo|decideType") << "ECAL, stored index, continue";
1210  }
1211  return 1; //continue
1212  case PFBlockElement::HCAL:
1213  if (active[iEle]) {
1214  if (elements[iEle].clusterRef()->flags() & reco::CaloCluster::badHcalMarker) {
1215  LogTrace("PFAlgo|decideType") << "HCAL DEAD AREA: remember and skip.";
1216  active[iEle] = false;
1217  deadArea[iEle] = true;
1218  return 1; //continue
1219  }
1220  inds.hcalIs.push_back(iEle);
1221  LogTrace("PFAlgo|decideType") << "HCAL, stored index, continue";
1222  }
1223  return 1; //continue
1224  case PFBlockElement::HO:
1225  if (useHO_) {
1226  if (active[iEle]) {
1227  inds.hoIs.push_back(iEle);
1228  LogTrace("PFAlgo|decideType") << "HO, stored index, continue";
1229  }
1230  }
1231  return 1; //continue
1232  case PFBlockElement::HFEM:
1233  if (active[iEle]) {
1234  inds.hfEmIs.push_back(iEle);
1235  LogTrace("PFAlgo|decideType") << "HFEM, stored index, continue";
1236  }
1237  return 1; //continue
1238  case PFBlockElement::HFHAD:
1239  if (active[iEle]) {
1240  inds.hfHadIs.push_back(iEle);
1241  LogTrace("PFAlgo|decideType") << "HFHAD, stored index, continue";
1242  }
1243  return 1; //continue
1244  default:
1245  return 1; //continue
1246  }
1247  LogTrace("PFAlgo|decideType") << "Did not match type to anything, return 0";
1248  return 0;
1249 }

References reco::CaloCluster::badHcalMarker, ECAL, ElementIndices::ecalIs, bookConverter::elements, HLT_FULL_cff::flags, HCAL, ElementIndices::hcalIs, ElementIndices::hfEmIs, ElementIndices::hfHadIs, DigiToRawDM_cff::HO, ElementIndices::hoIs, LogTrace, ElementIndices::trackIs, and useHO_.

Referenced by elementLoop().

◆ egammaFilters()

void PFAlgo::egammaFilters ( const reco::PFBlockRef &  blockref,
std::vector< bool > &  active,
PFEGammaFilters const *  pfegamma 
)
private

Definition at line 220 of file PFAlgo.cc.

222  {
223  // const edm::ValueMap<reco::GsfElectronRef> & myGedElectronValMap(*valueMapGedElectrons_);
224 
225  unsigned int negmcandidates = pfEgammaCandidates_->size();
226  LogTrace("PFAlgo|egammaFilters") << "start of function PFAlgo::egammaFilters(), negmcandidates=" << negmcandidates;
227 
228  for (unsigned int ieg = 0; ieg < negmcandidates; ++ieg) {
229  // const reco::PFCandidate & egmcand((*pfEgammaCandidates_)[ieg]);
230  reco::PFCandidateRef pfEgmRef = pfEgammaCandidates_->refAt(ieg).castTo<reco::PFCandidateRef>();
231 
232  const PFCandidate::ElementsInBlocks& theElements = (*pfEgmRef).elementsInBlocks();
233  PFCandidate::ElementsInBlocks::const_iterator iegfirst = theElements.begin();
234  bool sameBlock = false;
235  bool isGoodElectron = false;
236  bool isGoodPhoton = false;
237  bool isPrimaryElectron = false;
238  if (iegfirst->first == blockref)
239  sameBlock = true;
240  if (sameBlock) {
241  LogTrace("PFAlgo|egammaFilters") << " I am in looping on EGamma Candidates: pt " << (*pfEgmRef).pt()
242  << " eta,phi " << (*pfEgmRef).eta() << ", " << (*pfEgmRef).phi() << " charge "
243  << (*pfEgmRef).charge();
244 
245  if ((*pfEgmRef).gsfTrackRef().isNonnull()) {
246  reco::GsfElectronRef gedEleRef = (*valueMapGedElectrons_)[pfEgmRef];
247  if (gedEleRef.isNonnull()) {
248  isGoodElectron = pfegamma->passElectronSelection(*gedEleRef, *pfEgmRef, nVtx_);
249  isPrimaryElectron = pfegamma->isElectron(*gedEleRef);
250  if (isGoodElectron)
251  LogTrace("PFAlgo|egammaFilters")
252  << "** Good Electron, pt " << gedEleRef->pt() << " eta, phi " << gedEleRef->eta() << ", "
253  << gedEleRef->phi() << " charge " << gedEleRef->charge() << " isPrimary " << isPrimaryElectron
254  << " isoDr03 "
255  << (gedEleRef->dr03TkSumPt() + gedEleRef->dr03EcalRecHitSumEt() + gedEleRef->dr03HcalTowerSumEt())
256  << " mva_isolated " << gedEleRef->mva_Isolated() << " mva_e_pi " << gedEleRef->mva_e_pi();
257  }
258  }
259  if ((*pfEgmRef).superClusterRef().isNonnull()) {
260  reco::PhotonRef gedPhoRef = (*valueMapGedPhotons_)[pfEgmRef];
261  if (gedPhoRef.isNonnull()) {
262  isGoodPhoton = pfegamma->passPhotonSelection(*gedPhoRef);
263  if (isGoodPhoton)
264  LogTrace("PFAlgo|egammaFilters") << "** Good Photon, pt " << gedPhoRef->pt() << " eta, phi "
265  << gedPhoRef->eta() << ", " << gedPhoRef->phi() << endl;
266  }
267  }
268  } // end sameBlock
269 
270  if (isGoodElectron && isGoodPhoton) {
271  if (isPrimaryElectron)
272  isGoodPhoton = false;
273  else
274  isGoodElectron = false;
275  }
276 
277  // isElectron
278  if (isGoodElectron) {
279  reco::GsfElectronRef gedEleRef = (*valueMapGedElectrons_)[pfEgmRef];
280  reco::PFCandidate myPFElectron = *pfEgmRef;
281  // run protections
282  bool lockTracks = false;
283  bool isSafe = true;
284  if (useProtectionsForJetMET_) {
285  lockTracks = true;
286  isSafe = pfegamma->isElectronSafeForJetMET(*gedEleRef, myPFElectron, primaryVertex_, lockTracks);
287  }
288 
289  if (isSafe) {
290  reco::PFCandidate::ParticleType particleType = reco::PFCandidate::e;
291  myPFElectron.setParticleType(particleType);
292  myPFElectron.setCharge(gedEleRef->charge());
293  myPFElectron.setP4(gedEleRef->p4());
294  myPFElectron.set_mva_e_pi(gedEleRef->mva_e_pi());
295  myPFElectron.set_mva_Isolated(gedEleRef->mva_Isolated());
296 
297  LogTrace("PFAlgo|egammaFilters") << " PFAlgo: found an electron with NEW EGamma code ";
298  LogTrace("PFAlgo|egammaFilters") << " myPFElectron: pt " << myPFElectron.pt() << " eta,phi "
299  << myPFElectron.eta() << ", " << myPFElectron.phi() << " mva "
300  << myPFElectron.mva_e_pi() << " charge " << myPFElectron.charge();
301 
302  LogTrace("PFAlgo|egammaFilters") << " THE BLOCK " << *blockref;
303  for (auto const& eb : theElements) {
304  active[eb.second] = false;
305  LogTrace("PFAlgo|egammaFilters") << " Elements used " << eb.second;
306  }
307 
308  // The electron is considered safe for JetMET and the additional tracks pointing to it are locked
309  if (lockTracks) {
310  const PFCandidate::ElementsInBlocks& extraTracks = myPFElectron.egammaExtraRef()->extraNonConvTracks();
311  for (auto const& trk : extraTracks) {
312  LogTrace("PFAlgo|egammaFilters") << " Extra locked track " << trk.second;
313  active[trk.second] = false;
314  }
315  }
316 
317  LogTrace("PFAlgo|egammaFilters") << "Creating PF electron: pt=" << myPFElectron.pt()
318  << " eta=" << myPFElectron.eta() << " phi=" << myPFElectron.phi();
319  pfCandidates_->push_back(myPFElectron);
320 
321  } else {
322  LogTrace("PFAlgo|egammaFilters") << "PFAlgo: Electron DISCARDED, NOT SAFE FOR JETMET ";
323  }
324  } //isGoodElectron
325 
326  if (isGoodPhoton) {
327  reco::PhotonRef gedPhoRef = (*valueMapGedPhotons_)[pfEgmRef];
328  reco::PFCandidate myPFPhoton = *pfEgmRef;
329  bool isSafe = true;
330  if (useProtectionsForJetMET_) {
331  isSafe = pfegamma->isPhotonSafeForJetMET(*gedPhoRef, myPFPhoton);
332  }
333 
334  if (isSafe) {
335  reco::PFCandidate::ParticleType particleType = reco::PFCandidate::gamma;
336  myPFPhoton.setParticleType(particleType);
337  myPFPhoton.setCharge(0);
338  myPFPhoton.set_mva_nothing_gamma(1.);
339  ::math::XYZPoint v(primaryVertex_.x(), primaryVertex_.y(), primaryVertex_.z());
340  myPFPhoton.setVertex(v);
341  myPFPhoton.setP4(gedPhoRef->p4());
342  LogTrace("PFAlgo|egammaFilters") << " PFAlgo: found a photon with NEW EGamma code ";
343  LogTrace("PFAlgo|egammaFilters") << " myPFPhoton: pt " << myPFPhoton.pt() << " eta,phi " << myPFPhoton.eta()
344  << ", " << myPFPhoton.phi() << " charge " << myPFPhoton.charge();
345 
346  // Lock all the elements
347  LogTrace("PFAlgo|egammaFilters") << " THE BLOCK " << *blockref;
348  for (auto const& eb : theElements) {
349  active[eb.second] = false;
350  LogTrace("PFAlgo|egammaFilters") << " Elements used " << eb.second;
351  }
352  LogTrace("PFAlgo|egammaFilters") << "Creating PF photon: pt=" << myPFPhoton.pt() << " eta=" << myPFPhoton.eta()
353  << " phi=" << myPFPhoton.phi();
354  pfCandidates_->push_back(myPFPhoton);
355 
356  } // end isSafe
357  } // end isGoodPhoton
358  } // end loop on EGM candidates
359  LogTrace("PFAlgo|egammaFilters") << "end of function PFAlgo::egammaFilters";
360 }

References reco::LeafCandidate::charge(), reco::PFCandidate::e, reco::PFCandidate::egammaExtraRef(), reco::LeafCandidate::eta(), reco::PFCandidate::gamma, PFEGammaFilters::isElectron(), PFEGammaFilters::isElectronSafeForJetMET(), edm::Ref< C, T, F >::isNonnull(), PFEGammaFilters::isPhotonSafeForJetMET(), LogTrace, reco::PFCandidate::mva_e_pi(), nVtx_, PbPb_ZMuSkimMuonDPG_cff::particleType, PFEGammaFilters::passElectronSelection(), PFEGammaFilters::passPhotonSelection(), pfCandidates_, pfEgammaCandidates_, reco::LeafCandidate::phi(), primaryVertex_, reco::LeafCandidate::pt(), edm::View< T >::refAt(), reco::PFCandidate::set_mva_e_pi(), reco::PFCandidate::set_mva_Isolated(), reco::PFCandidate::set_mva_nothing_gamma(), reco::LeafCandidate::setCharge(), reco::LeafCandidate::setP4(), reco::PFCandidate::setParticleType(), reco::LeafCandidate::setVertex(), edm::View< T >::size(), useProtectionsForJetMET_, findQualityFiles::v, reco::Vertex::x(), reco::Vertex::y(), and reco::Vertex::z().

Referenced by processBlock().

◆ elementLoop()

void PFAlgo::elementLoop ( const reco::PFBlock &  block,
reco::PFBlock::LinkData &  linkData,
const edm::OwnVector< reco::PFBlockElement > &  elements,
std::vector< bool > &  active,
const reco::PFBlockRef &  blockref,
ElementIndices &  inds,
std::vector< bool > &  deadArea 
)
private

Definition at line 989 of file PFAlgo.cc.

995  {
996  LogTrace("PFAlgo|elementLoop") << "start of function PFAlgo::elementLoop, elements.size()" << elements.size();
997 
998  for (unsigned iEle = 0; iEle < elements.size(); iEle++) {
999  PFBlockElement::Type type = elements[iEle].type();
1000 
1001  LogTrace("PFAlgo|elementLoop") << "elements[iEle=" << iEle << "]=" << elements[iEle];
1002  //only process TRACK elements, but fill the ElementIndices vector with indices for all elements.
1003  //Mark the active & deadArea for bad HCAL
1004  auto ret_decideType = decideType(elements, type, active, inds, deadArea, iEle);
1005  if (ret_decideType == 1) {
1006  LogTrace("PFAlgo|elementLoop") << "ret_decideType==1, continuing";
1007  continue;
1008  }
1009  LogTrace("PFAlgo|elementLoop") << "ret_decideType=" << ret_decideType << " type=" << type;
1010 
1011  active[iEle] = checkAndReconstructSecondaryInteraction(blockref, elements, active[iEle], iEle);
1012 
1013  if (!active[iEle]) {
1014  LogTrace("PFAlgo|elementLoop") << "Already used by electrons, muons, conversions";
1015  continue;
1016  }
1017 
1018  reco::TrackRef trackRef = elements[iEle].trackRef();
1019  assert(!trackRef.isNull());
1020 
1021  LogTrace("PFAlgo|elementLoop") << "PFAlgo:processBlock"
1022  << " trackIs.size()=" << inds.trackIs.size()
1023  << " ecalIs.size()=" << inds.ecalIs.size() << " hcalIs.size()=" << inds.hcalIs.size()
1024  << " hoIs.size()=" << inds.hoIs.size() << " hfEmIs.size()=" << inds.hfEmIs.size()
1025  << " hfHadIs.size()=" << inds.hfHadIs.size();
1026 
1027  // look for associated elements of all types
1028  //COLINFEB16
1029  // all types of links are considered.
1030  // the elements are sorted by increasing distance
1031  std::multimap<double, unsigned> ecalElems;
1032  block.associatedElements(iEle, linkData, ecalElems, reco::PFBlockElement::ECAL, reco::PFBlock::LINKTEST_ALL);
1033 
1034  std::multimap<double, unsigned> hcalElems;
1035  block.associatedElements(iEle, linkData, hcalElems, reco::PFBlockElement::HCAL, reco::PFBlock::LINKTEST_ALL);
1036 
1037  std::multimap<double, unsigned> hfEmElems;
1038  std::multimap<double, unsigned> hfHadElems;
1039  block.associatedElements(iEle, linkData, hfEmElems, reco::PFBlockElement::HFEM, reco::PFBlock::LINKTEST_ALL);
1040  block.associatedElements(iEle, linkData, hfHadElems, reco::PFBlockElement::HFHAD, reco::PFBlock::LINKTEST_ALL);
1041 
1042  LogTrace("PFAlgo|elementLoop") << "\tTrack " << iEle << " is linked to " << ecalElems.size() << " ecal and "
1043  << hcalElems.size() << " hcal and " << hfEmElems.size() << " hfEm and "
1044  << hfHadElems.size() << " hfHad elements";
1045 
1046 #ifdef EDM_ML_DEBUG
1047  for (const auto& pair : ecalElems) {
1048  LogTrace("PFAlgo|elementLoop") << "ecal: dist " << pair.first << "\t elem " << pair.second;
1049  }
1050  for (const auto& pair : hcalElems) {
1051  LogTrace("PFAlgo|elementLoop") << "hcal: dist " << pair.first << "\t elem " << pair.second
1052  << (deadArea[pair.second] ? " DEAD AREA MARKER" : "");
1053  }
1054 #endif
1055 
1056  const bool hasDeadHcal = checkHasDeadHcal(hcalElems, deadArea);
1057  if (hasDeadHcal) {
1058  hcalElems.clear();
1059  }
1060  const bool goodTrackDeadHcal = checkGoodTrackDeadHcal(trackRef, hasDeadHcal);
1061 
1062  // When a track has no HCAL cluster linked, but another track is linked to the same
1063  // ECAL cluster and an HCAL cluster, link the track to the HCAL cluster for
1064  // later analysis
1065  if (hcalElems.empty() && !ecalElems.empty() && !hasDeadHcal) {
1066  relinkTrackToHcal(block, ecalElems, hcalElems, active, linkData, iEle);
1067  }
1068 
1069  //MICHELE
1070  //TEMPORARY SOLUTION FOR ELECTRON REJECTION IN PFTAU
1071  //COLINFEB16
1072  // in case particle flow electrons are not reconstructed,
1073  // the mva_e_pi of the charged hadron will be set to 1
1074  // if a GSF element is associated to the current TRACK element
1075  // This information will be used in the electron rejection for tau ID.
1076  std::multimap<double, unsigned> gsfElems;
1077  block.associatedElements(iEle, linkData, gsfElems, reco::PFBlockElement::GSF);
1078 
1079  if (hcalElems.empty()) {
1080  LogTrace("PFAlgo|elementLoop") << "no hcal element connected to track " << iEle;
1081  }
1082 
1083  // will now loop on associated elements ...
1084  bool hcalFound = false;
1085  bool hfhadFound = false;
1086 
1087  LogTrace("PFAlgo|elementLoop") << "now looping on elements associated to the track: ecalElems";
1088 
1089  // ... first on associated ECAL elements
1090  // Check if there is still a free ECAL for this track
1091  for (auto const& ecal : ecalElems) {
1092  unsigned index = ecal.second;
1093  // Sanity checks and optional printout
1094  PFBlockElement::Type type = elements[index].type();
1095 #ifdef EDM_ML_DEBUG
1096  double dist = ecal.first;
1097  LogTrace("PFAlgo|elementLoop") << "\telement " << elements[index] << " linked with distance = " << dist;
1098  if (!active[index])
1099  LogTrace("PFAlgo|elementLoop") << "This ECAL is already used - skip it";
1100 #endif
1101  assert(type == PFBlockElement::ECAL);
1102 
1103  // This ECAL is not free (taken by an electron?) - just skip it
1104  if (!active[index])
1105  continue;
1106 
1107  // Flag ECAL clusters for which the corresponding track is not linked to an HCAL cluster
1108 
1109  //reco::PFClusterRef clusterRef = elements[index].clusterRef();
1110  //assert( !clusterRef.isNull() );
1111  if (!hcalElems.empty())
1112  LogTrace("PFAlgo|elementLoop") << "\t\tat least one hcal element connected to the track."
1113  << " Sparing Ecal cluster for the hcal loop";
1114 
1115  } //loop print ecal elements
1116 
1117  // tracks which are not linked to an HCAL (or HFHAD)
1118  // are reconstructed now.
1119 
1120  if (hcalElems.empty() && hfHadElems.empty()) {
1121  auto ret_continue = recoTracksNotHCAL(
1122  block, linkData, elements, blockref, active, goodTrackDeadHcal, hasDeadHcal, iEle, ecalElems, trackRef);
1123  if (ret_continue) {
1124  continue;
1125  }
1126  } // end if( hcalElems.empty() && hfHadElems.empty() )
1127 
1128  // In case several HCAL elements are linked to this track,
1129  // unlinking all of them except the closest.
1130  for (auto const& hcal : hcalElems) {
1131  unsigned index = hcal.second;
1132 
1133  PFBlockElement::Type type = elements[index].type();
1134 
1135 #ifdef EDM_ML_DEBUG
1136  double dist = block.dist(iEle, index, linkData, reco::PFBlock::LINKTEST_ALL);
1137  LogTrace("PFAlgo|elementLoop") << "\telement " << elements[index] << " linked with distance " << dist;
1138 #endif
1139  assert(type == PFBlockElement::HCAL);
1140 
1141  // all hcal clusters except the closest
1142  // will be unlinked from the track
1143  if (!hcalFound) { // closest hcal
1144  LogTrace("PFAlgo|elementLoop") << "\t\tclosest hcal cluster, doing nothing";
1145 
1146  hcalFound = true;
1147 
1148  // active[index] = false;
1149  // hcalUsed.push_back( index );
1150  } else { // other associated hcal
1151  // unlink from the track
1152  LogTrace("PFAlgo|elementLoop") << "\t\tsecondary hcal cluster. unlinking";
1153  block.setLink(iEle, index, -1., linkData, PFBlock::LINKTEST_RECHIT);
1154  }
1155  } //loop hcal elements
1156 
1157  // ---Same for HFHAD---
1158  // In case several HFHAD elements are linked to this track,
1159  // unlinking all of them except the closest.
1160  for (auto const& hfhad : hfHadElems) {
1161  unsigned index = hfhad.second;
1162 
1163  PFBlockElement::Type type = elements[index].type();
1164 
1165 #ifdef EDM_ML_DEBUG
1166  double dist = block.dist(iEle, index, linkData, reco::PFBlock::LINKTEST_ALL);
1167  LogTrace("PFAlgo|elementLoop") << "\telement " << elements[index] << " linked with distance " << dist;
1168 #endif
1169  assert(type == PFBlockElement::HFHAD);
1170 
1171  // all hfhad clusters except the closest
1172  // will be unlinked from the track
1173  if (!hfhadFound) { // closest hfhad
1174  LogTrace("PFAlgo|elementLoop") << "\t\tclosest hfhad cluster, doing nothing";
1175 
1176  hfhadFound = true;
1177 
1178  } else { // other associated hfhad
1179  // unlink from the track
1180  LogTrace("PFAlgo|elementLoop") << "\t\tsecondary hfhad cluster. unlinking";
1181  block.setLink(iEle, index, -1., linkData, PFBlock::LINKTEST_RECHIT);
1182  }
1183  } //loop hfhad elements
1184 
1185  LogTrace("PFAlgo|elementLoop") << "end of loop over iEle";
1186  } // end of outer loop on elements iEle of any type
1187  LogTrace("PFAlgo|elementLoop") << "end of function PFAlgo::elementLoop";
1188 }

References cms::cuda::assert(), groupFilesInBlocks::block, checkAndReconstructSecondaryInteraction(), checkGoodTrackDeadHcal(), checkHasDeadHcal(), decideType(), ECAL, reco::PFBlockElement::ECAL, ElementIndices::ecalIs, bookConverter::elements, reco::PFBlockElement::GSF, HCAL, reco::PFBlockElement::HCAL, ElementIndices::hcalIs, reco::PFBlockElement::HFEM, ElementIndices::hfEmIs, reco::PFBlockElement::HFHAD, ElementIndices::hfHadIs, ElementIndices::hoIs, edm::Ref< C, T, F >::isNull(), reco::PFBlock::LINKTEST_ALL, LogTrace, recoTracksNotHCAL(), relinkTrackToHcal(), and ElementIndices::trackIs.

Referenced by processBlock().

◆ getCleanedCandidates()

reco::PFCandidateCollection& PFAlgo::getCleanedCandidates ( )
inline
Returns
collection of cleaned HF candidates

Definition at line 112 of file PFAlgo.h.

112 { return pfCleanedCandidates_; }

References pfCleanedCandidates_.

Referenced by PFProducer::produce().

◆ getPFMuonAlgo()

PFMuonAlgo * PFAlgo::getPFMuonAlgo ( )

Definition at line 66 of file PFAlgo.cc.

66 { return pfmu_.get(); }

References pfmu_.

Referenced by PFProducer::produce().

◆ hfEnergyResolution()

double PFAlgo::hfEnergyResolution ( double  clusterEnergy) const
private

Definition at line 3370 of file PFAlgo.cc.

3370  {
3371  // Add a protection
3372  clusterEnergyHF = std::max(clusterEnergyHF, 1.);
3373 
3374  double resol =
3375  sqrt(resolHF_square_[0] / clusterEnergyHF + resolHF_square_[1] + resolHF_square_[2] / pow(clusterEnergyHF, 2));
3376  // 0: stochastic term, 1: constant term, 2: noise term
3377  // Note: resolHF_square_[0,1,2] should be already squared
3378 
3379  return resol;
3380 }

References SiStripPI::max, funct::pow(), l1ParticleFlow_cff::resol, resolHF_square_, and mathSSE::sqrt().

Referenced by createCandidatesHF().

◆ isFromSecInt()

bool PFAlgo::isFromSecInt ( const reco::PFBlockElement &  eTrack,
std::string  order 
) const
private

Definition at line 3462 of file PFAlgo.cc.

3462  {
3463  reco::PFBlockElement::TrackType T_TO_DISP = reco::PFBlockElement::T_TO_DISP;
3464  reco::PFBlockElement::TrackType T_FROM_DISP = reco::PFBlockElement::T_FROM_DISP;
3465  // reco::PFBlockElement::TrackType T_FROM_GAMMACONV = reco::PFBlockElement::T_FROM_GAMMACONV;
3466  reco::PFBlockElement::TrackType T_FROM_V0 = reco::PFBlockElement::T_FROM_V0;
3467 
3468  bool bPrimary = (order.find("primary") != string::npos);
3469  bool bSecondary = (order.find("secondary") != string::npos);
3470  bool bAll = (order.find("all") != string::npos);
3471 
3472  bool isToDisp = usePFNuclearInteractions_ && eTrack.trackType(T_TO_DISP);
3473  bool isFromDisp = usePFNuclearInteractions_ && eTrack.trackType(T_FROM_DISP);
3474 
3475  if (bPrimary && isToDisp)
3476  return true;
3477  if (bSecondary && isFromDisp)
3478  return true;
3479  if (bAll && (isToDisp || isFromDisp))
3480  return true;
3481 
3482  // bool isFromConv = usePFConversions_ && eTrack.trackType(T_FROM_GAMMACONV);
3483 
3484  // if ((bAll || bSecondary)&& isFromConv) return true;
3485 
3486  bool isFromDecay = (bAll || bSecondary) && usePFDecays_ && eTrack.trackType(T_FROM_V0);
3487 
3488  return isFromDecay;
3489 }

References eventshapeDQM_cfi::order, reco::PFBlockElement::T_FROM_DISP, reco::PFBlockElement::T_FROM_V0, reco::PFBlockElement::T_TO_DISP, reco::PFBlockElement::trackType(), usePFDecays_, and usePFNuclearInteractions_.

Referenced by checkAndReconstructSecondaryInteraction(), createCandidatesHCAL(), reconstructTrack(), and recoTracksNotHCAL().

◆ makeConnectedCandidates()

reco::PFCandidateCollection PFAlgo::makeConnectedCandidates ( )
inline
Returns
the collection of candidates

Definition at line 115 of file PFAlgo.h.

115 { return connector_.connect(*pfCandidates_); }

References PFCandConnector::connect(), connector_, and pfCandidates_.

Referenced by PFProducer::produce().

◆ neutralHadronEnergyResolution()

double PFAlgo::neutralHadronEnergyResolution ( double  clusterEnergy,
double  clusterEta 
) const
private

todo: use PFClusterTools for this

Definition at line 3353 of file PFAlgo.cc.

3353  {
3354  // Add a protection
3355  clusterEnergyHCAL = std::max(clusterEnergyHCAL, 1.);
3356 
3357  double resol = fabs(eta) < 1.48 ? sqrt(1.02 * 1.02 / clusterEnergyHCAL + 0.065 * 0.065)
3358  : sqrt(1.20 * 1.20 / clusterEnergyHCAL + 0.028 * 0.028);
3359 
3360  return resol;
3361 }

References PVValHelper::eta, SiStripPI::max, l1ParticleFlow_cff::resol, and mathSSE::sqrt().

Referenced by createCandidatesHCAL(), and recoTracksNotHCAL().

◆ nSigmaHCAL()

double PFAlgo::nSigmaHCAL ( double  clusterEnergy,
double  clusterEta 
) const
private

Definition at line 3363 of file PFAlgo.cc.

3363  {
3364  double nS = fabs(eta) < 1.48 ? nSigmaHCAL_ * (1. + exp(-clusterEnergyHCAL / nSigmaEConstHCAL))
3365  : nSigmaHCAL_ * (1. + exp(-clusterEnergyHCAL / nSigmaEConstHCAL));
3366 
3367  return nS;
3368 }

References PVValHelper::eta, JetChargeProducer_cfi::exp, nSigmaEConstHCAL, and nSigmaHCAL_.

Referenced by createCandidatesHCAL().

◆ nSigmaHFEM()

double PFAlgo::nSigmaHFEM ( double  clusterEnergy) const
private

Definition at line 3382 of file PFAlgo.cc.

3382  {
3383  double nS = nSigmaHFEM_ * (1. + exp(-clusterEnergyHF / nSigmaEConstHFEM));
3384  return nS;
3385 }

References JetChargeProducer_cfi::exp, nSigmaEConstHFEM, and nSigmaHFEM_.

Referenced by createCandidatesHF().

◆ nSigmaHFHAD()

double PFAlgo::nSigmaHFHAD ( double  clusterEnergy) const
private

Definition at line 3387 of file PFAlgo.cc.

3387  {
3388  double nS = nSigmaHFHAD_ * (1. + exp(-clusterEnergyHF / nSigmaEConstHFHAD));
3389  return nS;
3390 }

References JetChargeProducer_cfi::exp, nSigmaEConstHFHAD, and nSigmaHFHAD_.

Referenced by createCandidatesHF().

◆ postCleaning()

void PFAlgo::postCleaning ( )
private

Definition at line 3491 of file PFAlgo.cc.

3491  {
3492  //Compute met and met significance (met/sqrt(SumEt))
3493  double metX = 0.;
3494  double metY = 0.;
3495  double sumet = 0;
3496  std::vector<unsigned int> pfCandidatesToBeRemoved;
3497  for (auto const& pfc : *pfCandidates_) {
3498  metX += pfc.px();
3499  metY += pfc.py();
3500  sumet += pfc.pt();
3501  }
3502  double met2 = metX * metX + metY * metY;
3503  // Select events with large MET significance.
3504  double significance = std::sqrt(met2 / sumet);
3505  double significanceCor = significance;
3506  if (significance > minSignificance_) {
3507  double metXCor = metX;
3508  double metYCor = metY;
3509  double sumetCor = sumet;
3510  double met2Cor = met2;
3511  double deltaPhi = 3.14159;
3512  double deltaPhiPt = 100.;
3513  bool next = true;
3514  unsigned iCor = 1E9;
3515 
3516  // Find the HF candidate with the largest effect on the MET
3517  while (next) {
3518  double metReduc = -1.;
3519  // Loop on the candidates
3520  for (unsigned i = 0; i < pfCandidates_->size(); ++i) {
3521  const PFCandidate& pfc = (*pfCandidates_)[i];
3522 
3523  // Check that the pfCandidate is in the HF
3524  if (pfc.particleId() != reco::PFCandidate::h_HF && pfc.particleId() != reco::PFCandidate::egamma_HF)
3525  continue;
3526 
3527  // Check if has meaningful pt
3528  if (pfc.pt() < minHFCleaningPt_)
3529  continue;
3530 
3531  // Check that it is not already scheculed to be cleaned
3532  const bool skip = std::any_of(
3533  pfCandidatesToBeRemoved.begin(), pfCandidatesToBeRemoved.end(), [&](unsigned int j) { return i == j; });
3534  if (skip)
3535  continue;
3536 
3537  // Check that the pt and the MET are aligned
3538  deltaPhi = std::acos((metX * pfc.px() + metY * pfc.py()) / (pfc.pt() * std::sqrt(met2)));
3539  deltaPhiPt = deltaPhi * pfc.pt();
3540  if (deltaPhiPt > maxDeltaPhiPt_)
3541  continue;
3542 
3543  // Now remove the candidate from the MET
3544  double metXInt = metX - pfc.px();
3545  double metYInt = metY - pfc.py();
3546  double sumetInt = sumet - pfc.pt();
3547  double met2Int = metXInt * metXInt + metYInt * metYInt;
3548  if (met2Int < met2Cor) {
3549  metXCor = metXInt;
3550  metYCor = metYInt;
3551  metReduc = (met2 - met2Int) / met2Int;
3552  met2Cor = met2Int;
3553  sumetCor = sumetInt;
3554  significanceCor = std::sqrt(met2Cor / sumetCor);
3555  iCor = i;
3556  }
3557  }
3558  //
3559  // If the MET must be significanly reduced, schedule the candidate to be cleaned
3560  if (metReduc > minDeltaMet_) {
3561  pfCandidatesToBeRemoved.push_back(iCor);
3562  metX = metXCor;
3563  metY = metYCor;
3564  sumet = sumetCor;
3565  met2 = met2Cor;
3566  } else {
3567  // Otherwise just stop the loop
3568  next = false;
3569  }
3570  }
3571  //
3572  // The significance must be significantly reduced to indeed clean the candidates
3573  if (significance - significanceCor > minSignificanceReduction_ && significanceCor < maxSignificance_) {
3574  edm::LogInfo("PFAlgo|postCleaning") << "Significance reduction = " << significance << " -> " << significanceCor
3575  << " = " << significanceCor - significance;
3576  for (unsigned int toRemove : pfCandidatesToBeRemoved) {
3577  edm::LogInfo("PFAlgo|postCleaning") << "Removed : " << (*pfCandidates_)[toRemove];
3578  pfCleanedCandidates_.push_back((*pfCandidates_)[toRemove]);
3579  (*pfCandidates_)[toRemove].rescaleMomentum(1E-6);
3580  //reco::PFCandidate::ParticleType unknown = reco::PFCandidate::X;
3581  //(*pfCandidates_)[toRemove].setParticleType(unknown);
3582  }
3583  }
3584  } //significance
3585 } //postCleaning

References SiPixelRawToDigiRegional_cfi::deltaPhi, reco::PFCandidate::egamma_HF, reco::PFCandidate::h_HF, mps_fire::i, dqmiolumiharvest::j, maxDeltaPhiPt_, maxSignificance_, minDeltaMet_, minHFCleaningPt_, minSignificance_, minSignificanceReduction_, GetRecoTauVFromDQM_MC_cff::next, reco::PFCandidate::particleId(), pfCandidates_, pfCleanedCandidates_, reco::LeafCandidate::pt(), reco::LeafCandidate::px(), reco::LeafCandidate::py(), met_cff::significance, optionsL1T::skip, and mathSSE::sqrt().

Referenced by reconstructParticles().

◆ processBlock()

void PFAlgo::processBlock ( const reco::PFBlockRef &  blockref,
std::list< reco::PFBlockRef > &  hcalBlockRefs,
std::list< reco::PFBlockRef > &  ecalBlockRefs,
PFEGammaFilters const *  pfegamma 
)
private

process one block. can be reimplemented in more sophisticated algorithms

Definition at line 3053 of file PFAlgo.cc.

3056  {
3057  assert(!blockref.isNull());
3058  const reco::PFBlock& block = *blockref;
3059 
3060  LogTrace("PFAlgo|processBlock") << "start of function PFAlgo::processBlock, block=" << block;
3061 
3062  const edm::OwnVector<reco::PFBlockElement>& elements = block.elements();
3063  LogTrace("PFAlgo|processBlock") << "elements.size()=" << elements.size();
3064  // make a copy of the link data, which will be edited.
3065  PFBlock::LinkData linkData = block.linkData();
3066 
3067  // keep track of the elements which are still active.
3068  vector<bool> active(elements.size(), true);
3069 
3070  // //PFElectrons:
3071  // usePFElectrons_ external configurable parameter to set the usage of pf electron
3072  std::vector<reco::PFCandidate> tempElectronCandidates;
3073  tempElectronCandidates.clear();
3074 
3075  // New EGamma Reconstruction 10/10/2013
3076  if (useEGammaFilters_) {
3077  egammaFilters(blockref, active, pfegamma);
3078  } // end if use EGammaFilters
3079 
3080  //Lock extra conversion tracks not used by Photon Algo
3081  if (usePFConversions_) {
3082  conversionAlgo(elements, active);
3083  }
3084 
3085  // In the following elementLoop() function, the primary goal is to deal with tracks that are:
3086  // - not associated to an HCAL cluster
3087  // - not identified as an electron.
3088  // Those tracks should be predominantly relatively low energy charged
3089  // hadrons which are not detected in the ECAL.
3090 
3091  // The secondary goal is to prepare for the next loops
3092  // - The ecal and hcal elements are sorted in separate vectors in `ElementIndices inds`
3093  // which will be used as a base for the corresponding loops.
3094  // - For tracks which are connected to more than one HCAL cluster,
3095  // the links between the track and the cluster are cut for all clusters
3096  // but the closest one.
3097  // - HF only blocks ( HFEM, HFHAD, HFEM+HFAD) are identified
3098 
3099  // obsolete comments?
3100  // loop1:
3101  // - sort ecal and hcal elements in separate vectors
3102  // - for tracks:
3103  // - lock closest ecal cluster
3104  // - cut link to farthest hcal cluster, if more than 1.
3105 
3106  vector<bool> deadArea(elements.size(), false);
3107 
3108  // vectors to store element indices to ho, hcal and ecal elements, will be filled by elementLoop()
3109  ElementIndices inds;
3110 
3111  elementLoop(block, linkData, elements, active, blockref, inds, deadArea);
3112 
3113  // Reconstruct pfCandidate from HF (either EM-only, Had-only or both)
3114  // For phase2, process also pfblocks containing HF clusters and linked tracks
3115  if (!(inds.hfEmIs.empty() && inds.hfHadIs.empty())) {
3116  createCandidatesHF(block, linkData, elements, active, blockref, inds);
3117  if (inds.hcalIs.empty() && inds.ecalIs.empty())
3118  return;
3119  LogDebug("PFAlgo::processBlock")
3120  << "Block contains HF clusters, and also contains ECAL or HCAL clusters. Continue.\n"
3121  << block;
3122  }
3123 
3124  createCandidatesHCAL(block, linkData, elements, active, blockref, inds, deadArea);
3125  // COLINFEB16: now dealing with the HCAL elements that are not linked to any track
3126  createCandidatesHCALUnlinked(block, linkData, elements, active, blockref, inds, deadArea);
3127  createCandidatesECAL(block, linkData, elements, active, blockref, inds, deadArea);
3128 
3129  LogTrace("PFAlgo|processBlock") << "end of function PFAlgo::processBlock";
3130 } // end processBlock

References cms::cuda::assert(), groupFilesInBlocks::block, conversionAlgo(), createCandidatesECAL(), createCandidatesHCAL(), createCandidatesHCALUnlinked(), createCandidatesHF(), egammaFilters(), elementLoop(), bookConverter::elements, edm::Ref< C, T, F >::isNull(), LogDebug, LogTrace, useEGammaFilters_, and usePFConversions_.

Referenced by reconstructParticles().

◆ reconstructCluster()

unsigned PFAlgo::reconstructCluster ( const reco::PFCluster &  cluster,
double  particleEnergy,
bool  useDirection = false,
double  particleX = 0.,
double  particleY = 0.,
double  particleZ = 0. 
)
private

Reconstruct a neutral particle from a cluster. If chargedEnergy is specified, the neutral particle is created only if the cluster energy is significantly larger than the chargedEnergy. In this case, the energy of the neutral particle is cluster energy - chargedEnergy

Definition at line 3212 of file PFAlgo.cc.

3217  {
3218  LogTrace("PFAlgo|reconstructCluster") << "start of function PFAlgo::reconstructCluster, cluster=" << cluster
3219  << "particleEnergy=" << particleEnergy << "useDirection=" << useDirection
3220  << "particleX=" << particleX << "particleY=" << particleY
3221  << "particleZ=" << particleZ;
3222 
3223  reco::PFCandidate::ParticleType particleType = reco::PFCandidate::X;
3224 
3225  // need to convert the ::math::XYZPoint data member of the PFCluster class=
3226  // to a displacement vector:
3227 
3228  // Transform particleX,Y,Z to a position at ECAL/HCAL entrance
3229  double factor = 1.;
3230  if (useDirection) {
3231  switch (cluster.layer()) {
3232  case PFLayer::ECAL_BARREL:
3233  case PFLayer::HCAL_BARREL1:
3234  factor = std::sqrt(cluster.position().Perp2() / (particleX * particleX + particleY * particleY));
3235  break;
3236  case PFLayer::ECAL_ENDCAP:
3237  case PFLayer::HCAL_ENDCAP:
3238  case PFLayer::HF_HAD:
3239  case PFLayer::HF_EM:
3240  factor = cluster.position().Z() / particleZ;
3241  break;
3242  default:
3243  assert(0);
3244  }
3245  }
3246  //MIKE First of all let's check if we have vertex.
3247  ::math::XYZPoint vertexPos;
3248  if (useVertices_)
3249  vertexPos = ::math::XYZPoint(primaryVertex_.x(), primaryVertex_.y(), primaryVertex_.z());
3250  else
3251  vertexPos = ::math::XYZPoint(0.0, 0.0, 0.0);
3252 
3253  ::math::XYZVector clusterPos(cluster.position().X() - vertexPos.X(),
3254  cluster.position().Y() - vertexPos.Y(),
3255  cluster.position().Z() - vertexPos.Z());
3256  ::math::XYZVector particleDirection(
3257  particleX * factor - vertexPos.X(), particleY * factor - vertexPos.Y(), particleZ * factor - vertexPos.Z());
3258 
3259  //::math::XYZVector clusterPos( cluster.position().X(), cluster.position().Y(),cluster.position().Z() );
3260  //::math::XYZVector particleDirection ( particleX, particleY, particleZ );
3261 
3262  clusterPos = useDirection ? particleDirection.Unit() : clusterPos.Unit();
3263  clusterPos *= particleEnergy;
3264 
3265  // clusterPos is now a vector along the cluster direction,
3266  // with a magnitude equal to the cluster energy.
3267 
3268  double mass = 0;
3269  ROOT::Math::LorentzVector<ROOT::Math::PxPyPzM4D<double>> momentum(
3270  clusterPos.X(), clusterPos.Y(), clusterPos.Z(), mass);
3271  // mathcore is a piece of #$%
3272  ::math::XYZTLorentzVector tmp;
3273  // implicit constructor not allowed
3274  tmp = momentum;
3275 
3276  // Charge
3277  int charge = 0;
3278 
3279  // Type
3280  switch (cluster.layer()) {
3281  case PFLayer::ECAL_BARREL:
3282  case PFLayer::ECAL_ENDCAP:
3283  particleType = PFCandidate::gamma;
3284  break;
3285  case PFLayer::HCAL_BARREL1:
3286  case PFLayer::HCAL_ENDCAP:
3287  particleType = PFCandidate::h0;
3288  break;
3289  case PFLayer::HF_HAD:
3290  particleType = PFCandidate::h_HF;
3291  break;
3292  case PFLayer::HF_EM:
3293  particleType = PFCandidate::egamma_HF;
3294  break;
3295  default:
3296  assert(0);
3297  }
3298 
3299  // The pf candidate
3300  LogTrace("PFAlgo|reconstructCluster") << "Creating PFCandidate charge=" << charge << ", type=" << particleType
3301  << ", pt=" << tmp.pt() << ", eta=" << tmp.eta() << ", phi=" << tmp.phi();
3302  pfCandidates_->push_back(PFCandidate(charge, tmp, particleType));
3303 
3304  // The position at ECAL entrance (well: watch out, it is not true
3305  // for HCAL clusters... to be fixed)
3306  pfCandidates_->back().setPositionAtECALEntrance(
3307  ::math::XYZPointF(cluster.position().X(), cluster.position().Y(), cluster.position().Z()));
3308 
3309  //Set the cnadidate Vertex
3310  pfCandidates_->back().setVertex(vertexPos);
3311 
3312  // depth info
3313  setHcalDepthInfo(pfCandidates_->back(), cluster);
3314 
3315  //*TODO* cluster time is not reliable at the moment, so only use track timing
3316 
3317  LogTrace("PFAlgo|reconstructCluster") << "** candidate: " << pfCandidates_->back();
3318 
3319  // returns index to the newly created PFCandidate
3320  return pfCandidates_->size() - 1;
3321 }

References cms::cuda::assert(), ALCARECOTkAlJpsiMuMu_cff::charge, PFLayer::ECAL_BARREL, PFLayer::ECAL_ENDCAP, DQMScaleToClient_cfi::factor, CustomPhysics_cfi::gamma, PFLayer::HCAL_BARREL1, PFLayer::HCAL_ENDCAP, PFLayer::HF_EM, PFLayer::HF_HAD, reco::PFCluster::layer(), LogTrace, EgHLTOffHistBins_cfi::mass, PbPb_ZMuSkimMuonDPG_cff::particleType, pfElectronTranslator_cfi::PFCandidate, pfCandidates_, reco::CaloCluster::position(), primaryVertex_, setHcalDepthInfo(), mathSSE::sqrt(), createJobs::tmp, useVertices_, reco::PFCandidate::X, reco::Vertex::x(), reco::Vertex::y(), and reco::Vertex::z().

Referenced by checkCleaning(), createCandidatesECAL(), createCandidatesHCAL(), createCandidatesHCALUnlinked(), createCandidatesHF(), and recoTracksNotHCAL().

◆ reconstructParticles()

void PFAlgo::reconstructParticles ( const reco::PFBlockHandle &  blockHandle,
PFEGammaFilters const *  pfegamma 
)

reconstruct particles

Definition at line 130 of file PFAlgo.cc.

130  {
131  auto const& blocks = *blockHandle;
132 
133  // reset output collection
134  pfCandidates_->clear();
135 
136  LogTrace("PFAlgo|reconstructParticles")
137  << "start of function PFAlgo::reconstructParticles, blocks.size()=" << blocks.size();
138 
139  // sort elements in three lists:
140  std::list<reco::PFBlockRef> hcalBlockRefs;
141  std::list<reco::PFBlockRef> ecalBlockRefs;
142  std::list<reco::PFBlockRef> hoBlockRefs;
143  std::list<reco::PFBlockRef> otherBlockRefs;
144 
145  for (unsigned i = 0; i < blocks.size(); ++i) {
146  reco::PFBlockRef blockref = reco::PFBlockRef(blockHandle, i);
147 
148  const reco::PFBlock& block = *blockref;
149  const edm::OwnVector<reco::PFBlockElement>& elements = block.elements();
150 
151  bool singleEcalOrHcal = false;
152  if (elements.size() == 1) {
153  if (elements[0].type() == reco::PFBlockElement::ECAL) {
154  ecalBlockRefs.push_back(blockref);
155  singleEcalOrHcal = true;
156  }
157  if (elements[0].type() == reco::PFBlockElement::HCAL) {
158  hcalBlockRefs.push_back(blockref);
159  singleEcalOrHcal = true;
160  }
161  if (elements[0].type() == reco::PFBlockElement::HO) {
162  // Single HO elements are likely to be noise. Not considered for now.
163  hoBlockRefs.push_back(blockref);
164  singleEcalOrHcal = true;
165  }
166  }
167 
168  if (!singleEcalOrHcal) {
169  otherBlockRefs.push_back(blockref);
170  }
171  } //loop blocks
172 
173  LogTrace("PFAlgo|reconstructParticles")
174  << "# Ecal blocks: " << ecalBlockRefs.size() << ", # Hcal blocks: " << hcalBlockRefs.size()
175  << ", # HO blocks: " << hoBlockRefs.size() << ", # Other blocks: " << otherBlockRefs.size();
176 
177  // loop on blocks that are not single ecal,
178  // and not single hcal.
179 
180  unsigned nblcks = 0;
181  for (auto const& other : otherBlockRefs) {
182  LogTrace("PFAlgo|reconstructParticles") << "processBlock, Block number " << nblcks++;
183  processBlock(other, hcalBlockRefs, ecalBlockRefs, pfegamma);
184  }
185 
186  std::list<reco::PFBlockRef> empty;
187 
188  unsigned hblcks = 0;
189  // process remaining single hcal blocks
190  for (auto const& hcal : hcalBlockRefs) {
191  LogTrace("PFAlgo|reconstructParticles") << "processBlock, HCAL block number " << hblcks++;
192  processBlock(hcal, empty, empty, pfegamma);
193  }
194 
195  unsigned eblcks = 0;
196  // process remaining single ecal blocks
197  for (auto const& ecal : ecalBlockRefs) {
198  LogTrace("PFAlgo|reconstructParticles") << "processBlock, ECAL block number " << eblcks++;
199  processBlock(ecal, empty, empty, pfegamma);
200  }
201 
202  // Post HF Cleaning
203  pfCleanedCandidates_.clear();
204  // Check if the post HF Cleaning was requested - if not, do nothing
205  if (postHFCleaning_) {
206  postCleaning();
207  }
208 
209  //Muon post cleaning
210  pfmu_->postClean(pfCandidates_.get());
211 
212  //Add Missing muons
213  if (muonHandle_.isValid())
214  pfmu_->addMissingMuons(muonHandle_, pfCandidates_.get());
215 
216  LogTrace("PFAlgo|reconstructParticles")
217  << "end of function PFAlgo::reconstructParticles, pfCandidates_->size()=" << pfCandidates_->size();
218 }

References groupFilesInBlocks::block, gather_cfg::blocks, reco::PFBlockElement::ECAL, bookConverter::elements, relativeConstraints::empty, reco::PFBlockElement::HCAL, reco::PFBlockElement::HO, mps_fire::i, edm::HandleBase::isValid(), LogTrace, muonHandle_, trackingPlots::other, pfCandidates_, pfCleanedCandidates_, pfmu_, postCleaning(), postHFCleaning_, and processBlock().

Referenced by PFProducer::produce().

◆ reconstructTrack()

unsigned PFAlgo::reconstructTrack ( const reco::PFBlockElement &  elt,
bool  allowLoose = false 
)
private

Reconstruct a charged particle from a track Returns the index of the newly created candidate in pfCandidates_ Michalis added a flag here to treat muons inside jets

Definition at line 3133 of file PFAlgo.cc.

3133  {
3134  const auto* eltTrack = dynamic_cast<const reco::PFBlockElementTrack*>(&elt);
3135 
3136  const reco::TrackRef& trackRef = eltTrack->trackRef();
3137  const reco::Track& track = *trackRef;
3138  const reco::MuonRef& muonRef = eltTrack->muonRef();
3139  int charge = track.charge() > 0 ? 1 : -1;
3140 
3141  // Assume this particle is a charged Hadron
3142  double px = track.px();
3143  double py = track.py();
3144  double pz = track.pz();
3145  double energy = sqrt(track.p() * track.p() + 0.13957 * 0.13957);
3146 
3147  LogTrace("PFAlgo|reconstructTrack") << "Reconstructing PF candidate from track of pT = " << track.pt()
3148  << " eta = " << track.eta() << " phi = " << track.phi() << " px = " << px
3149  << " py = " << py << " pz = " << pz << " energy = " << energy;
3150 
3151  // Create a PF Candidate
3152  ::math::XYZTLorentzVector momentum(px, py, pz, energy);
3153  reco::PFCandidate::ParticleType particleType = reco::PFCandidate::h;
3154 
3155  // Add it to the stack
3156  LogTrace("PFAlgo|reconstructTrack") << "Creating PFCandidate charge=" << charge << ", type=" << particleType
3157  << ", pt=" << momentum.pt() << ", eta=" << momentum.eta()
3158  << ", phi=" << momentum.phi();
3159  pfCandidates_->push_back(PFCandidate(charge, momentum, particleType));
3160  //Set vertex and stuff like this
3161  pfCandidates_->back().setVertex(trackRef->vertex());
3162  pfCandidates_->back().setTrackRef(trackRef);
3163  pfCandidates_->back().setPositionAtECALEntrance(eltTrack->positionAtECALEntrance());
3164  if (muonRef.isNonnull())
3165  pfCandidates_->back().setMuonRef(muonRef);
3166 
3167  //Set time
3168  if (elt.isTimeValid())
3169  pfCandidates_->back().setTime(elt.time(), elt.timeError());
3170 
3171  //OK Now try to reconstruct the particle as a muon
3172  bool isMuon = pfmu_->reconstructMuon(pfCandidates_->back(), muonRef, allowLoose);
3173  bool isFromDisp = isFromSecInt(elt, "secondary");
3174 
3175  if ((!isMuon) && isFromDisp) {
3176  double dpt = trackRef->ptError();
3177  double dptRel = dpt / trackRef->pt() * 100;
3178  //If the track is ill measured it is better to not refit it, since the track information probably would not be used.
3179  //In the PFAlgo we use the trackref information. If the track error is too big the refitted information might be very different
3180  // from the not refitted one.
3181  if (dptRel < dptRel_DispVtx_) {
3182  LogTrace("PFAlgo|reconstructTrack")
3183  << "Not refitted px = " << px << " py = " << py << " pz = " << pz << " energy = " << energy;
3184  //reco::TrackRef trackRef = eltTrack->trackRef();
3185  reco::PFDisplacedVertexRef vRef =
3186  eltTrack->displacedVertexRef(reco::PFBlockElement::T_FROM_DISP)->displacedVertexRef();
3187  reco::Track trackRefit = vRef->refittedTrack(trackRef);
3188  //change the momentum with the refitted track
3189  ::math::XYZTLorentzVector momentum(
3190  trackRefit.px(), trackRefit.py(), trackRefit.pz(), sqrt(trackRefit.p() * trackRefit.p() + 0.13957 * 0.13957));
3191  LogTrace("PFAlgo|reconstructTrack")
3192  << "Refitted px = " << px << " py = " << py << " pz = " << pz << " energy = " << energy;
3193  }
3194  pfCandidates_->back().setFlag(reco::PFCandidate::T_FROM_DISP, true);
3195  pfCandidates_->back().setDisplacedVertexRef(
3196  eltTrack->displacedVertexRef(reco::PFBlockElement::T_FROM_DISP)->displacedVertexRef(),
3197  reco::PFCandidate::T_FROM_DISP);
3198  }
3199 
3200  // do not label as primary a track which would be recognised as a muon. A muon cannot produce NI. It is with high probability a fake
3201  if (isFromSecInt(elt, "primary") && !isMuon) {
3202  pfCandidates_->back().setFlag(reco::PFCandidate::T_TO_DISP, true);
3203  pfCandidates_->back().setDisplacedVertexRef(
3204  eltTrack->displacedVertexRef(reco::PFBlockElement::T_TO_DISP)->displacedVertexRef(),
3205  reco::PFCandidate::T_TO_DISP);
3206  }
3207 
3208  // returns index to the newly created PFCandidate
3209  return pfCandidates_->size() - 1;
3210 }

References ALCARECOTkAlJpsiMuMu_cff::charge, dptRel_DispVtx_, HCALHighEnergyHPDFilter_cfi::energy, reco::PFCandidate::h, isFromSecInt(), reco::isMuon(), edm::Ref< C, T, F >::isNonnull(), reco::PFBlockElement::isTimeValid(), LogTrace, reco::TrackBase::p(), PbPb_ZMuSkimMuonDPG_cff::particleType, pfElectronTranslator_cfi::PFCandidate, pfCandidates_, pfmu_, multPhiCorr_741_25nsDY_cfi::px, reco::TrackBase::px(), multPhiCorr_741_25nsDY_cfi::py, reco::TrackBase::py(), reco::TrackBase::pz(), mathSSE::sqrt(), reco::PFBlockElement::T_FROM_DISP, reco::PFCandidate::T_FROM_DISP, reco::PFBlockElement::T_TO_DISP, reco::PFCandidate::T_TO_DISP, reco::PFBlockElement::time(), reco::PFBlockElement::timeError(), and HLT_FULL_cff::track.

Referenced by checkAndReconstructSecondaryInteraction(), createCandidatesHCAL(), createCandidatesHF(), and recoTracksNotHCAL().

◆ recoTracksNotHCAL()

bool PFAlgo::recoTracksNotHCAL ( const reco::PFBlock &  block,
reco::PFBlock::LinkData &  linkData,
const edm::OwnVector< reco::PFBlockElement > &  elements,
const reco::PFBlockRef &  blockref,
std::vector< bool > &  active,
bool  goodTrackDeadHcal,
bool  hasDeadHcal,
unsigned int  iTrack,
std::multimap< double, unsigned > &  ecalElems,
reco::TrackRef &  trackRef 
)
private

Definition at line 391 of file PFAlgo.cc.

400  {
401  LogTrace("PFAlgo|recoTracksNotHCAL") << "start of function PFAlgo::recoTracksNotHCAL(), now dealing with tracks "
402  "linked to no HCAL clusters. Was HCal active? "
403  << (!hasDeadHcal);
404  // vector<unsigned> elementIndices;
405  // elementIndices.push_back(iTrack);
406 
407  // The track momentum
408  double trackMomentum = elements[iTrack].trackRef()->p();
409  LogTrace("PFAlgo|recoTracksNotHCAL") << elements[iTrack];
410 
411  // Is it a "tight" muon ? In that case assume no
412  //Track momentum corresponds to the calorimeter
413  //energy for now
414  bool thisIsAMuon = PFMuonAlgo::isMuon(elements[iTrack]);
415  if (thisIsAMuon)
416  trackMomentum = 0.;
417 
418  // If it is not a muon check if Is it a fake track ?
419  //Michalis: I wonder if we should convert this to dpt/pt?
420  if (!thisIsAMuon && elements[iTrack].trackRef()->ptError() > ptError_) {
421  double deficit = trackMomentum;
422  double resol = neutralHadronEnergyResolution(trackMomentum, elements[iTrack].trackRef()->eta());
423  resol *= trackMomentum;
424 
425  if (!ecalElems.empty()) {
426  unsigned thisEcal = ecalElems.begin()->second;
427  reco::PFClusterRef clusterRef = elements[thisEcal].clusterRef();
428  bool useCluster = true;
429  if (hasDeadHcal) {
430  std::multimap<double, unsigned> sortedTracks;
431  block.associatedElements(
432  thisEcal, linkData, sortedTracks, reco::PFBlockElement::TRACK, reco::PFBlock::LINKTEST_ALL);
433  useCluster = (sortedTracks.begin()->second == iTrack);
434  }
435  if (useCluster) {
436  deficit -= clusterRef->energy();
437  resol = neutralHadronEnergyResolution(trackMomentum, clusterRef->positionREP().Eta());
438  resol *= trackMomentum;
439  }
440  }
441 
442  bool isPrimary = isFromSecInt(elements[iTrack], "primary");
443 
444  if (deficit > nSigmaTRACK_ * resol && !isPrimary && !goodTrackDeadHcal) {
445  active[iTrack] = false;
446  LogTrace("PFAlgo|recoTracksNotHCAL")
447  << elements[iTrack] << std::endl
448  << " deficit " << deficit << " > " << nSigmaTRACK_ << " x " << resol << " track pt " << trackRef->pt()
449  << " +- " << trackRef->ptError() << " layers valid " << trackRef->hitPattern().trackerLayersWithMeasurement()
450  << ", lost " << trackRef->hitPattern().trackerLayersWithoutMeasurement(HitPattern::TRACK_HITS)
451  << ", lost outer " << trackRef->hitPattern().trackerLayersWithoutMeasurement(HitPattern::MISSING_OUTER_HITS)
452  << ", lost inner " << trackRef->hitPattern().trackerLayersWithoutMeasurement(HitPattern::MISSING_INNER_HITS)
453  << "(valid fraction " << trackRef->validFraction() << ")"
454  << " chi2/ndf " << trackRef->normalizedChi2() << " |dxy| "
455  << std::abs(trackRef->dxy(primaryVertex_.position())) << " +- " << trackRef->dxyError() << " |dz| "
456  << std::abs(trackRef->dz(primaryVertex_.position())) << " +- " << trackRef->dzError()
457  << "is probably a fake (1) --> lock the track";
458  return true;
459  }
460  } // !thisIsaMuon
461 
462  // Create a track candidate
463  // unsigned tmpi = reconstructTrack( elements[iTrack] );
464  //active[iTrack] = false;
465  std::vector<unsigned> tmpi;
466  std::vector<unsigned> kTrack;
467 
468  // Some cleaning : secondary tracks without calo's and large momentum must be fake
469  double dpt = trackRef->ptError();
470 
471  if (rejectTracks_Step45_ && ecalElems.empty() && trackMomentum > 30. && dpt > 0.5 &&
472  (PFTrackAlgoTools::step45(trackRef->algo()) && !goodTrackDeadHcal)) {
473  double dptRel = dpt / trackRef->pt() * 100;
474  bool isPrimaryOrSecondary = isFromSecInt(elements[iTrack], "all");
475 
476  if (isPrimaryOrSecondary && dptRel < dptRel_DispVtx_) {
477  return true;
478  };
479  unsigned nHits = elements[iTrack].trackRef()->hitPattern().trackerLayersWithMeasurement();
480  unsigned int NLostHit = trackRef->hitPattern().trackerLayersWithoutMeasurement(HitPattern::TRACK_HITS);
481 
482  LogTrace("PFAlgo|recoTracksNotHCAL") << "A track (algo = " << trackRef->algo() << ") with momentum "
483  << trackMomentum << " / " << elements[iTrack].trackRef()->pt() << " +/- "
484  << dpt << " / " << elements[iTrack].trackRef()->eta()
485  << " without any link to ECAL/HCAL and with " << nHits << " (" << NLostHit
486  << ") hits (lost hits) has been cleaned";
487 
488  active[iTrack] = false;
489  return true;
490  } //rejectTracks_Step45_ && ...
491 
492  tmpi.push_back(reconstructTrack(elements[iTrack]));
493 
494  kTrack.push_back(iTrack);
495  active[iTrack] = false;
496 
497  // No ECAL cluster either ... continue...
498  if (ecalElems.empty()) {
499  (*pfCandidates_)[tmpi[0]].setEcalEnergy(0., 0.);
500  (*pfCandidates_)[tmpi[0]].setHcalEnergy(0., 0.);
501  (*pfCandidates_)[tmpi[0]].setHoEnergy(0., 0.);
502  (*pfCandidates_)[tmpi[0]].setPs1Energy(0);
503  (*pfCandidates_)[tmpi[0]].setPs2Energy(0);
504  (*pfCandidates_)[tmpi[0]].addElementInBlock(blockref, kTrack[0]);
505  return true;
506  }
507 
508  // Look for closest ECAL cluster
509  const unsigned int thisEcal = ecalElems.begin()->second;
510  reco::PFClusterRef clusterRef = elements[thisEcal].clusterRef();
511  LogTrace("PFAlgo|recoTracksNotHCAL") << " is associated to " << elements[thisEcal];
512 
513  // Set ECAL energy for muons
514  if (thisIsAMuon) {
515  (*pfCandidates_)[tmpi[0]].setEcalEnergy(clusterRef->energy(), std::min(clusterRef->energy(), muonECAL_[0]));
516  (*pfCandidates_)[tmpi[0]].setHcalEnergy(0., 0.);
517  (*pfCandidates_)[tmpi[0]].setHoEnergy(0., 0.);
518  (*pfCandidates_)[tmpi[0]].setPs1Energy(0);
519  (*pfCandidates_)[tmpi[0]].setPs2Energy(0);
520  (*pfCandidates_)[tmpi[0]].addElementInBlock(blockref, kTrack[0]);
521  }
522 
523  double slopeEcal = 1.;
524  bool connectedToEcal = false;
525  unsigned iEcal = 99999;
526  double calibEcal = 0.;
527  double calibHcal = 0.;
528  double totalEcal = thisIsAMuon ? -muonECAL_[0] : 0.;
529 
530  // Consider charged particles closest to the same ECAL cluster
531  std::multimap<double, unsigned> sortedTracks;
532  block.associatedElements(thisEcal, linkData, sortedTracks, reco::PFBlockElement::TRACK, reco::PFBlock::LINKTEST_ALL);
533  LogTrace("PFAlgo|recoTracksNotHCAL") << "the closest ECAL cluster, id " << thisEcal << ", has " << sortedTracks.size()
534  << " associated tracks, now processing them. ";
535 
536  if (hasDeadHcal && !sortedTracks.empty()) {
537  // GP: only allow the ecal cluster closest to this track to be considered
538  LogTrace("PFAlgo|recoTracksNotHCAL") << " the closest track to ECAL " << thisEcal << " is "
539  << sortedTracks.begin()->second << " (distance " << sortedTracks.begin()->first
540  << ")";
541  if (sortedTracks.begin()->second != iTrack) {
542  LogTrace("PFAlgo|recoTracksNotHCAL")
543  << " the closest track to ECAL " << thisEcal << " is " << sortedTracks.begin()->second
544  << " which is not the one being processed. Will skip ECAL linking for this track";
545  (*pfCandidates_)[tmpi[0]].setEcalEnergy(0., 0.);
546  (*pfCandidates_)[tmpi[0]].setHcalEnergy(0., 0.);
547  (*pfCandidates_)[tmpi[0]].setHoEnergy(0., 0.);
548  (*pfCandidates_)[tmpi[0]].setPs1Energy(0);
549  (*pfCandidates_)[tmpi[0]].setPs2Energy(0);
550  (*pfCandidates_)[tmpi[0]].addElementInBlock(blockref, kTrack[0]);
551  return true;
552  } else {
553  LogTrace("PFAlgo|recoTracksNotHCAL")
554  << " the closest track to ECAL " << thisEcal << " is " << sortedTracks.begin()->second
555  << " which is the one being processed. Will skip ECAL linking for all other track";
556  sortedTracks.clear();
557  }
558  }
559 
560  for (auto const& trk : sortedTracks) {
561  unsigned jTrack = trk.second;
562 
563  // Don't consider already used tracks
564  if (!active[jTrack])
565  continue;
566 
567  // The loop is on the other tracks !
568  if (jTrack == iTrack)
569  continue;
570 
571  // Check if the ECAL cluster closest to this track is
572  // indeed the current ECAL cluster. Only tracks not
573  // linked to an HCAL cluster are considered here
574  // (GP: this is because if there's a jTrack linked
575  // to the same Ecal cluster and that also has an Hcal link
576  // we would have also linked iTrack to that Hcal above)
577  std::multimap<double, unsigned> sortedECAL;
578  block.associatedElements(jTrack, linkData, sortedECAL, reco::PFBlockElement::ECAL, reco::PFBlock::LINKTEST_ALL);
579  if (sortedECAL.begin()->second != thisEcal)
580  continue;
581 
582  // Check if this track is a muon
583  bool thatIsAMuon = PFMuonAlgo::isMuon(elements[jTrack]);
584  LogTrace("PFAlgo|recoTracksNotHCAL") << " found track " << jTrack << (thatIsAMuon ? " (muon) " : " (non-muon)")
585  << ", with distance = " << sortedECAL.begin()->first;
586 
587  // Check if this track is not a fake
588  bool rejectFake = false;
589  reco::TrackRef trackRef = elements[jTrack].trackRef();
590  if (!thatIsAMuon && trackRef->ptError() > ptError_) {
591  // GP: FIXME this selection should be adapted depending on hasDeadHcal
592  // but we don't know if jTrack is linked to a dead Hcal or not
593  double deficit = trackMomentum + trackRef->p() - clusterRef->energy();
594  double resol =
595  nSigmaTRACK_ * neutralHadronEnergyResolution(trackMomentum + trackRef->p(), clusterRef->positionREP().Eta());
596  resol *= (trackMomentum + trackRef->p());
597  if (deficit > nSigmaTRACK_ * resol && !goodTrackDeadHcal) {
598  rejectFake = true;
599  kTrack.push_back(jTrack);
600  active[jTrack] = false;
601  LogTrace("PFAlgo|recoTracksNotHCAL")
602  << elements[jTrack] << std::endl
603  << "is probably a fake (2) --> lock the track"
604  << "(trackMomentum = " << trackMomentum << ", clusterEnergy = " << clusterRef->energy()
605  << ", deficit = " << deficit << " > " << nSigmaTRACK_ << " x " << resol
606  << " assuming neutralHadronEnergyResolution "
607  << neutralHadronEnergyResolution(trackMomentum + trackRef->p(), clusterRef->positionREP().Eta()) << ") ";
608  }
609  }
610  if (rejectFake)
611  continue;
612 
613  // Otherwise, add this track momentum to the total track momentum
614  /* */
615  // reco::TrackRef trackRef = elements[jTrack].trackRef();
616  if (!thatIsAMuon) {
617  LogTrace("PFAlgo|recoTracksNotHCAL") << "Track momentum increased from " << trackMomentum << " GeV ";
618  trackMomentum += trackRef->p();
619  LogTrace("PFAlgo|recoTracksNotHCAL") << "to " << trackMomentum << " GeV.";
620  LogTrace("PFAlgo|recoTracksNotHCAL") << "with " << elements[jTrack];
621  } else {
622  totalEcal -= muonECAL_[0];
623  totalEcal = std::max(totalEcal, 0.);
624  }
625 
626  // And create a charged particle candidate !
627 
628  tmpi.push_back(reconstructTrack(elements[jTrack]));
629 
630  kTrack.push_back(jTrack);
631  active[jTrack] = false;
632 
633  if (thatIsAMuon) {
634  (*pfCandidates_)[tmpi.back()].setEcalEnergy(clusterRef->energy(), std::min(clusterRef->energy(), muonECAL_[0]));
635  (*pfCandidates_)[tmpi.back()].setHcalEnergy(0., 0.);
636  (*pfCandidates_)[tmpi.back()].setHoEnergy(0., 0.);
637  (*pfCandidates_)[tmpi.back()].setPs1Energy(0);
638  (*pfCandidates_)[tmpi.back()].setPs2Energy(0);
639  (*pfCandidates_)[tmpi.back()].addElementInBlock(blockref, kTrack.back());
640  }
641  }
642 
643  LogTrace("PFAlgo|recoTracksNotHCAL") << "Loop over all associated ECAL clusters";
644  // Loop over all ECAL linked clusters ordered by increasing distance.
645  for (auto const& ecal : ecalElems) {
646  const unsigned index = ecal.second;
647  const PFBlockElement::Type type = elements[index].type();
648  assert(type == PFBlockElement::ECAL);
649  LogTrace("PFAlgo|recoTracksNotHCAL") << elements[index];
650 
651  // Just skip clusters already taken
652  if (!active[index]) {
653  LogTrace("PFAlgo|recoTracksNotHCAL") << "is not active - ignore ";
654  continue;
655  }
656 
657  // Just skip this cluster if it's closer to another track
658  block.associatedElements(index, linkData, sortedTracks, reco::PFBlockElement::TRACK, reco::PFBlock::LINKTEST_ALL);
659 
660  const bool skip = std::none_of(
661  kTrack.begin(), kTrack.end(), [&](unsigned iTrack) { return sortedTracks.begin()->second == iTrack; });
662 
663  if (skip) {
664  LogTrace("PFAlgo|recoTracksNotHCAL") << "is closer to another track - ignore ";
665  continue;
666  }
667 
668  // The corresponding PFCluster and energy
669  const reco::PFClusterRef clusterRef = elements[index].clusterRef();
670  assert(!clusterRef.isNull());
671 
672  LogTrace("PFAlgo|recoTracksNotHCAL") << "Ecal cluster with raw energy = " << clusterRef->energy()
673  << " linked with distance = " << ecal.first;
674 
675  // Check the presence of preshower clusters in the vicinity
676  // Preshower cluster closer to another ECAL cluster are ignored.
677  //JOSH: This should use the association map already produced by the cluster corrector for consistency,
678  //but should be ok for now
679  vector<double> ps1Ene{0.};
680  associatePSClusters(index, reco::PFBlockElement::PS1, block, elements, linkData, active, ps1Ene);
681  vector<double> ps2Ene{0.};
682  associatePSClusters(index, reco::PFBlockElement::PS2, block, elements, linkData, active, ps2Ene);
683 
684  // KH: use raw ECAL energy for PF hadron calibration. use calibrated ECAL energy when adding PF photons
685  const double ecalEnergy = clusterRef->energy();
686  const double ecalEnergyCalibrated = clusterRef->correctedEnergy(); // calibrated based on the egamma hypothesis
687  LogTrace("PFAlgo|recoTracksNotHCAL") << "Corrected ECAL(+PS) energy = " << ecalEnergy;
688 
689  // Since the electrons were found beforehand, this track must be a hadron. Calibrate
690  // the energy under the hadron hypothesis.
691  totalEcal += ecalEnergy;
692  const double previousCalibEcal = calibEcal;
693  const double previousSlopeEcal = slopeEcal;
694  calibEcal = std::max(totalEcal, 0.);
695  calibHcal = 0.;
696  calibration_.energyEmHad(
697  trackMomentum, calibEcal, calibHcal, clusterRef->positionREP().Eta(), clusterRef->positionREP().Phi());
698  if (totalEcal > 0.)
699  slopeEcal = calibEcal / totalEcal;
700 
701  LogTrace("PFAlgo|recoTracksNotHCAL") << "The total calibrated energy so far amounts to = " << calibEcal
702  << " (slope = " << slopeEcal << ")";
703 
704  // Stop the loop when adding more ECAL clusters ruins the compatibility
705  if (connectedToEcal && calibEcal - trackMomentum >= 0.) {
706  // if ( connectedToEcal && calibEcal - trackMomentum >=
707  // nSigmaECAL_*neutralHadronEnergyResolution(trackMomentum,clusterRef->positionREP().Eta()) ) {
708  calibEcal = previousCalibEcal;
709  slopeEcal = previousSlopeEcal;
710  totalEcal = calibEcal / slopeEcal;
711 
712  // Turn this last cluster in a photon
713  // (The PS clusters are already locked in "associatePSClusters")
714  active[index] = false;
715 
716  // Find the associated tracks
717  std::multimap<double, unsigned> assTracks;
718  block.associatedElements(index, linkData, assTracks, reco::PFBlockElement::TRACK, reco::PFBlock::LINKTEST_ALL);
719 
720  auto& ecalCand = (*pfCandidates_)[reconstructCluster(
721  *clusterRef, ecalEnergyCalibrated)]; // KH: use the PF ECAL cluster calibrated energy
722  ecalCand.setEcalEnergy(clusterRef->energy(), ecalEnergyCalibrated);
723  ecalCand.setHcalEnergy(0., 0.);
724  ecalCand.setHoEnergy(0., 0.);
725  ecalCand.setPs1Energy(ps1Ene[0]);
726  ecalCand.setPs2Energy(ps2Ene[0]);
727  ecalCand.addElementInBlock(blockref, index);
728  // Check that there is at least one track
729  if (!assTracks.empty()) {
730  ecalCand.addElementInBlock(blockref, assTracks.begin()->second);
731 
732  // Assign the position of the track at the ECAL entrance
733  const ::math::XYZPointF& chargedPosition =
734  dynamic_cast<const reco::PFBlockElementTrack*>(&elements[assTracks.begin()->second])
735  ->positionAtECALEntrance();
736  ecalCand.setPositionAtECALEntrance(chargedPosition);
737  }
738  break;
739  }
740 
741  // Lock used clusters.
742  connectedToEcal = true;
743  iEcal = index;
744  active[index] = false;
745  for (unsigned ic : tmpi)
746  (*pfCandidates_)[ic].addElementInBlock(blockref, iEcal);
747 
748  } // Loop ecal elements
749 
750  bool bNeutralProduced = false;
751 
752  // Create a photon if the ecal energy is too large
753  if (connectedToEcal) {
754  reco::PFClusterRef pivotalRef = elements[iEcal].clusterRef();
755 
756  double neutralEnergy = calibEcal - trackMomentum;
757 
758  /*
759  // Check if there are other tracks linked to that ECAL
760  for(IE ie = sortedTracks.begin(); ie != sortedTracks.end() && neutralEnergy > 0; ++ie ) {
761  unsigned jTrack = ie->second;
762 
763  // The loop is on the other tracks !
764  if ( jTrack == iTrack ) continue;
765 
766  // Check if the ECAL closest to this track is the current ECAL
767  // Otherwise ignore this track in the neutral energy determination
768  std::multimap<double, unsigned> sortedECAL;
769  block.associatedElements( jTrack, linkData,
770  sortedECAL,
771  reco::PFBlockElement::ECAL,
772  reco::PFBlock::LINKTEST_ALL );
773  if ( sortedECAL.begin()->second != iEcal ) continue;
774 
775  // Check if this track is also linked to an HCAL
776  // (in which case it goes to the next loop and is ignored here)
777  std::multimap<double, unsigned> sortedHCAL;
778  block.associatedElements( jTrack, linkData,
779  sortedHCAL,
780  reco::PFBlockElement::HCAL,
781  reco::PFBlock::LINKTEST_ALL );
782  if ( sortedHCAL.size() ) continue;
783 
784  // Otherwise, subtract this track momentum from the ECAL energy
785  reco::TrackRef trackRef = elements[jTrack].trackRef();
786  neutralEnergy -= trackRef->p();
787 
788  } // End other tracks
789  */
790 
791  // Add a photon if the energy excess is large enough
792  double resol = neutralHadronEnergyResolution(trackMomentum, pivotalRef->positionREP().Eta());
793  resol *= trackMomentum;
794  if (neutralEnergy > std::max(0.5, nSigmaECAL_ * resol)) {
795  neutralEnergy /= slopeEcal;
796  unsigned tmpj = reconstructCluster(*pivotalRef, neutralEnergy);
797  (*pfCandidates_)[tmpj].setEcalEnergy(pivotalRef->energy(), neutralEnergy);
798  (*pfCandidates_)[tmpj].setHcalEnergy(0., 0.);
799  (*pfCandidates_)[tmpj].setHoEnergy(0., 0.);
800  (*pfCandidates_)[tmpj].setPs1Energy(0.);
801  (*pfCandidates_)[tmpj].setPs2Energy(0.);
802  (*pfCandidates_)[tmpj].addElementInBlock(blockref, iEcal);
803  bNeutralProduced = true;
804  for (unsigned ic = 0; ic < kTrack.size(); ++ic)
805  (*pfCandidates_)[tmpj].addElementInBlock(blockref, kTrack[ic]);
806  } // End neutral energy
807 
808  // Set elements in blocks and ECAL energies to all tracks
809  for (unsigned ic = 0; ic < tmpi.size(); ++ic) {
810  // Skip muons
811  if ((*pfCandidates_)[tmpi[ic]].particleId() == reco::PFCandidate::mu)
812  continue;
813 
814  double fraction = trackMomentum > 0 ? (*pfCandidates_)[tmpi[ic]].trackRef()->p() / trackMomentum : 0;
815  double ecalCal = bNeutralProduced ? (calibEcal - neutralEnergy * slopeEcal) * fraction : calibEcal * fraction;
816  double ecalRaw = totalEcal * fraction;
817 
818  LogTrace("PFAlgo|recoTracksNotHCAL")
819  << "The fraction after photon supression is " << fraction << " calibrated ecal = " << ecalCal;
820 
821  (*pfCandidates_)[tmpi[ic]].setEcalEnergy(ecalRaw, ecalCal);
822  (*pfCandidates_)[tmpi[ic]].setHcalEnergy(0., 0.);
823  (*pfCandidates_)[tmpi[ic]].setHoEnergy(0., 0.);
824  (*pfCandidates_)[tmpi[ic]].setPs1Energy(0);
825  (*pfCandidates_)[tmpi[ic]].setPs2Energy(0);
826  (*pfCandidates_)[tmpi[ic]].addElementInBlock(blockref, kTrack[ic]);
827  }
828 
829  } // End connected ECAL
830 
831  // Fill the element_in_block for tracks that are eventually linked to no ECAL clusters at all.
832  for (unsigned ic = 0; ic < tmpi.size(); ++ic) {
833  const PFCandidate& pfc = (*pfCandidates_)[tmpi[ic]];
834  const PFCandidate::ElementsInBlocks& eleInBlocks = pfc.elementsInBlocks();
835  if (eleInBlocks.empty()) {
836  LogTrace("PFAlgo|recoTracksNotHCAL") << "Single track / Fill element in block! ";
837  (*pfCandidates_)[tmpi[ic]].addElementInBlock(blockref, kTrack[ic]);
838  }
839  }
840  LogTrace("PFAlgo|recoTracksNotHCAL") << "end of function PFAlgo::recoTracksNotHCAL";
841  return false;
842 }

References funct::abs(), cms::cuda::assert(), associatePSClusters(), groupFilesInBlocks::block, calibration_, dptRel_DispVtx_, ECAL, reco::PFBlockElement::ECAL, bookConverter::elements, reco::PFCandidate::elementsInBlocks(), PFEnergyCalibration::energyEmHad(), PVValHelper::eta, HLT_FULL_cff::fraction, isFromSecInt(), PFMuonAlgo::isMuon(), edm::Ref< C, T, F >::isNull(), phase2tkutil::isPrimary(), reco::PFBlock::LINKTEST_ALL, LogTrace, SiStripPI::max, min(), reco::PFCandidate::mu, muonECAL_, neutralHadronEnergyResolution(), nHits, nSigmaECAL_, nSigmaTRACK_, pfCandidates_, reco::Vertex::position(), primaryVertex_, reco::PFBlockElement::PS1, reco::PFBlockElement::PS2, ptError_, reconstructCluster(), reconstructTrack(), rejectTracks_Step45_, l1ParticleFlow_cff::resol, optionsL1T::skip, PFTrackAlgoTools::step45(), reco::PFBlockElement::TRACK, and reco::btau::trackMomentum.

Referenced by elementLoop().

◆ relinkTrackToHcal()

void PFAlgo::relinkTrackToHcal ( const reco::PFBlock &  block,
std::multimap< double, unsigned > &  ecalElems,
std::multimap< double, unsigned > &  hcalElems,
const std::vector< bool > &  active,
reco::PFBlock::LinkData &  linkData,
unsigned int  iTrack 
)
private

Definition at line 933 of file PFAlgo.cc.

938  {
939  unsigned ntt = 1;
940  unsigned index = ecalElems.begin()->second;
941  std::multimap<double, unsigned> sortedTracks;
942  block.associatedElements(index, linkData, sortedTracks, reco::PFBlockElement::TRACK, reco::PFBlock::LINKTEST_ALL);
943  LogTrace("PFAlgo|elementLoop") << "The closest ECAL cluster is linked to " << sortedTracks.size()
944  << " tracks, with distance = " << ecalElems.begin()->first;
945 
946  LogTrace("PFAlgo|elementLoop") << "Looping over sortedTracks";
947  // Loop over all tracks
948  for (auto const& trk : sortedTracks) {
949  unsigned jTrack = trk.second;
950  LogTrace("PFAlgo|elementLoop") << "jTrack=" << jTrack;
951  // Track must be active
952  if (!active[jTrack])
953  continue;
954  LogTrace("PFAlgo|elementLoop") << "active[jTrack]=" << active[jTrack];
955 
956  // The loop is on the other tracks !
957  if (jTrack == iTrack)
958  continue;
959  LogTrace("PFAlgo|elementLoop") << "skipping jTrack=" << jTrack << " for same iTrack";
960 
961  // Check if the ECAL closest to this track is the current ECAL
962  // Otherwise ignore this track in the neutral energy determination
963  std::multimap<double, unsigned> sortedECAL;
964  block.associatedElements(jTrack, linkData, sortedECAL, reco::PFBlockElement::ECAL, reco::PFBlock::LINKTEST_ALL);
965  if (sortedECAL.begin()->second != index)
966  continue;
967  LogTrace("PFAlgo|elementLoop") << " track " << jTrack << " with closest ECAL identical ";
968 
969  // Check if this track is also linked to an HCAL
970  std::multimap<double, unsigned> sortedHCAL;
971  block.associatedElements(jTrack, linkData, sortedHCAL, reco::PFBlockElement::HCAL, reco::PFBlock::LINKTEST_ALL);
972  if (sortedHCAL.empty())
973  continue;
974  LogTrace("PFAlgo|elementLoop") << " and with an HCAL cluster " << sortedHCAL.begin()->second;
975  ntt++;
976 
977  // In that case establish a link with the first track
978  block.setLink(iTrack, sortedHCAL.begin()->second, sortedECAL.begin()->first, linkData, PFBlock::LINKTEST_RECHIT);
979 
980  } // End other tracks
981 
982  // Redefine HCAL elements
983  block.associatedElements(iTrack, linkData, hcalElems, reco::PFBlockElement::HCAL, reco::PFBlock::LINKTEST_ALL);
984 
985  if (!hcalElems.empty())
986  LogTrace("PFAlgo|elementLoop") << "Track linked back to HCAL due to ECAL sharing with other tracks";
987 }

References groupFilesInBlocks::block, reco::PFBlockElement::ECAL, reco::PFBlockElement::HCAL, reco::PFBlock::LINKTEST_ALL, LogTrace, and reco::PFBlockElement::TRACK.

Referenced by elementLoop().

◆ setCandConnectorParameters() [1/2]

void PFAlgo::setCandConnectorParameters ( bool  bCorrect,
bool  bCalibPrimary,
double  dptRel_PrimaryTrack,
double  dptRel_MergedTrack,
double  ptErrorSecondary,
const std::vector< double > &  nuclCalibFactors 
)
inline

Definition at line 72 of file PFAlgo.h.

77  {
78  connector_.setParameters(
79  bCorrect, bCalibPrimary, dptRel_PrimaryTrack, dptRel_MergedTrack, ptErrorSecondary, nuclCalibFactors);
80  }

References HLT_FULL_cff::bCalibPrimary, HLT_FULL_cff::bCorrect, connector_, HLT_FULL_cff::nuclCalibFactors, and PFCandConnector::setParameters().

◆ setCandConnectorParameters() [2/2]

void PFAlgo::setCandConnectorParameters ( const edm::ParameterSet &  iCfgCandConnector)
inline

Definition at line 68 of file PFAlgo.h.

68  {
69  connector_.setParameters(iCfgCandConnector);
70  }

References connector_, HLT_FULL_cff::iCfgCandConnector, and PFCandConnector::setParameters().

◆ setDisplacedVerticesParameters()

void PFAlgo::setDisplacedVerticesParameters ( bool  rejectTracks_Bad,
bool  rejectTracks_Step45,
bool  usePFNuclearInteractions,
bool  usePFConversions,
bool  usePFDecays,
double  dptRel_DispVtx 
)

Definition at line 96 of file PFAlgo.cc.

101  {
102  rejectTracks_Bad_ = rejectTracks_Bad;
103  rejectTracks_Step45_ = rejectTracks_Step45;
104  usePFNuclearInteractions_ = usePFNuclearInteractions;
105  usePFConversions_ = usePFConversions;
106  usePFDecays_ = usePFDecays;
107  dptRel_DispVtx_ = dptRel_DispVtx;
108 }

References HLT_FULL_cff::dptRel_DispVtx, dptRel_DispVtx_, HLT_FULL_cff::rejectTracks_Bad, rejectTracks_Bad_, HLT_FULL_cff::rejectTracks_Step45, rejectTracks_Step45_, HLT_FULL_cff::usePFConversions, usePFConversions_, HLT_FULL_cff::usePFDecays, usePFDecays_, HLT_FULL_cff::usePFNuclearInteractions, and usePFNuclearInteractions_.

◆ setEGammaCollections()

void PFAlgo::setEGammaCollections ( const edm::View< reco::PFCandidate > &  pfEgammaCandidates,
const edm::ValueMap< reco::GsfElectronRef > &  valueMapGedElectrons,
const edm::ValueMap< reco::PhotonRef > &  valueMapGedPhotons 
)

Definition at line 76 of file PFAlgo.cc.

78  {
79  if (useEGammaFilters_) {
80  pfEgammaCandidates_ = &pfEgammaCandidates;
81  valueMapGedElectrons_ = &valueMapGedElectrons;
82  valueMapGedPhotons_ = &valueMapGedPhotons;
83  }
84 }

References gedPhotonCore_cfi::pfEgammaCandidates, pfEgammaCandidates_, useEGammaFilters_, valueMapGedElectrons_, and valueMapGedPhotons_.

Referenced by PFProducer::produce().

◆ setEGammaParameters()

void PFAlgo::setEGammaParameters ( bool  use_EGammaFilters,
bool  useProtectionsForJetMET 
)

Definition at line 68 of file PFAlgo.cc.

68  {
69  useEGammaFilters_ = use_EGammaFilters;
70 
71  if (!useEGammaFilters_)
72  return;
73 
74  useProtectionsForJetMET_ = useProtectionsForJetMET;
75 }

References useEGammaFilters_, HLT_FULL_cff::useProtectionsForJetMET, and useProtectionsForJetMET_.

◆ setEGElectronCollection()

void PFAlgo::setEGElectronCollection ( const reco::GsfElectronCollection &  egelectrons)

◆ setHcalDepthInfo()

void PFAlgo::setHcalDepthInfo ( reco::PFCandidate &  cand,
const reco::PFCluster &  cluster 
) const
private

Definition at line 3323 of file PFAlgo.cc.

3323  {
3324  std::array<double, 7> energyPerDepth;
3325  std::fill(energyPerDepth.begin(), energyPerDepth.end(), 0.0);
3326  for (auto& hitRefAndFrac : cluster.recHitFractions()) {
3327  const auto& hit = *hitRefAndFrac.recHitRef();
3328  if (DetId(hit.detId()).det() == DetId::Hcal) {
3329  if (hit.depth() == 0) {
3330  edm::LogWarning("setHcalDepthInfo") << "Depth zero found";
3331  continue;
3332  }
3333  if (hit.depth() < 1 || hit.depth() > 7) {
3334  throw cms::Exception("CorruptData") << "Bogus depth " << hit.depth() << " at detid " << hit.detId() << "\n";
3335  }
3336  energyPerDepth[hit.depth() - 1] += hitRefAndFrac.fraction() * hit.energy();
3337  }
3338  }
3339  double sum = std::accumulate(energyPerDepth.begin(), energyPerDepth.end(), 0.);
3340  std::array<float, 7> depthFractions;
3341  if (sum > 0) {
3342  for (unsigned int i = 0; i < depthFractions.size(); ++i) {
3343  depthFractions[i] = energyPerDepth[i] / sum;
3344  }
3345  } else {
3346  std::fill(depthFractions.begin(), depthFractions.end(), 0.f);
3347  }
3348  cand.setHcalDepthEnergyFractions(depthFractions);
3349 }

References Exception, ntuplemaker::fill, DetId::Hcal, mps_fire::i, and reco::PFCluster::recHitFractions().

Referenced by createCandidatesHCAL(), and reconstructCluster().

◆ setHOTag()

void PFAlgo::setHOTag ( bool  ho)
inline

Definition at line 65 of file PFAlgo.h.

65 { useHO_ = ho; }

References photonIsolationHIProducer_cfi::ho, and useHO_.

◆ setMuonHandle()

void PFAlgo::setMuonHandle ( const edm::Handle< reco::MuonCollection > &  muons)
inline

Definition at line 66 of file PFAlgo.h.

66 { muonHandle_ = muons; }

References muonHandle_, and PDWG_BPHSkim_cff::muons.

Referenced by PFProducer::produce().

◆ setPFVertexParameters()

void PFAlgo::setPFVertexParameters ( bool  useVertex,
reco::VertexCollection const &  primaryVertices 
)

Definition at line 110 of file PFAlgo.cc.

110  {
111  useVertices_ = useVertex;
112 
113  //Set the vertices for muon cleaning
114  pfmu_->setInputsForCleaning(primaryVertices);
115 
116  //Now find the primary vertex!
117  bool primaryVertexFound = false;
118  nVtx_ = primaryVertices.size();
119  for (auto const& vertex : primaryVertices) {
120  if (vertex.isValid() && (!vertex.isFake())) {
121  primaryVertex_ = vertex;
122  primaryVertexFound = true;
123  break;
124  }
125  }
126  //Use vertices if the user wants to but only if it exists a good vertex
127  useVertices_ = useVertex && primaryVertexFound;
128 }

References nVtx_, pfmu_, primaryVertex_, HLT_FULL_cff::primaryVertices, HLT_FULL_cff::useVertex, useVertices_, and bphysicsOniaDQM_cfi::vertex.

Referenced by PFProducer::produce().

◆ setPostHFCleaningParameters()

void PFAlgo::setPostHFCleaningParameters ( bool  postHFCleaning,
const edm::ParameterSet &  pfHFCleaningParams 
)

Definition at line 86 of file PFAlgo.cc.

86  {
87  postHFCleaning_ = postHFCleaning;
88  minHFCleaningPt_ = pfHFCleaningParams.getParameter<double>("minHFCleaningPt");
89  minSignificance_ = pfHFCleaningParams.getParameter<double>("minSignificance");
90  maxSignificance_ = pfHFCleaningParams.getParameter<double>("maxSignificance");
91  minSignificanceReduction_ = pfHFCleaningParams.getParameter<double>("minSignificanceReduction");
92  maxDeltaPhiPt_ = pfHFCleaningParams.getParameter<double>("maxDeltaPhiPt");
93  minDeltaMet_ = pfHFCleaningParams.getParameter<double>("minDeltaMet");
94 }

References edm::ParameterSet::getParameter(), maxDeltaPhiPt_, maxSignificance_, minDeltaMet_, minHFCleaningPt_, minSignificance_, minSignificanceReduction_, HLT_FULL_cff::postHFCleaning, and postHFCleaning_.

Friends And Related Function Documentation

◆ operator<<

std::ostream& operator<< ( std::ostream &  out,
const PFAlgo &  algo 
)
friend

Member Data Documentation

◆ calibration_

PFEnergyCalibration& PFAlgo::calibration_
private

Definition at line 257 of file PFAlgo.h.

Referenced by createCandidatesHCAL(), createCandidatesHCALUnlinked(), and recoTracksNotHCAL().

◆ connector_

PFCandConnector PFAlgo::connector_
private

A tool used for a postprocessing of displaced vertices based on reconstructed PFCandidates

Definition at line 290 of file PFAlgo.h.

Referenced by makeConnectedCandidates(), and setCandConnectorParameters().

◆ dptRel_DispVtx_

double PFAlgo::dptRel_DispVtx_
private

Maximal relative uncertainty on the tracks going to or incoming from the displcaed vertex to be used in the PFAlgo

Definition at line 285 of file PFAlgo.h.

Referenced by createCandidatesHCAL(), reconstructTrack(), recoTracksNotHCAL(), and setDisplacedVerticesParameters().

◆ factors45_

std::vector<double> PFAlgo::factors45_
private

Definition at line 298 of file PFAlgo.h.

Referenced by createCandidatesHCAL(), and PFAlgo().

◆ goodPixelTrackDeadHcal_chi2n_

float PFAlgo::goodPixelTrackDeadHcal_chi2n_
private

Definition at line 310 of file PFAlgo.h.

Referenced by checkGoodTrackDeadHcal(), and PFAlgo().

◆ goodPixelTrackDeadHcal_dxy_

float PFAlgo::goodPixelTrackDeadHcal_dxy_
private

Definition at line 313 of file PFAlgo.h.

Referenced by checkGoodTrackDeadHcal(), and PFAlgo().

◆ goodPixelTrackDeadHcal_dz_

float PFAlgo::goodPixelTrackDeadHcal_dz_
private

Definition at line 314 of file PFAlgo.h.

Referenced by checkGoodTrackDeadHcal(), and PFAlgo().

◆ goodPixelTrackDeadHcal_maxLost3Hit_

int PFAlgo::goodPixelTrackDeadHcal_maxLost3Hit_
private

Definition at line 311 of file PFAlgo.h.

Referenced by checkGoodTrackDeadHcal(), and PFAlgo().

◆ goodPixelTrackDeadHcal_maxLost4Hit_

int PFAlgo::goodPixelTrackDeadHcal_maxLost4Hit_
private

Definition at line 312 of file PFAlgo.h.

Referenced by checkGoodTrackDeadHcal(), and PFAlgo().

◆ goodPixelTrackDeadHcal_maxPt_

float PFAlgo::goodPixelTrackDeadHcal_maxPt_
private

Definition at line 308 of file PFAlgo.h.

Referenced by checkGoodTrackDeadHcal(), and PFAlgo().

◆ goodPixelTrackDeadHcal_minEta_

float PFAlgo::goodPixelTrackDeadHcal_minEta_
private

Definition at line 307 of file PFAlgo.h.

Referenced by checkGoodTrackDeadHcal(), and PFAlgo().

◆ goodPixelTrackDeadHcal_ptErrRel_

float PFAlgo::goodPixelTrackDeadHcal_ptErrRel_
private

Definition at line 309 of file PFAlgo.h.

Referenced by checkGoodTrackDeadHcal(), and PFAlgo().

◆ goodTrackDeadHcal_chi2n_

float PFAlgo::goodTrackDeadHcal_chi2n_
private

Definition at line 302 of file PFAlgo.h.

Referenced by checkGoodTrackDeadHcal(), and PFAlgo().

◆ goodTrackDeadHcal_dxy_

float PFAlgo::goodTrackDeadHcal_dxy_
private

Definition at line 305 of file PFAlgo.h.

Referenced by checkGoodTrackDeadHcal(), and PFAlgo().

◆ goodTrackDeadHcal_layers_

int PFAlgo::goodTrackDeadHcal_layers_
private

Definition at line 303 of file PFAlgo.h.

Referenced by checkGoodTrackDeadHcal(), and PFAlgo().

◆ goodTrackDeadHcal_ptErrRel_

float PFAlgo::goodTrackDeadHcal_ptErrRel_
private

Variables for track cleaning in bad HCal areas.

Definition at line 301 of file PFAlgo.h.

Referenced by checkGoodTrackDeadHcal(), and PFAlgo().

◆ goodTrackDeadHcal_validFr_

float PFAlgo::goodTrackDeadHcal_validFr_
private

Definition at line 304 of file PFAlgo.h.

Referenced by checkGoodTrackDeadHcal(), and PFAlgo().

◆ maxDeltaPhiPt_

double PFAlgo::maxDeltaPhiPt_
private

Definition at line 323 of file PFAlgo.h.

Referenced by postCleaning(), and setPostHFCleaningParameters().

◆ maxSignificance_

double PFAlgo::maxSignificance_
private

Definition at line 321 of file PFAlgo.h.

Referenced by postCleaning(), and setPostHFCleaningParameters().

◆ minDeltaMet_

double PFAlgo::minDeltaMet_
private

Definition at line 324 of file PFAlgo.h.

Referenced by checkCleaning(), postCleaning(), and setPostHFCleaningParameters().

◆ minHFCleaningPt_

double PFAlgo::minHFCleaningPt_
private

Definition at line 319 of file PFAlgo.h.

Referenced by postCleaning(), and setPostHFCleaningParameters().

◆ minSignificance_

double PFAlgo::minSignificance_
private

Definition at line 320 of file PFAlgo.h.

Referenced by postCleaning(), and setPostHFCleaningParameters().

◆ minSignificanceReduction_

double PFAlgo::minSignificanceReduction_
private

Definition at line 322 of file PFAlgo.h.

Referenced by postCleaning(), and setPostHFCleaningParameters().

◆ muonECAL_

std::vector<double> PFAlgo::muonECAL_
private

Definition at line 294 of file PFAlgo.h.

Referenced by createCandidatesHCAL(), PFAlgo(), and recoTracksNotHCAL().

◆ muonHandle_

edm::Handle<reco::MuonCollection> PFAlgo::muonHandle_
private

Definition at line 331 of file PFAlgo.h.

Referenced by reconstructParticles(), and setMuonHandle().

◆ muonHCAL_

std::vector<double> PFAlgo::muonHCAL_
private

Variables for muons and fakes.

Definition at line 293 of file PFAlgo.h.

Referenced by createCandidatesHCAL(), and PFAlgo().

◆ muonHO_

std::vector<double> PFAlgo::muonHO_
private

Definition at line 295 of file PFAlgo.h.

Referenced by createCandidatesHCAL(), and PFAlgo().

◆ nSigmaECAL_

const double PFAlgo::nSigmaECAL_
private

number of sigma to judge energy excess in ECAL

Definition at line 245 of file PFAlgo.h.

Referenced by recoTracksNotHCAL().

◆ nSigmaEConstHCAL

const double PFAlgo::nSigmaEConstHCAL = 100.
private

Definition at line 334 of file PFAlgo.h.

Referenced by nSigmaHCAL().

◆ nSigmaEConstHFEM

const double PFAlgo::nSigmaEConstHFEM = 100.
private

Definition at line 335 of file PFAlgo.h.

Referenced by nSigmaHFEM().

◆ nSigmaEConstHFHAD

const double PFAlgo::nSigmaEConstHFHAD = 100.
private

Definition at line 336 of file PFAlgo.h.

Referenced by nSigmaHFHAD().

◆ nSigmaHCAL_

const double PFAlgo::nSigmaHCAL_
private

number of sigma to judge energy excess in HCAL

Definition at line 248 of file PFAlgo.h.

Referenced by nSigmaHCAL().

◆ nSigmaHFEM_

const double PFAlgo::nSigmaHFEM_
private

number of sigma to judge energy excess in HF

Definition at line 251 of file PFAlgo.h.

Referenced by nSigmaHFEM().

◆ nSigmaHFHAD_

const double PFAlgo::nSigmaHFHAD_
private

Definition at line 252 of file PFAlgo.h.

Referenced by nSigmaHFHAD().

◆ nSigmaTRACK_

double PFAlgo::nSigmaTRACK_
private

Definition at line 296 of file PFAlgo.h.

Referenced by createCandidatesHCAL(), PFAlgo(), and recoTracksNotHCAL().

◆ nVtx_

int PFAlgo::nVtx_
private

Definition at line 286 of file PFAlgo.h.

Referenced by egammaFilters(), and setPFVertexParameters().

◆ pfCandidates_

std::unique_ptr<reco::PFCandidateCollection> PFAlgo::pfCandidates_
private

Definition at line 226 of file PFAlgo.h.

Referenced by checkCleaning(), createCandidatesHCAL(), egammaFilters(), makeConnectedCandidates(), postCleaning(), reconstructCluster(), reconstructParticles(), reconstructTrack(), and recoTracksNotHCAL().

◆ pfCleanedCandidates_

reco::PFCandidateCollection PFAlgo::pfCleanedCandidates_
private

Definition at line 228 of file PFAlgo.h.

Referenced by getCleanedCandidates(), postCleaning(), and reconstructParticles().

◆ pfEgammaCandidates_

const edm::View<reco::PFCandidate>* PFAlgo::pfEgammaCandidates_
private

Definition at line 267 of file PFAlgo.h.

Referenced by egammaFilters(), and setEGammaCollections().

◆ pfmu_

std::unique_ptr<PFMuonAlgo> PFAlgo::pfmu_
private

Definition at line 262 of file PFAlgo.h.

Referenced by getPFMuonAlgo(), PFAlgo(), reconstructParticles(), reconstructTrack(), and setPFVertexParameters().

◆ postHFCleaning_

bool PFAlgo::postHFCleaning_
private

Definition at line 317 of file PFAlgo.h.

Referenced by reconstructParticles(), and setPostHFCleaningParameters().

◆ postMuonCleaning_

bool PFAlgo::postMuonCleaning_
private

Definition at line 318 of file PFAlgo.h.

◆ primaryVertex_

reco::Vertex PFAlgo::primaryVertex_
private

Definition at line 328 of file PFAlgo.h.

Referenced by checkGoodTrackDeadHcal(), egammaFilters(), reconstructCluster(), recoTracksNotHCAL(), and setPFVertexParameters().

◆ ptError_

double PFAlgo::ptError_
private

Definition at line 297 of file PFAlgo.h.

Referenced by createCandidatesHCAL(), PFAlgo(), and recoTracksNotHCAL().

◆ rejectTracks_Bad_

bool PFAlgo::rejectTracks_Bad_
private

Flags to use the protection against fakes and not reconstructed displaced vertices

Definition at line 276 of file PFAlgo.h.

Referenced by createCandidatesHCAL(), and setDisplacedVerticesParameters().

◆ rejectTracks_Step45_

bool PFAlgo::rejectTracks_Step45_
private

Definition at line 277 of file PFAlgo.h.

Referenced by createCandidatesHCAL(), recoTracksNotHCAL(), and setDisplacedVerticesParameters().

◆ resolHF_square_

const std::vector<double> PFAlgo::resolHF_square_
private

Definition at line 255 of file PFAlgo.h.

Referenced by hfEnergyResolution(), and PFAlgo().

◆ thepfEnergyCalibrationHF_

PFEnergyCalibrationHF& PFAlgo::thepfEnergyCalibrationHF_
private

Definition at line 258 of file PFAlgo.h.

Referenced by createCandidatesHF().

◆ useBestMuonTrack_

double PFAlgo::useBestMuonTrack_
private

Definition at line 325 of file PFAlgo.h.

◆ useEGammaFilters_

bool PFAlgo::useEGammaFilters_
private

Variables for NEW EGAMMA selection.

Definition at line 265 of file PFAlgo.h.

Referenced by processBlock(), setEGammaCollections(), and setEGammaParameters().

◆ useHO_

bool PFAlgo::useHO_
private

Definition at line 260 of file PFAlgo.h.

Referenced by createCandidatesHCAL(), createCandidatesHCALUnlinked(), decideType(), and setHOTag().

◆ usePFConversions_

bool PFAlgo::usePFConversions_
private

Definition at line 280 of file PFAlgo.h.

Referenced by processBlock(), and setDisplacedVerticesParameters().

◆ usePFDecays_

bool PFAlgo::usePFDecays_
private

Definition at line 281 of file PFAlgo.h.

Referenced by isFromSecInt(), and setDisplacedVerticesParameters().

◆ usePFMuonMomAssign_

bool PFAlgo::usePFMuonMomAssign_
private

Definition at line 272 of file PFAlgo.h.

◆ usePFNuclearInteractions_

bool PFAlgo::usePFNuclearInteractions_
private

Definition at line 279 of file PFAlgo.h.

Referenced by isFromSecInt(), and setDisplacedVerticesParameters().

◆ useProtectionsForJetMET_

bool PFAlgo::useProtectionsForJetMET_
private

Definition at line 266 of file PFAlgo.h.

Referenced by egammaFilters(), and setEGammaParameters().

◆ useVertices_

bool PFAlgo::useVertices_ = false
private

Definition at line 329 of file PFAlgo.h.

Referenced by reconstructCluster(), and setPFVertexParameters().

◆ valueMapGedElectrons_

const edm::ValueMap<reco::GsfElectronRef>* PFAlgo::valueMapGedElectrons_
private

Definition at line 268 of file PFAlgo.h.

Referenced by setEGammaCollections().

◆ valueMapGedPhotons_

const edm::ValueMap<reco::PhotonRef>* PFAlgo::valueMapGedPhotons_
private

Definition at line 269 of file PFAlgo.h.

Referenced by setEGammaCollections().

runTheMatrix.ret
ret
prodAgent to be discontinued
Definition: runTheMatrix.py:542
PDWG_BPHSkim_cff.muons
muons
Definition: PDWG_BPHSkim_cff.py:47
reco::PFBlockElement::isTimeValid
bool isTimeValid() const
do we have a valid time information
Definition: PFBlockElement.h:128
PFAlgo::goodPixelTrackDeadHcal_maxLost3Hit_
int goodPixelTrackDeadHcal_maxLost3Hit_
Definition: PFAlgo.h:311
reco::PFBlockElement::trackType
virtual bool trackType(TrackType trType) const
Definition: PFBlockElement.h:72
reco::PFBlockElement::HO
Definition: PFBlockElement.h:42
HLT_FULL_cff.postHFCleaning
postHFCleaning
Definition: HLT_FULL_cff.py:13642
PFAlgo::nSigmaEConstHFHAD
const double nSigmaEConstHFHAD
Definition: PFAlgo.h:336
PFAlgo::muonECAL_
std::vector< double > muonECAL_
Definition: PFAlgo.h:294
mps_fire.i
i
Definition: mps_fire.py:428
PFAlgo::calibration_
PFEnergyCalibration & calibration_
Definition: PFAlgo.h:257
PFMuonAlgo::isIsolatedMuon
static bool isIsolatedMuon(const reco::PFBlockElement &elt)
Definition: PFMuonAlgo.cc:94
PFAlgo::useEGammaFilters_
bool useEGammaFilters_
Variables for NEW EGAMMA selection.
Definition: PFAlgo.h:265
HLT_FULL_cff.track
track
Definition: HLT_FULL_cff.py:11713
PFAlgo::nSigmaHFEM
double nSigmaHFEM(double clusterEnergy) const
Definition: PFAlgo.cc:3382
l1ParticleFlow_cff.resol
resol
Definition: l1ParticleFlow_cff.py:17
PFAlgo::pfmu_
std::unique_ptr< PFMuonAlgo > pfmu_
Definition: PFAlgo.h:262
PFAlgo::goodTrackDeadHcal_validFr_
float goodTrackDeadHcal_validFr_
Definition: PFAlgo.h:304
PFAlgo::goodPixelTrackDeadHcal_maxLost4Hit_
int goodPixelTrackDeadHcal_maxLost4Hit_
Definition: PFAlgo.h:312
hit::y
double y
Definition: SiStripHitEffFromCalibTree.cc:90
HLT_FULL_cff.postMuonCleaning
postMuonCleaning
Definition: HLT_FULL_cff.py:13639
PFAlgo::decideType
int decideType(const edm::OwnVector< reco::PFBlockElement > &elements, const reco::PFBlockElement::Type type, std::vector< bool > &active, ElementIndices &inds, std::vector< bool > &deadArea, unsigned int iEle)
Definition: PFAlgo.cc:1193
PFAlgo::goodPixelTrackDeadHcal_chi2n_
float goodPixelTrackDeadHcal_chi2n_
Definition: PFAlgo.h:310
reco::PFCandidate::e
Definition: PFCandidate.h:47
PFAlgo::nVtx_
int nVtx_
Definition: PFAlgo.h:286
PFAlgo::goodTrackDeadHcal_chi2n_
float goodTrackDeadHcal_chi2n_
Definition: PFAlgo.h:302
PFAlgo::useVertices_
bool useVertices_
Definition: PFAlgo.h:329
custom_jme_cff.nMuons
nMuons
Definition: custom_jme_cff.py:148
muon
Definition: MuonCocktails.h:17
PFAlgo::goodPixelTrackDeadHcal_dz_
float goodPixelTrackDeadHcal_dz_
Definition: PFAlgo.h:314
PFAlgo::postCleaning
void postCleaning()
Definition: PFAlgo.cc:3491
DiDispStaMuonMonitor_cfi.pt
pt
Definition: DiDispStaMuonMonitor_cfi.py:39
reco::TrackBase::p
double p() const
momentum vector magnitude
Definition: TrackBase.h:631
muons2muons_cfi.chargedHadron
chargedHadron
Definition: muons2muons_cfi.py:26
reco::PFCluster::layer
PFLayer::Layer layer() const
cluster layer, see PFLayer.h in this directory
Definition: PFCluster.cc:56
X
#define X(str)
Definition: MuonsGrabber.cc:38
PFTrackAlgoTools::step45
bool step45(const reco::TrackBase::TrackAlgorithm &)
Definition: PFTrackAlgoTools.cc:221
PFAlgo::isFromSecInt
bool isFromSecInt(const reco::PFBlockElement &eTrack, std::string order) const
Definition: PFAlgo.cc:3462
min
T min(T a, T b)
Definition: MathUtil.h:58
reco::Vertex::z
double z() const
z coordinate
Definition: Vertex.h:120
reco::PFCandidate::mva_e_pi
float mva_e_pi() const
mva for electron-pion discrimination
Definition: PFCandidate.h:317
edm::View::refAt
RefToBase< value_type > refAt(size_type i) const
HLT_FULL_cff.useProtectionsForJetMET
useProtectionsForJetMET
Definition: HLT_FULL_cff.py:13736
PFTrackAlgoTools::step5
bool step5(const reco::TrackBase::TrackAlgorithm &)
Definition: PFTrackAlgoTools.cc:232
reco::PFBlockElement::time
float time() const
Definition: PFBlockElement.h:130
PFAlgo::valueMapGedElectrons_
const edm::ValueMap< reco::GsfElectronRef > * valueMapGedElectrons_
Definition: PFAlgo.h:268
hcal
Definition: ConfigurationDatabase.cc:13
edm::Ref::isNull
bool isNull() const
Checks for null.
Definition: Ref.h:235
multPhiCorr_741_25nsDY_cfi.py
py
Definition: multPhiCorr_741_25nsDY_cfi.py:12
AlCaHLTBitMon_ParallelJobs.p
p
Definition: AlCaHLTBitMon_ParallelJobs.py:153
reco::PFBlock
Block of elements.
Definition: PFBlock.h:26
reco::PFBlockElement::T_TO_DISP
Definition: PFBlockElement.h:47
ElementIndices::hfHadIs
std::vector< unsigned > hfHadIs
Definition: PFAlgo.h:50
HLT_FULL_cff.nuclCalibFactors
nuclCalibFactors
Definition: HLT_FULL_cff.py:13646
PFAlgo::rejectTracks_Step45_
bool rejectTracks_Step45_
Definition: PFAlgo.h:277
ElementIndices::ecalIs
std::vector< unsigned > ecalIs
Definition: PFAlgo.h:44
PFEnergyCalibration::energyEmHad
void energyEmHad(double t, double &e, double &h, double eta, double phi) const
Definition: PFEnergyCalibration.cc:195
HLT_FULL_cff.usePFConversions
usePFConversions
Definition: HLT_FULL_cff.py:13633
DetId::Hcal
Definition: DetId.h:28
PFAlgo::createCandidatesHCAL
void createCandidatesHCAL(const reco::PFBlock &block, reco::PFBlock::LinkData &linkData, const edm::OwnVector< reco::PFBlockElement > &elements, std::vector< bool > &active, const reco::PFBlockRef &blockref, ElementIndices &inds, std::vector< bool > &deadArea)
Definition: PFAlgo.cc:1694
PFLayer::HCAL_ENDCAP
Definition: PFLayer.h:37
reco::PFCluster::recHitFractions
const std::vector< reco::PFRecHitFraction > & recHitFractions() const
vector of rechit fractions
Definition: PFCluster.h:65
reco::PFCandidate::h
Definition: PFCandidate.h:46
PFAlgo::nSigmaEConstHFEM
const double nSigmaEConstHFEM
Definition: PFAlgo.h:335
reco::PFCandidate::h_HF
Definition: PFCandidate.h:51
XYZVector
math::XYZVector XYZVector
Definition: RawParticle.h:26
cms::cuda::assert
assert(be >=bs)
CustomPhysics_cfi.gamma
gamma
Definition: CustomPhysics_cfi.py:17
edm::second
U second(std::pair< T, U > const &p)
Definition: ParameterSet.cc:222
PFAlgo::dptRel_DispVtx_
double dptRel_DispVtx_
Definition: PFAlgo.h:285
HLT_FULL_cff.resolHF_square
resolHF_square
Definition: HLT_FULL_cff.py:13636
reco::TrackBase::px
double px() const
x coordinate of momentum vector
Definition: TrackBase.h:640
reco::Vertex::position
const Point & position() const
position
Definition: Vertex.h:114
reco::PFBlockElement::HCAL
Definition: PFBlockElement.h:36
DDAxes::x
PFAlgo::pfEgammaCandidates_
const edm::View< reco::PFCandidate > * pfEgammaCandidates_
Definition: PFAlgo.h:267
reco::LeafCandidate::pt
double pt() const final
transverse momentum
Definition: LeafCandidate.h:146
PFAlgo::muonHO_
std::vector< double > muonHO_
Definition: PFAlgo.h:295
edm::LogInfo
Log< level::Info, false > LogInfo
Definition: MessageLogger.h:125
findQualityFiles.v
v
Definition: findQualityFiles.py:179
PFAlgo::goodPixelTrackDeadHcal_ptErrRel_
float goodPixelTrackDeadHcal_ptErrRel_
Definition: PFAlgo.h:309
reco::PFCandidate::elementsInBlocks
const ElementsInBlocks & elementsInBlocks() const
Definition: PFCandidate.cc:636
createJobs.tmp
tmp
align.sh
Definition: createJobs.py:716
ecal
Definition: ElectronicsMappingGPU.h:13
convertSQLiteXML.ok
bool ok
Definition: convertSQLiteXML.py:98
PFAlgo::minHFCleaningPt_
double minHFCleaningPt_
Definition: PFAlgo.h:319
reco::PFBlockRef
edm::Ref< PFBlockCollection > PFBlockRef
persistent reference to PFCluster objects
Definition: PFBlockFwd.h:19
edm::LogWarning
Log< level::Warning, false > LogWarning
Definition: MessageLogger.h:122
PFAlgo::createCandidatesHCALUnlinked
void createCandidatesHCALUnlinked(const reco::PFBlock &block, reco::PFBlock::LinkData &linkData, const edm::OwnVector< reco::PFBlockElement > &elements, std::vector< bool > &active, const reco::PFBlockRef &blockref, ElementIndices &inds, std::vector< bool > &deadArea)
Definition: PFAlgo.cc:2805
PFAlgo::nSigmaECAL_
const double nSigmaECAL_
number of sigma to judge energy excess in ECAL
Definition: PFAlgo.h:245
reco::PFCandidate::X
Definition: PFCandidate.h:45
edm::Ref< PFClusterCollection >
hltEgammaHLTExtra_cfi.hcal
hcal
Definition: hltEgammaHLTExtra_cfi.py:42
PFLayer::ECAL_BARREL
Definition: PFLayer.h:33
reco::PFBlockElement::Type
Type
Definition: PFBlockElement.h:30
quality
const uint32_t *__restrict__ Quality * quality
Definition: CAHitNtupletGeneratorKernelsImpl.h:109
hit::x
double x
Definition: SiStripHitEffFromCalibTree.cc:89
PFAlgo::nSigmaHCAL
double nSigmaHCAL(double clusterEnergy, double clusterEta) const
Definition: PFAlgo.cc:3363
reco::PFCandidate::egammaExtraRef
reco::PFCandidateEGammaExtraRef egammaExtraRef() const
return a reference to the EGamma extra
Definition: PFCandidate.cc:554
PFAlgo::minDeltaMet_
double minDeltaMet_
Definition: PFAlgo.h:324
DetId
Definition: DetId.h:17
cmsdt::algo
algo
Definition: constants.h:165
optionsL1T.skip
skip
Definition: optionsL1T.py:30
pfElectronTranslator_cfi.PFCandidate
PFCandidate
Definition: pfElectronTranslator_cfi.py:6
HLT_FULL_cff.rejectTracks_Bad
rejectTracks_Bad
Definition: HLT_FULL_cff.py:13650
PFAlgo::usePFDecays_
bool usePFDecays_
Definition: PFAlgo.h:281
ECAL
Definition: HCALResponse.h:21
RPCNoise_example.check
check
Definition: RPCNoise_example.py:71
DivergingColor.frac
float frac
Definition: DivergingColor.py:175
PFAlgo::elementLoop
void elementLoop(const reco::PFBlock &block, reco::PFBlock::LinkData &linkData, const edm::OwnVector< reco::PFBlockElement > &elements, std::vector< bool > &active, const reco::PFBlockRef &blockref, ElementIndices &inds, std::vector< bool > &deadArea)
Definition: PFAlgo.cc:989
reco::PFCandidate::mu
Definition: PFCandidate.h:48
reco::TrackBase::py
double py() const
y coordinate of momentum vector
Definition: TrackBase.h:643
PFLayer::HF_EM
Definition: PFLayer.h:38
SiPixelRawToDigiRegional_cfi.deltaPhi
deltaPhi
Definition: SiPixelRawToDigiRegional_cfi.py:9
nHits
const caConstants::TupleMultiplicity *__restrict__ const HitsOnGPU *__restrict__ double *__restrict__ float *__restrict__ double *__restrict__ uint32_t nHits
Definition: BrokenLineFitOnGPU.h:27
PFAlgo::conversionAlgo
void conversionAlgo(const edm::OwnVector< reco::PFBlockElement > &elements, std::vector< bool > &active)
Definition: PFAlgo.cc:362
Calorimetry_cff.dp
dp
Definition: Calorimetry_cff.py:158
HLT_FULL_cff.fraction
fraction
Definition: HLT_FULL_cff.py:52802
reco::LeafCandidate::py
double py() const final
y coordinate of momentum vector
Definition: LeafCandidate.h:142
PFAlgo::createCandidatesHF
void createCandidatesHF(const reco::PFBlock &block, reco::PFBlock::LinkData &linkData, const edm::OwnVector< reco::PFBlockElement > &elements, std::vector< bool > &active, const reco::PFBlockRef &blockref, ElementIndices &inds)
Definition: PFAlgo.cc:1251
reco::PFBlockElement::TRACK
Definition: PFBlockElement.h:32
PVValHelper::eta
Definition: PVValidationHelpers.h:70
mathSSE::sqrt
T sqrt(T t)
Definition: SSEVec.h:19
PFLayer::HCAL_BARREL1
Definition: PFLayer.h:35
PFEnergyCalibrationHF::energyHad
double energyHad(double uncalibratedEnergyHCAL, double eta, double phi)
Definition: PFEnergyCalibrationHF.cc:63
reco::Track
Definition: Track.h:27
reco::btau::trackMomentum
Definition: TaggingVariable.h:41
PFEnergyCalibrationHF::energyEmHad
double energyEmHad(double uncalibratedEnergyECAL, double uncalibratedEnergyHCAL, double eta, double phi)
Definition: PFEnergyCalibrationHF.cc:76
reco::CaloCluster::badHcalMarker
Definition: CaloCluster.h:50
PFAlgo::checkAndReconstructSecondaryInteraction
bool checkAndReconstructSecondaryInteraction(const reco::PFBlockRef &blockref, const edm::OwnVector< reco::PFBlockElement > &elements, bool isActive, int iElement)
Definition: PFAlgo.cc:847
DQMScaleToClient_cfi.factor
factor
Definition: DQMScaleToClient_cfi.py:8
trackingPlots.other
other
Definition: trackingPlots.py:1460
HCALHighEnergyHPDFilter_cfi.energy
energy
Definition: HCALHighEnergyHPDFilter_cfi.py:5
PFAlgo::reconstructTrack
unsigned reconstructTrack(const reco::PFBlockElement &elt, bool allowLoose=false)
Definition: PFAlgo.cc:3133
HCAL
Definition: HCALResponse.h:21
DigiToRawDM_cff.HO
HO
Definition: DigiToRawDM_cff.py:23
b
double b
Definition: hdecay.h:118
PFAlgo::muonHCAL_
std::vector< double > muonHCAL_
Variables for muons and fakes.
Definition: PFAlgo.h:293
PFAlgo::minSignificanceReduction_
double minSignificanceReduction_
Definition: PFAlgo.h:322
reco::LeafCandidate::setVertex
void setVertex(const Point &vertex) override
set vertex
Definition: LeafCandidate.h:173
first
auto first
Definition: CAHitNtupletGeneratorKernelsImpl.h:112
HLT_FULL_cff.usePFNuclearInteractions
usePFNuclearInteractions
Definition: HLT_FULL_cff.py:13729
hit::z
double z
Definition: SiStripHitEffFromCalibTree.cc:91
reco::PFCandidate::ElementsInBlocks
std::vector< ElementInBlock > ElementsInBlocks
Definition: PFCandidate.h:379
gedPhotonCore_cfi.pfEgammaCandidates
pfEgammaCandidates
Definition: gedPhotonCore_cfi.py:8
reco::PFBlockElement::T_FROM_V0
Definition: PFBlockElement.h:47
edm::View::size
size_type size() const
reco::PFCandidate::set_mva_Isolated
void set_mva_Isolated(float mvaI)
Definition: PFCandidate.h:311
PFAlgo::nSigmaHFHAD
double nSigmaHFHAD(double clusterEnergy) const
Definition: PFAlgo.cc:3387
reco::PFBlockElement::TrackType
TrackType
Definition: PFBlockElement.h:47
PFAlgo::maxDeltaPhiPt_
double maxDeltaPhiPt_
Definition: PFAlgo.h:323
PFAlgo::goodTrackDeadHcal_layers_
int goodTrackDeadHcal_layers_
Definition: PFAlgo.h:303
ALCARECOTkAlJpsiMuMu_cff.charge
charge
Definition: ALCARECOTkAlJpsiMuMu_cff.py:47
PFLayer::HF_HAD
Definition: PFLayer.h:39
reco::PFBlockElement::T_FROM_DISP
Definition: PFBlockElement.h:47
ntuplemaker.fill
fill
Definition: ntuplemaker.py:304
reco::PFCandidate::set_mva_e_pi
void set_mva_e_pi(float mvaNI)
Definition: PFCandidate.h:315
PFAlgo::nSigmaHFEM_
const double nSigmaHFEM_
number of sigma to judge energy excess in HF
Definition: PFAlgo.h:251
bphysicsOniaDQM_cfi.vertex
vertex
Definition: bphysicsOniaDQM_cfi.py:7
LogDebug
#define LogDebug(id)
Definition: MessageLogger.h:233
edm::ParameterSet
Definition: ParameterSet.h:47
math::XYZPoint
XYZPointD XYZPoint
point in space with cartesian internal representation
Definition: Point3D.h:12
reco::PFBlockElement::HFEM
Definition: PFBlockElement.h:39
a
double a
Definition: hdecay.h:119
RecoEcal_cff.ecalClusters
ecalClusters
Definition: RecoEcal_cff.py:26
Phase1L1TJetCalibrator_cfi.calibration
calibration
Definition: Phase1L1TJetCalibrator_cfi.py:2
reco::PFBlockElement::GSF
Definition: PFBlockElement.h:37
PFCandConnector::setParameters
void setParameters(const edm::ParameterSet &iCfgCandConnector)
Definition: PFCandConnector.h:30
reco::LeafCandidate::eta
double eta() const final
momentum pseudorapidity
Definition: LeafCandidate.h:152
SiStripPI::max
Definition: SiStripPayloadInspectorHelper.h:169
PFAlgo::nSigmaHCAL_
const double nSigmaHCAL_
number of sigma to judge energy excess in HCAL
Definition: PFAlgo.h:248
reco::PFCandidate::particleId
virtual ParticleType particleId() const
Definition: PFCandidate.h:367
phase2tkutil::isPrimary
bool isPrimary(const SimTrack &simTrk, const PSimHit *simHit)
Definition: TrackerPhase2ValidationUtil.cc:2
type
type
Definition: SiPixelVCal_PayloadInspector.cc:37
reco::Vertex::x
double x() const
x coordinate
Definition: Vertex.h:116
edm::Ref::isNonnull
bool isNonnull() const
Checks for non-null.
Definition: Ref.h:238
reco::PFCandidate::T_FROM_DISP
Definition: PFCandidate.h:69
gainCalibHelper::gainCalibPI::type
type
Definition: SiPixelGainCalibHelper.h:40
PFAlgo::checkHasDeadHcal
bool checkHasDeadHcal(const std::multimap< double, unsigned > &hcalElems, const std::vector< bool > &deadArea)
Definition: PFAlgo.cc:867
reco::PFCandidate::egamma_HF
Definition: PFCandidate.h:52
HLT_FULL_cff.dptRel_DispVtx
dptRel_DispVtx
Definition: HLT_FULL_cff.py:13654
math::XYZVector
XYZVectorD XYZVector
spatial vector with cartesian internal representation
Definition: Vector3D.h:31
PFAlgo::ptError_
double ptError_
Definition: PFAlgo.h:297
reco::PFBlockElement::ECAL
Definition: PFBlockElement.h:35
PFAlgo::usePFNuclearInteractions_
bool usePFNuclearInteractions_
Definition: PFAlgo.h:279
cand
Definition: decayParser.h:32
eventshapeDQM_cfi.order
order
Definition: eventshapeDQM_cfi.py:8
PFAlgo::goodPixelTrackDeadHcal_minEta_
float goodPixelTrackDeadHcal_minEta_
Definition: PFAlgo.h:307
HLT_FULL_cff.useVertex
useVertex
Definition: HLT_FULL_cff.py:26519
PFAlgo::goodTrackDeadHcal_dxy_
float goodTrackDeadHcal_dxy_
Definition: PFAlgo.h:305
PFAlgo::checkGoodTrackDeadHcal
bool checkGoodTrackDeadHcal(const reco::TrackRef &trackRef, bool hasDeadHcal)
Definition: PFAlgo.cc:894
alignmentValidation.c1
c1
do drawing
Definition: alignmentValidation.py:1025
PFAlgo::rejectTracks_Bad_
bool rejectTracks_Bad_
Definition: PFAlgo.h:276
reco::PFBlockElement::timeError
float timeError() const
Definition: PFBlockElement.h:132
reco::PFBlockElement::T_FROM_GAMMACONV
Definition: PFBlockElement.h:47
PFAlgo::setHcalDepthInfo
void setHcalDepthInfo(reco::PFCandidate &cand, const reco::PFCluster &cluster) const
Definition: PFAlgo.cc:3323
HLT_FULL_cff.iCfgCandConnector
iCfgCandConnector
Definition: HLT_FULL_cff.py:13645
reco::PFCandidate::gamma
Definition: PFCandidate.h:49
photonIsolationHIProducer_cfi.ho
ho
Definition: photonIsolationHIProducer_cfi.py:10
reco::LeafCandidate::charge
int charge() const final
electric charge
Definition: LeafCandidate.h:106
HcalResponse_cfi.energyHF
energyHF
Definition: HcalResponse_cfi.py:1084
groupFilesInBlocks.block
block
Definition: groupFilesInBlocks.py:150
PFAlgo::goodTrackDeadHcal_ptErrRel_
float goodTrackDeadHcal_ptErrRel_
Variables for track cleaning in bad HCal areas.
Definition: PFAlgo.h:301
PFAlgo::createCandidatesECAL
void createCandidatesECAL(const reco::PFBlock &block, reco::PFBlock::LinkData &linkData, const edm::OwnVector< reco::PFBlockElement > &elements, std::vector< bool > &active, const reco::PFBlockRef &blockref, ElementIndices &inds, std::vector< bool > &deadArea)
Definition: PFAlgo.cc:3004
PFAlgo::maxSignificance_
double maxSignificance_
Definition: PFAlgo.h:321
reco::PFBlockElement::HFHAD
Definition: PFBlockElement.h:40
edm::LogError
Log< level::Error, false > LogError
Definition: MessageLogger.h:123
PFAlgo::primaryVertex_
reco::Vertex primaryVertex_
Definition: PFAlgo.h:328
PFAlgo::reconstructCluster
unsigned reconstructCluster(const reco::PFCluster &cluster, double particleEnergy, bool useDirection=false, double particleX=0., double particleY=0., double particleZ=0.)
Definition: PFAlgo.cc:3212
reco::PFBlock::LinkData
std::map< unsigned int, Link > LinkData
Definition: PFBlock.h:39
PFAlgo::nSigmaHFHAD_
const double nSigmaHFHAD_
Definition: PFAlgo.h:252
PFAlgo::valueMapGedPhotons_
const edm::ValueMap< reco::PhotonRef > * valueMapGedPhotons_
Definition: PFAlgo.h:269
PFAlgo::muonHandle_
edm::Handle< reco::MuonCollection > muonHandle_
Definition: PFAlgo.h:331
reco::PFBlock::LINKTEST_ALL
Definition: PFBlock.h:41
PFAlgo::useProtectionsForJetMET_
bool useProtectionsForJetMET_
Definition: PFAlgo.h:266
reco::TrackBase::conversionStep
Definition: TrackBase.h:102
XYZPoint
math::XYZVector XYZPoint
Definition: CalorimetryManager.cc:69
PFAlgo::recoTracksNotHCAL
bool recoTracksNotHCAL(const reco::PFBlock &block, reco::PFBlock::LinkData &linkData, const edm::OwnVector< reco::PFBlockElement > &elements, const reco::PFBlockRef &blockref, std::vector< bool > &active, bool goodTrackDeadHcal, bool hasDeadHcal, unsigned int iTrack, std::multimap< double, unsigned > &ecalElems, reco::TrackRef &trackRef)
Definition: PFAlgo.cc:391
reco::CaloCluster::position
const math::XYZPoint & position() const
cluster centroid position
Definition: CaloCluster.h:154
ElementIndices
Definition: PFAlgo.h:40
PFAlgo::pfCandidates_
std::unique_ptr< reco::PFCandidateCollection > pfCandidates_
Definition: PFAlgo.h:226
PFEnergyCalibrationHF::energyEm
double energyEm(double uncalibratedEnergyECAL, double eta, double phi)
Definition: PFEnergyCalibrationHF.cc:48
multPhiCorr_741_25nsDY_cfi.px
px
Definition: multPhiCorr_741_25nsDY_cfi.py:10
PFAlgo::nSigmaTRACK_
double nSigmaTRACK_
Definition: PFAlgo.h:296
reco::LeafCandidate::setP4
void setP4(const LorentzVector &p4) final
set 4-momentum
Definition: LeafCandidate.h:158
bookConverter.elements
elements
Definition: bookConverter.py:147
HLT_FULL_cff.rejectTracks_Step45
rejectTracks_Step45
Definition: HLT_FULL_cff.py:13723
reco::Vertex::y
double y() const
y coordinate
Definition: Vertex.h:118
reco::LeafCandidate::phi
double phi() const final
momentum azimuthal angle
Definition: LeafCandidate.h:148
PFAlgo::factors45_
std::vector< double > factors45_
Definition: PFAlgo.h:298
PFMuonAlgo::isMuon
static bool isMuon(const reco::PFBlockElement &elt)
Definition: PFMuonAlgo.cc:48
PFAlgo::nSigmaEConstHCAL
const double nSigmaEConstHCAL
Definition: PFAlgo.h:334
PFAlgo::connector_
PFCandConnector connector_
Definition: PFAlgo.h:290
PFAlgo::goodPixelTrackDeadHcal_dxy_
float goodPixelTrackDeadHcal_dxy_
Definition: PFAlgo.h:313
math::XYZTLorentzVector
XYZTLorentzVectorD XYZTLorentzVector
Lorentz vector with cylindrical internal representation using pseudorapidity.
Definition: LorentzVector.h:29
reco::PFCandidate::T_TO_DISP
Definition: PFCandidate.h:68
relativeConstraints.empty
bool empty
Definition: relativeConstraints.py:46
PFAlgo::egammaFilters
void egammaFilters(const reco::PFBlockRef &blockref, std::vector< bool > &active, PFEGammaFilters const *pfegamma)
Definition: PFAlgo.cc:220
PFAlgo::processBlock
void processBlock(const reco::PFBlockRef &blockref, std::list< reco::PFBlockRef > &hcalBlockRefs, std::list< reco::PFBlockRef > &ecalBlockRefs, PFEGammaFilters const *pfegamma)
Definition: PFAlgo.cc:3053
Exception
Definition: hltDiff.cc:245
reco::PFCandidateCollection
std::vector< reco::PFCandidate > PFCandidateCollection
collection of PFCandidates
Definition: PFCandidateFwd.h:12
reco::PFCluster
Particle flow cluster, see clustering algorithm in PFClusterAlgo.
Definition: PFCluster.h:42
reco::LeafCandidate::setCharge
void setCharge(Charge q) final
set electric charge
Definition: LeafCandidate.h:108
EgHLTOffHistBins_cfi.mass
mass
Definition: EgHLTOffHistBins_cfi.py:34
PFAlgo::usePFConversions_
bool usePFConversions_
Definition: PFAlgo.h:280
edm::ParameterSet::getParameter
T getParameter(std::string const &) const
Definition: ParameterSet.h:303
PFAlgo::pfCleanedCandidates_
reco::PFCandidateCollection pfCleanedCandidates_
Definition: PFAlgo.h:228
BeamSpotPI::Y
Definition: BeamSpotPayloadInspectorHelper.h:32
fftjetpileupestimator_calo_uncalib_cfi.c0
c0
Definition: fftjetpileupestimator_calo_uncalib_cfi.py:8
AlignmentPI::index
index
Definition: AlignmentPayloadInspectorHelper.h:46
PFAlgo::goodPixelTrackDeadHcal_maxPt_
float goodPixelTrackDeadHcal_maxPt_
Definition: PFAlgo.h:308
PFAlgo::hfEnergyResolution
double hfEnergyResolution(double clusterEnergy) const
Definition: PFAlgo.cc:3370
reco::PFRecHit
Particle flow rechit (rechit + geometry and topology information). See clustering algorithm in PFClus...
Definition: PFRecHit.h:31
reco::PFCandidate
Particle reconstructed by the particle flow algorithm.
Definition: PFCandidate.h:41
funct::pow
Power< A, B >::type pow(const A &a, const B &b)
Definition: Power.h:29
PFTrackAlgoTools::errorScale
double errorScale(const reco::TrackBase::TrackAlgorithm &, const std::vector< double > &)
Definition: PFTrackAlgoTools.cc:86
reco::TrackBase::pz
double pz() const
z coordinate of momentum vector
Definition: TrackBase.h:646
PFEnergyCalibrationHF::getcalibHF_use
const bool & getcalibHF_use() const
Definition: PFEnergyCalibrationHF.h:47
ElementIndices::hoIs
std::vector< unsigned > hoIs
Definition: PFAlgo.h:43
funct::abs
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
reco::PFBlockElement::PS1
Definition: PFBlockElement.h:33
LogTrace
#define LogTrace(id)
Definition: MessageLogger.h:234
ElementIndices::hfEmIs
std::vector< unsigned > hfEmIs
Definition: PFAlgo.h:49
BeamSpotPI::Z
Definition: BeamSpotPayloadInspectorHelper.h:33
PFCandConnector::connect
reco::PFCandidateCollection connect(reco::PFCandidateCollection &pfCand) const
Definition: PFCandConnector.cc:63
PbPb_ZMuSkimMuonDPG_cff.particleType
particleType
Definition: PbPb_ZMuSkimMuonDPG_cff.py:27
reco::PFCandidate::ParticleType
ParticleType
particle types
Definition: PFCandidate.h:44
dqmiolumiharvest.j
j
Definition: dqmiolumiharvest.py:66
JetChargeProducer_cfi.exp
exp
Definition: JetChargeProducer_cfi.py:6
HLT_FULL_cff.bCorrect
bCorrect
Definition: HLT_FULL_cff.py:13647
HLT_FULL_cff.usePFDecays
usePFDecays
Definition: HLT_FULL_cff.py:13721
edm::HandleBase::isValid
bool isValid() const
Definition: HandleBase.h:70
PFAlgo::useHO_
bool useHO_
Definition: PFAlgo.h:260
math::XYZPointF
ROOT::Math::PositionVector3D< ROOT::Math::Cartesian3D< float > > XYZPointF
point in space with cartesian internal representation
Definition: Point3D.h:10
HLT_FULL_cff.flags
flags
Definition: HLT_FULL_cff.py:13150
PFLayer::ECAL_ENDCAP
Definition: PFLayer.h:32
bTagCommon_cff.ratioMax
ratioMax
Definition: bTagCommon_cff.py:34
PFAlgo::neutralHadronEnergyResolution
double neutralHadronEnergyResolution(double clusterEnergy, double clusterEta) const
todo: use PFClusterTools for this
Definition: PFAlgo.cc:3353
PFAlgo::minSignificance_
double minSignificance_
Definition: PFAlgo.h:320
ElementIndices::hcalIs
std::vector< unsigned > hcalIs
Definition: PFAlgo.h:42
reco::LeafCandidate::px
double px() const final
x coordinate of momentum vector
Definition: LeafCandidate.h:140
PFAlgo::postHFCleaning_
bool postHFCleaning_
Definition: PFAlgo.h:317
gather_cfg.blocks
blocks
Definition: gather_cfg.py:90
HLT_FULL_cff.primaryVertices
primaryVertices
Definition: HLT_FULL_cff.py:51939
cuy.ii
ii
Definition: cuy.py:590
ElementIndices::trackIs
std::vector< unsigned > trackIs
Definition: PFAlgo.h:45
PFAlgo::relinkTrackToHcal
void relinkTrackToHcal(const reco::PFBlock &block, std::multimap< double, unsigned > &ecalElems, std::multimap< double, unsigned > &hcalElems, const std::vector< bool > &active, reco::PFBlock::LinkData &linkData, unsigned int iTrack)
Definition: PFAlgo.cc:933
reco::PFCandidate::setParticleType
void setParticleType(ParticleType type)
set Particle Type
Definition: PFCandidate.cc:256
PFMuonAlgo::isLooseMuon
static bool isLooseMuon(const reco::PFBlockElement &elt)
Definition: PFMuonAlgo.cc:57
PFAlgo::resolHF_square_
const std::vector< double > resolHF_square_
Definition: PFAlgo.h:255
hit
Definition: SiStripHitEffFromCalibTree.cc:88
muonDTDigis_cfi.pset
pset
Definition: muonDTDigis_cfi.py:27
reco::PFBlockElement::PS2
Definition: PFBlockElement.h:34
edm::OwnVector< reco::PFBlockElement >
met_cff.significance
significance
Definition: met_cff.py:18
PFAlgo::associatePSClusters
void associatePSClusters(unsigned iEcal, reco::PFBlockElement::Type psElementType, const reco::PFBlock &block, const edm::OwnVector< reco::PFBlockElement > &elements, const reco::PFBlock::LinkData &linkData, std::vector< bool > &active, std::vector< double > &psEne)
Associate PS clusters to a given ECAL cluster, and return their energy.
Definition: PFAlgo.cc:3417
reco::PFCandidate::set_mva_nothing_gamma
void set_mva_nothing_gamma(float mva)
set mva for gamma detection
Definition: PFCandidate.h:332
GetRecoTauVFromDQM_MC_cff.next
next
Definition: GetRecoTauVFromDQM_MC_cff.py:31
MillePedeFileConverter_cfg.e
e
Definition: MillePedeFileConverter_cfg.py:37
HLT_FULL_cff.bCalibPrimary
bCalibPrimary
Definition: HLT_FULL_cff.py:13648
PFAlgo::thepfEnergyCalibrationHF_
PFEnergyCalibrationHF & thepfEnergyCalibrationHF_
Definition: PFAlgo.h:258
reco::isMuon
bool isMuon(const Candidate &part)
Definition: pdgIdUtils.h:9
reco::TrackBase::highPurity
Definition: TrackBase.h:154
Generated for CMSSW Reference Manual by doxygen 1.8.16