idx.h

Gaussian Mixture Algorithm

/***************************************************************************
 *   Copyright (C) 2008 by Yann LeCun   *
 *   [e-mail]   *
 *   All rights reserved.
 *
 * Redistribution and use in source and binary forms, with or without
 * modification, are permitted provided that the following conditions are met:
 *     * Redistributions of source code must retain the above copyright
 *       notice, this list of conditions and the following disclaimer.
 *     * Redistributions in binary form must reproduce the above copyright
 *       notice, this list of conditions and the following disclaimer in the
 *       documentation and/or other materials provided with the distribution.
 *     * Redistribution under a license not approved by the Open Source
 *       Initiative (http://www.opensource.org) must display the
 *       following acknowledgement in all advertising material:
 *        This product includes software developed at the Courant
 *        Institute of Mathematical Sciences (http://cims.nyu.edu).
 *     * The names of the authors may not be used to endorse or promote products
 *       derived from this software without specific prior written permission.
 *
 * THIS SOFTWARE IS PROVIDED ``AS IS'' AND ANY EXPRESS OR IMPLIED
 * WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
 * WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
 * DISCLAIMED. IN NO EVENT SHALL ThE AUTHORS BE LIABLE FOR ANY
 * DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
 * (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
 * LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND
 * ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
 * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
 * SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
 ***************************************************************************/

#ifndef Idx_H
#define Idx_H

#define USING_STL_ITERS 1

#include <stdio.h>
#include <iostream>
#include <sstream>
#include <vector>
#include <algorithm>
#include <numeric>
#include "LibidxDefines.h"
#include "Srg.h"

namespace ebl {

const int MAXDIMS=8;

// Forward declarations of iterators
template<typename T> class ScalarIter;
template<typename T> class ReverseScalarIter;
template<typename T> class DimIter;
template<typename T> class ReverseDimIter;

////////////////////////////////////////////////////////////////

class dummyt {  bool someunk; };

//! IdxSpec contains all the characteristics of an Idx,
//! except the storage. It includes the order (number of dimensions)
//! the offset, and the dim and mod arrays.
//! having an IdxSpec class separate from Idx allows us to write
//! generic manipulation functions that do not depend on the
//! type of the storage.
class IdxSpec {

 private:

  //! private method to set the order (number of dims) to n.
  //! Calling this with n = 0 deallocates the dim and mod arrays.
  //! The dim and mod arrays are reallocated as necessary when
  //! the order is increased.
  int setndim(int n);

  //! private method to set the order to n, and the dim and
  //! mod arrays to pre-allocated arrays ldim and lmod.
  int setndim(int n, intg *ldim, intg *lmod);

  //! number of dimensions
  int ndim;
  //! offset in the storage
  intg offset;
  //! array of sizes in each dimension
  intg *dim;
  //! array of strides in each dimension
  intg *mod;

  //! resize the spec and return the new footprint.
  //! this is private because only Idx can call this
  //! so that the storage gets properly resized.
  //! We do not allow the order to be changed with this,
  //! only the size of each dimension can be modified.
  //! The resized spec will be contiguous but will have
  //! the same offset as the original.
  intg resize(intg s0=-1, intg s1=-1, intg s2=-1, intg s3=-1,
	      intg s4=-1, intg s5=-1, intg s6=-1, intg s7=-1);

  template<typename SizeIter> intg resize( SizeIter& dimsBegin, SizeIter& dimsEnd ){

  	const int nArgDims = std::distance(dimsBegin, dimsEnd);

  	// Error-check the supplied number of dims.
  	if( ndim == 0 ){
  		ylerror("Cannot call resize on a 0-dimensional IdxSpec.");
  	}
  	else if( ndim != nArgDims ){
  		std::ostringstream oss;
  		oss<<"Number of supplied dimension sizes ("<<nArgDims<<") doesn't match IdxSpec's number of dims ("<<ndim<<")";
  		ylerror(oss.str().c_str());
  	}

  	// copy dimensions to dim
  	std::copy(dimsBegin, dimsEnd, dim);

  	// set mod to be the partial sum of the dim sequence, in reverse order.
  	typedef std::reverse_iterator<SizeIter> RIter;
  	typedef std::reverse_iterator<intg*> RintgIter;
  	std::partial_sum( RIter(dimsEnd-1),
  					  RIter(dimsBegin-1),
  					  RintgIter(mod+(nArgDims-1)),
  					  std::multiplies<intg>() );

  	// return the memory footprint
  	return mod[0] * dim[0] + offset;
  }

