› › ›

Reaction-Diffusion using GLSL works different from normal

Lord_of_the_GalaxyLord_of_the_Galaxy
edited May 2017 in GLSL / Shaders

I managed to get Reaction Diffusion Algorithm to run on the GPU, but now I'm having problems. The code runs and all, but the effect is much different from the effect I get with Reaction Diffusion Algorithm on the CPU. Why?

Code for shader based Reaction-Diffusion -

PDE file -

PShader algorithm, pass;
PGraphics pg;

public static final float DA = 1.0;
public static final float DB = 0.5;
public static final float feed = 0.055;
public static final float kill = 0.062;
public static final float dt = 1.0;

public static final float min = 0.48;
public static final float max = 0.52;

public static final int N_PASS = 10;

boolean first = true;

void setup() {
  size(500, 500, P2D);

  algorithm = loadShader("algorithm.frag");
  pass = loadShader("pass.frag");

  pg = createGraphics(width, height, P2D);
  pg.noSmooth();

  algorithm.set("resolution", float(pg.width), float(pg.height));
  algorithm.set("DA", DA);
  algorithm.set("DB", DB);
  algorithm.set("feed", feed);
  algorithm.set("kill", kill);
  algorithm.set("dt", dt);
  algorithm.set("min", min);
  algorithm.set("max", max);
}
void draw() {
  algorithm.set("first", first);
  first = false;
  for(int i = 0; i < N_PASS; i++){
    pg.beginDraw();
    pg.background(0);
    pg.shader(algorithm);
    pg.rect(0, 0, pg.width, pg.height);
    pg.endDraw();  
  }
  shader(pass);
  image(pg, 0, 0);
  if(frameCount%30 == 0)println(frameRate);//forget this
}

algortihm.frag -

#ifdef GL_ES
precision mediump float;
#endif

uniform vec2 resolution;
uniform sampler2D ppixels;

uniform float DA;
uniform float DB;
uniform float feed;
uniform float kill;
uniform float dt;
uniform float min;
uniform float max;

uniform bool first;

void main (void) {
    vec2 position = ( gl_FragCoord.xy / resolution.xy );
    vec2 pixel = 1.0/resolution;

    if(first){
        if(position.x > min && position.x < max && position.y > min && position.y < max){
            gl_FragColor = vec4(1.0, 1.0, 0.0, 1.0);
        }else{
            gl_FragColor = vec4(1.0, 0.0, 0.0, 1.0);
        }
    }else{
        float A = texture(ppixels, position).r;

        float sumA = -A;

        sumA += 0.05*texture(ppixels, position + pixel * vec2(-1.0, -1.0)).r;
        sumA += 0.05*texture(ppixels, position + pixel * vec2(-1.0, 1.0)).r;
        sumA += 0.05*texture(ppixels, position + pixel * vec2(1.0, -1.0)).r;
        sumA += 0.05*texture(ppixels, position + pixel * vec2(1.0, 1.0)).r;
        sumA += 0.2*texture(ppixels, position + pixel * vec2(1.0, 0.0)).r;
        sumA += 0.2*texture(ppixels, position + pixel * vec2(-1.0, 0.0)).r;
        sumA += 0.2*texture(ppixels, position + pixel * vec2(0.0, -1.0)).r;
        sumA += 0.2*texture(ppixels, position + pixel * vec2(0.0, 1.0)).r;


        float B = texture(ppixels, position).g;
        float sumB = -B;
        sumB += 0.05*texture(ppixels, position + pixel * vec2(-1.0, -1.0)).g;
        sumB += 0.05*texture(ppixels, position + pixel * vec2(-1.0, 1.0)).g;
        sumB += 0.05*texture(ppixels, position + pixel * vec2(1.0, -1.0)).g;
        sumB += 0.05*texture(ppixels, position + pixel * vec2(1.0, 1.0)).g;
        sumB += 0.2*texture(ppixels, position + pixel * vec2(1.0, 0.0)).g;
        sumB += 0.2*texture(ppixels, position + pixel * vec2(-1.0, 0.0)).g;
        sumB += 0.2*texture(ppixels, position + pixel * vec2(0.0, -1.0)).g;
        sumB += 0.2*texture(ppixels, position + pixel * vec2(0.0, 1.0)).g;

        float nA = A + ((DA*sumA) - (A*B)*B + (feed*(1.0 - A)))*dt;
        float nB = B + ((DB*sumB) + (A*B)*B - ((kill + feed)*B))*dt;


        gl_FragColor = vec4(clamp(nA, 0.0, 1.0), clamp(nB, 0.0, 1.0), 0.0, 1.0);
    }
}

pass.frag -

uniform sampler2D texture;
varying vec4 vertTexCoord;

void main(){
    vec4 inColor = texture2D(texture, vertTexCoord.st);
    gl_FragColor = vec4(1.0 - inColor.g, 1.0 - inColor.g, 1.0 - inColor.g, 1.0);
}
Tagged:

Best Answers

