diff --git a/stim/math/plane.h b/stim/math/plane.h
index 26f28cb..1b41f45 100644
--- a/stim/math/plane.h
+++ b/stim/math/plane.h
@@ -31,6 +31,12 @@ template <typename T> class plane
 			init();
 		}
 
+		CUDA_CALLABLE plane(vec<T> p)
+		{
+			init();
+			P = p;
+		}
+
 		CUDA_CALLABLE plane(vec<T> n, vec<T> p = vec<T>(0, 0, 0))
 		{
 			init();
@@ -125,6 +131,11 @@ template <typename T> class plane
 			return v - perpendicular(v);
 		}
 
+		CUDA_CALLABLE void setU(vec<T> v)
+		{
+			U = (parallel(v.norm())).norm();		
+		}
+
 		CUDA_CALLABLE void decompose(vec<T> v, vec<T>& para, vec<T>& perp){
 			perp = N * v.dot(N);
 			para = v - perp;
@@ -176,6 +187,18 @@ template <typename T> class plane
 
 		}
 
+		CUDA_CALLABLE void rotate(vec<T> n, vec<T> &X, vec<T> &Y)
+		{
+			quaternion<T> q;
+			q.CreateRotation(N, n);
+			
+			N = q.toMatrix3() * N;
+			U = q.toMatrix3() * U;
+			X = q.toMatrix3() * X;
+			Y = q.toMatrix3() * Y;
+
+		}
+
 };
 		
 		
diff --git a/stim/math/rect.h b/stim/math/rect.h
index 336a3ac..e8d5bec 100644
--- a/stim/math/rect.h
+++ b/stim/math/rect.h
@@ -1,11 +1,11 @@
-#ifndef RTS_RECT_H
-#define RTS_RECT_H
+#ifndef STIM_RECT_H
+#define STIM_RECT_H
 
 //enable CUDA_CALLABLE macro
 #include <stim/cuda/cudatools/callable.h>
+#include <stim/math/plane.h>
 #include <stim/math/vector.h>
 #include <stim/math/triangle.h>
-#include <stim/math/quaternion.h>
 #include <iostream>
 #include <iomanip>
 #include <algorithm>
