idx.hpp

Gaussian Mixture Algorithm

#define idx_bloop4(dst0,src0,type0,dst1,src1,type1,dst2,src2,type2,dst3,src3,type3) \
  idx_checkdim4_all(src0, src1, src2, src3, 0); \
  IdxLooper<type0> dst0(src0,0);					\
  IdxLooper<type1> dst1(src1,0);					\
  IdxLooper<type2> dst2(src2,0);					\
  IdxLooper<type3> dst3(src3,0);					\
  for ( ; dst0.notdone(); dst0.next(), dst1.next(), dst2.next(), dst3.next())

// eloop macros

#define idx_eloop1(dst0,src0,type0) \
  IdxLooper<type0> dst0(src0,src0.order()-1);	\
  for ( ; dst0.notdone(); dst0.next())

#define idx_eloop2(dst0,src0,type0,dst1,src1,type1)				\
  if ( (src0).dim((src0).order()) != (src1).dim((src1).order()) ) ylerror("incompatible Idxs for bloop\n"); \
  IdxLooper<type0> dst0(src0,(src0).order()-1); \
  IdxLooper<type1> dst1(src1,(src1).order()-1); \
  for ( ; dst0.notdone(); dst0.next(), dst1.next())


////////////////////////////////////////////////////////////////
// aloop macros: loop over all elements

// Loops over all elements of an idx. This takes a pointer to
// the data type of idx elements, and a blank IdxIter object:
// idx_aloop1(data_pointer,idxiter,&idx) { do_stuff(data_pointer); }
// Example of use: add 1 to all element of m:
//  Idx<double> m(3,4);
//  IdxIter<double> p;
//  idx_aloop1(p,&m) { *p += 1; }
#define idx_aloop1_on(itr0,src0)					\
  for ( itr0.init(src0); itr0.notdone(); itr0.next())


// this loops simultaneously over all elements of 2 Idxs.
// The two Idxs can have different structures as long as they have
// the same total number of elements.
#define idx_aloop2_on(itr0,src0,itr1,src1)			\
	idx_checknelems2_all(src0, src1); \
  for ( itr0.init(src0), itr1.init(src1); itr0.notdone(); itr0.next(), itr1.next())

#define idx_aloop3_on(itr0,src0,itr1,src1,itr2,src2)	\
	idx_checknelems3_all(src0, src1, src2); \
  for (itr0.init(src0), itr1.init(src1), itr2.init(src2); itr0.notdone(); itr0.next(), itr1.next(), itr2.next())

// high level aloop macros.
// These should be enclosed in braces, to avoid name clashes
#define idx_aloop1(itr0,src0,type0) \
  IdxIter<type0> itr0; \
  idx_aloop1_on(itr0,src0)

#define idx_aloop2(itr0,src0,type0,itr1,src1,type1)		\
  IdxIter<type0> itr0; \
  IdxIter<type0> itr1; \
	idx_checknelems2_all(src0, src1); \
  for (itr0.init(src0), itr1.init(src1); itr0.notdone(); itr0.next(), itr1.next())

#define idx_aloop3(itr0,src0,type0,itr1,src1,type1,itr2,src2,type2) \
  IdxIter<type0> itr0; \
  IdxIter<type0> itr1; \
  IdxIter<type0> itr2; \
	idx_checknelems3_all(src0, src1, src2); \
  for (itr0.init(src0), itr1.init(src1), itr2.init(src2); itr0.notdone(); itr0.next(), itr1.next(), itr2.next())

#endif // if USING_STL_ITERS, else

////////////////////////////////////////////////////////////////
// Idx methods

template <class T> void Idx<T>::growstorage() {
  if ( storage->growsize(spec.footprint()) < 0) {
    ylerror("cannot grow storage");
  }
}

template <class T> void Idx<T>::growstorage_chunk(intg s_chunk){
	if( storage->growsize_chunk(spec.footprint(), s_chunk) < 0) {
		ylerror("cannot grow storage");
	}
}

/* the constructor doesn't exist and I didn't find one to replace it
template<typename T>
void Idx<T>::printElems(FILE* filePtr){
	this->printElems(std::ofstream(filePtr));
}
*/

template <typename T>
void Idx<T>::printElems( std::ostream& out ){
	printElems_impl(0, out);
}

template <typename T>
void Idx<T>::printElems(){
	this->printElems( std::cout );
}

template<class T> inline T printElems_impl_cast(T val) {
	return val;
}

// specialization for ubyte to print as unsigned integers.
inline unsigned int printElems_impl_cast(ubyte val) {
	return (unsigned int) val;
}

template <typename T>
void Idx<T>::printElems_impl( int indent, std::ostream& out ){
	static const std::string lbrace = "[";
	static const std::string rbrace = "]";
	static const std::string sep = " ";
	std::ostringstream oss;
	for( unsigned int ii = 0; ii < lbrace.length(); ++ii ){
		oss<<" ";
	}
	const std::string tab(oss.str());

	// printing a 0-dimensional tensor
	if( order() == 0 ){
		out<<lbrace<<"@"<<sep<< printElems_impl_cast(get()) <<sep<<rbrace<<"\n";
	}
	// printing a 1-D tensor
	else if( order() == 1 ){
		out<<lbrace<<sep;
		for( int ii = 0; ii < dim(0); ++ii ){
			out<< printElems_impl_cast(get(ii)) <<sep;
		}
		out<<rbrace<<"\n";
	}
	// printing a multidimensional tensor
	else{

		// opening brace
		out<<lbrace;

		// print subtensors.
		Idx<T> subtensor(storage, spec.offset);
		for( int dimInd = 0; dimInd < dim(0); ++dimInd ){

			// only print indent if this isn't the first subtensor.
			if( dimInd > 0 ){
				for( int ii = 0; ii < indent+1; ++ii ){
					out<<(tab);
				}
			}

			// print subtensor
			spec.select_into(&subtensor.spec, 0, dimInd);
			subtensor.printElems_impl( indent+1, out );

			// only print the newline if this isn't the last subtensor.
			if( dimInd < (dim(0)-1)){
				out<<"\n";
			}

		}

		// closing brace
		out<<rbrace<<"\n";
	}
}

