Skip to content

Commit

Permalink
Merged commit for PR
Browse files Browse the repository at this point in the history
Changes to existing code...

Revert "Changes to existing code..."

This reverts commit b5da60b.

Updated gitignore

Working producer on random data

Added back in ECAL nTuplizer

Plugged in DRN output

Everything works, but doesn't take input

Moved everything into classes that look like the EGM classes

Looping over scs, but weird segfaults

Fully plugged in on the right now

Everything is in except rho and H/E

Pushing for Rajdeep

Added rho and electron objects..

It all works! :D

Removed some troubleshooting code

Renamed folder

Propegate changed path to #include statements

Moved things into RecoEcal

5000 events

Rename with camelCase convention

Added rho max; new nTuplizer

Removed code that doesn't belong in CMSSW

Added test code back into RecoEcal

Fixed model path

Cleaned up a little

Added cfi for DRN producer

Added test code for DRN producer

Added some junk to gitignore

Reverted gitignore changes

Deleted stupid extra file
  • Loading branch information
ssrothman committed Oct 26, 2021
1 parent bec73cc commit 45b7741
Show file tree
Hide file tree
Showing 9 changed files with 1,047 additions and 0 deletions.
1 change: 1 addition & 0 deletions RecoEcal/EgammaClusterAlgos/BuildFile.xml
Expand Up @@ -14,6 +14,7 @@
<use name="Geometry/CaloTopology"/>
<use name="RecoEgamma/EgammaTools"/>
<use name="RecoParticleFlow/PFClusterTools"/>
<use name="HeterogeneousCore/SonicTriton"/>
<use name="clhep"/>
<export>
<lib name="1"/>
Expand Down
94 changes: 94 additions & 0 deletions RecoEcal/EgammaClusterAlgos/interface/SCEnergyCorrectorDRN.h
@@ -0,0 +1,94 @@
//--------------------------------------------------------------------------------------------------
//
// SCEnergyCorrectorDRN
//
// Helper Class for applying regression-based energy corrections with DRN implimentation
//
// Based on RecoEcal/EgammaClusterAlgos/SCEnergyCorrectorSemiParm
//
// Author: Simon Rothman (MIT, UMN)
//
//--------------------------------------------------------------------------------------------------

#ifndef Progression_EGM_DRN_SCEnergyCorrectorDRN_h
#define Progression_EGM_DRN_SCEnergyCorrectorDRN_h

#include "HeterogeneousCore/SonicTriton/interface/TritonEDProducer.h"

#include "FWCore/Framework/interface/EventSetup.h"
#include "FWCore/Framework/interface/ESHandle.h"
#include "FWCore/Framework/interface/ConsumesCollector.h"
#include "FWCore/ParameterSet/interface/ParameterSet.h"

#include "Geometry/Records/interface/CaloTopologyRecord.h"
#include "Geometry/CaloTopology/interface/CaloTopology.h"
#include "Geometry/Records/interface/CaloGeometryRecord.h"
#include "Geometry/CaloGeometry/interface/CaloGeometry.h"

#include "DataFormats/EgammaReco/interface/SuperCluster.h"
#include "DataFormats/EcalRecHit/interface/EcalRecHit.h"
#include "DataFormats/EcalRecHit/interface/EcalRecHitCollections.h"
#include "DataFormats/EcalDetId/interface/EcalSubdetector.h"
#include "DataFormats/ParticleFlowReco/interface/PFRecHit.h"
#include "DataFormats/EgammaReco/interface/SuperClusterFwd.h"
#include "DataFormats/VertexReco/interface/VertexFwd.h"

#include "CondFormats/GBRForest/interface/GBRForestD.h"
#include "CondFormats/DataRecord/interface/GBRDWrapperRcd.h"

#include "RecoEgamma/EgammaTools/interface/EgammaBDTOutputTransformer.h"
#include "RecoEgamma/EgammaTools/interface/HGCalShowerShapeHelper.h"

#include <sstream>
#include <string>
#include <vector>
#include <random>

