moving the OpenGL examples to the core

java/examples/OpenGL/Advanced/Particles/Particle.pde

@@ -0,0 +1,58 @@
+class Particle {
+
+  PVector velocity;
+  float lifespan = 255;
+  
+  PShape part;
+  float partSize;
+  
+  PVector gravity = new PVector(0,0.1);
+
+
+  Particle() {
+    partSize = random(10,60);
+    part = createShape(QUAD);
+    part.noStroke();
+    part.texture(sprite);
+    part.normal(0, 0, 1);
+    part.vertex(-partSize/2, -partSize/2, 0, 0);
+    part.vertex(+partSize/2, -partSize/2, sprite.width, 0);
+    part.vertex(+partSize/2, +partSize/2, sprite.width, sprite.height);
+    part.vertex(-partSize/2, +partSize/2, 0, sprite.height);
+    part.end();
+    
+    rebirth(width/2,height/2);
+    lifespan = random(255);
+  }
+
+  PShape getShape() {
+    return part;
+  }
+  
+  void rebirth(float x, float y) {
+    float a = random(TWO_PI);
+    float speed = random(0.5,4);
+    velocity = new PVector(cos(a), sin(a));
+    velocity.mult(speed);
+    lifespan = 255;   
+    part.resetMatrix();
+    part.translate(x, y); 
+  }
+  
+  boolean isDead() {
+    if (lifespan < 0) {
+     return true;
+    } else {
+     return false;
+    } 
+  }
+  
+
+  public void update() {
+    lifespan = lifespan - 1;
+    velocity.add(gravity);
+    
+    part.tint(255,lifespan);
+    part.translate(velocity.x, velocity.y);
+  }
+}

java/examples/OpenGL/Advanced/Particles/ParticleSystem.pde

@@ -0,0 +1,36 @@
+class ParticleSystem {
+  ArrayList<Particle> particles;
+
+  PShape particleShape;
+
+  ParticleSystem(int n) {
+    particles = new ArrayList<Particle>();
+    particleShape = createShape(PShape.GROUP);
+
+    for (int i = 0; i < n; i++) {
+      Particle p = new Particle();
+      particles.add(p);
+      particleShape.addChild(p.getShape());
+    }
+  }
+
+  void update() {
+    for (Particle p : particles) {
+      p.update();
+    }
+  }
+
+  void setEmitter(float x, float y) {
+    for (Particle p : particles) {
+      if (p.isDead()) {
+        p.rebirth(x, y);
+      }
+    }
+  }
+
+  void display() {
+
+    shape(particleShape);
+  }
+}
+

java/examples/OpenGL/Advanced/Particles/Particles.pde

@@ -0,0 +1,31 @@
+// Particles, by Daniel Shiffman
+
+ParticleSystem ps;
+PImage sprite;  
+
+void setup() {
+  size(640, 400, P3D);
+  orientation(LANDSCAPE);
+  sprite = loadImage("sprite.png");
+  ps = new ParticleSystem(10000);
+
+  // Writing to the depth buffer is disabled to avoid rendering
+  // artifacts due to the fact that the particles are semi-transparent
+  // but not z-sorted.
+  hint(DISABLE_DEPTH_MASK);
+} 
+
+void draw () {
+  background(0);
+  ps.update();
+  ps.display();
+  
+  ps.setEmitter(mouseX,mouseY);
+  
+  fill(255);
+  textSize(16);
+  text("Frame rate: " + int(frameRate),10,20);
+  
+}
+
+

java/examples/OpenGL/Advanced/Patch/Patch.pde

@@ -0,0 +1,161 @@
+// Bezier patch By Maritus Watz:
+// http://www.openprocessing.org/sketch/57709
+// Normal calculation added by Andres Colubri
+// Direct port of sample code by Paul Bourke.
+// Original code: http://paulbourke.net/geometry/bezier/
+
+int ni=4, nj=5, RESI=ni*10, RESJ=nj*10;
+PVector outp[][], inp[][];
+PVector normp[][]; 
+boolean autoNormals = false; 
+ 
+void setup() {
+  size(600, 600, P3D);
+  build();
+}
+ 
+void draw() {
+  background(255);
+  translate(width/2,height/2);
+  lights();
+  scale(0.9);
+  rotateY(map(mouseX,0,width,-PI,PI));
+  rotateX(map(mouseY,0,height,-PI,PI));
+   
+  noStroke();
+  fill(255);
+  for(int i=0; i<RESI-1; i++) {
+    beginShape(QUAD_STRIP);
+    for(int j=0; j<RESJ; j++) {
+      if (!autoNormals) {
+        normal(normp[i][j].x, normp[i][j].y, normp[i][j].z);
+      }
+      vertex(outp[i][j].x,outp[i][j].y,outp[i][j].z);
+      vertex(outp[i+1][j].x,outp[i+1][j].y,outp[i+1][j].z);     
+    }
+    endShape();
+  }
+}
+ 
+void keyPressed() {
+  if(key==' ') build();
+  saveFrame("bezPatch.png");
+}
+ 
+void build() {
+  int i, j, ki, kj;
+  double mui, muj, bi, bj, dbi, dbj;
+  
+  outp=new PVector[RESI][RESJ];
+  normp=new PVector[RESI][RESJ];
+  inp=new PVector[ni+1][nj+1];
+  PVector uitang = new PVector();
+  PVector ujtang = new PVector();
+  
+  for (i=0;i<=ni;i++) {
+    for (j=0;j<=nj;j++) {
+      inp[i][j]=new PVector(i,j,random(-3,3));
+    }
+  }
+  
+  for (i=0;i<RESI;i++) {
+    mui = i / (double)(RESI-1);
+    for (j=0;j<RESJ;j++) {
+      muj = j / (double)(RESJ-1);
+      outp[i][j]=new PVector();
+      uitang.set(0, 0, 0);
+      ujtang.set(0, 0, 0);
+      for (ki=0;ki<=ni;ki++) {
+        bi = BezierBlend(ki, mui, ni);
+        dbi = DBezierBlend(ki, mui, ni);
+        for (kj=0;kj<=nj;kj++) {
+          bj = BezierBlend(kj, muj, nj);
+          dbj = DBezierBlend(kj, muj, nj);
+          outp[i][j].x += (inp[ki][kj].x * bi * bj);
+          outp[i][j].y += (inp[ki][kj].y * bi * bj);
+          outp[i][j].z += (inp[ki][kj].z * bi * bj);
+           
+          uitang.x += (inp[ki][kj].x * dbi * bj);
+          uitang.y += (inp[ki][kj].y * dbi * bj);
+          uitang.z += (inp[ki][kj].z * dbi * bj);
+           
+          ujtang.x += (inp[ki][kj].x * bi * dbj);
+          ujtang.y += (inp[ki][kj].y * bi * dbj);
+          ujtang.z += (inp[ki][kj].z * bi * dbj);
+        }
+      }
+      outp[i][j].add(new PVector(-ni/2,-nj/2,0));
+      outp[i][j].mult(100);
+       
+      uitang.normalize();
+      ujtang.normalize();
+      normp[i][j] = uitang.cross(ujtang);     
+    }
+  }
+}
+ 
+double BezierBlend(int k, double mu, int n) {
+  int nn, kn, nkn;
+  double blend=1;
+ 
+  nn = n;
+  kn = k;
+  nkn = n - k;
+ 
+  while (nn >= 1) {
+    blend *= nn;
+    nn--;
+    if (kn > 1) {
+      blend /= (double)kn;
+      kn--;
+    }
+    if (nkn > 1) {
+      blend /= (double)nkn;
+      nkn--;
+    }
+  }
+  if (k > 0)
+    blend *= Math.pow(mu, (double)k);
+  if (n-k > 0)
+    blend *= Math.pow(1-mu, (double)(n-k));
+ 
+  return(blend);
+}
+
+double DBezierBlend(int k, double mu, int n) {
+  int nn, kn, nkn;
+  double dblendf = 1;
+  
+  nn = n;
+  kn = k;
+  nkn = n - k;
+  
+  while (nn >= 1) {
+    dblendf *= nn;
+    nn--;
+    if (kn > 1) {
+      dblendf /= (double)kn;
+      kn--;
+    }
+    if (nkn > 1) {
+      dblendf /= (double)nkn;
+      nkn--;
+    }
+  }
+   
+  double fk = 1;
+  double dk = 0;
+  double fnk = 1;
+  double dnk = 0;   
+  if (k > 0) {
+    fk = Math.pow(mu, (double)k);
+    dk = k*Math.pow(mu, (double)k-1);
+  } 
+  if (n-k > 0) {
+    fnk = Math.pow(1-mu, (double)(n-k));
+    dnk = (k-n)*Math.pow(1-mu, (double)(n-k-1));
+  }
+  dblendf *= (dk * fnk + fk * dnk);
+  
+  return(dblendf); 
+}

java/examples/OpenGL/Advanced/Planets/Perlin.pde

