diff --git a/examples/exFit.cpp b/examples/exFit.cpp
index b1fff6a..4a3376f 100644
--- a/examples/exFit.cpp
+++ b/examples/exFit.cpp
@@ -8,10 +8,27 @@ int main(void)
 {
     FitInterface f;
     
-    f.addXDim("a", 3);
-    f.addXDim("b", 4);
-    f.addXDim("c", 3);
-    f.addYDim("y");
+    f.addYDim("q1");
+    f.addYDim("q2");
+    f.addXDim("x1", 6);
+    f.addXDim("x2", 5);
+    f.addXDim("x3", 5);
+    f.registerDataPoint(f.dataIndex(0,0,0), 0);
+    f.registerDataPoint(f.dataIndex(1,1,1), 0);
+    f.registerDataPoint(f.dataIndex(2,2,2), 0);
+    f.registerDataPoint(f.dataIndex(2,3,3), 0);
+    f.registerDataPoint(f.dataIndex(0,0,0), 1);
+    f.registerDataPoint(f.dataIndex(1,1,1), 1);
+    f.registerDataPoint(f.dataIndex(2,2,3), 1);
+    f.fitPoint(false, f.dataIndex(1,1,1), 1);
+    cout << f << endl;
+    DEBUG_VAR(f.getYFitSize());
+    DEBUG_VAR(f.getYFitSize(0));
+    DEBUG_VAR(f.getYFitSize(1));
+    DEBUG_VAR(f.getXFitSize());
+    DEBUG_VAR(f.getXFitSize(0));
+    DEBUG_VAR(f.getXFitSize(1));
+    DEBUG_VAR(f.getXFitSize(2));
     
     return EXIT_SUCCESS;
 }
diff --git a/lib/FitInterface.cpp b/lib/FitInterface.cpp
index 39c545e..069e8d7 100644
--- a/lib/FitInterface.cpp
+++ b/lib/FitInterface.cpp
@@ -23,6 +23,50 @@
 using namespace std;
 using namespace Latan;
 
+// error checks ////////////////////////////////////////////////////////////////
+#define checkXDim(i)\
+if ((i) >= getNXDim())\
+{\
+    LATAN_ERROR(Range, "X dimension " + strFrom(i) + " out of range");\
+}
+
+#define checkXIndex(vi, i)\
+if ((vi) >= getXSize(i))\
+{\
+    LATAN_ERROR(Range, "index " + strFrom(vi) + " in X dimension "\
+                + strFrom(i) + " out of range");\
+}
+
+#define checkYDim(j)\
+if ((j) >= getNYDim())\
+{\
+    LATAN_ERROR(Range, "Y dimension " + strFrom(j) + " out of range");\
+}
+
+#define checkDataIndex(k)\
+if ((k) >= maxDataIndex_)\
+{\
+    LATAN_ERROR(Range, "data point index " + strFrom(k) + " invalid");\
+}
+
+#define checkDataCoord(v)\
+if (static_cast<Index>((v).size()) != getNXDim())\
+{\
+    LATAN_ERROR(Size, "number of coordinates and number of X dimensions "\
+                "mismatch");\
+}\
+for (unsigned int i_ = 0; i_ < (v).size(); ++i_)\
+{\
+    checkXIndex((v)[i_], i_);\
+}
+
+#define checkPoint(k, j)\
+if (!isXUsed(k, j))\
+{\
+    LATAN_ERROR(Range, "no data point in Y dimension " + strFrom(j)\
+                + " with index " + strFrom(k));\
+}
+
 /******************************************************************************
  *                     FitInterface implementation                            *
  ******************************************************************************/
@@ -33,42 +77,241 @@ FitInterface::FitInterface(void)
 // add dimensions //////////////////////////////////////////////////////////////
 void FitInterface::addXDim(const string name, const Index nData)
 {
-    xDimName_.push_back(name);
-    xSize_.push_back(nData);
-    xIndex_[name] = xDimName_.size();
-    updateYSize();
+    if (getYSize() != 0)
+    {
+        LATAN_ERROR(Logic, "cannot add an X dimension if fit data is "
+                    "not empty");
+    }
+    else
+    {
+        xDimName_.push_back(name);
+        xSize_.push_back(nData);
+        xDimIndex_[name]  = xDimName_.size();
+        maxDataIndex_    *= nData;
+        updateDataSize();
+    }
 }
 
 void FitInterface::addYDim(const string name)
 {
     yDimName_.push_back(name);
-    yIndex_[name] = yDimName_.size();
+    yDataIndex_.push_back(map<Index, bool>());
+    yDimIndex_[name] = yDimName_.size();
 }
 
-unsigned int FitInterface::getNXDim(void)
+// size access /////////////////////////////////////////////////////////////////
+Index FitInterface::getNXDim(void) const
 {
     return xDimName_.size();
 }
 
-unsigned int FitInterface::getNYDim(void)
+Index FitInterface::getNYDim(void) const
 {
     return yDimName_.size();
 }
 