class SCEnergyCorrectorDRN {
public:
SCEnergyCorrectorDRN();
//if you want override the default on where conditions are consumed, you need to use
//the other constructor and then call setTokens approprately
SCEnergyCorrectorDRN(const edm::ParameterSet& iConfig, edm::ConsumesCollector cc);

static void fillPSetDescription(edm::ParameterSetDescription& desc);
static edm::ParameterSetDescription makePSetDescription();

template <edm::Transition tr = edm::Transition::BeginLuminosityBlock>
void setTokens(const edm::ParameterSet& iConfig, edm::ConsumesCollector cc);

void setEventSetup(const edm::EventSetup& es);
void setEvent(const edm::Event& e);

//std::pair<double, double> getCorrections(const reco::SuperCluster& sc) const;

void makeInput(const edm::Event& iEvent, TritonInputMap& iInput, const reco::SuperClusterCollection& inputSCs) const;
TritonOutput<float> getOutput(const TritonOutputMap& iOutput);

//void modifyObject(reco::SuperCluster& sc) const;

private:
const CaloTopology* caloTopo_;
const CaloGeometry* caloGeom_;
edm::ESGetToken<CaloTopology, CaloTopologyRecord> caloTopoToken_;
edm::ESGetToken<CaloGeometry, CaloGeometryRecord> caloGeomToken_;

edm::EDGetTokenT<EcalRecHitCollection> tokenEBRecHits_;
edm::EDGetTokenT<EcalRecHitCollection> tokenEERecHits_;
edm::EDGetTokenT<double> rhoToken_;

edm::Handle<EcalRecHitCollection> recHitsEB_;
edm::Handle<EcalRecHitCollection> recHitsEE_;

edm::Handle<double> rhoHandle_;
};

template <edm::Transition esTransition>
void SCEnergyCorrectorDRN::setTokens(const edm::ParameterSet& iConfig, edm::ConsumesCollector cc) {
tokenEBRecHits_ = cc.consumes<EcalRecHitCollection>(iConfig.getParameter<edm::InputTag>("ecalRecHitsEB"));
tokenEERecHits_ = cc.consumes<EcalRecHitCollection>(iConfig.getParameter<edm::InputTag>("ecalRecHitsEE"));
caloGeomToken_ = cc.esConsumes<CaloGeometry, CaloGeometryRecord, esTransition>();
caloTopoToken_ = cc.esConsumes<CaloTopology, CaloTopologyRecord, esTransition>();
rhoToken_ = cc.consumes<double>(iConfig.getParameter<edm::InputTag>("rhoFastJet"));
}
#endif
116 changes: 116 additions & 0 deletions RecoEcal/EgammaClusterAlgos/src/SCEnergyCorrectorDRN.cc
@@ -0,0 +1,116 @@
#include "RecoEcal/EgammaClusterAlgos/interface/SCEnergyCorrectorDRN.h"

#include "FWCore/Utilities/interface/isFinite.h"
#include "FWCore/Utilities/interface/Transition.h"
#include "DataFormats/EcalDetId/interface/EcalSubdetector.h"
#include "DataFormats/VertexReco/interface/Vertex.h"
#include "DataFormats/Math/interface/deltaPhi.h"
#include "RecoEcal/EgammaCoreTools/interface/EcalTools.h"
#include "RecoEcal/EgammaCoreTools/interface/EcalClusterTools.h"
#include "RecoEgamma/EgammaTools/interface/EgammaHGCALIDParamDefaults.h"

#include <vdt/vdtMath.h>

#include <iostream>

static const float RHO_MAX=15.0f;

SCEnergyCorrectorDRN::SCEnergyCorrectorDRN()
: caloTopo_(nullptr),
caloGeom_(nullptr) {}

SCEnergyCorrectorDRN::SCEnergyCorrectorDRN(const edm::ParameterSet& iConfig, edm::ConsumesCollector cc)
: SCEnergyCorrectorDRN() {
setTokens(iConfig, cc);
}

void SCEnergyCorrectorDRN::fillPSetDescription(edm::ParameterSetDescription& desc) {
desc.add<edm::InputTag>("ecalRecHitsEE", edm::InputTag("ecalRecHit", "EcalRecHitsEE"));
desc.add<edm::InputTag>("ecalRecHitsEB", edm::InputTag("ecalRecHit", "EcalRecHitsEB"));
desc.add<edm::InputTag>("rhoFastJet");
}

edm::ParameterSetDescription SCEnergyCorrectorDRN::makePSetDescription() {
edm::ParameterSetDescription desc;
fillPSetDescription(desc);
return desc;
}

void SCEnergyCorrectorDRN::setEventSetup(const edm::EventSetup& es) {
caloTopo_ = &es.getData(caloTopoToken_);
caloGeom_ = &es.getData(caloGeomToken_);
}

void SCEnergyCorrectorDRN::setEvent(const edm::Event& event) {
event.getByToken(tokenEBRecHits_, recHitsEB_);
event.getByToken(tokenEERecHits_, recHitsEE_);
event.getByToken(rhoToken_, rhoHandle_);
}

