PBO stuff
authorLuke Lau <luke_lau@icloud.com>
Mon, 13 Apr 2020 23:18:37 +0000 (00:18 +0100)
committerLuke Lau <luke_lau@icloud.com>
Mon, 13 Apr 2020 23:18:37 +0000 (00:18 +0100)
Makefile
billboardfrag.glsl
clouds.cpp
debug.hpp
simulation.cpp
simulation.hpp [moved from simulation.h with 87% similarity]

index cf0c407248271c0911ebad0dc52473cee2c102bf..61f80a70d3596d4712183bef627e681b815f620b 100644 (file)
--- a/Makefile
+++ b/Makefile
@@ -1,5 +1,5 @@
 clouds: clouds.cpp program.cpp simulation.cpp
-       clang++ -std=c++17 $^ -O3 -march=native -o $@ \
+       clang++ -std=c++17 -O3 $^ -march=native -o $@ \
                -Wall -g \
                -I/usr/local/include -L/usr/local/lib \
                -framework OpenGL -framework glut -lglew
index 39fd71d135d1d36fb9fe820fc9eef3a79082e9cc..b4cfd511464569ea23c577cda41a5b261c7ba288 100644 (file)
@@ -25,7 +25,9 @@ void main() {
     // GL_MODULATE: C
     // C = Cf * Ct
     // A = Af * At
-               FragColor = color * (f + 0.03); // the +0.03 is a hack to get lighter clouds!
+    // the +0.06 is a hack to get lighter clouds!
+    // can be thought of as ambient light
+    FragColor = color * (f + 0.02);
   } else {
     // "That is, the colors in the frame buffer are multiplied by the
     // attenuation ratio of the billboard texture and then the colors in
index f6c0083837e15611cc628e5316e3c55fd63cd4d5..f0f3a957576d17da6fdf5b60c521174e060eaeec 100644 (file)
@@ -1,13 +1,13 @@
+#include "GL/glew.h"
 #include "debug.hpp"
 #include "program.hpp"
-#include "simulation.h"
-#include <GL/glew.h>
+#include "simulation.hpp"
 #include <GLUT/glut.h>
 #include <array>
+#include <chrono>
 #include <cmath>
 #include <cstdio>
 #include <cstdlib>
-#include <chrono>
 #include <glm/ext.hpp>
 #include <glm/glm.hpp>
 #include <sys/stat.h>
@@ -63,7 +63,7 @@ void precalculateBillboardTextures() {
         float r = distance(vec2(i, j), vec2(16, 16)) / 16;
         float density = (float)d / NQ;
         data[i + j * 32] =
-            1 - (3 * density * 0.7 * (metaballField(r) / normalizationFactor));
+            1 - fmin(1, (3 * density * (metaballField(r) / normalizationFactor)));
       }
     }
 
