/* MCS 572 Spring 2006 simple program to plot a Mandelbrot set as a ps file */

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

#define XOFFSET 54
#define YOFFSET 288
#define SCALE 432

int iterate ( double x, double y );
/* Returns the number of iterations for z^2 + c to grow larger than 2,
 * for c = x + i*y, where i = sqrt(-1), start at z = 0. */

void write_postscript_plot
( FILE *f, int rows, int columns, double a, double b, double c, double d );
/*
 * Writes a plot of a Mandelbrot set to the file in postscript format.
 * The plot is generated as a matrix of grayscales with the specified
 * number of rows and columns and written directly into the file.
 * To save toner for the default plot, the grayscales are inverted.
 *
 * The input parameters have the following meaning:
 *   f        pointer to a file which must be opened for writing;
 *   rows     horizontal resolution of the plot;
 *   columns  vertical resolution of the plot as a matrix of gray scales;
 *   a        start of the range for the real part of z = x + i*y;
 *   b        end of the range for the real part of z = x + i*y;
 *   c        start of the range for the imaginary part of z = x + i*y;
 *   d        end of the range for the imaginary part of z = x + i*y. */

int main ( int argc, char* argv[] )
{
   int interactive = 0; /* set to 1 for interactive version */
   int rows,columns;
   double a,b,c,d;
   char filename[80];
   FILE *outputfile;

   if(interactive==1)
   {
      printf("Give the number of rows : "); scanf("%d",&rows);
      printf("Give the number of columns : "); scanf("%d",&columns);
      printf("Give the lower bound for x : "); scanf("%lf",&a);
      printf("Give the upper bound for x : "); scanf("%lf",&b);
      printf("Give the lower bound for y : "); scanf("%lf",&c);
      printf("Give the upper bound for y : "); scanf("%lf",&d);
      printf("Give the name of the output file : "); scanf("%s",filename);
      outputfile = fopen(filename,"w");
      printf("See the file %s for results...\n",filename);
   }
   else
   {
      rows = 1000; columns = 2000;
      a = -2.0; b = 2.0; c = -2.0; d = 2.0;
      outputfile = fopen("mandelbrot.ps","w");
   }
   write_postscript_plot(outputfile,rows,columns,a,b,c,d);
   fclose(outputfile);

   return 0;
}

void write_postscript_plot
( FILE *f, int rows, int columns, double a, double b, double c, double d )
{
   int i,j,n;
   double x,y;
   double dx = (b-a)/(columns-1);
   double dy = (d-c)/(rows-1);

   fprintf(f,"%%!PS\n");
   fprintf(f,"%%%BoundingBox: 0 300 500 720\n");
   fprintf(f,"/picstr %d string def\n",rows);
   fprintf(f,"%d %d translate\n",XOFFSET,YOFFSET);
   fprintf(f,"%d %d scale\n",SCALE,SCALE);
   fprintf(f,"%d %d %d\n",columns,rows,8);
   fprintf(f,"[%d 0 0 -%d 0 %d]\n",columns,rows,rows);
   fprintf(f,"{currentfile picstr readhexstring pop}\n");
   fprintf(f,"image\n");
   for(i=0,y=d; i<rows; i++,y-=dy,fprintf(f,"\n"))
      for(j=0,x=a; j<columns; j++,x+=dx)
         fprintf(f,"%.2x",255-iterate(x,y)); /* inverted grayscale */
   fprintf(f,"showpage\n");
}

int iterate ( double x, double y )
{
   double wx,wy,v,xx;
   int k = 0;

   wx = 0.0; wy = 0.0; v = 0.0;
   while ((v < 4) && (k++ < 254))
   {
      xx = wx*wx - wy*wy;
      wy = 2.0*wx*wy;
      wx = xx + x;
      wy = wy + y;
      v = wx*wx + wy*wy; 
   }

   return k;
}