void SCEnergyCorrectorDRN::makeInput(const edm::Event& iEvent, TritonInputMap& iInput, const reco::SuperClusterCollection& inputSCs ) const {

std::vector<unsigned> nHits = {};
unsigned totalHits = 0;
unsigned n;
for(const auto& inputSC : inputSCs){
n = inputSC.hitsAndFractions().size();
totalHits+=n;
nHits.push_back(n);
}

//set shapes
auto& input1 = iInput.at("x__0");
input1.setShape(0, totalHits);
auto data1 = input1.allocate<float>();
auto& vdata1 = (*data1)[0];

auto& input2 = iInput.at("batch__1");
input2.setShape(0, totalHits);
auto data2 = input2.allocate<int64_t>();
auto& vdata2 = (*data2)[0];

auto& input3 = iInput.at("graphx__2");
input3.setShape(0, 2*nHits.size());
auto data3 = input3.allocate<float>();
auto& vdata3 = (*data3)[0];

//fill
unsigned batchNum = 0;
float En, frac, x, y, z;
for(const auto& inputSC : inputSCs){
const auto& hits = inputSC.hitsAndFractions();
const bool isEB = hits[0].first.subdetId() == EcalBarrel;
const auto& recHitsProduct = isEB ? recHitsEB_.product() : recHitsEE_.product();
for(const auto& hit : hits){
En = EcalClusterTools::recHitEnergy(hit.first, recHitsProduct);
frac = hit.second;
GlobalPoint position = caloGeom_->getGeometry(hit.first)->getPosition();
x = (position.x()+150.0)/300.0;
y = (position.y()+150.0)/300.0;
z = (position.z()+330.0)/660.0;
vdata1.push_back(x);
vdata1.push_back(y);
vdata1.push_back(z);
vdata1.push_back(En*frac/250.0);
vdata2.push_back(batchNum);
}
vdata3.push_back(*rhoHandle_/ RHO_MAX );
vdata3.push_back(0.0);
++batchNum;
}

// convert to server format
input1.toServer(data1);
input2.toServer(data2);
input3.toServer(data3);
}

TritonOutput<float> SCEnergyCorrectorDRN::getOutput(const TritonOutputMap& iOutput) {
//check the results
const auto& output1 = iOutput.begin()->second;
// convert from server format
const auto& serverout = output1.fromServer<float>();

return serverout;
}

1 change: 1 addition & 0 deletions RecoEcal/EgammaClusterProducers/BuildFile.xml
@@ -1,4 +1,5 @@
<use name="FWCore/Framework"/>
<use name="HeterogeneousCore/SonicTriton"/>
<use name="FWCore/ParameterSet"/>
<use name="DataFormats/EcalDigi"/>
<use name="DataFormats/EcalRecHit"/>
Expand Down
@@ -0,0 +1,46 @@
import FWCore.ParameterSet.Config as cms

DRNProducerEB = cms.EDProducer('SCEnergyCorrectorDRNProducer',
correctorCfg = cms.PSet(
ecalRecHitsEE = cms.InputTag('reducedEcalRecHitsEE'),
ecalRecHitsEB = cms.InputTag('reducedEcalRecHitsEB'),
rhoFastJet = cms.InputTag("fixedGridRhoAll"),
),
inputSCs = cms.InputTag('particleFlowSuperClusterECAL','particleFlowSuperClusterECALBarrel'),
Client = cms.PSet(
mode = cms.string("Async"),
preferredServer = cms.untracked.string(""),
timeout = cms.untracked.uint32(10),
modelName = cms.string("MustacheEB"),
modelVersion = cms.string(""),
modelConfigPath = cms.FileInPath("DRNData/models/MustacheEB/config.pbtxt"),
verbose = cms.untracked.bool(False),
allowedTries = cms.untracked.uint32(1),
useSharedMemory = cms.untracked.bool(True),
compression = cms.untracked.string(""),
),
)


DRNProducerEE = cms.EDProducer('SCEnergyCorrectorDRNProducer',
correctorCfg = cms.PSet(
ecalRecHitsEE = cms.InputTag('reducedEcalRecHitsEE'),
ecalRecHitsEB = cms.InputTag('reducedEcalRecHitsEB'),
rhoFastJet = cms.InputTag("fixedGridRhoAll"),
),
inputSCs = cms.InputTag('particleFlowSuperClusterECAL','particleFlowSuperClusterECALEndcapWithPreshower'),
Client = cms.PSet(
mode = cms.string("Async"),
preferredServer = cms.untracked.string(""),
timeout = cms.untracked.uint32(10),
modelName = cms.string('MustacheEE'),
modelVersion = cms.string(""),
modelConfigPath = cms.FileInPath("DRNData/models/MustacheEE/config.pbtxt"),
verbose = cms.untracked.bool(False),
allowedTries = cms.untracked.uint32(1),
useSharedMemory = cms.untracked.bool(True),
compression = cms.untracked.string(""),
),
)


116 changes: 116 additions & 0 deletions RecoEcal/EgammaClusterProducers/src/SCEnergyCorrectorDRNProducer.cc
@@ -0,0 +1,116 @@
#include "HeterogeneousCore/SonicTriton/interface/TritonEDProducer.h"
#include "HeterogeneousCore/SonicTriton/interface/TritonData.h"

