File indexing completed on 2023-03-17 11:17:22
0001 #include "CommonTools/Utils/interface/StringToEnumValue.h"
0002 #include "DataFormats/CaloRecHit/interface/CaloClusterFwd.h"
0003 #include "DataFormats/CaloRecHit/interface/CaloID.h"
0004 #include "DataFormats/Common/interface/Handle.h"
0005 #include "DataFormats/EcalDetId/interface/EBDetId.h"
0006 #include "DataFormats/EcalRecHit/interface/EcalRecHit.h"
0007 #include "DataFormats/EcalRecHit/interface/EcalRecHitCollections.h"
0008 #include "DataFormats/EgammaReco/interface/BasicCluster.h"
0009 #include "DataFormats/EgammaReco/interface/BasicClusterFwd.h"
0010 #include "FWCore/Framework/interface/ESHandle.h"
0011 #include "FWCore/Framework/interface/Event.h"
0012 #include "FWCore/Framework/interface/EventSetup.h"
0013 #include "FWCore/Framework/interface/stream/EDProducer.h"
0014 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0015 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0016 #include "FWCore/Utilities/interface/ESGetToken.h"
0017 #include "FWCore/Utilities/interface/Exception.h"
0018 #include "Geometry/CaloGeometry/interface/CaloCellGeometry.h"
0019 #include "Geometry/CaloGeometry/interface/CaloGeometry.h"
0020 #include "Geometry/CaloGeometry/interface/CaloSubdetectorGeometry.h"
0021 #include "Geometry/CaloTopology/interface/CaloSubdetectorTopology.h"
0022 #include "Geometry/CaloTopology/interface/EcalBarrelTopology.h"
0023 #include "Geometry/CaloTopology/interface/EcalEndcapTopology.h"
0024 #include "Geometry/Records/interface/CaloGeometryRecord.h"
0025 #include "RecoEcal/EgammaClusterAlgos/interface/Multi5x5ClusterAlgo.h"
0026 #include "RecoEcal/EgammaCoreTools/interface/PositionCalc.h"
0027
0028 #include <ctime>
0029 #include <iostream>
0030 #include <memory>
0031 #include <vector>
0032
0033 class Multi5x5ClusterProducer : public edm::stream::EDProducer<> {
0034 public:
0035 Multi5x5ClusterProducer(const edm::ParameterSet& ps);
0036
0037 ~Multi5x5ClusterProducer() override;
0038
0039 void produce(edm::Event&, const edm::EventSetup&) override;
0040
0041 private:
0042 int nMaxPrintout_;
0043 int nEvt_;
0044
0045
0046 bool doBarrel_;
0047 bool doEndcap_;
0048
0049 edm::EDGetTokenT<EcalRecHitCollection> barrelHitToken_;
0050 edm::EDGetTokenT<EcalRecHitCollection> endcapHitToken_;
0051 edm::ESGetToken<CaloGeometry, CaloGeometryRecord> caloGeometryToken_;
0052
0053 std::string barrelClusterCollection_;
0054 std::string endcapClusterCollection_;
0055
0056 PositionCalc posCalculator_;
0057 Multi5x5ClusterAlgo* island_p;
0058
0059 bool counterExceeded() const { return ((nEvt_ > nMaxPrintout_) || (nMaxPrintout_ < 0)); }
0060
0061 const EcalRecHitCollection* getCollection(edm::Event& evt, const edm::EDGetTokenT<EcalRecHitCollection>& token);
0062
0063 void clusterizeECALPart(edm::Event& evt,
0064 const edm::EventSetup& es,
0065 const edm::EDGetTokenT<EcalRecHitCollection>& token,
0066 const std::string& clusterCollection,
0067 const reco::CaloID::Detectors detector);
0068
0069 void outputValidationInfo(reco::CaloClusterPtrVector& clusterPtrVector);
0070 };
0071
0072 #include "FWCore/Framework/interface/MakerMacros.h"
0073 DEFINE_FWK_MODULE(Multi5x5ClusterProducer);
0074
0075 Multi5x5ClusterProducer::Multi5x5ClusterProducer(const edm::ParameterSet& ps) {
0076
0077 barrelHitToken_ = consumes<EcalRecHitCollection>(ps.getParameter<edm::InputTag>("barrelHitTag"));
0078
0079 endcapHitToken_ = consumes<EcalRecHitCollection>(ps.getParameter<edm::InputTag>("endcapHitTag"));
0080
0081
0082 caloGeometryToken_ = esConsumes<CaloGeometry, CaloGeometryRecord>();
0083
0084
0085 doEndcap_ = ps.getParameter<bool>("doEndcap");
0086 doBarrel_ = ps.getParameter<bool>("doBarrel");
0087
0088
0089 barrelClusterCollection_ = ps.getParameter<std::string>("barrelClusterCollection");
0090 endcapClusterCollection_ = ps.getParameter<std::string>("endcapClusterCollection");
0091
0092
0093 double barrelSeedThreshold = ps.getParameter<double>("IslandBarrelSeedThr");
0094 double endcapSeedThreshold = ps.getParameter<double>("IslandEndcapSeedThr");
0095
0096 const std::vector<std::string> flagnames = ps.getParameter<std::vector<std::string> >("RecHitFlagToBeExcluded");
0097
0098 const std::vector<int> v_chstatus = StringToEnumValue<EcalRecHit::Flags>(flagnames);
0099
0100
0101 edm::ParameterSet posCalcParameters = ps.getParameter<edm::ParameterSet>("posCalcParameters");
0102 posCalculator_ = PositionCalc(posCalcParameters);
0103
0104
0105 produces<reco::BasicClusterCollection>(endcapClusterCollection_);
0106 produces<reco::BasicClusterCollection>(barrelClusterCollection_);
0107
0108 bool reassignSeedCrysToClusterItSeeds = false;
0109 if (ps.exists("reassignSeedCrysToClusterItSeeds"))
0110 reassignSeedCrysToClusterItSeeds = ps.getParameter<bool>("reassignSeedCrysToClusterItSeeds");
0111
0112 island_p = new Multi5x5ClusterAlgo(
0113 barrelSeedThreshold, endcapSeedThreshold, v_chstatus, posCalculator_, reassignSeedCrysToClusterItSeeds);
0114
0115 nEvt_ = 0;
0116 }
0117
0118 Multi5x5ClusterProducer::~Multi5x5ClusterProducer() { delete island_p; }
0119
0120 void Multi5x5ClusterProducer::produce(edm::Event& evt, const edm::EventSetup& es) {
0121 if (doEndcap_) {
0122 clusterizeECALPart(evt, es, endcapHitToken_, endcapClusterCollection_, reco::CaloID::DET_ECAL_ENDCAP);
0123 }
0124 if (doBarrel_) {
0125 clusterizeECALPart(evt, es, barrelHitToken_, barrelClusterCollection_, reco::CaloID::DET_ECAL_BARREL);
0126 }
0127
0128 nEvt_++;
0129 }
0130
0131 const EcalRecHitCollection* Multi5x5ClusterProducer::getCollection(
0132 edm::Event& evt, const edm::EDGetTokenT<EcalRecHitCollection>& token) {
0133 edm::Handle<EcalRecHitCollection> rhcHandle;
0134 evt.getByToken(token, rhcHandle);
0135 return rhcHandle.product();
0136 }
0137
0138 void Multi5x5ClusterProducer::clusterizeECALPart(edm::Event& evt,
0139 const edm::EventSetup& es,
0140 const edm::EDGetTokenT<EcalRecHitCollection>& token,
0141 const std::string& clusterCollection,
0142 const reco::CaloID::Detectors detector) {
0143
0144 const EcalRecHitCollection* hitCollection_p = getCollection(evt, token);
0145
0146
0147 edm::ESHandle<CaloGeometry> geoHandle = es.getHandle(caloGeometryToken_);
0148
0149 const CaloSubdetectorGeometry* geometry_p;
0150 std::unique_ptr<CaloSubdetectorTopology> topology_p;
0151
0152 if (detector == reco::CaloID::DET_ECAL_BARREL) {
0153 geometry_p = geoHandle->getSubdetectorGeometry(DetId::Ecal, EcalBarrel);
0154 topology_p = std::make_unique<EcalBarrelTopology>(*geoHandle);
0155 } else {
0156 geometry_p = geoHandle->getSubdetectorGeometry(DetId::Ecal, EcalEndcap);
0157 topology_p = std::make_unique<EcalEndcapTopology>(*geoHandle);
0158 }
0159
0160 const CaloSubdetectorGeometry* geometryES_p;
0161 geometryES_p = geoHandle->getSubdetectorGeometry(DetId::Ecal, EcalPreshower);
0162
0163
0164 reco::BasicClusterCollection clusters;
0165 clusters = island_p->makeClusters(hitCollection_p, geometry_p, topology_p.get(), geometryES_p, detector);
0166
0167
0168 auto clusters_p = std::make_unique<reco::BasicClusterCollection>();
0169 clusters_p->assign(clusters.begin(), clusters.end());
0170 edm::OrphanHandle<reco::BasicClusterCollection> bccHandle;
0171 if (detector == reco::CaloID::DET_ECAL_BARREL)
0172 bccHandle = evt.put(std::move(clusters_p), barrelClusterCollection_);
0173 else
0174 bccHandle = evt.put(std::move(clusters_p), endcapClusterCollection_);
0175 }