@@ -13,14 +13,14 @@
 namespace stim{
 
 //template for a rectangle class in ND space
-template <class T>
-struct rect
+template <typename T>
+class rect : public plane <T>
 {
 	/*
 		^                   O
 		|                   
 		|                   
-		Y         C         
+		Y         P         
 		|                   
 		|                   
 		O---------X--------->
@@ -28,7 +28,6 @@ struct rect
 
 private:
 
-	stim::vec<T> C;
 	stim::vec<T> X;
 	stim::vec<T> Y;
 
@@ -42,26 +41,30 @@ private:
 public:
 
 	///base constructor.
-	CUDA_CALLABLE rect(){
+	CUDA_CALLABLE rect()
+	 : plane<T>()
+	{
 		init();
 	}
 
 	///create a rectangle given a size and position in Z space.
 	///@param size: size of the rectangle in ND space.
 	///@param z_pos z coordinate of the rectangle.
-	CUDA_CALLABLE rect(T size, T z_pos = (T)0){
+	CUDA_CALLABLE rect(T size, T z_pos = (T)0)
+	 : plane<T>(z_pos)
+	{
 		init();			//use the default setup
 		scale(size);	//scale the rectangle
-		C[2] = z_pos;
 	}
 
 	
 	///create a rectangle from a center point, normal
 	///@param c: x,y,z location of the center.
 	///@param n: x,y,z direction of the normal.
-	CUDA_CALLABLE rect(vec<T> c, vec<T> n = vec<T>(0, 0, 1)){
+	CUDA_CALLABLE rect(vec<T> c, vec<T> n = vec<T>(0, 0, 1))
+		: plane<T>(n, c)
+	{
 		init();			//start with the default setting
-		C = c;
 		normal(n);		//orient
 	}
 
@@ -69,11 +72,11 @@ public:
 	///@param c: x,y,z location of the center.
 	///@param s: size of the rectangle.
 	///@param n: x,y,z direction of the normal.
-	CUDA_CALLABLE rect(vec<T> c, T s, vec<T> n = vec<T>(0, 0, 1)){
+	CUDA_CALLABLE rect(vec<T> c, T s, vec<T> n = vec<T>(0, 0, 1))
+		: plane<T>(n, c)
+	{
 		init();			//start with the default setting
-		C = c;
 		scale(s);
-		normal(n);		//orient
 	}
 
 	///creates a rectangle from a centerpoint and an X and Y direction vectors.
@@ -81,8 +84,8 @@ public:
 	///@param directionX: u,v,w direction of the X vector.
 	///@param directionY: u,v,w direction of the Y vector.
 	CUDA_CALLABLE rect(vec<T> center, vec<T> directionX, vec<T> directionY )
+		 : plane<T>(p)
 	{
-		C = center;
 		X = directionX;
 		Y = directionY;
 	}
@@ -93,8 +96,8 @@ public:
 	///@param directionX: u,v,w direction of the X vector.
 	///@param directionY: u,v,w direction of the Y vector.
 	CUDA_CALLABLE rect(T size, vec<T> center, vec<T> directionX, vec<T> directionY )
+		: plane<T>(p)
 	{	
-		C = center;
 		X = directionX;
 		Y = directionY;
 		scale(size);
@@ -106,8 +109,8 @@ public:
 	///@param directionX: u,v,w direction of the X vector.
 	///@param directionY: u,v,w direction of the Y vector.
 	CUDA_CALLABLE rect(vec<T> size, vec<T> center, vec<T> directionX, vec<T> directionY )
+		: plane<T>(p)
 	{	
-		C = center;
 		X = directionX;
 		Y = directionY;
 		scale(size[0], size[1]);
@@ -116,28 +119,24 @@ public:
 	///scales a rectangle in ND space.
 	///@param factor1: size of the scale in the X-direction.
 	///@param factor2: size of the scale in the Y-direction.	
-	CUDA_CALLABLE void scale(T factor1, T factor2){
+	CUDA_CALLABLE void scale(T factor1, T factor2)
+	{
 		X *= factor1;
 		Y *= factor2;
 	}
 
 	///@param n; vector with the normal.
 	///Orients the rectangle along the normal n.
-	CUDA_CALLABLE void normal(vec<T> n){		//orient the rectangle along the specified normal
-
-		n = n.norm();								//normalize, just in case
-		vec<T> n_current = X.cross(Y).norm();	//compute the current normal
-		quaternion<T> q;							//create a quaternion
-		q.CreateRotation(n_current, n);				//initialize a rotation from n_current to n
-
-		//apply the quaternion to the vectors and position
-		X = q.toMatrix3() * X;
-		Y = q.toMatrix3() * Y;
+	CUDA_CALLABLE void normal(vec<T> n)
+	{	
+		//orient the rectangle along the specified normal
+		rotate(n, X, Y);
 	}
 
 	///general init method that sets a general rectangle.
-	CUDA_CALLABLE void init(){
-		C = vec<T>(0, 0, 0);
+	CUDA_CALLABLE void init()
+	{
+		P = vec<T>(0, 0, 0);
 		X = vec<T>(1, 0, 0);
 		Y = vec<T>(0, 1, 0);
 	}
@@ -145,26 +144,19 @@ public:
 	//boolean comparison
 	bool operator==(const rect<T> & rhs)
 	{
-		if(C == rhs.C && X == rhs.X && Y == rhs.Y)
+		if(P == rhs.P && X == rhs.X && Y == rhs.Y)
 			return true;
 		else
 			return false;
 	}
 
-	/*******************************************
-	Return the normal for the rect
-	*******************************************/
-	CUDA_CALLABLE stim::vec<T> n()
-	{
-        return (X.cross(Y)).norm();
-	}
 
 	//get the world space value given the planar coordinates a, b in [0, 1]
 	CUDA_CALLABLE stim::vec<T> p(T a, T b)
 	{
 		stim::vec<T> result;
 		//given the two parameters a, b = [0 1], returns the position in world space
-		vec<T> A = C - X * (T)0.5 - Y * (T)0.5;
+		vec<T> A = P - X * (T)0.5 - Y * (T)0.5;
 		result = A + X * a + Y * b;
 
 		return result;
@@ -179,12 +171,12 @@ public:
 	std::string str()
 	{
 		std::stringstream ss;
-		vec<T> A = C - X * (T)0.5 - Y * (T)0.5;
+		vec<T> A = P - X * (T)0.5 - Y * (T)0.5;
 		ss<<std::left<<"B="<<std::setfill('-')<<std::setw(20)<<A + Y<<">"<<"C="<<A + Y + X<<std::endl;
 		ss<<std::setfill(' ')<<std::setw(23)<<"|"<<"|"<<std::endl<<std::setw(23)<<"|"<<"|"<<std::endl;
 		ss<<std::left<<"A="<<std::setfill('-')<<std::setw(20)<<A<<">"<<"D="<<A + X;
 
-        return ss.str();
+        	return ss.str();
 
 	}
 
@@ -207,20 +199,20 @@ public:
 	{
         //compute the distance between a point and this rect
 
-		vec<T> A = C - X * (T)0.5 - Y * (T)0.5;
+		vec<T> A = P - X * (T)0.5 - Y * (T)0.5;
 
-        //first break the rect up into two triangles
-        triangle<T> T0(A, A+X, A+Y);
-        triangle<T> T1(A+X+Y, A+X, A+Y);
+		//first break the rect up into two triangles
+		triangle<T> T0(A, A+X, A+Y);
+		triangle<T> T1(A+X+Y, A+X, A+Y);
 
 
-        T d0 = T0.dist(p);
-        T d1 = T1.dist(p);
+		T d0 = T0.dist(p);
+		T d1 = T1.dist(p);
 
-        if(d0 < d1)
-            return d0;
-        else
-            return d1;
+		if(d0 < d1)
+		    return d0;
+		else
+		    return d1;
 	}
 
 	CUDA_CALLABLE T center(vec<T> p)
@@ -232,13 +224,13 @@ public:
 	///@param p: x, y, z point.
 	CUDA_CALLABLE T dist_max(vec<T> p)
 	{
-		vec<T> A = C - X * (T)0.5 - Y * (T)0.5;
-        T da = (A - p).len();
-        T db = (A+X - p).len();
-        T dc = (A+Y - p).len();
-        T dd = (A+X+Y - p).len();
+		vec<T> A = P - X * (T)0.5 - Y * (T)0.5;
+		T da = (A - p).len();
+		T db = (A+X - p).len();
+		T dc = (A+Y - p).len();
+		T dd = (A+X+Y - p).len();
 
-        return std::max( da, std::max(db, std::max(dc, dd) ) );
+		return std::max( da, std::max(db, std::max(dc, dd) ) );
 	}
 };
 
diff --git a/stim/math/rect_old.h b/stim/math/rect_old.h
new file mode 100644
index 0000000..336a3ac
--- /dev/null
+++ b/stim/math/rect_old.h
@@ -0,0 +1,255 @@
+#ifndef RTS_RECT_H
+#define RTS_RECT_H
+
+//enable CUDA_CALLABLE macro
+#include <stim/cuda/cudatools/callable.h>
+#include <stim/math/vector.h>
+#include <stim/math/triangle.h>
+#include <stim/math/quaternion.h>
+#include <iostream>
+#include <iomanip>
+#include <algorithm>
+
+namespace stim{
+
+//template for a rectangle class in ND space
+template <class T>
+struct rect
+{
+	/*
+		^                   O
+		|                   
+		|                   
+		Y         C         
+		|                   
+		|                   
+		O---------X--------->
+	*/
+
+private:
+
+	stim::vec<T> C;
+	stim::vec<T> X;
+	stim::vec<T> Y;
+
+	CUDA_CALLABLE void scale(T factor){
+		X *= factor;
+		Y *= factor;
+	}
+	
+
+
+public:
+
+	///base constructor.
+	CUDA_CALLABLE rect(){
+		init();
+	}
+
+	///create a rectangle given a size and position in Z space.
+	///@param size: size of the rectangle in ND space.
+	///@param z_pos z coordinate of the rectangle.
+	CUDA_CALLABLE rect(T size, T z_pos = (T)0){
+		init();			//use the default setup
+		scale(size);	//scale the rectangle
+		C[2] = z_pos;
+	}
+
+	
+	///create a rectangle from a center point, normal
+	///@param c: x,y,z location of the center.
+	///@param n: x,y,z direction of the normal.
+	CUDA_CALLABLE rect(vec<T> c, vec<T> n = vec<T>(0, 0, 1)){
+		init();			//start with the default setting
+		C = c;
+		normal(n);		//orient
+	}
+
+	///create a rectangle from a center point, normal, and size
+	///@param c: x,y,z location of the center.
+	///@param s: size of the rectangle.
+	///@param n: x,y,z direction of the normal.
+	CUDA_CALLABLE rect(vec<T> c, T s, vec<T> n = vec<T>(0, 0, 1)){
+		init();			//start with the default setting
+		C = c;
+		scale(s);
+		normal(n);		//orient
+	}
+
+	///creates a rectangle from a centerpoint and an X and Y direction vectors.
+	///@param center: x,y,z location of the center.
+	///@param directionX: u,v,w direction of the X vector.
+	///@param directionY: u,v,w direction of the Y vector.
+	CUDA_CALLABLE rect(vec<T> center, vec<T> directionX, vec<T> directionY )
+	{
+		C = center;
+		X = directionX;
+		Y = directionY;
+	}
+
+	///creates a rectangle from a size, centerpoint, X, and Y direction vectors.
+	///@param size of the rectangle in ND space.
+	///@param center: x,y,z location of the center.
+	///@param directionX: u,v,w direction of the X vector.
+	///@param directionY: u,v,w direction of the Y vector.
+	CUDA_CALLABLE rect(T size, vec<T> center, vec<T> directionX, vec<T> directionY )
+	{	
+		C = center;
+		X = directionX;
+		Y = directionY;
+		scale(size);
+	}
+	
+	///creates a rectangle from a size, centerpoint, X, and Y direction vectors.
+	///@param size of the rectangle in ND space, size[0] = size in X, size[1] = size in Y.
+	///@param center: x,y,z location of the center.
+	///@param directionX: u,v,w direction of the X vector.
+	///@param directionY: u,v,w direction of the Y vector.
+	CUDA_CALLABLE rect(vec<T> size, vec<T> center, vec<T> directionX, vec<T> directionY )
+	{	
+		C = center;
+		X = directionX;
+		Y = directionY;
+		scale(size[0], size[1]);
+	}
+
+	///scales a rectangle in ND space.
+	///@param factor1: size of the scale in the X-direction.
+	///@param factor2: size of the scale in the Y-direction.	
+	CUDA_CALLABLE void scale(T factor1, T factor2){
+		X *= factor1;
+		Y *= factor2;
+	}
+
+	///@param n; vector with the normal.
+	///Orients the rectangle along the normal n.
+	CUDA_CALLABLE void normal(vec<T> n){		//orient the rectangle along the specified normal
+
+		n = n.norm();								//normalize, just in case
+		vec<T> n_current = X.cross(Y).norm();	//compute the current normal
+		quaternion<T> q;							//create a quaternion
+		q.CreateRotation(n_current, n);				//initialize a rotation from n_current to n
+
+		//apply the quaternion to the vectors and position
+		X = q.toMatrix3() * X;
+		Y = q.toMatrix3() * Y;
+	}
+
+	///general init method that sets a general rectangle.
+	CUDA_CALLABLE void init(){
+		C = vec<T>(0, 0, 0);
+		X = vec<T>(1, 0, 0);
+		Y = vec<T>(0, 1, 0);
+	}
+
+	//boolean comparison
+	bool operator==(const rect<T> & rhs)
+	{
+		if(C == rhs.C && X == rhs.X && Y == rhs.Y)
+			return true;
+		else
+			return false;
+	}
+
+	/*******************************************
+	Return the normal for the rect
+	*******************************************/
+	CUDA_CALLABLE stim::vec<T> n()
+	{
+        return (X.cross(Y)).norm();
+	}
+
+	//get the world space value given the planar coordinates a, b in [0, 1]
+	CUDA_CALLABLE stim::vec<T> p(T a, T b)
+	{
+		stim::vec<T> result;
+		//given the two parameters a, b = [0 1], returns the position in world space
+		vec<T> A = C - X * (T)0.5 - Y * (T)0.5;
+		result = A + X * a + Y * b;
+
+		return result;
+	}
+
+	//parenthesis operator returns the world space given rectangular coordinates a and b in [0 1]
+	CUDA_CALLABLE stim::vec<T> operator()(T a, T b)
+	{
+		return p(a, b);
+	}
+
+	std::string str()
+	{
+		std::stringstream ss;
+		vec<T> A = C - X * (T)0.5 - Y * (T)0.5;
+		ss<<std::left<<"B="<<std::setfill('-')<<std::setw(20)<<A + Y<<">"<<"C="<<A + Y + X<<std::endl;
+		ss<<std::setfill(' ')<<std::setw(23)<<"|"<<"|"<<std::endl<<std::setw(23)<<"|"<<"|"<<std::endl;
+		ss<<std::left<<"A="<<std::setfill('-')<<std::setw(20)<<A<<">"<<"D="<<A + X;
+
+        return ss.str();
+
+	}
+
+	///multiplication operator scales the rectangle by a value rhs.
+	CUDA_CALLABLE rect<T> operator*(T rhs)
+	{
+		//scales the plane by a scalar value
+
+		//create the new rectangle
+		rect<T> result = *this;
+		result.scale(rhs);
+
+		return result;
+
+	}
+
+	///computes the distance between the specified point and this rectangle.
+	///@param p: x, y, z coordinates of the point to calculate distance to.
+	CUDA_CALLABLE T dist(vec<T> p)
+	{
+        //compute the distance between a point and this rect
+
+		vec<T> A = C - X * (T)0.5 - Y * (T)0.5;
+
+        //first break the rect up into two triangles
+        triangle<T> T0(A, A+X, A+Y);
+        triangle<T> T1(A+X+Y, A+X, A+Y);
+
+
+        T d0 = T0.dist(p);
+        T d1 = T1.dist(p);
+
+        if(d0 < d1)
+            return d0;
+        else
+            return d1;
+	}
+
+	CUDA_CALLABLE T center(vec<T> p)
+	{
+		C = p;
+	}
+
+	///Returns the maximum distance of the rectangle from a point p to the sides of the rectangle.
+	///@param p: x, y, z point.
+	CUDA_CALLABLE T dist_max(vec<T> p)
+	{
+		vec<T> A = C - X * (T)0.5 - Y * (T)0.5;
+        T da = (A - p).len();
+        T db = (A+X - p).len();
+        T dc = (A+Y - p).len();
+        T dd = (A+X+Y - p).len();
+
+        return std::max( da, std::max(db, std::max(dc, dd) ) );
+	}
+};
+
+}	//end namespace rts
+
+template <typename T, int N>
+std::ostream& operator<<(std::ostream& os, stim::rect<T> R)
+{
+    os<<R.str();
+    return os;
+}
+
+
+#endif
--
libgit2 0.21.4