File indexing completed on 2022-05-21 03:43:06
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020
0021
0022
0023
0024
0025
0026
0027 #ifndef RecoJets_FFTJetProducers_FFTJetProducer_h
0028 #define RecoJets_FFTJetProducers_FFTJetProducer_h
0029
0030 #include <memory>
0031
0032
0033 #include "fftjet/AbsRecombinationAlg.hh"
0034 #include "fftjet/AbsVectorRecombinationAlg.hh"
0035 #include "fftjet/SparseClusteringTree.hh"
0036
0037
0038 #include "fftjet/VectorRecombinationAlgFactory.hh"
0039 #include "fftjet/RecombinationAlgFactory.hh"
0040
0041 #include "DataFormats/JetReco/interface/FFTCaloJetCollection.h"
0042 #include "DataFormats/JetReco/interface/FFTGenJetCollection.h"
0043 #include "DataFormats/JetReco/interface/FFTPFJetCollection.h"
0044 #include "DataFormats/JetReco/interface/FFTJPTJetCollection.h"
0045 #include "DataFormats/JetReco/interface/FFTBasicJetCollection.h"
0046 #include "DataFormats/JetReco/interface/FFTTrackJetCollection.h"
0047 #include "DataFormats/JetReco/interface/FFTJetProducerSummary.h"
0048 #include "RecoJets/FFTJetProducers/interface/FFTJetParameterParser.h"
0049
0050 #include "RecoJets/FFTJetAlgorithms/interface/clusteringTreeConverters.h"
0051 #include "RecoJets/FFTJetAlgorithms/interface/jetConverters.h"
0052 #include "RecoJets/FFTJetAlgorithms/interface/matchOneToOne.h"
0053 #include "RecoJets/FFTJetAlgorithms/interface/JetToPeakDistance.h"
0054 #include "RecoJets/FFTJetAlgorithms/interface/adjustForPileup.h"
0055
0056 #include "DataFormats/JetReco/interface/DiscretizedEnergyFlow.h"
0057
0058
0059 #include "FWCore/Framework/interface/Frameworkfwd.h"
0060 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0061 #include "FWCore/Framework/interface/Event.h"
0062
0063 #include "DataFormats/Candidate/interface/CandidateFwd.h"
0064
0065
0066 #include "RecoJets/FFTJetAlgorithms/interface/fftjetTypedefs.h"
0067 #include "RecoJets/FFTJetAlgorithms/interface/AbsPileupCalculator.h"
0068 #include "RecoJets/FFTJetProducers/interface/FFTJetInterface.h"
0069
0070 namespace fftjetcms {
0071 class DiscretizedEnergyFlow;
0072 }
0073
0074 class CaloGeometry;
0075 class CaloGeometryRecord;
0076 class HcalTopology;
0077 class HcalRecNumberingRecord;
0078
0079
0080
0081
0082 class FFTJetProducer : public fftjetcms::FFTJetInterface {
0083 public:
0084 typedef fftjet::RecombinedJet<fftjetcms::VectorLike> RecoFFTJet;
0085 typedef fftjet::SparseClusteringTree<fftjet::Peak, long> SparseTree;
0086
0087
0088
0089
0090 enum StatusBits {
0091 RESOLUTION = 0xff,
0092 CONSTITUENTS_RESUMMED = 0x100,
0093 PILEUP_CALCULATED = 0x200,
0094 PILEUP_SUBTRACTED_4VEC = 0x400,
0095 PILEUP_SUBTRACTED_PT = 0x800
0096 };
0097
0098 enum Resolution { FIXED = 0, MAXIMALLY_STABLE, GLOBALLY_ADAPTIVE, LOCALLY_ADAPTIVE, FROM_GENJETS };
0099
0100 explicit FFTJetProducer(const edm::ParameterSet&);
0101
0102 FFTJetProducer() = delete;
0103 FFTJetProducer(const FFTJetProducer&) = delete;
0104 FFTJetProducer& operator=(const FFTJetProducer&) = delete;
0105
0106 ~FFTJetProducer() override;
0107
0108
0109 static Resolution parse_resolution(const std::string& name);
0110
0111 protected:
0112
0113 void beginStream(edm::StreamID) override;
0114 void produce(edm::Event&, const edm::EventSetup&) override;
0115
0116
0117
0118
0119
0120
0121 virtual void selectPreclusters(const SparseTree& tree,
0122 const fftjet::Functor1<bool, fftjet::Peak>& peakSelector,
0123 std::vector<fftjet::Peak>* preclusters);
0124
0125
0126 virtual void genJetPreclusters(const SparseTree& tree,
0127 edm::Event&,
0128 const edm::EventSetup&,
0129 const fftjet::Functor1<bool, fftjet::Peak>& peakSelector,
0130 std::vector<fftjet::Peak>* preclusters);
0131
0132
0133
0134
0135
0136 virtual void assignMembershipFunctions(std::vector<fftjet::Peak>* preclusters);
0137
0138
0139 virtual std::unique_ptr<fftjet::Functor1<bool, fftjet::Peak> > parse_peakSelector(const edm::ParameterSet&);
0140
0141
0142 virtual std::unique_ptr<fftjet::ScaleSpaceKernel> parse_jetMembershipFunction(const edm::ParameterSet&);
0143
0144
0145 virtual std::unique_ptr<fftjetcms::AbsBgFunctor> parse_bgMembershipFunction(const edm::ParameterSet&);
0146
0147
0148 virtual std::unique_ptr<fftjet::Functor1<double, fftjet::Peak> > parse_recoScaleCalcPeak(const edm::ParameterSet&);
0149
0150
0151 virtual std::unique_ptr<fftjet::Functor1<double, fftjet::Peak> > parse_recoScaleRatioCalcPeak(
0152 const edm::ParameterSet&);
0153
0154
0155 virtual std::unique_ptr<fftjet::Functor1<double, fftjet::Peak> > parse_memberFactorCalcPeak(const edm::ParameterSet&);
0156
0157
0158 virtual std::unique_ptr<fftjet::Functor1<double, RecoFFTJet> > parse_recoScaleCalcJet(const edm::ParameterSet&);
0159
0160 virtual std::unique_ptr<fftjet::Functor1<double, RecoFFTJet> > parse_recoScaleRatioCalcJet(const edm::ParameterSet&);
0161
0162 virtual std::unique_ptr<fftjet::Functor1<double, RecoFFTJet> > parse_memberFactorCalcJet(const edm::ParameterSet&);
0163
0164
0165
0166 virtual std::unique_ptr<fftjet::Functor2<double, RecoFFTJet, RecoFFTJet> > parse_jetDistanceCalc(
0167 const edm::ParameterSet&);
0168
0169
0170 virtual std::unique_ptr<fftjetcms::AbsPileupCalculator> parse_pileupDensityCalc(const edm::ParameterSet& ps);
0171
0172
0173
0174
0175 void selectTreeNodes(const SparseTree& tree,
0176 const fftjet::Functor1<bool, fftjet::Peak>& peakSelect,
0177 std::vector<SparseTree::NodeId>* nodes);
0178
0179 private:
0180 typedef fftjet::AbsVectorRecombinationAlg<fftjetcms::VectorLike, fftjetcms::BgData> RecoAlg;
0181 typedef fftjet::AbsRecombinationAlg<fftjetcms::Real, fftjetcms::VectorLike, fftjetcms::BgData> GridAlg;
0182
0183
0184 template <class Real>
0185 void loadSparseTreeData(const edm::Event&,
0186 const edm::EDGetTokenT<reco::PattRecoTree<Real, reco::PattRecoPeak<Real> > >& tok);
0187
0188 void removeFakePreclusters();
0189
0190
0191
0192 bool loadEnergyFlow(const edm::Event& iEvent, std::unique_ptr<fftjet::Grid2d<fftjetcms::Real> >& flow);
0193 void buildGridAlg();
0194 void prepareRecombinationScales();
0195 bool checkConvergence(const std::vector<RecoFFTJet>& previousIterResult, std::vector<RecoFFTJet>& thisIterResult);
0196 void determineGriddedConstituents();
0197 void determineVectorConstituents();
0198 void saveResults(edm::Event& iEvent, const edm::EventSetup&, unsigned nPreclustersFound);
0199
0200 template <typename Jet>
0201 void writeJets(edm::Event& iEvent, const edm::EventSetup&);
0202
0203 template <typename Jet>
0204 void makeProduces(const std::string& alias, const std::string& tag);
0205
0206
0207
0208
0209 virtual void determinePileupDensityFromConfig(const edm::Event& iEvent,
0210 std::unique_ptr<fftjet::Grid2d<fftjetcms::Real> >& density);
0211
0212
0213 virtual void determinePileupDensityFromDB(const edm::Event& iEvent,
0214 const edm::EventSetup& iSetup,
0215 std::unique_ptr<fftjet::Grid2d<fftjetcms::Real> >& density);
0216
0217
0218
0219 void determinePileup();
0220
0221
0222
0223
0224
0225
0226
0227 unsigned iterateJetReconstruction();
0228
0229
0230 static void setJetStatusBit(RecoFFTJet* jet, int mask, bool value);
0231
0232
0233
0234
0235
0236
0237 const edm::ParameterSet myConfiguration;
0238
0239
0240 const edm::InputTag treeLabel;
0241
0242
0243
0244 const bool useGriddedAlgorithm;
0245
0246
0247
0248 const bool reuseExistingGrid;
0249
0250
0251 const unsigned maxIterations;
0252
0253
0254 const unsigned nJetsRequiredToConverge;
0255 const double convergenceDistance;
0256
0257
0258 const bool assignConstituents;
0259
0260
0261
0262
0263 const bool resumConstituents;
0264
0265
0266
0267 const bool calculatePileup;
0268 const bool subtractPileup;
0269 const bool subtractPileupAs4Vec;
0270
0271
0272
0273 const edm::InputTag pileupLabel;
0274
0275
0276 const double fixedScale;
0277
0278
0279 const double minStableScale;
0280 const double maxStableScale;
0281
0282
0283 const double stabilityAlpha;
0284
0285
0286
0287 const double noiseLevel;
0288
0289
0290 const unsigned nClustersRequested;
0291
0292
0293 const double gridScanMaxEta;
0294
0295
0296 const std::string recombinationAlgorithm;
0297 const bool isCrisp;
0298 const double unlikelyBgWeight;
0299 const double recombinationDataCutoff;
0300
0301
0302 const edm::InputTag genJetsLabel;
0303
0304
0305
0306
0307 const unsigned maxInitialPreclusters;
0308
0309
0310
0311
0312 Resolution resolution;
0313
0314
0315
0316 std::string pileupTableRecord;
0317 std::string pileupTableName;
0318 std::string pileupTableCategory;
0319 bool loadPileupFromDB;
0320
0321
0322 std::unique_ptr<std::vector<double> > iniScales;
0323
0324
0325 SparseTree sparseTree;
0326
0327
0328 std::unique_ptr<fftjet::Functor1<bool, fftjet::Peak> > peakSelector;
0329
0330
0331 std::unique_ptr<RecoAlg> recoAlg;
0332 std::unique_ptr<GridAlg> gridAlg;
0333 std::unique_ptr<fftjet::ScaleSpaceKernel> jetMembershipFunction;
0334 std::unique_ptr<fftjetcms::AbsBgFunctor> bgMembershipFunction;
0335
0336
0337 std::unique_ptr<fftjet::Functor1<double, fftjet::Peak> > recoScaleCalcPeak;
0338
0339
0340 std::unique_ptr<fftjet::Functor1<double, fftjet::Peak> > recoScaleRatioCalcPeak;
0341
0342
0343 std::unique_ptr<fftjet::Functor1<double, fftjet::Peak> > memberFactorCalcPeak;
0344
0345
0346 std::unique_ptr<fftjet::Functor1<double, RecoFFTJet> > recoScaleCalcJet;
0347 std::unique_ptr<fftjet::Functor1<double, RecoFFTJet> > recoScaleRatioCalcJet;
0348 std::unique_ptr<fftjet::Functor1<double, RecoFFTJet> > memberFactorCalcJet;
0349
0350
0351
0352 std::unique_ptr<fftjet::Functor2<double, RecoFFTJet, RecoFFTJet> > jetDistanceCalc;
0353
0354
0355 std::vector<SparseTree::NodeId> nodes;
0356
0357
0358 std::vector<fftjet::Peak> preclusters;
0359
0360
0361 std::vector<RecoFFTJet> recoJets;
0362
0363
0364 std::vector<unsigned> occupancy;
0365
0366
0367 std::vector<double> thresholds;
0368
0369
0370 unsigned minLevel, maxLevel, usedLevel;
0371
0372
0373 fftjetcms::VectorLike unclustered;
0374 double unused;
0375
0376
0377 std::vector<fftjet::Peak> iterPreclusters;
0378 std::vector<RecoFFTJet> iterJets;
0379 unsigned iterationsPerformed;
0380
0381
0382 std::vector<std::vector<reco::CandidatePtr> > constituents;
0383
0384
0385
0386 std::vector<fftjetcms::VectorLike> pileup;
0387
0388
0389
0390
0391 std::unique_ptr<fftjet::Grid2d<fftjetcms::Real> > pileupEnergyFlow;
0392
0393
0394 std::unique_ptr<fftjetcms::AbsPileupCalculator> pileupDensityCalc;
0395
0396
0397 std::vector<fftjet::AbsKernel2d*> memFcns2dVec;
0398 std::vector<double> doubleBuf;
0399 std::vector<unsigned> cellCountsVec;
0400
0401
0402 edm::EDGetTokenT<reco::PattRecoTree<double, reco::PattRecoPeak<double> > > input_recotree_token_d_;
0403 edm::EDGetTokenT<reco::PattRecoTree<float, reco::PattRecoPeak<float> > > input_recotree_token_f_;
0404 edm::EDGetTokenT<std::vector<reco::FFTAnyJet<reco::GenJet> > > input_genjet_token_;
0405 edm::EDGetTokenT<reco::DiscretizedEnergyFlow> input_energyflow_token_;
0406 edm::EDGetTokenT<reco::FFTJetPileupSummary> input_pusummary_token_;
0407
0408 edm::ESGetToken<CaloGeometry, CaloGeometryRecord> geometry_token_;
0409 edm::ESGetToken<HcalTopology, HcalRecNumberingRecord> topology_token_;
0410 };
0411
0412 #endif