/* mandel.c - Mandelbrot image generator
 * Written 6 Nov 2008 by John T. Wodder II
 * Last edited 6 Nov 2008 by John Wodder
 */

#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <unistd.h>

_Bool escapes(double a, double b);

int main(int argc, char** argv) {
 double minx = -2.0, maxx = 0.5, miny = -1.0, maxy = 1.0;
 int precision = 100, ch;
 while ((ch = getopt(argc, argv, "x:y:X:Y:p:")) != -1) {
  switch (ch) {
   case 'x': minx = strtod(optarg, NULL); break;
   case 'y': miny = strtod(optarg, NULL); break;
   case 'X': maxx = strtod(optarg, NULL); break;
   case 'Y': maxy = strtod(optarg, NULL); break;
   case 'p': precision = strtol(optarg, NULL, 0); break;
   default:
    fprintf(stderr, "Usage: mandel [-x minx] [-y miny] [-X maxx] [-Y maxy]"
     " [-p precision] [outfile]\n");
    return 2;
  }
 }
 char* filename = (optind < argc) ? argv[optind] : "mandelbrot.eps";
 FILE* eps = fopen(filename, "w");
 if (eps == NULL) {
  fprintf(stderr, "mandel: %s: ", filename);
  perror(NULL);
  return 5;
 }
 fprintf(eps, "%%!PS-Adobe-3.0 EPSF-3.0\n%%%%BoundingBox: %.1f %ld %.1f %ld\n",
  minx * precision - 0.5, lround(miny * precision), maxx * precision + 0.5,
  lround(maxy * precision + 1));
 for (double x = minx; x <= maxx; x += 1.0 / precision) {
  for (double y = miny; y <= maxy; y += 1.0 / precision) {
   if (!escapes(x, y)) {
    fprintf(eps, "%ld %ld moveto 0 1 rlineto\n", lround(x * precision),
     lround(y * precision));
   }
  }
 }
 fprintf(eps, "stroke showpage\n");
 fclose(eps);
 return 0;
}

_Bool escapes(double a, double b) {
 if (hypot(a, b) > 2.0) return 1;
 double x = 0.0, y = 0.0;
 for (int i=0; i<100; i++) {
  double xTmp = x;
  x = x*x - y*y + a;
  y = 2 * xTmp * y + b;
  if (hypot(x, y) > 2.0) return 1;
 }
 return 0;
}