Merge branch 'master' of git.stim.ee.uh.edu:codebase/stimlib

David Mayerich
2 parents 474c351a 9c5854ec
Showing 2 changed files with 151 additions and 66 deletions Show diff stats
stim/biomodels/cellset.h
stim/gl/gl_spider.h
@@ -3,6 +3,7 @@
 #include <stim/math/vec3.h>
 #include <vector>
+#include <list>
 #include <unordered_map>
 #include <fstream>
@@ -129,4 +130,4 @@ public:
 };		//end class cellset
 };		//end namespace stim
-#endif
 \ No newline at end of file
+#endif
@@ -68,6 +68,7 @@ class gl_spider // : public virtual gl_texture&lt;T&gt;
 		std::vector<stim::vec3<float> > dV;	//A list of all the direction vectors.
 		std::vector<stim::vec3<float> > pV;	//A list of all test positions (relative to p)
 		std::vector<float> mV;	//A list of all the size vectors.
+		std::vector<float> lV;	//A list of all the size vectors.
 		stim::matrix<float, 4> cT;			//current Transformation matrix (tissue)->(texture)
 		GLuint texID;						//OpenGL ID for the texture to be traced
@@ -93,6 +94,9 @@ class gl_spider // : public virtual gl_texture&lt;T&gt;
 		GLuint radius_buffID;				//framebuffer ID, radius templates
 		GLuint radius_texID;		//texture ID, radius templates
+		GLuint length_buffID;				//framebuffer ID, radius templates
+		GLuint length_texID;		//texture ID, radius templates
+
 		GLuint cylinder_buffID;				//framebuffer ID, cylinder (surrounding fiber)
 		GLuint cylinder_texID;		//texture ID, cylinder
@@ -113,6 +117,7 @@ class gl_spider // : public virtual gl_texture&lt;T&gt;
 		std::vector< stim::vec3<float> > cL;	//centerline up to the current point
 		std::vector< stim::vec3<float> > cD;	//directions up to the current point (debugging)
 		std::vector< float > cM;		//radius up to the current point
+		std::vector< float > cLen;		//radius up to the current point
 		stim::glnetwork<float> nt;				//network object holding the currently traced centerlines
 		stim::glObj<float> sk;					//OBJ file storing the network (identical to above)
@@ -132,6 +137,7 @@ class gl_spider // : public virtual gl_texture&lt;T&gt;
 		stim::cuda::cuda_texture t_dir;				//cuda_texture object used as an interface between OpenGL and cuda for direction vectors.
 		stim::cuda::cuda_texture t_pos;				//cuda_texture object used as an interface between OpenGL and cuda for position vectors.
 		stim::cuda::cuda_texture t_mag;				//cuda_texture object used as an interface between OpenGL and cuda for size vectors.
+		stim::cuda::cuda_texture t_len;				//cuda_texture object used as an interface between OpenGL and cuda for size vectors.
 		#ifdef DEBUG
@@ -222,7 +228,7 @@ class gl_spider // : public virtual gl_texture&lt;T&gt;
 		/// changes the x, y, z size of the spider to minimize the cost
 		/// function. 
 		void
-		findOptimalScale()
+		findOptimalRadius()
 		{
 			#ifdef TIMING
 				gpuStartTimer();
@@ -243,6 +249,31 @@ class gl_spider // : public virtual gl_texture&lt;T&gt;
 			setMagnitude(m*mV[best]);	//adjust the magnitude.
 		}
+		/// Method for finding the best length for the spider.
+		/// changes the x, y, z size of the spider to minimize the cost
+		/// function. 
+		void
+		findOptimalLength()
+		{
+			#ifdef TIMING
+				gpuStartTimer();
+			#endif
+			setMatrix();			//create the transformation.
+			glCallList(dList+3);		//move the templates to p, d, m.
+			glFinish();			//flush the drawing pipeline.
+			#ifdef TIMING
+				size_time += gpuStopTimer();
+			#endif
+			int best = getCost(t_len.getTexture(), t_len.getAuxArray(), numSamplesMag);		//get best cost.
+			#ifdef DEBUG
+//				name.str("");
+//				name << "Final_Cost_Size_" << iter << "_" <<  iter_siz << ".bmp";	
+//				test(t_mag.getTexture(), n_pixels*2.0, numSamplesMag*n_pixels, name.str());
+//				iter_siz++;
+			#endif
+			setLength(mV[best]);	//adjust the magnitude.
+		}
+
@@ -485,6 +516,48 @@ class gl_spider // : public virtual gl_texture&lt;T&gt;
 			glEndList();					///finilize the displaylist.
 		}
