#include <RecoJets/FFTJetProducer/plugins/FFTJetPatRecoProducer.cc>

Inheritance diagram for FFTJetPatRecoProducer:

Public Member Functions
	FFTJetPatRecoProducer (const edm::ParameterSet &)

	~FFTJetPatRecoProducer ()

Public Member Functions inherited from fftjetcms::FFTJetInterface
virtual	~FFTJetInterface ()

Public Member Functions inherited from edm::EDProducer
	EDProducer ()

ModuleDescription const &	moduleDescription () const

virtual	~EDProducer ()

Public Member Functions inherited from edm::ProducerBase
	ProducerBase ()

void	registerProducts (ProducerBase , ProductRegistry , ModuleDescription const &)

std::function< void(BranchDescription const &)>	registrationCallback () const
	used by the fwk to register list of products More...

virtual	~ProducerBase ()

Public Member Functions inherited from edm::EDConsumerBase
std::vector< ConsumesInfo >	consumesInfo () const

	EDConsumerBase ()

ProductHolderIndexAndSkipBit	indexFrom (EDGetToken, BranchType, TypeID const &) const

void	itemsMayGet (BranchType, std::vector< ProductHolderIndexAndSkipBit > &) const

void	itemsToGet (BranchType, std::vector< ProductHolderIndexAndSkipBit > &) const

std::vector < ProductHolderIndexAndSkipBit > const &	itemsToGetFromEvent () const

void	labelsForToken (EDGetToken iToken, Labels &oLabels) const

void	modulesDependentUpon (const std::string &iProcessName, std::vector< const char * > &oModuleLabels) const

void	modulesWhoseProductsAreConsumed (std::vector< ModuleDescription const * > &modules, ProductRegistry const &preg, std::map< std::string, ModuleDescription const * > const &labelsToDesc, std::string const &processName) const

bool	registeredToConsume (ProductHolderIndex, bool, BranchType) const

bool	registeredToConsumeMany (TypeID const &, BranchType) const

void	updateLookup (BranchType iBranchType, ProductHolderIndexHelper const &)

virtual	~EDConsumerBase ()

Protected Types
typedef fftjet::ProximityClusteringTree < fftjet::Peak, long >	ClusteringTree

typedef fftjet::ClusteringSequencer < Real >	Sequencer

typedef fftjet::SparseClusteringTree < fftjet::Peak, long >	SparseTree

typedef fftjet::ClusteringTreeSparsifier < fftjet::Peak, long >	Sparsifier

Protected Member Functions
void	beginJob () override

template<class Real >
void	buildDenseProduct (edm::Event &) const

void	buildKernelConvolver (const edm::ParameterSet &)

fftjet::PeakFinder	buildPeakFinder (const edm::ParameterSet &)

template<class Real >
void	buildSparseProduct (edm::Event &) const

void	endJob () override

void	produce (edm::Event &, const edm::EventSetup &) override

Protected Member Functions inherited from fftjetcms::FFTJetInterface
template<class Ptr >
void	checkConfig (const Ptr &ptr, const char *message)

void	discretizeEnergyFlow ()

	FFTJetInterface (const edm::ParameterSet &)

double	getEventScale () const

void	loadInputCollection (const edm::Event &)

bool	storeInSinglePrecision () const

const reco::Particle::Point &	vertexUsed () const

Protected Member Functions inherited from edm::ProducerBase
void	callWhenNewProductsRegistered (std::function< void(BranchDescription const &)> const &func)

Protected Member Functions inherited from edm::EDConsumerBase
template<typename ProductType , BranchType B = InEvent>
EDGetTokenT< ProductType >	consumes (edm::InputTag const &tag)

EDGetToken	consumes (const TypeToGet &id, edm::InputTag const &tag)

template<BranchType B>
EDGetToken	consumes (TypeToGet const &id, edm::InputTag const &tag)

ConsumesCollector	consumesCollector ()
	Use a ConsumesCollector to gather consumes information from helper functions. More...

template<typename ProductType , BranchType B = InEvent>
void	consumesMany ()

void	consumesMany (const TypeToGet &id)

template<BranchType B>
void	consumesMany (const TypeToGet &id)

template<typename ProductType , BranchType B = InEvent>
EDGetTokenT< ProductType >	mayConsume (edm::InputTag const &tag)

EDGetToken	mayConsume (const TypeToGet &id, edm::InputTag const &tag)

template<BranchType B>
EDGetToken	mayConsume (const TypeToGet &id, edm::InputTag const &tag)