@@ -0,0 +1,261 @@
+// Implementation of 1D, 2D, and 3D Perlin noise. Based on the 
+// C code by Paul Bourke:
+// http://local.wasp.uwa.edu.au/~pbourke/texture_colour/perlin/
+class Perlin {
+  int B = 0x100;
+  int BM = 0xff;
+  int N = 0x1000;
+  int NP = 12; 
+  int NM = 0xfff;
+
+  int p[];
+  float g3[][];
+  float g2[][];
+  float g1[];
+
+  void normalize2(float v[]) {
+    float s = sqrt(v[0] * v[0] + v[1] * v[1]);
+    v[0] = v[0] / s;
+    v[1] = v[1] / s;
+  }
+
+  void normalize3(float v[]) {
+    float s = sqrt(v[0] * v[0] + v[1] * v[1] + v[2] * v[2]);
+    v[0] = v[0] / s;
+    v[1] = v[1] / s;
+    v[2] = v[2] / s;
+  }
+
+  float sCurve(float t) {
+    return t * t * (3.0 - 2.0 * t);
+  }
+
+  float at2(float q[], float rx, float ry) { 
+    return rx * q[0] + ry * q[1];
+  }
+
+  float at3(float q[], float rx, float ry, float rz) { 
+    return rx * q[0] + ry * q[1] + rz * q[2];
+  }
+
+  Perlin() {
+    p = new int[B + B + 2];
+    g3 = new float[B + B + 2][3];
+    g2 = new float[B + B + 2][2];
+    g1 = new float[B + B + 2];  
+
+    init();
+  }  
+
+  void init() {
+    int i, j, k;
+
+    for (i = 0 ; i < B ; i++) {
+      p[i] = i;
+      g1[i] = (random(B + B) - B) / B;
+
+      for (j = 0 ; j < 2 ; j++)
+        g2[i][j] = (random(B + B) - B) / B;
+      normalize2(g2[i]);
+
+      for (j = 0 ; j < 3 ; j++)
+        g3[i][j] = (random(B + B) - B) / B;
+      normalize3(g3[i]);
+    }
+
+    while (0 < --i) {
+      k = p[i];
+      p[i] = p[j = int(random(B))];
+      p[j] = k;
+    }
+
+    for (i = 0 ; i < B + 2 ; i++) {
+      p[B + i] = p[i];
+      g1[B + i] = g1[i];
+      for (j = 0 ; j < 2 ; j++)
+        g2[B + i][j] = g2[i][j];
+      for (j = 0 ; j < 3 ; j++)
+        g3[B + i][j] = g3[i][j];
+    }
+  }
+
+  float noise1(float[] vec) {
+    int bx0, bx1;
+    float rx0, rx1, sx, t, u, v;
+
+    t = vec[0] + N;
+    bx0 = int(t) & BM;
+    bx1 = (bx0 + 1) & BM;
+    rx0 = t - int(t);
+    rx1 = rx0 - 1.0; 
+
+    sx = sCurve(rx0);
+    u = rx0 * g1[p[bx0]];
+    v = rx1 * g1[p[bx1]];
+
+    return lerp(u, v, sx);
+  }
+
+  float noise2(float[] vec) {
+    int bx0, bx1, by0, by1, b00, b10, b01, b11;
+    float rx0, rx1, ry0, ry1, sx, sy, a, b, t, u, v;
+    float[] q;    
+    int i, j;
+
+    t = vec[0] + N;
+    bx0 = int(t) & BM;
+    bx1 = (bx0 + 1) & BM;
+    rx0 = t - int(t);
+    rx1 = rx0 - 1.0; 
+
+    t = vec[1] + N;
+    by0 = int(t) & BM;
+    by1 = (by0 + 1) & BM;
+    ry0 = t - int(t);
+    ry1 = ry0 - 1.0;
+
+    i = p[bx0];
+    j = p[bx1];
+
+    b00 = p[i + by0];
+    b10 = p[j + by0];
+    b01 = p[i + by1];
+    b11 = p[j + by1];
+
+    sx = sCurve(rx0);
+    sy = sCurve(ry0);
+
+    q = g2[b00]; 
+    u = at2(q, rx0, ry0);
+    q = g2[b10]; 
+    v = at2(q, rx1, ry0);
+    a = lerp(u, v, sx);
+
+    q = g2[b01] ; 
+    u = at2(q, rx0, ry1);
+    q = g2[b11] ; 
+    v = at2(q, rx1, ry1);
+    b = lerp(u, v, sx);
+
+    return lerp(a, b, sy);
+  }
+
+  float noise3(float[] vec) {
+    int bx0, bx1, by0, by1, bz0, bz1, b00, b10, b01, b11;
+    float rx0, rx1, ry0, ry1, rz0, rz1, sy, sz, a, b, c, d, t, u, v;
+    float[] q;
+    int i, j;
+
+    t = vec[0] + N;
+    bx0 = int(t) & BM;
+    bx1 = (bx0 + 1) & BM;
+    rx0 = t - int(t);
+    rx1 = rx0 - 1.0;
+
+    t = vec[1] + N;
+    by0 = int(t) & BM;
+    by1 = (by0 + 1) & BM;
+    ry0 = t - int(t);
+    ry1 = ry0 - 1.0;
+
+    t = vec[2] + N;
+    bz0 = int(t) & BM;
+    bz1 = (bz0 + 1) & BM;
+    rz0 = t - int(t);
+    rz1 = rz0 - 1.0;
+
+    i = p[bx0];
+    j = p[bx1];
+
+    b00 = p[i + by0];
+    b10 = p[j + by0];
+    b01 = p[i + by1];
+    b11 = p[j + by1];
+
+    t  = sCurve(rx0);
+    sy = sCurve(ry0);
+    sz = sCurve(rz0);
+
+    q = g3[b00 + bz0]; 
+    u = at3(q, rx0, ry0, rz0);
+    q = g3[b10 + bz0]; 
+    v = at3(q, rx1, ry0, rz0);
+    a = lerp(u, v, t);
+
+    q = g3[b01 + bz0]; 
+    u = at3(q, rx0, ry1, rz0);
+    q = g3[b11 + bz0]; 
+    v = at3(q, rx1, ry1, rz0);
+    b = lerp(u, v, t);
+
+    c = lerp(a, b, sy);
+
+    q = g3[b00 + bz1]; 
+    u = at3(q, rx0, ry0, rz1);
+    q = g3[b10 + bz1]; 
+    v = at3(q, rx1, ry0, rz1);
+    a = lerp(u, v, t);
+
+    q = g3[b01 + bz1]; 
+    u = at3(q, rx0, ry1, rz1);
+    q = g3[b11 + bz1]; 
+    v = at3(q, rx1, ry1, rz1);
+    b = lerp(u, v, t);
+
+    d = lerp(a, b, sy);
+
+    return lerp(c, d, sz);
+  }
+
+  // In what follows "nalpha" is the weight when the sum is formed.
+  // Typically it is 2, as this approaches 1 the function is noisier.
+  // "nbeta" is the harmonic scaling/spacing, typically 2. n is the
+  // number of harmonics added up in the final result. Higher number 
+  // results in more detailed noise.
+
+  float noise1D(float x, float nalpha, float nbeta, int n) {
+    float val, sum = 0;
+    float v[] = {x};
+    float nscale = 1;
+
+    for (int i = 0; i < n; i++) {
+      val = noise1(v);
+      sum += val / nscale;
+      nscale *= nalpha;
+      v[0] *= nbeta;
+    }
+    return sum;
+  }
+
+  float noise2D(float x, float y, float nalpha, float nbeta, int n) {
+   float val,sum = 0;
+   float v[] = {x, y};
+   float nscale = 1;
+
+   for (int i = 0; i < n; i++) {
+      val = noise2(v);
+      sum += val / nscale;
+      nscale *= nalpha;
+      v[0] *= nbeta;
+      v[1] *= nbeta;
+   }
+   return sum;
+  }
+
+  float noise3D(float x, float y, float z, float nalpha, float nbeta, int n) {
+    float val, sum = 0;
+    float v[] = {x, y, z};
+    float nscale = 1;
+
+    for (int i = 0 ; i < n; i++) {
+      val = noise3(v);
+      sum += val / nscale;
+      nscale *= nalpha;
+      v[0] *= nbeta;
+      v[1] *= nbeta;
+      v[2] *= nbeta;
+    }
+    return sum;
+  }
+}
+

java/examples/OpenGL/Advanced/Planets/Planets.pde

@@ -0,0 +1,122 @@
+// Planets, by Andres Colubri
+//
+// Sun and mercury textures from http://planetpixelemporium.com
+// Star field picture from http://www.galacticimages.com/
+
+PImage starfield;
+
+PShape sun;
+PImage suntex;
+
+PShape planet1;
+PImage surftex1;
+PImage cloudtex;
+
+PShape planet2;
+PImage surftex2;
+
+void setup() {
+  size(800,  480,  P3D);
+  
+  starfield = loadImage("starfield.jpg");
+  suntex = loadImage("sun.jpg");  
+  surftex1 = loadImage("planet.jpg");  
+   
+  // We need trilinear sampling for this texture so it looks good
+  // even when rendered very small.
+  //PTexture.Parameters params1 = PTexture.newParameters(ARGB, TRILINEAR);  
+  surftex2 = loadImage("mercury.jpg");  
+  
+  /*
+  // The clouds texture will "move" having the values of its u
+  // texture coordinates displaced by adding a constant increment
+  // in each frame. This requires REPEAT wrapping mode so texture 
+  // coordinates can be larger than 1.
+  //PTexture.Parameters params2 = PTexture.newParameters();
+  //params2.wrapU = REPEAT;
+  cloudtex = createImage(512, 256);
+
+  // Using 3D Perlin noise to generate a clouds texture that is seamless on
+  // its edges so it can be applied on a sphere.
+  cloudtex.loadPixels();
+  Perlin perlin = new Perlin();
+  for (int j = 0; j < cloudtex.height; j++) {
+    for (int i = 0; i < cloudtex.width; i++) {
+      // The angle values corresponding to each u,v pair:
+      float u = float(i) / cloudtex.width;
+      float v = float(j) / cloudtex.height;
+      float phi = map(u, 0, 1, TWO_PI, 0); 
+      float theta = map(v, 0, 1, -HALF_PI, HALF_PI);
+      // The x, y, z point corresponding to these angles:
+      float x = cos(phi) * cos(theta);
+      float y = sin(theta);            
+      float z = sin(phi) * cos(theta);      
+      float n = perlin.noise3D(x, y, z, 1.2, 2, 8);
+      cloudtex.pixels[j * cloudtex.width + i] = color(255, 255,  255, 255 * n * n);
+    }
+  }  
+  cloudtex.updatePixels();
+  */
+
+  noStroke();
+  fill(255);
+  sphereDetail(40);
+
+  sun = createShape(SPHERE, 150);
+  sun.texture(suntex);  
+
+  planet1 = createShape(SPHERE, 150);
+  planet1.texture(surftex1);
+  
+  planet2 = createShape(SPHERE, 50);
+  planet2.texture(surftex2);
+}
+
+void draw() {
+  // Even we draw a full screen image after this, it is recommended to use
+  // background to clear the screen anyways, otherwise A3D will think
+  // you want to keep each drawn frame in the framebuffer, which results in 
+  // slower rendering.
+  background(0);
+  
+  // Disabling writing to the depth mask so the 
+  // background image doesn't occludes any 3D object.
+  hint(DISABLE_DEPTH_MASK);
+  image(starfield, 0, 0, width, height);
+  hint(ENABLE_DEPTH_MASK);
+  
+  /*
+  // Displacing the u texture coordinate of layer 1 in planet
+  // so it creates the effect of moving clouds.
+  PShape3D p = (PShape3D)planet1;
+  p.loadTexcoords(1);
+  for (int i = 0; i < p.getVertexCount(); i++) {
+      float u = p.texcoords[2 * i + 0];
+      u += 0.002;
+      p.texcoords[2 * i + 0] = u;
+  }
+  p.updateTexcoords();
+  */
+  
+  pushMatrix();
+  translate(width/2, height/2, -300);  
+  
+  pushMatrix();
+  rotateY(PI * frameCount / 500);
+  shape(sun);
+  popMatrix();
+
+  pointLight(255,  255,  255,  0,  0,  0);  
+  rotateY(PI * frameCount / 300);
+  translate(0, 0, 300);
+
+  shape(planet2);  
+  
+  popMatrix();
+  
+  noLights();
+  pointLight(255,  255,  255,  0,  0,  -150); 
+  
+  translate(0.75 * width,  0.6 * height,  50);
+  shape(planet1);
+}

