// marVector.cpp: implementation of the marVector class.
//
//////////////////////////////////////////////////////////////////////
#include "marVector.h"
#include <assert.h>
#include <math.h>


//////////////////////////////////////////////////////////////////////
// Constructeurs
//////////////////////////////////////////////////////////////////////
marVector::marVector(size_t s)
	:shallowCopy(false)
{
	if (s!=0)
	{
		_data=new double[s];
		_size=s;
	}
	else
	{
		_data=NULL;
		_size=0;
		assert(false);
	}
}

marVector::marVector(const marVector &v)
	:shallowCopy(false)
{
	_size=v._size;
	_data=new double[_size];
	int i;
	for (i=0;i<_size;i++)
	{
		_data[i]=v._data[i];
	}
}

marVector::marVector(double*v,size_t s)
	:shallowCopy(true)
{
	if (s!=0)
	{
		_size=s;
		_data=v;
// PS -> 		_data=new double[_size];
// PS -> 		for (int i=0;i<_size;i++)
// PS -> 		{
// PS -> 			_data[i]=v[i];
// PS -> 		}
		
	}
	else
	{
		_data=NULL;
		_size=0;
		assert(false);
	}
	
}

//////////////////////////////////////////////////////////////////////
// Destructeur
//////////////////////////////////////////////////////////////////////
marVector::~marVector()
{
	if (!shallowCopy) 
	{
		delete _data;
	}
	else
	{
		_data=NULL;
	}
}

//////////////////////////////////////////////////////////////////////
// Opérateurs
//////////////////////////////////////////////////////////////////////
// ---------------------------------------------------------------------
/*Affichage*/
std::ostream& operator<<(std::ostream& os,const marVector& v)
{
	int i;
    for(i=0;i < v.size();i++) 
	{
		os << " " << v(i);
	}
	os <<"\n";
    return(os);
}

// ---------------------------------------------------------------------
/*Casting*/
marVector::operator double*() const 
{ 
	return(_data);
}

// ---------------------------------------------------------------------
/*Indexation*/
double& marVector::operator()(size_t i) 
{ 
	return(_data[i]);
}

const double& marVector::operator()(size_t i) const
{
	return(_data[i]);
}

// ---------------------------------------------------------------------
/*Affectation*/
marVector& marVector::operator=(const marVector& o)
{
	_size=o._size;
	if (!(_data==NULL)&&(!shallowCopy))
	{
		delete _data;
	}
	shallowCopy=false;
	_data=new double[_size];
	int i;
	for (i=0;i<_size;i++)
	{
		_data[i]=o._data[i];
	}
	return (*this);
}

marVector& marVector::operator=(double o)
{
	if (_size!=0)
	{
		if (!(_data==NULL)&&(!shallowCopy))
		{
			delete _data;
		}
		shallowCopy=false;
		_data=new double[_size];
		int i;
		for (i=0;i<_size;i++)
		{
			_data[i]=o;
		}
	}
	else
	{
		assert(false);
	}
	return (*this);
}

marVector& marVector::operator=(double*o)
{
	if (_size!=0)
	{
		if (!(_data==NULL)&&(!shallowCopy))
		{
			delete _data;
		}
		shallowCopy=false;
		_data=new double[_size];
		int i;
		for (i=0;i<_size;i++)
		{
			_data[i]=o[i];
		}
	}
	else
	{
		assert(false);
	}
	return (*this);
}

// ---------------------------------------------------------------------
/*Comparaison*/
bool marVector::operator==(const marVector& o) const
{
    bool equal=true;
	
    if(_size!=o._size)
	{
		return(false);
	}
	
	int i;
	for(i=0;i<_size && equal;i++)
	{
		equal=equal && (_data[i]==o._data[i]); 
	}
    return(equal);
}

bool marVector::operator!=(const marVector& o) const
{
	return(!((*this)==o));
}

// ---------------------------------------------------------------------
/*Addition*/
marVector marVector::operator+(const marVector& o)
{
	marVector result(*this);
	size_t s=result._size;
	
	if ((o._size!=s) || (s==0))
	{
		assert(false);
	}
	else 
	{
		for (int i=0;i < s;i++)
		{
			result._data[i] +=o._data[i];
		}
	}
	return result;
}

marVector marVector::operator+(double o)
{
	marVector result(*this);
	size_t s=result._size;
	
	if (s==0)
	{
		assert(false);
	}
	else 
	{
		for (int i=0;i < s;i++)
		{
			result._data[i] +=o;
		}
	}
	return result;
}