template <class T> void Idx<T>::pretty(FILE *f) {
  fprintf(f,"Idx: at address %ld\n",(intg)this);
  fprintf(f,"  storage=%ld (size=%ld)\n",(intg)storage,storage->size());
  spec.pretty(f);
}

template <class T> void Idx<T>::pretty(std::ostream& out){
	out << "Idx: at address " << (intg)this << "\n";
	out << "  storage=" <<  (intg)storage << "(size=" << storage->size() << "\n";
	spec.pretty(out);
}

template <class T> Idx<T>::~Idx() {
  DEBUG("Idx::destructor %ld\n",long(this));
  storage->unlock();
}

// fake constructor called by IdxLooper constructor
template <class T> Idx<T>::Idx(dummyt *dummy) {
  spec.dim = NULL;
  spec.mod = NULL;
  storage = NULL;
}

//template <typename T>
//Idx<T>::Idx( Idx<T>& other )
//	:storage(other.storage),
//	 spec(other.spec)
// {
//	storage->lock();
// }
//
//template <typename T>
//Idx<T>::Idx( const Idx<T>& other )
//	:storage(other.storage),
//	 spec(other.spec)
// {
//	storage->lock();
// }

template <class T> Idx<T>::Idx(Srg<T> *srg, IdxSpec &s) {
  spec = s;
  storage = srg;
  growstorage();
  storage->lock();
}

template <class T> Idx<T>::Idx(Srg<T> *srg, intg o) : spec(o) {
  storage = srg;
  growstorage();
  storage->lock();
}

template <class T> Idx<T>::Idx() : spec(0) {
  storage = new Srg<T>();
  growstorage();
  storage->lock();
}

template <class T> Idx<T>::Idx(Srg<T> *srg, intg o, intg size0) : spec(o,size0) {
  storage = srg;
  growstorage();
  storage->lock();
}

template <class T> Idx<T>::Idx(intg size0) : spec(0,size0) {
  storage = new Srg<T>();
  growstorage();
  storage->lock();
}

template <class T> Idx<T>::Idx(Srg<T> *srg, intg o, intg size0, intg size1) : spec(o,size0,size1) {
  storage = srg;
  growstorage();
  storage->lock();
}

template <class T> Idx<T>::Idx(intg size0, intg size1) : spec(0,size0,size1) {
  storage = new Srg<T>();
  growstorage();
  storage->lock();
}

template <class T> Idx<T>::Idx(Srg<T> *srg, intg o, intg size0, intg size1, intg size2) : spec(o,size0,size1,size2) {
  storage = srg;
  growstorage();
  storage->lock();
}

template <class T> Idx<T>::Idx(intg size0, intg size1, intg size2) : spec(0,size0,size1,size2) {
  storage = new Srg<T>();
  growstorage();
  storage->lock();
}

template <class T> Idx<T>::Idx(Srg<T> *srg, intg o, intg s0, intg s1, intg s2, intg s3, intg s4, intg s5, intg s6, intg s7) : spec(o,s0,s1,s2,s3,s4,s5,s6,s7) {
  storage = srg;
  growstorage();
  storage->lock();
}

template <class T> Idx<T>::Idx(intg s0, intg s1, intg s2, intg s3, intg s4, intg s5, intg s6, intg s7) : spec(0,s0,s1,s2,s3,s4,s5,s6,s7) {
  storage = new Srg<T>();
  growstorage();
  storage->lock();
}


template <class T> intg Idx<T>::setoffset(intg o) {
  if (o<0) { ylerror("Idx::setoffset: offset must be positive"); }
  if (o > spec.offset) {
    spec.setoffset(o);
    growstorage();
    return o;
  } else {
    spec.setoffset(o);
    return o;
  }
}

template <class T> void Idx<T>::resize(intg s0, intg s1, intg s2, intg s3, intg s4, intg s5, intg s6, intg s7) {
  spec.resize(s0,s1,s2,s3,s4,s5,s6,s7);
  growstorage();
}

template <class T> void Idx<T>::resize_chunk(intg s_chunk, intg s0, intg s1, intg s2, intg s3, intg s4, intg s5, intg s6, intg s7) {
  spec.resize(s0,s1,s2,s3,s4,s5,s6,s7);
  growstorage_chunk(s_chunk);
}

//template<typename T, typename SizeIter>
//void Idx<T>::resize( SizeIter& sizesBegin, SizeIter& sizesEnd ){
//	const int ndims = std::distance( sizesBegin, sizesEnd );
//	if ( ndims > MAXDIMS ){
//		std::ostringstream oss;
//		oss<<"Number of dimensions ("<<ndims<<") exceeds MAX_DIMS ("<<MAXDIMS<<")";
//		ylerror(oss.str());
//	}
//
//	std::vector<T> sizes(ndims+1);
//	std::copy(sizesBegin, sizesEnd, sizes.begin());
//	sizes.back() = -1;
//	spec.resize( sizesBegin, sizesEnd );
//	growstorage();
//}


template <class T> Idx<T> Idx<T>::select(int d, intg i) {
  Idx<T> r(storage,spec.getoffset());