java/examples/OpenGL/Advanced/Ribbons/ArcBall.pde

@@ -0,0 +1,190 @@
+// Ariel and V3ga's arcball class with a couple tiny mods by Robert Hodgin
+
+class Arcball {
+  float center_x, center_y, radius;
+  Vec3 v_down, v_drag;
+  Quat q_now, q_down, q_drag;
+  Vec3[] axisSet;
+  int axis;
+  float mxv, myv;
+  float x, y;
+  
+  Arcball(float center_x, float center_y, float radius){
+    this.center_x = center_x;
+    this.center_y = center_y;
+    this.radius = radius;
+
+    v_down = new Vec3();
+    v_drag = new Vec3();
+
+    q_now = new Quat();
+    q_down = new Quat();
+    q_drag = new Quat();
+
+    axisSet = new Vec3[] {new Vec3(1.0f, 0.0f, 0.0f), new Vec3(0.0f, 1.0f, 0.0f), new Vec3(0.0f, 0.0f, 1.0f)};
+    axis = -1;  // no constraints...    
+  }
+
+  void mousePressed(){
+    v_down = mouse_to_sphere(mouseX, mouseY);
+    q_down.set(q_now);
+    q_drag.reset();
+  }
+
+  void mouseDragged(){
+    v_drag = mouse_to_sphere(mouseX, mouseY);
+    q_drag.set(Vec3.dot(v_down, v_drag), Vec3.cross(v_down, v_drag));
+  }
+
+  void run(){
+    q_now = Quat.mul(q_drag, q_down);
+    applyQuat2Matrix(q_now);
+    
+    x += mxv;
+    y += myv;
+    mxv -= mxv * .01;
+    myv -= myv * .01;
+  }
+  
+  Vec3 mouse_to_sphere(float x, float y){
+    Vec3 v = new Vec3();
+    v.x = (x - center_x) / radius;
+    v.y = (y - center_y) / radius;
+
+    float mag = v.x * v.x + v.y * v.y;
+    if (mag > 1.0f){
+      v.normalize();
+    } else {
+      v.z = sqrt(1.0f - mag);
+    }
+
+    return (axis == -1) ? v : constrain_vector(v, axisSet[axis]);
+  }
+
+  Vec3 constrain_vector(Vec3 vector, Vec3 axis){
+    Vec3 res = new Vec3();
+    res.sub(vector, Vec3.mul(axis, Vec3.dot(axis, vector)));
+    res.normalize();
+    return res;
+  }
+
+  void applyQuat2Matrix(Quat q){
+    // instead of transforming q into a matrix and applying it...
+
+    float[] aa = q.getValue();
+    rotate(aa[0], aa[1], aa[2], aa[3]);
+  }
+}
+
+static class Vec3{
+  float x, y, z;
+
+  Vec3(){
+  }
+
+  Vec3(float x, float y, float z){
+    this.x = x;
+    this.y = y;
+    this.z = z;
+  }
+
+  void normalize(){
+    float length = length();
+    x /= length;
+    y /= length;
+    z /= length;
+  }
+
+  float length(){
+    return (float) Math.sqrt(x * x + y * y + z * z);
+  }
+
+  static Vec3 cross(Vec3 v1, Vec3 v2){
+    Vec3 res = new Vec3();
+    res.x = v1.y * v2.z - v1.z * v2.y;
+    res.y = v1.z * v2.x - v1.x * v2.z;
+    res.z = v1.x * v2.y - v1.y * v2.x;
+    return res;
+  }
+
+  static float dot(Vec3 v1, Vec3 v2){
+    return v1.x * v2.x + v1.y * v2.y + v1.z * v2.z;
+  }
+  
+  static Vec3 mul(Vec3 v, float d){
+    Vec3 res = new Vec3();
+    res.x = v.x * d;
+    res.y = v.y * d;
+    res.z = v.z * d;
+    return res;
+  }
+
+  void sub(Vec3 v1, Vec3 v2){
+    x = v1.x - v2.x;
+    y = v1.y - v2.y;
+    z = v1.z - v2.z;
+  }
+}
+
+static class Quat{
+  float w, x, y, z;
+
+  Quat(){
+    reset();
+  }
+
+  Quat(float w, float x, float y, float z){
+    this.w = w;
+    this.x = x;
+    this.y = y;
+    this.z = z;
+  }
+
+  void reset(){
+    w = 1.0f;
+    x = 0.0f;
+    y = 0.0f;
+    z = 0.0f;
+  }
+
+  void set(float w, Vec3 v){
+    this.w = w;
+    x = v.x;
+    y = v.y;
+    z = v.z;
+  }
+
+  void set(Quat q){
+    w = q.w;
+    x = q.x;
+    y = q.y;
+    z = q.z;
+  }
+
+  static Quat mul(Quat q1, Quat q2){
+    Quat res = new Quat();
+    res.w = q1.w * q2.w - q1.x * q2.x - q1.y * q2.y - q1.z * q2.z;
+    res.x = q1.w * q2.x + q1.x * q2.w + q1.y * q2.z - q1.z * q2.y;
+    res.y = q1.w * q2.y + q1.y * q2.w + q1.z * q2.x - q1.x * q2.z;
+    res.z = q1.w * q2.z + q1.z * q2.w + q1.x * q2.y - q1.y * q2.x;
+    return res;
+  }
+  
+  float[] getValue(){
+    // transforming this quat into an angle and an axis vector...
+
+    float[] res = new float[4];
+
+    float sa = (float) Math.sqrt(1.0f - w * w);
+    if (sa < EPSILON){
+      sa = 1.0f;
+    }
+
+    res[0] = (float) Math.acos(w) * 2.0f;
+    res[1] = x / sa;
+    res[2] = y / sa;
+    res[3] = z / sa;
+
+    return res;
+  }
+}

java/examples/OpenGL/Advanced/Ribbons/BSpline.pde

