CMS 3D CMS Logo

CaloParticleValidation.cc
Go to the documentation of this file.
1 // -*- C++ -*-
2 //
3 // Class: CaloParticleValidation
4 // Original Author: Marco Rovere
5 // Created: Thu, 18 Jan 2018 15:54:55 GMT
6 //
7 //
8 
9 #include <string>
10 #include <unordered_map>
11 
12 // user include files
15 
18 
20 
28 
29 //
30 // class declaration
31 //
32 
42  dqm::reco::MonitorElement* pfcandidate_vect_sum_pt_; // This is indeed a cumulative istogram
43 };
44 
45 using Histograms_CaloParticleValidation = std::unordered_map<int, Histogram_CaloParticleSingle>;
46 
47 class CaloParticleValidation : public DQMGlobalEDAnalyzer<Histograms_CaloParticleValidation> {
48 public:
50  ~CaloParticleValidation() override;
51 
52  static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
53 
54 private:
56  edm::Run const&,
57  edm::EventSetup const&,
58  Histograms_CaloParticleValidation&) const override;
59 
60  void dqmAnalyze(edm::Event const&, edm::EventSetup const&, Histograms_CaloParticleValidation const&) const override;
61 
62  // ----------member data ---------------------------
66 };
67 
68 //
69 // constants, enums and typedefs
70 //
71 
72 //
73 // static data member definitions
74 //
75 
76 //
77 // constructors and destructor
78 //
80  : folder_(iConfig.getParameter<std::string>("folder")),
81  simPFClusters_(consumes<std::vector<reco::SuperCluster>>(iConfig.getParameter<edm::InputTag>("simPFClusters"))),
82  simPFCandidates_(consumes<reco::PFCandidateCollection>(iConfig.getParameter<edm::InputTag>("simPFCandidates"))) {
83  //now do what ever initialization is needed
84 }
85 
87  // do anything here that needs to be done at desctruction time
88  // (e.g. close files, deallocate resources etc.)
89 }
90 
91 //
92 // member functions
93 //
94 
95 // ------------ method called for each event ------------
96 
98  edm::EventSetup const& iSetup,
100  using namespace edm;
101 
102  Handle<std::vector<reco::SuperCluster>> simPFClustersHandle;
103  iEvent.getByToken(simPFClusters_, simPFClustersHandle);
104  std::vector<reco::SuperCluster> const& simPFClusters = *simPFClustersHandle;
105 
106  Handle<reco::PFCandidateCollection> simPFCandidatesHandle;
107  iEvent.getByToken(simPFCandidates_, simPFCandidatesHandle);
108  reco::PFCandidateCollection const& simPFCandidates = *simPFCandidatesHandle;
109 
110  // simPFSuperClusters
111  for (auto const& sc : simPFClusters) {
112  histos.at(0).simPFSuperClusterSize_->Fill((float)sc.clustersSize());
113  histos.at(0).simPFSuperClusterEnergy_->Fill(sc.rawEnergy());
114  }
115 
116  // simPFCandidates
117  int offset = 100000;
118  double ptx_tot = 0.;
119  double pty_tot = 0.;
120  for (auto const& pfc : simPFCandidates) {
121  size_t type = offset + pfc.particleId();
122  ptx_tot += pfc.px();
123  pty_tot += pfc.py();
124  histos.at(offset).pfcandidateType_->Fill(type - offset);
125  auto& histo = histos.at(type);
126  histo.pfcandidateEnergy_->Fill(pfc.energy());
127  histo.pfcandidatePt_->Fill(pfc.pt());
128  histo.pfcandidateEta_->Fill(pfc.eta());
129  histo.pfcandidatePhi_->Fill(pfc.phi());
130  histo.pfcandidateElementsInBlocks_->Fill(pfc.elementsInBlocks().size());
131  }
132  auto& histo = histos.at(offset);
133  histo.pfcandidate_vect_sum_pt_->Fill(std::sqrt(ptx_tot * ptx_tot + pty_tot * pty_tot));
134 }
135 
137  edm::Run const& run,
138  edm::EventSetup const& iSetup,
140  int offset = 100000;
141  ibook.setCurrentFolder(folder_ + "PFCandidates");
142  histos[offset].pfcandidateType_ = ibook.book1D("PFCandidateType", "PFCandidateType", 10, 0, 10);
143  histos[offset].pfcandidate_vect_sum_pt_ = ibook.book1D("PFCandidatePtVectSum", "PFCandidatePtVectSum", 200, 0., 200.);
145  ibook.setCurrentFolder(folder_ + "PFCandidates/" + std::to_string(type));
146  auto& histo = histos[offset + type];
147  histo.pfcandidateEnergy_ = ibook.book1D("PFCandidateEnergy", "PFCandidateEnergy", 250, 0., 250.);
148  histo.pfcandidatePt_ = ibook.book1D("PFCandidatePt", "PFCandidatePt", 250, 0., 250.);
149  histo.pfcandidateEta_ = ibook.book1D("PFCandidateEta", "PFCandidateEta", 100, -5., 5.);
150  histo.pfcandidatePhi_ = ibook.book1D("PFCandidatePhi", "PFCandidatePhi", 100, -4., 4.);
151  histo.pfcandidateElementsInBlocks_ = ibook.book1D("PFCandidateElements", "PFCandidateElements", 20, 0., 20.);
152  }
153  // Folder '0' is meant to be cumulative, with no connection to pdgId
155  histos[0].simPFSuperClusterSize_ = ibook.book1D("SimPFSuperClusterSize", "SimPFSuperClusterSize", 40, 0., 40.);
156  histos[0].simPFSuperClusterEnergy_ =
157  ibook.book1D("SimPFSuperClusterEnergy", "SimPFSuperClusterEnergy", 250, 0., 500.);
158 }
159 
160 // ------------ method fills 'descriptions' with the allowed parameters for the module ------------
162  //The following says we do not know what parameters are allowed so do no validation
163  // Please change this to state exactly what you do use, even if it is no parameters
165  desc.add<std::string>("folder", "HGCAL/"); // Please keep the trailing '/'
166  desc.add<edm::InputTag>("simPFClusters", edm::InputTag("simPFProducer", "perfect"));
167  desc.add<edm::InputTag>("simPFCandidates", edm::InputTag("simPFProducer"));
168  descriptions.add("caloparticlevalidationDefault", desc);
169 }
170 
171 //define this as a plug-in
std::unordered_map< int, Histogram_CaloParticleSingle > Histograms_CaloParticleValidation
dqm::reco::MonitorElement * simPFSuperClusterSize_
dqm::reco::MonitorElement * pfcandidateType_
CaloParticleValidation(const edm::ParameterSet &)
std::string folder_
virtual void setCurrentFolder(std::string const &fullpath)
Definition: DQMStore.cc:36
std::vector< l1t::PFCandidate > PFCandidateCollection
Definition: PFCandidate.h:86
edm::EDGetTokenT< std::vector< reco::SuperCluster > > simPFClusters_
dqm::reco::MonitorElement * pfcandidatePt_
dqm::reco::MonitorElement * pfcandidateElementsInBlocks_
static std::string to_string(const XMLCh *ch)
void dqmAnalyze(edm::Event const &, edm::EventSetup const &, Histograms_CaloParticleValidation const &) const override
dqm::reco::MonitorElement * pfcandidateEnergy_
int iEvent
Definition: GenABIO.cc:224
edm::EDGetTokenT< reco::PFCandidateCollection > simPFCandidates_
dqm::reco::MonitorElement * pfcandidate_vect_sum_pt_
T sqrt(T t)
Definition: SSEVec.h:23
dqm::reco::MonitorElement * pfcandidateEta_
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
dqm::reco::MonitorElement * simPFSuperClusterEnergy_
void bookHistograms(DQMStore::IBooker &, edm::Run const &, edm::EventSetup const &, Histograms_CaloParticleValidation &) const override
std::vector< reco::PFCandidate > PFCandidateCollection
collection of PFCandidates
void add(std::string const &label, ParameterSetDescription const &psetDescription)
fixed size matrix
HLT enums.
dqm::reco::MonitorElement * pfcandidatePhi_
MonitorElement * book1D(TString const &name, TString const &title, int const nchX, double const lowX, double const highX, FUNC onbooking=NOOP())
Definition: DQMStore.h:98
static void fillDescriptions(edm::ConfigurationDescriptions &descriptions)
Definition: Run.h:45