+		///@param float delta, How much the rectangles are allowed to expand.
+		///Method for populating the buffer with the sampled texture.
+		///Objects created are rectangles the with the created sizes.
+		///All points are sampled from a texture.
+		///Stored in a display list.
+		///uses the default m <1,1,0>
+		void
+		genLengthVectors(float delta = 0.70)
+		{
+			
+			//Set up the vectors necessary for Rectangle creation.
+			stim::vec3<float> Y(1.0, 0.0, 0.0);		//orthogonal vec.
+			stim::vec3<float> pos(0.0, 0.0, 0.0);		//center of the future rect.
+			float mag = 1.0;				///size of the rectangle
+			stim::vec3<float> dir(0.0, 0.0, 1.0);		///normal of the rectangle plane.
+			stim::vec<float> temp(0.0,0.0,0.0);
+
+			//Set up the variable necessary for vector creation.
+			float min 	= 1.0-delta;			///smallest size
+			float max 	= 1.0+delta;			///largers size.
+			float step	= (max-min)/(numSamplesMag-1);	///the size variation from one rect to the next.
+			float factor;
+			glNewList(dList+3, GL_COMPILE);
+			for(int i = 0; i < numSamplesMag; i++){		///for the number of position samples
+				//Create linear index
+				factor = (min+step*i)*mag;		///scaling factor
+				lV.push_back(factor);			///save the size factor for further use.
+				temp[0] = factor;
+				temp[1] = mag;
+				hor = stim::rect<float>(temp,		///generate a rectangle with the new vector as a normal.
+					 pos, dir, 
+       	       				((Y.cross(d)).cross(d))
+					.norm());
+               			ver = stim::rect<float>(temp,		///generate another rectangle that's perpendicular the first but parallel to the cart vector.
+					 pos, dir,
+                       			hor.n());
+				UpdateBuffer(0.0, 0.0+i*n_pixels);	///sample the necessary points and put them into a display list.
+			CHECK_OPENGL_ERROR
+			}
+			glEndList();					///finilize the displaylist.
+		}
+
 		///@param float v_x x-coordinate in buffer-space,
 		///@param float v_y y-coordinate in buffer-space.
 		///Samples the texture space.
@@ -654,21 +727,21 @@ class gl_spider // : public virtual gl_texture&lt;T&gt;
 				GL_COLOR_ATTACHMENT0,
 				GL_TEXTURE_2D,
 				direction_texID,
-				0);
-			glBindFramebuffer(GL_FRAMEBUFFER, direction_buffID);
-			GLenum DrawBuffers[1] = {GL_COLOR_ATTACHMENT0};
-			glDrawBuffers(1, DrawBuffers);
-			glBindTexture(GL_TEXTURE_2D, direction_texID);
-			glClearColor(1,1,1,1);
-			glClear(GL_COLOR_BUFFER_BIT);
-			glMatrixMode(GL_PROJECTION);
+				0);			///Bind the texture to the 0th color attachement of the framebuffer
+			glBindFramebuffer(GL_FRAMEBUFFER, direction_buffID);	///Bind the buffer again (safety operation).
+			GLenum DrawBuffers[1] = {GL_COLOR_ATTACHMENT0};		///Designate the texture to be the drawbuffer of the framebuffer
+			glDrawBuffers(1, DrawBuffers);				///Set the current drawbuffer to the texture.
+			glBindTexture(GL_TEXTURE_2D, direction_texID);		///Bind the Texture
+			glClearColor(1,1,1,1);					///Set clear color to white
+			glClear(GL_COLOR_BUFFER_BIT);				///Clear the texture
+			glMatrixMode(GL_PROJECTION);				
 			glLoadIdentity();
 			glMatrixMode(GL_MODELVIEW);
-			glLoadIdentity();
-			glViewport(0,0,2.0*n_pixels, numSamples*n_pixels);
+			glLoadIdentity();					///Load identity matrix into the projection and modelview
+			glViewport(0,0,2.0*n_pixels, numSamples*n_pixels);	///Designate the viewport and orth
 			gluOrtho2D(0.0,2.0*n_pixels,0.0,numSamples*n_pixels);
 			glEnable(GL_TEXTURE_3D);
-			glBindTexture(GL_TEXTURE_3D, texID);
+			glBindTexture(GL_TEXTURE_3D, texID);			///Bind the larger 3D texture.
 			CHECK_OPENGL_ERROR
 		}
