Draw the sun, debug the colour and variable size

[clouds.git] / simulation.cpp

diff --git a/simulation.cpp b/simulation.cpp
index 360fb07481a2702188c6afd63ea9d34b72c27846..4e107e72670a2d86452d0383fcdf8563c786ce4a 100644 (file)
--- a/simulation.cpp
+++ b/simulation.cpp
@@ -7,10 +7,10 @@ inline float randf() {
 }
 
 // Helper to account for bounds
-inline void set(float x[CLOUD_DIM][CLOUD_DIM][CLOUD_DIM], int i, int j, int k, float y) {
-  if (i < 0 || i >= CLOUD_DIM ||
-      j < 0 || j >= CLOUD_DIM ||
-      k < 0 || k >= CLOUD_DIM)
+inline void set(float x[CLOUD_DIM_X][CLOUD_DIM_Y][CLOUD_DIM_Z], int i, int j, int k, float y) {
+  if (i < 0 || i >= CLOUD_DIM_X ||
+      j < 0 || j >= CLOUD_DIM_Y ||
+      k < 0 || k >= CLOUD_DIM_Z)
     return;
   x[i][j][k] = y;
 }
@@ -20,12 +20,12 @@ inline void set(float x[CLOUD_DIM][CLOUD_DIM][CLOUD_DIM], int i, int j, int k, f
 #define P_ACT 0.001
 
 void initClouds(Clouds *cs) {
-  for (int i = 0; i < CLOUD_DIM; i++) {
-    for (int j = 0; j < CLOUD_DIM; j++) {
-      for (int k = 0; k < CLOUD_DIM; k++) {
+  for (int i = 0; i < CLOUD_DIM_X; i++) {
+    for (int j = 0; j < CLOUD_DIM_Y; j++) {
+      for (int k = 0; k < CLOUD_DIM_Z; k++) {
         cs->act[i][j][k] = randf() < 0.01;
         cs->cld[i][j][k] = false;
-        cs->hum[i][j][k] = randf() < 0.01;
+        cs->hum[i][j][k] = randf() < 0.1;
         cs->p_ext[i][j][k] = 0.f;
         cs->p_hum[i][j][k] = 0.f;
         cs->p_act[i][j][k] = 0.f;
@@ -33,22 +33,25 @@ void initClouds(Clouds *cs) {
     }
   }
 
-  for (int i = 0; i < CLOUD_DIM; i++) {
-    for (int j = 0; j < CLOUD_DIM; j++) {
-      for (int k = 0; k < CLOUD_DIM; k++) {
+  for (int i = 0; i < CLOUD_DIM_X; i++) {
+    for (int j = 0; j < CLOUD_DIM_Y; j++) {
+      for (int k = 0; k < CLOUD_DIM_Z; k++) {
         assert(cs->p_act[i][j][k] == 0.f);
         assert(cs->p_ext[i][j][k] == 0.f);
         assert(cs->p_hum[i][j][k] == 0.f);
       }
     }
   }
+  for (int k = 0; k < CLOUD_DIM_Z; k++)
+    cs->vz[k] = floor(randf() * 3);
 
   // generate ellipsoids of probability
-  for (int n = 0; n < 6; n++) {
+  const int numEllipsoids = CLOUD_DIM_X * CLOUD_DIM_Y * CLOUD_DIM_Z * 0.002;
+  for (int n = 0; n < numEllipsoids; n++) {
     const float maxSize = 8, minSize = 4;
     float delta = maxSize - minSize;
     int width = minSize + randf() * delta, height = minSize + randf() * delta, depth = minSize + randf() * delta;
-    int x = randf() * CLOUD_DIM, y = randf() * CLOUD_DIM, z = randf() * CLOUD_DIM;
+    int x = randf() * CLOUD_DIM_X, y = randf() * CLOUD_DIM_Y, z = randf() * CLOUD_DIM_Z;
     glm::vec3 center(x + width / 2, y + height / 2, z + depth / 2);
 
     for (int i = x; i < x + width; i++) {
@@ -65,10 +68,10 @@ void initClouds(Clouds *cs) {
 }
 
 // Helper to account for bounds
-inline bool get(bool x[CLOUD_DIM][CLOUD_DIM][CLOUD_DIM], int i, int j, int k) {
-  if (i < 0 || i >= CLOUD_DIM ||
-      j < 0 || j >= CLOUD_DIM ||
-      k < 0 || k >= CLOUD_DIM)
+inline bool get(bool x[CLOUD_DIM_X][CLOUD_DIM_Y][CLOUD_DIM_Z], int i, int j, int k) {
+  if (i < 0 || i >= CLOUD_DIM_X ||
+      j < 0 || j >= CLOUD_DIM_Y ||
+      k < 0 || k >= CLOUD_DIM_Z)
     return false;
   return x[i][j][k];
 }
@@ -83,9 +86,9 @@ inline bool f_act(Clouds *cs, int i, int j, int k) {
 void growth(Clouds *cs) {
   Clouds ncs = *cs;
 
-  for (int i = 0; i < CLOUD_DIM; i++) {
-    for (int j = 0; j < CLOUD_DIM; j++) {
-      for (int k = 0; k < CLOUD_DIM; k++) {
+  for (int i = 0; i < CLOUD_DIM_X; i++) {
+    for (int j = 0; j < CLOUD_DIM_Y; j++) {
+      for (int k = 0; k < CLOUD_DIM_Z; k++) {
         ncs.hum[i][j][k] = cs->hum[i][j][k] && !cs->act[i][j][k];
         ncs.cld[i][j][k] = cs->cld[i][j][k] || cs->act[i][j][k];
         ncs.act[i][j][k] = !cs->act[i][j][k] && cs->hum[i][j][k] && f_act(cs, i, j, k);
@@ -98,9 +101,9 @@ void growth(Clouds *cs) {
 
 void extinction(Clouds *cs) {
   Clouds ncs = *cs;
-  for (int i = 0; i < CLOUD_DIM; i++) {
-    for (int j = 0; j < CLOUD_DIM; j++) {
-      for (int k = 0; k < CLOUD_DIM; k++) {
+  for (int i = 0; i < CLOUD_DIM_X; i++) {
+    for (int j = 0; j < CLOUD_DIM_Y; j++) {
+      for (int k = 0; k < CLOUD_DIM_Z; k++) {
         ncs.cld[i][j][k] = cs->cld[i][j][k] && (randf() > cs->p_ext[i][j][k]);
         ncs.hum[i][j][k] = cs->hum[i][j][k] || (randf() < cs->p_hum[i][j][k]);
         ncs.act[i][j][k] = cs->act[i][j][k] || (randf() < cs->p_act[i][j][k]);
@@ -110,16 +113,34 @@ void extinction(Clouds *cs) {
   *cs = ncs;
 }
 
+void advection(Clouds *cs) {
+  Clouds ncs = *cs;
+
+  for (int i = 0; i < CLOUD_DIM_X; i++) {
+    for (int j = 0; j < CLOUD_DIM_Y; j++) {
+      for (int k = 0; k < CLOUD_DIM_Z; k++) {
+        int v = cs->vz[k];
+        ncs.hum[i][j][k] = i - v > 0 ? cs->hum[i - v][j][k] : 0;
+        ncs.cld[i][j][k] = i - v > 0 ? cs->cld[i - v][j][k] : 0;
+        ncs.act[i][j][k] = i - v > 0 ? cs->act[i - v][j][k] : 0;
+      }
+    }
+  }
+  
+  *cs = ncs;
+}
+
 /** Weighting function */
+// TODO: fill this out
 float w(int ip, int jp, int kp) { return 1; }
 
 void calcContDist(Clouds *cls) {
   const int i0 = 2, j0 = 2, k0 = 2, t0 = 0;
   const float divisor =
       1.f / ((2 * t0 + 1) * (2 * k0 + 1) * (2 * j0 + 1) * (2 * i0 + 1));
-  for (int i = 0; i < CLOUD_DIM; i++) {
-    for (int j = 0; j < CLOUD_DIM; j++) {
-      for (int k = 0; k < CLOUD_DIM; k++) {
+  for (int i = 0; i < CLOUD_DIM_X; i++) {
+    for (int j = 0; j < CLOUD_DIM_Y; j++) {
+      for (int k = 0; k < CLOUD_DIM_Z; k++) {
         float sum = 0;
 
         // sum
@@ -127,9 +148,9 @@ void calcContDist(Clouds *cls) {
           for (int ip = -i0; ip <= i0; ip++) {
             for (int jp = -j0; jp <= j0; jp++) {
               for (int kp = -k0; kp <= k0; kp++) {
-                if (i + ip < 0 || i + ip >= CLOUD_DIM ||
-                    j + jp < 0 || j + jp >= CLOUD_DIM ||
-                    k + kp < 0 || k + kp >= CLOUD_DIM)
+                if (i + ip < 0 || i + ip >= CLOUD_DIM_X ||
+                    j + jp < 0 || j + jp >= CLOUD_DIM_Y ||
+                    k + kp < 0 || k + kp >= CLOUD_DIM_Z)
                   continue;
 
                 sum += w(ip, jp, kp) * (float)cls->cld[i + ip][j + jp][k + kp];
@@ -148,5 +169,6 @@ void calcContDist(Clouds *cls) {
 void stepClouds(Clouds *cs) {
   growth(cs);
   extinction(cs);
+  advection(cs);
   calcContDist(cs);
 }