Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 12:02:16

0001 #include <algorithm>
0002 #include "FWCore/Utilities/interface/Exception.h"
0003 
0004 #include "CondFormats/HcalObjects/interface/HcalPiecewiseLinearFunctor.h"
0005 
0006 inline static double interpolateSimple(
0007     const double x0, const double x1, const double y0, const double y1, const double x) {
0008   return y0 + (y1 - y0) * ((x - x0) / (x1 - x0));
0009 }
0010 
0011 HcalPiecewiseLinearFunctor::HcalPiecewiseLinearFunctor()
0012     : leftExtrapolationLinear_(false), rightExtrapolationLinear_(false) {}
0013 
0014 HcalPiecewiseLinearFunctor::HcalPiecewiseLinearFunctor(const std::vector<std::pair<double, double> >& points,
0015                                                        const bool leftExtrapolationLinear,
0016                                                        const bool rightExtrapolationLinear)
0017     : leftExtrapolationLinear_(leftExtrapolationLinear), rightExtrapolationLinear_(rightExtrapolationLinear) {
0018   const std::size_t sz = points.size();
0019   if (sz) {
0020     std::vector<std::pair<double, double> > tmp(points);
0021     std::sort(tmp.begin(), tmp.end());
0022     abscissae_.reserve(sz);
0023     values_.reserve(sz);
0024     for (std::size_t i = 0; i < sz; ++i) {
0025       abscissae_.push_back(tmp[i].first);
0026       values_.push_back(tmp[i].second);
0027     }
0028     if (!isStrictlyIncreasing(abscissae_.begin(), abscissae_.end()))
0029       throw cms::Exception(
0030           "In HcalPiecewiseLinearFunctor constructor:"
0031           " abscissae must not coincide");
0032   }
0033 }
0034 
0035 double HcalPiecewiseLinearFunctor::operator()(const double x) const {
0036   double result = 0.0;
0037   const std::size_t sz = abscissae_.size();
0038   if (sz) {
0039     if (sz > 1) {
0040       const std::size_t szm1 = sz - 1;
0041       if (x >= abscissae_[szm1]) {
0042         if (rightExtrapolationLinear_)
0043           result = interpolateSimple(abscissae_[sz - 2], abscissae_[szm1], values_[sz - 2], values_[szm1], x);
0044         else
0045           result = values_[szm1];
0046       } else if (x <= abscissae_[0]) {
0047         if (leftExtrapolationLinear_)
0048           result = interpolateSimple(abscissae_[0], abscissae_[1], values_[0], values_[1], x);
0049         else
0050           result = values_[0];
0051       } else {
0052         const std::size_t cell = std::upper_bound(abscissae_.begin(), abscissae_.end(), x) - abscissae_.begin() - 1;
0053         result = interpolateSimple(abscissae_[cell], abscissae_[cell + 1], values_[cell], values_[cell + 1], x);
0054       }
0055     } else
0056       result = values_[0];
0057   }
0058   return result;
0059 }
0060 
0061 bool HcalPiecewiseLinearFunctor::isStrictlyMonotonous() const {
0062   return isStrictlyIncreasing(values_.begin(), values_.end()) || isStrictlyDecreasing(values_.begin(), values_.end());
0063 }
0064 
0065 HcalPiecewiseLinearFunctor HcalPiecewiseLinearFunctor::inverse() const {
0066   if (!isStrictlyMonotonous())
0067     throw cms::Exception(
0068         "In HcalPiecewiseLinearFunctor::inverse:"
0069         " can't invert non-monotonous functor");
0070   const std::size_t sz = abscissae_.size();
0071   std::vector<std::pair<double, double> > points;
0072   points.reserve(sz);
0073   for (std::size_t i = 0; i < sz; ++i)
0074     points.push_back(std::make_pair(values_[i], abscissae_[i]));
0075   bool l = leftExtrapolationLinear_;
0076   bool r = rightExtrapolationLinear_;
0077   if (values_[0] > values_[sz - 1])
0078     std::swap(l, r);
0079   return HcalPiecewiseLinearFunctor(points, l, r);
0080 }
0081 
0082 double HcalPiecewiseLinearFunctor::xmin() const {
0083   double result = 0.0;
0084   if (!abscissae_.empty())
0085     result = abscissae_[0];
0086   return result;
0087 }
0088 
0089 double HcalPiecewiseLinearFunctor::xmax() const {
0090   double result = 0.0;
0091   if (!abscissae_.empty())
0092     result = abscissae_.back();
0093   return result;
0094 }
0095 
0096 BOOST_CLASS_EXPORT_IMPLEMENT(HcalPiecewiseLinearFunctor)