/* Mandelbrot Generator */
/* Paul David Buchan, December 20, 2004 */

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

int main (void)
{
int n, ncolor, index, i, j, dum, x, y, xres, yres, maxit, temp[4], pad;
int *tr, *tg, *tb, *c, buf, r, g, b, filesize, iwidth, iheight, bitmapsize;
float rlo, rhi, ilo, ihi, u, v, stepu, stepv;
float r1, r2, i1, i2;
FILE *fp, *fo;

/* Determine how many colors in color palette */
fp=fopen("palette", "r");
if (fp==NULL) {
  printf ("Can't open palette file.\n");
  exit (EXIT_FAILURE);
}
i=0;
while ((n=fgetc(fp)) !=EOF) {
  if (n=='\n') {
    i=i+1;
  }
}
fclose(fp);
ncolor=i;

/* Allocate memory for color palette */
tr=(int *) calloc(ncolor,sizeof(int));
tg=(int *) calloc(ncolor,sizeof(int));
tb=(int *) calloc(ncolor,sizeof(int));

/* Read in color palette */
fp=fopen("palette","r");
for (i=0; i<=ncolor-1; i=i+1) {
  fscanf (fp, "%i %i %i %i",&dum, &tr[i], &tg[i], &tb[i]);
}
fclose(fp);

/* Ask for image dimensions (px) */
printf ("\nWhat is the horizontal dimension (px) ? ");
scanf ("%u", &xres);
printf ("\nWhat is the vertical dimension (px) ? ");
scanf ("%u", &yres);

/* Allocate memory for the array containing iterations */
c=(int *)calloc(xres*yres,sizeof(int));

/* Ask for range on real scale */
printf ("\nWhat is the lowest real value ? ");
scanf ("%f", &rlo);
printf ("\nWhat is the highest real value ? ");
scanf ("%f", &rhi);

/* Ask for range on imaginary scale */
printf ("\nWhat is the lowest imaginary value ? ");
scanf ("%f", &ilo);
printf ("\nWhat is the highest imaginary value ? ");
scanf ("%f", &ihi);

/* Ask for maximum allowable number of iterations */
printf ("\nWhat is the maximum allowable number of iterations ? ");
scanf ("%i", &maxit);

/* Open output file */
fo = fopen ("output.bmp", "wb");
if (fo==NULL) {
  printf ("Can't open new bitmap file.\n");
  exit (EXIT_FAILURE);
}

/* Determine step-sizes in (u,v) plane */
stepu=(rhi-rlo)/xres;
stepv=(ihi-ilo)/yres;

/* File is for Windows O/S */
buf=66;
fputc (buf, fo);
buf=77;
fputc (buf, fo);

/* File size */
filesize=(xres*yres*3)+54;
pad=0;
if ((filesize%4)>0) {
  filesize=filesize+1;
  pad=pad+1;
}
if ((filesize%4)>0) {
  filesize=filesize+1;
  pad=pad+1;
}
if ((filesize%4)>0) {
  filesize=filesize+1;
  pad=pad+1;
}
temp[3]=(filesize/(256*256*256));
filesize=filesize-(temp[3]*256*256*256);
temp[2]=(filesize/(256*256));
filesize=filesize-(temp[2]*256*256);
temp[1]=(filesize/256);
temp[0]=(filesize-(temp[1]*256));
fputc (temp[0], fo);
fputc (temp[1], fo);
fputc (temp[2], fo);
fputc (temp[3], fo);

/* Reserved for future use */
buf=0;
fputc (buf, fo);
fputc (buf, fo);
fputc (buf, fo);
fputc (buf, fo);

/* Offset to BMP data */
buf=54;
fputc (buf, fo);
buf=0;
fputc (buf, fo);
fputc (buf, fo);
fputc (buf, fo);

/* Header is for Windows O/S */
buf=40;
fputc (buf, fo);
buf=0;
fputc (buf, fo);
fputc (buf, fo);
fputc (buf, fo);

/* Image width (px) */
iwidth=xres;
temp[3]=(iwidth/(256*256*256));
iwidth=iwidth-(temp[3]*256*256*256);
temp[2]=(iwidth/(256*256));
iwidth=iwidth-(temp[2]*256*256);
temp[1]=(iwidth/256);
temp[0]=(iwidth-(temp[1]*256));
fputc (temp[0], fo);
fputc (temp[1], fo);
fputc (temp[2], fo);
fputc (temp[3], fo);

/* Image height (px) */
iheight=yres;
temp[3]=(iheight/(256*256*256));
iheight=iheight-(temp[3]*256*256*256);
temp[2]=(iheight/(256*256));
iheight=iheight-(temp[2]*256*256);
temp[1]=(iheight/256);
temp[0]=(iheight-(temp[1]*256));
fputc (temp[0], fo);
fputc (temp[1], fo);
fputc (temp[2], fo);
fputc (temp[3], fo);

/* Number of planes */
buf=1;
fputc (buf, fo);
buf=0;
fputc (buf, fo);

/* Bit depth of image */
buf=24;
fputc (buf, fo);
buf=0;
fputc (buf, fo);

/* No compression */
buf=0;
fputc (buf, fo);
fputc (buf, fo);
fputc (buf, fo);
fputc (buf, fo);

/* BMP data size */
bitmapsize=filesize-54;
temp[3]=(bitmapsize/(256*256*256));
bitmapsize=bitmapsize-(temp[3]*256*256*256);
temp[2]=(bitmapsize/(256*256));
bitmapsize=bitmapsize-(temp[2]*256*256);
temp[1]=(bitmapsize/256);
temp[0]=(bitmapsize-(temp[1]*256));
fputc (temp[0], fo);
fputc (temp[1], fo);
fputc (temp[2], fo);
fputc (temp[3], fo);

/* Horizontal resolution 72 dpi */
buf=18;
fputc (buf, fo);
buf=11;
fputc (buf, fo);
buf=0;
fputc (buf, fo);
fputc (buf, fo);

/* Vertical resolution 72 dpi */
buf=18;
fputc (buf, fo);
buf=11;
fputc (buf, fo);
buf=0;
fputc (buf, fo);
fputc (buf, fo);

/* Use maximum number of colors */
buf=0;
fputc (buf, fo);
fputc (buf, fo);
fputc (buf, fo);
fputc (buf, fo);

/* All colors are important */
buf=0;
fputc (buf, fo);
fputc (buf, fo);
fputc (buf, fo);
fputc (buf, fo);

/* Plot selected area, iterating on each point */
/* Pixels are in (x,y) plane */
i=0;
v=ilo;
for (y=1; y<=yres; y=y+1) {
  u=rlo;
  for (x=1; x<=xres; x=x+1) {
    r1=u;
    i1=v;

/* Iterate until either: maxit is reached, or abs value > 2.0) */
/* Variable c counts iterations */
    c[i]=0;
    r2=0.0;
    i2=0.0;
    while ((((r2*r2)+(i2*i2))<4.0)&&(c[i]<=maxit)) {
      r2=(r1*r1)-(i1*i1)+u;
      i2=(2.0*i1*r1)+v;
      c[i]=c[i]+1;
      r1=r2;
      i1=i2;
    }
    i=i+1;
    u=u+stepu;
  }
  v=v+stepv;
}

/* Find minimum number of iterations taken */
/* in order to scale color range */
i=0;
j=maxit;
for (y=1; y<=yres; y=y+1) {
  for (x=1; x<=xres; x=x+1) {
    if (c[i]<j) {
      j=c[i];
    }
    i=i+1;
  }
}

/* Assign a color to each pixel */
i=0;
for (y=1; y<=yres; y=y+1) {
  for (x=1; x<=xres; x=x+1) {
    if (c[i]>maxit) {  
      r=0;
      g=0;
      b=0; }
    else {
      index=(int)c[i]%ncolor;
      r=tr[index];
      g=tg[index];
      b=tb[index];
    }
    i=i+1;

/* Write pixel to file (RGB) */
    buf=b;
    fputc (buf, fo);
    buf=g;
    fputc (buf, fo);
    buf=r;
    fputc (buf, fo);
  }
}

/* Add file padding to reach 4-byte boundary */
buf=0;
for (i=0; i<pad; i=i+1) {
  fputc (buf, fo);
}
fclose (fo);
free (tr);
free (tg);
free (tb);
free (c);
return (EXIT_SUCCESS);
}