/* ImLib3D
 * Copyright (c) 2001, ULP-IPB Strasbourg.
 *
 * This program is free software; you can redistribute it and/or modify
 * it under the terms of the GNU General Public License as published by
 * the Free Software Foundation; either version 2 of the License, or (at
 * your option) any later version.
 * 
 * This program is distributed in the hope that it will be useful, but
 * WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
 * General Public License for more details.
 * 
 * You should have received a copy of the GNU General Public License
 * along with this program; if not, write to the Free Software
 * Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
 */
/*!  \file   Container3DIterators.hpp
*/

#ifndef _Container3D_internal
#error Image3DIterators should only be included by Container3D.hpp
#endif

// **************************************************
// ****************** iterators *********************
// **************************************************

// FIXME!! this is a temporary fix for gcc>=3.1
#if __GNUC__==3 && __GNUC_MINOR__>=1

template<typename _Tp, typename _Distance>
  struct random_access_iterator {
    typedef random_access_iterator_tag iterator_category;
    typedef _Tp                        value_type;
    typedef _Distance                  difference_type;
    typedef _Tp*                       pointer;
    typedef _Tp&                       reference;
  };

#endif


//! basic iterator for walking through the WHOLE container
template<class Container3D,class _Ref,class _Ptr>
class CTR3D_iteratorFast : public random_access_iterator<typename Container3D::value_type,size_t>
{
  	typedef CTR3D_iteratorFast<Container3D,_Ref,_Ptr> _Self;
public:
	Container3D *ima;
	_Ptr pos;
public:
	void AssertBounds()
	{
		if(pos<ima->begin().pos || pos>=ima->end().pos)
		{
			ThrowError("iteratorFast: out of bounds:pos:%x begin:%x end:%x\n",pos,ima->begin().pos,ima->end().pos);
		}
	}
	inline _Ref operator*(){return *pos;}
	double Progress() const {if(!this){return 0;}return (pos-ima->GetData())/(float)ima->GetNVoxels();}
	inline _Self  operator++(int)     {_Self __tmp=*this;pos++;return __tmp;}
	inline _Self &operator++()        {pos++;return *this;}
	inline void   operator+=(size_t m){pos+=m;}
	inline _Self  operator+ (size_t m){_Self __tmp=*this;__tmp.pos+=m;return __tmp;}
	inline _Self  operator- (size_t m){_Self __tmp=*this;__tmp.pos-=m;return __tmp;}
	inline _Self  operator--(int)     {_Self __tmp=*this;pos--;return __tmp;}
	inline _Self &operator--()        {pos--;return *this;}
	inline void   operator-=(size_t m){pos-=m;}
	inline int    operator- (const _Self &other){return pos-other.pos;}
	inline bool   operator< (const _Self &other){return pos< other.pos;}
	inline bool   operator<=(const _Self &other){return pos<=other.pos;}
	inline bool   operator> (const _Self &other){return pos> other.pos;}
	inline bool   operator>=(const _Self &other){return pos>=other.pos;}

	inline Vect3Di ComputePosition() const
	{
		int offset=pos-ima->GetData();
		Vect3Di res;
		res.z=offset/(ima->Height()*ima->Width());
		offset=offset%(ima->Height()*ima->Width());
		res.y=offset/ima->Width();
		res.x=offset%ima->Width();
		return res;
	}
	void operator+=(const Vect3Di &dV)
	{
		pos+=dV.x+ima->Width()*dV.y+ima->Width()*ima->Height()*dV.z;
	}
	inline void operator=(const Vect3Di &P)
	{
		pos=ima->GetData()+P.x+ima->Width()*P.y+ima->Width()*ima->Height()*P.z;
	}
	inline bool  operator!=(const _Self  &other){return pos!=other.pos; }
	inline bool  operator==(const _Self  &other){return pos==other.pos; }
	CTR3D_iteratorFast():ima(NULL),pos(NULL){}
	CTR3D_iteratorFast(Container3D *_ima,_Ptr _pos):ima(_ima),pos(_pos){}
};

