CMS 3D CMS Logo
PdfSystematicsAnalyzer.cc
Go to the documentation of this file.
1 #include "FWCore/Framework/interface/EDFilter.h"
3 #include "FWCore/Utilities/interface/InputTag.h"
4 #include "DataFormats/Common/interface/TriggerResults.h"
5 
6 class PdfSystematicsAnalyzer : public edm::EDFilter {
7 public:
8  PdfSystematicsAnalyzer(const edm::ParameterSet& pset);
9  ~PdfSystematicsAnalyzer() override;
10  bool filter(edm::Event&, const edm::EventSetup&) override;
11  void beginJob() override;
12  void endJob() override;
13 
14 private:
15  std::string selectorPath_;
16  std::vector<edm::InputTag> pdfWeightTags_;
17  std::vector<edm::EDGetTokenT<std::vector<double> > > pdfWeightTokens_;
18  edm::EDGetTokenT<edm::TriggerResults> triggerResultsToken_;
19  unsigned int originalEvents_;
20  unsigned int selectedEvents_;
21  std::vector<int> pdfStart_;
22  std::vector<double> weightedSelectedEvents_;
23  std::vector<double> weighted2SelectedEvents_;
24  std::vector<double> weightedEvents_;
25 };
26 
28 #include "FWCore/Framework/interface/MakerMacros.h"
29 #include "FWCore/Framework/interface/Frameworkfwd.h"
30 #include "FWCore/Framework/interface/Event.h"
31 #include "FWCore/Framework/interface/EventSetup.h"
32 #include "FWCore/Framework/interface/ESHandle.h"
33 #include "DataFormats/Common/interface/Handle.h"
34 
35 #include "FWCore/Common/interface/TriggerNames.h"
36 
37 #include "FWCore/Utilities/interface/transform.h"
38 
40 PdfSystematicsAnalyzer::PdfSystematicsAnalyzer(const edm::ParameterSet& pset)
41  : selectorPath_(pset.getUntrackedParameter<std::string>("SelectorPath", "")),
42  pdfWeightTags_(pset.getUntrackedParameter<std::vector<edm::InputTag> >("PdfWeightTags")),
43  pdfWeightTokens_(edm::vector_transform(
44  pdfWeightTags_, [this](edm::InputTag const& tag) { return consumes<std::vector<double> >(tag); })),
45  triggerResultsToken_(consumes<edm::TriggerResults>(edm::InputTag("TriggerResults"))) {}
46 
48 PdfSystematicsAnalyzer::~PdfSystematicsAnalyzer() {}
49 
51 void PdfSystematicsAnalyzer::beginJob() {
52  originalEvents_ = 0;
53  selectedEvents_ = 0;
54  edm::LogVerbatim("PDFAnalysis") << "PDF uncertainties will be determined for the following sets: ";
55  for (unsigned int i = 0; i < pdfWeightTags_.size(); ++i) {
56  edm::LogVerbatim("PDFAnalysis") << "\t" << pdfWeightTags_[i].instance();
57  pdfStart_.push_back(-1);
58  }
59 }
60 
62 void PdfSystematicsAnalyzer::endJob() {
63  if (originalEvents_ == 0) {
64  edm::LogVerbatim("PDFAnalysis") << "NO EVENTS => NO RESULTS";
65  return;
66  }
67  if (selectedEvents_ == 0) {
68  edm::LogVerbatim("PDFAnalysis") << "NO SELECTED EVENTS => NO RESULTS";
69  return;
70  }
71 
72  edm::LogVerbatim("PDFAnalysis") << "\n>>>> Begin of PDF weight systematics summary >>>>";
73  edm::LogVerbatim("PDFAnalysis") << "Total number of analyzed data: " << originalEvents_ << " [events]";
74  double originalAcceptance = double(selectedEvents_) / originalEvents_;
75  edm::LogVerbatim("PDFAnalysis") << "Total number of selected data: " << selectedEvents_
76  << " [events], corresponding to acceptance: [" << originalAcceptance * 100 << " +- "
77  << 100 * sqrt(originalAcceptance * (1. - originalAcceptance) / originalEvents_)
78  << "] %";
79 
80  edm::LogVerbatim("PDFAnalysis") << "\n>>>>> PDF UNCERTAINTIES ON RATE >>>>>>";
81  for (unsigned int i = 0; i < pdfWeightTags_.size(); ++i) {
82  bool nnpdfFlag = (pdfWeightTags_[i].instance().substr(0, 5) == "NNPDF");
83  unsigned int nmembers = weightedSelectedEvents_.size() - pdfStart_[i];
84  if (i < pdfWeightTags_.size() - 1)
85  nmembers = pdfStart_[i + 1] - pdfStart_[i];
86  unsigned int npairs = (nmembers - 1) / 2;
87  edm::LogVerbatim("PDFAnalysis") << "RATE Results for PDF set " << pdfWeightTags_[i].instance() << " ---->";
88 
89  double events_central = weightedSelectedEvents_[pdfStart_[i]];
90  edm::LogVerbatim("PDFAnalysis") << "\tEstimate for central PDF member: " << int(events_central) << " [events]";
91  double events2_central = weighted2SelectedEvents_[pdfStart_[i]];
92  edm::LogVerbatim("PDFAnalysis")
93  << "\ti.e. [" << std::setprecision(4) << 100 * (events_central - selectedEvents_) / selectedEvents_ << " +- "
94  << 100 * sqrt(events2_central - events_central + selectedEvents_ * (1 - originalAcceptance)) / selectedEvents_
95  << "] % relative variation with respect to original PDF";
96 
97  if (npairs > 0) {
98  edm::LogVerbatim("PDFAnalysis") << "\tNumber of eigenvectors for uncertainty estimation: " << npairs;
99  double wplus = 0.;
100  double wminus = 0.;
101  unsigned int nplus = 0;
102  unsigned int nminus = 0;
103  for (unsigned int j = 0; j < npairs; ++j) {
104  double wa = weightedSelectedEvents_[pdfStart_[i] + 2 * j + 1] / events_central - 1.;
105  double wb = weightedSelectedEvents_[pdfStart_[i] + 2 * j + 2] / events_central - 1.;
106  if (nnpdfFlag) {
107  if (wa > 0.) {
108  wplus += wa * wa;
109  nplus++;
110  } else {
111  wminus += wa * wa;
112  nminus++;
113  }
114  if (wb > 0.) {
115  wplus += wb * wb;
116  nplus++;
117  } else {
118  wminus += wb * wb;
119  nminus++;
120  }
121  } else {
122  if (wa > wb) {
123  if (wa < 0.)
124  wa = 0.;
125  if (wb > 0.)
126  wb = 0.;
127  wplus += wa * wa;
128  wminus += wb * wb;
129  } else {
130  if (wb < 0.)
131  wb = 0.;
132  if (wa > 0.)
133  wa = 0.;
134  wplus += wb * wb;
135  wminus += wa * wa;
136  }
137  }
138  }
139  if (wplus > 0)
140  wplus = sqrt(wplus);
141  if (wminus > 0)
142  wminus = sqrt(wminus);
143  if (nnpdfFlag) {
144  if (nplus > 0)
145  wplus /= sqrt(nplus);
146  if (nminus > 0)
147  wminus /= sqrt(nminus);
148  }
149  edm::LogVerbatim("PDFAnalysis") << "\tRelative uncertainty with respect to central member: +"
150  << std::setprecision(4) << 100. * wplus << " / -" << std::setprecision(4)
151  << 100. * wminus << " [%]";
152  } else {
153  edm::LogVerbatim("PDFAnalysis") << "\tNO eigenvectors for uncertainty estimation";
154  }
155  }
156 
157  edm::LogVerbatim("PDFAnalysis") << "\n>>>>> PDF UNCERTAINTIES ON ACCEPTANCE >>>>>>";
158  for (unsigned int i = 0; i < pdfWeightTags_.size(); ++i) {
159  bool nnpdfFlag = (pdfWeightTags_[i].instance().substr(0, 5) == "NNPDF");
160  unsigned int nmembers = weightedEvents_.size() - pdfStart_[i];
161  if (i < pdfWeightTags_.size() - 1)
162  nmembers = pdfStart_[i + 1] - pdfStart_[i];
163  unsigned int npairs = (nmembers - 1) / 2;
164  edm::LogVerbatim("PDFAnalysis") << "ACCEPTANCE Results for PDF set " << pdfWeightTags_[i].instance() << " ---->";
165 
166  double acc_central = 0.;
167  double acc2_central = 0.;
168  if (weightedEvents_[pdfStart_[i]] > 0) {
169  acc_central = weightedSelectedEvents_[pdfStart_[i]] / weightedEvents_[pdfStart_[i]];
170  acc2_central = weighted2SelectedEvents_[pdfStart_[i]] / weightedEvents_[pdfStart_[i]];
171  }
172  double waverage = weightedEvents_[pdfStart_[i]] / originalEvents_;
173  edm::LogVerbatim("PDFAnalysis") << "\tEstimate for central PDF member acceptance: [" << acc_central * 100 << " +- "
174  << 100 *
175  sqrt((acc2_central / waverage - acc_central * acc_central) / originalEvents_)
176  << "] %";
177  double xi = acc_central - originalAcceptance;
178  double deltaxi = (acc2_central - (originalAcceptance + 2 * xi + xi * xi)) / originalEvents_;
179  if (deltaxi > 0)
180  deltaxi = sqrt(deltaxi); //else deltaxi = 0.;
181  edm::LogVerbatim("PDFAnalysis") << "\ti.e. [" << std::setprecision(4) << 100 * xi / originalAcceptance << " +- "
182  << std::setprecision(4) << 100 * deltaxi / originalAcceptance
183  << "] % relative variation with respect to the original PDF";
184 
185  if (npairs > 0) {
186  edm::LogVerbatim("PDFAnalysis") << "\tNumber of eigenvectors for uncertainty estimation: " << npairs;
187  double wplus = 0.;
188  double wminus = 0.;
189  unsigned int nplus = 0;
190  unsigned int nminus = 0;
191  for (unsigned int j = 0; j < npairs; ++j) {
192  double wa = 0.;
193  if (weightedEvents_[pdfStart_[i] + 2 * j + 1] > 0)
194  wa = (weightedSelectedEvents_[pdfStart_[i] + 2 * j + 1] / weightedEvents_[pdfStart_[i] + 2 * j + 1]) /
195  acc_central -
196  1.;
197  double wb = 0.;
198  if (weightedEvents_[pdfStart_[i] + 2 * j + 2] > 0)
199  wb = (weightedSelectedEvents_[pdfStart_[i] + 2 * j + 2] / weightedEvents_[pdfStart_[i] + 2 * j + 2]) /
200  acc_central -
201  1.;
202  if (nnpdfFlag) {
203  if (wa > 0.) {
204  wplus += wa * wa;
205  nplus++;
206  } else {
207  wminus += wa * wa;
208  nminus++;
209  }
210  if (wb > 0.) {
211  wplus += wb * wb;
212  nplus++;
213  } else {
214  wminus += wb * wb;
215  nminus++;
216  }
217  } else {
218  if (wa > wb) {
219  if (wa < 0.)
220  wa = 0.;
221  if (wb > 0.)
222  wb = 0.;
223  wplus += wa * wa;
224  wminus += wb * wb;
225  } else {
226  if (wb < 0.)
227  wb = 0.;
228  if (wa > 0.)
229  wa = 0.;
230  wplus += wb * wb;
231  wminus += wa * wa;
232  }
233  }
234  }
235  if (wplus > 0)
236  wplus = sqrt(wplus);
237  if (wminus > 0)
238  wminus = sqrt(wminus);
239  if (nnpdfFlag) {
240  if (nplus > 0)
241  wplus /= sqrt(nplus);
242  if (nminus > 0)
243  wminus /= sqrt(nminus);
244  }
245  edm::LogVerbatim("PDFAnalysis") << "\tRelative uncertainty with respect to central member: +"
246  << std::setprecision(4) << 100. * wplus << " / -" << std::setprecision(4)
247  << 100. * wminus << " [%]";
248  } else {
249  edm::LogVerbatim("PDFAnalysis") << "\tNO eigenvectors for uncertainty estimation";
250  }
251  }
252  edm::LogVerbatim("PDFAnalysis") << ">>>> End of PDF weight systematics summary >>>>";
253 }
254 
256 bool PdfSystematicsAnalyzer::filter(edm::Event& ev, const edm::EventSetup&) {
257  edm::Handle<std::vector<double> > weightHandle;
258  for (unsigned int i = 0; i < pdfWeightTags_.size(); ++i) {
259  if (!ev.getByToken(pdfWeightTokens_[i], weightHandle)) {
260  if (originalEvents_ == 0) {
261  edm::LogError("PDFAnalysis") << ">>> WARNING: some weights not found!";
262  edm::LogError("PDFAnalysis") << ">>> But maybe OK, if you are prefiltering!";
263  edm::LogError("PDFAnalysis") << ">>> If things are OK, this warning should disappear after a while!";
264  }
265  return false;
266  }
267  }
268 
269  originalEvents_++;
270 
271  bool selectedEvent = false;
272  edm::Handle<edm::TriggerResults> triggerResults;
273  if (!ev.getByToken(triggerResultsToken_, triggerResults)) {
274  edm::LogError("PDFAnalysis") << ">>> TRIGGER collection does not exist !!!";
275  return false;
276  }
277 
278  const edm::TriggerNames& trigNames = ev.triggerNames(*triggerResults);
279  unsigned int pathIndex = trigNames.triggerIndex(selectorPath_);
280  bool pathFound = (pathIndex < trigNames.size()); // pathIndex >= 0, since pathIndex is unsigned
281  if (pathFound) {
282  if (triggerResults->accept(pathIndex))
283  selectedEvent = true;
284  }
285  //edm::LogVerbatim("PDFAnalysis") << ">>>> Path Name: " << selectorPath_ << ", selected? " << selectedEvent;
286 
287  if (selectedEvent)
288  selectedEvents_++;
289 
290  for (unsigned int i = 0; i < pdfWeightTags_.size(); ++i) {
291  if (!ev.getByToken(pdfWeightTokens_[i], weightHandle))
292  return false;
293  std::vector<double> weights = (*weightHandle);
294  unsigned int nmembers = weights.size();
295  // Set up arrays the first time wieghts are read
296  if (pdfStart_[i] < 0) {
297  pdfStart_[i] = weightedEvents_.size();
298  for (unsigned int j = 0; j < nmembers; ++j) {
299  weightedEvents_.push_back(0.);
300  weightedSelectedEvents_.push_back(0.);
301  weighted2SelectedEvents_.push_back(0.);
302  }
303  }
304 
305  for (unsigned int j = 0; j < nmembers; ++j) {
306  weightedEvents_[pdfStart_[i] + j] += weights[j];
307  if (selectedEvent) {
308  weightedSelectedEvents_[pdfStart_[i] + j] += weights[j];
309  weighted2SelectedEvents_[pdfStart_[i] + j] += weights[j] * weights[j];
310  }
311  }
312 
313  /*
314  printf("\n>>>>>>>>> Run %8d Event %d, members %3d PDF set %s : Weights >>>> \n", ev.id().run(), ev.id().event(), nmembers, pdfWeightTags_[i].instance().data());
315  for (unsigned int i=0; i<nmembers; i+=5) {
316  for (unsigned int j=0; ((j<5)&&(i+j<nmembers)); ++j) {
317  printf(" %2d: %7.4f", i+j, weights[i+j]);
318  }
319  safe_printf("\n");
320  }
321  */
322  }
323 
324  return true;
325 }
326 
327 DEFINE_FWK_MODULE(PdfSystematicsAnalyzer);
Handle.h
HLT_2018_cff.weights
weights
Definition: HLT_2018_cff.py:87167
mps_fire.i
i
Definition: mps_fire.py:355
PdfSystematicsAnalyzer::beginJob
void beginJob() override
Definition: PdfSystematicsAnalyzer.cc:51
PdfSystematicsAnalyzer::filter
bool filter(edm::Event &, const edm::EventSetup &) override
Definition: PdfSystematicsAnalyzer.cc:256
trigNames
static const char *const trigNames[]
Definition: EcalDumpRaw.cc:57
ESHandle.h
TriggerResults.h
edm::EDGetTokenT< edm::TriggerResults >
edm
HLT enums.
Definition: AlignableModifier.h:19
PdfSystematicsAnalyzer::weighted2SelectedEvents_
std::vector< double > weighted2SelectedEvents_
Definition: PdfSystematicsAnalyzer.cc:23
hybridSuperClusters_cfi.xi
xi
Definition: hybridSuperClusters_cfi.py:10
triggerResults
static const std::string triggerResults
Definition: EdmProvDump.cc:45
EDFilter.h
watchdog.const
const
Definition: watchdog.py:83
edm::Handle
Definition: AssociativeIterator.h:50
PdfSystematicsAnalyzer::pdfWeightTokens_
std::vector< edm::EDGetTokenT< std::vector< double > > > pdfWeightTokens_
Definition: PdfSystematicsAnalyzer.cc:17
MakerMacros.h
DEFINE_FWK_MODULE
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
GlobalPosition_Frontier_DevDB_cff.tag
tag
Definition: GlobalPosition_Frontier_DevDB_cff.py:11
PdfSystematicsAnalyzer::pdfStart_
std::vector< int > pdfStart_
Definition: PdfSystematicsAnalyzer.cc:21
mathSSE::sqrt
T sqrt(T t)
Definition: SSEVec.h:19
PdfSystematicsAnalyzer::originalEvents_
unsigned int originalEvents_
Definition: PdfSystematicsAnalyzer.cc:19
PdfSystematicsAnalyzer::selectorPath_
std::string selectorPath_
Definition: PdfSystematicsAnalyzer.cc:15
PdfSystematicsAnalyzer::~PdfSystematicsAnalyzer
~PdfSystematicsAnalyzer() override
Definition: PdfSystematicsAnalyzer.cc:48
AlCaHLTBitMon_QueryRunRegistry.string
string
Definition: AlCaHLTBitMon_QueryRunRegistry.py:256
edm::vector_transform
auto vector_transform(std::vector< InputType > const &input, Function predicate) -> std::vector< typename std::remove_cv< typename std::remove_reference< decltype(predicate(input.front()))>::type >::type >
Definition: transform.h:11
HLT_2018_cff.InputTag
InputTag
Definition: HLT_2018_cff.py:79016
edm::ParameterSet
Definition: ParameterSet.h:36
edm::LogError
Definition: MessageLogger.h:183
Event.h
createfilelist.int
int
Definition: createfilelist.py:10
TriggerNames.h
PdfSystematicsAnalyzer::endJob
void endJob() override
Definition: PdfSystematicsAnalyzer.cc:62
edm::LogVerbatim
Definition: MessageLogger.h:297
edm::EDFilter
Definition: EDFilter.h:39
PdfSystematicsAnalyzer::weightedEvents_
std::vector< double > weightedEvents_
Definition: PdfSystematicsAnalyzer.cc:24
edm::EventSetup
Definition: EventSetup.h:57
InputTag.h
PdfSystematicsAnalyzer::pdfWeightTags_
std::vector< edm::InputTag > pdfWeightTags_
Definition: PdfSystematicsAnalyzer.cc:16
std
Definition: JetResolutionObject.h:76
Frameworkfwd.h
transform.h
ev
bool ev
Definition: Hydjet2Hadronizer.cc:95
edm::TriggerNames
Definition: TriggerNames.h:55
EventSetup.h
PdfSystematicsAnalyzer::triggerResultsToken_
edm::EDGetTokenT< edm::TriggerResults > triggerResultsToken_
Definition: PdfSystematicsAnalyzer.cc:18
dqmiolumiharvest.j
j
Definition: dqmiolumiharvest.py:66
edm::Event
Definition: Event.h:73
edm::InputTag
Definition: InputTag.h:15
PdfSystematicsAnalyzer::PdfSystematicsAnalyzer
PdfSystematicsAnalyzer(const edm::ParameterSet &pset)
Definition: PdfSystematicsAnalyzer.cc:40
PdfSystematicsAnalyzer::weightedSelectedEvents_
std::vector< double > weightedSelectedEvents_
Definition: PdfSystematicsAnalyzer.cc:22
muonDTDigis_cfi.pset
pset
Definition: muonDTDigis_cfi.py:27
PdfSystematicsAnalyzer::selectedEvents_
unsigned int selectedEvents_
Definition: PdfSystematicsAnalyzer.cc:20
PdfSystematicsAnalyzer
Definition: PdfSystematicsAnalyzer.cc:6
Generated for CMSSW Reference Manual by doxygen 1.8.16