CMS 3D CMS Logo

PFMETFilter.cc
Go to the documentation of this file.
2 
11 
13  collections_ = iConfig.getParameter<std::vector<std::string>>("Collections");
14  variables_ = iConfig.getParameter<std::vector<std::string>>("Variables");
15  min_ = iConfig.getParameter<std::vector<double>>("Mins");
16  max_ = iConfig.getParameter<std::vector<double>>("Maxs");
17  doMin_ = iConfig.getParameter<std::vector<int>>("DoMin");
18  doMax_ = iConfig.getParameter<std::vector<int>>("DoMax");
19  TrueMET_ = iConfig.getParameter<std::string>("TrueMET");
20  DeltaMEXsigma_ = iConfig.getParameter<double>("DeltaMEXsigma");
21  sigma_a_ = iConfig.getParameter<double>("sigma_a");
22  sigma_b_ = iConfig.getParameter<double>("sigma_b");
23  sigma_c_ = iConfig.getParameter<double>("sigma_c");
24  verbose_ = iConfig.getParameter<bool>("verbose");
25 }
26 
28 
30  if (collections_.size() != min_.size()) {
31  std::cout << "Error: in PFMETFilter: collections_.size()!=min_.size()" << std::endl;
32  std::cout << "collections_.size() = " << collections_.size() << std::endl;
33  std::cout << "min_.size() = " << min_.size() << std::endl;
34  return false;
35  }
36  if (collections_.size() != max_.size()) {
37  std::cout << "Error: in PFMETFilter: collections_.size()!=max_.size()" << std::endl;
38  std::cout << "collections_.size() = " << collections_.size() << std::endl;
39  std::cout << "max_.size() = " << max_.size() << std::endl;
40  return false;
41  }
42  if (collections_.size() != doMin_.size()) {
43  std::cout << "Error: in PFMETFilter: collections_.size()!=min_.size()" << std::endl;
44  std::cout << "collections_.size() = " << collections_.size() << std::endl;
45  std::cout << "doMin_.size() = " << doMin_.size() << std::endl;
46  return false;
47  }
48  if (collections_.size() != doMax_.size()) {
49  std::cout << "Error: in PFMETFilter: collections_.size()!=min_.size()" << std::endl;
50  std::cout << "collections_.size() = " << collections_.size() << std::endl;
51  std::cout << "doMax_.size() = " << doMax_.size() << std::endl;
52  return false;
53  }
54  if (collections_.size() != variables_.size()) {
55  std::cout << "Error: in PFMETFilter: collections_.size()!=variables_.size()" << std::endl;
56  std::cout << "collections_.size() = " << collections_.size() << std::endl;
57  std::cout << "variables_.size() = " << variables_.size() << std::endl;
58  return false;
59  }
60  return true;
61 }
62 
64  // std::cout << "FL: beginJob" << std::endl;
65 }
66 
68  // std::cout << "FL: filter" << std::endl;
69  // std::cout << "FL: Mins = " << min_ << std::endl;
70 
71  if (!checkInput())
72  return true; // no filtering !
73 
74  bool skip = false;
75 
76  for (unsigned int varc = 0; varc < collections_.size(); ++varc) {
77  // std::cout << "FL: var[" << varc << "] = " << collections_[varc] <<
78  // std::endl; std::cout << "FL: var[0] = " << collections_[0] << std::endl;
79 
80  // if the collection is collection1-collection2:
81  const unsigned int minuspos = collections_[varc].find("-");
82  if (minuspos < collections_[varc].size()) {
83  std::string collection1;
84  collection1.assign(collections_[varc], 0, minuspos);
85  // std::cout << "collection1 = " << collection1 << std::endl;
86  std::string collection2;
87  collection2.assign(collections_[varc], minuspos + 1, collections_[varc].size());
88  // std::cout << "collection2 = " << collection2 << std::endl;
89 
90  const edm::View<reco::Candidate> *var1;
92  bool isVar1 = iEvent.getByLabel(collection1, var1_hnd);
93  if (!isVar1) {
94  std::cout << "Warning : no " << collection1 << " in input !" << std::endl;
95  return false;
96  }
97  var1 = var1_hnd.product();
98  const reco::Candidate *var10 = &(*var1)[0];
99  // std::cout << "FL: var10.pt = " << var10->et() << std::endl;
100  double coll_var1;
101  if (variables_[varc] == "et")
102  coll_var1 = var10->et();
103  else if (variables_[varc] == "phi")
104  coll_var1 = var10->phi();
105  else if (variables_[varc] == "eta")
106  coll_var1 = var10->eta();
107  else {
108  std::cout << "Error: PFMETFilter: variable unknown: " << variables_[varc] << std::endl;
109  return true;
110  }
111  // std::cout << "FL: coll_var1[" << variables_[varc] << "] = " <<
112  // coll_var1 << std::endl;
113 
114  const edm::View<reco::Candidate> *var2;
116  bool isVar2 = iEvent.getByLabel(collection2, var2_hnd);
117  if (!isVar2) {
118  std::cout << "Warning : no " << collection2 << " in input !" << std::endl;
119  return false;
120  }
121  var2 = var2_hnd.product();
122  const reco::Candidate *var20 = &(*var2)[0];
123  // std::cout << "FL: var20.pt = " << var20->et() << std::endl;
124  double coll_var2;
125  if (variables_[varc] == "et")
126  coll_var2 = var20->et();
127  else if (variables_[varc] == "phi")
128  coll_var2 = var20->phi();
129  else if (variables_[varc] == "eta")
130  coll_var2 = var20->eta();
131  else {
132  std::cout << "Error: PFMETFilter: variable unknown: " << variables_[varc] << std::endl;
133  return true;
134  }
135  // std::cout << "FL: coll_var2[" << variables_[varc] << "] = " <<
136  // coll_var2 << std::endl; std::cout << "FL: max_[varc] = " << max_[varc]
137  // << std::endl; std::cout << "FL: min_[varc] = " << min_[varc] <<
138  // std::endl;
139 
140  // Delta computation
141  double delta = coll_var1 - coll_var2;
142  if (variables_[varc] == "phi") {
143  if (coll_var1 > M_PI)
144  coll_var1 -= ceil((coll_var1 - M_PI) / (2 * M_PI)) * 2 * M_PI;
145  if (coll_var1 <= -M_PI)
146  coll_var1 += ceil((coll_var1 + M_PI) / (-2. * M_PI)) * 2. * M_PI;
147  if (coll_var2 > M_PI)
148  coll_var2 -= ceil((coll_var2 - M_PI) / (2 * M_PI)) * 2 * M_PI;
149  if (coll_var2 <= -M_PI)
150  coll_var2 += ceil((coll_var2 + M_PI) / (-2. * M_PI)) * 2 * M_PI;
151 
152  double deltaphi = -999.0;
153  if (fabs(coll_var1 - coll_var2) < M_PI) {
154  deltaphi = (coll_var1 - coll_var2);
155  } else {
156  if ((coll_var1 - coll_var2) > 0.0) {
157  deltaphi = (2 * M_PI - fabs(coll_var1 - coll_var2));
158  } else {
159  deltaphi = -(2 * M_PI - fabs(coll_var1 - coll_var2));
160  }
161  }
162  delta = deltaphi;
163  }
164 
165  // cuts
166  if (doMin_[varc] && doMax_[varc] && max_[varc] < min_[varc]) {
167  if (delta > max_[varc] && delta < min_[varc])
168  skip = true;
169  } else {
170  if (doMin_[varc] && delta < min_[varc])
171  skip = true;
172  if (doMax_[varc] && delta > max_[varc])
173  skip = true;
174  }
175  // std::cout << "skip = " << skip << std::endl;
176  } else {
177  // get the variable:
178  const edm::View<reco::Candidate> *var0;
180  bool isVar0 = iEvent.getByLabel(collections_[varc], var0_hnd);
181  if (!isVar0) {
182  std::cout << "Warning : no " << collections_[varc] << " in input !" << std::endl;
183  return false;
184  }
185  var0 = var0_hnd.product();
186  const reco::Candidate *var00 = &(*var0)[0];
187  // std::cout << "FL: var00.pt = " << var00->et() << std::endl;
188  double coll_var;
189  if (variables_[varc] == "et")
190  coll_var = var00->et();
191  else if (variables_[varc] == "phi")
192  coll_var = var00->phi();
193  else if (variables_[varc] == "eta")
194  coll_var = var00->eta();
195  else if (variables_[varc] == "DeltaMEXcut") {
196  const edm::View<reco::Candidate> *truevar0;
198  bool istrueVar0 = iEvent.getByLabel(TrueMET_, truevar0_hnd);
199  if (!istrueVar0) {
200  std::cout << "Warning : no " << TrueMET_ << " in input !" << std::endl;
201  return false;
202  }
203  truevar0 = truevar0_hnd.product();
204  const reco::Candidate *truevar00 = &(*truevar0)[0];
205 
206  const double DeltaMEX = var00->px() - truevar00->px();
207  const double DeltaMEY = var00->py() - truevar00->py();
208  const double cutvalc = sqrt(DeltaMEX * DeltaMEX + DeltaMEY * DeltaMEY);
209  const reco::MET *met = static_cast<const reco::MET *>(truevar00);
210  const double SETc = met->sumEt();
211  // std::cout << "FL: SETc = " << SETc << std::endl;
212  const double sigmac = sigma_a_ + sigma_b_ * sqrt(SETc) + sigma_c_ * SETc;
213  if (cutvalc > DeltaMEXsigma_ * sigmac) {
214  if (verbose_) {
215  std::cout << "DeltaMET = " << var00->et() - truevar00->et() << std::endl;
216  std::cout << "trueSET = " << SETc << std::endl;
217  std::cout << "pfMET = " << var00->et() << std::endl;
218  std::cout << "trueMET = " << truevar00->et() << std::endl;
219  std::cout << "DeltaMEX = " << DeltaMEX << std::endl;
220  std::cout << "DeltaMEY = " << DeltaMEY << std::endl;
221  std::cout << "cutvalc = " << cutvalc << std::endl;
222  std::cout << "sigmac = " << sigmac << std::endl;
223  std::cout << "cutvalc/sigmac = " << cutvalc / sigmac << std::endl;
224  }
225  return true;
226  } else {
227  if (verbose_ && (var00->et() - truevar00->et()) > 300.0) {
228  std::cout << "EVENT NOT KEPT:" << std::endl;
229  std::cout << "DeltaMET = " << var00->et() - truevar00->et() << std::endl;
230  std::cout << "SETc = " << SETc << std::endl;
231  std::cout << "pfMET = " << var00->et() << std::endl;
232  std::cout << "trueMET = " << truevar00->et() << std::endl;
233  std::cout << "DeltaMEX = " << DeltaMEX << std::endl;
234  std::cout << "DeltaMEY = " << DeltaMEY << std::endl;
235  std::cout << "cutvalc = " << cutvalc << std::endl;
236  std::cout << "sigmac = " << sigmac << std::endl;
237  std::cout << "cutvalc/sigmac = " << cutvalc / sigmac << std::endl;
238  }
239  return false;
240  }
241  } else {
242  std::cout << "Error: PFMETFilter: variable unknown: " << variables_[varc] << std::endl;
243  return true;
244  }
245  // std::cout << "FL: coll_var[" << variables_[varc] << "] = " << coll_var
246  // << std::endl; std::cout << "FL: max_[varc] = " << max_[varc] <<
247  // std::endl; std::cout << "FL: min_[varc] = " << min_[varc] << std::endl;
248 
249  // cuts
250  if (doMin_[varc] && doMax_[varc] && max_[varc] < min_[varc]) {
251  if (coll_var > max_[varc] && coll_var < min_[varc])
252  skip = true;
253  } else {
254  if (doMin_[varc] && coll_var < min_[varc])
255  skip = true;
256  if (doMax_[varc] && coll_var > max_[varc])
257  skip = true;
258  }
259  // std::cout << "skip = " << skip << std::endl;
260  }
261  }
262  // std::cout << "final skip = " << skip << std::endl;
263  return !skip;
264 }
265 
267  // std::cout << "FL: endJob" << std::endl;
268 }
size
Write out results.
dbl * delta
Definition: mlp_gen.cc:36
T getParameter(std::string const &) const
std::vector< double > max_
Definition: PFMETFilter.h:25
double sigma_c_
Definition: PFMETFilter.h:35
std::vector< double > min_
Definition: PFMETFilter.h:24
PFMETFilter(const edm::ParameterSet &)
Definition: PFMETFilter.cc:12
std::vector< int > doMax_
Definition: PFMETFilter.h:27
bool filter(edm::Event &, const edm::EventSetup &) override
Definition: PFMETFilter.cc:67
virtual double et() const =0
transverse energy
virtual double py() const =0
y coordinate of momentum vector
int iEvent
Definition: GenABIO.cc:224
double sumEt() const
Definition: MET.h:56
double sigma_a_
Definition: PFMETFilter.h:33
Definition: MET.h:42
T sqrt(T t)
Definition: SSEVec.h:18
bool checkInput()
Definition: PFMETFilter.cc:29
std::string TrueMET_
Definition: PFMETFilter.h:31
void beginJob() override
Definition: PFMETFilter.cc:63
~PFMETFilter() override
Definition: PFMETFilter.cc:27
double DeltaMEXsigma_
Definition: PFMETFilter.h:32
bool getByLabel(InputTag const &tag, Handle< PROD > &result) const
Definition: Event.h:480
#define M_PI
virtual double eta() const =0
momentum pseudorapidity
std::vector< std::string > collections_
Definition: PFMETFilter.h:22
T const * product() const
Definition: Handle.h:74
std::vector< int > doMin_
Definition: PFMETFilter.h:26
met
===> hadronic RAZOR
void endJob() override
Definition: PFMETFilter.cc:266
double sigma_b_
Definition: PFMETFilter.h:34
virtual double px() const =0
x coordinate of momentum vector
virtual double phi() const =0
momentum azimuthal angle
std::vector< std::string > variables_
Definition: PFMETFilter.h:23
bool verbose_
Definition: PFMETFilter.h:36