Index: trunk/npstat/nm/FejerQuadrature.icc
===================================================================
--- trunk/npstat/nm/FejerQuadrature.icc	(revision 576)
+++ trunk/npstat/nm/FejerQuadrature.icc	(revision 577)
@@ -1,51 +1,51 @@
 #include <cmath>
 #include <algorithm>
 
 namespace npstat {
     template <typename FcnResult, typename FcnArg>
     inline long double FejerQuadrature::integrate(
         const Functor1<FcnResult,FcnArg>& function,
         const long double left, const long double right) const
     {
         if (buf_.size() != npoints_)
             buf_.resize(npoints_);
         std::pair<long double, long double>* results = &buf_[0];
         const long double midpoint = (left + right)/2.0L;
         const long double unit = (right - left)/2.0L;
         FcnArg a;
         long double v;
         for (unsigned i=0; i<npoints_; ++i)
         {
             a = static_cast<FcnArg>(midpoint + unit*a_[i]);
             v = w_[i]*static_cast<long double>(function(a));
             results[i].first = std::abs(v);
             results[i].second = v;
         }
         std::sort(buf_.begin(), buf_.end());
         v = 0.0L;
         for (unsigned i=0; i<npoints_; ++i)
             v += results[i].second;
         return v*unit;
     }
 
     template <class Functor>
     inline std::vector<std::pair<double,double> >
     FejerQuadrature::weightedIntegrationPoints(
         const Functor& weight,
         const long double left, const long double right) const
     {
         typedef std::pair<double,double> DDPair;
 
         std::vector<DDPair> result;
         result.reserve(npoints_);
         const long double midpoint = (left + right)/2.0L;
         const long double unit = (right - left)/2.0L;
         for (unsigned i=0; i<npoints_; ++i)
         {
             const long double x = midpoint + unit*a_[i];
-            const double w = w_[i]*weight(x);
+            const double w = unit*w_[i]*weight(x);
             result.push_back(DDPair(x, w));
         }
         return result;
     }
 }
Index: trunk/npstat/nm/GaussLegendreQuadrature.icc
===================================================================
--- trunk/npstat/nm/GaussLegendreQuadrature.icc	(revision 576)
+++ trunk/npstat/nm/GaussLegendreQuadrature.icc	(revision 577)
@@ -1,81 +1,81 @@
 #include <cmath>
 #include <stdexcept>
 #include <algorithm>
 
 namespace npstat {
     template <typename FcnResult, typename FcnArg>
     inline long double GaussLegendreQuadrature::integrate(
         const Functor1<FcnResult,FcnArg>& function,
         const long double left, const long double right) const
     {
         if (buf_.size() != npoints_)
             buf_.resize(npoints_);
         std::pair<long double, long double>* results = &buf_[0];
         const unsigned halfpoints = npoints_/2;
         const long double midpoint = (left + right)/2.0L;
         const long double unit = (right - left)/2.0L;
         FcnArg a;
         long double v;
         for (unsigned i=0; i<halfpoints; ++i)
         {
             a = static_cast<FcnArg>(midpoint + unit*a_[i]);
             v = w_[i]*static_cast<long double>(function(a));
             results[2*i].first = std::abs(v);
             results[2*i].second = v;
             a = static_cast<FcnArg>(midpoint - unit*a_[i]);
             v = w_[i]*static_cast<long double>(function(a));
             results[2*i+1].first = std::abs(v);
             results[2*i+1].second = v;
         }
         std::sort(buf_.begin(), buf_.end());
         v = 0.0L;
         for (unsigned i=0; i<npoints_; ++i)
             v += results[i].second;
         return v*unit;
     }
 
     template <typename FcnResult, typename FcnArg>
     inline long double GaussLegendreQuadrature::integrate(
         const Functor1<FcnResult,FcnArg>& function,
         const long double left, const long double right,
         const unsigned nsplit) const
     {
         if (!nsplit) throw std::invalid_argument(
             "In npstat::GaussLegendreQuadrature::integrate: "
             "number of subintervals must be positive");
         const long double step = (right - left)/nsplit;
         long double acc = 0.0L;
         for (unsigned i=0; i<nsplit; ++i)
         {
             const long double a = left + step*i;
             const long double b = (i == nsplit - 1U ? right : a + step);
             acc += integrate(function, a, b);
         }
         return acc;
     }
 
     template <class Functor>
     inline std::vector<std::pair<double,double> >
     GaussLegendreQuadrature::weightedIntegrationPoints(
         const Functor& weight,
         const long double left, const long double right) const
     {
         typedef std::pair<double,double> DDPair;
 
         std::vector<DDPair> result;
         result.reserve(npoints_);
         const unsigned halfpoints = npoints_/2;
         const long double midpoint = (left + right)/2.0L;
         const long double unit = (right - left)/2.0L;
         for (unsigned i=0; i<halfpoints; ++i)
         {
             long double x = midpoint + unit*a_[i];
-            double w = w_[i]*weight(x);
+            double w = unit*w_[i]*weight(x);
             result.push_back(DDPair(x, w));
             x = midpoint - unit*a_[i];
-            w = w_[i]*weight(x);
+            w = unit*w_[i]*weight(x);
             result.push_back(DDPair(x, w));
         }
         return result;
     }
 }