@@ -0,0 +1,307 @@
+final int MAX_BEZIER_ORDER = 10; // Maximum curve order.
+
+final float[][] BSplineMatrix = {
+  {-1.0/6.0,  1.0/2.0, -1.0/2.0, 1.0/6.0},
+  { 1.0/2.0,     -1.0,  1.0/2.0,     0.0},
+  {-1.0/2.0,      0.0,  1.0/2.0,     0.0},
+  { 1.0/6.0,  2.0/3.0,  1.0/6.0,     0.0}
+};
+
+// The element(i, n) of this array contains the binomial coefficient
+// C(i, n) = n!/(i!(n-i)!)
+final int[][] BinomialCoefTable = {
+  {1, 1, 1, 1,  1,  1,  1,  1,   1,   1},
+  {1, 2, 3, 4,  5,  6,  7,  8,   9,  10},
+  {0, 1, 3, 6, 10, 15, 21, 28,  36,  45},
+  {0, 0, 1, 4, 10, 20, 35, 56,  84, 120},
+  {0, 0, 0, 1,  5, 15, 35, 70, 126, 210},
+  {0, 0, 0, 0,  1,  6, 21, 56, 126, 252},
+  {0, 0, 0, 0,  0,  1,  7, 28,  84, 210},
+  {0, 0, 0, 0,  0,  0,  1,  8,  36, 120},
+  {0, 0, 0, 0,  0,  0,  0,  1,   9,  45},
+  {0, 0, 0, 0,  0,  0,  0,  0,   1,  10},
+  {0, 0, 0, 0,  0,  0,  0,  0,   0,   1}
+};
+
+// The element of this(i, j) of this table contains(i/10)^(3-j).
+final float[][] TVectorTable = {  
+//   t^3,  t^2, t^1, t^0
+  {    0,    0,   0,   1}, // t = 0.0
+  {0.001, 0.01, 0.1,   1}, // t = 0.1
+  {0.008, 0.04, 0.2,   1}, // t = 0.2
+  {0.027, 0.09, 0.3,   1}, // t = 0.3
+  {0.064, 0.16, 0.4,   1}, // t = 0.4
+  {0.125, 0.25, 0.5,   1}, // t = 0.5
+  {0.216, 0.36, 0.6,   1}, // t = 0.6
+  {0.343, 0.49, 0.7,   1}, // t = 0.7
+  {0.512, 0.64, 0.8,   1}, // u = 0.8
+  {0.729, 0.81, 0.9,   1}, // t = 0.9
+  {    1,    1,   1,   1}  // t = 1.0
+};
+
+// The element of this(i, j) of this table contains(3-j)*(i/10)^(2-j) if
+// j < 3, 0 otherwise.
+final float[][] DTVectorTable = { 
+// 3t^2,  2t^1, t^0
+  {   0,     0,   1, 0}, // t = 0.0
+  {0.03,   0.2,   1, 0}, // t = 0.1
+  {0.12,   0.4,   1, 0}, // t = 0.2
+  {0.27,   0.6,   1, 0}, // t = 0.3
+  {0.48,   0.8,   1, 0}, // t = 0.4
+  {0.75,   1.0,   1, 0}, // t = 0.5
+  {1.08,   1.2,   1, 0}, // t = 0.6
+  {1.47,   1.4,   1, 0}, // t = 0.7
+  {1.92,   1.6,   1, 0}, // t = 0.8
+  {2.43,   1.8,   1, 0}, // t = 0.9
+  {   3,     2,   1, 0}  // t = 1.0
+};
+
+abstract class Curve3D {
+  abstract void feval(float t, PVector p);
+  abstract void deval(float t, PVector d);
+  abstract float fevalX(float t);
+  abstract float fevalY(float t);
+  abstract float fevalZ(float t);
+  abstract float devalX(float t);
+  abstract float devalY(float t);
+  abstract float devalZ(float t);
+}
+
+abstract class Spline extends Curve3D {
+  // The factorial of n.
+  int factorial(int n) { 
+    return n <= 0 ? 1 : n * factorial(n - 1); 
+  }
+  
+  // Gives n!/(i!(n-i)!).
+  int binomialCoef(int i, int n) {
+    if ((i <= MAX_BEZIER_ORDER) && (n <= MAX_BEZIER_ORDER)) return BinomialCoefTable[i][n - 1];
+    else return int(factorial(n) / (factorial(i) * factorial(n - i)));
+  }
+  
+  // Evaluates the Berstein polinomial(i, n) at u.
+  float bersteinPol(int i, int n, float u) {
+    return binomialCoef(i, n) * pow(u, i) * pow(1 - u, n - i);
+  }
+  
+  // The derivative of the Berstein polinomial.
+  float dbersteinPol(int i, int n, float u) {
+    float s1, s2; 
+    if (i == 0) s1 = 0; 
+    else s1 = i * pow(u, i-1) * pow(1 - u, n - i);
+    if (n == i) s2 = 0; 
+    else s2 = -(n - i) * pow(u, i) * pow(1 - u, n - i - 1);
+    return binomialCoef(i, n) *(s1 + s2);
+  }
+}
+
+class BSpline extends Spline {
+  // Control points.
+  float[][] bsplineCPoints;
+
+  // Parameters.
+  boolean lookup;
+
+  // Auxiliary arrays used in the calculations.
+  float[][] M3;
+  float[] TVector, DTVector;
+
+  // Point and tangent vectors.
+  float[] pt, tg;
+  
+  BSpline() { 
+    initParameters(true); 
+  }
+  
+  BSpline(boolean t) { 
+    initParameters(t); 
+  }
+
+  // Sets lookup table use.
+  void initParameters(boolean t) { 
+    bsplineCPoints = new float[4][3];
+    TVector = new float[4];
+    DTVector = new float[4]; 
+    M3 = new float[4][3];
+    pt = new float[3];
+    tg = new float[3];        
+    lookup = t;
+  }
+  
+  // Sets n-th control point.
+  void setCPoint(int n, PVector P) {
+    bsplineCPoints[n][0] = P.x;
+    bsplineCPoints[n][1] = P.y;
+    bsplineCPoints[n][2] = P.z;        
+    updateMatrix3();
+  }
+
+  // Gets n-th control point.
+  void getCPoint(int n, PVector P) {
+    P.set(bsplineCPoints[n]);
+  }
+
+  // Replaces the current B-spline control points(0, 1, 2) with(1, 2, 3). This
+  // is used when a new spline is to be joined to the recently drawn.
+  void shiftBSplineCPoints() {
+    for (int i = 0; i < 3; i++) {
+      bsplineCPoints[0][i] = bsplineCPoints[1][i];
+      bsplineCPoints[1][i] = bsplineCPoints[2][i];
+      bsplineCPoints[2][i] = bsplineCPoints[3][i];
+    }
+    updateMatrix3();
+  }
+
+  void copyCPoints(int n_source, int n_dest) {
+    for (int i = 0; i < 3; i++) {
+      bsplineCPoints[n_dest][i] = bsplineCPoints[n_source][i];
+    }
+  }
+
+  // Updates the temporal matrix used in order 3 calculations.
+  void updateMatrix3() {
+    float s; 
+    int i, j, k;
+    for(i = 0; i < 4; i++) {
+      for(j = 0; j < 3; j++) {
+        s = 0;
+        for(k = 0; k < 4; k++) s += BSplineMatrix[i][k] * bsplineCPoints[k][j];
+        M3[i][j] = s;
+      }
+    }
+  }    
+
+  void feval(float t, PVector p) { 
+    evalPoint(t); 
+    p.set(pt); 
+  }
+  
+  void deval(float t, PVector d) { 
+    evalTangent(t); 
+    d.set(tg); 
+  }
+
+  float fevalX(float t) { 
+    evalPoint(t); 
+    return pt[0]; 
+  }
+  
+  float fevalY(float t) { 
+    evalPoint(t); 
+    return pt[1]; 
+  }
+  
+  float fevalZ(float t) { 
+    evalPoint(t); 
+    return pt[2]; 
+  }
+
+  float devalX(float t) { 
+    evalTangent(t); 
+    return tg[0]; 
+  }
+  
+  float devalY(float t) { 
+    evalTangent(t); 
+    return tg[1]; 
+  }
+  
+  float devalZ(float t) { 
+    evalTangent(t); 
+    return tg[2]; 
+  }
+
+  // Point evaluation.
+  void evalPoint(float t) {
+    if (lookup) {
+      bsplinePointI(int(10 * t));
+    } else {
+      bsplinePoint(t);
+    }
+  }    
+
+  // Tangent evaluation.
+  void evalTangent(float t) {
+    if (lookup) {
+      bsplineTangentI(int(10 * t));
+    } else {
+      bsplineTangent(t);
+    }
+  }    
+
+  // Calculates the point on the cubic spline corresponding to the parameter value t in [0, 1].
+  void bsplinePoint(float t) {
+    // Q(u) = UVector * BSplineMatrix * BSplineCPoints
+
+    float s;
+    int i, j, k;
+
+    for(i = 0; i < 4; i++) {
+      TVector[i] = pow(t, 3 - i);
+    }
+
+    for(j = 0; j < 3; j++) {
+      s = 0;
+      for(k = 0; k < 4; k++) {
+        s += TVector[k] * M3[k][j];
+      }
+      pt[j] = s;
+    }
+  }
+
+  // Calculates the tangent vector of the spline at t.
+  void bsplineTangent(float t) {
+    // Q(u) = DTVector * BSplineMatrix * BSplineCPoints
+
+    float s;
+    int i, j, k;
+
+    for(i = 0; i < 4; i++) {
+      if (i < 3) {
+        DTVector[i] = (3 - i) * pow(t, 2 - i);
+      } else {
+        DTVector[i] = 0;
+      }
+    }
+
+    for(j = 0; j < 3; j++) {
+      s = 0;
+      for(k = 0; k < 4; k++) {
+        s += DTVector[k] * M3[k][j];
+      }
+      tg[j] = s;
+    }
+  }
+
+  // Gives the point on the cubic spline corresponding to t/10(using the lookup table).
+  void bsplinePointI(int t) {
+    // Q(u) = TVectorTable[u] * BSplineMatrix * BSplineCPoints
+
+    float s;
+    int j, k;
+
+    for(j = 0; j < 3; j++) {
+      s = 0;
+      for(k = 0; k < 4; k++) {
+        s += TVectorTable[t][k] * M3[k][j];
+      }
+      pt[j] = s;
+    }
+  }
+
+  // Calulates the tangent vector of the spline at t/10.
+  void bsplineTangentI(int t) {
+    // Q(u) = DTVectorTable[u] * BSplineMatrix * BSplineCPoints
+
+    float s;
+    int j, k;
+
+    for(j = 0; j < 3; j++) {
+      s = 0;
+      for(k = 0; k < 4; k++) {
+        s += DTVectorTable[t][k] * M3[k][j];
+      }
+      tg[j] = s;
+    }
+  }    
+}

java/examples/OpenGL/Advanced/Ribbons/Geometry.pde