marVector marVector::operator+(double*o)
{
	marVector result(*this);
	size_t s=result._size;
	
	if (s==0)
	{
		assert(false);
	}
	else 
	{
		for (int i=0;i < s;i++)
		{
			result._data[i] +=o[i];
		}
	}
	return result;
}

// ---------------------------------------------------------------------
/*Addition + Affectation*/
marVector& marVector::operator+=(const marVector& o)
{
	if ((o._size!=_size) || (_size==0))
	{
		assert(false);  	
	}
	else 
	{
		for (int i=0;i < _size;i++)
		{
			_data[i] +=o._data[i];
		}
	}
	return (*this);
}

marVector& marVector::operator+=(double o)
{
	if (_size==0)
	{
		assert(false);  
	}
	else 
	{
		for (int i=0;i < _size;i++)
		{
			_data[i] +=o;
		}
	}
    return(*this);
}

marVector& marVector::operator+=(double*o)
{
	if (_size==0)
	{
		assert(false);  	
	}
	else 
	{
		for (int i=0;i < _size;i++)
		{
			_data[i] +=o[i];
		}
	}
    return(*this);
}

// ---------------------------------------------------------------------
/*Soustraction*/
marVector marVector::operator-(const marVector& o)
{
	marVector result(*this);
	size_t s=result._size;
	
	if ((o._size !=s) || (s==0))
	{
		assert(false);  	
	}
	else 
	{
		for (int i=0;i < s;i++)
		{
			result._data[i] -=o._data[i];
		}
	}
	return result;
}

marVector marVector::operator-(double o)
{
	marVector result(*this);
	size_t s=result._size;
	
	if (s==0)
	{
		assert(false);  
	}
	else 
	{
		for (int i=0;i < s;i++)
		{
			result._data[i] -=o;
		}
	}
	return result;
}

marVector marVector::operator-(double*o)
{
	marVector result(*this);
	size_t s=result._size;
	
	if (s==0)
	{
		assert(false);  	
	}
	else 
	{
		for (int i=0;i < s;i++)
		{
			result._data[i] -=o[i];
		}
	}
	return result;
}

// ---------------------------------------------------------------------
/*Sooustraction + Affection*/
marVector& marVector::operator-=(const marVector& o)
{
	if ((o._size!=_size) || (_size==0))
	{
		assert(false);  	
	}
	else 
	{
		for (int i=0;i < _size;i++)
		{
			_data[i] -=o._data[i];
		}
	}
	return (*this);
}

marVector& marVector::operator-=(double o)
{
	if (_size==0)
	{
		assert(false);  
	}
	else 
	{
		for (int i=0;i < _size;i++)
		{
			_data[i] -=o;
		}
	}
    return(*this);
}

marVector& marVector::operator-=(double*o)
{
	if (_size==0)
	{
		assert(false);  	
	}
	else 
	{
		for (int i=0;i < _size;i++)
		{
			_data[i] -=o[i];
		}
	}
    return(*this);
}

// ---------------------------------------------------------------------
/*Multiplication (produit scalaire)*/
marVector marVector::operator*(double o)
{
	marVector result(*this);
	size_t s=result._size;
	
	if (s==0)
	{
		assert(false);  
	}
	else 
	{
		for (int i=0;i < s;i++)
		{
			result._data[i]*=o;
		}
	}
	return result;
}

// ---------------------------------------------------------------------
/*Multiplication (produit scalaire) + Affectation*/
marVector& marVector::operator*=(double o)
{
	if (_size==0)
	{
		assert(false);  
	}
	else 
	{
		for (int i=0;i < _size;i++)
		{
			_data[i]*=o;
		}
	}
    return(*this);
}

// ---------------------------------------------------------------------
/*Division (division scalaire)*/
marVector marVector::operator/(double o)
{
	marVector result(*this);
	size_t s=result._size;
	
	if (s==0)
	{
		assert(false);  
	}
	else 
	{
		for (int i=0;i < s;i++)
		{
			result._data[i] /=o;
		}
	}
	return result;
}

// ---------------------------------------------------------------------
/*Division (division scalaire) + Affectation*/
marVector& marVector::operator/=(double o)
{
	if (_size==0)
	{
		assert(false);  
	}
	else 
	{
		for (int i=0;i < _size;i++)
		{
			_data[i] /=o;
		}
	}
    return(*this);
}

