#ifndef STIM_AABB
#define STIM_AABB

namespace stim{


/// This class describes a structure for an axis-aligned bounding box
template< typename T >
class aaboundingbox{

public:
	bool set;				//has the bounding box been set to include any points?
	stim::vec3<T> A;			//minimum point in the bounding box
	stim::vec3<T> B;			//maximum point in the bounding box

	aaboundingbox(){					//constructor generates an empty bounding box
		set = false;
	}


	/// Test if a point is inside of the bounding box and returns true if it is.

	/// @param p is the point to be tested
	bool test(stim::vec3<T> p){

		for(unsigned d = 0; d < p.size(); p++){		//for each dimension
			if(p[d] < A[d]) return false;			//if the point is less than the minimum bound, return false
			if(p[d] > B[d]) return false;			//if the point is greater than the max bound, return false
		}
		return true;
	}

	/// Expand the bounding box to include the specified point.

	/// @param p is the point to be included
	void expand(stim::vec3<T> p){

		if(!set){							//if the bounding box is empty, fill it with the current point
			A = B = p;
			set = true;
		}

		for(unsigned d = 0; d < p.size(); d++){		//for each dimension
			if(p[d] < A[d]) A[d] = p[d];			//expand the bounding box as necessary
			if(p[d] > B[d]) B[d] = p[d];
		}
	}

	/// Return the center point of the bounding box as a stim::vec
	stim::vec3<T> center(){
		return (B + A) * 0.5;
	}

	/// Return the size of the bounding box as a stim::vec
	stim::vec3<T> size(){
		return (B - A);
	}

	/// Generate a string for the bounding box
	std::string str(){
		std::stringstream ss;
		ss<<A.str()<<"----->"<<B.str();
		return ss.str();
	}

        ///resamples the boundingbox by a factor of N
        std::vector<stim::aaboundingbox<T> >
        resample(unsigned int n)
        {   
                stim::vec3<T> step_size = (B-A)/n;
                std::vector<stim::aaboundingbox<T> > boundingboxes(n*n*n);
                for(int x = 0; x < n; x++)
                {   
                        for(int y = 0; y < n; y++)
                        {   
                                for(int z = 0; z < n; z++)
                                {
                                        ///create the bounding box
                                        stim::aaboundingbox<T> box;
                                        box.A[0] = A[0]+step_size[0]*x;
                                        box.A[1] = A[1]+step_size[1]*y;
                                        box.A[2] = A[2]+step_size[2]*z;

                                        box.B[0] = A[0]+step_size[0]*(x+1);
                                        box.B[1] = A[1]+step_size[1]*(y+1);
                                        box.B[2] = A[2]+step_size[2]*(z+1);
                                        ///put the new bounding box into the right place
                                        boundingboxes[x*n*n + y*n + z] = box;
                                }
                        }
                }
                return boundingboxes;
        }


};		//end stim::aabb


};		//end namespace stim

#endif