//! iterator that keeps track of x,y,z position
/*! for walking through the WHOLE container */
template<class Container3D,class _Ref,class _Ptr>
class CTR3D_iteratorXYZ :    public CTR3D_iteratorFast<Container3D,_Ref,_Ptr>
{
	typedef CTR3D_iteratorFast<Container3D,_Ref,_Ptr> _Base;
  	typedef CTR3D_iteratorXYZ <Container3D,_Ref,_Ptr> _Self;
	void Reposition(){Vect3Di v=ComputePosition();x=v.x;y=v.y;z=v.z;}
public:
	int x,y,z;

	void AssertBounds()
	{		
		if(!ima->IsInside(x,y,z) || pos<ima->GetData() || pos>=ima->GetData()+ima->GetNVoxels())
		{
			ThrowError("iteratorXYZ: out of bounds %d %d %d: image %d %d %d :pos:%x begin:%x end:%x\n",
					 x,y,z,
					 ima->Width(),ima->Height(),ima->Depth(),
					 pos,ima->GetData(),ima->GetData()+ima->GetNVoxels());
		}
	}

	inline _Self  operator++(int){_Self __tmp=*this;++(*this);return __tmp;}
	inline _Self &operator++()
	{
		pos++;
		x++;
		if(x>=ima->Width())
		{
			x=0;
			y++;
			if(y>=ima->Height())
			{
				y=0;
				z++;
			}
		}
		return *this;
	}
	inline _Self  operator--(int){_Self __tmp=*this;--(*this);return __tmp;}
	inline _Self &operator--()
	{
		pos--;
		x--;
		if(x<0)
		{
			x=ima->Width()-1;
			y--;
			if(y<0)
			{
				y=ima->Height()-1;
				z--;
			}
		}
		return *this;
	}
	//! deprecated!!!
	inline Vect3Di V() const {return Vect3Di(x,y,z);} 
	//! use this instead
	inline Vect3Di Pos() const {return Vect3Di(x,y,z);}
	// jump to new location
	inline void operator=(const Vect3Di &P)
	{
		x=P.x;
		y=P.y;
		z=P.z;
		pos=&(*ima)(P);
	}

	inline void   operator+=(size_t m){pos+=m;Reposition();}
	inline void   operator-=(size_t m){pos-=m;Reposition();}
	inline _Self  operator+ (size_t m){_Self __tmp=*this;__tmp.pos+=m;__tmp.Reposition();return __tmp;}
	inline _Self  operator- (size_t m){_Self __tmp=*this;__tmp.pos-=m;__tmp.Reposition();return __tmp;}

	void operator= (const CTR3D_iteratorFast<Container3D,_Ref,_Ptr> &other)
	{
		ima=other.ima;
		pos=other.pos;		
		if(pos==ima->GetData()){x=0;y=0;z=0;}
		else{Reposition();}
	}
	CTR3D_iteratorXYZ():_Base(){;}
	CTR3D_iteratorXYZ(const CTR3D_iteratorFast<Container3D,_Ref,_Ptr> &other) : _Base(other)
	{
		operator=(other);
	}
	CTR3D_iteratorXYZ(Container3D *_ima,_Ptr _pos) : _Base(_ima,_pos)
	{
		// optionally move xyz if pos is not at begining of ima
		if(pos!=ima->GetData())
		{_Base::operator=(CTR3D_iteratorFast<Container3D,_Ref,_Ptr>(ima,pos).ComputePosition());}
	}

};