-void FitInterface::updateYSize(void)
+Index FitInterface::getXSize(void) const
 {
-    ySize_ = 1;
-    for (Index size: xSize_)
+    Index size = 0;
+    
+    for (Index i = 0; i < getNXDim(); ++i)
     {
-        ySize_ *= size;
+        size += getXSize(i);
     }
-    isFitPoint_.resize(ySize_);
-    isFitPoint_.fill(1);
+    
+    return size;
+}
+
+Index FitInterface::getXSize(const Index i) const
+{
+    checkXDim(i);
+    
+    return xSize_[i];
+}
+
+Index FitInterface::getYSize(void) const
+{
+    Index size = 0;
+    
+    for (Index j = 0; j < getNYDim(); ++j)
+    {
+        size += getYSize(j);
+    }
+    
+    return size;
+}
+
+Index FitInterface::getYSize(const Index j) const
+{
+    checkYDim(j);
+    
+    return static_cast<Index>(yDataIndex_[j].size());
+}
+
+Index FitInterface::getXFitSize(void) const
+{
+    set<Index> fitIndex;
+    
+    for (Index j = 0; j < getNYDim(); ++j)
+    {
+        for (auto &p: yDataIndex_[j])
+        {
+            if (p.second)
+            {
+                fitIndex.insert(p.first);
+            }
+        }
+    }
+    
+    return static_cast<Index>(fitIndex.size());
+}
+
+Index FitInterface::getXFitSize(const Index i) const
+{
+    set<Index>    fitCoord;
+    vector<Index> v;
+    
+    checkXDim(i);
+    for (Index j = 0; j < getNYDim(); ++j)
+    {
+        for (auto &p: yDataIndex_[j])
+        {
+            if (p.second)
+            {
+                v = dataCoord(p.first);
+                fitCoord.insert(v[i]);
+            }
+        }
+    }
+    
+    return fitCoord.size();
+}
+
+Index FitInterface::getYFitSize(void) const
+{
+    Index size = 0;
+    
+    for (Index j = 0; j < getNYDim(); ++j)
+    {
+        size += getYFitSize(j);
+    }
+    
+    return size;
+}
+
+Index FitInterface::getYFitSize(const Index j) const
+{
+    Index size;
+    auto  pred = [](const pair<Index, bool> &p)
+    {
+        return p.second;
+    };
+
+    checkYDim(j);
+    size = count_if(yDataIndex_[j].begin(), yDataIndex_[j].end(), pred);
+    
+    return size;
+}
+
+// Y dimension index helper ////////////////////////////////////////////////////
+Index FitInterface::dataIndex(const vector<Index> &v) const
+{
+    Index k, n = static_cast<Index>(v.size());
+    
+    checkDataCoord(v);
+    k = xSize_[1]*v[0];
+    for (unsigned int d = 1; d < n-1; ++d)
+    {
+        k = xSize_[d+1]*(v[d] + k);
+    }
+    k += v[n-1];
+    
+    return k;
+}
+
+vector<Index> FitInterface::dataCoord(const Index k) const
+{
+    vector<Index> v(getNXDim());
+    Index         buf, dimProd;
+    
+    checkDataIndex(k);
+    buf     = k;
+    dimProd = 1;
+    for (Index d = getNXDim() - 1; d >= 0; --d)
+    {
+        v[d]     = (buf/dimProd)%xSize_[d];
+        buf     -= dimProd*v[d];
+        dimProd *= xSize_[d];
+    }
+    
+    return v;
 }
 
 // enable fit points ///////////////////////////////////////////////////////////
-void FitInterface::fitPoint(const bool isFitPoint, const Index j)
+void FitInterface::fitPoint(const bool isFitPoint, const Index k, const Index j)
 {
-    DEBUG_VAR(j);
-    isFitPoint_(j) = isFitPoint ? 1 : 0;
+    checkPoint(k, j);
+    yDataIndex_[j][k] = isFitPoint;
+}
+
+// tests ///////////////////////////////////////////////////////////////////////
+bool FitInterface::isXUsed(const Index k) const
+{
+    bool isUsed = false;
+    
+    for (Index j = 0; j < getNYDim(); ++j)
+    {
+        isUsed = isUsed or isXUsed(k, j);
+    }
+    
+    return isUsed;
+}
+
+bool FitInterface::isXUsed(const Index k, const Index j) const
+{
+    checkYDim(j);
+    
+    return !(yDataIndex_[j].find(k) == yDataIndex_[j].end());
+}
+
+bool FitInterface::isFitPoint(const Index k, const Index j) const
+{
+    checkPoint(k, j);
+    
+    return yDataIndex_[j].at(k);
+}
+
+// register a data point ///////////////////////////////////////////////////////
+void FitInterface::registerDataPoint(const Index k, const Index j)
+{
+    checkYDim(j);
+    yDataIndex_[j][k] = true;
+}
+
+// IO //////////////////////////////////////////////////////////////////////////
+ostream & Latan::operator<<(ostream &out, FitInterface &f)
+{
+    out << "X dimensions: " << f.getNXDim() << endl;
+    for (Index i = 0; i < f.getNXDim(); ++i)
+    {
+        out << "  * " << i << " \"" << f.xDimName_[i] << "\": ";
+        out << f.getXSize(i) << " value(s)" << endl;
+    }
+    out << "Y dimensions: " << f.getNYDim() << endl;
+    for (Index j = 0; j < f.getNYDim(); ++j)
+    {
+        out << "  * " << j << " \"" << f.yDimName_[j] << "\": ";
+        out << f.getYSize(j) << " value(s)" << endl;
+        for (auto &p: f.yDataIndex_[j])
+        {
+            out << "    " << setw(3) << p.first << " (";
+            for (auto vi: f.dataCoord(p.first))
+            {
+                out << vi << ",";
+            }
+            out << "\b) fit: " << (p.second ? "true" : "false") << endl;
+        }
+    }
+    
+    return out;
 }
