1 files changed, 107 insertions, 21 deletions

diff --git a/Smoke/fluids.c b/Smoke/fluids.c
index 9ae7d7f..4c96636 100644
--- a/Smoke/fluids.c
+++ b/Smoke/fluids.c
@@ -27,7 +27,8 @@ fftw_real *vx, *vy;             //(vx,vy)   = velocity field at the current mome
 fftw_real *vx0, *vy0;           //(vx0,vy0) = velocity field at the previous moment
 fftw_real *fx, *fy;	            //(fx,fy)   = user-controlled simulation forces, steered with the mouse 
 fftw_real *rho, *rho0;		    	//smoke density at the current (rho) and previous (rho0) moment 
-fftw_real *height_array;         //used for height plot
+fftw_real *height_array;        //used for height plot
+struct point *normal_array;     //used for normal vectors
 int *frame_hist;
 int frame_index = 0;
 rfftwnd_plan plan_rc, plan_cr;  //simulation domain discretization
@@ -81,7 +82,8 @@ fftw_real *init_simulation(int n)
 	plan_cr = rfftw2d_create_plan(n, n, FFTW_COMPLEX_TO_REAL, FFTW_IN_PLACE);
 
   height_array = (fftw_real*) malloc(dim);
-  frame_hist  = (fftw_real*) malloc(dim * n);
+  normal_array = (struct point*) malloc(dim);
+  frame_hist   = (fftw_real*) malloc(dim * n);
 
   for (i = 0; i <= n; i++)
   {
@@ -212,7 +214,107 @@ void set_forces(void)
 }
 
 
-void calculate_hight_plot(void)
+
+void vector_normal(struct point vert1, struct point vert2, struct point vert3, struct point *normal)
+{
+  struct point v1, v2;
+  float length;
+
+  // vector v1
+  v1.x = vert1.x - vert2.x;
+  v1.y = vert1.y - vert2.y;
+  v1.z = vert1.z - vert2.z;
+
+  // vector v2
+  v2.x = vert2.x - vert3.x;
+	v2.y = vert2.y - vert3.y;
+	v2.z = vert2.z - vert3.z;
+
+  // calculate cross produkt
+  normal->x = v1.y * v2.z - v1.z * v2.y;
+  normal->y = v1.z * v2.x - v1.x * v2.z;
+  normal->z = v1.x * v2.y - v1.y * v2.x;
+
+  length = vec_len3f(normal->x, normal->y, normal->z);
+
+  if(length == 0.0f)
+  {
+    length = 1.0f;
+  }
+
+  // normalize
+  normal->x /= length;
+  normal->y /= length;
+  normal->z /= length;
+}
+
+
+void calculate_normal_vectors(void)
+{
+  int i;
+  struct point p1, p2, p3;
+
+	for (i = 0; i < DIM * DIM; i++) 
+	{
+    switch (vis_dataset)
+    {
+		  case DATASET_FORCE:
+        p1.x = fx[i];
+        p1.y = fy[i];
+        p1.z = height_array[i];
+
+        p2.x = fx[i + 1];
+        p2.y = fy[i + 1];
+        p2.z = height_array[i + 1];
+
+        p3.x = fx[i + DIM];
+        p3.y = fy[i + DIM];
+        p3.z = height_array[i + DIM];
+			  break;
+
+		  case DATASET_VEL:
+        p1.x = vx[i];
+        p1.y = vy[i];
+        p1.z = height_array[i];
+
+        p2.x = vx[i + 1];
+        p2.y = vy[i + 1];
+        p2.z = height_array[i + 1];
+
+        p3.x = vx[i + DIM];
+        p3.y = vy[i + DIM];
+        p3.z = height_array[i + DIM];
+			  break;
+
+		  case DATASET_RHO:
+        p1.x = rho[i];
+        p1.y = rho[i];
+        p1.z = height_array[i];
+
+        p2.x = rho[i + 1];
+        p2.y = rho[i + 1];
+        p2.z = height_array[i + 1];
+
+        p3.x = rho[i + DIM];
+        p3.y = rho[i + DIM];
+        p3.z = height_array[i + DIM];
+			  break;
+
+      default:
+		  case DATASET_DIVV:
+		  case DATASET_DIVF:
+        p1.x = 0.0f; p1.y = 0.0f; p1.z = 0.0f;
+        p2.x = 0.0f; p2.y = 0.0f; p2.z = 0.0f;
+        p3.x = 0.0f; p3.y = 0.0f; p3.z = 0.0f;
+			  break;
+	  }
+
+    vector_normal(p1, p2, p3, &normal_array[i]);
+  }
+}
+
+
+void calculate_height_plot(void)
 {
   int i;
 	for (i = 0; i < DIM * DIM; i++) 
@@ -240,23 +342,6 @@ void calculate_hight_plot(void)
 void copy_frame(fftw_real *field)
 {
   memcpy(field, rho, DIM * 2 * (DIM / 2 + 1) * sizeof(fftw_real));
-#if 0
-  switch (vis_dataset)
-  {
-	  case DATASET_FORCE:
-		  break;
-	  case DATASET_VEL:
-		  break;
-	  case DATASET_RHO:
-      //memcpy(frame_hist[frame_index], rho, sizeof(DIM * 2 * (DIM / 2 + 1) * sizeof(fftw_real)));
-      memcpy(frame_hist[frame_index], rho, DIM * 2 * (DIM / 2 + 1) * sizeof(fftw_real));
-      break;
-    default:
-	  case DATASET_DIVV:
-	  case DATASET_DIVF:
-		  break;
-  }
-#endif
 }
 
 //do_one_simulation_step: Do one complete cycle of the simulation:
@@ -272,7 +357,8 @@ void calculate_one_simulation_step(fftw_real *field)
       set_forces();
       solve(DIM, vx, vy, vx0, vy0, visc, dt);
       diffuse_matter(DIM, vx, vy, rho, rho0, dt);
-      calculate_hight_plot();
+      calculate_height_plot();
+      calculate_normal_vectors();
       copy_frame(field);
   }
 }