Answers

  • Lord_of_the_GalaxyLord_of_the_Galaxy

    And the other code (warning - it is really slow):

    float[][][] grid;
float[][][] next;

float cW = 0.05;//weight for diagonal place
float sW = 0.2;//weight for adjacent place

float dA = 1.0;
float dB = 0.5;

float feed = 0.055;
float kill = 0.062;

float dt = 1;

void setup(){
  size(250, 250);
  grid = new float[width][height][2];
  next = new float[width][height][2];
  for(int x = 0; x < width; x++){
    for(int y = 0; y < height; y++){
      grid[x][y][0] = 1.0;
      grid[x][y][1] = 0.0;
      next[x][y][0] = 1.0;
      next[x][y][1] = 0.0;
    }
  }
  for(int x = width/2 - 5; x < width/2 + 5; x++){
    for(int y = height/2 - 5; y < height/2 + 5; y++){
      grid[x][y][0] = 0.0;
      grid[x][y][1] = 1.0;
      next[x][y][0] = 0.0;
      next[x][y][1] = 1.0;
    }
  }
}

void draw(){
    calc();
    loadPixels();
    for(int y = 0; y < height; y++){
      for(int x = 0; x < width; x++){
        int i = x + y*width;
        float r = 255 - grid[x][y][0]*255;
        float g = 255 - grid[x][y][1]*255;
        float b = 0;
        pixels[i] = color(g, g, g);
      }
    }
    updatePixels();
}

void calc(){
  for(int x = 1; x < width-1; x++){
    for(int y = 1; y < height-1; y++){
      float a = grid[x][y][0];
      float b = grid[x][y][1];
      next[x][y][0] = a + (dA*laplace(x, y, 0) - a*b*b + feed*(1 - a))*dt;
      next[x][y][1] = b + (dB*laplace(x, y, 1) + a*b*b - (kill + feed)*b)*dt;
    }
  }
  swap();
}

float laplace(int x, int y, int i){
  float r = -grid[x][y][i];
  r += sW * grid[x+1][y][i];
  r += sW * grid[x-1][y][i];
  r += sW * grid[x][y+1][i];
  r += sW * grid[x][y-1][i];
  r += cW * grid[x+1][y+1][i];
  r += cW * grid[x+1][y-1][i];
  r += cW * grid[x-1][y+1][i];
  r += cW * grid[x-1][y-1][i];
  return r;
}

void swap(){
  float[][][] temp = grid;
  grid = next;
  next = temp;
}
  • nabrnabr
    edited May 2017

    looks interesting! are you sure you are multiplying the same way as in glsl shader? the growing speed is just slower ... did also some code cleanup

    float[][][] grid;
float[][][] next;

float cW = 0.05;//weight for diagonal place
float sW = 0.2;//weight for adjacent place
float dA = 1.0;
float dB = 0.5;
float feed = 0.055;
float kill = 0.062;
float dt = 1;