  //! set the offset and return the new value
  intg setoffset(intg o) { return offset = o; }

 public:

  //! the destructor of IdxSpec deallocates the dim and mod arrays.
  ~IdxSpec();

  //! assignment operator overloading.
  const IdxSpec& operator=( const IdxSpec& src);

  //! copy IdxSpec src into current IdxSpec
  void copy( const IdxSpec& src);

  //! copy constructor from IdxSpec src
  IdxSpec( const IdxSpec& src);

  //! This creates an IdxSpec0 with offset 0.
  //! This can be used to build an empty/blank IdxSpec.
  IdxSpec();

  //! Creates an IdxSpec0 with offset o.
  IdxSpec(intg o);

  //! Creates an IdxSpec1 with offset o.
  IdxSpec(intg o, intg size0);

  //! Creates an IdxSpec2 with offset o.
  IdxSpec(intg o, intg size0, intg size1);

  //! Creates an IdxSpec3 with offset o.
  IdxSpec(intg o, intg size0, intg size1, intg size2);

  //! Generic constructor with offset.
  IdxSpec(intg o, intg s0, intg s1, intg s2, intg s3,
	  intg s4=-1, intg s5=-1, intg s6=-1, intg s7=-1);

  //! Creates an IdxSpec of order n from arrays of dims and mods.
  //! The arrays are copied.
  IdxSpec(intg o, int n, intg *ldim, intg *lmod);

  //! return the offset of IdxSpec
  intg getoffset() { return offset; }

  //! return the order (number of dimensions).
  int getndim() { return ndim; }

  //! return the memory footprint, including the offset.
  //! The storage of an Idx containing this IdxSpec must
  //! be have at least this size.
  intg footprint()  {
    intg r = offset + 1;
    for(int i=0; i<ndim; i++){ r += mod[i]*(dim[i]-1); }
    return r;
  }

  //! total number of elements accessed by IdxSpec
  intg nelements() {
    intg r = 1;
    for(int i=0; i<ndim; i++){ r *= dim[i]; }
    return r;
  }

  //! returns true if the IdxSpec elements are
  //! in a continuous chunk of memory. This is useful
  //! to optimize iterators over the data.
  bool contiguousp() {
    intg size = 1; bool r = true;
    for(int i=ndim-1; i>=0; i--){
      if (size != mod[i]) r = false;
      size *= dim[i];
    }
    return r;
  }


  //! pretty-prints the IdxSpec on the specified file.
  void pretty(FILE *f);
  void pretty(std::ostream& out);

  //! select: return a new IdxSpec corresponding to
  //! a slice of the current IdxSpec with slice i
  //! in dimension d. In other words, if m is an
  //! IdxSpec of order 2 of size (10,4), the call
  //! IdxSpec p = m.select(0,3) will set p to
  //! an IdxSpec or order 1 containing the 4-th
  //! row of m.
  IdxSpec select(int d, intg i);
  //! select_into: same as select, but modifies an existing
  //! IdxSpec instead of returning a new one.
  intg select_into(IdxSpec *dst, int d, intg n);
  //! select_inplace: same as select, but modifies the
  //! current IdxSpec.
  intg select_inplace(int d, intg i);

  //! narrow: return a new IdxSpec in which the d-th
  //! dimension has been reduced to size s, starting
  //! at item o. In other words, if m is an
  //! IdxSpec of order 2 of size (10,4), the call
  //! IdxSpec p = m.narrow(0,6,2) will set p to
  //! an IdxSpec or order 2 of size (6,4) whose rows
  //! are rows 2 to 7 of m.
  IdxSpec narrow(int d, intg s, intg o);
  //! narrow_into: same as narrow, but modifies an existing
  //! IdxSpec instead of returning a new one.
  intg narrow_into(IdxSpec *dst, int d, intg s, intg o);
  //! narrow_inplace: same as narrow, but modifies the
  //! current IdxSpec.
  intg narrow_inplace(int d, intg s, intg o);

