stable quantization for Raw version

FMKGPHyperparameterOptimization.cpp

@@ -671,7 +671,8 @@ void FMKGPHyperparameterOptimization::setFastMinKernel ( FastMinKernel * _fmk )
   //
   if ( this->q != NULL )
   {  
-    this->q->computeParametersFromData ( &(this->fmk->featureMatrix()) );
+    NICE::Vector _maxValuesPerDimension = this->fmk->featureMatrix().getLargestValuePerDimension();
+    this->q->computeParametersFromData ( _maxValuesPerDimension );
   }
 }
 

GMHIKernelRaw.cpp

@@ -17,52 +17,71 @@ using namespace std;
 
 
 GMHIKernelRaw::GMHIKernelRaw( const std::vector< const NICE::SparseVector *> &_examples,
-                              const double _d_noise
-                              const NICE::Quantization * _q
+                              const double _d_noise,
+                              NICE::Quantization * _q
                             )
 {
     this->examples_raw = NULL;
     this->nnz_per_dimension = NULL;
     this->table_A = NULL;
     this->table_B = NULL;
+    this->table_T = NULL;
     this->d_noise = _d_noise;
     this->q       = _q;
 
-    initData(_examples);
+    this->initData(_examples);
 }
 
 GMHIKernelRaw::~GMHIKernelRaw()
 {
-    cleanupData();
+    this->cleanupData();
 }
 
 void GMHIKernelRaw::cleanupData()
 {
-    if ( this->examples_raw != NULL ) {
+    // data structure of examples
+    if ( this->examples_raw != NULL )
+    {
         for ( uint d = 0; d < this->num_dimension; d++ )
             if (examples_raw[d] != NULL)
                 delete [] examples_raw[d];
         delete [] this->examples_raw;
         this->examples_raw = NULL;
     }
-    if ( this->nnz_per_dimension != NULL ) {
+
+    // counter of non-zero examples in each dimension
+    if ( this->nnz_per_dimension != NULL )
+    {
         delete [] this->nnz_per_dimension;
         this->nnz_per_dimension = NULL;
     }
-    if ( this->table_A != NULL ) {
+
+    // LUT A for classification without quantization
+    if ( this->table_A != NULL )
+    {
         for ( uint d = 0; d < this->num_dimension; d++ )
             if (table_A[d] != NULL)
                 delete [] table_A[d];
         delete [] this->table_A;
         this->table_A = NULL;
     }
-    if ( this->table_B != NULL ) {
+
+    // LUT B for classification without quantization
+    if ( this->table_B != NULL )
+    {
         for ( uint d = 0; d < this->num_dimension; d++ )
             if (table_B[d] != NULL)
                 delete [] table_B[d];
         delete [] this->table_B;
         this->table_B = NULL;
     }
+
+    // LUT T for classification with quantization
+    if ( this->table_T != NULL )
+    {
+        delete [] this->table_T;
+        this->table_T = NULL;
+    }
 }
 
 void GMHIKernelRaw::initData ( const std::vector< const NICE::SparseVector *> &_examples )
@@ -72,10 +91,10 @@ void GMHIKernelRaw::initData ( const std::vector< const NICE::SparseVector *> &_
 
     cleanupData();
 
-    this->num_dimension = _examples[0]->getDim();
-    this->examples_raw = new sparseVectorElement *[num_dimension];
+    this->num_dimension     = _examples[0]->getDim();
+    this->examples_raw      = new sparseVectorElement *[num_dimension];
     this->nnz_per_dimension = new uint [num_dimension];
-    this->num_examples = _examples.size();
+    this->num_examples      = _examples.size();
 
     // waste memory and allocate a non-sparse data block
     sparseVectorElement **examples_raw_increment = new sparseVectorElement *[num_dimension];
@@ -97,10 +116,11 @@ void GMHIKernelRaw::initData ( const std::vector< const NICE::SparseVector *> &_
     NICE::Vector::iterator itDiagEl = this->diagonalElements.begin();
 
     // minor pre-allocation
-    uint index;
+    uint i_dimNonZero;
     double value;
     double l1norm;
 
+    // iterate over all provided training examples to process their data
     for ( std::vector< const NICE::SparseVector * >::const_iterator i = _examples.begin();
           i != _examples.end();
           i++, example_index++, itDiagEl++
@@ -108,15 +128,18 @@ void GMHIKernelRaw::initData ( const std::vector< const NICE::SparseVector *> &_
     {
         l1norm = 0.0;
         const NICE::SparseVector *x = *i;
+        // loop over all non-zero dimensions, copy dimension and value into our data structure, and compute the L1 norm
         for ( NICE::SparseVector::const_iterator j = x->begin(); j != x->end(); j++ )
         {
-            index = j->first;
-            value = j->second;
-            examples_raw_increment[index]->value = value;
-            examples_raw_increment[index]->example_index = example_index;
-            // move to the next element
-            examples_raw_increment[index]++;
-            this->nnz_per_dimension[index]++;
+            i_dimNonZero = j->first;
+            value        = j->second;
+
+            examples_raw_increment[i_dimNonZero]->value = value;
+            examples_raw_increment[i_dimNonZero]->example_index = example_index;
+
+            // move data pointer to the next element in the current dimension
+            examples_raw_increment[i_dimNonZero]++;
+            this->nnz_per_dimension[i_dimNonZero]++;
 
             l1norm = l1norm + value;
         }
@@ -141,7 +164,8 @@ void GMHIKernelRaw::initData ( const std::vector< const NICE::SparseVector *> &_
     if ( this->q != NULL )
     {
       // (1) if yes, setup the parameters of the quantization object
-      this->q->computeParametersFromData ( this );
+      NICE::Vector _maxValuesPerDimension = this->getLargestValuePerDimension();
+      this->q->computeParametersFromData ( _maxValuesPerDimension );
       this->table_T = this->allocateTableT();
     }
 }
@@ -185,43 +209,16 @@ void GMHIKernelRaw::copyTableAorB(double **src, double **dst) const
 
 void GMHIKernelRaw::copyTableT(double *_src, double *_dst) const
 {
-  double p_src = _src;
-  double p_dst = _dst;
-  for ( int i = 0; i < this->num_dimension * this->q->getNumberOfBins(); i++ )
+  double * p_src = _src;
+  double * p_dst = _dst;
+  for ( int i = 0; i < this->num_dimension * this->q->getNumberOfBins(); i++, p_src++, p_dst++ )
   {
     *p_dst = *p_src;
-    p_src++;
-    p_dst++;    
   }
 }
 
-void GMHIKernelRaw::updateTables ( const NICE::Vector _x ) const
+void GMHIKernelRaw::updateTablesAandB ( const NICE::Vector _x ) const
 {
-    // pre-computions if quantization is activated
-    double * prototypes;
-    double * p_prototypes;
-    uint hmax;
-
-    // store prototypes
-    if ( this->q != NULL)
-    {
-        // number of quantization bins
-         hmax = this->q->getNumberOfBins();
-
-
-        double * prototypes   = new double [ hmax * this->num_dimension ];
-        double * p_prototypes = prototypes;
-
-        for (uint dim = 0; dim < this->num_dimension; dim++)
-        {
-          for ( uint i = 0 ; i < hmax ; i++ )
-          {
-            *p_prototypes = this->q->getPrototype( i, dim );
-             p_prototypes++;
-          }
-        }
-    }
-
     // start the actual computations of A, B, and optionally T
     for (uint dim = 0; dim < this->num_dimension; dim++)
     {
@@ -245,35 +242,65 @@ void GMHIKernelRaw::updateTables ( const NICE::Vector _x ) const
 
         alpha_sum += _x[index];
         this->table_B[dim][cntNonzeroFeat] = alpha_sum;
-      }
+      }      
+    }
+}
 
-      if ( this->q != NULL)
+void GMHIKernelRaw::updateTableT ( const NICE::Vector _x ) const
+{
+    // sanity check
+    if ( this->q == NULL)
+    {
+        return;
+    }
+
+
+
+    // number of quantization bins
+    uint hmax = this->q->getNumberOfBins();
+
+
+    double * prototypes;
+    prototypes   = new double [ hmax * this->num_dimension ];
+
+    double * p_prototypes;
+    p_prototypes = prototypes;
+
+    // compute all prototypes to compare against lateron
+    for (uint dim = 0; dim < this->num_dimension; dim++)
+    {
+      for ( uint i = 0 ; i < hmax ; i++ )
       {
-        //////////
-        // variables which are needed for computing T
-        uint idxProto ( 0 ); // previously j
-        double t;
+        *p_prototypes = this->q->getPrototype( i, dim );
+         p_prototypes++;
+      }
+    }
+
+    // start the actual computation of  T
+    for (uint dim = 0; dim < this->num_dimension; dim++)
+    {
+      uint nnz = nnz_per_dimension[dim];
+
+        uint idxProtoElem; // denotes the bin number in dim i of a quantized example, previously termed qBin
 
-        uint idxProtoElem; // previously qBin
-        
         sparseVectorElement * i            = examples_raw[dim];
         sparseVectorElement * iPredecessor = examples_raw[dim];
-        
+
         // index of the element, which is always bigger than the current value fval
-        int indexElem = 0;//training_values_in_dim->example_index;
+        int indexElem = 0;
         // element of the feature
-        double elem = i->value;//training_values_in_dim->value;        
-           
+        double elem = i->value;
+
         for (uint idxProto = 0; idxProto < hmax; idxProto++) // previously j
         {
           double fvalProto = prototypes[ dim*hmax + idxProto ];
           double t;
-          
-          
-          idxProtoElem = this->q->quantize ( elem, dim ); 
-            
 
-          if (  (indexElem == 0) && (idxProto < idxProtoElem) ) 
+
+          idxProtoElem = this->q->quantize ( elem, dim );
+
+
+          if (  (indexElem == 0) && (idxProto < idxProtoElem) )
           {
             // current prototype is smaller than everything else
             // resulting value = fval * sum_l=1^n alpha_l
@@ -296,9 +323,11 @@ void GMHIKernelRaw::updateTables ( const NICE::Vector _x ) const
                }
                // compute current element in the lookup table and keep in mind that
                // indexElem is the next element and not the previous one
+
+
                if ( (idxProto >= idxProtoElem) && ( indexElem==( nnz-1 ) ) )
                {
-                 // the current prototype is equal or bigger to the largest training example in this dimension
+                 // the current prototype is equal to or larger than the largest training example in this dimension
                  // -> the term B[ dim ][ nnz-1 ] - B[ dim ][ indexElem ] is equal to zero and vanishes, which is logical, since all elements are smaller than j!
                  t = table_A[ dim ][ indexElem ];
                }
@@ -307,13 +336,12 @@ void GMHIKernelRaw::updateTables ( const NICE::Vector _x ) const
                  // standard case
                  t = table_A[ dim ][ indexElem-1 ] + fvalProto*( table_B[ dim ][ nnz-1 ] - table_B[ dim ][ indexElem-1 ] );
                }
-             }
 
-             this->table_T[ dim*hmax + idxProto ] = t;
            }
-        }
-      }
-    }
+
+           this->table_T[ dim*hmax + idxProto ] = t;
+        }//for-loop over prototypes
+    }//for-loop over dimensions
 
 
     // clean-up prototypes
@@ -327,7 +355,7 @@ void GMHIKernelRaw::updateTables ( const NICE::Vector _x ) const
 void GMHIKernelRaw::multiply (NICE::Vector & _y, const NICE::Vector & _x) const
 {
   // STEP 1: initialize tables A and B
-  updateTables(_x);
+  this->updateTablesAandB(_x);
 
   _y.resize( this->num_examples );
   _y.set(0.0);
@@ -390,11 +418,11 @@ double **GMHIKernelRaw::getTableB() const
     return t;
 }
 
-double **GMHIKernelRaw::getTableT() const
+double * GMHIKernelRaw::getTableT() const
 {
-    double **t = allocateTableT();
-    copyTableT(this->table_T, t);
-    return t;
+    double * T = this->allocateTableT();
+    copyTableT(this->table_T, T);
+    return T;
 }
 
 uint *GMHIKernelRaw::getNNZPerDimension() const
@@ -416,47 +444,20 @@ void NICE::GMHIKernelRaw::getDiagonalElements( NICE::Vector & _diagonalElements)
     _diagonalElements = this->diagonalElements;
 }
 
-double NICE::GMHIKernelRaw::getLargestValue ( ) const
-{
-  double vmax (0.0);
-  double vtmp (0.0);
-
-  uint tmpIdx ( 0 );
-  // compare largest elements of all dimensions
-  for (uint d = 0; d < this->num_dimension; d++)
-  {
-      uint nnz = this->nnz_per_dimension[d];
-
-      if ( nnz > 1 )
-      {
-          tmpIdx = tmpIdx + nnz;
-          vtmp   = this->examples_raw[tmpIdx];
-          if ( vtmp > vmax )
-          {
-            vmax = vtmp;
-          }
-      }
-  }
-
-  return vmax;
-}
-
 
 NICE::Vector NICE::GMHIKernelRaw::getLargestValuePerDimension ( ) const
 {
-  NICE::Vector vmax ( this->get_d() );
+  NICE::Vector vmax ( this->num_dimension );
 
   NICE::Vector::iterator vmaxIt = vmax.begin();
 
-  uint tmpIdx ( 0 );
   for (uint d = 0; d < this->num_dimension; d++, vmaxIt++)
   {
       uint nnz = this->nnz_per_dimension[d];
 
-      if ( nnz > 1 )
+      if ( nnz > 0 )
       {
-          tmpIdx  = tmpIdx + nnz;
-          *vmaxIt = this->examples_raw[tmpIdx];
+          *vmaxIt = this->examples_raw[ d ][ nnz-1 ].value;
       }
       else
       {

GMHIKernelRaw.h

@@ -69,6 +69,9 @@ class GMHIKernelRaw : public GenericMatrix
     void copyTableAorB(double **src, double **dst) const;
     void copyTableT(double *src, double *dst) const;
 
+    void clearTablesAandB();
+    void clearTablesT();
+
 
     double * computeTableT ( const NICE::Vector & _alpha
                            );
@@ -82,9 +85,9 @@ class GMHIKernelRaw : public GenericMatrix
   public:
 
     /** simple constructor */
-    GMHIKernelRaw( const std::vector< const NICE::SparseVector *> & examples,
-                   const double d_noise = 0.1,
-                   const NICE::Quantization * _q = NULL
+    GMHIKernelRaw( const std::vector< const NICE::SparseVector *> & _examples,
+                   const double _d_noise = 0.1,
+                   NICE::Quantization * _q = NULL
                  );
 
     /** multiply with a vector: A*x = y; this is not really const anymore!! */
@@ -109,14 +112,12 @@ class GMHIKernelRaw : public GenericMatrix
     virtual ~GMHIKernelRaw();
 
     sparseVectorElement **getDataMatrix() const { return examples_raw; };
-    void updateTables ( const NICE::Vector _x ) const;
+    void updateTablesAandB ( const NICE::Vector _x ) const;
+    void updateTableT ( const NICE::Vector _x ) const;
 
     /** get the diagonal elements of the current matrix */
     void getDiagonalElements ( NICE::Vector & _diagonalElements ) const;
 
-
-    double getLargestValue ( ) const;
-
     NICE::Vector getLargestValuePerDimension ( ) const;
 
 };

GPHIKRawClassifier.cpp

@@ -35,64 +35,53 @@ using namespace NICE;
 /////////////////////////////////////////////////////
 
 
-double * GPHIKRawClassifier::computeTableT ( const NICE::Vector & _alpha
-                                           )
+void GPHIKRawClassifier::clearSetsOfTablesAandB( )
 {
-    if (this->q == NULL )
-    {
-        return NULL;
-    }
-
-    //
-    // number of quantization bins
-    uint hmax = _q->getNumberOfBins();
 
-    // store (transformed) prototypes
-    double * prototypes   = new double [ hmax * this->num_dimension ];
-    double * p_prototypes = prototypes;
-
-    for (uint dim = 0; dim < this->num_dimension; dim++)
+    // delete all LUTs A which are needed when no quantization is activated
+    for ( std::map< uint,PrecomputedType >::iterator itA = this->precomputedA.begin();
+          itA != this->precomputedA.end();
+          itA++
+        )
     {
-      for ( uint i = 0 ; i < hmax ; i++ )
-      {
-        if ( _pf != NULL )
-        {
-          *p_prototypes = _pf->f ( dim, _q->getPrototype( i, dim ) );
-        } else
+        for ( uint idxDim = 0 ; idxDim < this->num_dimension; idxDim++ )
         {
-          *p_prototypes = _q->getPrototype( i, dim );
+            if ( (itA->second)[idxDim] != NULL )
+                delete [] (itA->second)[idxDim];
         }
-
-        p_prototypes++;
-      }
+        delete [] itA->second;
     }
+    this->precomputedA.clear();
 
-    //allocate memory for the LUT
-    double *Tlookup = new double [ hmax * this->num_dimension ];
 
-    // loop through all dimensions
-    for (uint dim = 0; dim < this->ui_d; dim++)
+    // delete all LUTs B which are needed when no quantization is activated
+    for ( std::map< uint,PrecomputedType >::iterator itB = this->precomputedB.begin();
+          itB != this->precomputedB.end();
+          itB++
+        )
     {
-        if ( nnz_per_dimension[dim] == 0 )
-            continue;
-
-        double alphaSumTotalInDim(0.0);
-        double alphaTimesXSumTotalInDim(0.0);
-
-        for ( SortedVectorSparse<double>::const_elementpointer i = nonzeroElements.begin(); i != nonzeroElements.end(); i++ )
+        for ( uint idxDim = 0 ; idxDim < this->num_dimension; idxDim++ )
         {
-          alphaSumTotalInDim += _alpha[i->second.first];
-          alphaTimesXSumTotalInDim += _alpha[i->second.first] * i->second.second;
+            if ( (itB->second)[idxDim] != NULL )
+                delete [] (itB->second)[idxDim];
         }
+        delete [] itB->second;
     }
-
-    //don't waste memory
-    delete [] prototypes;
-
-    return Tlookup;
-
+    this->precomputedB.clear();
 }
 
+void GPHIKRawClassifier::clearSetsOfTablesT( )
+{
+    // delete all LUTs used for quantization
+    for ( std::map< uint, double * >::iterator itT = this->precomputedT.begin();
+          itT != this->precomputedT.end();
+          itT++
+         )
+    {
+        delete [] itT->second;
+    }
+    this->precomputedT.clear();
+}
 
 /////////////////////////////////////////////////////
 /////////////////////////////////////////////////////
@@ -101,11 +90,16 @@ double * GPHIKRawClassifier::computeTableT ( const NICE::Vector & _alpha
 /////////////////////////////////////////////////////
 GPHIKRawClassifier::GPHIKRawClassifier( )
 {
-  this->b_isTrained = false;
-  this->confSection = "";
+  this->b_isTrained   = false;
+  this->confSection   = "";
+
   this->nnz_per_dimension = NULL;
-  this->q = NULL;
-  this->gm = NULL;
+  this->num_examples  = 0;
+  this->num_dimension = 0;
+
+  this->q             = NULL;
+  this->gm            = NULL;
+
 
   // in order to be sure about all necessary variables be setup with default values, we
   // run initFromConfig with an empty config
@@ -124,7 +118,11 @@ GPHIKRawClassifier::GPHIKRawClassifier( const Config *_conf,
 
   this->b_isTrained = false;
   this->confSection = "";
+
   this->nnz_per_dimension = NULL;
+  this->num_examples  = 0;
+  this->num_dimension = 0;
+
   this->q = NULL;
   this->gm = NULL;
 
@@ -151,11 +149,26 @@ GPHIKRawClassifier::GPHIKRawClassifier( const Config *_conf,
 
 GPHIKRawClassifier::~GPHIKRawClassifier()
 {
-  delete this->solver;
-  this->solver = NULL;
+  if ( this->solver != NULL )
+  {
+    delete this->solver;
+    this->solver = NULL;
+  }
 
-  if (gm != NULL)
-    delete gm;
+  if ( this->gm != NULL)
+  {
+    delete this->gm;
+    this->gm = NULL;
+  }
+
+  this->clearSetsOfTablesAandB();
+  this->clearSetsOfTablesT();
+
+  if ( this->q != NULL )
+  {
+      delete this->q;
+      this->q = NULL;
+  }
 }
 
 void GPHIKRawClassifier::initFromConfig(const Config *_conf,
@@ -262,81 +275,101 @@ void GPHIKRawClassifier::classify ( const NICE::SparseVector * _xstar,
      fthrow(Exception, "Classifier not trained yet -- aborting!" );
   _scores.clear();
 
-  GMHIKernelRaw::sparseVectorElement **dataMatrix = gm->getDataMatrix();
 
-  uint maxClassNo = 0;
-  for ( std::map<uint, PrecomputedType>::const_iterator i = this->precomputedA.begin() ; i != this->precomputedA.end(); i++ )
-  {
-    uint classno = i->first;
-    maxClassNo = std::max ( maxClassNo, classno );
-    double beta = 0;
-
-    if ( this->q != NULL ) {
-      std::map<uint, double *>::const_iterator j = this->precomputedT.find ( classno );
-      double *T = j->second;
-      for (SparseVector::const_iterator i = _xstar->begin(); i != _xstar->end(); i++ )
-      {
-        uint dim = i->first;
-        double v = i->second;
-        uint qBin = q->quantize( v, dim );
-
-        beta += T[dim * q->getNumberOfBins() + qBin];
-      }
-    } else {
-      const PrecomputedType & A = i->second;
-      std::map<uint, PrecomputedType>::const_iterator j = this->precomputedB.find ( classno );
-      const PrecomputedType & B = j->second;
-
-      for (SparseVector::const_iterator i = _xstar->begin(); i != _xstar->end(); i++)
-      {
-        uint dim = i->first;
-        double fval = i->second;
-
-        uint nnz = this->nnz_per_dimension[dim];
-        uint nz = this->num_examples - nnz;
-
-        if ( nnz == 0 ) continue;
-        // useful
-        //if ( fval < this->f_tolerance ) continue;
-
-        uint position = 0;
-
-        //this->X_sorted.findFirstLargerInDimension(dim, fval, position);
-        GMHIKernelRaw::sparseVectorElement fval_element;
-        fval_element.value = fval;
-
-        //std::cerr << "value to search for " << fval << endl;
-        //std::cerr << "data matrix in dimension " << dim << endl;
-        //for (int j = 0; j < nnz; j++)
-        //    std::cerr << dataMatrix[dim][j].value << std::endl;
-
-        GMHIKernelRaw::sparseVectorElement *it = upper_bound ( dataMatrix[dim], dataMatrix[dim] + nnz, fval_element );
-        position = distance ( dataMatrix[dim], it );
-        // add zero elements
-        if ( fval_element.value > 0.0 )
-            position += nz;
-
-
-        bool posIsZero ( position == 0 );
-        if ( !posIsZero )
-            position--;
-
-
-        double firstPart = 0.0;
-        if ( !posIsZero && ((position-nz) < this->num_examples) )
-          firstPart = (A[dim][position-nz]);
-
-        double secondPart( B[dim][this->num_examples-1-nz]);
-        if ( !posIsZero && (position >= nz) )
-            secondPart -= B[dim][position-nz];
-
-        // but apply using the transformed one
-        beta += firstPart + secondPart* fval;
-      }
+    // classification with quantization of test inputs
+    if ( this->q != NULL )
+    {
+        uint maxClassNo = 0;
+        for ( std::map< uint, double * >::const_iterator itT = this->precomputedT.begin() ;
+              itT != this->precomputedT.end();
+              itT++
+            )
+        {
+          uint classno = itT->first;
+          maxClassNo   = std::max ( maxClassNo, classno );
+          double beta  = 0;
+          double *T    = itT->second;
+
+          for (SparseVector::const_iterator i = _xstar->begin(); i != _xstar->end(); i++ )
+          {
+            uint dim = i->first;
+            double v = i->second;
+            uint qBin = this->q->quantize( v, dim );
+
+            beta += T[dim * this->q->getNumberOfBins() + qBin];
+          }//for-loop over dimensions of test input
+
+          _scores[ classno ] = beta;
+
+        }//for-loop over 1-vs-all models
     }
+    // classification with exact test inputs, i.e., no quantization involved
+    else
+    {
+        uint maxClassNo = 0;
+        for ( std::map<uint, PrecomputedType>::const_iterator i = this->precomputedA.begin() ; i != this->precomputedA.end(); i++ )
+        {
+          uint classno = i->first;
+          maxClassNo = std::max ( maxClassNo, classno );
+          double beta = 0;
+          GMHIKernelRaw::sparseVectorElement **dataMatrix = gm->getDataMatrix();
 
-    _scores[ classno ] = beta;
-  }
+          const PrecomputedType & A = i->second;
+          std::map<uint, PrecomputedType>::const_iterator j = this->precomputedB.find ( classno );
+          const PrecomputedType & B = j->second;
+
+          for (SparseVector::const_iterator i = _xstar->begin(); i != _xstar->end(); i++)
+          {
+            uint dim = i->first;
+            double fval = i->second;
+
+            uint nnz = this->nnz_per_dimension[dim];
+            uint nz = this->num_examples - nnz;
+
+            if ( nnz == 0 ) continue;
+            // useful
+            //if ( fval < this->f_tolerance ) continue;
+
+            uint position = 0;
+
+            //this->X_sorted.findFirstLargerInDimension(dim, fval, position);
+            GMHIKernelRaw::sparseVectorElement fval_element;
+            fval_element.value = fval;
+
+            //std::cerr << "value to search for " << fval << endl;
+            //std::cerr << "data matrix in dimension " << dim << endl;
+            //for (int j = 0; j < nnz; j++)
+            //    std::cerr << dataMatrix[dim][j].value << std::endl;
+
+            GMHIKernelRaw::sparseVectorElement *it = upper_bound ( dataMatrix[dim], dataMatrix[dim] + nnz, fval_element );
+            position = distance ( dataMatrix[dim], it );
+            // add zero elements
+            if ( fval_element.value > 0.0 )
+                position += nz;
+
+
+            bool posIsZero ( position == 0 );
+            if ( !posIsZero )
+                position--;
+
+
+            double firstPart = 0.0;
+            if ( !posIsZero && ((position-nz) < this->num_examples) )
+              firstPart = (A[dim][position-nz]);
+
+            double secondPart( B[dim][this->num_examples-1-nz]);
+            if ( !posIsZero && (position >= nz) )
+                secondPart -= B[dim][position-nz];
+
+            // but apply using the transformed one
+            beta += firstPart + secondPart* fval;
+          }//for-loop over dimensions of test input
+
+          _scores[ classno ] = beta;
+
+        }//for-loop over 1-vs-all models
+
+    } // if-condition wrt quantization
   _scores.setDim ( *this->knownClasses.rbegin() + 1 );
 
 
@@ -419,9 +452,8 @@ void GPHIKRawClassifier::train ( const std::vector< const NICE::SparseVector *>
   Timer t;
   t.start();
 
-  precomputedA.clear();
-  precomputedB.clear();
-  precomputedT.clear();
+  this->clearSetsOfTablesAandB();
+  this->clearSetsOfTablesT();
 
 
   // sort examples in each dimension and "transpose" the feature matrix
@@ -482,8 +514,8 @@ void GPHIKRawClassifier::train ( const std::vector< const NICE::SparseVector *>
 
     solver->solveLin( *gm, y, alpha );
 
-    // TODO: get lookup tables, A, B, etc. and store them
-    gm->updateTables(alpha);
+    // get lookup tables, A, B, etc. and store them
+    gm->updateTablesAandB( alpha );
     double **A = gm->getTableA();
     double **B = gm->getTableB();
 
@@ -493,14 +525,19 @@ void GPHIKRawClassifier::train ( const std::vector< const NICE::SparseVector *>
     // Quantization for classification?
     if ( this->q != NULL )
     {
-      // (2) then compute the corresponding look-up table T
-        double **B = gm->getTableT();
+      gm->updateTableT( alpha );
+      double *T = gm->getTableT ( );
+      precomputedT.insert( pair<uint, double * > ( classno, T ) );
 
-      double *T = this->computeTableT ( alpha );
-      precomputedT.insert( pair<uint, PrecomputedType> ( classno, T ) );
     }
   }
 
+  // NOTE if quantization is turned on, we do not need LUTs A and B anymore
+  if ( this->q != NULL )
+  {
+    this->clearSetsOfTablesAandB();
+  }
+
 
   t.stop();
   if ( this->b_verbose )

GPHIKRawClassifier.h

@@ -96,6 +96,9 @@ class GPHIKRawClassifier
     /////////////////////////
     /////////////////////////
 
+    void clearSetsOfTablesAandB();
+    void clearSetsOfTablesT();
+
 
     /////////////////////////
     /////////////////////////

matlab/GPHIKRawClassifierMex.cpp

@@ -56,6 +56,7 @@ NICE::Config parseParametersGPHIKRawClassifier(const mxArray *prhs[], int nrhs)
     /////////////////////////////////////////
     if( (variable == "verbose") ||
         (variable == "debug") ||
+        (variable == "use_quantization") ||
         (variable == "ils_verbose")
       )
     {
@@ -92,7 +93,9 @@ NICE::Config parseParametersGPHIKRawClassifier(const mxArray *prhs[], int nrhs)
     /////////////////////////////////////////
     // READ STRICT POSITIVE INT VARIABLES
     /////////////////////////////////////////
-    if ( variable == "ils_max_iterations" )
+    if ( (variable == "num_bins") ||
+         (variable == "ils_max_iterations" )
+       )
     {
       if ( mxIsDouble( prhs[i+1] ) )
       {
@@ -164,6 +167,14 @@ NICE::Config parseParametersGPHIKRawClassifier(const mxArray *prhs[], int nrhs)
         conf.sS("GPHIKRawClassifier", variable, value);
     }
 
+    if(variable == "s_quantType")
+    {
+      string value = MatlabConversion::convertMatlabToString( prhs[i+1] );
+      if( value != "1d-aequi-0-1" && value != "1d-aequi-0-max" && value != "nd-aequi-0-max" )
+        mexErrMsgIdAndTxt("mexnice:error","Unexpected parameter value for \'s_quantType\'. \'1d-aequi-0-1\' , \'1d-aequi-0-max\' or \'nd-aequi-0-max\' expected.");
+        conf.sS("GPHIKRawClassifier", variable, value);
+    }
+
   }
 
 

quantization/Quantization.h

@@ -13,11 +13,6 @@
 // 
 #include <core/vector/VectorT.h>
 
-// gp-hik-core includes
-#include "gp-hik-core/FeatureMatrixT.h"
-// #include "gp-hik-core/GMHIKernelRaw.h"
-//FIXME cyclic dependency!
-
 namespace NICE {
   
  /** 
@@ -99,8 +94,9 @@ class Quantization  : public NICE::Persistent
                         
                         
   //FIXME should the argument _fm be templated?
-  virtual void computeParametersFromData ( const NICE::FeatureMatrix *  _fm ) = 0;
-  virtual void computeParametersFromData ( const NICE::GMHIKernelRaw *  _gm ) = 0;
+  virtual void computeParametersFromData ( const NICE::Vector & _maxValuesPerDimension ) = 0;
+//  FeatureMatrix *  _fm
+//  virtual void computeParametersFromData ( const NICE::GMHIKernelRaw *  _gm ) = 0;
   
   ///////////////////// INTERFACE PERSISTENT /////////////////////
   // interface specific methods for store and restore

quantization/Quantization1DAequiDist0To1.cpp

@@ -52,15 +52,11 @@ uint Quantization1DAequiDist0To1::quantize ( double _value,
     return static_cast<uint> ( _value * (this->ui_numBins-1) + 0.5 );
 }
 
-void Quantization1DAequiDist0To1::computeParametersFromData ( const NICE::FeatureMatrix *  _fm )
+void Quantization1DAequiDist0To1::computeParametersFromData ( const NICE::Vector & _maxValuesPerDimension )
 {
   // nothing to do here...
 }
 
-void Quantization1DAequiDist0To1::computeParametersFromData ( const NICE::GMHIKernelRaw *  _gm )
-{
-  // nothing to do here...
-}
 // ---------------------- STORE AND RESTORE FUNCTIONS ----------------------
 
 void Quantization1DAequiDist0To1::restore ( std::istream & _is, 

quantization/Quantization1DAequiDist0To1.h

@@ -77,8 +77,7 @@ class Quantization1DAequiDist0To1  : public NICE::Quantization
                           const uint & _dim = 0
                         ) const;
                         
-  virtual void computeParametersFromData ( const NICE::FeatureMatrix *  _fm ) ;
-  virtual void computeParametersFromData ( const NICE::GMHIKernelRaw *  _gm ) ;
+  virtual void computeParametersFromData ( const NICE::Vector & _maxValuesPerDimension );
   
   ///////////////////// INTERFACE PERSISTENT /////////////////////
   // interface specific methods for store and restore

quantization/Quantization1DAequiDist0ToMax.cpp

@@ -58,19 +58,12 @@ uint Quantization1DAequiDist0ToMax::quantize ( double _value,
 
 
 
-void Quantization1DAequiDist0ToMax::computeParametersFromData ( const NICE::FeatureMatrix *  _fm )
+void Quantization1DAequiDist0ToMax::computeParametersFromData ( const NICE::Vector & _maxValuesPerDimension )
 {
-      double vmax = ( _fm->getLargestValue( ) );       
       this->v_upperBounds.resize ( 1 );
-      this->v_upperBounds ( 0 ) = vmax;
+      this->v_upperBounds ( 0 ) = _maxValuesPerDimension.Max();
 }
 
-void Quantization1DAequiDist0ToMax::computeParametersFromData ( const NICE::GMHIKernelRaw *  _gm )
-{
-      double vmax = ( _gm->getLargestValue( ) );
-      this->v_upperBounds.resize ( 1 );
-      this->v_upperBounds ( 0 ) = vmax;
-}
 // ---------------------- STORE AND RESTORE FUNCTIONS ----------------------
 
 void Quantization1DAequiDist0ToMax::restore ( std::istream & _is, 

quantization/Quantization1DAequiDist0ToMax.h

@@ -78,8 +78,7 @@ class Quantization1DAequiDist0ToMax  : public NICE::Quantization
                         ) const;
                         
                         
-  virtual void computeParametersFromData ( const NICE::FeatureMatrix *  _fm ) ;
-  virtual void computeParametersFromData ( const NICE::GMHIKernelRaw *  _gm ) ;
+  virtual void computeParametersFromData ( const NICE::Vector & _maxValuesPerDimension );
   
   ///////////////////// INTERFACE PERSISTENT /////////////////////
   // interface specific methods for store and restore

quantization/QuantizationNDAequiDist0ToMax.cpp

@@ -53,18 +53,9 @@ uint QuantizationNDAequiDist0ToMax::quantize ( double _value,
 }
 
 
-void QuantizationNDAequiDist0ToMax::computeParametersFromData ( const NICE::FeatureMatrix *  _fm )
-{     
-  // 100% quantile...
-  double d_quantile ( 1.00 ); 
-  this->v_upperBounds = _fm->getLargestValuePerDimension( d_quantile );  
-}
-
-void QuantizationNDAequiDist0ToMax::computeParametersFromData ( const NICE::GMHIKernelRaw *  _gm )
+void QuantizationNDAequiDist0ToMax::computeParametersFromData ( const NICE::Vector & _maxValuesPerDimension )
 {
-  // 100% quantile...
-  double d_quantile ( 1.00 );
-  this->v_upperBounds = _gm->getLargestValuePerDimension( d_quantile );
+  this->v_upperBounds = _maxValuesPerDimension;
 }
 
 

quantization/QuantizationNDAequiDist0ToMax.h

@@ -78,8 +78,7 @@ class QuantizationNDAequiDist0ToMax  : public NICE::Quantization
                         ) const;
                         
                         
-  virtual void computeParametersFromData ( const NICE::FeatureMatrix *  _fm ) ;
-  virtual void computeParametersFromData ( const NICE::GMHIKernelRaw *  _gm ) ;
+  virtual void computeParametersFromData ( const NICE::Vector & _maxValuesPerDimension );
                           
   ///////////////////// INTERFACE PERSISTENT /////////////////////
   // interface specific methods for store and restore