Protected Attributes
std::auto_ptr< MyFFTEngine >	anotherEngine

ClusteringTree *	clusteringTree

const double	completeEventDataCutoff

std::auto_ptr < fftjet::AbsConvolverBase < Real > >	convolver

std::auto_ptr < fftjet::AbsDistanceCalculator < fftjet::Peak > >	distanceCalc

std::auto_ptr< MyFFTEngine >	engine

std::auto_ptr < fftjet::AbsFrequencyKernel1d >	etaKernel

std::auto_ptr < fftjet::AbsFrequencyKernel >	kernel2d

const bool	makeClusteringTree

std::auto_ptr < fftjet::Functor1< bool, fftjet::Peak > >	peakSelector

std::auto_ptr < fftjet::AbsFrequencyKernel1d >	phiKernel

std::auto_ptr< Sequencer >	sequencer

SparseTree	sparseTree

std::auto_ptr< Sparsifier >	sparsifier

const bool	sparsify

const bool	storeDiscretizationGrid

const bool	verifyDataConversion

Protected Attributes inherited from fftjetcms::FFTJetInterface
const AnomalousTower	anomalous

std::vector< unsigned >	candidateIndex

const bool	doPVCorrection

std::auto_ptr< fftjet::Grid2d < fftjetcms::Real > >	energyFlow

const std::vector< double >	etaDependentMagnutideFactors

std::vector < fftjetcms::VectorLike >	eventData

edm::Handle< reco::CandidateView >	inputCollection

const edm::InputTag	inputLabel

const JetType	jetType

const std::string	outputLabel

const edm::InputTag	srcPVs

Private Member Functions
	FFTJetPatRecoProducer ()

	FFTJetPatRecoProducer (const FFTJetPatRecoProducer &)

FFTJetPatRecoProducer &	operator= (const FFTJetPatRecoProducer &)

Private Attributes
std::ofstream	externalGridStream

fftjet::Grid2d< float > *	extGrid

bool	storeGridsExternally

Additional Inherited Members
Public Types inherited from edm::EDProducer
typedef EDProducer	ModuleType

Public Types inherited from edm::ProducerBase
typedef ProductRegistryHelper::TypeLabelList	TypeLabelList

Static Public Member Functions inherited from edm::EDProducer
static const std::string &	baseType ()

static void	fillDescriptions (ConfigurationDescriptions &descriptions)

static void	prevalidate (ConfigurationDescriptions &descriptions)

Detailed Description

Description: Runs FFTJet pattern recognition stage and saves the results

Implementation: [Notes on implementation]

Definition at line 63 of file FFTJetPatRecoProducer.cc.

Member Typedef Documentation

typedef fftjet::ProximityClusteringTree<fftjet::Peak,long> FFTJetPatRecoProducer::ClusteringTree

protected

Definition at line 71 of file FFTJetPatRecoProducer.cc.

typedef fftjet::ClusteringSequencer<Real> FFTJetPatRecoProducer::Sequencer

protected

Definition at line 73 of file FFTJetPatRecoProducer.cc.

typedef fftjet::SparseClusteringTree<fftjet::Peak,long> FFTJetPatRecoProducer::SparseTree

protected

Definition at line 72 of file FFTJetPatRecoProducer.cc.

typedef fftjet::ClusteringTreeSparsifier<fftjet::Peak,long> FFTJetPatRecoProducer::Sparsifier

protected

Definition at line 74 of file FFTJetPatRecoProducer.cc.

Constructor & Destructor Documentation

FFTJetPatRecoProducer::FFTJetPatRecoProducer ( const edm::ParameterSet & ps )

explicit

Definition at line 152 of file FFTJetPatRecoProducer.cc.