@@ -0,0 +1,495 @@
+BSpline splineSide1;
+BSpline splineCenter;
+BSpline splineSide2;
+PVector flipTestV;
+int uspacing;
+
+int HELIX = 0;
+int STRAND = 1;
+int COIL = 2;
+int LHANDED = -1;
+int RHANDED = 1;
+
+void createRibbonModel(ArrayList residues, PShape model, ArrayList trj) {
+  // For line ribbons
+  ArrayList vertices0 = new ArrayList();
+  ArrayList vertices1 = new ArrayList();
+  ArrayList vertices2 = new ArrayList();
+  
+  // For flat ribbons
+  ArrayList vertices = new ArrayList();
+  ArrayList normals = new ArrayList();
+    
+  if (ribbonDetail == 1) uspacing = 10;
+  else if (ribbonDetail == 2) uspacing = 5;
+  else if (ribbonDetail == 3) uspacing = 2;
+  else uspacing = 1;
+
+  flipTestV = new PVector();
+  splineSide1 = new BSpline(false);
+  splineCenter = new BSpline(false);
+  splineSide2 = new BSpline(false);
+
+  int[] ss = new int[residues.size()];
+  int[] handness = new int[residues.size()];
+
+  calculateSecStr(residues, ss, handness);   
+
+  for (int i = 0; i < residues.size(); i++) {
+    constructControlPoints(residues, i, ss[i], handness[i]);
+
+    if (renderMode == 0) {
+      generateSpline(0, vertices0);
+      generateSpline(1, vertices1);        
+      generateSpline(2, vertices2);        
+    } 
+    else {
+      generateFlatRibbon(vertices, normals);
+    }
+  }  
+
+  if (renderMode == 0) {
+    model = createShape();
+    model.stroke(ribbonColor);
+    model.noFill();    
+    model.beginContour();
+    for (int i = 0; i < vertices0.size(); i++) {
+      PVector posVec = (PVector)vertices0.get(i);
+      model.vertex(posVec.x, posVec.y, posVec.z);    
+    }
+    model.endContour();
+    model.beginContour();
+    for (int i = 0; i < vertices1.size(); i++) {
+      PVector posVec = (PVector)vertices1.get(i);
+      model.vertex(posVec.x, posVec.y, posVec.z);    
+    }
+    model.endContour();    
+    model.beginContour();
+    for (int i = 0; i < vertices2.size(); i++) {
+      PVector posVec = (PVector)vertices2.get(i);
+      model.vertex(posVec.x, posVec.y, posVec.z);    
+    }
+    model.endContour();    
+    model.end(OPEN);    
+  } else {
+    // The ribbon construction is fairly inneficient here, since
+    // it could use triangle strips instead to avoid duplicating
+    // shared vertices...
+    model = createShape(TRIANGLES);
+    model.noStroke();
+    model.fill(ribbonColor);
+    for (int i = 0; i < vertices.size(); i++) {
+      PVector posVec = (PVector)vertices.get(i);
+      PVector normVec = (PVector)normals.get(i);
+      model.normal(-normVec.x, -normVec.y, -normVec.z);
+      model.vertex(posVec.x, posVec.y, posVec.z);
+    }
+    model.end();
+  }
+  
+  trj.add(model);
+
+  if (renderMode == 0) {
+    int totCount = vertices0.size() + vertices1.size() + vertices2.size();
+    println("Adding new model with " + totCount + " vertices.");
+  } else {
+    println("Adding new model with " + vertices.size() + " vertices.");
+  }
+}
+
+float calculateGyrRadius(ArrayList atoms) {
+  PVector ati, atj;
+  float dx, dy, dz;    
+  float r = 0;
+  for (int i = 0; i < atoms.size(); i++) {
+    ati = (PVector)atoms.get(i);
+    for (int j = i + 1; j < atoms.size(); j++) {  
+      atj = (PVector)atoms.get(j);
+
+      dx = ati.x - atj.x;
+      dy = ati.y - atj.y;
+      dz = ati.z - atj.z;
+      r +=  dx * dx + dy * dy + dz * dz;
+    }
+  }
+  return sqrt(r) / (atoms.size() + 1);
+}
+
+// Does a cheap and dirty secondary structure assignment to the protein
+// residues given in the array.
+void calculateSecStr(ArrayList residues, int[] ss, int[] handness) {
+  PVector c0, n1, ca1, c1, n2;
+  HashMap res0, res1, res2;
+  int n = residues.size();
+
+  float[] phi = new float[n];
+  float[] psi = new float[n];
+
+  for (int i = 0; i < n; i++) {
+    if (i == 0 || i == n - 1) {
+      phi[i] = 90;
+      psi[i] = 90;              
+    } else {
+      res0 = (HashMap)residues.get(i - 1);
+      res1 = (HashMap)residues.get(i);
+      res2 = (HashMap)residues.get(i + 1);
+
+      c0 = (PVector)res0.get("C");
+      n1 = (PVector)res1.get("N");
+      ca1 = (PVector)res1.get("CA"); 
+      c1 = (PVector)res1.get("C");
+      n2 = (PVector)res2.get("N");
+
+      phi[i] = calculateTorsionalAngle(c0, n1, ca1, c1);
+      psi[i] = calculateTorsionalAngle(n1, ca1, c1, n2);
+    }
+  }  
+
+  int firstHelix = 0;
+  int nconsRHelix = 0;
+  int nconsLHelix = 0;
+  int firstStrand = 0;
+  int nconsStrand = 0;
+  for (int i = 0; i < n; i++) {
+    // Right-handed helix      
+    if ((dist(phi[i], psi[i], -60, -45) < 30) && (i < n - 1)) {
+      if (nconsRHelix == 0) firstHelix = i;
+      nconsRHelix++;
+    } 
+    else {
+      if (3 <= nconsRHelix) {
+        for (int k = firstHelix; k < i; k++) {
+          ss[k] = HELIX;
+          handness[k] = RHANDED;                  
+        }
+      }
+      nconsRHelix = 0;
+    }
+
+    // Left-handed helix
+    if ((dist(phi[i], psi[i], +60, +45) < 30) && (i < n - 1)) {
+      if (nconsLHelix == 0) firstHelix = i;
+      nconsLHelix++;
+
+    } else {
+      if (3 <= nconsLHelix) {
+        for (int k = firstHelix; k < i; k++) {
+          ss[k] = HELIX;
+          handness[k] = LHANDED;
+        }
+      }
+      nconsLHelix = 0;
+    }
+
+    // Strand
+    if ((dist(phi[i], psi[i], -110, +130) < 30) && (i < n - 1)) {
+      if (nconsStrand == 0) firstStrand = i;
+      nconsStrand++;
+    } else {
+      if (2 <= nconsStrand) {
+        for (int k = firstStrand; k < i; k++) {
+          ss[k] = STRAND;
+          handness[k] = RHANDED;
+
+        }
+      }
+      nconsStrand = 0;
+    }        
+
+    ss[i] = COIL;
+    handness[i] = RHANDED;
+  }
+}
+
+// Calculates the torsional angle defined by four atoms with positions at0, at1, at2 and at3.
+float calculateTorsionalAngle(PVector at0, PVector at1, PVector at2, PVector at3) {
+  PVector r01 = PVector.sub(at0, at1);
+  PVector r32 = PVector.sub(at3, at2);
+  PVector r12 = PVector.sub(at1, at2);
+
+  PVector p = r12.cross(r01);
+  PVector q = r12.cross(r32);
+  PVector r = r12.cross(q);
+
+  float u = q.dot(q);
+  float v = r.dot(r);    
+
+  float a;
+  if (u <= 0.0 || v <= 0.0) {
+    a = 360.0;
+  } else {
+    float u1 = p.dot(q); // u1 = p * q
+    float v1 = p.dot(r); // v1 = p * r
+
+      u = u1 / sqrt(u);
+    v = v1 / sqrt(v);
+
+    if (abs(u) > 0.01 || abs(v) > 0.01) a = degrees(atan2(v, u));
+    else a = 360.0;
+  }    
+  return a;
+}
+
+void generateSpline(int n, ArrayList vertices) {
+  int ui;
+  float u;
+  PVector v0, v1;
+
+  v1 = new PVector();
+
+  if (n == 0) splineSide1.feval(0, v1); 
+  else if (n == 1) splineCenter.feval(0, v1);
+  else splineSide2.feval(0, v1);
+  vertices.add(new PVector(v1.x, v1.y, v1.z));
+
+  for (ui = 1; ui <= 10; ui ++) {
+    if (ui % uspacing == 0) {
+      u = 0.1 * ui; 
+
+      if (n == 0) splineSide1.feval(u, v1); 
+      else if (n == 1) splineCenter.feval(u, v1); 
+      else splineSide2.feval(u, v1);
+
+      vertices.add(new PVector(v1.x, v1.y, v1.z));            
+    }
+  }
+}
+
+void generateFlatRibbon(ArrayList vertices, ArrayList normals) {
+  PVector CentPoint0, CentPoint1;
+  PVector Sid1Point0, Sid1Point1;
+  PVector Sid2Point0, Sid2Point1;
+  PVector Transversal, Tangent;
+  PVector Normal0, Normal1;
+  int ui;
+  float u;
+
+  CentPoint0 = new PVector();
+  CentPoint1 = new PVector();
+  Sid1Point0 = new PVector();
+  Sid1Point1 = new PVector();
+  Sid2Point0 = new PVector();
+  Sid2Point1 = new PVector();
+  Transversal = new PVector();
+  Tangent = new PVector();
+  Normal0 = new PVector();
+  Normal1 = new PVector();
+
+  // The initial geometry is generated.
+  splineSide1.feval(0, Sid1Point1);
+  splineCenter.feval(0, CentPoint1);
+  splineSide2.feval(0, Sid2Point1);
+
+  // The tangents at the three previous points are the same.
+  splineSide2.deval(0, Tangent);
+
+  // Vector transversal to the ribbon.    
+  Transversal = PVector.sub(Sid1Point1, Sid2Point1);
+
+  // The normal is calculated.
+  Normal1 = Transversal.cross(Tangent);
+  Normal1.normalize();
+
+  for (ui = 1; ui <= 10; ui ++) {
+    if (ui % uspacing == 0) {
+      u = 0.1 * ui;
+
+      // The geometry of the previous iteration is saved.
+      Sid1Point0.set(Sid1Point1);
+      CentPoint0.set(CentPoint1);
+      Sid2Point0.set(Sid2Point1);
+      Normal0.set(Normal1);
+
+      // The new geometry is generated.
+      splineSide1.feval(u, Sid1Point1);
+      splineCenter.feval(u, CentPoint1);
+      splineSide2.feval(u, Sid2Point1);
+
+      // The tangents at the three previous points are the same.
+      splineSide2.deval(u, Tangent);
+      // Vector transversal to the ribbon.
+      Transversal = PVector.sub(Sid1Point1, Sid2Point1);
+      // The normal is calculated.
+      Normal1 = Transversal.cross(Tangent);
+      Normal1.normalize();
+
+      // The (Sid1Point0, Sid1Point1, CentPoint1) triangle is added.
+      vertices.add(new PVector(Sid1Point0.x, Sid1Point0.y, Sid1Point0.z));
+      normals.add(new PVector(Normal0.x, Normal0.y, Normal0.z));
+
+      vertices.add(new PVector(Sid1Point1.x, Sid1Point1.y, Sid1Point1.z));
+      normals.add(new PVector(Normal1.x, Normal1.y, Normal1.z));
+
+      vertices.add(new PVector(CentPoint1.x, CentPoint1.y, CentPoint1.z));
+      normals.add(new PVector(Normal1.x, Normal1.y, Normal1.z));      
+
+      // The (Sid1Point0, CentPoint1, CentPoint0) triangle is added.
+      vertices.add(new PVector(Sid1Point0.x, Sid1Point0.y, Sid1Point0.z));
+      normals.add(new PVector(Normal0.x, Normal0.y, Normal0.z));
+
+      vertices.add(new PVector(CentPoint1.x, CentPoint1.y, CentPoint1.z));
+      normals.add(new PVector(Normal1.x, Normal1.y, Normal1.z));      
+
+      vertices.add(new PVector(CentPoint0.x, CentPoint0.y, CentPoint0.z));
+      normals.add(new PVector(Normal0.x, Normal0.y, Normal0.z));
+      
+      // (Sid2Point0, Sid2Point1, CentPoint1) triangle is added.
+      vertices.add(new PVector(Sid2Point0.x, Sid2Point0.y, Sid2Point0.z));
+      normals.add(new PVector(Normal0.x, Normal0.y, Normal0.z));
+
+      vertices.add(new PVector(Sid2Point1.x, Sid2Point1.y, Sid2Point1.z));
+      normals.add(new PVector(Normal1.x, Normal1.y, Normal1.z));                        
+
+      vertices.add(new PVector(CentPoint1.x, CentPoint1.y, CentPoint1.z));
+      normals.add(new PVector(Normal1.x, Normal1.y, Normal1.z));
+
+      // (Sid2Point0, CentPoint1, CentPoint0) triangle is added.
+      vertices.add(new PVector(Sid2Point0.x, Sid2Point0.y, Sid2Point0.z));
+      normals.add(new PVector(Normal0.x, Normal0.y, Normal0.z));                        
+     
+      vertices.add(new PVector(CentPoint1.x, CentPoint1.y, CentPoint1.z));
+      normals.add(new PVector(Normal1.x, Normal1.y, Normal1.z));
+
+      vertices.add(new PVector(CentPoint0.x, CentPoint0.y, CentPoint0.z));
+      normals.add(new PVector(Normal0.x, Normal0.y, Normal0.z));
+    }
+  }
+}
+
+/******************************************************************************
+ * The code in the following three functions is based in the method introduced 
+ * in this paper:
+ * "Algorithm for ribbon models of proteins."
+ * Authors: Mike Carson and Charles E. Bugg
+ * Published in: J.Mol.Graphics 4, pp. 121-122 (1986)
+ ******************************************************************************/
+
+// Shifts the control points one place to the left.
+void shiftControlPoints() {
+  splineSide1.shiftBSplineCPoints();
+  splineCenter.shiftBSplineCPoints();
+  splineSide2.shiftBSplineCPoints();
+}
+
+// Adds a new control point to the arrays CPCenter, CPRight and CPLeft
+void addControlPoints(PVector ca0, PVector ox0, PVector ca1, int ss, int handness) {
+  PVector A, B, C, D, p0, cpt0, cpt1, cpt2;
+
+  A = PVector.sub(ca1, ca0);
+  B = PVector.sub(ox0, ca0);
+
+  // Vector normal to the peptide plane (pointing outside in the case of the
+  // alpha helix).
+  C = A.cross(B);
+
+  // Vector contained in the peptide plane (perpendicular to its direction).
+  D = C.cross(A);
+
+  // Normalizing vectors.
+  C.normalize();
+  D.normalize();
+
+  // Flipping test (to avoid self crossing in the strands).
+  if ((ss != HELIX) && (90.0 < degrees(PVector.angleBetween(flipTestV, D)))) {
+    // Flip detected. The plane vector is inverted.
+    D.mult(-1.0);
+  }
+
+  // The central control point is constructed.
+  cpt0 = linearComb(0.5, ca0, 0.5, ca1);
+  splineCenter.setCPoint(3, cpt0);
+
+  if (ss == HELIX) {
+    // When residue i is contained in a helix, the control point is moved away
+    // from the helix axis, along the C direction. 
+    p0 = new PVector();
+    splineCenter.getCPoint(3, p0);
+    cpt0 = linearComb(1.0, p0, handness * helixDiam, C);
+    splineCenter.setCPoint(3, cpt0);
+  }
+
+  // The control points for the side ribbons are constructed.
+  cpt1 = linearComb(1.0, cpt0, +ribbonWidth[ss], D);
+  splineSide1.setCPoint(3, cpt1);
+
+  cpt2 = linearComb(1.0, cpt0, -ribbonWidth[ss], D);
+  splineSide2.setCPoint(3, cpt2);
+
+  // Saving the plane vector (for the flipping test in the next call).
+  flipTestV.set(D);
+}
+
+void constructControlPoints(ArrayList residues, int res, int ss, int handness) {
+  PVector ca0, ox0, ca1;
+  PVector p0, p1, p2, p3;
+
+  p1 = new PVector();        
+  p2 = new PVector();
+  p3 = new PVector();
+
+  HashMap res0, res1;
+
+  res0 = res1 = null;
+  if (res == 0) {
+    // The control points 2 and 3 are created.
+    flipTestV.set(0, 0, 0);
+
+    res0 = (HashMap)residues.get(res);
+    res1 = (HashMap)residues.get(res + 1);
+    ca0 = (PVector)res0.get("CA");
+    ox0 = (PVector)res0.get("O");
+    ca1 = (PVector)res1.get("CA");        
+    addControlPoints(ca0, ox0, ca1, ss, handness);
+    splineSide1.copyCPoints(3, 2);
+    splineCenter.copyCPoints(3, 2);
+    splineSide2.copyCPoints(3, 2);        
+
+    res0 = (HashMap)residues.get(res + 1);
+    res1 = (HashMap)residues.get(res + 2);
+    ca0 = (PVector)res0.get("CA");
+    ox0 = (PVector)res0.get("O");
+    ca1 = (PVector)res1.get("CA"); 
+    addControlPoints(ca0, ox0, ca1, ss, handness);
+
+    // We still need the two first control points.
+    // Moving backwards along the cp_center[2] - cp_center[3] direction.
+    splineCenter.getCPoint(2, p2);
+    splineCenter.getCPoint(3, p3);
+
+    p1 = linearComb(2.0, p2, -1, p3);
+    splineCenter.setCPoint(1, p1);        
+    splineSide1.setCPoint(1, linearComb(1.0, p1, +ribbonWidth[ss], flipTestV)); 
+    splineSide2.setCPoint(1, linearComb(1.0, p1, -ribbonWidth[ss], flipTestV));        
+
+    p0 = linearComb(2.0, p1, -1, p2);
+    splineCenter.setCPoint(0, p0);
+    splineSide1.setCPoint(0, linearComb(1.0, p0, +ribbonWidth[ss], flipTestV));
+    splineSide2.setCPoint(0, linearComb(1.0, p0, -ribbonWidth[ss], flipTestV));
+  } else {
+    shiftControlPoints();
+    if ((residues.size() - 1 == res) || (residues.size() - 2 == res)) { 
+      // Moving forward along the cp_center[1] - cp_center[2] direction.
+      splineCenter.getCPoint(1, p1);             
+      splineCenter.getCPoint(2, p2);
+
+      p3 = linearComb(2.0, p2, -1, p1);
+      splineCenter.setCPoint(3, p3);
+      splineSide1.setCPoint(3, linearComb(1.0, p3, +ribbonWidth[ss], flipTestV));
+      splineSide2.setCPoint(3, linearComb(1.0, p3, -ribbonWidth[ss], flipTestV));
+    } else {
+      res0 = (HashMap)residues.get(res + 1);
+      res1 = (HashMap)residues.get(res + 2);
+      ca0 = (PVector)res0.get("CA");
+      ox0 = (PVector)res0.get("O");
+      ca1 = (PVector)res1.get("CA");        
+      addControlPoints(ca0, ox0, ca1, ss, handness);
+    }
+  }
+  splineSide1.updateMatrix3();
+  splineCenter.updateMatrix3();
+  splineSide2.updateMatrix3();
+}
+
+PVector linearComb(float scalar0, PVector vector0, float scalar1, PVector vector1) {
+  return PVector.add(PVector.mult(vector0, scalar0), PVector.mult(vector1, scalar1));
+}