@@ -93,7 +93,7 @@ struct Metaball {
 array<Metaball, CLOUD_DIM_X * CLOUD_DIM_Y * CLOUD_DIM_Z> metaballs;
 
 const float cloudScale = metaballR;
-const float metaballScale = metaballR * 3.f;
+const float metaballScale = metaballR * 1.5f;
 Clouds cs;
 
 void calculateMetaballs() {
@@ -123,7 +123,7 @@ void calculateMetaballs() {
    */
 }
 
-vec3 sunPos = {0, 10, 0}, sunDir = {0, -1, 0};
+vec3 sunPos = {0, 5, 0}, sunDir = {0, -1, 0};
 size_t envColorIdx = 0;
 // First color is sun color, second is sky color
 std::array<std::array<vec4, 2>, 3> envColors{
@@ -134,8 +134,8 @@ vec3 camPos = {0, 0, -3}, viewPos = {0, 0, 0};
 mat4 proj; // projection matrix
 mat4 view; // view matrix
 float znear = 0.001, zfar = 1000;
-float width = 512, height = 512;
-float aspect = width / height;
+// for performance with glReadPixels these should be powers of 2!
+float width = 1200, height = 800;
 
 void setProjectionAndViewUniforms(GLuint progId) {
   GLuint projLoc = glGetUniformLocation(progId, "projection");
@@ -160,13 +160,36 @@ mat4 faceView(mat4 m) {
   return m;
 }
 
+#define PBO
+
+/* const int shadeWidth = 256, shadeHeight = 256; */
+const int shadeWidth = 256, shadeHeight = 256;
+
 #ifdef PBO
 const int numPbos = 64;
 GLuint pboBufs[numPbos];
-GLbyte sink[1000 * 1000 * 4];
+GLbyte sink[shadeWidth * shadeHeight * 4];
+
+void inline mapPixelRead(int pboBuf, int metaball) {
+      glBindBuffer(GL_PIXEL_PACK_BUFFER, pboBufs[pboBuf]);
+      GLubyte *src = (GLubyte *)glMapBufferRange(GL_PIXEL_PACK_BUFFER, 0, 4 * sizeof(GLubyte),
+          GL_MAP_READ_BIT);
+      vec4 pixel = vec4(src[0], src[1], src[2], src[3]) / vec4(255.f);
+      glUnmapBuffer(GL_PIXEL_PACK_BUFFER);
+
+      glBindBuffer(GL_PIXEL_PACK_BUFFER, 0);
+
+      // Multiply the pixel value by the sunlight color.
+      pixel *= envColors[envColorIdx][0];
+
+      // Store the color for the previous set of pixels
+      metaballs[metaball].col = pixel;
+}
+
 #endif
 
 void shadeClouds() {
+
   glDisable(GL_DEPTH_TEST);
   // shaderOutput * 0 + buffer * shader alpha
   glBlendFunc(GL_ZERO, GL_SRC_ALPHA);
@@ -183,7 +206,7 @@ void shadeClouds() {
   GLuint modelLoc = glGetUniformLocation(bbProg, "model");
   glUniform1i(glGetUniformLocation(bbProg, "debug"), 0);
 
-  /* glViewport(0, 0, shadeWidth, shadeHeight); */
+  glViewport(0, 0, shadeWidth, shadeHeight);
 
 #ifdef PBO
   int pboIdx = 0;
@@ -192,7 +215,7 @@ void shadeClouds() {
   auto begin_time = std::chrono::system_clock::now();
   size_t i = 0;
   for (auto &k : metaballs) {
-    fprintf(stderr, "\rShading metaball %lu/%lu...", i++, metaballs.size());
+    /* fprintf(stderr, "\rShading metaball %lu/%lu...", i, metaballs.size()); */
     // place the billboard at the center of k
     mat4 model = translate(mat4(1), k.pos);
 
@@ -225,48 +248,55 @@ void shadeClouds() {
     // 3. Multiply by (width * height)
     ivec2 screenPos =
         ((vec2(proj * view * model * vec4(0, 0, 0, 1)) + vec2(1)) / vec2(2)) *
-        vec2(width, height);
-    vec4 pixel;
+        vec2(shadeWidth, shadeHeight);
 
 #ifndef PBO
+    vec4 pixel;
     // TODO: This is a huge bottleneck
     glReadPixels(screenPos.x, screenPos.y, 1, 1, GL_RGBA, GL_FLOAT, &pixel);
+
+    // Multiply the pixel value by the sunlight color.
+    pixel *= envColors[envColorIdx][0];
+
+    // Store the color into an array C[k] as the color of the billboard.
+    k.col = pixel;
 #else
 
     glBindBuffer(GL_PIXEL_PACK_BUFFER, pboBufs[pboIdx]);
-    glReadPixels(0, 0, width, height, GL_RGBA, GL_UNSIGNED_BYTE,
+    glReadPixels(screenPos.x, screenPos.y, 1, 1, GL_RGBA, GL_UNSIGNED_BYTE,
                  NULL);
+    // TODO: use this
+    /* glReadPixels(0, 0, shadeWidth, shadeHeight, GL_RGBA, GL_UNSIGNED_BYTE, */
+    /*              NULL); */
 
-    glBindBuffer(GL_PIXEL_PACK_BUFFER, pboBufs[(pboIdx + 1) % numPbos]);
-    GLuint offset = screenPos.x * 4 + screenPos.y * (int)width * 4;
-    ivec4 *src = (ivec4 *)glMapBufferRange(GL_PIXEL_PACK_BUFFER, offset, 4, GL_MAP_READ_BIT);
-    checkError();
-    if (src) {
-      ivec4 t = src[(int)screenPos.x + (int)screenPos.y * (int)width];
-      pixel = vec4(t) / vec4(255.f);
+    int nextPbo = (pboIdx + 1) % numPbos;
+    if (i >= numPbos - 1) {
+      // start mapping the read values back
+      mapPixelRead(nextPbo, i - numPbos + 1);
     }
-    glUnmapBuffer(GL_PIXEL_PACK_BUFFER);
-
-    pboIdx = (pboIdx + 1) % numPbos;
-
-    glBindBuffer(GL_PIXEL_PACK_BUFFER, 0);
+    pboIdx = nextPbo;
 #endif
     
-    // Multiply the pixel value by the sunlight color.
-    pixel *= envColors[envColorIdx][0];
+    i++;
+  }
+  /* fprintf(stderr, "\n"); */
 
-    // Store the color into an array C[k] as the color of the billboard.
-    k.col = pixel;
+#ifdef PBO
+  // sink remaining reads
+  for (int i = 0; i < numPbos; i++) {
+    mapPixelRead(i, metaballs.size() - numPbos + i);
   }
-  fprintf(stderr, "\n");
+#endif
 
   auto elapsed = std::chrono::system_clock::now() - begin_time;
-  double elapsed_seconds = std::chrono::duration_cast<std::chrono::duration<double> >(elapsed).count();
-  fprintf(stderr, "time with pbo: %f\n", elapsed_seconds);
+  double elapsed_seconds =
+      std::chrono::duration_cast<std::chrono::duration<double>>(elapsed)
+          .count();
+  fprintf(stderr, "Time taken to shade: %fs\n", elapsed_seconds);
 
   saveFBO();
   checkError();
-  /* glViewport(0, 0, width, height); */
+  glViewport(0, 0, width, height);
 }
 
 void renderObject() {
@@ -274,7 +304,7 @@ void renderObject() {
   // render the sun
   glUseProgram(sunProg);
   mat4 model = translate(mat4(1), sunPos);
-  /* model = lookAt(sunPos, sunPos + sunDir, {0, 1, 0}) * model; */
+  model = lookAt(sunPos, sunPos + sunDir, {0, 1, 0}) * model;
   model = translate(scale(translate(model, -sunPos), vec3(0.3)), sunPos);
   glUniformMatrix4fv(glGetUniformLocation(sunProg, "model"), 1, GL_FALSE,
                      glm::value_ptr(model));
@@ -341,7 +371,7 @@ void renderClouds() {
         debugVal = cs.p_ext[k.coords.x][k.coords.y][k.coords.z] / P_EXT;
       glUniform1f(glGetUniformLocation(bbProg, "debugVal"), debugVal);
       glBlendFunc(GL_SRC_ALPHA, GL_ONE_MINUS_SRC_ALPHA);
-      model = scale(model, vec3(0.2));
+      model = scale(model, vec3(0.1));
       glUniformMatrix4fv(modelLoc, 1, GL_FALSE, glm::value_ptr(model));
     }
 
@@ -350,9 +380,10 @@ void renderClouds() {
   }
 }
 
-bool needsReshading = true;
+bool curBeenShaded = false;
+
 void display() {
-  if (needsReshading) {
+  if (!curBeenShaded && (curMode == render || curMode == debugColor)) {
     // TODO: find a way to make sure there's no clipping
     view = glm::lookAt(sunPos + sunDir * vec3(100.f), sunPos, {0, 0, 1});
     // TODO: calculate bounds so everything is covered
@@ -363,10 +394,11 @@ void display() {
     glClearColor(1, 1, 1, 1);
     glClear(GL_COLOR_BUFFER_BIT);
     shadeClouds();
-    needsReshading = false;
+    curBeenShaded = true;
   }
 
   view = glm::lookAt(camPos, viewPos, {0, 1, 0});
+  const float aspect = width / height;
   proj = glm::perspective(45.f, aspect, znear, zfar);
   glUseProgram(sunProg);
   setProjectionAndViewUniforms(sunProg);
@@ -385,7 +417,6 @@ void display() {
 
 bool needsRedisplay = false;
 void timer(int _) {
-  /* calculateMetaballs(); */
   if (needsRedisplay) {
     glutPostRedisplay();
   }
@@ -397,10 +428,9 @@ void keyboard(unsigned char key, int x, int y) {
   if (key == ' ') {
     calculateMetaballs();
     needsRedisplay = true;
-    needsReshading = curMode == render;
+    curBeenShaded = false;
   }
   if (key == '0') {
-    needsReshading = curMode != render;
     curMode = render;
     needsRedisplay = true;
   }
@@ -423,7 +453,7 @@ void keyboard(unsigned char key, int x, int y) {
   if (key == 's') {
     envColorIdx = (envColorIdx + 1) % envColors.size();
     needsRedisplay = true;
-    needsReshading = true;
+    curBeenShaded = false;
   }
 }
 
@@ -453,6 +483,11 @@ void passiveMotion(int x, int y) {
   prevMouseY = y;
 }
 
+void reshape(int w, int h) {
+  width = w;
+  height = h;
+}
+
 int main(int argc, char **argv) {
   glutInit(&argc, argv);
   glutInitDisplayMode(GLUT_DEPTH | GLUT_DOUBLE | GLUT_RGBA | GLUT_ALPHA |
@@ -460,6 +495,11 @@ int main(int argc, char **argv) {
   glutInitWindowSize(width, height);
   glutCreateWindow("Clouds");
   glutDisplayFunc(display);
+  glutReshapeFunc(reshape);
+  glutKeyboardFunc(keyboard);
+  glutMotionFunc(motion);
+  glutPassiveMotionFunc(passiveMotion);
+  glutTimerFunc(16, timer, 0);
 
   glewInit();
 
@@ -497,17 +537,12 @@ int main(int argc, char **argv) {
   initClouds(&cs);
   calculateMetaballs();
 
-  glutKeyboardFunc(keyboard);
-  glutMotionFunc(motion);
-  glutPassiveMotionFunc(passiveMotion);
-  glutTimerFunc(16, timer, 0);
-
 #ifdef PBO
   // setup PBOs for buffering readPixels
   glGenBuffers(numPbos, pboBufs);
   for (int i = 0; i < numPbos; i++) {
     glBindBuffer(GL_PIXEL_PACK_BUFFER, pboBufs[i]);
-    glBufferData(GL_PIXEL_PACK_BUFFER, width * height * 4, NULL,
+    glBufferData(GL_PIXEL_PACK_BUFFER, shadeWidth * shadeHeight * 4, NULL,
                  GL_DYNAMIC_READ);
   }
   glBindBuffer(GL_PIXEL_PACK_BUFFER, 0);
index 87281f994995011a81baa14dd1a9c621e3431b18..daa756a2df82d4bb0549bfac628a3ca0b1cba8a7 100644 (file)
--- a/debug.hpp
+++ b/debug.hpp
@@ -4,6 +4,8 @@
 #include <cstring>
 #include <glm/glm.hpp>
 
+#pragma clang diagnostic ignored "-Wdeprecated-declarations"
+
 inline void dump(glm::mat4 m) {
   for (int i = 0; i < 4; i++)
     fprintf(stderr, "%f,%f,%f,%f\n", m[i][0], m[i][1], m[i][2], m[i][3]);
index 33b58de92c21f03d3e39a7310e33ad20a386dcc7..2268650cc4250f13be0386a6e0b1ca33b264aa81 100644 (file)
@@ -1,4 +1,4 @@
-#include "simulation.h"
+#include "simulation.hpp"
 #include <cstdlib>
 #include <glm/glm.hpp>
 
@@ -19,7 +19,7 @@ void initClouds(Clouds *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++) {
-        cs->act[i][j][k] = randf() < 0.01;
+        cs->act[i][j][k] = randf() < 0.005;
         cs->cld[i][j][k] = false;
         cs->hum[i][j][k] = randf() < 0.01;
         cs->p_ext[i][j][k] = 0.f;
@@ -29,17 +29,8 @@ void initClouds(Clouds *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++) {
-        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);
+  for (int k = 0; k < CLOUD_DIM_Y; k++)
+    cs->vy[k] = floor(randf() * 2);
 
   // generate ellipsoids of probability
   const int numEllipsoids = CLOUD_DIM_X * CLOUD_DIM_Y * CLOUD_DIM_Z * 0.001;
@@ -48,12 +39,15 @@ void initClouds(Clouds *cs) {
     float delta = maxSize - minSize;
     int width = minSize + randf() * delta, height = minSize + randf() * delta, depth = minSize + randf() * delta;
     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);
 
+    const float maxDist = glm::distance(glm::vec3(0), glm::vec3(0.5,0.5,0.5));
     for (int i = x; i < x + width; i++) {
       for (int j = y; j < y + height; j++) {
         for (int k = z; k < z + depth; k++) {
-          float dist = glm::distance(glm::vec3(i,j,k), center) / maxSize;
+          // between [0,0,0] and [1,1,1]
+          glm::vec3 uniformPos = glm::vec3(i - x, j - y, k - z) / glm::vec3(width, height, depth);
+
+          float dist = glm::distance(uniformPos, glm::vec3(0.5,0.5,0.5)) / maxDist;
           set(cs->p_ext, i, j, k, P_EXT * dist);
           set(cs->p_hum, i, j, k, P_HUM * (1.f - dist));
           set(cs->p_act, i, j, k, P_ACT * (1.f - dist));
@@ -79,8 +73,11 @@ inline bool f_act(Clouds *cs, int i, int j, int k) {
     || get(cs->act, i, j - 2, k) || get(cs->act, i , j + 2, k) || get(cs->act, i, j, k - 2);
 }
 
+// scratch for updates on the heap
+Clouds ncs;
+
 void growth(Clouds *cs) {
-  Clouds ncs = *cs;
+  ncs = *cs;
 
   for (int i = 0; i < CLOUD_DIM_X; i++) {
     for (int j = 0; j < CLOUD_DIM_Y; j++) {
@@ -96,7 +93,7 @@ void growth(Clouds *cs) {
 }
 
 void extinction(Clouds *cs) {
-  Clouds ncs = *cs;
+  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++) {
@@ -110,12 +107,12 @@ void extinction(Clouds *cs) {
 }
 
 void advection(Clouds *cs) {
-  Clouds ncs = *cs;
+  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];
+        int v = cs->vy[j];
         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;
similarity index 87%
rename from simulation.h
rename to simulation.hpp
index 1985ee3662fb6af8eec728d2a21d48ba6fb9dcb7..2a81c873c34c57262369093ff98bcb9b0b7043fb 100644 (file)
@@ -18,8 +18,10 @@ struct Clouds {
   float p_ext[CLOUD_DIM_X][CLOUD_DIM_Y][CLOUD_DIM_Z];
   float p_hum[CLOUD_DIM_X][CLOUD_DIM_Y][CLOUD_DIM_Z];
   float p_act[CLOUD_DIM_X][CLOUD_DIM_Y][CLOUD_DIM_Z];
-  /** Wind velocity */
-  float vz[CLOUD_DIM_Z];
+  /** Wind velocity 
+   * Called vz in the paper, but here our z axis is y
+   */
+  float vy[CLOUD_DIM_Y];
   /** continuous distribution */
   float q[CLOUD_DIM_X][CLOUD_DIM_Y][CLOUD_DIM_Z];
 };