//! iterator for walking through a rectangular portion
/*! keeps track of x,y,z position */
template<class Container3D,class _Ref,class _Ptr>
class CTR3D_iteratorZone : public CTR3D_iteratorXYZ<Container3D,_Ref,_Ptr>
{
	typedef CTR3D_iteratorXYZ <Container3D,_Ref,_Ptr> _Base;
  	typedef CTR3D_iteratorZone<Container3D,_Ref,_Ptr> _Self;
public:
  	RectZone3Di zone;
	double Progress(){return (pos-ima->GetData())/(float)zone.Volume();}
	void AssertBounds()
	{		
		if(!ima->IsInside(x,y,z)  || !zone.IsInside(x,y,z) || 
		   pos<ima->GetData() || pos>=ima->GetData()+ima->GetNVoxels())
		{
			ThrowError("iteratorZone: out of bounds %d %d %d:zone:%d %d %d - %d %d %d: image %d %d %d : pos:%x begin:%x end:%x\n",
					 x,y,z,
					 zone.x0,zone.y0,zone.z0,
					 zone.x1,zone.y1,zone.z1,
					 ima->Width(),ima->Height(),ima->Depth(),
					 pos,ima->GetData(),ima->GetData()+ima->GetNVoxels());
		}
	}
	inline _Self  operator++(int){_Self __tmp=*this;++(*this);return __tmp;}
	inline _Self &operator++()
	{
  		pos++;
  		if(zone.ToNext(x,y,z))
  		{
  			if(zone.IsInside(x,y,z))
  			{
  				pos=&(*ima)(x,y,z);
  			}
  			else{pos=ima->GetDataEnd();}
  		}
		return *this;
	}
	inline void operator--(int)
	{
		pos--;
  		if(zone.ToPrev(x,y,z))
  		{
  			if(zone.IsInside(x,y,z))
  			{
  				pos=&(*ima)(x,y,z);
  			}
  			else{pos=NULL;}
  		}
	}

	inline void  operator+=(size_t m){while (m--) operator++(0);}
	inline void  operator-=(size_t m){while (m--) operator--(0);}
	inline _Self operator+ (size_t m){_Self __tmp=*this;__tmp+=m;return __tmp;}
	inline _Self operator- (size_t m){_Self __tmp=*this;__tmp-=m;return __tmp;}

	void operator= (const CTR3D_iteratorFast<Container3D,_Ref,_Ptr> &other)
	{
		ima=other.ima;		
		pos=other.pos;		
		if(!ima->IsInside(zone.GetP0()) || !ima->IsInside(zone.GetP1()))
		{
			ThrowError("iteratorZone zone does not fit in image");
		}		
		// if iter fast is begining of image put it at begining of zone
		if(pos==ima->GetData()){zone.ToStart(x,y,z);}
		else{_Base::operator=(other.ComputePosition());}
		pos=&(*ima)(x,y,z);
	}
	CTR3D_iteratorZone():_Base(){}
//  	CTR3D_iteratorZone(const CTR3D_iteratorFast<Container3D,_Ref,_Ptr> &other) : _Base()
//  	{
//  		operator=(other);
//  	}
  	CTR3D_iteratorZone(const RectZone3Di &_zone):
  		_Base(),zone(_zone)
  	{
	}
};



//! iterator for walking through the WHOLE container, but only on unmasked values
template<class Container3D,class _Ref,class _Ptr>
class CTR3D_iteratorFastMasked : public CTR3D_iteratorFast<Container3D,_Ref,_Ptr>
{
	typedef CTR3D_iteratorFast      <Container3D,_Ref,_Ptr> _Base;
  	typedef CTR3D_iteratorFastMasked<Container3D,_Ref,_Ptr> _Self;
	const byte *pmask;
public:
	void AssertMasked()
	{
		if(*pmask==0)
		{
			ThrowError("iteratorFastMasked: mask at pos==0 :pos:%x begin:%x end:%x\n",_Base::pos,_Base::begin,_Base::end);
		}
	}
	inline _Self  operator++(int){_Self __tmp=*this;++(*this);return __tmp;}
	inline _Self &operator++()
	{
		// for lazy programers that don't want to write masked and unmasked versions
		// of their algorithms... 
		// this works on unmasked images too but warning, this is a litle bit slower
		if(!pmask){_Base::operator++(0);return *this;}
		do
		{
			pmask++;
			_Base::operator++(0);
		}while(!(*pmask) && _Base::pos!=ima->GetDataEnd());
		return *this;
	}	
	inline _Self  operator--(int){_Self __tmp=*this;--(*this);return __tmp;}
	inline _Self &operator--()
	{
		// waring: this is dangerous beacause there is no 100% safe way of checking if we are under begining
		// for lazy programers that don't want to write masked and unmasked versions
		// of their algorithms... 
		// this works on unmasked images too but warning, this is a litle bit slower
		if(!pmask){_Base::operator--(0);return *this;}
		do
		{
			pmask--;
			_Base::operator--(0);
		}while(!(*pmask) && _Base::pos!=ima->GetData());
		return *this;
	}	


