FFTJetProducers/src/FFTJetInterface.cc

0001 #include <cassert>
0002
0003 #include "RecoJets/FFTJetProducers/interface/FFTJetInterface.h"
0004
0005 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0006 #include "FWCore/Utilities/interface/isFinite.h"
0007
0008 #include "DataFormats/JetReco/interface/Jet.h"
0009 #include "DataFormats/CaloTowers/interface/CaloTower.h"
0010
0011 #define init_param(type, varname) varname(ps.getParameter<type>(#varname))
0012
0013 namespace fftjetcms {
0014
0015   bool FFTJetInterface::storeInSinglePrecision() const { return true; }
0016
0017   FFTJetInterface::FFTJetInterface(const edm::ParameterSet& ps)
0018       : inputLabel(ps.getParameter<edm::InputTag>("src")),
0019         init_param(std::string, outputLabel),
0020         jetType(parseJetType(ps.getParameter<std::string>("jetType"))),
0021         init_param(bool, doPVCorrection),
0022         init_param(edm::InputTag, srcPVs),
0023         etaDependentMagnutideFactors(ps.getParameter<std::vector<double> >("etaDependentMagnutideFactors")),
0024         init_param(edm::ParameterSet, anomalous),
0025         init_param(bool, insertCompleteEvent),
0026         init_param(double, completeEventScale),
0027         vertex_(0.0, 0.0, 0.0) {
0028     if (insertCompleteEvent && completeEventScale <= 0.0)
0029       throw cms::Exception("FFTJetBadConfig") << "Bad scale for the complete event : must be positive" << std::endl;
0030
0031     inputToken = consumes<reco::CandidateView>(inputLabel);
0032
0033     if (doPVCorrection && jetType == CALOJET)
0034       srcPVsToken = consumes<reco::VertexCollection>(srcPVs);
0035   }
0036
0037   double FFTJetInterface::getEventScale() const { return insertCompleteEvent ? completeEventScale : 0.0; }
0038
0039   void FFTJetInterface::loadInputCollection(const edm::Event& iEvent) {
0040     // Figure out if we are going to perform the vertex adjustment
0041     const bool adjustForVertex = doPVCorrection && jetType == CALOJET;
0042
0043     // Figure out the vertex
0044     if (adjustForVertex) {
0045       edm::Handle<reco::VertexCollection> pvCollection;
0046       iEvent.getByToken(srcPVsToken, pvCollection);
0047       if (pvCollection->empty())
0048         vertex_ = reco::Particle::Point(0.0, 0.0, 0.0);
0049       else
0050         vertex_ = pvCollection->begin()->position();
0051     }
0052
0053     // Get the input collection
0054     iEvent.getByToken(inputToken, inputCollection);
0055
0056     // Create the set of 4-vectors needed by the algorithm
0057     eventData.clear();
0058     candidateIndex.clear();
0059     unsigned index = 0;
0060     const reco::CandidateView::const_iterator end = inputCollection->end();
0061     for (reco::CandidateView::const_iterator it = inputCollection->begin(); it != end; ++it, ++index) {
0062       const reco::Candidate& item(*it);
0063       if (anomalous(item))
0064         continue;
0065       if (edm::isNotFinite(item.pt()))
0066         continue;
0067
0068       if (adjustForVertex) {
0069         const CaloTower& tower(dynamic_cast<const CaloTower&>(item));
0070         eventData.push_back(VectorLike(tower.p4(vertex_)));
0071       } else
0072         eventData.push_back(item.p4());
0073       candidateIndex.push_back(index);
0074     }
0075     assert(eventData.size() == candidateIndex.size());
0076   }
0077
0078   void FFTJetInterface::discretizeEnergyFlow() {
0079     // It is a bug to call this function before defining the energy flow grid
0080     assert(energyFlow.get());
0081
0082     fftjet::Grid2d<Real>& g(*energyFlow);
0083     g.reset();
0084
0085     const unsigned nInputs = eventData.size();
0086     const VectorLike* inp = nInputs ? &eventData[0] : nullptr;
0087     for (unsigned i = 0; i < nInputs; ++i) {
0088       const VectorLike& item(inp[i]);
0089       g.fill(item.Eta(), item.Phi(), item.Et());
0090     }
0091
0092     if (!etaDependentMagnutideFactors.empty()) {
0093       if (etaDependentMagnutideFactors.size() != g.nEta())
0094         throw cms::Exception("FFTJetBadConfig") << "ERROR in FFTJetInterface::discretizeEnergyFlow() :"
0095                                                    " number of elements in the \"etaDependentMagnutideFactors\""
0096                                                    " vector is inconsistent with the grid binning"
0097                                                 << std::endl;
0098       g.scaleData(&etaDependentMagnutideFactors[0], etaDependentMagnutideFactors.size());
0099     }
0100   }
0101
0102 }  // namespace fftjetcms