#define setbit(bitmap, n)	(bitmap[n/8] |= (0x80 >> (n % 8)))
#define clearbit(bitmap, n)	(bitmap[n/8] &= (~(0x80 >> (n % 8))))
/***************************************************************************
 * Func: error_diffusion                                                   *
 *                                                                         *
 * Desc: halftones an image via error diffusion, implements the zig-zag    *
 *       processing pattern to eliminate avalanche patterns                *
 *                                                                         *
 * Params: image - pointer to image in memory                              *
 *         rows - number of rows in the image                              *
 *         cols - number of columns in the image                           *
 *         threshold - threshold                                           *
 *         filename - name of output file to create                        *
 * Returns: nothing                                                        *
 ***************************************************************************/

void error_diffusion(image_ptr image, int rows, int cols,
		     int threshold, char *filename)
    {
    float error;                  /* error to propagate */
    FILE *fp;                     /* file pointer for output file */
    unsigned long index;          /* index within input image */
    int row_size;                 /* size of row in bytes */
    unsigned long i, j;           /* loop variables */
    unsigned char line_buff[1024];/* storage for line to write out */
    unsigned char scanline;       /* indicates even or odd scanline */

    /* open new output file */
    if((fp=fopen(filename, "wb")) == NULL)
	{
	printf("Unable to open %s for output\n",filename);
	exit(1);
	}

    /* print out the portable bitmap header */
    fprintf(fp, "P4\n%d %d\n", cols-2, rows-1);
    row_size = (cols-2)/8 + (1 * (((cols-2) % 8)) != 0) ;

    scanline = 1;
    for(i=0; i<(rows-1); i++)
	{
	if(scanline)
	    index = i * cols;           /* start at beginning of scanline */
	else
	    index = i * cols + cols - 1; /* end of scanline */

	for(j=1; j<(cols-1); j++)
	    {
	   if(image[index] < threshold)
		{
		setbit(line_buff, j);   /* this function can be swapped */
		error = image[index];
		}
	    else
		{
		clearbit(line_buff, j); /* with this one if necessary */
		error = image[index] - 255;
		}
	    if(scanline)
		{
		image[index+1] = CLAMP(image[index+1] + (int) (error * 0.4375), 0, 255);      /* 7/16 */
		image[index+cols-1] = CLAMP(image[index+cols-1] + (int) (error * 0.1875), 0, 255); /* 3/16 */
		image[index+cols] = CLAMP(image[index+cols] + (int) (error * 0.3125), 0, 255);   /* 5/16 */
		image[index+cols+1] = CLAMP(image[index+cols+1] + (int)(error * 0.0625), 0, 255);  /* 1/16 */
		}
	    else
		{
		image[index-1] = CLAMP(image[index-1] + (int) (error * 0.4375), 0, 255);      /* 7/16 */
		image[index+cols+1] = CLAMP(image[index+cols+1] + (int) (error * 0.1875), 0, 255); /* 3/16 */
		image[index+cols] = CLAMP(image[index+cols] + (int) (error * 0.3125), 0, 255);   /* 5/16 */
		image[index+cols-1] = CLAMP(image[index+cols-1] + (int)(error * 0.0625), 0, 255);  /* 1/16 */
		}
	    if(scanline)
		index++;
	    else
		index--;
	    }
	fwrite(line_buff, 1, row_size, fp);    /* write out line */
	scanline = ~scanline;                  /* toggle even/odd indicator */
	}
    fclose(fp);
    }