// sketch:   PG_TuringTunnel.pde
// purpose:  combination of a simple tunnel sketch with a turing pattern simulation
// version:  v1.0  2012-01-18
//
// for standard tunnel example see...  
//   http://processing.org/learning/topics/tunnel.html
// or for a tiny tunnel sketch
//   http://www.openprocessing.org/visuals/?visualID=3621
/*
keys:
[b] toggle bluring
[f] toggle fps
[h] toggle help
[m] toggle map
[o] toggle rotation
[p] change pattern
[r] reset
[s] save picture
[t] toggle turning simulation
[ ] blur once
[+],[-] change limits
[up],[down] change resolution
[left],[right] change iteration
*/

int ts = 512;  // 800   
PImage wallpaper;
int[] distances;
int[] angles;

void setup() 
{
  size(800, 600);
  frameRate (30);
  float w2 = width / 2.0;
  distances = new int[width*height];
  angles = new int[width*height];
  wallpaper = createImage(ts,ts,RGB);
  
  // prepare distances and angles
  for (int y=0, ni=0; y < height; y++)
  { float y1 = y - height * 0.5;
    float y2 = y1*y1;
    for (int x=0; x < width; x++, ni++) 
    { distances[ni] = (int) ((33 * ts / sqrt(sq(x - w2) + y2)) % ts);
      angles[ni] = (int) (0.5 * ts * atan2(y1, x - w2) / PI);
    }
  }
  resetTuring();
}

int shiftX = ts;
int shiftY = ts;

void tunnel(boolean animateTunnel)
{
  if (animateTunnel) 
  { shiftX = ts + frameCount;  
    shiftY = ts + frameCount*3;
  }
  loadPixels();
  for (int y = 0, ni = 0; y < height; y++) 
  { for (int x = 0; x < width; x++, ni++) 
    { int c_x = (distances[ni] + shiftX) % ts;
      int c_y = (angles[ni] + shiftY) % ts;
      pixels[ni] = wallpaper.pixels[c_x + c_y * ts];
    }
  }
  updatePixels();
}

boolean animateTunnel = true;
boolean simulateTuring = true;
boolean showMap = false;
boolean showHelp = false;
boolean showFPS = false;

void keyPressed()
{ 
  switch(key) 
  { case ' ': if (!simulateTuring) doTuring();  break;
    case 'f': showFPS = !showFPS;    break;
    case 'h': showHelp = !showHelp;  break;
    case 'm': showMap = !showMap;    break;
    case 'o': animateTunnel = !animateTunnel;    break;
    case 't': simulateTuring = !simulateTuring;  break;
    case 's': save("TuringTunnel.png");   break;
    default: keyPressed2();
  }
}

void draw()
{
  // draw img
  if (simulateTuring) doTuring();
  // draw wallpaper
  wallpaper.copy(img, 0,0, img.width, img.height, 0,0, wallpaper.width, wallpaper.height);
  // draw image
  tunnel(animateTunnel);
  // draw turing map   
  if (showMap) image(img, 22, 22);
  // draw help
  if (showHelp) text ("keys: +,-,b,f,h,m,o,p,r,s,t,cursor keys",20,14);
  // draw fps
  if (showFPS) text (round(frameRate) + " fps", width-60,20);
}

/////////////////////////////////////////////////
//                                             //
//    The Secret Life of Turing Patterns       //
//                                             //
/////////////////////////////////////////////////

// Inspired by the work of Jonathan McCabe
// (c) Martin Schneider 2010

// for original code of "Turning Patterns" see...
//   http://openprocessing.org/visuals/?visualID=17043

int scl = 4, dirs = 19, lim = 128;
int res = 5, patternId = 2, bluring = 0;
int dx, dy, w, h, imgSize;
float[] pat;
PImage img;
  
void keyPressed2() 
{ 
  switch(key) {
    case 'b': bluring = (bluring + 1) % 2;  doTuring();  break;
    case 'p': patternId = (patternId + 1) % 3;  break;
    case 'r': res = 3;  resetTuring();  doTuring(); break;
    case '+': lim = min(lim+8, 255);  break;
    case '-': lim = max(lim-8, 0);  break;
    case CODED:
      switch(keyCode) {
        case LEFT:  scl = max(scl-1, 2); break;
        case RIGHT: scl = min(scl+1, 6); break;
        case UP:    res = max(res-1, 2); resetTuring(); doTuring(); break;
        case DOWN:  res = min(res+1, 4); resetTuring(); doTuring(); break;
      }
      break;
  }
}
 
// moving the canvas
void mouseDragged()
{
  if(mousePressed) 
  {
    dy = mod(dx - mouseX + pmouseX, width);
    dx = mod(dy + mouseY - pmouseY, height);
  }
}
 
