#include "AffineTransformation.h"
#include "SimplexPoint.h"
#include "../Index/MatrixIdx.h"

namespace Domain
{

FILE* AffineTransformation::datafp; 



AffineTransformation* AffineTransformation::getExample()
{
	AffineTransformation *t;
	
	double *dv1,*dv2,*dv3;
	dv1 = new double[3];
	dv1[0]= 1.0;
	dv1[1]= 0.0;
	dv1[2]= 0.0;

	dv2 = new double[3];
	dv2[0]= 0.0;
	dv2[1]= 1.0;
	dv2[2]= 0.0;

	dv3 = new double[3];
	dv3[0]= 1;
	dv3[1]= 1;
	dv3[2]= -1.0;

	Domain::SimplexPoint* spx1 = Domain::SimplexPoint::fromDoubleVector(dv1,2);
	//spx1->print(stderr);

	Domain::SimplexPoint* spx2 = Domain::SimplexPoint::fromDoubleVector(dv2,2);
	//spx2->print(stderr);

	Domain::SimplexPoint* spx3 = Domain::SimplexPoint::fromDoubleVector(dv3,2);
	//spx3->print(stderr);
	

	
	std::vector<Domain::SimplexPoint*> v;
	std::vector<Domain::SimplexPoint*>::iterator iter;
	
	Domain::SimplexPoint** spxvec;
	spxvec = new Domain::SimplexPoint*[3];
	spxvec[0] = spx1;
	spxvec[1] = spx2;
	spxvec[2] = spx3;
	
	
	//T = Domain::AffineTransformation::getIdentity(2,2);
	//T->v->coord[0] = 
	
	t = Domain::AffineTransformation::getFromImage(spxvec,3);
	
	return t;

}


AffineTransformation* AffineTransformation::getExample2()
{
	AffineTransformation *t;
	
	double *dv1,*dv2,*dv3;
	dv1 = new double[2];
	dv1[0]= 0.2;
	dv1[1]= 0.2;


	dv2 = new double[2];
	dv2[0]= 0.8;
	dv2[1]= 0;


	dv3 = new double[2];
	dv3[0]= 0;
	dv3[1]= 0.8;


	Domain::SimplexPoint* spx1 = Domain::SimplexPoint::fromDoubleVector(dv1,1);
	//spx1->print(stderr);

	Domain::SimplexPoint* spx2 = Domain::SimplexPoint::fromDoubleVector(dv2,1);
	//spx2->print(stderr);

	Domain::SimplexPoint* spx3 = Domain::SimplexPoint::fromDoubleVector(dv3,1);
	//spx3->print(stderr);
	

	
	std::vector<Domain::SimplexPoint*> v;
	std::vector<Domain::SimplexPoint*>::iterator iter;
	
	Domain::SimplexPoint** spxvec;
	spxvec = new Domain::SimplexPoint*[3];
	spxvec[0] = spx1;
	spxvec[1] = spx2;
	spxvec[2] = spx3;
	
	
	//T = Domain::AffineTransformation::getIdentity(2,2);
	//T->v->coord[0] = 
	
	t = Domain::AffineTransformation::getFromImage(spxvec,3);
	
	return t;

}



AffineTransformation** AffineTransformation::getExample3()
{
	AffineTransformation **t;
	double *dv1,*dv2;
	Domain::SimplexPoint  *spx1,*spx2;
	Domain::SimplexPoint** spxvec;
	
	
	t = new AffineTransformation*[2];
	
	
	
	dv1 = new double[3];
	dv1[0]= 1;
	dv1[1]= 0.5;
	dv1[2]= 0;


	dv2 = new double[3];
	dv2[0]= 0;
	dv2[1]= 0.5;
	dv2[2]= 1;

	spx1 = Domain::SimplexPoint::fromDoubleVector(dv1,2);
	//spx1->print(stderr);

	spx2 = Domain::SimplexPoint::fromDoubleVector(dv2,2);
	

	
	
	spxvec = new Domain::SimplexPoint*[2];
	spxvec[0] = spx1;
	spxvec[1] = spx2;
	
	
	
	t[0] = Domain::AffineTransformation::getFromColumns(spxvec,1);
	//---------------------------------------------------------------
	

	dv1[0]= 0;
	dv1[1]= 1;
	dv1[2] = 0;
	
	dv2[0]= 1;
	dv2[1]= 0;
	dv2[2] = 1;
	
	spx1 = Domain::SimplexPoint::fromDoubleVector(dv1,2);
	//spx1->print(stderr);

	spx2 = Domain::SimplexPoint::fromDoubleVector(dv2,2);
	

	spxvec[0] = spx1;
	spxvec[1] = spx2;
	
	
	
	t[1] = Domain::AffineTransformation::getFromColumns(spxvec,1);
	
	
	
	
	
	
	
	
	
	
	
	
	
	return t;

}

AffineTransformation* AffineTransformation::getFromColumns(SimplexPoint** cols, int destDim)
{
	AffineTransformation *t;
	
	t = new  AffineTransformation(cols[0]->dim,destDim);
		
//	t->v = new SimplexPoint*[origDim+1];

	int i;

	for(i=0;i<=destDim;i++)
	{
		t->v[i] = cols[i]->clone();
		//t->v[i]->setCoord(i,1);
	}
	
	return t;

}

AffineTransformation* AffineTransformation::getFromImage(SimplexPoint** img, int size)
{
	AffineTransformation *t;
	
	t = new  AffineTransformation(img[0]->dim,size-1);
		
//	t->v = new SimplexPoint*[origDim+1];

	int i;

	for(i=0;i<size;i++)
	{
		t->v[i] = img[i]->clone();
		//t->v[i]->setCoord(i,1);
	}
	
	return t;

}

AffineTransformation::AffineTransformation(int origDim,int destDim )
{
	this->origDim = origDim;
	this->destDim = destDim;
	
	if(AffineTransformation::datafp==NULL)
		AffineTransformation::datafp = stderr;
	
	
	this->v = new SimplexPoint*[destDim];
	//int i;
	//for(i=0;i<origDim;i++)
	//{
	//	this->v[i] = new SimplexPoint(origDim);
	//}
}


double AffineTransformation::powerMatrixIdx2(Index::MatrixIdx *e)
{
	int i,j;
	
	double res=1;
	double exp,base;
	
	//this->print(stderr);
	
	for(i=0;i<= this->destDim;i++)
	{
		for(j=0;j<=this->origDim;j++)
		{
			exp = (e->line[j])->idx[i];
			base = (this->v[i]->coord[j]);
			
			if(base==0)
			{
				if(exp!=0)
					return 0;
			}
			else
				res = res* pow(base , exp);
		}
			
	}
	
	
	return res;
	
}


double AffineTransformation::powerMatrixIdx(Index::MatrixIdx *e)
{
	int i,j;
	
	double res=1;
	double exp,base;
	
	//this->print(stderr);
	
	for(i=0;i<= this->destDim;i++)
	{
		for(j=0;j<=this->origDim;j++)
		{
			exp = (e->line[i])->idx[j];
			base = (this->v[j]->coord[i]);
			
			res = res* pow(base , exp);
		}
			
	}
	
	
	return res;
	
}


Float* AffineTransformation::mppowerMatrixIdx(Index::MatrixIdx *e)
{
	int i,j;
	
	Float *res;
	Float *Faux;
	
	int exp;
	double base;
	
	
	res = new Float((double)1.f);
	
	
	//this->print(stderr);
	
	for(i=0;i<= this->destDim;i++)
	{
		for(j=0;j<=this->origDim;j++)
		{
			exp = (e->line[i])->idx[j];
			
			base = (this->v[j]->coord[i]);
			
			Faux = Float::ipow(base,exp);
			
			res->Float::mul_(Faux);
			
			delete Faux;
		}
			
	}
	
	
	return res;
	
}


AffineTransformation* AffineTransformation::getIdentity(int origDim,int destDim )
{
	AffineTransformation *t;
	
	t = new  AffineTransformation(origDim,destDim);
		
	//t->v = new SimplexPoint*[origDim+1];
	int i;
	for(i=0;i<= destDim;i++)
	{
		t->v[i] = new SimplexPoint(origDim,0);
		t->v[i]->setCoord(i,1);
	}
	
	return t;
}


SimplexPoint* AffineTransformation::transform(SimplexPoint* u)
{
	SimplexPoint* v;
	
	int i;
	
	if(u->dim!=this->origDim)
		return NULL;
	
	v = new SimplexPoint(this->destDim);
	
	for(i=0;i<= this->destDim;i++)
		v->coord[i]= u->dot( this->v[i] );		

	
	
	return v;
	
}




void AffineTransformation::print()
{
	this->print(AffineTransformation::datafp);
}

void AffineTransformation::print(FILE* datafp)
{
	int i,j;
	
	for(i=0; i<= this->origDim; i++)
	{
		for(j=0;j <= this->destDim; j++)
			fprintf(datafp," %g",this->v[j]->coord[i]);
		
		fprintf(datafp,"\n");
		
	}
	
}


AffineTransformation::~AffineTransformation()
{
}

}