diff --git a/lib/FitInterface.hpp b/lib/FitInterface.hpp
index d9334e5..2e48b34 100644
--- a/lib/FitInterface.hpp
+++ b/lib/FitInterface.hpp
@@ -27,7 +27,6 @@ BEGIN_LATAN_NAMESPACE
 /******************************************************************************
  *                             FitInterface                                   *
  ******************************************************************************/
-
 class FitInterface
 {
 public:
@@ -36,63 +35,61 @@ public:
     // destructor
     virtual ~FitInterface(void) = default;
     // add dimensions
-    void         addXDim(const std::string name, const Index nData);
-    void         addYDim(const std::string name);
-    unsigned int getNXDim(void);
-    unsigned int getNYDim(void);
+    void  addXDim(const std::string name, const Index nData);
+    void  addYDim(const std::string name);
+    // size access
+    Index getNXDim(void) const;
+    Index getNYDim(void) const;
+    Index getXSize(void) const;
+    Index getXSize(const Index i) const;
+    Index getYSize(void) const;
+    Index getYSize(const Index j) const;
+    Index getXFitSize(void) const;
+    Index getXFitSize(const Index i) const;
+    Index getYFitSize(void) const;
+    Index getYFitSize(const Index j) const;
     // Y dimension index helper
     template <typename... Ts>
-    Index yIndex(const Ts... is);
+    Index              dataIndex(const Ts... is) const;
+    Index              dataIndex(const std::vector<Index> &v) const;
+    std::vector<Index> dataCoord(const Index k) const;
     // enable fit points
-    void fitPoint(const bool isFitPoint, const Index j);
-private:
-    void updateYSize(void);
-private:
-    std::vector<std::string>     xDimName_, yDimName_;
-    std::map<std::string, Index> xIndex_, yIndex_;
-    std::vector<Index>           xSize_;
-    Index                        ySize_;
+    void fitPoint(const bool isFitPoint, const Index k, const Index j = 0);
+    // tests
+    bool isXUsed(const Index k) const;
+    bool isXUsed(const Index k, const Index j) const;
+    bool isFitPoint(const Index k, const Index j) const;
+    // IO
+    friend std::ostream & operator<<(std::ostream &out, FitInterface &f);
+protected:
 public:
-    IVec                         isFitPoint_;
+    // register a data point
+    void registerDataPoint(const Index k, const Index j = 0);
+    // abstract method to update data container size
+    virtual void updateDataSize(void) {};
+private:
+    std::vector<std::string>           xDimName_, yDimName_;
+    std::map<std::string, Index>       xDimIndex_, yDimIndex_;
+    std::vector<Index>                 xSize_;
+    std::vector<std::map<Index, bool>> yDataIndex_;
+    Index                              maxDataIndex_{1};
 };
 
+std::ostream & operator<<(std::ostream &out, FitInterface &f);
+
 /******************************************************************************
  *                  FitInterface template implementation                      *
  ******************************************************************************/
 // Y dimension index helper ////////////////////////////////////////////////////
 template <typename... Ts>
-Index FitInterface::yIndex(const Ts... is)
+Index FitInterface::dataIndex(const Ts... coords) const
 {
     static_assert(static_or<std::is_convertible<Index, Ts>::value...>::value,
                   "fitPoint arguments are not compatible with Index");
     
-    constexpr size_t n = sizeof...(is);
-    
-    if (n != getNXDim())
-    {
-        LATAN_ERROR(Size, "number of arguments and number of X dimensions "
-                          "mismatch");
-    }
-    
-    constexpr std::array<Index, sizeof...(is)> i = {is...};
-    Index                                      j;
-    
-    for (unsigned int d = 0; d < n; ++d)
-    {
-        if (i[d] >= xSize_[d])
-        {
-            LATAN_ERROR(Range, "index in X dimension " + strFrom(d)
-                               + " out of range");
-        }
-    }
-    j = xSize_[1]*i[0];
-    for (unsigned int d = 1; d < n-1; ++d)
-    {
-        j = xSize_[d+1]*(i[d] + j);
-    }
-    j += i[n-1];
+    const std::vector<Index> coord = {coords...};
 
-    return j;
+    return dataIndex(coord);
 }
 
 END_LATAN_NAMESPACE