void resetTuring() 
{
  colorMode(HSB);
  w = 512/res;
  h = 512/res;
  imgSize = w*h;
  img = createImage(w, h, RGB);
  pat = new float[imgSize];
  // random init
  for(int i=0; i<imgSize; i++) 
    pat[i] = floor(random(256));
}
 
void doTuring() 
{
  // calculate a single pattern step
  pattern();
   
  // draw chemicals to the canvas
  img.loadPixels();
  for(int y=0; y<h; y++) 
    for(int x=0; x<w; x++)
    { int c = (x+dx/res)%w + ((y+dy/res)%h)*w;
      float val = pat[x+y*w];
      img.pixels[c] = color(val, 255-val, 100 + val / 2); 
//      img.pixels[c] = color(255-val, 255-val, 100+val / 2); 
    }
  img.updatePixels();
   
  if (bluring == 1) img.filter(BLUR);
}
 
// floor modulo
final int mod(int a, int n) 
{
  return a>=0 ? a%n : (n-1)-(-a-1)%n;
}


//--------------------------------------------------------
// this is where the magic happens ...
 
void pattern() 
{
  // random angular offset
  float R = random(TWO_PI);
 
  // copy chemicals
  float[] pnew = new float[imgSize];
  for(int i=0; i<imgSize; i++) pnew[i] = pat[i];
 
  // create matrices
  float[][] pmedian = new float[imgSize][scl];
  float[][] prange = new float[imgSize][scl];
  float[][] pvar = new float[imgSize][scl];
 
  // iterate over increasing distances
  for(int i=0; i<scl; i++) 
  {
    float d = (2<<i) ;
     
    // update median matrix
    for(int j=0; j<dirs; j++) 
    {
      float dir = j*TWO_PI/dirs + R;
      int dx = int (d * cos(dir));
      int dy = int (d * sin(dir));
      for(int l=0; l<imgSize; l++) 
      { 
        // coordinates of the connected cell
        int x1 = l%w + dx, y1 = l/w + dy;
        // skip if the cell is beyond the border or wrap around
        if(x1<0) x1 = w-1-(-x1-1)% w; else if(x1>=w) x1 = x1%w;
        if(y1<0) y1 = h-1-(-y1-1)% h; else if(y1>=h) y1 = y1%h;
        // update median
        pmedian[l][i] += pat[x1+y1*w] / dirs;
      }
    }
     
    // update range and variance matrix
    for(int j=0; j<dirs; j++) 
    {
      float dir = j*TWO_PI/dirs + R;
      int dx = int (d * cos(dir));
      int dy = int (d * sin(dir));
      for(int l=0; l<imgSize; l++) 
      { 
        // coordinates of the connected cell
        int x1 = l%w + dx, y1 = l/w + dy;
        // skip if the cell is beyond the border or wrap around
        if(x1<0) x1 = w-1-(-x1-1)% w; else if(x1>=w) x1 = x1%w;
        if(y1<0) y1 = h-1-(-y1-1)% h; else if(y1>=h) y1 = y1%h;
        // update variance
        pvar[l][i] += abs( pat[x1+y1*w]  - pmedian[l][i] ) / dirs;
        // update range
        prange[l][i] += pat[x1+y1*w] > (lim + i*10) ? +1 : -1;   
      }
    }    
  }
 
  for(int l=0; l<imgSize; l++) 
  {    
    // find min and max variation
    int imin=0, imax=scl;
    float vmin = MAX_FLOAT;
    float vmax = -MAX_FLOAT;
    for(int i=0; i<scl; i++) 
    {
      if (pvar[l][i] <= vmin) { vmin = pvar[l][i]; imin = i; }
      if (pvar[l][i] >= vmax) { vmax = pvar[l][i]; imax = i; }
    }
     
    // turing pattern variants
    switch(patternId) 
    { case 0: for(int i=0;    i<=imin; i++) pnew[l] += prange[l][i]; break;
      case 1: for(int i=imin; i<=imax; i++) pnew[l] += prange[l][i]; break;
      case 2: for(int i=imin; i<=imax; i++) pnew[l] += prange[l][i] + pvar[l][i]/2; break;
    }
  }
 
  // rescale values
  float vmin = MAX_FLOAT;
  float vmax = -MAX_FLOAT;
  for(int i=0; i<imgSize; i++)  
  {
    vmin = min(vmin, pnew[i]);
    vmax = max(vmax, pnew[i]);
  }      
  float dv = vmax - vmin;
  for(int i=0; i<imgSize; i++)
    pat[i] = (pnew[i] - vmin) * 255 / dv;
}