	inline void  operator+=(size_t m){while (m--) operator++(0);}
	inline void  operator-=(size_t m){while (m--) operator--(0);}
	inline _Self operator+ (size_t m){_Self __tmp=*this;__tmp+=m;return __tmp;}
	inline _Self operator- (size_t m){_Self __tmp=*this;__tmp-=m;return __tmp;}

	void operator=(const CTR3D_iteratorFast<Container3D,_Ref,_Ptr> &other)
	{ 
		_Base::operator=(other);
		if(ima->HasMask()){pmask=ima->Mask().GetData();}
		else
		{pmask=NULL;}
	}

	CTR3D_iteratorFastMasked():_Base(),pmask(NULL){}
	CTR3D_iteratorFastMasked(const CTR3D_iteratorFast<Container3D,_Ref,_Ptr> &other):
		_Base()
	{
		operator=(other);
	}
};


//! iterator that keeps track of x,y,z position
/*! for walking through the WHOLE container */
template<class Container3D,class _Ref,class _Ptr>
class CTR3D_iteratorXYZMasked :    public CTR3D_iteratorXYZ<Container3D,_Ref,_Ptr>
{
	typedef CTR3D_iteratorXYZ      <Container3D,_Ref,_Ptr> _Base;
  	typedef CTR3D_iteratorXYZMasked<Container3D,_Ref,_Ptr> _Self;
	const byte *pmask;
public:
	inline _Self  operator++(int){_Self __tmp=*this;++(*this);return __tmp;}
	inline _Self &operator++()
	{
		// for lazy programers that don't want to write masked and unmasked versions
		// of their algorithms... 
		// this works on unmasked images too but warning, this is a litle bit slower
		if(!pmask){_Base::operator++(0);return *this;}
		do
		{
			pmask++;
			_Base::operator++(0);
		}while(!(*pmask) && _Base::pos!=ima->GetDataEnd());
		return *this;
	}	
	inline _Self  operator--(int){_Self __tmp=*this;--(*this);return __tmp;}
	inline _Self &operator--()
	{
		// waring: this is dangerous beacause there is no 100% safe way of checking if we are under begining
		// for lazy programers that don't want to write masked and unmasked versions
		// of their algorithms... 
		// this works on unmasked images too but warning, this is a litle bit slower
		if(!pmask){_Base::operator--(0);return *this;}
		do
		{
			pmask--;
			_Base::operator--(0);
		}while(!(*pmask) && _Base::pos!=ima->GetData());
		return *this;
	}	
	// jump to new location
	inline void operator=(const Vect3Di &P)
	{
		ThrowError("iteratorXYZMasked operator=V unimplemented");
	}

	inline void  operator+=(size_t m){while (m--) operator++(0);}
	inline void  operator-=(size_t m){while (m--) operator--(0);}
	inline _Self operator+ (size_t m){_Self __tmp=*this;__tmp+=m;return __tmp;}
	inline _Self operator- (size_t m){_Self __tmp=*this;__tmp-=m;return __tmp;}

	void operator= (const CTR3D_iteratorFast<Container3D,_Ref,_Ptr> &other)
	{
		_Base::operator=(other);
		if(ima->HasMask()){pmask=ima->Mask().GetData();}
		else
		{pmask=NULL;}
	}
	CTR3D_iteratorXYZMasked():_Base(),pmask(NULL){;}
	CTR3D_iteratorXYZMasked(const CTR3D_iteratorFast<Container3D,_Ref,_Ptr> &other):
		_Base()
	{
		operator=(other);
	}
};