References buildKernelConvolver(), buildPeakFinder(), fftjetcms::FFTJetInterface::checkConfig(), clusteringTree, convolver, fwrapper::cs, distanceCalc, fftjetcms::FFTJetInterface::energyFlow, edm::hlt::Exception, externalGridStream, fftjetcms::fftjet_ClusteringTreeSparsifier_parser(), fftjetcms::fftjet_DistanceCalculator_parser(), fftjetcms::fftjet_Grid2d_parser(), fftjetcms::fftjet_PeakSelector_parser(), fftjetcms::fftjet_ScaleSet_parser(), fftjetcms::FFTJetInterface::getEventScale(), edm::ParameterSet::getParameter(), i, makeClusteringTree, dbtoconf::out, fftjetcms::FFTJetInterface::outputLabel, peakSelector, sequencer, sparsifier, storeDiscretizationGrid, storeGridsExternally, fftjetcms::FFTJetInterface::storeInSinglePrecision(), and AlCaHLTBitMon_QueryRunRegistry::string.

     : FFTJetInterface(ps),
       clusteringTree(0),
       completeEventDataCutoff(ps.getParameter<double>("completeEventDataCutoff")),
       makeClusteringTree(ps.getParameter<bool>("makeClusteringTree")),
       verifyDataConversion(ps.getUntrackedParameter<bool>("verifyDataConversion",false)),
       storeDiscretizationGrid(ps.getParameter<bool>("storeDiscretizationGrid")),
       sparsify(ps.getParameter<bool>("sparsify")),
       extGrid(0)
 {
     // register your products
     if (makeClusteringTree)
     {
         if (storeInSinglePrecision())
             produces<reco::PattRecoTree<float,reco::PattRecoPeak<float> > >(outputLabel);
         else
             produces<reco::PattRecoTree<double,reco::PattRecoPeak<double> > >(outputLabel);
     }
     if (storeDiscretizationGrid)
         produces<reco::DiscretizedEnergyFlow>(outputLabel);
 
     // Check if we want to write the grids into an external file
     const std::string externalGridFile(ps.getParameter<std::string>("externalGridFile"));
     storeGridsExternally = externalGridFile.size() > 0;
     if (storeGridsExternally)
     {
         externalGridStream.open(externalGridFile.c_str(), std::ios_base::out | 
                                                           std::ios_base::binary);
         if (!externalGridStream.is_open())
             throw cms::Exception("FFTJetBadConfig")
                 << "FFTJetPatRecoProducer failed to open file "
                 << externalGridFile << std::endl;
     }
 
     if (!makeClusteringTree && !storeDiscretizationGrid && !storeGridsExternally)
     {
         throw cms::Exception("FFTJetBadConfig")
             << "FFTJetPatRecoProducer is not configured to produce anything"
             << std::endl;
     }
 
     // now do whatever other initialization is needed
 
     // Build the discretization grid
     energyFlow = fftjet_Grid2d_parser(
         ps.getParameter<edm::ParameterSet>("GridConfiguration"));
     checkConfig(energyFlow, "invalid discretization grid");
 
     // Build the FFT engine(s), pattern recognition kernel(s),
     // and the kernel convolver
     buildKernelConvolver(ps);
 
     // Build the peak selector
     peakSelector = fftjet_PeakSelector_parser(
         ps.getParameter<edm::ParameterSet>("PeakSelectorConfiguration"));
     checkConfig(peakSelector, "invalid peak selector");
 
     // Build the initial set of pattern recognition scales
     std::auto_ptr<std::vector<double> > iniScales = fftjet_ScaleSet_parser(
         ps.getParameter<edm::ParameterSet>("InitialScales"));
     checkConfig(iniScales, "invalid set of scales");
 
     // Do we want to use the adaptive clustering tree algorithm?
     const unsigned maxAdaptiveScales = 
         ps.getParameter<unsigned>("maxAdaptiveScales");
     const double minAdaptiveRatioLog = 
         ps.getParameter<double>("minAdaptiveRatioLog");
     if (minAdaptiveRatioLog <= 0.0)
         throw cms::Exception("FFTJetBadConfig")
             << "bad adaptive ratio logarithm limit" << std::endl;
 
     // Make sure that all standard scales are larger than the
     // complete event scale
     if (getEventScale() > 0.0)
     {
       const double cs = getEventScale();
       const unsigned nscales = iniScales->size();
       for (unsigned i=0; i<nscales; ++i)
         if (cs >= (*iniScales)[i])
       throw cms::Exception("FFTJetBadConfig")
         << "incompatible scale for complete event" << std::endl;
     }
 
     // At this point we are ready to build the clustering sequencer
     sequencer = std::auto_ptr<Sequencer>(new Sequencer(
         convolver.get(), peakSelector.get(), buildPeakFinder(ps),
         *iniScales, maxAdaptiveScales, minAdaptiveRatioLog));
 
     // Build the clustering tree sparsifier
     const edm::ParameterSet& SparsifierConfiguration(
         ps.getParameter<edm::ParameterSet>("SparsifierConfiguration"));
     sparsifier = fftjet_ClusteringTreeSparsifier_parser(
         SparsifierConfiguration);
     checkConfig(sparsifier, "invalid sparsifier parameters");
 
     // Build the distance calculator for the clustering tree
     const edm::ParameterSet& TreeDistanceCalculator(
         ps.getParameter<edm::ParameterSet>("TreeDistanceCalculator"));
     distanceCalc = fftjet_DistanceCalculator_parser(TreeDistanceCalculator);
     checkConfig(distanceCalc, "invalid tree distance calculator");
 
     // Build the clustering tree itself
     clusteringTree = new ClusteringTree(distanceCalc.get());
 }

