/* built with Studio Sketchpad: 
 *   https://sketchpad.cc
 * 
 * observe the evolution of this sketch: 
 *   https://studio.sketchpad.cc/sp/pad/view/ro.I1AGekjngta/rev.867
 * 
 * authors: 
 *   Nick Bennett

 * license (unless otherwise specified): 
 *   creative commons attribution-share alike 3.0 license.
 *   https://creativecommons.org/licenses/by-sa/3.0/ 
 */ 



float min_y = -2;
float max_y = 2;
float min_x = -2;
float max_x = 2;

int box_size = 4;

float max_iterations = 30;
// How many actual shades of grey are possible to display on the Kindle
float grey_depth = 16;

void setup() {
  size(400, 400);
  noLoop();
}

void draw() {
// This is very unoptimized code; making this depend in real time on the mouse
// position would crash your browser
//  Complex c = new Complex(map(mouseX, 0, width, min_x, max_x), 
//                          map(mouseY, 0, height, min_y, max_y));
  Complex c = new Complex(random(1), random(1));
  draw_juliaset(c);
}

void mouseClicked() {
  Complex c = new Complex(map(mouseX, 0, width, min_x, max_x), 
                          map(mouseY, 0, height, min_y, max_y));
  draw_juliaset(c);
}

void draw_juliaset(Complex c) {
  for(int y_screen = 0; y_screen < height; y_screen += box_size) {
    float y = map(y_screen, 0, height, min_y, max_y);
    for(int x_screen = 0; x_screen < width; x_screen += box_size) {
      float x = map(x_screen, 0, width, min_x, max_x);
      // For the julia set, c is pre-set for all points and z_0 is the original point in the complex plane
      // For the mandelbrot set, c is set as the original complex point and z_0 = 0 (0 + 0i)
      Complex z = new Complex(x, y);
      int iterations_passed = 0;
      while((z.magnitude_squared() < 4) && (iterations_passed < max_iterations)) {
        // I wish I could overload arithmetic operators for the Complex class
        // Then we could have a simple and pretty z = z*z + c
        // Change the exponent of z here for interesting results
        z = z.pow(2);
        z = z.plus(c);
        iterations_passed += 1;
      }
      int grey_shade = int(map(iterations_passed, 0, max_iterations, 0, grey_depth));
      color pixel_color = color(map(grey_shade, 0, grey_depth, 0, 255));
      fill(pixel_color);
      stroke(pixel_color);
//      noStroke();
      rect(x_screen, y_screen, box_size, box_size);
    }
  }
}

class Complex {
  // We are representing a complex number
  // a = x + yi where i is the square root
  // of -1
  float x;
  float y;
  
  Complex(float new_x, float new_y) {
    x = new_x;
    y = new_y;
  }
  
  Complex(Complex target) {
    x = target.x;
    y = target.y;
  }
  
  float magnitude() {
    // Taking a complex number as a vector in the plane, this is the magnitude of that vector
    float mag_val = sqrt(x*x + y*y);
    return mag_val;
  }
  
  float magnitude_squared() {
    float mag_squared = x*x + y*y;
    return mag_squared;    
  }
  
  float distance(Complex target) {
    // Calculate the distance between this
    // complex point and the target
    Complex t = new Complex(target.x-x, target.y-y);
    return t.abs();
  }

  Complex pow(int n) {
      // This only handles n >= 0
      if(n==0) return 1;
      if(n==1) return this;
      if(n==2) return this.times(this);
      if(n%2==0) {
          Complex c = new Complex(this.pow(n/2));
          return c.pow(2);
      }
      if(n%2==1) {
          return this.times(this.pow(n-1));
      }
  }

  Complex times(Complex a) {
    // (x+yi)*(a.x+a.yi) = x*a.x-y*a.y + (a.x*y+x*a.y)i
    float new_x = x*a.x - y*a.y;
    float new_y = x*a.y + y*a.x;
    Complex product = new Complex(new_x, new_y);
    return product;
  }
  
  Complex plus(Complex a) {
    return new Complex(x+a.x, y+a.y);
  }
  
  Complex minus(Complex b) {
    return new Complex(x-b.x, y-b.y);
  }
  
  Complex div_scalar(float n) {
    // not self-modifying
    Complex ret_val = new Complex(x/n, y/n);
    return ret_val;
  }
  
  void div_scalar_self(float n) {
    // IS self-modifying
    x /= n;
    y /= n;
  }
  
  Complex inverse() {
    // (x+yi)*(x-yi) = x*x+y*y
    // 1/(x+yi) = (x-yi)/(x*x+y*y)
    float denominator = x*x + y*y;
    Complex inverted = new Complex(x, -y);
    return inverted.div_scalar(denominator);
  }  
}