#include "FWCore/Framework/interface/stream/EDProducer.h"
#include "FWCore/Framework/interface/Event.h"
#include "FWCore/Framework/interface/MakerMacros.h"
#include "FWCore/ParameterSet/interface/ParameterSet.h"

#include "DataFormats/EgammaReco/interface/SuperCluster.h"
#include "DataFormats/EgammaReco/interface/SuperClusterFwd.h"
#include "DataFormats/Common/interface/ValueMap.h"

#include "DataFormats/EgammaCandidates/interface/GsfElectron.h"

#include "RecoEcal/EgammaClusterAlgos/interface/SCEnergyCorrectorDRN.h"

#include <sstream>
#include <string>
#include <vector>
#include <random>

/*
* SCEnergyCorrectorDRNProducer
*
* Simple producer to generate a set of corrected superclusters with the DRN regression
* Based on RecoEcal/EgammaClusterProducers/SCEnergyCorrectorProducer by S. Harper (RAL/CERN)
*
* Author: Simon Rothman (UMN, MIT)
*
*/

static float sigmoid(float x){
return 1.0f/(1.0f + exp(-x));
}
static float ln2 = log(2);

static float logcorrection(float x){
return ln2 * 2 * (sigmoid(x) - 0.5);
}

static float correction(float x){
return exp(-logcorrection(x));
}

class SCEnergyCorrectorDRNProducer : public TritonEDProducer<> {
public:
explicit SCEnergyCorrectorDRNProducer(const edm::ParameterSet& iConfig);

void beginLuminosityBlock(const edm::LuminosityBlock& iLumi, const edm::EventSetup& iSetup) override;

void acquire(edm::Event const& iEvent, edm::EventSetup const& iSetup, Input& input) override;
void produce(edm::Event& iEvent, const edm::EventSetup& iSetup, Output const& iOutput) override;

static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);

private:
SCEnergyCorrectorDRN energyCorrector_;
edm::EDGetTokenT<reco::SuperClusterCollection> inputSCToken_;
};

SCEnergyCorrectorDRNProducer::SCEnergyCorrectorDRNProducer(const edm::ParameterSet& iConfig)
: TritonEDProducer<>(iConfig, "SCEnergyCorrectorDRNProducer"),
energyCorrector_(iConfig.getParameterSet("correctorCfg"), consumesCollector()),
inputSCToken_(consumes<reco::SuperClusterCollection>(iConfig.getParameter<edm::InputTag>("inputSCs")))
{
produces<reco::SuperClusterCollection>();
}

void SCEnergyCorrectorDRNProducer::beginLuminosityBlock(const edm::LuminosityBlock& iLumi, const edm::EventSetup& iSetup) {
energyCorrector_.setEventSetup(iSetup);
}

void SCEnergyCorrectorDRNProducer::acquire(edm::Event const& iEvent, edm::EventSetup const& iSetup, Input& iInput){
auto inputSCs = iEvent.get(inputSCToken_);

if(inputSCs.size()==0){
client_->setBatchSize(0);
return;
} else{
client_->setBatchSize(1);
}

energyCorrector_.setEvent(iEvent);
energyCorrector_.makeInput(iEvent, iInput, inputSCs);
}

void SCEnergyCorrectorDRNProducer::produce(edm::Event& iEvent, const edm::EventSetup& iSetup, Output const& iOutput){
auto inputSCs = iEvent.get(inputSCToken_);
if (inputSCs.size()==0)
return;

const auto& serverout = energyCorrector_.getOutput(iOutput);

auto corrSCs = std::make_unique<reco::SuperClusterCollection>();
unsigned i=0;
for (const auto& inputSC : inputSCs) {
float corrEn = correction(serverout[0][0+6*i]) * inputSC.rawEnergy();
corrSCs->push_back(inputSC);
corrSCs->back().setEnergy(corrEn);
corrSCs->back().setCorrectedEnergy(corrEn);
//corrSCs->back().setCorrectedEnergyUncertainty(-1.0);
++i;
}

auto scHandle = iEvent.put(std::move(corrSCs));
}

void SCEnergyCorrectorDRNProducer::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
edm::ParameterSetDescription desc;
desc.add<edm::ParameterSetDescription>("correctorCfg", SCEnergyCorrectorDRN::makePSetDescription());
TritonClient::fillPSetDescription(desc);
desc.add<edm::InputTag>("inputSCs", edm::InputTag("particleFlowSuperClusterECAL"));
descriptions.add("scEnergyCorrectorDRNProducer", desc);
}

DEFINE_FWK_MODULE(SCEnergyCorrectorDRNProducer);

0 comments on commit 45b7741

Please sign in to comment.