FFTJetPatRecoProducer::~FFTJetPatRecoProducer ( )

Definition at line 356 of file FFTJetPatRecoProducer.cc.

References clusteringTree, and extGrid.

 {
     // do anything here that needs to be done at desctruction time
     // (e.g. close files, deallocate resources etc.)
     delete clusteringTree;
     delete extGrid;
 }

FFTJetPatRecoProducer::FFTJetPatRecoProducer ( )

private

FFTJetPatRecoProducer::FFTJetPatRecoProducer ( const FFTJetPatRecoProducer & )

private

Member Function Documentation

void FFTJetPatRecoProducer::beginJob ( void )

overrideprotectedvirtual

Reimplemented from edm::EDProducer.

Definition at line 497 of file FFTJetPatRecoProducer.cc.

498 {

499 }

template<class Real >

void FFTJetPatRecoProducer::buildDenseProduct ( edm::Event & ev ) const

protected

Definition at line 396 of file FFTJetPatRecoProducer.cc.

References DDVectorGetter::check(), clusteringTree, fftjetcms::densePeakTreeFromStorable(), fftjetcms::densePeakTreeToStorable(), distanceCalc, edm::hlt::Exception, fftjetcms::FFTJetInterface::getEventScale(), fftjetcms::FFTJetInterface::outputLabel, edm::Event::put(), sequencer, fftjetcms::FFTJetInterface::storeInSinglePrecision(), MainPageGenerator::tree, and verifyDataConversion.

 {
     typedef reco::PattRecoTree<Real,reco::PattRecoPeak<Real> > StoredTree;
 
     std::auto_ptr<StoredTree> tree(new StoredTree());
 
     densePeakTreeToStorable(*clusteringTree,
                             sequencer->maxAdaptiveScales(),
                             tree.get());
 
     // Check that we can restore the tree
     if (verifyDataConversion && !storeInSinglePrecision())
     {
         ClusteringTree check(distanceCalc.get());
         const std::vector<double>& scalesUsed(sequencer->getInitialScales());
         densePeakTreeFromStorable(*tree, &scalesUsed, getEventScale(), &check);
         if (*clusteringTree != check)
             throw cms::Exception("FFTJetInterface")
                 << "Data conversion failed for dense clustering tree"
                 << std::endl;
     }
 
     ev.put(tree, outputLabel);
 }

void FFTJetPatRecoProducer::buildKernelConvolver ( const edm::ParameterSet & ps )

protected

Definition at line 258 of file FFTJetPatRecoProducer.cc.

References anotherEngine, convolver, fftjetcms::FFTJetInterface::energyFlow, engine, etaKernel, edm::hlt::Exception, edm::ParameterSet::getParameter(), kernel2d, M_PI, and phiKernel.

