core: vertical_segment: Confidence intervals are a first-class notion.

sources/core/data_storage.cpp

@@ -21,6 +21,7 @@
 # include <endian.h>
 #endif
 
+using namespace alta;
 using namespace Eigen;
 
 // A deleter for arrays, to work around the lack of array support in C++11's
@@ -66,17 +67,10 @@ static bool within_bounds(const vecref v,
         && (v.array() > max.array()).all();
 }
 
-// The type of confidence interval.
-enum ci_kind
-{
-    NO_CONFIDENCE_INTERVAL = 0,
-    SYMMETRICAL_CONFIDENCE_INTERVAL,
-    ASYMMETRICAL_CONFIDENCE_INTERVAL
-};
-
 // Read a confidence interval on the output parameters from INPUT into V.
 static void read_confidence_interval(std::istream& input,
-                                     vecref v, ci_kind kind,
+                                     vecref v,
+                                     vertical_segment::ci_kind kind,
                                      unsigned int dimX,
                                      unsigned int dimY,
                                      const alta::arguments& args)
@@ -87,14 +81,14 @@ static void read_confidence_interval(std::istream& input,
     {
         double min_dt = 0.0, max_dt = 0.0;
 
-        if(i == 0 && kind == ASYMMETRICAL_CONFIDENCE_INTERVAL)
+        if(i == 0 && kind == vertical_segment::ASYMMETRICAL_CONFIDENCE_INTERVAL)
         {
             input >> min_dt ;
             input >> max_dt ;
             min_dt = min_dt-v(dimX + i);
             max_dt = max_dt-v(dimX + i);
         }
-        else if(i == 0 && kind == SYMMETRICAL_CONFIDENCE_INTERVAL)
+        else if(i == 0 && kind == vertical_segment::SYMMETRICAL_CONFIDENCE_INTERVAL)
         {
             double dt ;
             input >> dt ;
@@ -176,10 +170,10 @@ alta::data* alta::load_data_from_text(std::istream& input,
 #endif
 
   int vs_value = header.get_int("VS");
-  ci_kind kind =
-      vs_value == 2 ? ASYMMETRICAL_CONFIDENCE_INTERVAL
-      : (vs_value == 1 ? SYMMETRICAL_CONFIDENCE_INTERVAL
-         : NO_CONFIDENCE_INTERVAL);
+  vertical_segment::ci_kind kind =
+      vs_value == 2 ? vertical_segment::ASYMMETRICAL_CONFIDENCE_INTERVAL
+      : (vs_value == 1 ? vertical_segment::SYMMETRICAL_CONFIDENCE_INTERVAL
+         : vertical_segment::NO_CONFIDENCE_INTERVAL);
   std::vector<double> content;
 
   // Now read the body.

sources/core/vertical_segment.cpp

@@ -64,20 +64,26 @@ static vec data_max(Ref<MatrixXd> data)
 // Return a matrix view of PTR showing only the 'dimX' first rows of that
 // matrix.
 static Ref<MatrixXd>
-x_view(double *ptr, size_t cols, const parameters& params)
+x_view(double *ptr, size_t cols, const parameters& params,
+       vertical_segment::ci_kind kind)
 {
+    auto stride = params.dimX() + params.dimY()
+        + vertical_segment::confidence_interval_columns(kind, params);
+
     return Map<MatrixXd, 0, OuterStride<> >(ptr, params.dimX(), cols,
-                                            OuterStride<>(params.dimX()
-                                                          + 3 * params.dimY()));
+                                            OuterStride<>(stride));
 }
 
 vertical_segment::vertical_segment(const parameters& params,
                                    size_t size,
-                                   std::shared_ptr<double> input_data)
+                                   std::shared_ptr<double> input_data,
+                                   ci_kind kind)
     : data(params, size,
-           data_min(x_view(input_data.get(), size, params)),
-           data_max(x_view(input_data.get(), size, params))),
-      _data(input_data), _is_absolute(true), _dt(0.1)
+           data_min(x_view(input_data.get(), size, params, kind)),
+           data_max(x_view(input_data.get(), size, params, kind))),
+      _data(input_data),
+      _ci_kind(kind),
+      _is_absolute(true), _dt(0.1)
 {
 }
 
@@ -99,6 +105,9 @@ vertical_segment::vertical_segment(const parameters& params, unsigned int rows):
 
 void vertical_segment::get(int i, vec& x, vec& yl, vec& yu) const
 {
+    // Make sure we have the lower and upper bounds of Y.
+    assert(confidence_interval_kind() == ASYMMETRICAL_CONFIDENCE_INTERVAL);
+
     auto matrix = matrix_view();
 
 #ifdef DEBUG
@@ -118,6 +127,9 @@ void vertical_segment::get(int i, vec& x, vec& yl, vec& yu) const
 
 void vertical_segment::get(int i, vec& yl, vec& yu) const
 {
+    // Make sure we have the lower and upper bounds of Y.
+    assert(confidence_interval_kind() == ASYMMETRICAL_CONFIDENCE_INTERVAL);
+
     auto matrix = matrix_view();
 
     yl.resize(_parameters.dimY()) ; yu.resize(_parameters.dimY()) ;

sources/core/vertical_segment.h

@@ -78,11 +78,26 @@ namespace alta {
  */
 class vertical_segment : public data
 {
+   public: // types
+
+      // The type of confidence interval.
+      enum ci_kind
+      {
+          NO_CONFIDENCE_INTERVAL = 0,
+          SYMMETRICAL_CONFIDENCE_INTERVAL,
+          ASYMMETRICAL_CONFIDENCE_INTERVAL
+      };
+
+
    public: // methods
 
       vertical_segment(const parameters& params,
                        size_t size,
-                       std::shared_ptr<double> data);
+                       std::shared_ptr<double> data,
+
+                       // Currently we default to always storing asymmetrical
+                       // CI data for Y.
+                       ci_kind kind = ASYMMETRICAL_CONFIDENCE_INTERVAL);
 
       vertical_segment(const parameters& params, unsigned int size)
           ALTA_DEPRECATED;
@@ -105,25 +120,60 @@ class vertical_segment : public data
       //! ordinate segment.
       virtual void get(int i, vec& yl, vec& yu) const ;
 
-   private: // method
+      //! \brief Return the type of CI data provided by this object.
+      ci_kind confidence_interval_kind() const
+      {
+          return _ci_kind;
+      }
+
+      //! \brief Return the number of columns used to store confidence
+      //! interval data.
+      size_t confidence_interval_columns() const
+      {
+          return confidence_interval_columns(_ci_kind, _parameters);
+      }
+
+      //! \brief Return the number of columns used to store confidence
+      //! interval data for KIND.
+      static size_t confidence_interval_columns(ci_kind kind,
+                                                const parameters& params)
+      {
+          // Currently we're only storing CI data for Y, not for X.
+          switch (kind)
+          {
+          case NO_CONFIDENCE_INTERVAL:
+              return 0;
+          case SYMMETRICAL_CONFIDENCE_INTERVAL:
+              return params.dimY();
+          case ASYMMETRICAL_CONFIDENCE_INTERVAL:
+              return 2 * params.dimY();
+          default:
+              abort();
+          }
+      }
 
       //! \brief Return the number of columns in each row.
       size_t column_number() const
       {
-          return _parameters.dimX() + 3 * _parameters.dimY();
+          return _parameters.dimX() + _parameters.dimY()
+              + confidence_interval_columns();
       }
 
-      // Return a matrix view of the data: all the Xi and Yi followed by
-      // confidence interval data (lower and upper bound of the Yi).  Thus,
-      // it has (dimX + 3 * dimY) columns and SIZE rows.
+      //! \brief Return a matrix view of the data: all the Xi and Yi followed
+      // by confidence interval data (lower and upper bound of the Yi).
+      // Thus, it has (dimX + dimY + N * dimY) columns, where N is between 0
+      // and 2 depending on the confidence interval data available, and SIZE
+      // rows.
       Eigen::Map<Eigen::MatrixXd> matrix_view() const
       {
           return Eigen::Map<Eigen::MatrixXd>(_data.get(), size(),
                                              column_number());
       }
 
-      // Return a matrix view of DATA that excludes confidence interval data.
-      // It has (dimX + dimY) columns and SIZE rows.
+   private: // method
+
+      //! \brief Return a matrix view of DATA that excludes confidence
+      // interval data.  It has (dimX + dimY) columns and SIZE rows.
       Eigen::Map<Eigen::MatrixXd, 0, Eigen::OuterStride<> > data_view() const
       {
           return Eigen::Map<Eigen::MatrixXd, 0, Eigen::OuterStride<> >
@@ -143,6 +193,9 @@ class vertical_segment : public data
       // Store for each point of data, the upper and lower value.
       std::shared_ptr<double> _data;
 
+      // Type of confidence interval data available.
+      const ci_kind _ci_kind;
+
       // Store the different arguments for the vertical segment: is it using
       // relative or absolute intervals? What is the dt used ?
       bool   _is_absolute;

sources/plugins/rational_fitter_cgal/rational_fitter_cgal.cpp

 /* ALTA --- Analysis of Bidirectional Reflectance Distribution Functions
 
-   Copyright (C) 2013, 2014 Inria
+   Copyright (C) 2013, 2014, 2016 Inria
 
    This file is part of ALTA.
 
@@ -46,7 +46,8 @@ bool rational_fitter_cgal::fit_data(const ptr<data>& dat, ptr<function>& fit, co
 {
 	ptr<rational_function> r = dynamic_pointer_cast<rational_function>(fit) ;
 	const ptr<vertical_segment> d = dynamic_pointer_cast<vertical_segment>(dat) ;
-	if(!r || !d)
+	if(!r || !d
+     || d->confidence_interval_kind() != vertical_segment::ASYMMETRICAL_CONFIDENCE_INTERVAL)
 	{
 		std::cerr << "<<ERROR>> not passing the correct class to the fitter" << std::endl ;
 		return false ;

sources/plugins/rational_fitter_dca/rational_fitter.cpp

 /* ALTA --- Analysis of Bidirectional Reflectance Distribution Functions
 
-   Copyright (C) 2013, 2014 Inria
+   Copyright (C) 2013, 2014, 2016 Inria
 
    This file is part of ALTA.
 
@@ -42,7 +42,8 @@ bool rational_fitter_dca::fit_data(const ptr<data>& dat, ptr<function>& fit, con
 {
 	ptr<rational_function> r = dynamic_pointer_cast<rational_function>(fit) ;
 	const ptr<vertical_segment> d = dynamic_pointer_cast<vertical_segment>(dat) ;
-	if(!r || !d)
+	if(!r || !d
+     || d->confidence_interval_kind() != vertical_segment::ASYMMETRICAL_CONFIDENCE_INTERVAL)
 	{
 		std::cerr << "<<ERROR>> not passing the correct class to the fitter" << std::endl ;
 		return false ;

sources/plugins/rational_fitter_eigen/rational_fitter.cpp

@@ -46,7 +46,8 @@ bool rational_fitter_eigen::fit_data(const ptr<data>& dat, ptr<function>& fit, c
 	} 
 	
 	ptr<vertical_segment> d = dynamic_pointer_cast<vertical_segment>(dat) ;
-	if(!d)
+	if(!d
+     || d->confidence_interval_kind() != vertical_segment::ASYMMETRICAL_CONFIDENCE_INTERVAL)
 	{
 		std::cerr << "<<ERROR>> not passing the correct data object to the fitter" << std::endl ;
 		return false ;

sources/plugins/rational_fitter_leastsquare/rational_fitter.cpp

@@ -40,7 +40,8 @@ bool rational_fitter_leastsquare::fit_data(const ptr<data>& dat, ptr<function>&
 {
 	ptr<rational_function> r = dynamic_pointer_cast<rational_function>(fit) ;
 	const ptr<vertical_segment> d = dynamic_pointer_cast<vertical_segment>(dat) ;
-	if(!r || !d)
+	if(!r || !d
+     || d->confidence_interval_kind() != vertical_segment::ASYMMETRICAL_CONFIDENCE_INTERVAL)
 	{
 		std::cerr << "<<ERROR>> not passing the correct class to the fitter" << std::endl ;
 		return false ;

sources/plugins/rational_fitter_matlab/rational_fitter.cpp

 /* ALTA --- Analysis of Bidirectional Reflectance Distribution Functions
 
-   Copyright (C) 2013, 2014 Inria
+   Copyright (C) 2013, 2014, 2016 Inria
 
    This file is part of ALTA.
 
@@ -50,7 +50,8 @@ bool rational_fitter_matlab::fit_data(const ptr<data>& dat, ptr<function>& fit,
 {
 	ptr<rational_function> r = dynamic_pointer_cast<rational_function>(fit) ;
 	const ptr<vertical_segment> d = dynamic_pointer_cast<vertical_segment>(dat) ;
-	if(!r || !d || ep == NULL)
+	if(!r || !d || ep == NULL
+     || d->confidence_interval_kind() != vertical_segment::ASYMMETRICAL_CONFIDENCE_INTERVAL)
 	{
 		std::cerr << "<<ERROR>> not passing the correct class to the fitter" << std::endl ;
 		return false ;

sources/plugins/rational_fitter_parallel/rational_fitter.cpp

@@ -54,7 +54,8 @@ bool rational_fitter_parallel::fit_data(const ptr<data>& dat, ptr<function>& fit
   }
 
   ptr<vertical_segment> d = dynamic_pointer_cast<vertical_segment>(dat) ;
-  if(!d)
+  if(!d
+     || d->confidence_interval_kind() != vertical_segment::ASYMMETRICAL_CONFIDENCE_INTERVAL)
   {
     std::cerr << "<<WARNING>> automatic convertion of the data object to vertical_segment," << std::endl;
     std::cerr << "<<WARNING>> we advise you to perform convertion with a separate command." << std::endl;

sources/plugins/rational_fitter_parsec_multi/rational_fitter.cpp

@@ -72,7 +72,8 @@ bool rational_fitter_parsec_multi::fit_data(const ptr<data>& dat, ptr<function>&
 
     ptr<rational_function> r = dynamic_pointer_cast<rational_function>(fit);
     const ptr<vertical_segment>& d = dynamic_pointer_cast<vertical_segment>(dat);
-    if(!r || !d)
+    if(!r || !d
+       || d->confidence_interval_kind() != vertical_segment::ASYMMETRICAL_CONFIDENCE_INTERVAL)
 	{
 	    std::cerr << "<<ERROR>> not passing the correct class to the fitter" << std::endl;
 	    return false;

sources/plugins/rational_fitter_parsec_one/rational_fitter.cpp

@@ -70,7 +70,8 @@ bool rational_fitter_parsec_multi::fit_data(const ptr<data>& dat, ptr<function>&
 
     ptr<rational_function> r = dynamic_pointer_cast<rational_function>(fit);
     const ptr<vertical_segment>& d = dynamic_pointer_cast<vertical_segment>(dat);
-    if(!r || !d)
+    if(!r || !d
+       || d->confidence_interval_kind() != vertical_segment::ASYMMETRICAL_CONFIDENCE_INTERVAL)
 	{
 	    std::cerr << "<<ERROR>> not passing the correct class to the fitter" << std::endl;
 	    return false;

sources/plugins/rational_fitter_quadprog/rational_fitter.cpp

@@ -46,7 +46,8 @@ bool rational_fitter_quadprog::fit_data(const ptr<data>& dat, ptr<function>& fit
 {
   ptr<rational_function> r = dynamic_pointer_cast<rational_function>(fit) ;
   const ptr<vertical_segment>& d = dynamic_pointer_cast<vertical_segment>(dat) ;
-  if(!r || !d)
+  if(!r || !d
+     || d->confidence_interval_kind() != vertical_segment::ASYMMETRICAL_CONFIDENCE_INTERVAL)
   {
     std::cerr << "<<ERROR>> not passing the correct class to the fitter" << std::endl ;
     return false ;