PFProducer/interface/PFEGammaAlgo.h

0001 #ifndef PFProducer_PFEGammaAlgo_H
0002 #define PFProducer_PFEGammaAlgo_H
0003
0004 //
0005 // Rewrite for GED integration:  Lindsey Gray (FNAL): lagray@fnal.gov
0006 //
0007 // Original Authors: Fabian Stoeckli: fabian.stoeckli@cern.ch
0008 //                   Nicholas Wardle: nckw@cern.ch
0009 //                   Rishi Patel: rpatel@cern.ch
0010 //                   Josh Bendavid : Josh.Bendavid@cern.ch
0011 //
0012
0013 #include "DataFormats/ParticleFlowReco/interface/PFBlockFwd.h"
0014 #include "DataFormats/ParticleFlowReco/interface/PFBlock.h"
0015 #include "DataFormats/ParticleFlowReco/interface/PFCluster.h"
0016 #include "DataFormats/ParticleFlowCandidate/interface/PFCandidate.h"
0017 #include "DataFormats/ParticleFlowCandidate/interface/PFCandidatePhotonExtra.h"
0018 #include "DataFormats/ParticleFlowCandidate/interface/PFCandidatePhotonExtraFwd.h"
0019 #include "DataFormats/Common/interface/Handle.h"
0020 #include "DataFormats/Common/interface/OrphanHandle.h"
0021 #include "DataFormats/CaloRecHit/interface/CaloCluster.h"
0022 #include "DataFormats/CaloRecHit/interface/CaloClusterFwd.h"
0023 #include "DataFormats/GsfTrackReco/interface/GsfTrack.h"
0024 #include "DataFormats/GsfTrackReco/interface/GsfTrackExtra.h"
0025 #include "DataFormats/ParticleFlowReco/interface/PFBlockElementTrack.h"
0026 #include "DataFormats/ParticleFlowReco/interface/PFBlockElementSuperCluster.h"
0027 #include "DataFormats/ParticleFlowReco/interface/PFBlockElementBrem.h"
0028 #include "DataFormats/ParticleFlowReco/interface/PFBlockElementCluster.h"
0029 #include "DataFormats/ParticleFlowReco/interface/PFBlockElementGsfTrack.h"
0030 #include "DataFormats/VertexReco/interface/Vertex.h"
0031 #include "DataFormats/EgammaCandidates/interface/GsfElectronFwd.h"
0032 #include "DataFormats/EgammaCandidates/interface/GsfElectron.h"
0033 #include "DataFormats/ParticleFlowCandidate/interface/PFCandidateFwd.h"
0034 #include "DataFormats/ParticleFlowCandidate/interface/PFCandidateElectronExtraFwd.h"
0035 #include "DataFormats/ParticleFlowCandidate/interface/PFCandidateElectronExtra.h"
0036 #include "DataFormats/ParticleFlowCandidate/interface/PFCandidateEGammaExtraFwd.h"
0037 #include "DataFormats/ParticleFlowCandidate/interface/PFCandidateEGammaExtra.h"
0038 #include "DataFormats/EgammaReco/interface/ElectronSeed.h"
0039 #include "CondFormats/GBRForest/interface/GBRForest.h"
0040
0041 #include "CondFormats/ESObjects/interface/ESEEIntercalibConstants.h"
0042 #include "CondFormats/ESObjects/interface/ESChannelStatus.h"
0043
0044 #include <iostream>
0045 #include <TH2D.h>
0046
0047 #include <list>
0048 #include <forward_list>
0049 #include <unordered_map>
0050
0051 #include "RecoParticleFlow/PFProducer/interface/FlaggedPtr.h"
0052 #include "RecoParticleFlow/PFProducer/interface/CommutativePairs.h"
0053 #include "RecoParticleFlow/PFClusterTools/interface/PFEnergyCalibration.h"
0054
0055 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0056 #include "CommonTools/MVAUtils/interface/GBRForestTools.h"
0057
0058 #include <memory>
0059
0060 class PFEnergyCalibration;
0061
0062 class PFEGammaAlgo {
0063 public:
0064   typedef reco::PFCluster::EEtoPSAssociation EEtoPSAssociation;
0065   typedef reco::PFBlockElementSuperCluster PFSCElement;
0066   typedef reco::PFBlockElementBrem PFBremElement;
0067   typedef reco::PFBlockElementGsfTrack PFGSFElement;
0068   typedef reco::PFBlockElementTrack PFKFElement;
0069   typedef reco::PFBlockElementCluster PFClusterElement;
0070   typedef std::unordered_map<const PFKFElement*, float> KFValMap;
0071
0072   using ClusterMap = std::unordered_map<PFClusterElement const*, std::vector<PFClusterElement const*>>;
0073
0074   class GBRForests {
0075   public:
0076     GBRForests(const edm::ParameterSet& conf)
0077         : ele_(createGBRForest(conf.getParameter<edm::FileInPath>("pf_electronID_mvaWeightFile"))),
0078           singleLeg_(createGBRForest(conf.getParameter<edm::FileInPath>("pf_convID_mvaWeightFile"))) {}
0079
0080     const std::unique_ptr<const GBRForest> ele_;
0081     const std::unique_ptr<const GBRForest> singleLeg_;
0082   };
0083
0084   struct ProtoEGObject {
0085     reco::PFBlockRef parentBlock;
0086     const PFSCElement* parentSC = nullptr;  // if ECAL driven
0087     reco::ElectronSeedRef electronSeed;     // if there is one
0088     // this is a mutable list of clusters
0089     // if ECAL driven we take the PF SC and refine it
0090     // if Tracker driven we add things to it as we discover more valid clusters
0091     std::vector<FlaggedPtr<const PFClusterElement>> ecalclusters;
0092     ClusterMap ecal2ps;
0093     // associations to tracks of various sorts
0094     std::vector<PFGSFElement const*> primaryGSFs;
0095     std::vector<PFKFElement const*> primaryKFs;
0096     std::vector<PFBremElement const*> brems;  // these are tangent based brems
0097     // for manual brem recovery
0098     std::vector<PFGSFElement const*> secondaryGSFs;
0099     std::vector<PFKFElement const*> secondaryKFs;
0100     KFValMap singleLegConversionMvaMap;
0101     // for track-HCAL cluster linking
0102     std::vector<PFClusterElement const*> hcalClusters;
0103     CommutativePairs<const reco::PFBlockElement*> localMap;
0104     // cluster closest to the gsf track(s), primary kf if none for gsf
0105     // last brem tangent cluster if neither of those work
0106     std::vector<const PFClusterElement*> electronClusters;
0107     int firstBrem, lateBrem, nBremsWithClusters;
0108   };
0109
0110   struct PFEGConfigInfo {
0111     double mvaEleCut;
0112     bool applyCrackCorrections;
0113     bool produceEGCandsWithNoSuperCluster;
0114     double mvaConvCut;
0115   };
0116
0117   struct EgammaObjects {
0118     reco::PFCandidateCollection candidates;
0119     reco::PFCandidateEGammaExtraCollection candidateExtras;
0120     reco::SuperClusterCollection refinedSuperClusters;
0121   };
0122
0123   //constructor
0124   PFEGammaAlgo(const PFEGConfigInfo&,
0125                GBRForests const& gbrForests,
0126                EEtoPSAssociation const& eetops,
0127                ESEEIntercalibConstants const& esEEInterCalib,
0128                ESChannelStatus const& channelStatus,
0129                reco::Vertex const& primaryVertex);
0130
0131   // this runs the functions below
0132   EgammaObjects operator()(const reco::PFBlockRef& block);
0133
0134 private:
0135   GBRForests const& gbrForests_;
0136
0137   PFEnergyCalibration thePFEnergyCalibration_;
0138
0139   // ------ rewritten basic processing pieces and cleaning algorithms
0140
0141   // useful pre-cached mappings:
0142   // hopefully we get an enum that lets us just make an array in the future
0143   reco::PFCluster::EEtoPSAssociation const& eetops_;
0144   reco::PFBlockRef _currentblock;
0145   reco::PFBlock::LinkData _currentlinks;
0146   // keep a map of pf indices to the splayed block for convenience
0147   // sadly we're mashing together two ways of thinking about the block
0148   std::vector<std::vector<FlaggedPtr<const reco::PFBlockElement>>> _splayedblock;
0149
0150   // pre-cleaning for the splayed block
0151   bool isMuon(const reco::PFBlockElement&);
0152   // pre-processing of ECAL clusters near non-primary KF tracks
0153   void removeOrLinkECALClustersToKFTracks();
0154
0155   // functions:
0156
0157   // build proto eg object using all available unflagged resources in block.
0158   // this will be kind of like the old 'SetLinks' but with simplified and
0159   // maximally inclusive logic that builds a list of 'refinable' objects
0160   // that we will perform operations on to clean/remove as needed
0161   void initializeProtoCands(std::list<ProtoEGObject>&);
0162
0163   // turn a supercluster into a map of ECAL cluster elements
0164   // related to PS cluster elements
0165   bool unwrapSuperCluster(const reco::PFBlockElementSuperCluster*,
0166                           std::vector<FlaggedPtr<const PFClusterElement>>&,
0167                           ClusterMap&);
0168
0169   int attachPSClusters(const PFClusterElement*, ClusterMap::mapped_type&);
0170
0171   void dumpCurrentRefinableObjects() const;
0172
0173   // wax on
0174
0175   // the key merging operation, done after building up links
0176   void mergeROsByAnyLink(std::list<ProtoEGObject>&);
0177
0178   // refining steps you can do with tracks
0179   void linkRefinableObjectGSFTracksToKFs(ProtoEGObject&);
0180   void linkRefinableObjectPrimaryKFsToSecondaryKFs(ProtoEGObject&);
0181   void linkRefinableObjectPrimaryGSFTrackToECAL(ProtoEGObject&);
0182   void linkRefinableObjectPrimaryGSFTrackToHCAL(ProtoEGObject&);
0183   void linkRefinableObjectKFTracksToECAL(ProtoEGObject&);
0184   void linkRefinableObjectBremTangentsToECAL(ProtoEGObject&);
0185   // WARNING! this should be ONLY used after doing the ECAL->track
0186   // reverse lookup after the primary linking!
0187   void linkRefinableObjectConvSecondaryKFsToSecondaryKFs(ProtoEGObject&);
0188   void linkRefinableObjectSecondaryKFsToECAL(ProtoEGObject&);
0189   // helper function for above
0190   void linkKFTrackToECAL(PFKFElement const*, ProtoEGObject&);
0191
0192   // refining steps doing the ECAL -> track piece
0193   // this is the factorization of the old PF photon algo stuff
0194   // which through arcane means I came to understand was conversion matching
0195   void linkRefinableObjectECALToSingleLegConv(ProtoEGObject&);
0196
0197   // wax off
0198
0199   // refining steps remove things from the built-up objects
0200   // original bits were for removing bad KF tracks
0201   // new (experimental) piece to remove clusters associated to these tracks
0202   // behavior determined by bools passed to unlink_KFandECALMatchedToHCAL
0203   void unlinkRefinableObjectKFandECALWithBadEoverP(ProtoEGObject&);
0204   void unlinkRefinableObjectKFandECALMatchedToHCAL(ProtoEGObject&,
0205                                                    bool removeFreeECAL = false,
0206                                                    bool removeSCECAL = false);
0207
0208   // things for building the final candidate and refined SC collections
0209   EgammaObjects fillPFCandidates(const std::list<ProtoEGObject>&);
0210   reco::SuperCluster buildRefinedSuperCluster(const ProtoEGObject&);
0211
0212   // helper functions for that
0213
0214   float calculateEleMVA(const ProtoEGObject&, reco::PFCandidateEGammaExtra&) const;
0215   void fillExtraInfo(const ProtoEGObject&, reco::PFCandidateEGammaExtra&);
0216
0217   // ------ end of new stuff
0218
0219   bool isPrimaryTrack(const reco::PFBlockElementTrack& KfEl, const reco::PFBlockElementGsfTrack& GsfEl);
0220
0221   PFEGConfigInfo const& cfg_;
0222   reco::Vertex const& primaryVertex_;
0223
0224   ESChannelStatus const& channelStatus_;
0225
0226   float evaluateSingleLegMVA(const reco::PFBlockRef& blockref, const reco::Vertex& primaryVtx, unsigned int trackIndex);
0227 };
0228
0229 #endif