|  | @@ -91,30 +91,25 @@ void CubicHermiteSpline(const Eigen::Matrix<double, kDataDimension, 1>& p0,
 | 
	
		
			
				|  |  |    }
 | 
	
		
			
				|  |  |  }
 | 
	
		
			
				|  |  |  
 | 
	
		
			
				|  |  | -// Given as input a one dimensional array like object, which provides
 | 
	
		
			
				|  |  | -// the following interface.
 | 
	
		
			
				|  |  | +// Given as input an infinite one dimensional grid, which provides the
 | 
	
		
			
				|  |  | +// following interface.
 | 
	
		
			
				|  |  |  //
 | 
	
		
			
				|  |  | -//   struct Array {
 | 
	
		
			
				|  |  | +//   class Grid {
 | 
	
		
			
				|  |  | +//    public:
 | 
	
		
			
				|  |  |  //     enum { DATA_DIMENSION = 2; };
 | 
	
		
			
				|  |  |  //     void GetValue(int n, double* f) const;
 | 
	
		
			
				|  |  | -//     int NumValues() const;
 | 
	
		
			
				|  |  |  //   };
 | 
	
		
			
				|  |  |  //
 | 
	
		
			
				|  |  | -// Where, GetValue gives us the value of a function f (possibly vector
 | 
	
		
			
				|  |  | -// valued) on the integers:
 | 
	
		
			
				|  |  | -//
 | 
	
		
			
				|  |  | -//   [0, ..., NumValues() - 1].
 | 
	
		
			
				|  |  | +// Here, GetValue gives the value of a function f (possibly vector
 | 
	
		
			
				|  |  | +// valued) for any integer n.
 | 
	
		
			
				|  |  |  //
 | 
	
		
			
				|  |  | -// and the enum DATA_DIMENSION indicates the dimensionality of the
 | 
	
		
			
				|  |  | -// function being interpolated. For example if you are interpolating a
 | 
	
		
			
				|  |  | -// color image with three channels (Red, Green & Blue), then
 | 
	
		
			
				|  |  | -// DATA_DIMENSION = 3.
 | 
	
		
			
				|  |  | +// The enum DATA_DIMENSION indicates the dimensionality of the
 | 
	
		
			
				|  |  | +// function being interpolated. For example if you are interpolating
 | 
	
		
			
				|  |  | +// rotations in axis-angle format over time, then DATA_DIMENSION = 3.
 | 
	
		
			
				|  |  |  //
 | 
	
		
			
				|  |  |  // CubicInterpolator uses cubic Hermite splines to produce a smooth
 | 
	
		
			
				|  |  |  // approximation to it that can be used to evaluate the f(x) and f'(x)
 | 
	
		
			
				|  |  | -// at any real valued point in the interval:
 | 
	
		
			
				|  |  | -//
 | 
	
		
			
				|  |  | -//   [0, NumValues() - 1].
 | 
	
		
			
				|  |  | +// at any point on the real number line.
 | 
	
		
			
				|  |  |  //
 | 
	
		
			
				|  |  |  // For more details on cubic interpolation see
 | 
	
		
			
				|  |  |  //
 | 
	
	
		
			
				|  | @@ -123,117 +118,122 @@ void CubicHermiteSpline(const Eigen::Matrix<double, kDataDimension, 1>& p0,
 | 
	
		
			
				|  |  |  // Example usage:
 | 
	
		
			
				|  |  |  //
 | 
	
		
			
				|  |  |  //  const double data[] = {1.0, 2.0, 5.0, 6.0};
 | 
	
		
			
				|  |  | -//  Array1D<double, 1> array(x, 4);
 | 
	
		
			
				|  |  | -//  CubicInterpolator<Array1D<double, 1> > interpolator(array);
 | 
	
		
			
				|  |  | +//  Grid1D<double, 1> grid(x, 0, 4);
 | 
	
		
			
				|  |  | +//  CubicInterpolator<Grid1D<double, 1> > interpolator(grid);
 | 
	
		
			
				|  |  |  //  double f, dfdx;
 | 
	
		
			
				|  |  | -//  CHECK(interpolator.Evaluator(1.5, &f, &dfdx));
 | 
	
		
			
				|  |  | -template<typename Array>
 | 
	
		
			
				|  |  | +//  interpolator.Evaluator(1.5, &f, &dfdx);
 | 
	
		
			
				|  |  | +template<typename Grid>
 | 
	
		
			
				|  |  |  class CERES_EXPORT CubicInterpolator {
 | 
	
		
			
				|  |  |   public:
 | 
	
		
			
				|  |  | -  explicit CubicInterpolator(const Array& array)
 | 
	
		
			
				|  |  | -      : array_(array) {
 | 
	
		
			
				|  |  | -    CHECK_GT(array.NumValues(), 1);
 | 
	
		
			
				|  |  | +  explicit CubicInterpolator(const Grid& grid)
 | 
	
		
			
				|  |  | +      : grid_(grid) {
 | 
	
		
			
				|  |  |      // The + casts the enum into an int before doing the
 | 
	
		
			
				|  |  |      // comparison. It is needed to prevent
 | 
	
		
			
				|  |  |      // "-Wunnamed-type-template-args" related errors.
 | 
	
		
			
				|  |  | -    CHECK_GE(+Array::DATA_DIMENSION, 1);
 | 
	
		
			
				|  |  | +    CHECK_GE(+Grid::DATA_DIMENSION, 1);
 | 
	
		
			
				|  |  |    }
 | 
	
		
			
				|  |  |  
 | 
	
		
			
				|  |  | -  bool Evaluate(double x, double* f, double* dfdx) const {
 | 
	
		
			
				|  |  | -    const int num_values = array_.NumValues();
 | 
	
		
			
				|  |  | -    if (x < 0 || x > num_values - 1) {
 | 
	
		
			
				|  |  | -      LOG(ERROR) << "x =  " << x
 | 
	
		
			
				|  |  | -                 << " is not in the interval [0, " << num_values - 1 << "].";
 | 
	
		
			
				|  |  | -      return false;
 | 
	
		
			
				|  |  | -    }
 | 
	
		
			
				|  |  | -
 | 
	
		
			
				|  |  | -    int n = floor(x);
 | 
	
		
			
				|  |  | -    // Deal with the case where the point sits exactly on the right
 | 
	
		
			
				|  |  | -    // boundary.
 | 
	
		
			
				|  |  | -    if (n == num_values - 1) {
 | 
	
		
			
				|  |  | -      n -= 1;
 | 
	
		
			
				|  |  | -    }
 | 
	
		
			
				|  |  | -
 | 
	
		
			
				|  |  | -    Eigen::Matrix<double, Array::DATA_DIMENSION, 1> p0, p1, p2, p3;
 | 
	
		
			
				|  |  | -
 | 
	
		
			
				|  |  | -    // The point being evaluated is now expected to lie in the
 | 
	
		
			
				|  |  | -    // internal corresponding to p1 and p2.
 | 
	
		
			
				|  |  | -    array_.GetValue(n, p1.data());
 | 
	
		
			
				|  |  | -    array_.GetValue(n + 1, p2.data());
 | 
	
		
			
				|  |  | -
 | 
	
		
			
				|  |  | -    // If we are at n >=1, the choose the element at n - 1, otherwise
 | 
	
		
			
				|  |  | -    // linearly interpolate from p1 and p2.
 | 
	
		
			
				|  |  | -    if (n > 0) {
 | 
	
		
			
				|  |  | -      array_.GetValue(n - 1, p0.data());
 | 
	
		
			
				|  |  | -    } else {
 | 
	
		
			
				|  |  | -      p0 = 2 * p1 - p2;
 | 
	
		
			
				|  |  | -    }
 | 
	
		
			
				|  |  | -
 | 
	
		
			
				|  |  | -    // If we are at n < num_values_ - 2, then choose the element n +
 | 
	
		
			
				|  |  | -    // 2, otherwise linearly interpolate from p1 and p2.
 | 
	
		
			
				|  |  | -    if (n < num_values - 2) {
 | 
	
		
			
				|  |  | -      array_.GetValue(n + 2, p3.data());
 | 
	
		
			
				|  |  | -    } else {
 | 
	
		
			
				|  |  | -      p3 = 2 * p2 - p1;
 | 
	
		
			
				|  |  | -    }
 | 
	
		
			
				|  |  | -
 | 
	
		
			
				|  |  | -    CubicHermiteSpline<Array::DATA_DIMENSION>(p0, p1, p2, p3, x - n, f, dfdx);
 | 
	
		
			
				|  |  | -
 | 
	
		
			
				|  |  | -    return true;
 | 
	
		
			
				|  |  | +  void Evaluate(double x, double* f, double* dfdx) const {
 | 
	
		
			
				|  |  | +    const int n = std::floor(x);
 | 
	
		
			
				|  |  | +    Eigen::Matrix<double, Grid::DATA_DIMENSION, 1> p0, p1, p2, p3;
 | 
	
		
			
				|  |  | +    grid_.GetValue(n - 1, p0.data());
 | 
	
		
			
				|  |  | +    grid_.GetValue(n,     p1.data());
 | 
	
		
			
				|  |  | +    grid_.GetValue(n + 1, p2.data());
 | 
	
		
			
				|  |  | +    grid_.GetValue(n + 2, p3.data());
 | 
	
		
			
				|  |  | +    CubicHermiteSpline<Grid::DATA_DIMENSION>(p0, p1, p2, p3, x - n, f, dfdx);
 | 
	
		
			
				|  |  |    }
 | 
	
		
			
				|  |  |  
 | 
	
		
			
				|  |  |    // The following two Evaluate overloads are needed for interfacing
 | 
	
		
			
				|  |  |    // with automatic differentiation. The first is for when a scalar
 | 
	
		
			
				|  |  |    // evaluation is done, and the second one is for when Jets are used.
 | 
	
		
			
				|  |  | -  bool Evaluate(const double& x, double* f) const {
 | 
	
		
			
				|  |  | -    return Evaluate(x, f, NULL);
 | 
	
		
			
				|  |  | +  void Evaluate(const double& x, double* f) const {
 | 
	
		
			
				|  |  | +    Evaluate(x, f, NULL);
 | 
	
		
			
				|  |  |    }
 | 
	
		
			
				|  |  |  
 | 
	
		
			
				|  |  | -  template<typename JetT> bool Evaluate(const JetT& x, JetT* f) const {
 | 
	
		
			
				|  |  | -    double fx[Array::DATA_DIMENSION], dfdx[Array::DATA_DIMENSION];
 | 
	
		
			
				|  |  | -    if (!Evaluate(x.a, fx, dfdx)) {
 | 
	
		
			
				|  |  | -      return false;
 | 
	
		
			
				|  |  | -    }
 | 
	
		
			
				|  |  | -
 | 
	
		
			
				|  |  | -    for (int i = 0; i < Array::DATA_DIMENSION; ++i) {
 | 
	
		
			
				|  |  | +  template<typename JetT> void Evaluate(const JetT& x, JetT* f) const {
 | 
	
		
			
				|  |  | +    double fx[Grid::DATA_DIMENSION], dfdx[Grid::DATA_DIMENSION];
 | 
	
		
			
				|  |  | +    Evaluate(x.a, fx, dfdx);
 | 
	
		
			
				|  |  | +    for (int i = 0; i < Grid::DATA_DIMENSION; ++i) {
 | 
	
		
			
				|  |  |        f[i].a = fx[i];
 | 
	
		
			
				|  |  |        f[i].v = dfdx[i] * x.v;
 | 
	
		
			
				|  |  |      }
 | 
	
		
			
				|  |  | -    return true;
 | 
	
		
			
				|  |  |    }
 | 
	
		
			
				|  |  |  
 | 
	
		
			
				|  |  | -  int NumValues() const { return array_.NumValues(); }
 | 
	
		
			
				|  |  | + private:
 | 
	
		
			
				|  |  | +  const Grid& grid_;
 | 
	
		
			
				|  |  | +};
 | 
	
		
			
				|  |  | +
 | 
	
		
			
				|  |  | +// An object that implements an infinite one dimensional grid needed
 | 
	
		
			
				|  |  | +// by the CubicInterpolator where the source of the function values is
 | 
	
		
			
				|  |  | +// an array of type T on the interval
 | 
	
		
			
				|  |  | +//
 | 
	
		
			
				|  |  | +//   [begin, ..., end - 1]
 | 
	
		
			
				|  |  | +//
 | 
	
		
			
				|  |  | +// Since the input array is finite and the grid is infinite, values
 | 
	
		
			
				|  |  | +// outside this interval needs to be computed. Grid1D uses the value
 | 
	
		
			
				|  |  | +// from the nearest edge.
 | 
	
		
			
				|  |  | +//
 | 
	
		
			
				|  |  | +// The function being provided can be vector valued, in which case
 | 
	
		
			
				|  |  | +// kDataDimension > 1. The dimensional slices of the function maybe
 | 
	
		
			
				|  |  | +// interleaved, or they maybe stacked, i.e, if the function has
 | 
	
		
			
				|  |  | +// kDataDimension = 2, if kInterleaved = true, then it is stored as
 | 
	
		
			
				|  |  | +//
 | 
	
		
			
				|  |  | +//   f01, f02, f11, f12 ....
 | 
	
		
			
				|  |  | +//
 | 
	
		
			
				|  |  | +// and if kInterleaved = false, then it is stored as
 | 
	
		
			
				|  |  | +//
 | 
	
		
			
				|  |  | +//  f01, f11, .. fn1, f02, f12, .. , fn2
 | 
	
		
			
				|  |  | +//
 | 
	
		
			
				|  |  | +template <typename T,
 | 
	
		
			
				|  |  | +          int kDataDimension = 1,
 | 
	
		
			
				|  |  | +          bool kInterleaved = true>
 | 
	
		
			
				|  |  | +struct Grid1D {
 | 
	
		
			
				|  |  | + public:
 | 
	
		
			
				|  |  | +  enum { DATA_DIMENSION = kDataDimension };
 | 
	
		
			
				|  |  | +
 | 
	
		
			
				|  |  | +  Grid1D(const T* data, const int begin, const int end)
 | 
	
		
			
				|  |  | +      : data_(data), begin_(begin), end_(end), num_values_(end - begin) {
 | 
	
		
			
				|  |  | +    CHECK_LT(begin, end);
 | 
	
		
			
				|  |  | +  }
 | 
	
		
			
				|  |  | +
 | 
	
		
			
				|  |  | +  EIGEN_STRONG_INLINE void GetValue(const int n, double* f) const {
 | 
	
		
			
				|  |  | +    const int idx = std::min(std::max(begin_, n), end_ - 1) - begin_;
 | 
	
		
			
				|  |  | +    if (kInterleaved) {
 | 
	
		
			
				|  |  | +      for (int i = 0; i < kDataDimension; ++i) {
 | 
	
		
			
				|  |  | +        f[i] = static_cast<double>(data_[kDataDimension * idx + i]);
 | 
	
		
			
				|  |  | +      }
 | 
	
		
			
				|  |  | +    } else {
 | 
	
		
			
				|  |  | +      for (int i = 0; i < kDataDimension; ++i) {
 | 
	
		
			
				|  |  | +        f[i] = static_cast<double>(data_[i * num_values_ + idx]);
 | 
	
		
			
				|  |  | +      }
 | 
	
		
			
				|  |  | +    }
 | 
	
		
			
				|  |  | +  }
 | 
	
		
			
				|  |  |  
 | 
	
		
			
				|  |  | -private:
 | 
	
		
			
				|  |  | -  const Array& array_;
 | 
	
		
			
				|  |  | + private:
 | 
	
		
			
				|  |  | +  const T* data_;
 | 
	
		
			
				|  |  | +  const int begin_;
 | 
	
		
			
				|  |  | +  const int end_;
 | 
	
		
			
				|  |  | +  const int num_values_;
 | 
	
		
			
				|  |  |  };
 | 
	
		
			
				|  |  |  
 | 
	
		
			
				|  |  | -// Given as input a two dimensional array like object, which provides
 | 
	
		
			
				|  |  | -// the following interface:
 | 
	
		
			
				|  |  | +// Given as input an infinite two dimensional grid like object, which
 | 
	
		
			
				|  |  | +// provides the following interface:
 | 
	
		
			
				|  |  |  //
 | 
	
		
			
				|  |  | -//   struct Array {
 | 
	
		
			
				|  |  | +//   struct Grid {
 | 
	
		
			
				|  |  |  //     enum { DATA_DIMENSION = 1 };
 | 
	
		
			
				|  |  |  //     void GetValue(int row, int col, double* f) const;
 | 
	
		
			
				|  |  | -//     int NumRows() const;
 | 
	
		
			
				|  |  | -//     int NumCols() const;
 | 
	
		
			
				|  |  |  //   };
 | 
	
		
			
				|  |  |  //
 | 
	
		
			
				|  |  |  // Where, GetValue gives us the value of a function f (possibly vector
 | 
	
		
			
				|  |  | -// valued) on the integer grid:
 | 
	
		
			
				|  |  | -//
 | 
	
		
			
				|  |  | -//   [0, ..., NumRows() - 1] x [0, ..., NumCols() - 1]
 | 
	
		
			
				|  |  | -//
 | 
	
		
			
				|  |  | -// and the enum DATA_DIMENSION indicates the dimensionality of the
 | 
	
		
			
				|  |  | -// function being interpolated. For example if you are interpolating a
 | 
	
		
			
				|  |  | -// color image with three channels (Red, Green & Blue), then
 | 
	
		
			
				|  |  | -// DATA_DIMENSION = 3.
 | 
	
		
			
				|  |  | +// valued) for any pairs of integers (row, col), and the enum
 | 
	
		
			
				|  |  | +// DATA_DIMENSION indicates the dimensionality of the function being
 | 
	
		
			
				|  |  | +// interpolated. For example if you are interpolating a color image
 | 
	
		
			
				|  |  | +// with three channels (Red, Green & Blue), then DATA_DIMENSION = 3.
 | 
	
		
			
				|  |  |  //
 | 
	
		
			
				|  |  |  // BiCubicInterpolator uses the cubic convolution interpolation
 | 
	
		
			
				|  |  |  // algorithm of R. Keys, to produce a smooth approximation to it that
 | 
	
		
			
				|  |  |  // can be used to evaluate the f(r,c), df(r, c)/dr and df(r,c)/dc at
 | 
	
		
			
				|  |  | -// any real valued point in the quad:
 | 
	
		
			
				|  |  | -//
 | 
	
		
			
				|  |  | -//   [0, NumRows() - 1] x [0, NumCols() - 1]
 | 
	
		
			
				|  |  | +// any point in the real plane.
 | 
	
		
			
				|  |  |  //
 | 
	
		
			
				|  |  |  // For more details on the algorithm used here see:
 | 
	
		
			
				|  |  |  //
 | 
	
	
		
			
				|  | @@ -249,55 +249,29 @@ private:
 | 
	
		
			
				|  |  |  // const double data[] = {1.0, 3.0, -1.0, 4.0,
 | 
	
		
			
				|  |  |  //                         3.6, 2.1,  4.2, 2.0,
 | 
	
		
			
				|  |  |  //                        2.0, 1.0,  3.1, 5.2};
 | 
	
		
			
				|  |  | -//  Array2D<double, 1>  array(data, 3, 4);
 | 
	
		
			
				|  |  | -//  BiCubicInterpolator<Array2D<double, 1> > interpolator(array);
 | 
	
		
			
				|  |  | +//  Grid2D<double, 1>  grid(data, 3, 4);
 | 
	
		
			
				|  |  | +//  BiCubicInterpolator<Grid2D<double, 1> > interpolator(grid);
 | 
	
		
			
				|  |  |  //  double f, dfdr, dfdc;
 | 
	
		
			
				|  |  | -//  CHECK(interpolator.Evaluate(1.2, 2.5, &f, &dfdr, &dfdc));
 | 
	
		
			
				|  |  | +//  interpolator.Evaluate(1.2, 2.5, &f, &dfdr, &dfdc);
 | 
	
		
			
				|  |  |  
 | 
	
		
			
				|  |  | -template<typename Array>
 | 
	
		
			
				|  |  | +template<typename Grid>
 | 
	
		
			
				|  |  |  class CERES_EXPORT BiCubicInterpolator {
 | 
	
		
			
				|  |  |   public:
 | 
	
		
			
				|  |  | -  explicit BiCubicInterpolator(const Array& array)
 | 
	
		
			
				|  |  | -      : array_(array) {
 | 
	
		
			
				|  |  | -    CHECK_GT(array.NumRows(), 1);
 | 
	
		
			
				|  |  | -    CHECK_GT(array.NumCols(), 1);
 | 
	
		
			
				|  |  | +  explicit BiCubicInterpolator(const Grid& grid)
 | 
	
		
			
				|  |  | +      : grid_(grid) {
 | 
	
		
			
				|  |  |      // The + casts the enum into an int before doing the
 | 
	
		
			
				|  |  |      // comparison. It is needed to prevent
 | 
	
		
			
				|  |  |      // "-Wunnamed-type-template-args" related errors.
 | 
	
		
			
				|  |  | -    CHECK_GE(+Array::DATA_DIMENSION, 1);
 | 
	
		
			
				|  |  | +    CHECK_GE(+Grid::DATA_DIMENSION, 1);
 | 
	
		
			
				|  |  |    }
 | 
	
		
			
				|  |  |  
 | 
	
		
			
				|  |  |    // Evaluate the interpolated function value and/or its
 | 
	
		
			
				|  |  |    // derivative. Returns false if r or c is out of bounds.
 | 
	
		
			
				|  |  | -  bool Evaluate(double r, double c,
 | 
	
		
			
				|  |  | +  void Evaluate(double r, double c,
 | 
	
		
			
				|  |  |                  double* f, double* dfdr, double* dfdc) const {
 | 
	
		
			
				|  |  | -    const int num_rows = array_.NumRows();
 | 
	
		
			
				|  |  | -    const int num_cols = array_.NumCols();
 | 
	
		
			
				|  |  | -
 | 
	
		
			
				|  |  | -    if (r < 0 || r > num_rows - 1 || c < 0 || c > num_cols - 1) {
 | 
	
		
			
				|  |  | -      LOG(ERROR) << "(r, c) =  (" << r << ", " << c << ")"
 | 
	
		
			
				|  |  | -                 << " is not in the square defined by [0, 0] "
 | 
	
		
			
				|  |  | -                 << " and [" << num_rows - 1 << ", " << num_cols - 1 << "]";
 | 
	
		
			
				|  |  | -      return false;
 | 
	
		
			
				|  |  | -    }
 | 
	
		
			
				|  |  | -
 | 
	
		
			
				|  |  | -    int row = floor(r);
 | 
	
		
			
				|  |  | -    // Handle the case where the point sits exactly on the bottom
 | 
	
		
			
				|  |  | -    // boundary.
 | 
	
		
			
				|  |  | -    if (row == num_rows - 1) {
 | 
	
		
			
				|  |  | -      row -= 1;
 | 
	
		
			
				|  |  | -    }
 | 
	
		
			
				|  |  | -
 | 
	
		
			
				|  |  | -    int col = floor(c);
 | 
	
		
			
				|  |  | -    // Handle the case where the point sits exactly on the right
 | 
	
		
			
				|  |  | -    // boundary.
 | 
	
		
			
				|  |  | -    if (col == num_cols - 1) {
 | 
	
		
			
				|  |  | -      col -= 1;
 | 
	
		
			
				|  |  | -    }
 | 
	
		
			
				|  |  | -
 | 
	
		
			
				|  |  |      // BiCubic interpolation requires 16 values around the point being
 | 
	
		
			
				|  |  |      // evaluated.  We will use pij, to indicate the elements of the
 | 
	
		
			
				|  |  | -    // 4x4 array of values.
 | 
	
		
			
				|  |  | +    // 4x4 grid of values.
 | 
	
		
			
				|  |  |      //
 | 
	
		
			
				|  |  |      //          col
 | 
	
		
			
				|  |  |      //      p00 p01 p02 p03
 | 
	
	
		
			
				|  | @@ -308,200 +282,89 @@ class CERES_EXPORT BiCubicInterpolator {
 | 
	
		
			
				|  |  |      // The point (r,c) being evaluated is assumed to lie in the square
 | 
	
		
			
				|  |  |      // defined by p11, p12, p22 and p21.
 | 
	
		
			
				|  |  |  
 | 
	
		
			
				|  |  | -    Eigen::Matrix<double, Array::DATA_DIMENSION, 1> p00, p01, p02, p03;
 | 
	
		
			
				|  |  | -    Eigen::Matrix<double, Array::DATA_DIMENSION, 1> p10, p11, p12, p13;
 | 
	
		
			
				|  |  | -    Eigen::Matrix<double, Array::DATA_DIMENSION, 1> p20, p21, p22, p23;
 | 
	
		
			
				|  |  | -    Eigen::Matrix<double, Array::DATA_DIMENSION, 1> p30, p31, p32, p33;
 | 
	
		
			
				|  |  | -
 | 
	
		
			
				|  |  | -    array_.GetValue(row,     col,     p11.data());
 | 
	
		
			
				|  |  | -    array_.GetValue(row,     col + 1, p12.data());
 | 
	
		
			
				|  |  | -    array_.GetValue(row + 1, col,     p21.data());
 | 
	
		
			
				|  |  | -    array_.GetValue(row + 1, col + 1, p22.data());
 | 
	
		
			
				|  |  | -
 | 
	
		
			
				|  |  | -    // If we are in rows >= 1, then choose the element from the row - 1,
 | 
	
		
			
				|  |  | -    // otherwise linearly interpolate from row and row + 1.
 | 
	
		
			
				|  |  | -    if (row > 0) {
 | 
	
		
			
				|  |  | -      array_.GetValue(row - 1, col,     p01.data());
 | 
	
		
			
				|  |  | -      array_.GetValue(row - 1, col + 1, p02.data());
 | 
	
		
			
				|  |  | -    } else {
 | 
	
		
			
				|  |  | -      p01 = 2 * p11 - p21;
 | 
	
		
			
				|  |  | -      p02 = 2 * p12 - p22;
 | 
	
		
			
				|  |  | -    }
 | 
	
		
			
				|  |  | +    const int row = std::floor(r);
 | 
	
		
			
				|  |  | +    const int col = std::floor(c);
 | 
	
		
			
				|  |  |  
 | 
	
		
			
				|  |  | -    // If we are in row < num_rows - 2, then pick the element from the
 | 
	
		
			
				|  |  | -    // row + 2, otherwise linearly interpolate from row and row + 1.
 | 
	
		
			
				|  |  | -    if (row < num_rows - 2) {
 | 
	
		
			
				|  |  | -      array_.GetValue(row + 2, col,     p31.data());
 | 
	
		
			
				|  |  | -      array_.GetValue(row + 2, col + 1, p32.data());
 | 
	
		
			
				|  |  | -    } else {
 | 
	
		
			
				|  |  | -      p31 = 2 * p21 - p11;
 | 
	
		
			
				|  |  | -      p32 = 2 * p22 - p12;
 | 
	
		
			
				|  |  | -    }
 | 
	
		
			
				|  |  | -
 | 
	
		
			
				|  |  | -    // Same logic as above, applies to the columns instead of rows.
 | 
	
		
			
				|  |  | -    if (col > 0) {
 | 
	
		
			
				|  |  | -      array_.GetValue(row,     col - 1, p10.data());
 | 
	
		
			
				|  |  | -      array_.GetValue(row + 1, col - 1, p20.data());
 | 
	
		
			
				|  |  | -    } else {
 | 
	
		
			
				|  |  | -      p10 = 2 * p11 - p12;
 | 
	
		
			
				|  |  | -      p20 = 2 * p21 - p22;
 | 
	
		
			
				|  |  | -    }
 | 
	
		
			
				|  |  | -
 | 
	
		
			
				|  |  | -    if (col < num_cols - 2) {
 | 
	
		
			
				|  |  | -      array_.GetValue(row,     col + 2, p13.data());
 | 
	
		
			
				|  |  | -      array_.GetValue(row + 1, col + 2, p23.data());
 | 
	
		
			
				|  |  | -    } else {
 | 
	
		
			
				|  |  | -      p13 = 2 * p12 - p11;
 | 
	
		
			
				|  |  | -      p23 = 2 * p22 - p21;
 | 
	
		
			
				|  |  | -    }
 | 
	
		
			
				|  |  | -
 | 
	
		
			
				|  |  | -    // The four corners of the block require a bit more care.  Let us
 | 
	
		
			
				|  |  | -    // consider the evaluation of p00, the other three corners follow
 | 
	
		
			
				|  |  | -    // in the same manner.
 | 
	
		
			
				|  |  | -    //
 | 
	
		
			
				|  |  | -    // There are four cases in which we need to evaluate p00.
 | 
	
		
			
				|  |  | -    //
 | 
	
		
			
				|  |  | -    // row > 0, col > 0 : v(row, col)
 | 
	
		
			
				|  |  | -    // row = 0, col > 0 : Interpolate p10 & p20
 | 
	
		
			
				|  |  | -    // row > 0, col = 0 : Interpolate p01 & p02
 | 
	
		
			
				|  |  | -    // row = 0, col = 0 : Interpolate p10 & p20, or p01 & p02.
 | 
	
		
			
				|  |  | -    if (row > 0) {
 | 
	
		
			
				|  |  | -      if (col > 0) {
 | 
	
		
			
				|  |  | -        array_.GetValue(row - 1, col - 1, p00.data());
 | 
	
		
			
				|  |  | -      } else {
 | 
	
		
			
				|  |  | -        p00 = 2 * p01 - p02;
 | 
	
		
			
				|  |  | -      }
 | 
	
		
			
				|  |  | -
 | 
	
		
			
				|  |  | -      if (col < num_cols - 2) {
 | 
	
		
			
				|  |  | -        array_.GetValue(row - 1, col + 2, p03.data());
 | 
	
		
			
				|  |  | -      } else {
 | 
	
		
			
				|  |  | -        p03 = 2 * p02 - p01;
 | 
	
		
			
				|  |  | -      }
 | 
	
		
			
				|  |  | -    } else {
 | 
	
		
			
				|  |  | -      p00 = 2 * p10 - p20;
 | 
	
		
			
				|  |  | -      p03 = 2 * p13 - p23;
 | 
	
		
			
				|  |  | -    }
 | 
	
		
			
				|  |  | -
 | 
	
		
			
				|  |  | -    if (row < num_rows - 2) {
 | 
	
		
			
				|  |  | -      if (col > 0) {
 | 
	
		
			
				|  |  | -        array_.GetValue(row + 2, col - 1, p30.data());
 | 
	
		
			
				|  |  | -      } else {
 | 
	
		
			
				|  |  | -        p30 = 2 * p31 - p32;
 | 
	
		
			
				|  |  | -      }
 | 
	
		
			
				|  |  | -
 | 
	
		
			
				|  |  | -      if (col < num_cols - 2) {
 | 
	
		
			
				|  |  | -        array_.GetValue(row + 2, col + 2, p33.data());
 | 
	
		
			
				|  |  | -      } else {
 | 
	
		
			
				|  |  | -        p33 = 2 * p32 - p31;
 | 
	
		
			
				|  |  | -      }
 | 
	
		
			
				|  |  | -    } else {
 | 
	
		
			
				|  |  | -      p30 = 2 * p20 - p10;
 | 
	
		
			
				|  |  | -      p33 = 2 * p23 - p13;
 | 
	
		
			
				|  |  | -    }
 | 
	
		
			
				|  |  | +    Eigen::Matrix<double, Grid::DATA_DIMENSION, 1> p0, p1, p2, p3;
 | 
	
		
			
				|  |  |  
 | 
	
		
			
				|  |  |      // Interpolate along each of the four rows, evaluating the function
 | 
	
		
			
				|  |  |      // value and the horizontal derivative in each row.
 | 
	
		
			
				|  |  | -    Eigen::Matrix<double, Array::DATA_DIMENSION, 1> f0, f1, f2, f3;
 | 
	
		
			
				|  |  | -    Eigen::Matrix<double, Array::DATA_DIMENSION, 1> df0dc, df1dc, df2dc, df3dc;
 | 
	
		
			
				|  |  | -
 | 
	
		
			
				|  |  | -    CubicHermiteSpline<Array::DATA_DIMENSION>(p00, p01, p02, p03, c - col,
 | 
	
		
			
				|  |  | -                                              f0.data(), df0dc.data());
 | 
	
		
			
				|  |  | -    CubicHermiteSpline<Array::DATA_DIMENSION>(p10, p11, p12, p13, c - col,
 | 
	
		
			
				|  |  | -                                              f1.data(), df1dc.data());
 | 
	
		
			
				|  |  | -    CubicHermiteSpline<Array::DATA_DIMENSION>(p20, p21, p22, p23, c - col,
 | 
	
		
			
				|  |  | -                                              f2.data(), df2dc.data());
 | 
	
		
			
				|  |  | -    CubicHermiteSpline<Array::DATA_DIMENSION>(p30, p31, p32, p33, c - col,
 | 
	
		
			
				|  |  | -                                              f3.data(), df3dc.data());
 | 
	
		
			
				|  |  | +    Eigen::Matrix<double, Grid::DATA_DIMENSION, 1> f0, f1, f2, f3;
 | 
	
		
			
				|  |  | +    Eigen::Matrix<double, Grid::DATA_DIMENSION, 1> df0dc, df1dc, df2dc, df3dc;
 | 
	
		
			
				|  |  | +
 | 
	
		
			
				|  |  | +    grid_.GetValue(row - 1, col - 1, p0.data());
 | 
	
		
			
				|  |  | +    grid_.GetValue(row - 1, col    , p1.data());
 | 
	
		
			
				|  |  | +    grid_.GetValue(row - 1, col + 1, p2.data());
 | 
	
		
			
				|  |  | +    grid_.GetValue(row - 1, col + 2, p3.data());
 | 
	
		
			
				|  |  | +    CubicHermiteSpline<Grid::DATA_DIMENSION>(p0, p1, p2, p3, c - col,
 | 
	
		
			
				|  |  | +                                             f0.data(), df0dc.data());
 | 
	
		
			
				|  |  | +
 | 
	
		
			
				|  |  | +    grid_.GetValue(row, col - 1, p0.data());
 | 
	
		
			
				|  |  | +    grid_.GetValue(row, col    , p1.data());
 | 
	
		
			
				|  |  | +    grid_.GetValue(row, col + 1, p2.data());
 | 
	
		
			
				|  |  | +    grid_.GetValue(row, col + 2, p3.data());
 | 
	
		
			
				|  |  | +    CubicHermiteSpline<Grid::DATA_DIMENSION>(p0, p1, p2, p3, c - col,
 | 
	
		
			
				|  |  | +                                             f1.data(), df1dc.data());
 | 
	
		
			
				|  |  | +
 | 
	
		
			
				|  |  | +    grid_.GetValue(row + 1, col - 1, p0.data());
 | 
	
		
			
				|  |  | +    grid_.GetValue(row + 1, col    , p1.data());
 | 
	
		
			
				|  |  | +    grid_.GetValue(row + 1, col + 1, p2.data());
 | 
	
		
			
				|  |  | +    grid_.GetValue(row + 1, col + 2, p3.data());
 | 
	
		
			
				|  |  | +    CubicHermiteSpline<Grid::DATA_DIMENSION>(p0, p1, p2, p3, c - col,
 | 
	
		
			
				|  |  | +                                             f2.data(), df2dc.data());
 | 
	
		
			
				|  |  | +
 | 
	
		
			
				|  |  | +    grid_.GetValue(row + 2, col - 1, p0.data());
 | 
	
		
			
				|  |  | +    grid_.GetValue(row + 2, col    , p1.data());
 | 
	
		
			
				|  |  | +    grid_.GetValue(row + 2, col + 1, p2.data());
 | 
	
		
			
				|  |  | +    grid_.GetValue(row + 2, col + 2, p3.data());
 | 
	
		
			
				|  |  | +    CubicHermiteSpline<Grid::DATA_DIMENSION>(p0, p1, p2, p3, c - col,
 | 
	
		
			
				|  |  | +                                             f3.data(), df3dc.data());
 | 
	
		
			
				|  |  |  
 | 
	
		
			
				|  |  |      // Interpolate vertically the interpolated value from each row and
 | 
	
		
			
				|  |  |      // compute the derivative along the columns.
 | 
	
		
			
				|  |  | -    CubicHermiteSpline<Array::DATA_DIMENSION>(f0, f1, f2, f3, r - row, f, dfdr);
 | 
	
		
			
				|  |  | +    CubicHermiteSpline<Grid::DATA_DIMENSION>(f0, f1, f2, f3, r - row, f, dfdr);
 | 
	
		
			
				|  |  |      if (dfdc != NULL) {
 | 
	
		
			
				|  |  |        // Interpolate vertically the derivative along the columns.
 | 
	
		
			
				|  |  | -      CubicHermiteSpline<Array::DATA_DIMENSION>(df0dc, df1dc, df2dc, df3dc,
 | 
	
		
			
				|  |  | -                                                r - row, dfdc, NULL);
 | 
	
		
			
				|  |  | +      CubicHermiteSpline<Grid::DATA_DIMENSION>(df0dc, df1dc, df2dc, df3dc,
 | 
	
		
			
				|  |  | +                                               r - row, dfdc, NULL);
 | 
	
		
			
				|  |  |      }
 | 
	
		
			
				|  |  | -
 | 
	
		
			
				|  |  | -    return true;
 | 
	
		
			
				|  |  |    }
 | 
	
		
			
				|  |  |  
 | 
	
		
			
				|  |  |    // The following two Evaluate overloads are needed for interfacing
 | 
	
		
			
				|  |  |    // with automatic differentiation. The first is for when a scalar
 | 
	
		
			
				|  |  |    // evaluation is done, and the second one is for when Jets are used.
 | 
	
		
			
				|  |  | -  bool Evaluate(const double& r, const double& c, double* f) const {
 | 
	
		
			
				|  |  | -    return Evaluate(r, c, f, NULL, NULL);
 | 
	
		
			
				|  |  | +  void Evaluate(const double& r, const double& c, double* f) const {
 | 
	
		
			
				|  |  | +    Evaluate(r, c, f, NULL, NULL);
 | 
	
		
			
				|  |  |    }
 | 
	
		
			
				|  |  |  
 | 
	
		
			
				|  |  | -  template<typename JetT> bool Evaluate(const JetT& r,
 | 
	
		
			
				|  |  | +  template<typename JetT> void Evaluate(const JetT& r,
 | 
	
		
			
				|  |  |                                          const JetT& c,
 | 
	
		
			
				|  |  |                                          JetT* f) const {
 | 
	
		
			
				|  |  | -    double frc[Array::DATA_DIMENSION];
 | 
	
		
			
				|  |  | -    double dfdr[Array::DATA_DIMENSION];
 | 
	
		
			
				|  |  | -    double dfdc[Array::DATA_DIMENSION];
 | 
	
		
			
				|  |  | -    if (!Evaluate(r.a, c.a, frc, dfdr, dfdc)) {
 | 
	
		
			
				|  |  | -      return false;
 | 
	
		
			
				|  |  | -    }
 | 
	
		
			
				|  |  | -
 | 
	
		
			
				|  |  | -    for (int i = 0; i < Array::DATA_DIMENSION; ++i) {
 | 
	
		
			
				|  |  | +    double frc[Grid::DATA_DIMENSION];
 | 
	
		
			
				|  |  | +    double dfdr[Grid::DATA_DIMENSION];
 | 
	
		
			
				|  |  | +    double dfdc[Grid::DATA_DIMENSION];
 | 
	
		
			
				|  |  | +    Evaluate(r.a, c.a, frc, dfdr, dfdc);
 | 
	
		
			
				|  |  | +    for (int i = 0; i < Grid::DATA_DIMENSION; ++i) {
 | 
	
		
			
				|  |  |        f[i].a = frc[i];
 | 
	
		
			
				|  |  |        f[i].v = dfdr[i] * r.v + dfdc[i] * c.v;
 | 
	
		
			
				|  |  |      }
 | 
	
		
			
				|  |  | -
 | 
	
		
			
				|  |  | -    return true;
 | 
	
		
			
				|  |  |    }
 | 
	
		
			
				|  |  |  
 | 
	
		
			
				|  |  | -  int NumRows() const { return array_.NumRows(); }
 | 
	
		
			
				|  |  | -  int NumCols() const { return array_.NumCols(); }
 | 
	
		
			
				|  |  | -
 | 
	
		
			
				|  |  |   private:
 | 
	
		
			
				|  |  | -  const Array& array_;
 | 
	
		
			
				|  |  | +  const Grid& grid_;
 | 
	
		
			
				|  |  |  };
 | 
	
		
			
				|  |  |  
 | 
	
		
			
				|  |  | -// An object that implements the one dimensional array like object
 | 
	
		
			
				|  |  | -// needed by the CubicInterpolator where the source of the function
 | 
	
		
			
				|  |  | -// values is an array of type T.
 | 
	
		
			
				|  |  | +// An object that implements an infinite two dimensional grid needed
 | 
	
		
			
				|  |  | +// by the BiCubicInterpolator where the source of the function values
 | 
	
		
			
				|  |  | +// is an grid of type T on the grid
 | 
	
		
			
				|  |  |  //
 | 
	
		
			
				|  |  | -// The function being provided can be vector valued, in which case
 | 
	
		
			
				|  |  | -// kDataDimension > 1. The dimensional slices of the function maybe
 | 
	
		
			
				|  |  | -// interleaved, or they maybe stacked, i.e, if the function has
 | 
	
		
			
				|  |  | -// kDataDimension = 2, if kInterleaved = true, then it is stored as
 | 
	
		
			
				|  |  | +//   [(row_start,   col_start), ..., (row_start,   col_end - 1)]
 | 
	
		
			
				|  |  | +//   [                          ...                            ]
 | 
	
		
			
				|  |  | +//   [(row_end - 1, col_start), ..., (row_end - 1, col_end - 1)]
 | 
	
		
			
				|  |  |  //
 | 
	
		
			
				|  |  | -//   f01, f02, f11, f12 ....
 | 
	
		
			
				|  |  | -//
 | 
	
		
			
				|  |  | -// and if kInterleaved = false, then it is stored as
 | 
	
		
			
				|  |  | -//
 | 
	
		
			
				|  |  | -//  f01, f11, .. fn1, f02, f12, .. , fn2
 | 
	
		
			
				|  |  | -template <typename T, int kDataDimension = 1, bool kInterleaved = true>
 | 
	
		
			
				|  |  | -struct Array1D {
 | 
	
		
			
				|  |  | -  enum { DATA_DIMENSION = kDataDimension };
 | 
	
		
			
				|  |  | -
 | 
	
		
			
				|  |  | -  Array1D(const T* data, const int num_values)
 | 
	
		
			
				|  |  | -      : data_(data), num_values_(num_values) {
 | 
	
		
			
				|  |  | -  }
 | 
	
		
			
				|  |  | -
 | 
	
		
			
				|  |  | -  void GetValue(const int n, double* f) const {
 | 
	
		
			
				|  |  | -    DCHECK_GE(n, 0);
 | 
	
		
			
				|  |  | -    DCHECK_LT(n, num_values_);
 | 
	
		
			
				|  |  | -
 | 
	
		
			
				|  |  | -    for (int i = 0; i < kDataDimension; ++i) {
 | 
	
		
			
				|  |  | -      if (kInterleaved) {
 | 
	
		
			
				|  |  | -        f[i] = static_cast<double>(data_[kDataDimension * n + i]);
 | 
	
		
			
				|  |  | -      } else {
 | 
	
		
			
				|  |  | -        f[i] = static_cast<double>(data_[i * num_values_ + n]);
 | 
	
		
			
				|  |  | -      }
 | 
	
		
			
				|  |  | -    }
 | 
	
		
			
				|  |  | -  }
 | 
	
		
			
				|  |  | -
 | 
	
		
			
				|  |  | -  int NumValues() const { return num_values_; }
 | 
	
		
			
				|  |  | -
 | 
	
		
			
				|  |  | - private:
 | 
	
		
			
				|  |  | -  const T* data_;
 | 
	
		
			
				|  |  | -  const int num_values_;
 | 
	
		
			
				|  |  | -};
 | 
	
		
			
				|  |  | -
 | 
	
		
			
				|  |  | -// An object that implements the two dimensional array like object
 | 
	
		
			
				|  |  | -// needed by the BiCubicInterpolator where the source of the function
 | 
	
		
			
				|  |  | -// values is an array of type T.
 | 
	
		
			
				|  |  | +// Since the input grid is finite and the grid is infinite, values
 | 
	
		
			
				|  |  | +// outside this interval needs to be computed. Grid2D uses the value
 | 
	
		
			
				|  |  | +// from the nearest edge.
 | 
	
		
			
				|  |  |  //
 | 
	
		
			
				|  |  |  // The function being provided can be vector valued, in which case
 | 
	
		
			
				|  |  |  // kDataDimension > 1. The data maybe stored in row or column major
 | 
	
	
		
			
				|  | @@ -522,37 +385,55 @@ template <typename T,
 | 
	
		
			
				|  |  |            int kDataDimension = 1,
 | 
	
		
			
				|  |  |            bool kRowMajor = true,
 | 
	
		
			
				|  |  |            bool kInterleaved = true>
 | 
	
		
			
				|  |  | -struct Array2D {
 | 
	
		
			
				|  |  | +struct Grid2D {
 | 
	
		
			
				|  |  | + public:
 | 
	
		
			
				|  |  |    enum { DATA_DIMENSION = kDataDimension };
 | 
	
		
			
				|  |  |  
 | 
	
		
			
				|  |  | -  Array2D(const T* data, const int num_rows, const int num_cols)
 | 
	
		
			
				|  |  | -      : data_(data), num_rows_(num_rows), num_cols_(num_cols) {
 | 
	
		
			
				|  |  | +  Grid2D(const T* data,
 | 
	
		
			
				|  |  | +         const int row_begin, const int row_end,
 | 
	
		
			
				|  |  | +         const int col_begin, const int col_end)
 | 
	
		
			
				|  |  | +      : data_(data),
 | 
	
		
			
				|  |  | +        row_begin_(row_begin), row_end_(row_end),
 | 
	
		
			
				|  |  | +        col_begin_(col_begin), col_end_(col_end),
 | 
	
		
			
				|  |  | +        num_rows_(row_end - row_begin), num_cols_(col_end - col_begin),
 | 
	
		
			
				|  |  | +        num_values_(num_rows_ * num_cols_) {
 | 
	
		
			
				|  |  |      CHECK_GE(kDataDimension, 1);
 | 
	
		
			
				|  |  | +    CHECK_LT(row_begin, row_end);
 | 
	
		
			
				|  |  | +    CHECK_LT(col_begin, col_end);
 | 
	
		
			
				|  |  |    }
 | 
	
		
			
				|  |  |  
 | 
	
		
			
				|  |  | -  void GetValue(const int r, const int c, double* f) const {
 | 
	
		
			
				|  |  | -    DCHECK_GE(r, 0);
 | 
	
		
			
				|  |  | -    DCHECK_LT(r, num_rows_);
 | 
	
		
			
				|  |  | -    DCHECK_GE(c, 0);
 | 
	
		
			
				|  |  | -    DCHECK_LT(c, num_cols_);
 | 
	
		
			
				|  |  | +  EIGEN_STRONG_INLINE void GetValue(const int r, const int c, double* f) const {
 | 
	
		
			
				|  |  | +    const int row_idx =
 | 
	
		
			
				|  |  | +        std::min(std::max(row_begin_, r), row_end_ - 1) - row_begin_;
 | 
	
		
			
				|  |  | +    const int col_idx =
 | 
	
		
			
				|  |  | +        std::min(std::max(col_begin_, c), col_end_ - 1) - col_begin_;
 | 
	
		
			
				|  |  |  
 | 
	
		
			
				|  |  | -    const int n = (kRowMajor) ? num_cols_ * r + c : num_rows_ * c + r;
 | 
	
		
			
				|  |  | -    for (int i = 0; i < kDataDimension; ++i) {
 | 
	
		
			
				|  |  | -      if (kInterleaved) {
 | 
	
		
			
				|  |  | +    const int n =
 | 
	
		
			
				|  |  | +        (kRowMajor)
 | 
	
		
			
				|  |  | +        ? num_cols_ * row_idx + col_idx
 | 
	
		
			
				|  |  | +        : num_rows_ * col_idx + row_idx;
 | 
	
		
			
				|  |  | +
 | 
	
		
			
				|  |  | +
 | 
	
		
			
				|  |  | +    if (kInterleaved) {
 | 
	
		
			
				|  |  | +      for (int i = 0; i < kDataDimension; ++i) {
 | 
	
		
			
				|  |  |          f[i] = static_cast<double>(data_[kDataDimension * n + i]);
 | 
	
		
			
				|  |  | -      } else {
 | 
	
		
			
				|  |  | -        f[i] = static_cast<double>(data_[i * (num_rows_ * num_cols_) + n]);
 | 
	
		
			
				|  |  | +      }
 | 
	
		
			
				|  |  | +    } else {
 | 
	
		
			
				|  |  | +      for (int i = 0; i < kDataDimension; ++i) {
 | 
	
		
			
				|  |  | +        f[i] = static_cast<double>(data_[i * num_values_ + n]);
 | 
	
		
			
				|  |  |        }
 | 
	
		
			
				|  |  |      }
 | 
	
		
			
				|  |  |    }
 | 
	
		
			
				|  |  |  
 | 
	
		
			
				|  |  | -  int NumRows() const { return num_rows_; }
 | 
	
		
			
				|  |  | -  int NumCols() const { return num_cols_; }
 | 
	
		
			
				|  |  | -
 | 
	
		
			
				|  |  |   private:
 | 
	
		
			
				|  |  |    const T* data_;
 | 
	
		
			
				|  |  | +  const int row_begin_;
 | 
	
		
			
				|  |  | +  const int row_end_;
 | 
	
		
			
				|  |  | +  const int col_begin_;
 | 
	
		
			
				|  |  | +  const int col_end_;
 | 
	
		
			
				|  |  |    const int num_rows_;
 | 
	
		
			
				|  |  |    const int num_cols_;
 | 
	
		
			
				|  |  | +  const int num_values_;
 | 
	
		
			
				|  |  |  };
 | 
	
		
			
				|  |  |  
 | 
	
		
			
				|  |  |  }  // namespace ceres
 |