Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 12:32:17

0001 #include "CommonTools/UtilAlgos/interface/TFileService.h"
0002 #include "DataFormats/Math/interface/GeantUnits.h"
0003 
0004 #include "FWCore/Framework/interface/Event.h"
0005 #include "FWCore/Framework/interface/EventSetup.h"
0006 #include "FWCore/Framework/interface/Frameworkfwd.h"
0007 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0008 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0009 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0010 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0011 #include "FWCore/ServiceRegistry/interface/Service.h"
0012 #include "FWCore/Utilities/interface/InputTag.h"
0013 
0014 #include "SimDataFormats/CaloHit/interface/MaterialInformation.h"
0015 
0016 #include "TProfile.h"
0017 #include "TProfile2D.h"
0018 
0019 #include <iostream>
0020 #include <string>
0021 #include <vector>
0022 
0023 //#define EDM_ML_DEBUG
0024 
0025 using namespace geant_units::operators;
0026 
0027 class MaterialBudgetVolumeAnalysis : public edm::one::EDAnalyzer<edm::one::SharedResources> {
0028 public:
0029   explicit MaterialBudgetVolumeAnalysis(edm::ParameterSet const&);
0030 
0031   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0032 
0033 private:
0034   void analyze(edm::Event const&, edm::EventSetup const&) override;
0035   void bookHisto();
0036 
0037   const std::vector<std::string> names_;
0038   const edm::InputTag tag_;
0039   const int binEta_, binPhi_;
0040   const double etaLow_, etaHigh_, phiLow_, phiHigh_;
0041   edm::EDGetTokenT<edm::MaterialInformationContainer> tok_info_;
0042   std::vector<TProfile*> meStepEta_, meStepPhi_;
0043   std::vector<TProfile*> meRadLEta_, meRadLPhi_;
0044   std::vector<TProfile*> meIntLEta_, meIntLPhi_;
0045   std::vector<TProfile2D*> meStepEtaPhi_, meRadLEtaPhi_, meIntLEtaPhi_;
0046 };
0047 
0048 MaterialBudgetVolumeAnalysis::MaterialBudgetVolumeAnalysis(const edm::ParameterSet& p)
0049     : names_(p.getParameter<std::vector<std::string> >("names")),
0050       tag_(p.getParameter<edm::InputTag>("inputTag")),
0051       binEta_(p.getParameter<int>("nBinEta")),
0052       binPhi_(p.getParameter<int>("nBinPhi")),
0053       etaLow_(p.getParameter<double>("etaLow")),
0054       etaHigh_(p.getParameter<double>("etaHigh")),
0055       phiLow_(-1._pi),
0056       phiHigh_(1._pi) {
0057   usesResource(TFileService::kSharedResource);
0058   tok_info_ = consumes<edm::MaterialInformationContainer>(tag_);
0059 
0060   edm::LogVerbatim("MaterialBudget") << "MaterialBudgetVolumeAnalysis: Eta plot: NX " << binEta_ << " Range "
0061                                      << -etaLow_ << ":" << etaHigh_ << " Phi plot: NX " << binPhi_ << " Range "
0062                                      << -1._pi << ":" << 1._pi << " for " << names_.size() << " detectors from "
0063                                      << tag_;
0064   std::ostringstream st1;
0065   for (unsigned int k = 0; k < names_.size(); ++k)
0066     st1 << " [" << k << "] " << names_[k];
0067   edm::LogVerbatim("MaterialBudget") << "MaterialBudgetVolume: " << st1.str();
0068   bookHisto();
0069 }
0070 
0071 void MaterialBudgetVolumeAnalysis::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0072   edm::ParameterSetDescription desc;
0073   std::vector<std::string> names = {"BEAM",
0074                                     "BEAM1",
0075                                     "BEAM2",
0076                                     "BEAM3",
0077                                     "BEAM4",
0078                                     "Tracker",
0079                                     "ECAL",
0080                                     "HCal",
0081                                     "MUON",
0082                                     "VCAL",
0083                                     "MGNT",
0084                                     "OQUA",
0085                                     "CALOEC",
0086                                     "HFNoseVol"};
0087   desc.add<std::vector<std::string> >("names", names);
0088   desc.add<edm::InputTag>("inputTag", edm::InputTag("g4SimHits", "MaterialInformation"));
0089   desc.add<int>("nBinEta", 300);
0090   desc.add<int>("nBinPhi", 180);
0091   desc.add<double>("etaLow", -6.0);
0092   desc.add<double>("etaHigh", 6.0);
0093   descriptions.add("materialBudgetVolumeAnalysis", desc);
0094 }
0095 
0096 void MaterialBudgetVolumeAnalysis::analyze(const edm::Event& iEvent, const edm::EventSetup&) {
0097   edm::Handle<edm::MaterialInformationContainer> materialInformationContainer;
0098   iEvent.getByToken(tok_info_, materialInformationContainer);
0099 #ifdef EDM_ML_DEBUG
0100   unsigned int nsize(0), ntot(0), nused(0);
0101 #endif
0102   if (materialInformationContainer.isValid()) {
0103 #ifdef EDM_ML_DEBUG
0104     nsize = materialInformationContainer->size();
0105 #endif
0106     for (const auto& it : *(materialInformationContainer.product())) {
0107 #ifdef EDM_ML_DEBUG
0108       ntot++;
0109 #endif
0110       if (std::find(names_.begin(), names_.end(), it.vname()) != names_.end()) {
0111 #ifdef EDM_ML_DEBUG
0112         nused++;
0113 #endif
0114         unsigned int k =
0115             static_cast<unsigned int>(std::find(names_.begin(), names_.end(), it.vname()) - names_.begin());
0116         meStepEta_[k]->Fill(it.trackEta(), it.stepLength());
0117         meRadLEta_[k]->Fill(it.trackEta(), it.radiationLength());
0118         meIntLEta_[k]->Fill(it.trackEta(), it.interactionLength());
0119         meStepPhi_[k]->Fill(it.trackPhi(), it.stepLength());
0120         meRadLPhi_[k]->Fill(it.trackPhi(), it.radiationLength());
0121         meIntLPhi_[k]->Fill(it.trackPhi(), it.interactionLength());
0122         meStepEtaPhi_[k]->Fill(it.trackEta(), it.trackPhi(), it.stepLength());
0123         meRadLEtaPhi_[k]->Fill(it.trackEta(), it.trackPhi(), it.radiationLength());
0124         meIntLEtaPhi_[k]->Fill(it.trackEta(), it.trackPhi(), it.interactionLength());
0125       }
0126     }
0127   }
0128 #ifdef EDM_ML_DEBUG
0129   edm::LogVerbatim("MaterialBudget") << "MaterialInformation with " << nsize << ":" << ntot << " elements of which "
0130                                      << nused << " are used";
0131 #endif
0132 }
0133 
0134 void MaterialBudgetVolumeAnalysis::bookHisto() {
0135   edm::Service<TFileService> fs;
0136   char name[40], title[100];
0137   for (unsigned int k = 0; k < names_.size(); ++k) {
0138     sprintf(name, "stepEta%s", names_[k].c_str());
0139     sprintf(title, "MB(Step) vs #eta for %s", names_[k].c_str());
0140     meStepEta_.emplace_back(fs->make<TProfile>(name, title, binEta_, etaLow_, etaHigh_));
0141     sprintf(name, "radlEta%s", names_[k].c_str());
0142     sprintf(title, "MB(X0) vs #eta for %s", names_[k].c_str());
0143     meRadLEta_.emplace_back(fs->make<TProfile>(name, title, binEta_, etaLow_, etaHigh_));
0144     sprintf(name, "intlEta%s", names_[k].c_str());
0145     sprintf(title, "MB(L0) vs #eta for %s", names_[k].c_str());
0146     meIntLEta_.emplace_back(fs->make<TProfile>(name, title, binEta_, etaLow_, etaHigh_));
0147     sprintf(name, "stepPhi%s", names_[k].c_str());
0148     sprintf(title, "MB(Step) vs #phi for %s", names_[k].c_str());
0149     meStepPhi_.emplace_back(fs->make<TProfile>(name, title, binPhi_, phiLow_, phiHigh_));
0150     sprintf(name, "radlPhi%s", names_[k].c_str());
0151     sprintf(title, "MB(X0) vs #phi for %s", names_[k].c_str());
0152     meRadLPhi_.emplace_back(fs->make<TProfile>(name, title, binPhi_, phiLow_, phiHigh_));
0153     sprintf(name, "intlPhi%s", names_[k].c_str());
0154     sprintf(title, "MB(L0) vs #phi for %s", names_[k].c_str());
0155     meIntLPhi_.emplace_back(fs->make<TProfile>(name, title, binPhi_, phiLow_, phiHigh_));
0156     sprintf(name, "stepEtaPhi%s", names_[k].c_str());
0157     sprintf(title, "MB(Step) vs #eta and #phi for %s", names_[k].c_str());
0158     meStepEtaPhi_.emplace_back(
0159         fs->make<TProfile2D>(name, title, binEta_ / 2, etaLow_, etaHigh_, binPhi_ / 2, phiLow_, phiHigh_));
0160     sprintf(name, "radlEtaPhi%s", names_[k].c_str());
0161     sprintf(title, "MB(X0) vs #eta and #phi for %s", names_[k].c_str());
0162     meRadLEtaPhi_.emplace_back(
0163         fs->make<TProfile2D>(name, title, binEta_ / 2, etaLow_, etaHigh_, binPhi_ / 2, phiLow_, phiHigh_));
0164     sprintf(name, "intlEtaPhi%s", names_[k].c_str());
0165     sprintf(title, "MB(L0) vs #eta and #phi for %s", names_[k].c_str());
0166     meIntLEtaPhi_.emplace_back(
0167         fs->make<TProfile2D>(name, title, binEta_ / 2, etaLow_, etaHigh_, binPhi_ / 2, phiLow_, phiHigh_));
0168   }
0169 }
0170 
0171 // define this as a plug-in
0172 #include "FWCore/Framework/interface/MakerMacros.h"
0173 DEFINE_FWK_MODULE(MaterialBudgetVolumeAnalysis);