//////////////////////////////////////////////////////////////////////
// Opérations
//////////////////////////////////////////////////////////////////////
// ---------------------------------------------------------------------
/*Produit scalaire*/
double marVector::dot(const marVector& o)
{
    double result;
	if ((_size!=o._size) || (_size==0))
	{
		assert(false); 
	}
	else
	{
		result=0.0;
		for (int i=0;i<_size;i++)
		{
			result+=(_data[i]*o._data[i]);
		}
	}
	return result;
}

double marVector::dot(double*o)
{
    double result;
	if (_size==0)
	{
		assert(false); 
	}
	else
	{
		result=0.0;
		for (int i=0;i<_size;i++)
		{
			result+=(_data[i]*o[i]);
		}
	}
	return result;
}

// ---------------------------------------------------------------------
/*Produit vectoriel*/
marVector marVector::cross(const marVector& o)
{
    marVector result(*this);
	if ((o._size!=_size) || (_size!=3))
	{
		assert(false); 
	}
	else
	{
		double S,s;
		
		// 1st element
		S=_data[1]*o._data[2];
		s=_data[2]*o._data[1];
		result._data[0]=S - s;
		
		// 2nd element
		S=_data[2]*o._data[0];
		s=_data[0]*o._data[2];
		result._data[1]=S - s;
		
		// 3rd element
		S=_data[0]*o._data[1];
		s=_data[1]*o._data[0];
		result._data[2]=S - s;
	}
    return(result);    
}

marVector marVector::cross(double*o)
{
    marVector result(*this);
	if (_size!=3)
	{
		assert(false); 
	}
	else
	{
		double S,s;
		
		// 1st element
		S=_data[1]*o[2];
		s=_data[2]*o[1];
		result._data[0]=S - s;
		
		// 2nd element
		S=_data[2]*o[0];
		s=_data[0]*o[2];
		result._data[1]=S - s;
		
		// 3rd element
		S=_data[0]*o[1];
		s=_data[1]*o[0];
		result._data[2]=S - s;
	}
    return(result);
}

int marVector::scross(const marVector& o)
{
	if ((o._size!=_size) && (_size!=3))
	{
		return (1);
	}
	else
	{
		*this=cross(o);
		return(0);
	} 
}

int marVector::scross(double*o)
{
	if (_size!=3)
	{
		return (1);
	}
	else
	{
		*this=cross(o);
		return(0);
	} 
}

// ---------------------------------------------------------------------
/*Norme euclidienne*/
double marVector::norm2()
// Calcule la norme euclidienne du vecteur 
// en utilisant l'agorithme de Blue
{
	
	double norme2=0.0;
	if (_size==0)
	{
		assert(false);
	}
	else
	{
		const double seuil1=pow( 2., 154. ),seuil2=1/seuil1;
		double n1=0; /* les sommes partielles */
		double n2=0;
		double n3=0;
		double n4;
		int i;
		for ( i=0; i<_size; i++ ) 
		{
			double x=_data[i];	/* x=abs(vecteur[i]) */
			
			if (x<0) 
			{ 
				x=-x;
			}
			/* sommation par classe */
			if ( x > seuil1 ) 
			{
				x *=seuil2;
				n1 +=x*x;
			} 
			else if ( x < seuil2 ) 
			{
				x  *=seuil1;
				n2 +=x*x;
			} 
			else
			{
				n3 +=x*x;
			}
		}
		
		n3=sqrt(n3);
		if (n1>0)
		{
			n4=seuil1 * sqrt(n1);
		}
		else
		{
			n4=seuil2 * sqrt(n2);
		}
		
// EED 
/*
		if (n3 < n4) 
		{
			n3 /=n4;
			norme2=n4 * sqrt( 1 + n3*n3 );
		} 
		else 
		{
			n4 /=n3;
			norme2=n3 * sqrt( 1 + n4*n4 );
		}
*/
		norme2=n3;
	}
	return norme2;
}

// ---------------------------------------------------------------------
/*Normalisation*/
marVector marVector::normalize( )
{
	marVector result=marVector(*this);
	result *= (1/result.norm2());
	return result;
}

int marVector::snormalize( )
{
	(*this)*=(1/this->norm2());
	return(0);
}// renvoie 0 si OK, 1 sinon

//////////////////////////////////////////////////////////////////////
// Méthodes
//////////////////////////////////////////////////////////////////////

size_t marVector::size() const
{
	return _size;
}