Calorimetry/plugins/TreeWriterForEcalCorrection.cc

0001 // system include files
0002 #include <memory>
0003
0004 // user include files
0005 #include "FWCore/Framework/interface/Frameworkfwd.h"
0006 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0007
0008 #include "FWCore/Framework/interface/Event.h"
0009 #include "FWCore/Framework/interface/MakerMacros.h"
0010
0011 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0012
0013 #include "SimDataFormats/CaloHit/interface/PCaloHitContainer.h"
0014 #include "DataFormats/HepMCCandidate/interface/GenParticle.h"
0015
0016 #include "FWCore/Framework/interface/ESHandle.h"
0017
0018 #include "FWCore/ServiceRegistry/interface/Service.h"
0019 #include "CommonTools/UtilAlgos/interface/TFileService.h"
0020
0021 #include "TTree.h"
0022
0023 class TreeWriterForEcalCorrection : public edm::one::EDAnalyzer<edm::one::SharedResources> {
0024 public:
0025   explicit TreeWriterForEcalCorrection(const edm::ParameterSet&);
0026   ~TreeWriterForEcalCorrection() override = default;
0027
0028   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0029
0030 private:
0031   void analyze(const edm::Event&, const edm::EventSetup&) override;
0032
0033   const edm::EDGetTokenT<reco::GenParticleCollectio> tok_gen_;
0034   const edm::EDGetTokenT<edm::PCaloHitContainer> tok_ebf_;
0035   const edm::EDGetTokenT<edm::PCaloHitContainer> tok_eef_;
0036   const edm::EDGetTokenT<edm::PCaloHitContainer> tok_esf_;
0037   const edm::EDGetTokenT<edm::PCaloHitContainer> tok_ebs_;
0038   const edm::EDGetTokenT<edm::PCaloHitContainer> tok_ees_;
0039   const edm::EDGetTokenT<edm::PCaloHitContainer> tok_ess_;
0040   TTree* tree_;
0041   float tree_e_, tree_eta_, tree_response_;
0042 };
0043
0044 TreeWriterForEcalCorrection::TreeWriterForEcalCorrection(const edm::ParameterSet& iConfig)
0045     : tok_gen_(consumes<reco::GenParticleCollectio>(edm::InputTag("genParticles", ""))),
0046       tok_ebf_(consumes<edm::PCaloHitContainer>(edm::InputTag("fastSimProducer", "EcalHitsEB"))),
0047       tok_eef_(consumes<edm::PCaloHitContainer>(edm::InputTag("fastSimProducer", "EcalHitsEE"))),
0048       tok_esf_(consumes<edm::PCaloHitContainer>(edm::InputTag("fastSimProducer", "EcalHitsES"))),
0049       tok_ebs_(consumes<edm::PCaloHitContainer>(edm::InputTag("g4SimHits", "EcalHitsEB"))),
0050       tok_ees_(consumes<edm::PCaloHitContainer>(edm::InputTag("g4SimHits", "EcalHitsEE"))),
0051       tok_ess_(consumes<edm::PCaloHitContainer>(edm::InputTag("g4SimHits", "EcalHitsES"))) {
0052   usesResource(TFileService::kSharedResource);
0053   edm::Service<TFileService> file;
0054   tree_ = file->make<TTree>("responseTree", "same info as 3dhisto");
0055   tree_->Branch("e", &tree_e_, "e/F");
0056   tree_->Branch("eta", &tree_eta_, "eta/F");
0057   tree_->Branch("r", &tree_response_, "r/F");
0058 }
0059
0060 void TreeWriterForEcalCorrection::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0061   // get generated particles
0062   const edm::Handle<reco::GenParticleCollection>& GenParticles = iEvent.getHandle(tok_gen_);
0063
0064   // As this module is intended for single particle guns, there should be
0065   // exactly one generated particle.
0066   if (GenParticles->size() != 1) {
0067     throw cms::Exception("MismatchedInputFiles") << "Intended for particle guns only\n";
0068   }
0069
0070   // I assume here that the tracker simulation is disabled and no vertex
0071   // smearing is done, so that the generated particles is exactly the same
0072   // particle as the particle hitting the entrace of the ECAL.
0073   reco::GenParticle gen = GenParticles->at(0);
0074   float genE = gen.energy();
0075   float genEta = gen.eta();
0076   // genEta is positive for my photongun sample per definition.
0077   // If not, you could take the absolute value here.
0078
0079   // get sim hits
0080   edm::Handle<edm::PCaloHitContainer> SimHitsEB = iEvent.getHandle(tok_ebf_);
0081   // Finds out automatically, if this is fullsim or fastsim
0082   bool isFastSim = SimHitsEB.isValid();
0083   if (!isFastSim)
0084     SimHitsEB = iEvent.getHandle(tok_ebs_);
0085   const edm::Handle<edm::PCaloHitContainer> SimHitsEE =
0086       isFastSim ? iEvent.getHandle(tok_eef_) : iEvent.getHandle(tok_ees_);
0087   const edm::Handle<edm::PCaloHitContainer> SimHitsES =
0088       isFastSim ? iEvent.getHandle(tok_esf_) : iEvent.getHandle(tok_ess_);
0089
0090   // merge them into one single vector
0091   auto SimHits = *(SimHitsEB.product());
0092   SimHits.insert(SimHits.end(), SimHitsEE->begin(), SimHitsEE->end());
0093   SimHits.insert(SimHits.end(), SimHitsES->begin(), SimHitsES->end());
0094
0095   // As we only had one generated particle (and hopefully no pileup),
0096   // the total energy is due to the generated particle only
0097   float energyTotal = 0;
0098   for (auto const& Hit : SimHits) {
0099     energyTotal += Hit.energy();
0100   }
0101
0102   tree_e_ = genE;
0103   tree_eta_ = genEta;
0104   tree_response_ = energyTotal / genE;
0105   tree_->Fill();
0106 }
0107
0108 void TreeWriterForEcalCorrection::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0109   edm::ParameterSetDescription desc;
0110   descriptions.add("ecalScaleFactorCalculator", desc);
0111 }
0112
0113 DEFINE_FWK_MODULE(TreeWriterForEcalCorrection);