Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-11 23:28:11

0001 #include <Eigen/Core>
0002 #include <Eigen/Dense>
0003 #include <algorithm>
0004 #include <array>
0005 #include <cmath>
0006 #include <iostream>
0007 #include <memory>
0008 #include <string>
0009 #include <type_traits>
0010 #include <unordered_map>
0011 #include <utility>
0012 #include <vector>
0013 
0014 #include "FWCore/Framework/interface/stream/EDProducer.h"
0015 #include "FWCore/Utilities/interface/StreamID.h"
0016 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0017 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0018 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0019 #include "FWCore/Utilities/interface/EDGetToken.h"
0020 #include "FWCore/Utilities/interface/InputTag.h"
0021 #include "FWCore/Framework/interface/Event.h"
0022 #include "FWCore/Framework/interface/EventSetup.h"
0023 #include "FWCore/Utilities/interface/EDPutToken.h"
0024 
0025 #include "CondFormats/DataRecord/interface/HcalPFCutsRcd.h"
0026 #include "CondTools/Hcal/interface/HcalPFCutsHandler.h"
0027 
0028 #include "DataFormats/Common/interface/Handle.h"
0029 #include "DataFormats/DetId/interface/DetId.h"
0030 #include "DataFormats/ParticleFlowReco/interface/PFRecHitHostCollection.h"
0031 #include "DataFormats/ParticleFlowReco/interface/PFClusterHostCollection.h"
0032 #include "DataFormats/ParticleFlowReco/interface/PFRecHitFractionHostCollection.h"
0033 #include "HeterogeneousCore/CUDACore/interface/JobConfigurationGPURecord.h"
0034 #include "RecoParticleFlow/PFClusterProducer/interface/PFCPositionCalculatorBase.h"
0035 
0036 class LegacyPFClusterProducer : public edm::stream::EDProducer<> {
0037 public:
0038   LegacyPFClusterProducer(edm::ParameterSet const& config)
0039       : pfClusterSoAToken_(consumes(config.getParameter<edm::InputTag>("src"))),
0040         pfRecHitFractionSoAToken_(consumes(config.getParameter<edm::InputTag>("src"))),
0041         InputPFRecHitSoA_Token_{consumes(config.getParameter<edm::InputTag>("PFRecHitsLabelIn"))},
0042         legacyPfClustersToken_(produces()),
0043         recHitsLabel_(consumes(config.getParameter<edm::InputTag>("recHitsSource"))),
0044         hcalCutsToken_(esConsumes<HcalPFCuts, HcalPFCutsRcd>(edm::ESInputTag("", "withTopo"))),
0045         cutsFromDB_(config.getParameter<bool>("usePFThresholdsFromDB")) {
0046     edm::ConsumesCollector cc = consumesCollector();
0047 
0048     //setup pf cluster builder if requested
0049     const edm::ParameterSet& pfcConf = config.getParameterSet("pfClusterBuilder");
0050     if (!pfcConf.empty()) {
0051       if (pfcConf.exists("positionCalc")) {
0052         const edm::ParameterSet& acConf = pfcConf.getParameterSet("positionCalc");
0053         const std::string& algoac = acConf.getParameter<std::string>("algoName");
0054         positionCalc_ = PFCPositionCalculatorFactory::get()->create(algoac, acConf, cc);
0055       }
0056 
0057       if (pfcConf.exists("allCellsPositionCalc")) {
0058         const edm::ParameterSet& acConf = pfcConf.getParameterSet("allCellsPositionCalc");
0059         const std::string& algoac = acConf.getParameter<std::string>("algoName");
0060         allCellsPositionCalc_ = PFCPositionCalculatorFactory::get()->create(algoac, acConf, cc);
0061       }
0062     }
0063   }
0064 
0065   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0066     edm::ParameterSetDescription desc;
0067     desc.add<edm::InputTag>("src", edm::InputTag("pfClusterSoAProducer"));
0068     desc.add<edm::InputTag>("PFRecHitsLabelIn", edm::InputTag("pfRecHitSoAProducerHCAL"));
0069     desc.add<edm::InputTag>("recHitsSource", edm::InputTag("legacyPFRecHitProducer"));
0070     desc.add<bool>("usePFThresholdsFromDB", true);
0071     {
0072       edm::ParameterSetDescription pfClusterBuilder;
0073       pfClusterBuilder.add<unsigned int>("maxIterations", 5);
0074       pfClusterBuilder.add<double>("minFracTot", 1e-20);
0075       pfClusterBuilder.add<double>("minFractionToKeep", 1e-7);
0076       pfClusterBuilder.add<bool>("excludeOtherSeeds", true);
0077       pfClusterBuilder.add<double>("showerSigma", 10.);
0078       pfClusterBuilder.add<double>("stoppingTolerance", 1e-8);
0079       pfClusterBuilder.add<double>("timeSigmaEB", 10.);
0080       pfClusterBuilder.add<double>("timeSigmaEE", 10.);
0081       pfClusterBuilder.add<double>("maxNSigmaTime", 10.);
0082       pfClusterBuilder.add<double>("minChi2Prob", 0.);
0083       pfClusterBuilder.add<bool>("clusterTimeResFromSeed", false);
0084       pfClusterBuilder.add<std::string>("algoName", "");
0085       {
0086         edm::ParameterSetDescription validator;
0087         validator.add<std::string>("detector", "");
0088         validator.add<std::vector<int>>("depths", {});
0089         validator.add<std::vector<double>>("recHitEnergyNorm", {});
0090         std::vector<edm::ParameterSet> vDefaults(2);
0091         vDefaults[0].addParameter<std::string>("detector", "HCAL_BARREL1");
0092         vDefaults[0].addParameter<std::vector<int>>("depths", {1, 2, 3, 4});
0093         vDefaults[0].addParameter<std::vector<double>>("recHitEnergyNorm", {0.1, 0.2, 0.3, 0.3});
0094         vDefaults[1].addParameter<std::string>("detector", "HCAL_ENDCAP");
0095         vDefaults[1].addParameter<std::vector<int>>("depths", {1, 2, 3, 4, 5, 6, 7});
0096         vDefaults[1].addParameter<std::vector<double>>("recHitEnergyNorm", {0.1, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2});
0097         pfClusterBuilder.addVPSet("recHitEnergyNorms", validator, vDefaults);
0098       }
0099       {
0100         edm::ParameterSetDescription bar;
0101         bar.add<std::string>("algoName", "Basic2DGenericPFlowPositionCalc");
0102         bar.add<double>("minFractionInCalc", 1e-9);
0103         bar.add<int>("posCalcNCrystals", 5);
0104         {
0105           edm::ParameterSetDescription validator;
0106           validator.add<std::string>("detector", "");
0107           validator.add<std::vector<int>>("depths", {});
0108           validator.add<std::vector<double>>("logWeightDenominator", {});
0109           std::vector<edm::ParameterSet> vDefaults(2);
0110           vDefaults[0].addParameter<std::string>("detector", "HCAL_BARREL1");
0111           vDefaults[0].addParameter<std::vector<int>>("depths", {1, 2, 3, 4});
0112           vDefaults[0].addParameter<std::vector<double>>("logWeightDenominator", {0.1, 0.2, 0.3, 0.3});
0113           vDefaults[1].addParameter<std::string>("detector", "HCAL_ENDCAP");
0114           vDefaults[1].addParameter<std::vector<int>>("depths", {1, 2, 3, 4, 5, 6, 7});
0115           vDefaults[1].addParameter<std::vector<double>>("logWeightDenominator", {0.1, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2});
0116           bar.addVPSet("logWeightDenominatorByDetector", validator, vDefaults);
0117         }
0118         bar.add<double>("minAllowedNormalization", 1e-9);
0119         pfClusterBuilder.add("positionCalc", bar);
0120       }
0121       {
0122         edm::ParameterSetDescription bar;
0123         bar.add<std::string>("algoName", "Basic2DGenericPFlowPositionCalc");
0124         bar.add<double>("minFractionInCalc", 1e-9);
0125         bar.add<int>("posCalcNCrystals", -1);
0126         {
0127           edm::ParameterSetDescription validator;
0128           validator.add<std::string>("detector", "");
0129           validator.add<std::vector<int>>("depths", {});
0130           validator.add<std::vector<double>>("logWeightDenominator", {});
0131           std::vector<edm::ParameterSet> vDefaults(2);
0132           vDefaults[0].addParameter<std::string>("detector", "HCAL_BARREL1");
0133           vDefaults[0].addParameter<std::vector<int>>("depths", {1, 2, 3, 4});
0134           vDefaults[0].addParameter<std::vector<double>>("logWeightDenominator", {0.1, 0.2, 0.3, 0.3});
0135           vDefaults[1].addParameter<std::string>("detector", "HCAL_ENDCAP");
0136           vDefaults[1].addParameter<std::vector<int>>("depths", {1, 2, 3, 4, 5, 6, 7});
0137           vDefaults[1].addParameter<std::vector<double>>("logWeightDenominator", {0.1, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2});
0138           bar.addVPSet("logWeightDenominatorByDetector", validator, vDefaults);
0139         }
0140         bar.add<double>("minAllowedNormalization", 1e-9);
0141         pfClusterBuilder.add("allCellsPositionCalc", bar);
0142       }
0143       {
0144         edm::ParameterSetDescription bar;
0145         bar.add<double>("corrTermLowE", 0.);
0146         bar.add<double>("threshLowE", 6.);
0147         bar.add<double>("noiseTerm", 21.86);
0148         bar.add<double>("constantTermLowE", 4.24);
0149         bar.add<double>("noiseTermLowE", 8.);
0150         bar.add<double>("threshHighE", 15.);
0151         bar.add<double>("constantTerm", 2.82);
0152         pfClusterBuilder.add("timeResolutionCalcBarrel", bar);
0153       }
0154       {
0155         edm::ParameterSetDescription bar;
0156         bar.add<double>("corrTermLowE", 0.);
0157         bar.add<double>("threshLowE", 6.);
0158         bar.add<double>("noiseTerm", 21.86);
0159         bar.add<double>("constantTermLowE", 4.24);
0160         bar.add<double>("noiseTermLowE", 8.);
0161         bar.add<double>("threshHighE", 15.);
0162         bar.add<double>("constantTerm", 2.82);
0163         pfClusterBuilder.add("timeResolutionCalcEndcap", bar);
0164       }
0165       {
0166         edm::ParameterSetDescription bar;
0167         pfClusterBuilder.add("positionReCalc", bar);
0168       }
0169       {
0170         edm::ParameterSetDescription bar;
0171         pfClusterBuilder.add("energyCorrector", bar);
0172       }
0173       desc.add("pfClusterBuilder", pfClusterBuilder);
0174     }
0175     descriptions.addWithDefaultLabel(desc);
0176   }
0177 
0178 private:
0179   void produce(edm::Event&, const edm::EventSetup&) override;
0180   const edm::EDGetTokenT<reco::PFClusterHostCollection> pfClusterSoAToken_;
0181   const edm::EDGetTokenT<reco::PFRecHitFractionHostCollection> pfRecHitFractionSoAToken_;
0182   const edm::EDGetTokenT<reco::PFRecHitHostCollection> InputPFRecHitSoA_Token_;
0183   const edm::EDPutTokenT<reco::PFClusterCollection> legacyPfClustersToken_;
0184   const edm::EDGetTokenT<reco::PFRecHitCollection> recHitsLabel_;
0185   const edm::ESGetToken<HcalPFCuts, HcalPFCutsRcd> hcalCutsToken_;
0186   const bool cutsFromDB_;
0187   // the actual algorithm
0188   std::unique_ptr<PFCPositionCalculatorBase> positionCalc_;
0189   std::unique_ptr<PFCPositionCalculatorBase> allCellsPositionCalc_;
0190 };
0191 
0192 void LegacyPFClusterProducer::produce(edm::Event& event, const edm::EventSetup& setup) {
0193   const reco::PFRecHitHostCollection& pfRecHits = event.get(InputPFRecHitSoA_Token_);
0194 
0195   HcalPFCuts const* paramPF = cutsFromDB_ ? &setup.getData(hcalCutsToken_) : nullptr;
0196 
0197   auto const& pfClusterSoA = event.get(pfClusterSoAToken_).const_view();
0198   auto const& pfRecHitFractionSoA = event.get(pfRecHitFractionSoAToken_).const_view();
0199 
0200   int nRH = 0;
0201   if (pfRecHits->metadata().size() != 0)
0202     nRH = pfRecHits.view().size();
0203   reco::PFClusterCollection out;
0204   out.reserve(nRH);
0205 
0206   auto const rechitsHandle = event.getHandle(recHitsLabel_);
0207 
0208   if (nRH != 0) {
0209     // Build PFClusters in legacy format
0210     std::vector<int> nTopoSeeds(nRH, 0);
0211 
0212     for (int i = 0; i < pfClusterSoA.nSeeds(); i++) {
0213       nTopoSeeds[pfClusterSoA[i].topoId()]++;
0214     }
0215 
0216     // Looping over SoA PFClusters to produce legacy PFCluster collection
0217     for (int i = 0; i < pfClusterSoA.nSeeds(); i++) {
0218       unsigned int n = pfClusterSoA[i].seedRHIdx();
0219       reco::PFCluster temp;
0220       temp.setSeed((*rechitsHandle)[n].detId());  // Pulling the detId of this PFRecHit from the legacy format input
0221       int offset = pfClusterSoA[i].rhfracOffset();
0222       for (int k = offset; k < (offset + pfClusterSoA[i].rhfracSize()) && k >= 0;
0223            k++) {  // Looping over PFRecHits in the same topo cluster
0224         if (pfRecHitFractionSoA[k].pfrhIdx() < nRH && pfRecHitFractionSoA[k].pfrhIdx() > -1 &&
0225             pfRecHitFractionSoA[k].frac() > 0.0) {
0226           const reco::PFRecHitRef& refhit = reco::PFRecHitRef(rechitsHandle, pfRecHitFractionSoA[k].pfrhIdx());
0227           temp.addRecHitFraction(reco::PFRecHitFraction(refhit, pfRecHitFractionSoA[k].frac()));
0228         }
0229       }
0230 
0231       // Now PFRecHitFraction of this PFCluster is set. Now compute calculateAndSetPosition (energy, position etc)
0232       if (nTopoSeeds[pfClusterSoA[i].topoId()] == 1 && allCellsPositionCalc_) {
0233         allCellsPositionCalc_->calculateAndSetPosition(temp, paramPF);
0234       } else {
0235         positionCalc_->calculateAndSetPosition(temp, paramPF);
0236       }
0237       out.emplace_back(std::move(temp));
0238     }
0239   }
0240 
0241   event.emplace(legacyPfClustersToken_, std::move(out));
0242 }
0243 
0244 #include "FWCore/Framework/interface/MakerMacros.h"
0245 DEFINE_FWK_MODULE(LegacyPFClusterProducer);