void setup(){
  size(250, 250);
  grid = new float[width][height][2];
  next = new float[width][height][2];
  for(int x = 0; x < width; x++){
    for(int y = 0; y < height; y++){
      grid[x][y][0] = 1.0;
      grid[x][y][1] = 0.0;
      next[x][y][0] = 1.0;
      next[x][y][1] = 0.0; 
      if(( x >width/2-5) &&(x< width/2+5)
       &&( y >height/2-5)&&(y< height/2+5)){
         grid[x][y][0] = 0.0;
         grid[x][y][1] = 1.0;
         next[x][y][0] = 0.0;
         next[x][y][1] = 1.0;
      }
    }
  }
}
void draw(){
   for(int x = 1; x < width-1; x++){
    for(int y = 1; y < height-1; y++){
      float a = grid[x][y][0];
      float b = grid[x][y][1];
      next[x][y][0] = a + (dA*laplace(x, y, 0) - a*b*b + feed*(1 - a))*dt;
      next[x][y][1] = b + (dB*laplace(x, y, 1) + a*b*b - (kill + feed)*b)*dt;
    }
  }
 swap();
  loadPixels();
    for(int y = 0; y < height; y++){
      for(int x = 0; x < width; x++){
        int i = x + y*width;
        float r = 255 - grid[x][y][0]*255;
        float g = 255 - grid[x][y][1]*255;
        float b = 0;
        pixels[i] = color(g, g,g);
      }
    }; updatePixels();
}
float laplace(int x, int y, int i){
  float r = -grid[x][y][i];
  r += sW * grid[x+1][y][i];
  r += sW * grid[x-1][y][i];
  r += sW * grid[x][y+1][i];
  r += sW * grid[x][y-1][i];
  r += cW * grid[x+1][y+1][i];
  r += cW * grid[x+1][y-1][i];
  r += cW * grid[x-1][y+1][i];
  r += cW * grid[x-1][y-1][i];
  return r;
}
void swap(){
  float[][][] temp = grid;
  grid = next;
  next = temp;
}
  • Lord_of_the_GalaxyLord_of_the_Galaxy

    @nabr The growing speed is slower for a reason - in the GLSL version, I'm running the shader ten times each frame, whereas in the CPU version, I'm only doing it once per frame. If I try to do it ten times per frame, the frame rate drops to 10-12 fps, however, the GLSL version remains at 60 fps even with 10 times.

    But if you let the CPU version run for some time, it will create a nice pattern even with the exact same variable values. Why??

  • nabrnabr
    edited May 2017

    @Lord_of_the_Galaxy (if their is no bug in your code) generally speaking your GPU is faster then the CPU.

  • Lord_of_the_GalaxyLord_of_the_Galaxy
    edited May 2017

    @nabr Indeed, but that is not the point. The point is that the pattern doesn't form with the GPU version.

  • nabrnabr
    edited May 2017

    @Lord_of_the_Galaxy

    First, sorry, i kind of chatting here, back to professionalism

    can you share picures how RDA Pattern looks on your system?

    the CPU version you provided (copy +paste from your post) works fine on my side, but the GPU stay in the same state without growing. so it's not working. interessting.

    CPU version lap

  • Lord_of_the_GalaxyLord_of_the_Galaxy
    edited May 2017

    @nabr Exactly. That's exactly what's happening, the GPU version gets stuck, and the CPU version gives the same pattern as you posted. Strange, right?
    And yes, we can chat about this after it is solved.

  • Lord_of_the_GalaxyLord_of_the_Galaxy

    So no one knows why this is happening? I thought it must be an error on my side, but I think it may be some sort of bug or something. Should I submit a bug report?

  • nabrnabr
    edited May 2017

    @Lord_of_the_Galaxy

    Sorry i can't help you. Becourse i have to learn more about Reaction-Diffusion and i currently busy with my oft stuff.

    What i can say: If i start your shader it stops, - my idea is, it stops becourse the algorithm done his job to quickly, you have to "feed" it with more data.

    for(int i = 0; i < N_PASS; i++){
    pg.beginDraw();
    pg.background(0);
    pg.shader(algorithm);
    pg.rect(0, 0, pg.width, pg.height);
    pg.endDraw();  
  }

    The funtion call chain is:
    =>setup
    =>compile shader
    =>run draw loop

    if you run a for loop in the draw over a shader method ,
    "pg.shader(algorithm);" The shader does not recompile. I think, in this case, - you just say, make a copy N_PASS x times.     So you place the same image over and over again.

    shader(pass);
image(pg, 0, 0);

    And for the algorithm to work you have to say, place the image with more and more data in the "second" pass.

    Try to get more data with shader.set("uniform", data) with the set function your shader forced to recompile.

    Edit. @Lord_of_the_Galaxy

    here is an working example: https://forum.processing.org/two/discussion/9100/unknow-error-for-simple-fragment-shader-reaction-diffusion

    Edit 2:
    @Lord_of_Galaxy vec2 position = ( gl_FragCoord.xy / resolution.xy ); will give you what ?
    .5 ? with min is .48 and max .52. not much space for Reaction–diffusion

    you have to chage sumA += 0.05*texture(ppixels, position + pixel * vec2(-1.0, -1.0)).r;

  • nabrnabr
    edited May 2017

    your loop N_PASS x times works perfectly, nice!

    not sure ppixels has probably pre caltucated texture coords

    This is not a completely fix of Reaction-Diffusion, just something to get started

    algortihm.frag -

    #ifdef GL_ES
precision mediump float;
#endif

uniform vec2 resolution;
uniform sampler2D ppixels;
uniform float DA;
uniform float DB;
uniform float feed;
uniform float kill;
uniform float dt;
uniform float min;
uniform float max;

