CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
FFTJetParameterParser.cc
Go to the documentation of this file.
1 #include <cstdlib>
2 #include <cstring>
3 #include <string>
4 #include <fstream>
5 
7 
14 
17 
18 #include "fftjet/PeakSelectors.hh"
19 #include "fftjet/Kernel2dFactory.hh"
20 #include "fftjet/GaussianNoiseMembershipFcn.hh"
21 #include "fftjet/EquidistantSequence.hh"
22 #include "fftjet/PeakEtaPhiDistance.hh"
23 #include "fftjet/PeakEtaDependentDistance.hh"
24 #include "fftjet/JetProperty.hh"
25 #include "fftjet/InterpolatedMembershipFcn.hh"
26 #include "fftjet/CompositeKernel.hh"
27 #include "fftjet/InterpolatedKernel.hh"
28 #include "fftjet/InterpolatedKernel3d.hh"
29 #include "fftjet/MagneticSmearingKernel.hh"
30 
31 #define make_param(type, varname) const \
32  type & varname (ps.getParameter< type >( #varname ))
33 
34 using namespace fftjetcms;
35 
36 typedef fftjet::RecombinedJet<VectorLike> RecoFFTJet;
37 
38 static bool parse_peak_member_function(const char* fname,
39  fftjet::JetProperty<fftjet::Peak>::JetMemberFunction *f)
40 {
41  if (strcmp(fname, "eta") == 0)
42  *f = &fftjet::Peak::eta;
43  else if (strcmp(fname, "phi") == 0)
44  *f = &fftjet::Peak::phi;
45  else if (strcmp(fname, "magnitude") == 0)
46  *f = &fftjet::Peak::magnitude;
47  else if (strcmp(fname, "driftSpeed") == 0)
48  *f = &fftjet::Peak::driftSpeed;
49  else if (strcmp(fname, "magSpeed") == 0)
50  *f = &fftjet::Peak::magSpeed;
51  else if (strcmp(fname, "lifetime") == 0)
52  *f = &fftjet::Peak::lifetime;
53  else if (strcmp(fname, "mergeTime") == 0)
54  *f = &fftjet::Peak::mergeTime;
55  else if (strcmp(fname, "splitTime") == 0)
56  *f = &fftjet::Peak::splitTime;
57  else if (strcmp(fname, "scale") == 0)
58  *f = &fftjet::Peak::scale;
59  else if (strcmp(fname, "nearestNeighborDistance") == 0)
60  *f = &fftjet::Peak::nearestNeighborDistance;
61  else if (strcmp(fname, "membershipFactor") == 0)
62  *f = &fftjet::Peak::membershipFactor;
63  else if (strcmp(fname, "recoScale") == 0)
64  *f = &fftjet::Peak::recoScale;
65  else if (strcmp(fname, "recoScaleRatio") == 0)
66  *f = &fftjet::Peak::recoScaleRatio;
67  else if (strcmp(fname, "laplacian") == 0)
68  *f = &fftjet::Peak::laplacian;
69  else if (strcmp(fname, "hessianDeterminant") == 0)
70  *f = &fftjet::Peak::hessianDeterminant;
71  else if (strcmp(fname, "clusterRadius") == 0)
72  *f = &fftjet::Peak::clusterRadius;
73  else if (strcmp(fname, "clusterSeparation") == 0)
74  *f = &fftjet::Peak::clusterSeparation;
75  else
76  {
77  return false;
78  }
79  return true;
80 }
81 
82 static bool parse_jet_member_function(const char* fname,
83  fftjet::JetProperty<RecoFFTJet>::JetMemberFunction *f)
84 {
85  if (strcmp(fname, "peakEta") == 0)
86  *f = &RecoFFTJet::peakEta;
87  else if (strcmp(fname, "peakPhi") == 0)
88  *f = &RecoFFTJet::peakPhi;
89  else if (strcmp(fname, "peakMagnitude") == 0)
90  *f = &RecoFFTJet::peakMagnitude;
91  else if (strcmp(fname, "magnitude") == 0)
92  *f = &RecoFFTJet::magnitude;
93  else if (strcmp(fname, "driftSpeed") == 0)
94  *f = &RecoFFTJet::driftSpeed;
95  else if (strcmp(fname, "magSpeed") == 0)
96  *f = &RecoFFTJet::magSpeed;
97  else if (strcmp(fname, "lifetime") == 0)
98  *f = &RecoFFTJet::lifetime;
99  else if (strcmp(fname, "mergeTime") == 0)
100  *f = &RecoFFTJet::mergeTime;
101  else if (strcmp(fname, "splitTime") == 0)
102  *f = &RecoFFTJet::splitTime;
103  else if (strcmp(fname, "scale") == 0)
104  *f = &RecoFFTJet::scale;
105  else if (strcmp(fname, "nearestNeighborDistance") == 0)
106  *f = &RecoFFTJet::nearestNeighborDistance;
107  else if (strcmp(fname, "membershipFactor") == 0)
108  *f = &RecoFFTJet::membershipFactor;
109  else if (strcmp(fname, "recoScale") == 0)
110  *f = &RecoFFTJet::recoScale;
111  else if (strcmp(fname, "recoScaleRatio") == 0)
112  *f = &RecoFFTJet::recoScaleRatio;
113  else if (strcmp(fname, "laplacian") == 0)
114  *f = &RecoFFTJet::laplacian;
115  else if (strcmp(fname, "hessianDeterminant") == 0)
116  *f = &RecoFFTJet::hessianDeterminant;
117  else if (strcmp(fname, "clusterRadius") == 0)
118  *f = &RecoFFTJet::clusterRadius;
119  else if (strcmp(fname, "clusterSeparation") == 0)
120  *f = &RecoFFTJet::clusterSeparation;
121  else
122  {
123  return false;
124  }
125  return true;
126 }
127 
128 
129 namespace fftjetcms {
130 
131 std::auto_ptr<fftjet::Grid2d<Real> >
133 {
134  typedef std::auto_ptr<fftjet::Grid2d<Real> > return_type;
135  fftjet::Grid2d<Real> *g = 0;
136 
137  // Check if the grid should be read from file
138  if (ps.exists("file"))
139  {
140  const std::string file = ps.getParameter<std::string>("file");
141  std::ifstream in(file.c_str(),
142  std::ios_base::in | std::ios_base::binary);
143  if (!in.is_open())
144  throw cms::Exception("FFTJetBadConfig")
145  << "Failed to open file " << file << std::endl;
147  if (g == 0)
148  throw cms::Exception("FFTJetBadConfig")
149  << "Failed to read file " << file << std::endl;
150  }
151  else
152  {
153  const unsigned nEtaBins = ps.getParameter<unsigned>("nEtaBins");
154  const Real etaMin = ps.getParameter<Real>("etaMin");
155  const Real etaMax = ps.getParameter<Real>("etaMax");
156  const unsigned nPhiBins = ps.getParameter<unsigned>("nPhiBins");
157  const Real phiBin0Edge = ps.getParameter<Real>("phiBin0Edge");
159  "title", "");
160 
161  if (nEtaBins == 0 || nPhiBins == 0 || etaMin >= etaMax)
162  return return_type(NULL);
163 
164  g = new fftjet::Grid2d<Real>(
165  nEtaBins,
166  etaMin,
167  etaMax,
168  nPhiBins,
169  phiBin0Edge,
170  title.c_str()
171  );
172 
173  // Check if the grid data is provided
174  if (ps.exists("data"))
175  {
176  const std::vector<Real>& data =
177  ps.getParameter<std::vector<Real> >("data");
178  if (data.size() == nEtaBins*nPhiBins)
179  g->blockSet(&data[0], nEtaBins, nPhiBins);
180  else
181  {
182  delete g;
183  g = 0;
184  }
185  }
186  }
187 
188  return return_type(g);
189 }
190 
191 
192 std::auto_ptr<fftjet::Functor1<bool,fftjet::Peak> >
194 {
195  typedef std::auto_ptr<fftjet::Functor1<bool,fftjet::Peak> > return_type;
196 
197  const std::string peakselector_type = ps.getParameter<std::string>(
198  "Class");
199 
200  if (!peakselector_type.compare("AllPeaksPass"))
201  {
202  return return_type(new fftjet::AllPeaksPass());
203  }
204 
205  if (!peakselector_type.compare("EtaAndPtDependentPeakSelector"))
206  {
207  const std::string file = ps.getParameter<std::string>("file");
208  std::ifstream in(file.c_str(),
209  std::ios_base::in | std::ios_base::binary);
210  if (!in.is_open())
211  throw cms::Exception("FFTJetBadConfig")
212  << "Failed to open file " << file << std::endl;
215  if (!ptr->isValid())
216  throw cms::Exception("FFTJetBadConfig")
217  << "Failed to read file " << file << std::endl;
218  return return_type(ptr);
219  }
220 
221  if (!peakselector_type.compare("EtaAndPtLookupPeakSelector"))
222  {
223  make_param(unsigned, nx);
224  make_param(unsigned, ny);
225  make_param(double, xmin);
226  make_param(double, xmax);
227  make_param(double, ymin);
228  make_param(double, ymax);
229  make_param(std::vector<double>, data);
230 
231  if (xmin >= xmax || ymin >= ymax || !nx || !ny || data.size() != nx*ny)
232  throw cms::Exception("FFTJetBadConfig")
233  << "Failed to configure EtaAndPtLookupPeakSelector" << std::endl;
234 
235  return return_type(new EtaAndPtLookupPeakSelector(
236  nx, xmin, xmax, ny, ymin, ymax, data));
237  }
238 
239  if (!peakselector_type.compare("SimplePeakSelector"))
240  {
241  const double magCut = ps.getParameter<double>("magCut");
242  const double driftSpeedCut = ps.getParameter<double>("driftSpeedCut");
243  const double magSpeedCut = ps.getParameter<double>("magSpeedCut");
244  const double lifeTimeCut = ps.getParameter<double>("lifeTimeCut");
245  const double NNDCut = ps.getParameter<double>("NNDCut");
246  const double etaCut = ps.getParameter<double>("etaCut");
247  const double splitTimeCut = ps.getParameter<double>("splitTimeCut");
248  const double mergeTimeCut = ps.getParameter<double>("mergeTimeCut");
249 
250  return return_type(new fftjet::SimplePeakSelector(
251  magCut, driftSpeedCut, magSpeedCut, lifeTimeCut, NNDCut,
252  etaCut, splitTimeCut, mergeTimeCut));
253  }
254 
255  if (!peakselector_type.compare("ScalePowerPeakSelector"))
256  {
257  const double a = ps.getParameter<double>("a");
258  const double p = ps.getParameter<double>("p");
259  const double b = ps.getParameter<double>("b");
260  const double etaCut = ps.getParameter<double>("etaCut");
261 
262  return return_type(new fftjet::ScalePowerPeakSelector(
263  a, p, b, etaCut));
264  }
265 
266  return return_type(NULL);
267 }
268 
269 
270 std::auto_ptr<fftjet::ScaleSpaceKernel>
272 {
273  typedef std::auto_ptr<fftjet::ScaleSpaceKernel> return_type;
274 
275  const std::string MembershipFunction_type = ps.getParameter<std::string>(
276  "Class");
277 
278  // Parse special cases first
279  if (!MembershipFunction_type.compare("InterpolatedMembershipFcn"))
280  {
281  // This is a kernel defined by a 4d (sparsified) lookup table.
282  // Here, it is simply loaded from a file using a built-in
283  // method from fftjet. Note that the table representation
284  // must be native binary (this will not work on platforms with
285  // different endianity of floating point standard).
286  const std::string file = ps.getParameter<std::string>("file");
287  std::ifstream in(file.c_str(),
288  std::ios_base::in | std::ios_base::binary);
289  if (!in.is_open())
290  throw cms::Exception("FFTJetBadConfig")
291  << "Failed to open file " << file << std::endl;
293  }
294 
295  if (!MembershipFunction_type.compare("Composite"))
296  {
297  throw cms::Exception("FFTJetBadConfig")
298  << "Parsing of CompositeKernel objects is not implemented yet"
299  << std::endl;
300  }
301 
302  if (!MembershipFunction_type.compare("MagneticSmearing"))
303  {
304  // This kernel represents smearing of a jet in phi
305  // in a magnetic field. The meaning of the parameters
306  // is explained in the comments in the MagneticSmearingKernel.hh
307  // header file of the fftjet package.
308  make_param(std::vector<double>, fragmentationData);
309  make_param(double, numeratorConst);
310  make_param(double, charge1Fraction);
311  make_param(double, charge2Fraction);
312  make_param(unsigned, samplesPerBin);
313 
314  if (fragmentationData.empty())
315  throw cms::Exception("FFTJetBadConfig")
316  << "Fragmentation function data not defined for "
317  "MagneticSmearingKernel" << std::endl;
318  if (samplesPerBin < 1U)
319  throw cms::Exception("FFTJetBadConfig")
320  << "Bad number of samples per bin in "
321  "MagneticSmearingKernel" << std::endl;
322 
323  fftjet::LinearInterpolator1d* fragmentationFunction =
324  new fftjet::LinearInterpolator1d(
325  &fragmentationData[0], fragmentationData.size(), 0.0, 1.0);
326 
327  return return_type(
328  new fftjet::MagneticSmearingKernel<fftjet::LinearInterpolator1d>(
329  fragmentationFunction, numeratorConst,
330  charge1Fraction, charge2Fraction,
331  samplesPerBin, true));
332  }
333 
334  if (!MembershipFunction_type.compare("Interpolated"))
335  {
336  // This is a kernel defined by a histogram-like 2d lookup table
337  make_param(double, sx);
338  make_param(double, sy);
339  make_param(int, scalePower);
340  make_param(unsigned, nxbins);
341  make_param(double, xmin);
342  make_param(double, xmax);
343  make_param(unsigned, nybins);
344  make_param(double, ymin);
345  make_param(double, ymax);
346  make_param(std::vector<double>, data);
347 
348  if (data.size() != nxbins*nybins)
349  throw cms::Exception("FFTJetBadConfig")
350  << "Bad number of data points for Interpolated kernel"
351  << std::endl;
352 
353  return return_type(new fftjet::InterpolatedKernel(
354  sx, sy, scalePower,
355  &data[0], nxbins, xmin, xmax,
356  nybins, ymin, ymax));
357  }
358 
359  if (!MembershipFunction_type.compare("Interpolated3d"))
360  {
361  // This is a kernel defined by a histogram-like 3d lookup table
362  make_param(std::vector<double>, data);
363  make_param(std::vector<double>, scales);
364  make_param(bool, useLogSpaceForScale);
365  make_param(unsigned, nxbins);
366  make_param(double, xmin);
367  make_param(double, xmax);
368  make_param(unsigned, nybins);
369  make_param(double, ymin);
370  make_param(double, ymax);
371 
372  if (data.size() != nxbins*nybins*scales.size())
373  throw cms::Exception("FFTJetBadConfig")
374  << "Bad number of data points for Interpolated3d kernel"
375  << std::endl;
376 
377  return return_type(new fftjet::InterpolatedKernel3d(
378  &data[0], scales, useLogSpaceForScale,
379  nxbins, xmin, xmax, nybins, ymin, ymax));
380  }
381 
382  // This is not a special kernel. Try one of the classes
383  // in the kernel factory provided by FFTJet.
384  fftjet::DefaultKernel2dFactory factory;
385  if (factory[MembershipFunction_type] == NULL) {
386  return return_type(NULL);
387  }
388 
389  make_param(double, sx);
390  make_param(double, sy);
391  make_param(int, scalePower);
392  make_param(std::vector<double>, kernelParameters);
393 
394  const int n_expected = factory[MembershipFunction_type]->nParameters();
395  if (n_expected >= 0)
396  if (static_cast<unsigned>(n_expected) != kernelParameters.size())
397  throw cms::Exception("FFTJetBadConfig")
398  << "Bad number of kernel parameters" << std::endl;
399 
400  return std::auto_ptr<fftjet::ScaleSpaceKernel>(
401  factory[MembershipFunction_type]->create(
402  sx, sy, scalePower, kernelParameters));
403 }
404 
405 
406 std::auto_ptr<AbsBgFunctor>
408 {
409  const std::string bg_Membership_type = ps.getParameter<std::string>(
410  "Class");
411 
412  if (!bg_Membership_type.compare("GaussianNoiseMembershipFcn"))
413  {
414  const double minWeight = ps.getParameter<double>("minWeight");
415  const double prior = ps.getParameter<double>("prior");
416  return std::auto_ptr<AbsBgFunctor>(
417  new fftjet::GaussianNoiseMembershipFcn(minWeight,prior));
418  }
419 
420  return std::auto_ptr<AbsBgFunctor>(NULL);
421 }
422 
423 
424 std::auto_ptr<std::vector<double> >
426 {
427  typedef std::auto_ptr<std::vector<double> > return_type;
428 
429  const std::string className = ps.getParameter<std::string>("Class");
430 
431  if (!className.compare("EquidistantInLinearSpace") ||
432  !className.compare("EquidistantInLogSpace"))
433  {
434  const double minScale = ps.getParameter<double>("minScale");
435  const double maxScale = ps.getParameter<double>("maxScale");
436  const unsigned nScales = ps.getParameter<unsigned>("nScales");
437 
438  if (minScale <= 0.0 || maxScale <= 0.0 ||
439  nScales == 0 || minScale == maxScale)
440  return return_type(NULL);
441 
442  // Can't return pointers to EquidistantInLinearSpace
443  // or EquidistantInLogSpace directly because std::vector
444  // destructor is not virtual.
445  if (!className.compare("EquidistantInLinearSpace"))
446  return return_type(new std::vector<double>(
447  fftjet::EquidistantInLinearSpace(
448  minScale, maxScale, nScales)));
449  else
450  return return_type(new std::vector<double>(
451  fftjet::EquidistantInLogSpace(
452  minScale, maxScale, nScales)));
453  }
454 
455  if (!className.compare("UserSet"))
456  {
457  return_type scales(new std::vector<double>(
458  ps.getParameter<std::vector<double> >("scales")));
459 
460  // Verify that all scales are positive and unique
461  const unsigned nscales = scales->size();
462  for (unsigned i=0; i<nscales; ++i)
463  if ((*scales)[i] <= 0.0)
464  return return_type(NULL);
465 
466  for (unsigned i=1; i<nscales; ++i)
467  for (unsigned j=0; j<i; ++j)
468  if ((*scales)[i] == (*scales)[j])
469  return return_type(NULL);
470 
471  return scales;
472  }
473 
474  return return_type(NULL);
475 }
476 
477 
478 std::auto_ptr<fftjet::ClusteringTreeSparsifier<fftjet::Peak,long> >
480 {
481  typedef std::auto_ptr<fftjet::ClusteringTreeSparsifier<fftjet::Peak,long> > return_type;
482 
483  const int maxLevelNumber = ps.getParameter<int>("maxLevelNumber");
484  const unsigned filterMask = ps.getParameter<unsigned>("filterMask");
485  const std::vector<double> userScalesV =
486  ps.getParameter<std::vector<double> >("userScales");
487  const unsigned nUserScales = userScalesV.size();
488  const double* userScales = nUserScales ? &userScalesV[0] : NULL;
489 
490  return return_type(
491  new fftjet::ClusteringTreeSparsifier<fftjet::Peak,long>(
492  maxLevelNumber,
493  filterMask,
494  userScales,
495  nUserScales
496  )
497  );
498 }
499 
500 
501 std::auto_ptr<fftjet::AbsDistanceCalculator<fftjet::Peak> >
503 {
504  typedef std::auto_ptr<fftjet::AbsDistanceCalculator<fftjet::Peak> > return_type;
505 
506  const std::string calc_type = ps.getParameter<std::string>("Class");
507 
508  if (!calc_type.compare("PeakEtaPhiDistance"))
509  {
510  const double etaToPhiBandwidthRatio = ps.getParameter<double>(
511  "etaToPhiBandwidthRatio");
512  return return_type(new fftjet::PeakEtaPhiDistance(etaToPhiBandwidthRatio));
513  }
514 
515  if (!calc_type.compare("PeakEtaDependentDistance"))
516  {
517  std::auto_ptr<fftjet::LinearInterpolator1d> interp =
519  ps.getParameter<edm::ParameterSet>("Interpolator"));
520  const fftjet::LinearInterpolator1d* ip = interp.get();
521  if (ip == NULL)
522  return return_type(NULL);
523 
524  // Check that the interpolator is always positive
525  const unsigned n = ip->nx();
526  const double* data = ip->getData();
527  for (unsigned i=0; i<n; ++i)
528  if (data[i] <= 0.0)
529  return return_type(NULL);
530  if (ip->fLow() <= 0.0 || ip->fHigh() <= 0.0)
531  return return_type(NULL);
532 
533  return return_type(new fftjet::PeakEtaDependentDistance(*ip));
534  }
535 
536  return return_type(NULL);
537 }
538 
539 
540 std::auto_ptr<fftjetcms::LinInterpolatedTable1D>
542 {
543  const double xmin = ps.getParameter<double>("xmin");
544  const double xmax = ps.getParameter<double>("xmax");
545  const bool leftExtrapolationLinear =
546  ps.getParameter<bool>("leftExtrapolationLinear");
547  const bool rightExtrapolationLinear =
548  ps.getParameter<bool>("rightExtrapolationLinear");
549  const std::vector<double> data(
550  ps.getParameter<std::vector<double> >("data"));
551  if (data.empty())
552  return std::auto_ptr<fftjetcms::LinInterpolatedTable1D>(NULL);
553  else
554  return std::auto_ptr<fftjetcms::LinInterpolatedTable1D>(
556  &data[0], data.size(), xmin, xmax,
557  leftExtrapolationLinear, rightExtrapolationLinear));
558 }
559 
560 
561 std::auto_ptr<fftjet::LinearInterpolator1d>
563 {
564  const double xmin = ps.getParameter<double>("xmin");
565  const double xmax = ps.getParameter<double>("xmax");
566  const double flow = ps.getParameter<double>("flow");
567  const double fhigh = ps.getParameter<double>("fhigh");
568  const std::vector<double> data(
569  ps.getParameter<std::vector<double> >("data"));
570  if (data.empty())
571  return std::auto_ptr<fftjet::LinearInterpolator1d>(NULL);
572  else
573  return std::auto_ptr<fftjet::LinearInterpolator1d>(
574  new fftjet::LinearInterpolator1d(
575  &data[0], data.size(), xmin, xmax, flow, fhigh));
576 }
577 
578 
579 std::auto_ptr<fftjet::LinearInterpolator2d>
581 {
582  const std::string file = ps.getParameter<std::string>("file");
583  std::ifstream in(file.c_str(),
584  std::ios_base::in | std::ios_base::binary);
585  if (!in.is_open())
586  throw cms::Exception("FFTJetBadConfig")
587  << "Failed to open file " << file << std::endl;
588  fftjet::LinearInterpolator2d* ip = fftjet::LinearInterpolator2d::read(in);
589  if (!ip)
590  throw cms::Exception("FFTJetBadConfig")
591  << "Failed to read file " << file << std::endl;
592  return std::auto_ptr<fftjet::LinearInterpolator2d>(ip);
593 }
594 
595 
596 std::auto_ptr<fftjet::Functor1<double,fftjet::Peak> >
598 {
599  typedef fftjet::Functor1<double,fftjet::Peak> ptr_type;
600  typedef std::auto_ptr<ptr_type> return_type;
601 
602  const std::string property_type = ps.getParameter<std::string>("Class");
603 
604  if (!property_type.compare("Log"))
605  {
606  return_type wrapped = fftjet_PeakFunctor_parser(
607  ps.getParameter<edm::ParameterSet>("function"));
608  ptr_type* wr = wrapped.get();
609  if (wr)
610  {
611  return_type result = return_type(
612  new fftjet::LogProperty<fftjet::Peak>(wr, true));
613  wrapped.release();
614  return result;
615  }
616  }
617 
618  if (!property_type.compare("PeakProperty"))
619  {
620  const std::string member = ps.getParameter<std::string>("member");
621  fftjet::JetProperty<fftjet::Peak>::JetMemberFunction fcn;
622  if (parse_peak_member_function(member.c_str(), &fcn))
623  return return_type(
624  new fftjet::JetProperty<fftjet::Peak>(fcn));
625  else
626  return return_type(NULL);
627  }
628 
629  if (!property_type.compare("MinusScaledLaplacian"))
630  {
631  const double sx = ps.getParameter<double>("sx");
632  const double sy = ps.getParameter<double>("sx");
633  return return_type(
634  new fftjet::MinusScaledLaplacian<fftjet::Peak>(sx, sy));
635  }
636 
637  if (!property_type.compare("ScaledHessianDet"))
638  {
639  return return_type(
640  new fftjet::ScaledHessianDet<fftjet::Peak>());
641  }
642 
643  if (!property_type.compare("ScaledMagnitude"))
644  {
645  return return_type(
646  new fftjet::ScaledMagnitude<fftjet::Peak>());
647  }
648 
649  if (!property_type.compare("ScaledMagnitude2"))
650  {
651  return return_type(
652  new fftjet::ScaledMagnitude2<fftjet::Peak>());
653  }
654 
655  if (!property_type.compare("ConstDouble"))
656  {
657  const double value = ps.getParameter<double>("value");
658  return return_type(new ConstDouble<fftjet::Peak>(value));
659  }
660 
661  if (!property_type.compare("ProportionalToScale"))
662  {
663  const double value = ps.getParameter<double>("value");
664  return return_type(
666  }
667 
668  if (!property_type.compare("MultiplyByConst"))
669  {
670  const double factor = ps.getParameter<double>("factor");
671  return_type function = fftjet_PeakFunctor_parser(
672  ps.getParameter<edm::ParameterSet>("function"));
673  ptr_type* ptr = function.get();
674  if (ptr)
675  {
676  return_type result = return_type(
677  new MultiplyByConst<fftjet::Peak>(factor, ptr, true));
678  function.release();
679  return result;
680  }
681  }
682 
683  if (!property_type.compare("CompositeFunctor"))
684  {
685  std::auto_ptr<fftjet::Functor1<double,double> > fcn1 =
687  ps.getParameter<edm::ParameterSet>("function1"));
688  return_type fcn2 = fftjet_PeakFunctor_parser(
689  ps.getParameter<edm::ParameterSet>("function2"));
690  fftjet::Functor1<double,double>* f1 = fcn1.get();
691  ptr_type* f2 = fcn2.get();
692  if (f1 && f2)
693  {
694  return_type result = return_type(
695  new CompositeFunctor<fftjet::Peak>(f1, f2, true));
696  fcn1.release();
697  fcn2.release();
698  return result;
699  }
700  }
701 
702  if (!property_type.compare("ProductFunctor"))
703  {
704  return_type fcn1 = fftjet_PeakFunctor_parser(
705  ps.getParameter<edm::ParameterSet>("function1"));
706  return_type fcn2 = fftjet_PeakFunctor_parser(
707  ps.getParameter<edm::ParameterSet>("function2"));
708  ptr_type* f1 = fcn1.get();
709  ptr_type* f2 = fcn2.get();
710  if (f1 && f2)
711  {
712  return_type result = return_type(
713  new ProductFunctor<fftjet::Peak>(f1, f2, true));
714  fcn1.release();
715  fcn2.release();
716  return result;
717  }
718  }
719 
720  if (!property_type.compare("MagnitudeDependent"))
721  {
722  std::auto_ptr<fftjet::Functor1<double,double> > fcn1 =
724  ps.getParameter<edm::ParameterSet>("function"));
725  fftjet::Functor1<double,double>* f1 = fcn1.get();
726  if (f1)
727  {
728  return_type result = return_type(
729  new MagnitudeDependent<fftjet::Peak>(f1, true));
730  fcn1.release();
731  return result;
732  }
733  }
734 
735  if (!property_type.compare("PeakEtaDependent"))
736  {
737  std::auto_ptr<fftjet::Functor1<double,double> > fcn1 =
739  ps.getParameter<edm::ParameterSet>("function"));
740  fftjet::Functor1<double,double>* f1 = fcn1.get();
741  if (f1)
742  {
743  return_type result = return_type(new PeakEtaDependent(f1, true));
744  fcn1.release();
745  return result;
746  }
747  }
748 
749  return return_type(NULL);
750 }
751 
752 
753 std::auto_ptr<fftjet::Functor1<double,fftjet::RecombinedJet<VectorLike> > >
755 {
756  typedef fftjet::Functor1<double,RecoFFTJet> ptr_type;
757  typedef std::auto_ptr<ptr_type> return_type;
758 
759  const std::string property_type = ps.getParameter<std::string>("Class");
760 
761  if (!property_type.compare("Log"))
762  {
763  return_type wrapped = fftjet_JetFunctor_parser(
764  ps.getParameter<edm::ParameterSet>("function"));
765  fftjet::Functor1<double,RecoFFTJet>* wr = wrapped.get();
766  if (wr)
767  {
768  return_type result = return_type(
769  new fftjet::LogProperty<RecoFFTJet>(wr, true));
770  wrapped.release();
771  return result;
772  }
773  }
774 
775  if (!property_type.compare("JetEtaDependent"))
776  {
777  std::auto_ptr<fftjet::Functor1<double,double> > fcn1 =
779  ps.getParameter<edm::ParameterSet>("function"));
780  fftjet::Functor1<double,double>* f1 = fcn1.get();
781  if (f1)
782  {
783  return_type result = return_type(new JetEtaDependent(f1, true));
784  fcn1.release();
785  return result;
786  }
787  }
788 
789  if (!property_type.compare("JetProperty"))
790  {
791  const std::string member = ps.getParameter<std::string>("member");
792  fftjet::JetProperty<RecoFFTJet>::JetMemberFunction fcn;
793  if (parse_jet_member_function(member.c_str(), &fcn))
794  return return_type(
795  new fftjet::JetProperty<RecoFFTJet>(fcn));
796  else
797  return return_type(NULL);
798  }
799 
800  if (!property_type.compare("ConstDouble"))
801  {
802  const double value = ps.getParameter<double>("value");
803  return return_type(new ConstDouble<RecoFFTJet>(value));
804  }
805 
806  if (!property_type.compare("ProportionalToScale"))
807  {
808  const double value = ps.getParameter<double>("value");
809  return return_type(
811  }
812 
813  if (!property_type.compare("MultiplyByConst"))
814  {
815  const double factor = ps.getParameter<double>("factor");
816  return_type function = fftjet_JetFunctor_parser(
817  ps.getParameter<edm::ParameterSet>("function"));
818  ptr_type* ptr = function.get();
819  if (ptr)
820  {
821  return_type result = return_type(
822  new MultiplyByConst<RecoFFTJet>(factor, ptr, true));
823  function.release();
824  return result;
825  }
826  }
827 
828  if (!property_type.compare("CompositeFunctor"))
829  {
830  std::auto_ptr<fftjet::Functor1<double,double> > fcn1 =
832  ps.getParameter<edm::ParameterSet>("function1"));
833  return_type fcn2 = fftjet_JetFunctor_parser(
834  ps.getParameter<edm::ParameterSet>("function2"));
835  fftjet::Functor1<double,double>* f1 = fcn1.get();
836  ptr_type* f2 = fcn2.get();
837  if (f1 && f2)
838  {
839  return_type result = return_type(
840  new CompositeFunctor<RecoFFTJet>(f1, f2, true));
841  fcn1.release();
842  fcn2.release();
843  return result;
844  }
845  }
846 
847  if (!property_type.compare("ProductFunctor"))
848  {
849  return_type fcn1 = fftjet_JetFunctor_parser(
850  ps.getParameter<edm::ParameterSet>("function1"));
851  return_type fcn2 = fftjet_JetFunctor_parser(
852  ps.getParameter<edm::ParameterSet>("function2"));
853  ptr_type* f1 = fcn1.get();
854  ptr_type* f2 = fcn2.get();
855  if (f1 && f2)
856  {
857  return_type result = return_type(
858  new ProductFunctor<RecoFFTJet>(f1, f2, true));
859  fcn1.release();
860  fcn2.release();
861  return result;
862  }
863  }
864 
865  if (!property_type.compare("MagnitudeDependent"))
866  {
867  std::auto_ptr<fftjet::Functor1<double,double> > fcn1 =
869  ps.getParameter<edm::ParameterSet>("function"));
870  fftjet::Functor1<double,double>* f1 = fcn1.get();
871  if (f1)
872  {
873  return_type result = return_type(
874  new MagnitudeDependent<RecoFFTJet>(f1, true));
875  fcn1.release();
876  return result;
877  }
878  }
879 
880  return return_type(NULL);
881 }
882 
883 
884 std::auto_ptr<fftjet::Functor2<double,
885  fftjet::RecombinedJet<VectorLike>,
886  fftjet::RecombinedJet<VectorLike> > >
888 {
889  typedef std::auto_ptr<fftjet::Functor2<
890  double,
891  fftjet::RecombinedJet<VectorLike>,
892  fftjet::RecombinedJet<VectorLike> > > return_type;
893 
894  const std::string distance_type = ps.getParameter<std::string>(
895  "Class");
896 
897  if (!distance_type.compare("JetConvergenceDistance"))
898  {
899  make_param(double, etaToPhiBandwidthRatio);
900  make_param(double, relativePtBandwidth);
901 
902  if (etaToPhiBandwidthRatio > 0.0 && relativePtBandwidth > 0.0)
903  return return_type(new JetConvergenceDistance(
904  etaToPhiBandwidthRatio, relativePtBandwidth));
905  }
906 
907  return return_type(NULL);
908 }
909 
910 
911 std::auto_ptr<fftjet::Functor1<double,double> >
913 {
914  typedef std::auto_ptr<fftjet::Functor1<double,double> > return_type;
915 
916  const std::string fcn_type = ps.getParameter<std::string>("Class");
917 
918  if (!fcn_type.compare("LinearInterpolator1d"))
919  {
920  std::auto_ptr<fftjet::LinearInterpolator1d> p =
922  fftjet::LinearInterpolator1d* ptr = p.get();
923  if (ptr)
924  {
925  p.release();
926  return return_type(ptr);
927  }
928  }
929 
930  if (!fcn_type.compare("LinInterpolatedTable1D"))
931  {
932  std::auto_ptr<fftjetcms::LinInterpolatedTable1D> p =
934  fftjetcms::LinInterpolatedTable1D* ptr = p.get();
935  if (ptr)
936  {
937  p.release();
938  return return_type(ptr);
939  }
940  }
941 
942  if (!fcn_type.compare("Polynomial"))
943  {
944  std::vector<double> coeffs;
945  for (unsigned i=0; ; ++i)
946  {
947  std::ostringstream s;
948  s << 'c' << i;
949  if (ps.exists(s.str()))
950  coeffs.push_back(ps.getParameter<double>(s.str()));
951  else
952  break;
953  }
954  return return_type(new Polynomial(coeffs));
955  }
956 
957  return return_type(NULL);
958 }
959 
960 
961 std::auto_ptr<AbsPileupCalculator>
963 {
964  typedef std::auto_ptr<AbsPileupCalculator> return_type;
965 
966  const std::string fcn_type = ps.getParameter<std::string>("Class");
967 
968  if (!fcn_type.compare("EtaDependentPileup"))
969  {
970  std::auto_ptr<fftjet::LinearInterpolator2d> interp =
972  ps.getParameter<edm::ParameterSet>("Interpolator2d"));
973  const double inputRhoFactor = ps.getParameter<double>("inputRhoFactor");
974  const double outputRhoFactor = ps.getParameter<double>("outputRhoFactor");
975 
976  const fftjet::LinearInterpolator2d* ip = interp.get();
977  if (ip)
978  return return_type(new EtaDependentPileup(
979  *ip, inputRhoFactor, outputRhoFactor));
980  else
981  return return_type(NULL);
982  }
983 
984  if (!fcn_type.compare("PileupGrid2d"))
985  {
986  std::auto_ptr<fftjet::Grid2d<Real> > grid =
988  ps.getParameter<edm::ParameterSet>("Grid2d"));
989  const double rhoFactor = ps.getParameter<double>("rhoFactor");
990 
991  const fftjet::Grid2d<Real>* g = grid.get();
992  if (g)
993  return return_type(new PileupGrid2d(*g, rhoFactor));
994  else
995  return return_type(NULL);
996  }
997 
998  return return_type(NULL);
999 }
1000 
1001 
1002 std::auto_ptr<fftjet::JetMagnitudeMapper2d <fftjet::Peak> >
1004 {
1005  std::auto_ptr<fftjet::LinearInterpolator2d> responseCurve =
1007  ps.getParameter<edm::ParameterSet>("responseCurve"));
1008 
1009  const double minPredictor = ps.getParameter<double>("minPredictor");
1010  const double maxPredictor = ps.getParameter<double>("maxPredictor");
1011  const unsigned nPredPoints = ps.getParameter<unsigned>("nPredPoints");
1012  const double maxMagnitude = ps.getParameter<double>("maxMagnitude");
1013  const unsigned nMagPoints = ps.getParameter<unsigned>("nMagPoints");
1014 
1015  return (std::auto_ptr<fftjet::JetMagnitudeMapper2d<fftjet::Peak> >
1016  (new fftjet::JetMagnitudeMapper2d<fftjet::Peak>(
1017  *responseCurve,
1019  true,minPredictor,maxPredictor,nPredPoints,
1020  maxMagnitude,nMagPoints)));
1021 }
1022 
1023 
1024 std::auto_ptr<fftjet::JetMagnitudeMapper2d <fftjet::RecombinedJet<VectorLike> > >
1026 {
1027  std::auto_ptr<fftjet::LinearInterpolator2d> responseCurve =
1029  ps.getParameter<edm::ParameterSet>("responseCurve"));
1030 
1031  const double minPredictor = ps.getParameter<double>("minPredictor");
1032  const double maxPredictor = ps.getParameter<double>("maxPredictor");
1033  const unsigned nPredPoints = ps.getParameter<unsigned>("nPredPoints");
1034  const double maxMagnitude = ps.getParameter<double>("maxMagnitude");
1035  const unsigned nMagPoints = ps.getParameter<unsigned>("nMagPoints");
1036 
1037  return (std::auto_ptr<fftjet::JetMagnitudeMapper2d<RecoFFTJet> >
1038  (new fftjet::JetMagnitudeMapper2d<RecoFFTJet>(
1039  *responseCurve,
1041  true,minPredictor,maxPredictor,nPredPoints,
1042  maxMagnitude,nMagPoints)));
1043 }
1044 
1045 }
T getParameter(std::string const &) const
T getUntrackedParameter(std::string const &, T const &) const
int i
Definition: DBlmapReader.cc:9
std::auto_ptr< fftjetcms::LinInterpolatedTable1D > fftjet_LinInterpolatedTable1D_parser(const edm::ParameterSet &ps)
std::auto_ptr< fftjet::Grid2d< Real > > fftjet_Grid2d_parser(const edm::ParameterSet &ps)
std::auto_ptr< fftjet::JetMagnitudeMapper2d< fftjet::Peak > > fftjet_PeakMagnitudeMapper2d_parser(const edm::ParameterSet &ps)
bool exists(std::string const &parameterName) const
checks if a parameter exists
#define NULL
Definition: scimark2.h:8
std::auto_ptr< fftjet::ScaleSpaceKernel > fftjet_MembershipFunction_parser(const edm::ParameterSet &ps)
std::auto_ptr< AbsPileupCalculator > fftjet_PileupCalculator_parser(const edm::ParameterSet &ps)
The Signals That Services Can Subscribe To This is based on ActivityRegistry and is current per Services can connect to the signals distributed by the ActivityRegistry in order to monitor the activity of the application Each possible callback has some defined which we here list in angle e g
Definition: Activities.doc:4
tuple result
Definition: mps_fire.py:95
fftjet::RecombinedJet< VectorLike > RecoFFTJet
std::auto_ptr< fftjet::Functor1< bool, fftjet::Peak > > fftjet_PeakSelector_parser(const edm::ParameterSet &ps)
std::auto_ptr< fftjet::Functor2< double, fftjet::RecombinedJet< VectorLike >, fftjet::RecombinedJet< VectorLike > > > fftjet_JetDistance_parser(const edm::ParameterSet &ps)
std::auto_ptr< fftjet::Functor1< double, double > > fftjet_Function_parser(const edm::ParameterSet &ps)
static bool parse_peak_member_function(const char *fname, fftjet::JetProperty< fftjet::Peak >::JetMemberFunction *f)
int j
Definition: DBlmapReader.cc:9
double f[11][100]
#define make_param(type, varname)
std::auto_ptr< fftjet::Functor1< double, fftjet::Peak > > fftjet_PeakFunctor_parser(const edm::ParameterSet &ps)
double Real
std::auto_ptr< fftjet::JetMagnitudeMapper2d< fftjet::RecombinedJet< VectorLike > > > fftjet_JetMagnitudeMapper2d_parser(const edm::ParameterSet &ps)
void fcn(int &, double *, double &, double *, int)
std::auto_ptr< fftjet::LinearInterpolator1d > fftjet_LinearInterpolator1d_parser(const edm::ParameterSet &ps)
double b
Definition: hdecay.h:120
Geom::Phi< T > phi() const
dictionary prior
string fname
main script
std::auto_ptr< fftjet::AbsDistanceCalculator< fftjet::Peak > > fftjet_DistanceCalculator_parser(const edm::ParameterSet &ps)
double a
Definition: hdecay.h:121
std::auto_ptr< fftjet::ClusteringTreeSparsifier< fftjet::Peak, long > > fftjet_ClusteringTreeSparsifier_parser(const edm::ParameterSet &ps)
static bool parse_jet_member_function(const char *fname, fftjet::JetProperty< RecoFFTJet >::JetMemberFunction *f)
std::auto_ptr< AbsBgFunctor > fftjet_BgFunctor_parser(const edm::ParameterSet &ps)
std::auto_ptr< std::vector< double > > fftjet_ScaleSet_parser(const edm::ParameterSet &ps)
std::string className(const T &t)
Definition: ClassName.h:30
std::auto_ptr< fftjet::Functor1< double, fftjet::RecombinedJet< VectorLike > > > fftjet_JetFunctor_parser(const edm::ParameterSet &ps)
std::auto_ptr< fftjet::LinearInterpolator2d > fftjet_LinearInterpolator2d_parser(const edm::ParameterSet &ps)