Referenced by FFTJetPatRecoProducer().

 {
     // Check the parameter named "etaDependentScaleFactors". If the vector
     // of scales is empty we will use 2d kernel, otherwise use 1d kernels
     const std::vector<double> etaDependentScaleFactors(
         ps.getParameter<std::vector<double> >("etaDependentScaleFactors"));
 
     // Make sure that the number of scale factors provided is correct
     const bool use2dKernel = etaDependentScaleFactors.empty();
     if (!use2dKernel)
         if (etaDependentScaleFactors.size() != energyFlow->nEta())
             throw cms::Exception("FFTJetBadConfig")
                 << "invalid number of eta-dependent scale factors"
                 << std::endl;
 
     // Get the eta and phi scales for the kernel(s)
     double kernelEtaScale = ps.getParameter<double>("kernelEtaScale");
     const double kernelPhiScale = ps.getParameter<double>("kernelPhiScale");
     if (kernelEtaScale <= 0.0 || kernelPhiScale <= 0.0)
         throw cms::Exception("FFTJetBadConfig")
             << "invalid kernel scale" << std::endl;
 
     // FFT assumes that the grid extent in eta is 2*Pi. Adjust the
     // kernel scale in eta to compensate.
     kernelEtaScale *= (2.0*M_PI/(energyFlow->etaMax() - energyFlow->etaMin()));
 
     // Are we going to try to fix the efficiency near detector edges?
     const bool fixEfficiency = ps.getParameter<bool>("fixEfficiency");
 
     // Minimum and maximum eta bin for the convolver
     unsigned convolverMinBin = 0, convolverMaxBin = 0;
     if (fixEfficiency || !use2dKernel)
     {
         convolverMinBin = ps.getParameter<unsigned>("convolverMinBin");
         convolverMaxBin = ps.getParameter<unsigned>("convolverMaxBin");
     }
 
     if (use2dKernel)
     {
         // Build the FFT engine
         engine = std::auto_ptr<MyFFTEngine>(
             new MyFFTEngine(energyFlow->nEta(), energyFlow->nPhi()));
 
         // 2d kernel
         kernel2d = std::auto_ptr<fftjet::AbsFrequencyKernel>(
             new fftjet::DiscreteGauss2d(
                 kernelEtaScale, kernelPhiScale,
                 energyFlow->nEta(), energyFlow->nPhi()));
 
         // 2d convolver
         convolver = std::auto_ptr<fftjet::AbsConvolverBase<Real> >(
             new fftjet::FrequencyKernelConvolver<Real,Complex>(
                 engine.get(), kernel2d.get(),
                 convolverMinBin, convolverMaxBin));
     }
     else
     {
         // Need separate FFT engines for eta and phi
         engine = std::auto_ptr<MyFFTEngine>(
             new MyFFTEngine(1, energyFlow->nEta()));
         anotherEngine = std::auto_ptr<MyFFTEngine>(
             new MyFFTEngine(1, energyFlow->nPhi()));
 
         // 1d kernels
         etaKernel = std::auto_ptr<fftjet::AbsFrequencyKernel1d>(
             new fftjet::DiscreteGauss1d(kernelEtaScale, energyFlow->nEta()));
 
         phiKernel = std::auto_ptr<fftjet::AbsFrequencyKernel1d>(
             new fftjet::DiscreteGauss1d(kernelPhiScale, energyFlow->nPhi()));
 
         // Sequential convolver
         convolver = std::auto_ptr<fftjet::AbsConvolverBase<Real> >(
             new fftjet::FrequencySequentialConvolver<Real,Complex>(
                 engine.get(), anotherEngine.get(),
                 etaKernel.get(), phiKernel.get(),
                 etaDependentScaleFactors, convolverMinBin,
                 convolverMaxBin, fixEfficiency));
     }
 }

fftjet::PeakFinder FFTJetPatRecoProducer::buildPeakFinder ( const edm::ParameterSet & ps )

protected

Definition at line 339 of file FFTJetPatRecoProducer.cc.

References fftjetcms::FFTJetInterface::energyFlow, edm::hlt::Exception, and edm::ParameterSet::getParameter().

Referenced by FFTJetPatRecoProducer().

 {
     const double peakFinderMaxEta = ps.getParameter<double>("peakFinderMaxEta");
     if (peakFinderMaxEta <= 0.0)
         throw cms::Exception("FFTJetBadConfig")
             << "invalid peak finder eta cut" << std::endl;
     const double maxMagnitude = ps.getParameter<double>("peakFinderMaxMagnitude");
     int minBin = energyFlow->getEtaBin(-peakFinderMaxEta);
     if (minBin < 0)
         minBin = 0;
     int maxBin = energyFlow->getEtaBin(peakFinderMaxEta) + 1;
     if (maxBin < 0)
         maxBin = 0;
     return fftjet::PeakFinder(maxMagnitude, true, minBin, maxBin);
 }

template<class Real >

void FFTJetPatRecoProducer::buildSparseProduct ( edm::Event & ev ) const

protected

Definition at line 369 of file FFTJetPatRecoProducer.cc.