  //! transpose: transpose dimensions d1 and d2.
  IdxSpec transpose(int d1, int d2);
  //! transpose all dimensions through permutation matrix p.
  IdxSpec transpose(int *p);
  //! same as transpose, but modifies an existing
  //! IdxSpec instead of returning a new one.
  int transpose_into(IdxSpec *dst, int d1, int d2);
  //! same as transpose, but modifies an existing
  //! IdxSpec instead of returning a new one.
  int transpose_into(IdxSpec *dst, int *p);
  //! transpose_inplace: same as transpose, but modifies the
  //! current IdxSpec.
  int transpose_inplace(int d1, int d2);
  //! transpose_inplace: same as transpose, but modifies the
  //! current IdxSpec.
  int transpose_inplace(int *p);

  //! unfold: prepare an IdxSpec for a convolution.
  //! Returns an idx on the same storage as m (pointing to the
  //! same data) with an added dimension at the end obtained by
  //! "unfolding" the n -th dimension. The size of the new dimension
  //! is k. This essentially manipulates the mod array to make
  //! convolutions look like matrix-vector multiplies.
  IdxSpec unfold(int d, intg k, intg s);
  //! same as unfold, but modifies an existing
  //! IdxSpec instead of returning a new one.
  intg unfold_into(IdxSpec *dst, int d, intg k, intg s);
  //! unfold_into: same as unfold, but modifies the
  //! current IdxSpec.
  intg unfold_inplace(int d, intg k, intg s);

  // horrible syntax to declare a template class friend.
  // isn't C++ wonderful?
  template <class T> friend class IdxIter;
  template <class T> friend class IdxLooper;
  template <class T> friend class Idx;
  template <class T> friend class ScalarIter_Base;
  template <class T> friend class ScalarIter;
  template <class T> friend class ReverseScalarIter;
  template <class T> friend class DimIter_Base;
  template <class T> friend class DimIter;
  template <class T> friend class ReverseDimIter;


  friend bool same_dim(IdxSpec &s1, IdxSpec &s2);
};

//! return true if two idxspec have the same dimensions,
//! i.e. if all their dimensions are equal (regardless
//! of strides).
bool same_dim(IdxSpec &s1, IdxSpec &s2);


////////////////////////////////////////////////////////////////

//! Idx: main tensor class. This can represent vectors,
//! matrices, and tensors up to 8 dimensions.
//! An Idx is merely an access structure that
//! points to the data. Several Idx can point to the same
//! data.
template <class T> class Idx {

	//! Pretty-prints elements to a stream.
	//friend std::ostream& operator<<( std::ostream& out, Idx<T>& tensor );

 private:

  //! pointer to the Srg structure that contains the data.
  Srg<T> *storage;

  //! increase size of storage to fit the current dimension
  void growstorage();

  //! increase size of storage to fit the current dimension + a given size
  void growstorage_chunk(intg s_chunk);

  //! Implementation of public printElems() method.
  void printElems_impl( int indent, std::ostream& );

 protected:
  //! fake constructor that does nothing.
  //! This is called by the IdxLooper constructor
  Idx(dummyt *dummy);

 public:

  // STL-like typedefs
  typedef T value_type;
  typedef ScalarIter<value_type> scalar_iterator;
  typedef ReverseScalarIter<value_type> reverse_scalar_iterator;
  typedef DimIter<value_type> dimension_iterator;
  typedef ReverseDimIter<value_type> reverse_dimension_iterator;

  /* -----  STL-style iterator creators  ----- */

  //! Returns an iterator over all the elements.
  scalar_iterator scalars_begin();

  //! Returns a terminated iterator over all the elements.
  scalar_iterator scalars_end();

  //! Returns a reversed-order iterator over all the elements.
  reverse_scalar_iterator scalars_rbegin();

  //! Returns a terminated reversed-order iterator over all the elements.
  reverse_scalar_iterator scalars_rend();

  //! Returns an iterator for the given dimension.
  dimension_iterator dim_begin( int dim );

  //! Returns a terminated iterator for the given dimension.
  dimension_iterator dim_end( int dim );

  //! Returns an iterator for the given dimension.
  reverse_dimension_iterator dim_rbegin( int dim );

  //! Returns a terminated iterator for the given dimension.
  reverse_dimension_iterator dim_rend( int dim );


  //! IdxSpec that contains the order,
  //! offset, dimensions and strides.
  IdxSpec spec;