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/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::EDAnalyzer {
0024 public:
0025   explicit TreeWriterForEcalCorrection(const edm::ParameterSet&);
0026   ~TreeWriterForEcalCorrection() override{};
0027
0028   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0029
0030 private:
0031   void analyze(const edm::Event&, const edm::EventSetup&) override;
0032
0033   edm::Service<TFileService> file;
0034   TTree* tree;
0035   float tree_e, tree_eta, tree_response;
0036 };
0037
0038 TreeWriterForEcalCorrection::TreeWriterForEcalCorrection(const edm::ParameterSet& iConfig) {
0039   tree = file->make<TTree>("responseTree", "same info as 3dhisto");
0040   tree->Branch("e", &tree_e, "e/F");
0041   tree->Branch("eta", &tree_eta, "eta/F");
0042   tree->Branch("r", &tree_response, "r/F");
0043 }
0044
0045 void TreeWriterForEcalCorrection::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0046   // get generated particles
0047   edm::Handle<reco::GenParticleCollection> GenParticles;
0048   iEvent.getByLabel("genParticles", "", GenParticles);
0049
0050   // As this module is intended for single particle guns, there should be
0051   // exactly one generated particle.
0052   if (GenParticles->size() != 1) {
0053     throw cms::Exception("MismatchedInputFiles") << "Intended for particle guns only\n";
0054   }
0055
0056   // I assume here that the tracker simulation is disabled and no vertex
0057   // smearing is done, so that the generated particles is exactly the same
0058   // particle as the particle hitting the entrace of the ECAL.
0059   reco::GenParticle gen = GenParticles->at(0);
0060   float genE = gen.energy();
0061   float genEta = gen.eta();
0062   // genEta is positive for my photongun sample per definition.
0063   // If not, you could take the absolute value here.
0064
0065   // get sim hits
0066   edm::Handle<edm::PCaloHitContainer> SimHitsEB;
0067   edm::Handle<edm::PCaloHitContainer> SimHitsEE;
0068   edm::Handle<edm::PCaloHitContainer> SimHitsES;
0069
0070   // Finds out automatically, if this is fullsim or fastsim
0071   bool isFastSim = iEvent.getByLabel("fastSimProducer", "EcalHitsEB", SimHitsEB);
0072   if (isFastSim) {
0073     iEvent.getByLabel("fastSimProducer", "EcalHitsEE", SimHitsEE);
0074     iEvent.getByLabel("fastSimProducer", "EcalHitsES", SimHitsES);
0075   } else {
0076     iEvent.getByLabel("g4SimHits", "EcalHitsEB", SimHitsEB);
0077     iEvent.getByLabel("g4SimHits", "EcalHitsEE", SimHitsEE);
0078     iEvent.getByLabel("g4SimHits", "EcalHitsES", SimHitsES);
0079   }
0080
0081   // merge them into one single vector
0082   auto SimHits = *SimHitsEB;
0083   SimHits.insert(SimHits.end(), SimHitsEE->begin(), SimHitsEE->end());
0084   SimHits.insert(SimHits.end(), SimHitsES->begin(), SimHitsES->end());
0085
0086   // As we only had one generated particle (and hopefully no pileup),
0087   // the total energy is due to the generated particle only
0088   float energyTotal = 0;
0089   for (auto const& Hit : SimHits) {
0090     energyTotal += Hit.energy();
0091   }
0092
0093   tree_e = genE;
0094   tree_eta = genEta;
0095   tree_response = energyTotal / genE;
0096   tree->Fill();
0097 }
0098
0099 void TreeWriterForEcalCorrection::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0100   edm::ParameterSetDescription desc;
0101   descriptions.add("ecalScaleFactorCalculator", desc);
0102 }
0103
0104 DEFINE_FWK_MODULE(TreeWriterForEcalCorrection);