@@ -735,12 +808,12 @@ class gl_spider // : public virtual gl_texture&lt;T&gt;
 		getCost(cudaTextureObject_t tObj, float* result, int n)
 		{
 			#ifdef TIMING
-				gpuStartTimer();
+				gpuStartTimer();			///Add timing variables
 			#endif
 			stim::vec<int> cost = 
 				stim::cuda::get_cost(tObj, result, n, 2*n_pixels, n_pixels);	///call the cuda function with the appropriate texture buffer.
 			#ifdef TIMING
-				cost_time += gpuStopTimer();
+				cost_time += gpuStopTimer();		
 			#endif
 			current_cost = cost[1];							///current cost.
 			return cost[0];
@@ -761,9 +834,6 @@ class gl_spider // : public virtual gl_texture&lt;T&gt;
 		stim::rect<float> ver;	
 		//Timing variables.
-		#ifdef TESTING
-			std::clock_t start;
-		#endif
 //--------------------------------------------------------------------------//
 //-----------------------------CONSTRUCTORS---------------------------------//
@@ -845,7 +915,7 @@ class gl_spider // : public virtual gl_texture&lt;T&gt;
 		#endif
 		}
-		///destructor	
+		///destructor deletes all the texture and buffer objects.
 		~gl_spider
 		(void)
 		{
@@ -856,6 +926,8 @@ class gl_spider // : public virtual gl_texture&lt;T&gt;
 			glDeleteBuffers(1, &position_buffID);
 			glDeleteTextures(1, &radius_texID);
 			glDeleteBuffers(1, &radius_buffID);     
+			glDeleteTextures(1, &length_texID);
+			glDeleteBuffers(1, &length_buffID);     
 			glDeleteTextures(1, &cylinder_texID);
 			glDeleteBuffers(1, &cylinder_buffID);
 		}
@@ -882,7 +954,7 @@ class gl_spider // : public virtual gl_texture&lt;T&gt;
 			iter_dir = 0;
 			iter_siz = 0;
 #endif
-			stepsize = 10.0;
+			stepsize = 3.0;
 			n_pixels = 16.0;
 			srand(100);	
@@ -893,6 +965,8 @@ class gl_spider // : public virtual gl_texture&lt;T&gt;
 			CHECK_OPENGL_ERROR
 			GenerateFBO(n_pixels*2, numSamplesMag*n_pixels, radius_texID, radius_buffID); 
 			CHECK_OPENGL_ERROR
+			GenerateFBO(n_pixels*2, numSamplesMag*n_pixels, length_texID, length_buffID); 
+			CHECK_OPENGL_ERROR
 			GenerateFBO(16, 216, cylinder_texID, cylinder_buffID); 
 			CHECK_OPENGL_ERROR
 			t_dir.MapCudaTexture(direction_texID, GL_TEXTURE_2D);
@@ -901,8 +975,10 @@ class gl_spider // : public virtual gl_texture&lt;T&gt;
 			t_pos.Alloc(numSamplesPos);
 			t_mag.MapCudaTexture(radius_texID, GL_TEXTURE_2D);
 			t_mag.Alloc(numSamplesMag);
+			t_len.MapCudaTexture(length_texID, GL_TEXTURE_2D);
+			t_len.Alloc(numSamplesMag);
 			setMatrix();
-			dList = glGenLists(3);
+			dList = glGenLists(4);
 			glListBase(dList);
 			Bind(direction_texID, direction_buffID, numSamples, n_pixels);
 				genDirectionVectors(5*stim::PI/4);
@@ -913,8 +989,8 @@ class gl_spider // : public virtual gl_texture&lt;T&gt;
 			Bind(radius_texID, radius_buffID, numSamplesMag, n_pixels);
 				genMagnitudeVectors();
 			Unbind();
-			Bind(cylinder_texID, cylinder_buffID, 27);
-//				DrawCylinder();
+			Bind(length_texID, length_buffID, numSamplesMag, n_pixels);
+				genLengthVectors();
 			Unbind();
 		}
@@ -984,12 +1060,20 @@ class gl_spider // : public virtual gl_texture&lt;T&gt;
 		///@param float mag, size of the sampled region.
-		///Sets the m vector to the input mag for both templates.
+		///Sets the m value to the input mag for both templates.
 		void
 		setMagnitude(float mag)
 		{
 			m = mag;
 		}
+
+		///@param float len, length of the sampled region.
+		///Sets the length value to the input len for both templates.
+		void
+		setLength(float len)
+		{
+			length = len;
+		}
 		///@param float x, voxel size in the x direction.
 		///@param float y, voxel size in the y direction.