java/examples/OpenGL/Advanced/Ribbons/PDB.pde

@@ -0,0 +1,115 @@
+void readPDB(String filename) {
+  String strLines[];
+
+  float xmin, xmax, ymin, ymax, zmin, zmax;
+
+  String xstr, ystr, zstr;
+  float x, y, z;
+  int res, res0;
+  int nmdl;
+  String atstr, resstr;
+
+  PShape model;
+  ArrayList atoms;  
+  ArrayList residues;
+  HashMap residue;
+  PVector v;
+  String s;
+  strLines = loadStrings(filename);
+
+  models = new ArrayList();
+
+  xmin = ymin = zmin = 10000;
+  xmax = ymax = zmax = -10000; 
+
+  atoms = null;
+  residues = null;
+  residue = null;
+  model = null;
+  res0 = -1;    
+  nmdl = -1;
+  for (int i = 0; i < strLines.length; i++) {
+    s = strLines[i];
+
+    if (s.startsWith("MODEL") || (s.startsWith("ATOM") && res0 == -1)) {
+      nmdl++;
+
+      res0 = -1;
+
+      atoms = new ArrayList();
+      residues = new ArrayList();
+    }
+
+    if (s.startsWith("ATOM")) {
+      atstr = s.substring(12, 15);
+      atstr = atstr.trim();
+      resstr = s.substring(22, 26);
+      resstr = resstr.trim();
+      res = parseInt(resstr);
+
+      xstr = s.substring(30, 37);
+      xstr = xstr.trim();
+      ystr = s.substring(38, 45);
+      ystr = ystr.trim();            
+      zstr = s.substring(46, 53);
+      zstr = zstr.trim();
+
+      x = scaleFactor * parseFloat(xstr);
+      y = scaleFactor * parseFloat(ystr);
+      z = scaleFactor * parseFloat(zstr);            
+      v = new PVector(x, y, z);
+
+      xmin = min(xmin, x);
+      xmax = max(xmax, x);
+    
+      ymin = min(ymin, y);
+      ymax = max(ymax, y);
+
+      zmin = min(zmin, z);
+      zmax = max(zmax, z);
+
+      atoms.add(v);
+
+      if (res0 != res) {
+        if (residue != null) residues.add(residue);
+        residue = new HashMap();
+      }
+      residue.put(atstr, v);
+
+      res0 = res;
+    }
+
+    if (s.startsWith("ENDMDL") || s.startsWith("TER")) {
+      if (residue != null) residues.add(residue);
+
+      createRibbonModel(residues, model, models);
+      float rgyr = calculateGyrRadius(atoms);
+
+      res0 = -1;
+      residue = null;  
+      atoms = null;
+      residues = null;
+    }
+  }
+
+  if (residue != null) {
+    if (residue != null) residues.add(residue);
+
+    createRibbonModel(residues, model, models);
+    float rgyr = calculateGyrRadius(atoms);
+
+    atoms = null;
+    residues = null;
+  }
+
+  // Centering models at (0, 0, 0).
+  float dx = -0.5f * (xmin + xmax);
+  float dy = -0.5f * (ymin + ymax);
+  float dz = -0.5f * (zmin + zmax);
+  for (int n = 0; n < models.size(); n++) {
+    model = (PShape) models.get(n);
+    model.translate(dx, dy, dz);
+  }
+
+  println("Loaded PDB file with " + models.size() + " models.");
+}