//! iterator for walking through image by concentric boxes centered on start point
/*! keeps track of x,y,z position */
template<class Container3D,class _Ref,class _Ptr>
class CTR3D_iteratorBoxConcentric : public CTR3D_iteratorFast<Container3D,_Ref,_Ptr>
{
	typedef CTR3D_iteratorFast         <Container3D,_Ref,_Ptr> _Base;
  	typedef CTR3D_iteratorBoxConcentric<Container3D,_Ref,_Ptr> _Self;
public:
	// intersection of current box and image box
  	RectZone3Di currentBox;
	RectZone3Di imazone;
	int x,y,z;
	void AssertBounds()
	{		
//  		if(!ima->IsInside(x,y,z)  || !zone.IsInside(x,y,z) || 
//  		   pos<ima->GetData() || pos>=ima->GetData()+ima->GetNVoxels())
//  		{
//  			ThrowError("iteratorZone: out of bounds %d %d %d:zone:%d %d %d - %d %d %d: image %d %d %d : pos:%x begin:%x end:%x\n",
//  					 x,y,z,
//  					 zone.x0,zone.y0,zone.z0,
//  					 zone.x1,zone.y1,zone.z1,
//  					 ima->Width(),ima->Height(),ima->Depth(),
//  					 pos,ima->GetData(),ima->GetData()+ima->GetNVoxels());
//  		}
	}
	inline _Self  operator++(int){_Self __tmp=*this;++(*this);return __tmp;}
	inline _Self &operator++()
	{
#ifdef SAFE_DEBUG
		AssertBounds();
#endif
		// go to next point on current box that's inside image
		do{currentBox.ToNextOnSurface(x,y,z);}
		while(!imazone.IsInside(x,y,z) && currentBox.IsInside(x,y,z));

		if(currentBox.IsInside(x,y,z))
		{
			pos=&(*ima)(x,y,z);
		}
		else
		{
			// finished this box, now find next concentric box
			currentBox=RectZone3Di(currentBox.GetP0()-Vect3Di(1,1,1),currentBox.GetP1()+Vect3Di(1,1,1));
			// check if this was last box
			if(currentBox.IsInside(imazone.GetP0()) && !currentBox.IsOnSurface(imazone.GetP0()) &&
			   currentBox.IsInside(imazone.GetP1()) && !currentBox.IsOnSurface(imazone.GetP1())    )
			{pos=ima->GetDataEnd();return *this;}
			// place current pos at first point in box that's ok
			currentBox.ToStart(x,y,z);
			while(!imazone.IsInside(x,y,z) ){currentBox.ToNextOnSurface(x,y,z);}
			pos=&(*ima)(x,y,z);
		}
		return *this;
	}
	inline Vect3Di V() const {return Vect3Di(x,y,z);}
	inline Vect3Di Pos() const {return Vect3Di(x,y,z);}
	inline _Ref operator*()
	{
#ifdef SAFE_DEBUG
		AssertBounds();
#endif
		return *pos;
	}
	void operator= (const CTR3D_iteratorFast<Container3D,_Ref,_Ptr> &other)
	{
		_Base::operator=(other);
		if(pos!=ima->GetData()){ThrowError("BoxConcentric::operator= jump to rand pos unimplemented");}
		imazone=ima->GetZone();
		currentBox.ToStart(x,y,z);
		if(!imazone.IsInside(x,y,z))
		{
			Vect3Di p=Vect3Di(x,y,z);
			Vect3Di cp=imazone.ClosestPoint(p);
			Vect3Di dp=p-cp;
			int mx=max(abs(dp.x),max(abs(dp.y),abs(dp.z)));
			dp=Vect3Di(mx,mx,mx);
			currentBox=RectZone3Di(p-dp,p+dp);
			currentBox.ToStart(x,y,z);
			while(!imazone.IsInside(x,y,z) ){currentBox.ToNextOnSurface(x,y,z);}
		}
		pos=&(*ima)(x,y,z);
	}
	CTR3D_iteratorBoxConcentric():_Base(){}
  	CTR3D_iteratorBoxConcentric(const Vect3Di &p0):
  		_Base(),currentBox(p0,p0),imazone(0,0,0,0,0,0)
  	{
	}
};