uniform bool first;
void main (void) {
    vec2 position = ( gl_FragCoord.xy / resolution.xy );
    vec2 pixel = 1.0/resolution;

    if(first){
        if(position.x > min && position.x < max && position.y > min && position.y < max){
            gl_FragColor = vec4(1.0, 1.0, 0.0, 1.0);
        }else{
            gl_FragColor = vec4(1.0, 0.0, 0.0, 1.0);
        }
    }else{
        float A = texture(ppixels, position).r;

        float sumA = -A;

        sumA += 0.05*texture(ppixels, position + pixel * vec2( .5*-1.0, .5*-1.0)).r;
        sumA += 0.05*texture(ppixels, position + pixel * vec2(.5*-1.0, .5*1.0)).r;
        sumA += 0.05*texture(ppixels, position + pixel * vec2(.5*1.0, .5*-1.0)).r;
        sumA += 0.05*texture(ppixels, position + pixel * vec2(.5*1.0, .5*1.0)).r;
        sumA += 0.2*texture(ppixels, position + pixel * vec2(1.0, 0.0)).r;
        sumA += 0.2*texture(ppixels, position + pixel * vec2(-1.0, 0.0)).r;
        sumA += 0.2*texture(ppixels, position + pixel * vec2(0.0, -1.0)).r;
        sumA += 0.2*texture(ppixels, position + pixel * vec2(0.0, 1.0)).r;


        float B = texture(ppixels, position).g;
        float sumB = -B;
        sumB += 0.05*texture(ppixels, position + pixel * vec2(.05*-1.0, .05*-1.0)).g;
        sumB += 0.05*texture(ppixels, position + pixel * vec2(.05*-1.0, .05*1.0)).g;
        sumB += 0.05*texture(ppixels, position + pixel * vec2(.05*1.0, .05*-1.0)).g;
        sumB += 0.05*texture(ppixels, position + pixel * vec2(.05*1.0, .05*1.0)).g;
        sumB += 0.2*texture(ppixels, position + pixel * vec2(1.0, 0.0)).g;
        sumB += 0.2*texture(ppixels, position + pixel * vec2(-1.0, 0.0)).g;
        sumB += 0.2*texture(ppixels, position + pixel * vec2(0.0, -1.0)).g;
        sumB += 0.2*texture(ppixels, position + pixel * vec2(0.0, 1.0)).g;

        float nA = A + ((DA*sumA) - (A*B)*B + (feed*(1.0 - A)))*dt;
        float nB = B + ((DB*sumB) + (A*B)*B - ((kill + feed)*B))*dt;


        gl_FragColor = vec4(clamp(nA, 0.0, 1.0), clamp(nB, 0.0, 1.0), 0.0, 1.0);
    }
}
  • nabrnabr
    edited May 2017 Answer ✓

    @Lord_of_the_Galaxy:

    Still unsolved ? :((


    It's working for me, - but i have i idea of Reaction-Diffusion finally you have to fix the output.

    The Problem you are runing into:

    pass.frag -

    algortihm.frag -

    • sumA += 0.05*texture(ppixels, position + pixel * vec2(-1.0, -1.0)).r;
      etc ...

    When you run your original scetch, the algorithm stops becourse their is nothing more to do.

    .pde -

    • public static final float min = 0.48;
    • public static final float max =0.52;

    algortihm.frag -

    • if(position.x > min && position.x < max && position.y > min && position.y < max)

    you feed the algorithm e.g float min has only few (pixels) of data to move in the x direction 'til it get's to -1.

    reactiondiff

  • Lord_of_the_GalaxyLord_of_the_Galaxy

    @nabr I see. I thought due to the parallel nature of GLSL, it may just skip those problem pixels and go to next one. I try setting it this weekend and see if it works. Thanks!

  • Lord_of_the_GalaxyLord_of_the_Galaxy

    @nabr Unfortunately, it still won't work.

  • Lord_of_the_GalaxyLord_of_the_Galaxy

    I think, it is happening because of some sort of higher float inaccuracies in GLSL than in JVM (or whatever).

  • Lord_of_the_GalaxyLord_of_the_Galaxy

    Perhaps @cansik knows?

  • nabrnabr
    edited May 2017

    @Lord_of_the_Galaxy

    we were positing as the same time,I think, it is happening because of some sort of higher float inaccuracies in GLSL than in JVM (or whatever).
    I think i has more to do with the default texcoords*0.5+0.5; in processing. Idk.

    PShader algorithm, pass;
PGraphics pg;

public static final float DA = 1.0;
public static final float DB = 0.5;
public static final float feed = 0.055;
public static final float kill = 0.062;
public static final float dt = 1.0;

public static final float min = 0.48;
public static final float max = 0.52;

public static final int N_PASS = 10;

boolean first = true;

void setup() {
  size(500, 500, P2D);

  algorithm = loadShader("https://"+"pastebin.com/raw/zVQtz4pB");
  pass = loadShader("https://"+"pastebin.com/raw/irTCrBYa");

  pg = createGraphics(width, height, P2D);
  pg.noSmooth();

  algorithm.set("resolution", float(pg.width), float(pg.height));
  algorithm.set("DA", DA);
  algorithm.set("DB", DB);
  algorithm.set("feed", feed);
  algorithm.set("kill", kill);
  algorithm.set("dt", dt);
  algorithm.set("min", min);
  algorithm.set("max", max);
}
void draw() {
  algorithm.set("first", first);
  first = false;
  for(int i = 0; i < N_PASS; i++){
    pg.beginDraw();
    pg.background(0);
    pg.shader(algorithm);
    pg.rect(0, 0, pg.width, pg.height);
    pg.endDraw();  
  }
  shader(pass);
  image(pg, 0, 0);
  if(frameCount%30 == 0)println(frameRate);//forget this
}

                   alt text

  • Lord_of_the_GalaxyLord_of_the_Galaxy

    Yes, I get that. But it's still different from the original.

  • Lord_of_the_GalaxyLord_of_the_Galaxy

    @nabr Actually, I think it's due to an error on your part - for sumA, you multiply by 0.5, for sumB you multiply by 0.05. That's why you get a pattern and I don't. But the pattern is still not correct, I'm afraid.

  • nabrnabr
    edited May 2017

    @Lord_of_the_Galaxy

    Complicated math in a predefined environment. Things are gonna break.

    Normaly if you running in a bug, a common practice can be to comment out, turn the lights one by an other off.

    I did just a small offset value each step of interation. Yes for SumB is 0.05 but here is the point, without knowledge how Reaction -Diffusion works, i don't know what the right offset values can be. You can also trow 0.5 in SumB it will look again slighly different.

    For Debug the Problem -

    Start with hard coded nummers whereever it is possible:

    vec2 position = ( gl_FragCoord.xy / vec2(500,500));
    vec2 pixel = vec2(0.002); //1./resolution;

    gl_FragCoord can be a bit a of challenge
    khronos.org/opengl/wiki/Compute_eye_space_from_window_space

    Now we will compute my offset nummer
    Step 1. (position+pixel) * 0.5-1.0 =~-1.5
    Step 2 =Step1+ (position+pixel)*0.5-1.0=???

    Than substract from your origanl version and you get an avenge value to work with Reaction-Diffusion

    This might be become time consuming. What i would do, i would turn everything off, that i don't really have control over. In this case Processing.

    P2D has predefinded: Viewport, OrthoCam, TextureCoords etc.

    Try to define your vertexShader by your self. You own textureCoords.

    Status

    what is so special about 10px ?

                  offset1

  • Lord_of_the_GalaxyLord_of_the_Galaxy

    what is so special about 10px?

    Nothing really, except that it works in the CPU version.

  • Lord_of_the_GalaxyLord_of_the_Galaxy

    So basically, my problem can be simplified into a simple question - how to get the pixels right around the current fragment pixel? One would think that my original code does that, but something there is breaking.

  • Lord_of_the_GalaxyLord_of_the_Galaxy

    @nabr I think I know the problem now - in the shader version, the R and G values are automatically clamped in between 0 and 1, whereas that is not done in the CPU version. And my removing clamp from it does nothing because it is automatic - nothing I can do about it.

    So my new question - how to allow color values beyond 1 and below 0 in the shader?

  • nabrnabr

    @Lord_of_the_Galaxy I currently busy. But when, i will return to this. Can take a while.

    clamp means interpolation between values, like map function in processing. https://thebookofshaders.com/05/

  • Lord_of_the_GalaxyLord_of_the_Galaxy

    @nabr

    y = clamp(x,0.0,1.0); // constrain x to lie between 0.0 and 1.0.

    So clamp is like constrain, and not map.

  • nabrnabr

    @Lord_of_the_Galaxy

    yes, sure, your are right. https://processing.org/reference/map_.html

  • nabrnabr

    @Lord_of_the_Galaxy

    this is a working example in processing. https://github.com/abigpotostew/reaction_diffusion/blob/master/reaction_diffusion_shader/data/grey_scott.glsl

    can you adapt some of his technique ?

  • nabrnabr
    edited May 2017

    @Lord_of_the_Galaxy

    did some code clean up

    alt text

  • Lord_of_the_GalaxyLord_of_the_Galaxy

    @nabr Thanks, I'll look at it evening/tomorrow.

  • nabrnabr
    edited May 2017

    @Lord_of_the_Galaxy Take your time.

    Note the gif is renderend using my Geforce Card the result may varying (if i using Intel it want work probably)

  • nabrnabr
    edited May 2017

    Geforce v.s Intel

    same framecount

    alt text alt text

  • Lord_of_the_GalaxyLord_of_the_Galaxy
    edited May 2017

    @nabr How????????????????? This is strange.
    I think I'll shift to OpenCL (Aparapi probably). More stable.

  • nabrnabr
    edited May 2017 Answer ✓

    @Lord_of_the_Galaxy

    I'm not sure if it is a shaderbug - right now.

    I think, i happens becorse my Geforce is much more powerful.
    The NPASS Loop gets much faster unrolled, (the "time" is running faster ) the algorithm is very sensitive to small changes, the result looks totally different.

  • Lord_of_the_GalaxyLord_of_the_Galaxy

    @nabr And that's why I'll shift to OpenCL (Aparapi probably).

  • nabrnabr

    Sound great!

  • nabrnabr

    In case you already downloaded the files. - I updated the source file folder just to have the same conditions.

    Have a nice day. I realy need to work now. I spend already 2 hours with debugging
    Later.

  • Lord_of_the_GalaxyLord_of_the_Galaxy

    Thanks a lot for you effort!
    But like I said, I'll just use OpenCL, so sorry for not continuing your work.

  • nabrnabr
    edited May 2017

    @Lord_of_the_Galaxy

    No Problem. I did more source code optimization, just for fun. again with a slighly different result. i like the idea to map that thing to a mesh, and do some lighting, maybe i will return to it at a later point in time. I'm also sorry, i can't help you. Yes, learn OpenCL and get a programmer job!

    Best

  • Lord_of_the_GalaxyLord_of_the_Galaxy

    @nabr If I could really get a job simply by learning OpenCL (and programming itself), I would love that. But I'm probably not even old enough to qualify :)

  • nabrnabr
    edited May 2017 Answer ✓

    @T_D Interessting Offscreen Render Technique. You have to scroll up to the very first post, that is what matters. If you find any time for debug, any help would be much appreciated. The bug: it stops at 10 pixels. The question: why?

    @Lord_of_the_Galaxy: Visit his blog a Master of OpenCL

  • Lord_of_the_GalaxyLord_of_the_Galaxy

    @nabr Could you link to his blog please?

  • nabrnabr
    Answer ✓

    @Lord_of_the_Galaxy

    Hope you finde something useful http://thomasdiewald.com/blog/?p=2896

    Have a nice day.

  • Lord_of_the_GalaxyLord_of_the_Galaxy

    I didn't know PixelFlow uses OpenCL

  • T_DT_D
    edited June 2017 Answer ✓

    The main issue (but not the only one) here is that PGraphics textures use the internalformat RGBA8 which means you can't do a lot of multipass float operations because there is only 8 bit per channel ... just enough precission for RGB values (from 0 to 255).

    One solution is to use two 8bit-channels (2x8 = 16 bit) for encoding a float value. A RGBA8 texel can then hold 2 float values, encoded in the channels RG and BA, which is enough for the above example.

    The following demo uses some simple encoding, ... and fixes some of the other issues too. 

    resolution 800px /800px
75fps
100 passes / frame



    PShader shader_grayscott;
PShader shader_render;
PGraphics2D pg_src;
PGraphics2D pg_dst;

int pass = 0;

public void settings() {
  size(800, 800, P2D);
  smooth(0);
}

public void setup() {

  shader_grayscott = loadShader("grayscott2.frag");
  shader_render = loadShader("render2.frag");
  pg_src = createTexture(width, height);
  pg_dst = createTexture(width, height);
  
  // init
  pg_src.beginDraw();
  pg_src.background(0xFFFF0000);
  pg_src.fill(0x0000FFFF);
  pg_src.noStroke();
  pg_src.rectMode(CENTER);
  pg_src.rect(width/2, height/2, 20, 20);
  pg_src.endDraw();

  frameRate(1000);
}

public PGraphics2D createTexture(int w, int h){
  PGraphics2D pg = (PGraphics2D) createGraphics(w, h, P2D);
  pg.smooth(0);
  pg.beginDraw();
  pg.textureSampling(2);
  pg.blendMode(REPLACE);
  pg.clear();
  pg.noStroke();
  pg.endDraw();
  return pg;
}

public void swap(){
  PGraphics2D pg_tmp = pg_src;
  pg_src = pg_dst;
  pg_dst = pg_tmp;
}

public void reactionDiffusionPass(){
  pg_dst.beginDraw();
  shader_grayscott.set("dA"    , 1.0f  );
  shader_grayscott.set("dB"    , 0.5f  );
  shader_grayscott.set("feed"  , 0.055f);
  shader_grayscott.set("kill"  , 0.062f);
  shader_grayscott.set("dt"    , 1f    );
  shader_grayscott.set("wh_rcp", 1f/width, 1f/height);
  shader_grayscott.set("tex"   , pg_src);
  pg_dst.shader(shader_grayscott);
  pg_dst.rectMode(CORNER);
  pg_dst.rect(0, 0, width, height);
  pg_dst.endDraw(); 
  swap();
  pass++;
}

public void draw() {
  
  // multipass rendering, ping-pong 
  for(int i = 0; i < 100; i++){
    reactionDiffusionPass();
  }

  // display result
  shader_render.set("wh_rcp", 1f/width, 1f/height);
  shader_render.set("tex"   , pg_src);
  shader(shader_render);
  rect(0, 0, width, height);

  String txt_fps = String.format(getClass().getSimpleName()+ "   [size %d/%d]   [frame %d]   [fps %6.2f]", width, height, pass, frameRate);
  surface.setTitle(txt_fps);
}

public void keyReleased(){
  if(key == 's') saveFrame("grayscott.jpg");
}

    shader: grayscott2.frag

    
// file: grayscott2.frag
// author: diewald

#version 150

out vec4 fragColor;

uniform vec2 wh_rcp;
uniform sampler2D tex;

uniform float dA   = 1.0;
uniform float dB   = 0.5;
uniform float feed = 0.055;
uniform float kill = 0.062;
uniform float dt   = 1.0;

vec2 decode(vec4 dataf){
  ivec4 datai = ivec4(dataf * 255.0);
  float rg = (datai.r << 8 | datai.g) / 65535.0;
  float ba = (datai.b << 8 | datai.a) / 65535.0;
  return vec2(rg, ba);
}

vec4 encode(vec2 dataf){
  ivec2 datai = ivec2(dataf * 65535.0);
  int r = (datai.r >> 8) & 0xFF;
  int g = (datai.r     ) & 0xFF;
  int b = (datai.g >> 8) & 0xFF;
  int a = (datai.g     ) & 0xFF;
  return vec4(r,g,b,a) / 255.0;
}

void main () {

  vec2 posn = gl_FragCoord.xy * wh_rcp;
  
  vec2 val = decode(texture(tex, posn));
  
  vec2 laplace = -val;
  laplace += decode(textureOffset(tex, posn, ivec2(-1, 0))) * + 0.20;
  laplace += decode(textureOffset(tex, posn, ivec2(+1, 0))) * + 0.20;
  laplace += decode(textureOffset(tex, posn, ivec2( 0,-1))) * + 0.20;
  laplace += decode(textureOffset(tex, posn, ivec2( 0,+1))) * + 0.20;
  laplace += decode(textureOffset(tex, posn, ivec2(-1,-1))) * + 0.05;
  laplace += decode(textureOffset(tex, posn, ivec2(+1,-1))) * + 0.05;
  laplace += decode(textureOffset(tex, posn, ivec2(-1,+1))) * + 0.05;
  laplace += decode(textureOffset(tex, posn, ivec2(+1,+1))) * + 0.05;
  
  float nA = val.r + (dA * laplace.r - val.r * val.g * val.g + (feed * (1.0 - val.r))) * dt;
  float nB = val.g + (dB * laplace.g + val.r * val.g * val.g - ((kill + feed) * val.g)) * dt;

  fragColor = encode(clamp(vec2(nA, nB), vec2(0), vec2(1)));
}

    shader: render2.frag

    
// file: render2.frag
// author: diewald

#version 150

out vec4 fragColor;

uniform vec2 wh_rcp;
uniform sampler2D tex;


vec2 decode(vec4 dataf){
  ivec4 datai = ivec4(dataf * 255.0);
  float rg = (datai.r << 8 | datai.g) / 65535.0;
  float ba = (datai.b << 8 | datai.a) / 65535.0;
  return vec2(rg, ba);
}

vec4 encode(vec2 dataf){
  ivec2 datai = ivec2(dataf * 65535.0);
  int r = (datai.r >> 8) & 0xFF;
  int g = (datai.r     ) & 0xFF;
  int b = (datai.g >> 8) & 0xFF;
  int a = (datai.g     ) & 0xFF;
  return vec4(r,g,b,a) / 255.0;
}


void main () {
  vec2 val = decode(texture(tex, gl_FragCoord.xy * wh_rcp));
  fragColor = vec4(vec3(1.0 - val.g), 1);
}

    here is the output:
    grayscott2

  • nabrnabr

    impressed

  • Lord_of_the_GalaxyLord_of_the_Galaxy

    @T_D I should've realised that :) Thanks! I guess I've got a lot more to learn.

  • nabrnabr
    edited June 2017

    i need a gif for this

    (scaled to save loading time for .gif)

  • T_DT_D
    edited June 2017 Answer ✓

    I did a quick test to avoid float encoding/decoding using PixelFlow and custom texture formats GL_RG32F and GL_RG16F.
    It runs overall a bit faster now and heavy multipass rendering on float-texels is more straight forward.

    Notes:

    When GL_RG32F is used (32 bits per channel) we get the correct result, since we have enough precission. It runs about 85 fps (10fps more).

    when GL_RG16F, the outcome is noticeably different. There is a result, but it starts to freeze at some point, which indicates that 16 bit is clearly not enough. unless some better encoding handled by ourself. However, it is really fast, around 130 fps.

    the code is basically the same, the shaders got bit shorter.
    Download Sketch: ReactionDiffusion_PixelFlow.zip 

    
// Reaction Diffusion
// author: diewald

import com.jogamp.opengl.GL2;
import com.thomasdiewald.pixelflow.java.DwPixelFlow;
import com.thomasdiewald.pixelflow.java.dwgl.DwGLSLProgram;
import com.thomasdiewald.pixelflow.java.dwgl.DwGLTexture;
import com.thomasdiewald.pixelflow.java.imageprocessing.filter.DwFilter;

import processing.core.PApplet;
import processing.opengl.PGraphics2D;


DwGLSLProgram shader_grayscott;
DwGLSLProgram shader_render;

// multipass rendering texture
DwGLTexture.TexturePingPong tex_grayscott = new DwGLTexture.TexturePingPong();

// final render target for display
PGraphics2D tex_render;

DwPixelFlow context;
int pass = 0;

public void settings() {
  size(800, 800, P2D);
  smooth(0);
}

public void setup() {
  
  // pixelflow context
  context = new DwPixelFlow(this);

  // 1) 32 bit per channel
  tex_grayscott.resize(context, GL2.GL_RG32F, width, height, GL2.GL_RG, GL2.GL_FLOAT, GL2.GL_NEAREST, 2, 4);
  
  // 2) 16 bit per channel, lack of precision is obvious in the result, its fast though
  // tex_grayscott.resize(context, GL2.GL_RG16F, width, height, GL2.GL_RG, GL2.GL_FLOAT, GL2.GL_NEAREST, 2, 2);
  
  
  // glsl shader
  shader_grayscott = context.createShader("data/grayscott.frag");
  shader_render    = context.createShader("data/render.frag");
      
  // init
  tex_render = (PGraphics2D) createGraphics(width, height, P2D);
  tex_render.smooth(0);
  tex_render.beginDraw();
  tex_render.textureSampling(2);
  tex_render.blendMode(REPLACE);
  tex_render.clear();
  tex_render.noStroke();
  tex_render.background(0x00FF0000);
  tex_render.fill      (0x0000FF00);
  tex_render.noStroke();
  tex_render.rectMode(CENTER);
  tex_render.rect(width/2, height/2, 20, 20);
  tex_render.endDraw();

  // copy initial data to source texture
  DwFilter.get(context).copy.apply(tex_render, tex_grayscott.src);

  frameRate(1000);
}


public void reactionDiffusionPass(){
  context.beginDraw(tex_grayscott.dst);
  shader_grayscott.begin();
  shader_grayscott.uniform1f     ("dA"    , 1.0f  );
  shader_grayscott.uniform1f     ("dB"    , 0.5f  );
  shader_grayscott.uniform1f     ("feed"  , 0.055f);
  shader_grayscott.uniform1f     ("kill"  , 0.062f);
  shader_grayscott.uniform1f     ("dt"    , 1f    );
  shader_grayscott.uniform2f     ("wh_rcp", 1f/width, 1f/height);
  shader_grayscott.uniformTexture("tex"   , tex_grayscott.src);
  shader_grayscott.drawFullScreenQuad();
  shader_grayscott.end();
  context.endDraw("reactionDiffusionPass()"); 
  tex_grayscott.swap();
  pass++;
}


public void draw() {
  
  // multipass rendering, ping-pong 
  context.begin();
  for(int i = 0; i < 100; i++){
    reactionDiffusionPass();
  }

  // create display texture
  context.beginDraw(tex_render);
  shader_render.begin();
  shader_render.uniform2f     ("wh_rcp", 1f/width, 1f/height);
  shader_render.uniformTexture("tex"   , tex_grayscott.src);
  shader_render.drawFullScreenQuad();
  shader_render.end();
  context.endDraw("render()"); 
  context.end();
  

  // put it on the screen
  blendMode(REPLACE);
  image(tex_render, 0, 0);

  if(frameCount == 1000) saveFrame("data/grayscott.jpg");
  
  String txt_fps = String.format(getClass().getSimpleName()+ "   [size %d/%d]   [frame %d]   [fps %6.2f]", width, height, pass, frameRate);
  surface.setTitle(txt_fps);
}




    // file: grayscott.frag
// author: diewald

#version 150

out vec2 fragColor;

uniform vec2 wh_rcp;
uniform sampler2D tex;

uniform float dA   = 1.0;
uniform float dB   = 0.5;
uniform float feed = 0.055;
uniform float kill = 0.062;
uniform float dt   = 1.0;

void main () {

  vec2 posn = gl_FragCoord.xy * wh_rcp;
  
  vec4 val = texture(tex, posn);
  
  vec4 laplace = -val;
  laplace += textureOffset(tex, posn, ivec2(-1, 0)) * + 0.20;
  laplace += textureOffset(tex, posn, ivec2(+1, 0)) * + 0.20;
  laplace += textureOffset(tex, posn, ivec2( 0,-1)) * + 0.20;
  laplace += textureOffset(tex, posn, ivec2( 0,+1)) * + 0.20;
  laplace += textureOffset(tex, posn, ivec2(-1,-1)) * + 0.05;
  laplace += textureOffset(tex, posn, ivec2(+1,-1)) * + 0.05;
  laplace += textureOffset(tex, posn, ivec2(-1,+1)) * + 0.05;
  laplace += textureOffset(tex, posn, ivec2(+1,+1)) * + 0.05;
  
  float nA = val.r + (dA * laplace.r - val.r * val.g * val.g + (feed * (1.0 - val.r))) * dt;
  float nB = val.g + (dB * laplace.g + val.r * val.g * val.g - ((kill + feed) * val.g)) * dt;
  
  fragColor = clamp(vec2(nA, nB), vec2(0), vec2(1));
}



    // file: render.frag
// author: diewald

#version 150

out vec4 fragColor;

uniform vec2 wh_rcp;
uniform sampler2D tex;

void main () {
  vec2 val = texture(tex, gl_FragCoord.xy * wh_rcp).rg;
  fragColor = vec4(vec3(1.0 - val.g), 1);
}



    32 bit result



    16 bit result

  • T_DT_D

    I just figured, using the internalformat GL_RG16_SNORM instead of GL_RG16F gives roughly enough precision at around 120 fps.

    Using SNORM makes sense if the values are known to be in range of [-1, 1].

  • Lord_of_the_GalaxyLord_of_the_Galaxy

    @T_D Great Job!

    P.S. The imports on line 10-11 shouldn't be there. It's like its some sort of funny mixture of PDE and pure Java code.

  • BGADIIBGADII
    edited April 2018

    @T_D, your website is unavailable. It has been broken for about a month now. Also, does pixelflow work with android and shaders? Last question, when is your radiosity demo coming out? Haven't seen any updates in months and website is down, how do we know if you're still alive?

Sign In or Register to comment.