@@ -1038,32 +1122,32 @@ class gl_spider // : public virtual gl_texture&lt;T&gt;
 		stim::vec<float>
 		getRotation(stim::vec3<float> dir)
 		{
-			stim::vec<float> out(0.0,0.0,0.0,0.0);
-			stim::vec3<float> from(0.0,0.0,1.0);
-			out[0] = acos(dir.dot(from))*180/stim::PI;
-			#ifdef DEBUG
-				std::cout << "out is " << out << std::endl;
-				std::cout << "when rotating from " << from << " to " << dir << std::endl;
-			#endif
-			if(out[0] < 0.01){
+			stim::vec<float> out(0.0,0.0,0.0,0.0);	///The 4D rotation matrix for GL rotation
+			stim::vec3<float> from(0.0,0.0,1.0); 	///Converting from template always involves 0,0,1 as a starting vector
+			out[0] = acos(dir.dot(from))*180/stim::PI; ///angle of rotation
+			if(out[0] < 1.0){
 				out[0] = 0.0;
 				out[1] = 0.0;
 				out[2] = 0.0;
 				out[3] = 1.0;
-			}
-			else if(out[0] < -180.0+0.01)
+			}					///if we rotate less what one degree don't rotate
+			else if(out[0] < -179.0)		///if we rotate more than -179 degrees rotate 180.
 			{
 				out[0] = 180.0;
 				out[1] = 1.0;
 				out[2] = 0.0;
 				out[3] = 0.0;
-			} else {
+			} else {				///the rotational axis is the cross fromxdir.
 				stim::vec3<float> temp(0.0, 0.0, 0.0);;
 				temp = (from.cross(dir)).norm();
 				out[1] = temp[0];
 				out[2] = temp[1];
 				out[3] = temp[2];
 			}
+			#ifdef DEBUG
+				std::cout << "out is " << out << std::endl;
+				std::cout << "when rotating from " << from << " to " << dir << std::endl;
+			#endif
 			return out;
 		}
