spidx.hpp

来自「Gaussian Mixture Algorithm」· HPP 代码 · 共 626 行 · 第 1/2 页

/*
 * spIdx.hpp
 *
 *      Author: cyril Poulet
 */

#ifndef SPIDX_HPP_
#define SPIDX_HPP_

namespace ebl {

template<class T> spIdx<T>::spIdx(){
	ndim = 0;
	nelem = 0;
	mydims = NULL;
	myindex = NULL;
	myvalues = NULL;
}


template<class T> spIdx<T>::spIdx(intg &dims, intg nelems, Idx<T> &ind, Idx<T> &val):
	mydims(dims),
	myindex(ind),
	myvalues(val)
{
	nelem = nelems;
	ndim = myindex->dim(1);
}


template<class T> spIdx<T>::spIdx(intg Nelemmax, int Ndim, const intg* dims){
	intg N = Nelemmax;
	if(Nelemmax == 0) N = 1;
	ndim = Ndim;
	mydims = new intg[ndim];
	for(int i = 0; i < ndim; i++) mydims[i] = dims[i];
	myvalues = new Idx<T>(N);
	myindex = new Idx<intg>(N, Ndim);
	nelem = 0;
}

template<class T> spIdx<T>::spIdx(intg Nelemmax, int Ndim, const intg* dims, Srg<T> *srg, intg o){
	intg N = Nelemmax;
	if(Nelemmax == 0) N = 1;
	ndim = Ndim;
	mydims = new intg[ndim];
	for(int i = 0; i < ndim; i++) mydims[i] = dims[i];
	myvalues = new Idx<T>(srg, o, N);
	myindex = new Idx<intg>(N, Ndim);
	nelem = 0;
}

template<class T> spIdx<T>::spIdx(intg Nelemmax, int Ndim, const intg* dims, Srg<T> *srg):
	spIdx(Nelemmax, Ndim, dims, srg, srg->footprint())
	{}

template<class T> spIdx<T>::spIdx(intg Nelemmax, intg s0, intg s1, intg s2, intg s3, intg s4, intg s5, intg s6, intg s7){
	intg N = Nelemmax;
	if(Nelemmax == 0) N = 1;
	bool ndimset = false;
	mydims = NULL;
	try {
		if (s7>=0) {
			if (!ndimset) { ndim = 8; ndimset = true; mydims = new intg[8];}
			mydims[7] = s7;
		} else { if (ndimset) { throw(-8); } }
		if (s6>=0) {
			if (!ndimset) { ndim = 7; ndimset = true; mydims = new intg[7];}
			mydims[6] = s6;
		} else { if (ndimset) { throw(-7); } }
		if (s5>=0) {
			if (!ndimset) { ndim = 6; ndimset = true; mydims = new intg[6];}
			mydims[5] = s5;
		} else { if (ndimset) { throw(-6); } }
		if (s4>=0) {
			if (!ndimset) { ndim = 5; ndimset = true; mydims = new intg[5];}
			mydims[4] = s4;
		} else { if (ndimset) { throw(-5); } }
		if (s3>=0) {
			if (!ndimset) {ndim = 4; ndimset = true; mydims = new intg[4];}
			mydims[3] = s3;
		} else { if (ndimset) { throw(-4); } }
		if (s2>=0) {
			if (!ndimset) { ndim = 3; ndimset = true; mydims = new intg[3];}
			mydims[2] = s2;
		} else { if (ndimset) { throw(-3); } }
		if (s1>=0) {
			if (!ndimset) { ndim = 2; ndimset = true; mydims = new intg[2];}
			mydims[1] = s1;
		} else { if (ndimset) { throw(-2); } }
		if (s0>=0) {
			if (!ndimset) { ndim = 1; ndimset = true; mydims = new intg[1];}
			mydims[0] = s0;
		} else { if (ndimset) { throw(-1); } }
		if (!ndimset) { throw(0); }
	}
	catch(int v) { ylerror("spIdx: bad dimensions in constructor"); }
	myvalues = new Idx<T>(N);
	myindex = new Idx<intg>(N, ndim);
	nelem = 0;
}

template<class T> spIdx<T>::spIdx(intg Nelemmax, intg o, Srg<T> *srg, intg s0, intg s1, intg s2, intg s3, intg s4, intg s5, intg s6, intg s7){
	intg N = Nelemmax;
	if(Nelemmax == 0) N = 1;
	bool ndimset = false;
	mydims = NULL;
	try {
		if (s7>=0) {
			if (!ndimset) { ndim = 8; ndimset = true; mydims = new intg[8];}
			mydims[7] = s7;
		} else { if (ndimset) { throw(-8); } }
		if (s6>=0) {
			if (!ndimset) { ndim = 7; ndimset = true; mydims = new intg[7];}
			mydims[6] = s6;
		} else { if (ndimset) { throw(-7); } }
		if (s5>=0) {
			if (!ndimset) { ndim = 6; ndimset = true; mydims = new intg[6];}
			mydims[5] = s5;
		} else { if (ndimset) { throw(-6); } }
		if (s4>=0) {
			if (!ndimset) { ndim = 5; ndimset = true; mydims = new intg[5];}
			mydims[4] = s4;
		} else { if (ndimset) { throw(-5); } }
		if (s3>=0) {
			if (!ndimset) {ndim = 4; ndimset = true; mydims = new intg[4];}
			mydims[3] = s3;
		} else { if (ndimset) { throw(-4); } }
		if (s2>=0) {
			if (!ndimset) { ndim = 3; ndimset = true; mydims = new intg[3];}
			mydims[2] = s2;
		} else { if (ndimset) { throw(-3); } }
		if (s1>=0) {
			if (!ndimset) { ndim = 2; ndimset = true; mydims = new intg[2];}
			mydims[1] = s1;
		} else { if (ndimset) { throw(-2); } }
		if (s0>=0) {
			if (!ndimset) { ndim = 1; ndimset = true; mydims = new intg[1];}
			mydims[0] = s0;
		} else { if (ndimset) { throw(-1); } }
		if (!ndimset) { throw(0); }
	}
	catch(int v) { ylerror("spIdx: bad dimensions in constructor"); }
	myvalues = new Idx<T>(srg, o, N);
	myindex = new Idx<intg>(N, ndim);
	nelem = 0;
}

template<class T> spIdx<T>::spIdx(intg Nelemmax, Srg<T> *srg, intg s0, intg s1, intg s2, intg s3, intg s4, intg s5, intg s6, intg s7):
	spIdx(Nelemmax, srg->footprint(), srg, s0, s1, s2, s3, s4, s5, s6, s7)
	{}

template<class T> spIdx<T>::spIdx( const spIdx<T>& other ){
	this->myindex = new Idx<intg>(*(other.myindex));
	this->myvalues = new Idx<T>(*(other.myvalues));
	this->ndim = other.ndim;
	this->nelem = other.nelem;
	this->mydims = new intg[ndim];
	for(int i = 0; i < ndim; i++) this->mydims[i] = other.mydims[i];
	}

template<class T> spIdx<T>::~spIdx(){
	delete mydims;
	delete myindex;
	delete myvalues;
}


template<class T> spIdx<T> spIdx<T>::select(int d, intg i){
	if(ndim == 1) {ylerror("select : cannot use select on a 1-dim spIdx"); return *this;}
	if((d >= ndim) || (d < 0)) { ylerror("select : illegal dimension index"); return *this;}
	if((i > mydims[d]) || (i < 0 )) { ylerror("select : illegal offset"); return *this;}

	intg newdim[this->ndim-1];
	for(int j =0; j< d; j++) newdim[j] = mydims[j];
	for(int j = d; j< this->ndim -1; j++) newdim[j] = mydims[j+1];
	spIdx<T> res(0, this->ndim-1, newdim);
	if(nelem == 0) return res;
	int indx = 0;
	{ idx_bloop1(ind, *myindex, intg){
		if((ind.get(d) == i) && (indx < nelem)){
			int i0 = (d==0)? ind.get(1) : ind.get(0);
			int i1 = (ndim-1 < 2)? -1 : ((d<=1)? ind.get(2) : ind.get(1));
			int i2 = (ndim-1 < 3)? -1 : ((d<=2)? ind.get(3) : ind.get(2));
			int i3 = (ndim-1 < 4)? -1 : ((d<=3)? ind.get(4) : ind.get(3));
			int i4 = (ndim-1 < 5)? -1 : ((d<=4)? ind.get(5) : ind.get(4));
			int i5 = (ndim-1 < 6)? -1 : ((d<=5)? ind.get(6) : ind.get(5));
			int i6 = (ndim-1 < 7)? -1 : ((d<=6)? ind.get(7) : ind.get(6));
			res.set(myvalues->get(indx), i0, i1, i2, i3, i4, i5, i6);
		}
		indx++;
	}}
	return res;
}

template<class T> spIdx<T> spIdx<T>::narrow(int d, intg s, intg o){
	try {
		if ((d < 0) || (d>=ndim)) throw("narrow: illegal dimension index");
		if ((o < 0)||(s < 1)||(s+o > mydims[d])) throw("narrow: illegal size/offset");
	}
	catch(const char *s) { ylerror(s); return *this;}

	intg newdim[this->ndim];
	for(int j =0; j< this->ndim; j++){
		if(j == d) newdim[j] = s;
		else newdim[j] = mydims[j];
	}
	spIdx<T> res(0, this->ndim, newdim);
	if(nelem == 0) return res;
	int indx = 0;
	{ idx_bloop1(ind, *myindex, intg){
		int bla = ind.get(d);
		if((bla >= o) && (bla < o+s) && (indx < nelem)){
			int i0 = ind.get(0) - (d==0 ? o : 0);
			int i1 = (ndim < 2)? -1 : ind.get(1) - (d==1 ? o : 0);
			int i2 = (ndim < 3)? -1 : ind.get(2) - (d==2 ? o : 0);
			int i3 = (ndim < 4)? -1 : ind.get(3) - (d==3 ? o : 0);
			int i4 = (ndim < 5)? -1 : ind.get(4) - (d==4 ? o : 0);
			int i5 = (ndim < 6)? -1 : ind.get(5) - (d==5 ? o : 0);
			int i6 = (ndim < 7)? -1 : ind.get(6) - (d==6 ? o : 0);
			int i7 = (ndim < 8)? -1 : ind.get(7) - (d==7 ? o : 0);
			res.set(myvalues->get(indx), i0, i1, i2, i3, i4, i5, i6, i7);
		}
		indx++;
	}}
	return res;
}

template<class T> void spIdx<T>::resize(intg s0, intg s1, intg s2, intg s3, intg s4, intg s5, intg s6, intg s7){
	if(!has_right_dimension(s0, s1, s2, s3, s4, s5, s6, s7)) { ylerror("resize : wrong number of parameters"); return;}
	intg dims[8] = {s0, s1, s2, s3, s4, s5, s6, s7};
	for(int i = 0; i < ndim; i++) if(dims[i] <=0){ ylerror("resize : wrong parameters"); return;}
	if(!((s0 >= mydims[0])
			|| ((ndim >= 2) && (s1 >= mydims[1]))
			|| ((ndim >= 3) && (s2 >= mydims[2]))
			|| ((ndim >= 4) && (s3 >= mydims[3]))
			|| ((ndim >= 5) && (s4 >= mydims[4]))
			|| ((ndim >= 6) && (s5 >= mydims[5]))
			|| ((ndim >= 7) && (s6 >= mydims[6]))
			|| ((ndim >= 8) && (s7 >= mydims[7])))){
		for(intg i = 0; i< nelem; i++) {
			intg i0 = myindex->get(i, 0);
			intg i1 = (ndim >= 2)? myindex->get(i, 1) : -1;
			intg i2 = (ndim >= 3)? myindex->get(i, 2) : -1;
			intg i3 = (ndim >= 4)? myindex->get(i, 3) : -1;
			intg i4 = (ndim >= 5)? myindex->get(i, 4) : -1;
			intg i5 = (ndim >= 6)? myindex->get(i, 5) : -1;
			intg i6 = (ndim >= 7)? myindex->get(i, 6) : -1;
			intg i7 = (ndim >= 8)? myindex->get(i, 7) : -1;
			if((i0 >= s0)
					|| ((ndim >= 2) && (i1 >= s1))
					|| ((ndim >= 3) && (i2 >= s2))
					|| ((ndim >= 4) && (i3 >= s3))
					|| ((ndim >= 5) && (i4 >= s4))
					|| ((ndim >= 6) && (i5 >= s5))
					|| ((ndim >= 7) && (i6 >= s6))
					|| ((ndim >= 8) && (i7 >= s7))) {
				set( BACKGROUND, i0, i1, i2, i3, i4, i5, i6, i7);
				i--;
			}
		}
	}
	if(ndim >= 1) mydims[0] = s0;
	if(ndim >= 2) mydims[1] = s1;
	if(ndim >= 3) mydims[2] = s2;
	if(ndim >= 4) mydims[3] = s3;
	if(ndim >= 5) mydims[4] = s4;
	if(ndim >= 6) mydims[5] = s5;
	if(ndim >= 7) mydims[6] = s6;
	if(ndim >= 8) mydims[7] = s7;
}

template<class T> spIdx<T> spIdx<T>::transpose(int d1, int d2){
	if((d1 < 0) || (d1 > ndim-1) || (d2 < 0) || (d2 > ndim-1)) { ylerror("transpose : selected dimensions do not exist in spIdx\n"); return *this;}
	spIdx<T> res(myvalues->dim(0), ndim, mydims);
	idx_copy(*this, res);
	res.mydims[d1] = this->mydims[d2];
	res.mydims[d2] = this->mydims[d1];
	if(nelem == 0) return res;
	{ idx_bloop2(indt, *(this->myindex), intg, indr, *(res.myindex), intg){
		indr.set(indt.get(d2), d1);
		indr.set(indt.get(d1), d2);
	}}

	return res;
}

template<class T> spIdx<T> spIdx<T>::transpose(int *p){
	for (int i=0; i<ndim; i++) {
		if ((p[i] < 0) || (p[i] >= ndim)) ylerror("tranpose: illegal dimension index");
	}
	spIdx<T> res(myvalues->dim(0), this->ndim, this->mydims);
	res.set_nelements(nelem);
	for(int i = 0; i<ndim; i++){
		res.mydims[p[i]] = this->mydims[i];
	}
	idx_copy(*myvalues, *(res.myvalues));
	if(nelem == 0) return res;
	{ idx_bloop2(th, *(this->myindex), intg, re, *(res.myindex), intg){
		for(int i = 0; i<ndim; i++){
			re.set(th.get(p[i]), i);
		}
	}}

	return res;
}

template<class T> bool spIdx<T>::has_right_dimension(intg i0, intg i1, intg i2, intg i3, intg i4, intg i5, intg i6, intg i7){
	bool ndimset = false;
	int dim;
	try {
		if (i7>=0) {
			if (!ndimset) { dim = 8; ndimset = true; }