java/examples/OpenGL/Advanced/Ribbons/Ribbons.pde

@@ -0,0 +1,50 @@
+// Ribbons, by Andres Colubri
+// ArcBall class by Ariel, V3ga and Robert Hodgin (flight404)
+// This sketch loads 3D atomic coordinates of a protein molecule
+// from a file in PDB format (http://www.pdb.org/) and displays 
+// the structure using a ribbon representation.
+
+String pdbFile = "4HHB.pdb"; // PDB file to read
+//String pdbFile = "1CBS.pdb";
+//String pdbFile = "2POR.pdb";
+
+// Some parameters to control the visual appearance:
+float scaleFactor = 5;          // Size factor
+int renderMode = 1;             // 0 = lines, 1 = flat ribbons
+int ribbonDetail = 4;           // Ribbon detail: from 1 (lowest) to 4 (highest)
+float helixDiam = 10;           // Helix diameter.
+int[] ribbonWidth = {10, 7, 2}; // Ribbon widths for helix, strand and coil
+color ribbonColor = color(20, 30, 200, 255); // Ribbon color
+
+// All the molecular models read from the PDB file (it could contain more than one)
+ArrayList models;
+
+Arcball arcball;
+
+void setup() {
+  size(800, 600, P3D);
+  
+  arcball = new Arcball(width/2, height/2, 600);  
+  readPDB(pdbFile);
+}
+
+void draw() {
+  background(0);
+   
+  lights();
+   
+  translate(width/2, height/2, 200);
+  arcball.run();
+
+  for (int i = 0; i < models.size(); i++) {
+    shape((PShape)models.get(i));
+  }
+}
+
+void mousePressed(){
+  arcball.mousePressed();
+}
+
+void mouseDragged(){
+  arcball.mouseDragged();
+}

java/examples/OpenGL/Advanced/Trefoil/Surface.pde

@@ -0,0 +1,104 @@
+// Code to draw a trefoil knot surface, with normals and texture 
+// coordinates.
+// Adapted from the parametric equations example by Philip Rideout:
+// http://iphone-3d-programming.labs.oreilly.com/ch03.html
+
+// This function draws a trefoil knot surface as a triangle mesh derived
+// from its parametric equation.
+PShape createTrefoil(float s, int ny, int nx, PImage tex) {
+  PVector p0, p1, p2;
+  PVector n0, n1, n2;
+  float u0, u1, v0, v1;
+  
+  PShape obj = createShape(TRIANGLES);
+  obj.texture(tex);
+    
+  for (int j = 0; j < nx; j++) {
+    u0 = float(j) / nx;
+    u1 = float(j + 1) / nx;
+    for (int i = 0; i < ny; i++) {
+      v0 = float(i) / ny;
+      v1 = float(i + 1) / ny;
+      
+      p0 = evalPoint(u0, v0);
+      n0 = evalNormal(u0, v0);
+      
+      p1 = evalPoint(u0, v1);
+      n1 = evalNormal(u0, v1);
+      
+      p2 = evalPoint(u1, v1);
+      n2 = evalNormal(u1, v1);
+
+      // Triangle p0-p1-p2      
+      obj.normal(n0.x, n0.y, n0.z);
+      obj.vertex(s * p0.x, s * p0.y, s * p0.z, u0, v0);      
+      obj.normal(n1.x, n1.y, n1.z);
+      obj.vertex(s * p1.x, s * p1.y, s * p1.z, u0, v1);
+      obj.normal(n2.x, n2.y, n2.z);
+      obj.vertex(s * p2.x, s * p2.y, s * p2.z, u1, v1);
+
+      p1 = evalPoint(u1, v0);
+      n1 = evalNormal(u1, v0);
+
+      // Triangle p0-p2-p1      
+      obj.normal(n0.x, n0.y, n0.z);
+      obj.vertex(s * p0.x, s * p0.y, s * p0.z, u0, v0);      
+      obj.normal(n2.x, n2.y, n2.z);
+      obj.vertex(s * p2.x, s * p2.y, s * p2.z, u1, v1);
+      obj.normal(n1.x, n1.y, n1.z);
+      obj.vertex(s * p1.x, s * p1.y, s * p1.z, u1, v0);      
+    }
+  }
+  obj.end();
+  return obj;
+}
+
+// Evaluates the surface normal corresponding to normalized 
+// parameters (u, v)
+PVector evalNormal(float u, float v) {
+  // Compute the tangents and their cross product.
+  PVector p = evalPoint(u, v);
+  PVector tangU = evalPoint(u + 0.01, v);
+  PVector tangV = evalPoint(u, v + 0.01);
+  tangU.sub(p);
+  tangV.sub(p);
+  
+  PVector normUV = tangV.cross(tangU);
+  normUV.normalize();
+  return normUV;
+}
+
+// Evaluates the surface point corresponding to normalized 
+// parameters (u, v)
+PVector evalPoint(float u, float v) {
+  float a = 0.5;
+  float b = 0.3;
+  float c = 0.5;
+  float d = 0.1;
+  float s = TWO_PI * u;
+  float t = (TWO_PI * (1 - v)) * 2;  
+        
+  float r = a + b * cos(1.5 * t);
+  float x = r * cos(t);
+  float y = r * sin(t);
+  float z = c * sin(1.5 * t);
+        
+  PVector dv = new PVector();
+  dv.x = -1.5 * b * sin(1.5 * t) * cos(t) -
+         (a + b * cos(1.5 * t)) * sin(t);
+  dv.y = -1.5 * b * sin(1.5 * t) * sin(t) +
+         (a + b * cos(1.5 * t)) * cos(t);
+  dv.z = 1.5 * c * cos(1.5 * t);
+        
+  PVector q = dv;      
+  q.normalize();
+  PVector qvn = new PVector(q.y, -q.x, 0);
+  qvn.normalize();
+  PVector ww = q.cross(qvn);
+        
+  PVector pt = new PVector();
+  pt.x = x + d * (qvn.x * cos(s) + ww.x * sin(s));
+  pt.y = y + d * (qvn.y * cos(s) + ww.y * sin(s));
+  pt.z = z + d * ww.z * sin(s);
+  return pt;
+}

java/examples/OpenGL/Advanced/Trefoil/Trefoil.pde

@@ -0,0 +1,42 @@
+// Trefoil, by Andres Colubri
+// A parametric surface is textured procedurally
+// by drawing on an offscreen PGraphics surface.
+
+PGraphics pg;
+PShape trefoil;
+
+void setup() {
+  size(280, 400, P3D);
+  
+  textureMode(NORMAL);
+  noStroke();
+
+  // Creating offscreen surface for 3D rendering.
+  pg = createGraphics(32, 512, P3D);
+  pg.beginDraw();
+  pg.background(0, 0);
+  pg.noStroke();
+  pg.fill(255, 0, 0, 200);
+  pg.endDraw(); 
+
+  // Saving trefoil surface into a PShape3D object
+  trefoil = createTrefoil(250, 60, 15, pg);
+}
+
+void draw() {
+  background(0);
+  
+  pg.beginDraw();    
+  pg.ellipse(random(pg.width), random(pg.height), 4, 4);
+  pg.endDraw(); 
+
+  ambient(250, 250, 250);
+  pointLight(255, 255, 255, 0, 0, 200);
+     
+  pushMatrix();
+  translate(width/2, height/2, -200);
+  rotateX(frameCount * PI / 500);
+  rotateY(frameCount * PI / 500);      
+  shape(trefoil);
+  popMatrix();
+}

java/examples/OpenGL/Advanced/Wiggling/Wiggling.pde