References DDVectorGetter::check(), edm::hlt::Exception, fftjetcms::FFTJetInterface::getEventScale(), fftjetcms::FFTJetInterface::outputLabel, edm::Event::put(), sequencer, fftjetcms::sparsePeakTreeFromStorable(), fftjetcms::sparsePeakTreeToStorable(), sparseTree, fftjetcms::FFTJetInterface::storeInSinglePrecision(), MainPageGenerator::tree, and verifyDataConversion.

 {
     typedef reco::PattRecoTree<Real,reco::PattRecoPeak<Real> > StoredTree;
 
     std::auto_ptr<StoredTree> tree(new StoredTree());
 
     sparsePeakTreeToStorable(sparseTree,
                              sequencer->maxAdaptiveScales(),
                              tree.get());
 
     // Check that we can restore the tree
     if (verifyDataConversion && !storeInSinglePrecision())
     {
         SparseTree check;
         const std::vector<double>& scalesUsed(sequencer->getInitialScales());
         sparsePeakTreeFromStorable(*tree, &scalesUsed, getEventScale(), &check);
         if (sparseTree != check)
             throw cms::Exception("FFTJetInterface")
                 << "Data conversion failed for sparse clustering tree"
                 << std::endl;
     }
 
     ev.put(tree, outputLabel);
 }

void FFTJetPatRecoProducer::endJob ( void )

overrideprotectedvirtual

Reimplemented from edm::EDProducer.

Definition at line 503 of file FFTJetPatRecoProducer.cc.

References externalGridStream, and storeGridsExternally.

 {
     if (storeGridsExternally)
         externalGridStream.close();
 }

FFTJetPatRecoProducer& FFTJetPatRecoProducer::operator= ( const FFTJetPatRecoProducer & )

private

void FFTJetPatRecoProducer::produce	(	edm::Event &	iEvent,
		const edm::EventSetup &	iSetup
	)

overrideprotectedvirtual

Implements edm::EDProducer.

Definition at line 423 of file FFTJetPatRecoProducer.cc.

References assert(), DDVectorGetter::check(), clusteringTree, completeEventDataCutoff, fftjetcms::convert_Grid2d_to_float(), fftjetcms::copy_Grid2d_data(), fftjetcms::FFTJetInterface::discretizeEnergyFlow(), fftjetcms::FFTJetInterface::energyFlow, edm::hlt::Exception, externalGridStream, extGrid, g, fftjetcms::FFTJetInterface::getEventScale(), iEvent, fftjetcms::FFTJetInterface::loadInputCollection(), makeClusteringTree, fftjetcms::FFTJetInterface::outputLabel, edm::Event::put(), sequencer, sparseTree, sparsifier, sparsify, storeDiscretizationGrid, storeGridsExternally, fftjetcms::FFTJetInterface::storeInSinglePrecision(), and verifyDataConversion.

 {
     loadInputCollection(iEvent);
     discretizeEnergyFlow();
 
     if (makeClusteringTree)
     {
         sequencer->run(*energyFlow, clusteringTree);
         if (getEventScale() > 0.0)
         sequencer->insertCompleteEvent(getEventScale(), *energyFlow,
                                            clusteringTree, completeEventDataCutoff);
 
         if (sparsify)
         {
             sparsifier->sparsify(*clusteringTree, &sparseTree);
 
             // Do not call the "sortNodes" method of the sparse tree here.
             // Currently, the nodes are sorted by daughter number.
             // This is the way we want it in storage because the stored
             // tree does not include daughter ordering info explicitly.
 
             if (storeInSinglePrecision())
                 buildSparseProduct<float>(iEvent);
             else
                 buildSparseProduct<double>(iEvent);
         }
         else
         {
         if (storeInSinglePrecision())
                 buildDenseProduct<float>(iEvent);
             else
                 buildDenseProduct<double>(iEvent);
         }
     }
 
     if (storeDiscretizationGrid)
     {
         const fftjet::Grid2d<Real>& g(*energyFlow);
 
         std::auto_ptr<reco::DiscretizedEnergyFlow> flow(
             new reco::DiscretizedEnergyFlow(
                 g.data(), g.title(), g.etaMin(), g.etaMax(),
                 g.phiBin0Edge(), g.nEta(), g.nPhi()));
 
         if (verifyDataConversion)
         {
             fftjet::Grid2d<Real> check(
                 flow->nEtaBins(), flow->etaMin(), flow->etaMax(),
                 flow->nPhiBins(), flow->phiBin0Edge(), flow->title());
             check.blockSet(flow->data(), flow->nEtaBins(), flow->nPhiBins());
             assert(g == check);
         }
 
         iEvent.put(flow, outputLabel);
     }
 
     if (storeGridsExternally)
     {
         if (extGrid)
             copy_Grid2d_data(extGrid, *energyFlow);
         else
             extGrid = convert_Grid2d_to_float(*energyFlow);
         if (!extGrid->write(externalGridStream))
         {
             throw cms::Exception("FFTJetPatRecoProducer::produce")
                 << "Failed to write grid data into an external file"
                 << std::endl;
         }
     }
 }