@@ -1210,13 +1294,13 @@ class gl_spider // : public virtual gl_texture&lt;T&gt;
 		void
 		setSeeds(std::string file)
 		{
-			std::ifstream myfile(file.c_str());
+			std::ifstream myfile(file.c_str()); 	///open a stream
 			string line;
 			if(myfile.is_open())
            		{
                    		while (getline(myfile, line))
                    		{
-                           		float x, y, z, u, v, w, m;
+                           		float x, y, z, u, v, w, m;	///read the xyz uvw and m coordinates.
                            		myfile >> x >> y >> z >> u >> v >> w >> m;
 					setSeed(x, y, z);
 					setSeedVec(u, v, w);
@@ -1303,7 +1387,7 @@ class gl_spider // : public virtual gl_texture&lt;T&gt;
 				findOptimalPosition();
 			Unbind();
 			Bind(radius_texID, radius_buffID, numSamplesMag, n_pixels);
-				findOptimalScale();
+				findOptimalRadius();
 			Unbind();
 			CHECK_OPENGL_ERROR
@@ -1351,12 +1435,13 @@ class gl_spider // : public virtual gl_texture&lt;T&gt;
 			outFile.open("New_Pos_Vectors.txt");
 			outFile << name.str().c_str();
 		}
+/*
 		void
 		DrawCylinder()
 		{	 
 			 glNewList(dList+3, GL_COMPILE);
 			 float z0 = -0.5; float z1 = 0.5; float r0 = 0.5;
-			 float x,y;
+      			 float x,y;
 			 float xold = 0.5; float yold = 0.0;
 			 float step = 360.0/numSamples*32;
 			 //float step = 360.0/8.0;
@@ -1387,7 +1472,7 @@ class gl_spider // : public virtual gl_texture&lt;T&gt;
 		      	 glEnd();  
 			 glEndList();
 		}
-
+*/
 ///need to return the cylinder.
 ///SOMETHING MIGHT BE GOING ON HERE IN GENERATE BUFFER.
 		void
@@ -1455,7 +1540,7 @@ class gl_spider // : public virtual gl_texture&lt;T&gt;
 //			Unbind();
 //THIS IS EXPERIMENTAL
 			Bind(radius_texID, radius_buffID, numSamplesMag, n_pixels);
-				findOptimalScale();
+				findOptimalRadius();
 			Unbind();
 //THIS IS EXPERIMENTAL
@@ -1475,31 +1560,30 @@ class gl_spider // : public virtual gl_texture&lt;T&gt;
 				gpuStartTimer();
 			#endif
-//			float s = 3.0;
-			GLuint selectBuf[2048];
-			GLint hits;
-			glSelectBuffer(2048, selectBuf);
-			glDisable(GL_CULL_FACE);
-			(void) glRenderMode(GL_SELECT);
+			GLuint selectBuf[2048];			///size of the selection buffer in bytes.
+			GLint hits;				///hit fibers
+			glSelectBuffer(2048, selectBuf);	///bind the selection mode to the selection buffer.
+			glDisable(GL_CULL_FACE);		///Disable cullFace
+			(void) glRenderMode(GL_SELECT);		///initialize GL select mode.
 		//Init Names stack
-			glInitNames();
-			glPushName(1);
+			glInitNames();				///Initialize the naming array.
+			glPushName(1);				///Push a single name to the stack.
 			CHECK_OPENGL_ERROR
 		//What would that vessel see in front of it.
-				camSel.setPosition(loc);
-				camSel.setFocalDistance(mag/stepsize);
+				camSel.setPosition(loc);		///Set the viewing camera
+				camSel.setFocalDistance(mag/stepsize);	///Set how far the fiber looks forward.
 				camSel.LookAt((loc[0]+dir[0]*mag/stepsize),
 					 (loc[1]+dir[1]*mag/stepsize),
-					 (loc[2]+dir[2]*mag/stepsize));
-				ps  = camSel.getPosition();
+					 (loc[2]+dir[2]*mag/stepsize)); ///Set the look direction
+				ps  = camSel.getPosition();		///get all the necessary rotation variable for openGL
 				ups = camSel.getUp();
 				ds  = camSel.getLookAt();
-				glMatrixMode(GL_PROJECTION);
-				glPushMatrix();
+				glMatrixMode(GL_PROJECTION);		///Push the projection matrix.
+				glPushMatrix();				///Reset the current projection matrix
 				glLoadIdentity();
-				glOrtho(-mag/stepsize/2.0, mag/stepsize/2.0, -mag/stepsize/2.0, mag/stepsize/2.0, 0.0, mag/stepsize/2.0);
+				glOrtho(-mag/stepsize/2.0, mag/stepsize/2.0, -mag/stepsize/2.0, mag/stepsize/2.0, 0.0, mag/stepsize/2.0);			///Finalize the look paramenters
 				glMatrixMode(GL_MODELVIEW);
 				glPushMatrix();
 				glLoadIdentity();
@@ -1508,36 +1592,36 @@ class gl_spider // : public virtual gl_texture&lt;T&gt;
 				gluLookAt(ps[0], ps[1], ps[2],
 					 ds[0], ds[1], ds[2],
 					 ups[0], ups[1], ups[2]);
-
-				sk.Render();
+				///Set the look at distance
+				sk.Render();	///Render the network
 //				nt.Render();								
 				CHECK_OPENGL_ERROR
-				glLoadName((int) sk.numL());
+				glLoadName((int) sk.numL());		///Load all the names
 //				glLoadName(nt.sizeE());
-				sk.RenderLine(cL);
+				sk.RenderLine(cL);			///Render the current line.
 //				nt.RenderLine(cL);	
 //				glPopName();
-				glFlush();
+				glFlush();				///Flush the buffer
 				glMatrixMode(GL_PROJECTION);
 				glPopMatrix();
 				glMatrixMode(GL_MODELVIEW);
 				CHECK_OPENGL_ERROR
-				glPopMatrix();
+				glPopMatrix();				///clear the vis matrices and pop the matrix
 		//	glEnable(GL_CULL_FACE);
-			hits = glRenderMode(GL_RENDER);
-			int found_hits = processHits(hits, selectBuf);
+			hits = glRenderMode(GL_RENDER);			///Check for hits.
+			int found_hits = processHits(hits, selectBuf);	///Process the hits.
 			#ifdef TIMING
 				hit_time += gpuStopTimer();
 			#endif
-			return found_hits;
+			return found_hits;				///return whether we hit something or not.
 		}
 		//Given a size of the array (hits) and the memory holding it (buffer)
@@ -1545,7 +1629,7 @@ class gl_spider // : public virtual gl_texture&lt;T&gt;
 		int
 		processHits(GLint hits, GLuint buffer[])
 		{
-			GLuint *ptr;
+			GLuint *ptr;		///pointer to the detection buffer
 			ptr = (GLuint *) buffer;
 			ptr++;
 			ptr++; //Skip the minimum depth value.