@@ -0,0 +1,283 @@
+// Press 'w' to start wiggling, space to restore
+// original positions.
+
+PShape cube;
+float cubeSize = 200;
+float circleRad = 50;
+int circleRes = 20;
+float noiseMag = 1;
+
+boolean wiggling = false;
+
+void setup() {
+  size(400, 400, P3D);    
+
+  createCube();
+}
+
+void draw() {
+  background(0);
+
+  translate(width/2, height/2);
+  rotateX(frameCount * 0.01f);
+  rotateY(frameCount * 0.01f);
+
+  shape(cube);
+
+  if (wiggling) {
+    PVector pos = null;
+    for (int i = 0; i < cube.getChildCount(); i++) {
+      PShape face = cube.getChild(i);
+      for (int j = 0; j < face.getVertexCount(); j++) {
+        pos = face.getVertex(j, pos);
+        pos.x += random(-noiseMag/2, +noiseMag/2);
+        pos.y += random(-noiseMag/2, +noiseMag/2);
+        pos.z += random(-noiseMag/2, +noiseMag/2);
+        face.setVertex(j, pos.x, pos.y, pos.z);
+      }
+    }
+  }
+
+  println(frameRate);
+}
+
+public void keyPressed() {
+  if (key == 'w') {
+    wiggling = !wiggling;
+  } else if (key == ' ') {
+    restoreCube();
+  } else if (key == '1') {
+    cube.strokeWeight(1);
+  } else if (key == '2') {
+    cube.strokeWeight(5);
+  } else if (key == '3') {
+    cube.strokeWeight(10);
+  }
+}
+
+void createCube() {
+  cube = createShape(GROUP);  
+
+  PShape face;
+
+  // Front face         
+  face = createShape(POLYGON);
+  face.stroke(255, 0, 0);
+  face.fill(255);
+  face.beginContour();
+  face.vertex(-cubeSize/2, -cubeSize/2, +cubeSize/2);
+  face.vertex(+cubeSize/2, -cubeSize/2, +cubeSize/2);
+  face.vertex(+cubeSize/2, +cubeSize/2, +cubeSize/2);
+  face.vertex(-cubeSize/2, +cubeSize/2, +cubeSize/2);
+  face.endContour();
+  face.beginContour();
+  for (int i = 0; i < circleRes; i++) {
+    float angle = TWO_PI * i / circleRes;
+    float x = circleRad * sin(angle);
+    float y = circleRad * cos(angle);
+    float z = +cubeSize/2;
+    face.vertex(x, y, z);
+  }
+  face.endContour();
+  face.end(CLOSE);
+  cube.addChild(face);
+
+  // Back face
+  face = createShape(POLYGON);
+  face.stroke(255, 0, 0);
+  face.fill(255);
+  face.beginContour();
+  face.vertex(+cubeSize/2, -cubeSize/2, -cubeSize/2);
+  face.vertex(-cubeSize/2, -cubeSize/2, -cubeSize/2);
+  face.vertex(-cubeSize/2, +cubeSize/2, -cubeSize/2);
+  face.vertex(+cubeSize/2, +cubeSize/2, -cubeSize/2);
+  face.endContour();
+  face.beginContour();
+  for (int i = 0; i < circleRes; i++) {
+    float angle = TWO_PI * i / circleRes;
+    float x = circleRad * sin(angle);
+    float y = circleRad * cos(angle);
+    float z = -cubeSize/2;
+    face.vertex(x, y, z);
+  }
+  face.endContour();
+  face.end(CLOSE);
+  cube.addChild(face);
+
+  // Right face
+  face = createShape(POLYGON);
+  face.stroke(255, 0, 0);
+  face.fill(255);
+  face.beginContour();
+  face.vertex(+cubeSize/2, -cubeSize/2, +cubeSize/2);
+  face.vertex(+cubeSize/2, -cubeSize/2, -cubeSize/2);
+  face.vertex(+cubeSize/2, +cubeSize/2, -cubeSize/2);
+  face.vertex(+cubeSize/2, +cubeSize/2, +cubeSize/2);
+  face.endContour();
+  face.beginContour();
+  for (int i = 0; i < circleRes; i++) {
+    float angle = TWO_PI * i / circleRes;
+    float x = +cubeSize/2;
+    float y = circleRad * sin(angle);
+    float z = circleRad * cos(angle);
+    face.vertex(x, y, z);
+  }
+  face.endContour();
+  face.end(CLOSE);
+  cube.addChild(face);
+
+  // Left face
+  face = createShape(POLYGON);
+  face.stroke(255, 0, 0);
+  face.fill(255);
+  face.beginContour();
+  face.vertex(-cubeSize/2, -cubeSize/2, -cubeSize/2);
+  face.vertex(-cubeSize/2, -cubeSize/2, +cubeSize/2);
+  face.vertex(-cubeSize/2, +cubeSize/2, +cubeSize/2);
+  face.vertex(-cubeSize/2, +cubeSize/2, -cubeSize/2);
+  face.endContour();
+  face.beginContour();
+  for (int i = 0; i < circleRes; i++) {
+    float angle = TWO_PI * i / circleRes;
+    float x = -cubeSize/2;
+    float y = circleRad * sin(angle);
+    float z = circleRad * cos(angle);
+    face.vertex(x, y, z);
+  }
+  face.endContour();
+  face.end(CLOSE);
+  cube.addChild(face);
+
+  // Top face
+  face = createShape(POLYGON);
+  face.stroke(255, 0, 0);
+  face.fill(255);
+  face.beginContour();
+  face.vertex(-cubeSize/2, +cubeSize/2, +cubeSize/2);
+  face.vertex(+cubeSize/2, +cubeSize/2, +cubeSize/2);
+  face.vertex(+cubeSize/2, +cubeSize/2, -cubeSize/2);
+  face.vertex(-cubeSize/2, +cubeSize/2, -cubeSize/2);
+  face.endContour();
+  face.beginContour();
+  for (int i = 0; i < circleRes; i++) {
+    float angle = TWO_PI * i / circleRes;
+    float x = circleRad * sin(angle);
+    float y = +cubeSize/2;
+    float z = circleRad * cos(angle);
+    face.vertex(x, y, z);
+  }
+  face.endContour();
+  face.end(CLOSE);
+  cube.addChild(face);     
+
+  // Bottom face
+  face = createShape(POLYGON);
+  face.stroke(255, 0, 0);
+  face.fill(255);
+  face.beginContour();
+  face.vertex(+cubeSize/2, -cubeSize/2, +cubeSize/2);
+  face.vertex(-cubeSize/2, -cubeSize/2, +cubeSize/2);
+  face.vertex(-cubeSize/2, -cubeSize/2, -cubeSize/2);
+  face.vertex(+cubeSize/2, -cubeSize/2, -cubeSize/2);
+  face.endContour();
+  face.beginContour();
+  for (int i = 0; i < circleRes; i++) {
+    float angle = TWO_PI * i / circleRes;
+    float x = circleRad * sin(angle);
+    float y = -cubeSize/2;
+    float z = circleRad * cos(angle);
+    face.vertex(x, y, z);
+  }
+  face.endContour();
+  face.end(CLOSE);
+  cube.addChild(face);
+}
+
+void restoreCube() {
+  PShape face;
+
+  // Front face
+  face = cube.getChild(0);
+  face.setVertex(0, -cubeSize/2, -cubeSize/2, +cubeSize/2);
+  face.setVertex(1, +cubeSize/2, -cubeSize/2, +cubeSize/2);
+  face.setVertex(2, +cubeSize/2, +cubeSize/2, +cubeSize/2);
+  face.setVertex(3, -cubeSize/2, +cubeSize/2, +cubeSize/2);
+  for (int i = 0; i < circleRes; i++) {
+    float angle = TWO_PI * i / circleRes;
+    float x = circleRad * sin(angle);
+    float y = circleRad * cos(angle);
+    float z = +cubeSize/2;
+    face.setVertex(4 + i, x, y, z);
+  }
+
+  // Back face
+  face = cube.getChild(1);
+  face.setVertex(0, +cubeSize/2, -cubeSize/2, -cubeSize/2);
+  face.setVertex(1, -cubeSize/2, -cubeSize/2, -cubeSize/2);
+  face.setVertex(2, -cubeSize/2, +cubeSize/2, -cubeSize/2);
+  face.setVertex(3, +cubeSize/2, +cubeSize/2, -cubeSize/2);
+  for (int i = 0; i < circleRes; i++) {
+    float angle = TWO_PI * i / circleRes;
+    float x = circleRad * sin(angle);
+    float y = circleRad * cos(angle);
+    float z = -cubeSize/2;
+    face.setVertex(4 + i, x, y, z);
+  }
+
+  // Right face
+  face = cube.getChild(2);
+  face.setVertex(0, +cubeSize/2, -cubeSize/2, +cubeSize/2);
+  face.setVertex(1, +cubeSize/2, -cubeSize/2, -cubeSize/2);
+  face.setVertex(2, +cubeSize/2, +cubeSize/2, -cubeSize/2);
+  face.setVertex(3, +cubeSize/2, +cubeSize/2, +cubeSize/2);
+  for (int i = 0; i < circleRes; i++) {
+    float angle = TWO_PI * i / circleRes;
+    float x = +cubeSize/2;
+    float y = circleRad * sin(angle);
+    float z = circleRad * cos(angle);
+    face.setVertex(4 + i, x, y, z);
+  }
+
+  // Left face
+  face = cube.getChild(3);
+  face.setVertex(0, -cubeSize/2, -cubeSize/2, -cubeSize/2);
+  face.setVertex(1, -cubeSize/2, -cubeSize/2, +cubeSize/2);
+  face.setVertex(2, -cubeSize/2, +cubeSize/2, +cubeSize/2);
+  face.setVertex(3, -cubeSize/2, +cubeSize/2, -cubeSize/2);
+  for (int i = 0; i < circleRes; i++) {
+    float angle = TWO_PI * i / circleRes;
+    float x = -cubeSize/2;
+    float y = circleRad * sin(angle);
+    float z = circleRad * cos(angle);
+    face.setVertex(4 + i, x, y, z);
+  }    
+
+  // Top face
+  face = cube.getChild(4);
+  face.setVertex(0, -cubeSize/2, +cubeSize/2, +cubeSize/2);
+  face.setVertex(1, +cubeSize/2, +cubeSize/2, +cubeSize/2);
+  face.setVertex(2, +cubeSize/2, +cubeSize/2, -cubeSize/2);
+  face.setVertex(3, -cubeSize/2, +cubeSize/2, -cubeSize/2);
+  for (int i = 0; i < circleRes; i++) {
+    float angle = TWO_PI * i / circleRes;
+    float x = circleRad * sin(angle);
+    float y = +cubeSize/2;
+    float z = circleRad * cos(angle);
+    face.setVertex(4 + i, x, y, z);
+  }    
+
+  // Bottom face
+  face = cube.getChild(5);
+  face.setVertex(0, +cubeSize/2, -cubeSize/2, +cubeSize/2);
+  face.setVertex(1, -cubeSize/2, -cubeSize/2, +cubeSize/2);
+  face.setVertex(2, -cubeSize/2, -cubeSize/2, -cubeSize/2);
+  face.setVertex(3, +cubeSize/2, -cubeSize/2, -cubeSize/2);
+  for (int i = 0; i < circleRes; i++) {
+    float angle = TWO_PI * i / circleRes;
+    float x = circleRad * sin(angle);
+    float y = -cubeSize/2;
+    float z = circleRad * cos(angle);
+    face.setVertex(4 + i, x, y, z);
+  }
+}
+