/* +-------------------------------------------------------------------+ */
/* | Copyright 1993, David Koblas (koblas@netcom.com)                  | */
/* |                                                                   | */
/* | Permission to use, copy, modify, and to distribute this software  | */
/* | and its documentation for any purpose is hereby granted without   | */
/* | fee, provided that the above copyright notice appear in all       | */
/* | copies and that both that copyright notice and this permission    | */
/* | notice appear in supporting documentation.  There is no           | */
/* | representations about the suitability of this software for        | */
/* | any purpose.  this software is provided "as is" without express   | */
/* | or implied warranty.                                              | */
/* |                                                                   | */
/* +-------------------------------------------------------------------+ */

#include <X11/Intrinsic.h>
#include <X11/StringDefs.h>
#include <X11/Xatom.h>
#include <stdio.h>
#include "image.h"

/*
**  Some real image processing
*/

#define CLAMP(low, value, high) \
		if (value < low) value = low; else if (value > high) value = high

typedef float	ConvMatrix[3][3];

static Image *convolve(Image *input, ConvMatrix mat, 
			unsigned char *basePixel, Boolean absFlag)
{
	int		x,y,xx,yy, xv,yv;
	float		sum;
	unsigned char	*p;
	unsigned char	*op;
	float		r, g, b;
	int		ir, ig, ib;
	Image		*output;

	sum = 0;
	for (y = 0; y < 3; y++)
		for (x = 0; x < 3; x++)
			sum += mat[x][y];
	if (sum <= 0)
		sum = 0.5;

	output = ImageNew(input->width, input->height);
	op     = output->data;

	for (y = 0; y < input->height; y++) {
		for (x = 0; x < input->width; x++) {
			r = g = b = 0;
			sum = 0;
			for (yy = 0; yy < 3; yy++) {
				for (xx = 0; xx < 3; xx++) {
					xv = x + xx - 1;
					yv = y + yy - 1;
					if (xv < 0 || yv < 0)
						continue;
					if (xv >= input->width || yv >= input->height)
						continue;
					p = ImagePixel(input, xv, yv);
					r += (float)p[0] * mat[xx][yy];
					g += (float)p[1] * mat[xx][yy];
					b += (float)p[2] * mat[xx][yy];
					sum += mat[xx][yy];
				}
			}
			if (sum <= 0)
				sum = 1;
			if (absFlag) {
				if (r < 0) r = -r;
				if (g < 0) g = -g;
				if (b < 0) b = -b;
			}
			ir = r / sum; 
			if (basePixel) ir += basePixel[0];
			CLAMP(0, ir, 255);

			ig = g / sum; 
			if (basePixel) ig += basePixel[1];
			CLAMP(0, ig, 255);

			ib = b / sum; 
			if (basePixel) ib += basePixel[2];
			CLAMP(0, ib, 255);

			*op++ = ir;
			*op++ = ig;
			*op++ = ib;
		}

		if (y % 16 == 0)
			StateTimeStep();
	}

	return output;
}

/*
**  rescale values from 0..255
*/
static void normalize(Image *image)
{
	int		i, count;
	unsigned char	*sp;
	unsigned char	*ip;
	int		v, maxval = 0;

	if (image->cmapSize != 0) {
		sp = image->cmapData;
		count = image->cmapSize * 3;
	} else {
		sp = image->data;
		count = image->width * image->height * 3;
	}

	for (ip = sp, i = 0; i < count; i++, ip++)
		if (*ip > maxval) 
			maxval = *ip;
	if (maxval == 0)
		return;
	for (ip = sp, i = 0; i < count; i++, ip++) {
		v = ((int)*ip * 255) / maxval;
		CLAMP(0,v,255);
		*ip = v;
	}
}

/*
**  Convert image into a monochrome image
*/
static void monochrome(Image *image)
{
	int		i, count;
	unsigned char	*sp;
	unsigned char	*ip;
	int		v;

	if (image->cmapSize != 0) {
		sp = image->cmapData;
		count = image->cmapSize;
	} else {
		sp = image->data;
		count = image->width * image->height;
	}

	for (ip = sp, i = 0; i < count; i++, ip += 3) {
		v = (ip[0]*11 + ip[1]*16 + ip[2]*5) >> 5;  /* pp=.33R+.5G+.17B */
		ip[0] = v;
		ip[1] = v;
		ip[2] = v;
	}
}


Image *ImageSmooth(Image *input)
{
	static ConvMatrix	mat = {
			{ 1, 1, 1 },
			{ 1, 1, 1 },
			{ 1, 1, 1 }
		};
	return convolve(input, mat, NULL, False);
}

Image *ImageSharpen(Image *input)
{
	static ConvMatrix	mat = {
			{ -1, -2, -1 },
			{ -2, 20, -2 },
			{ -1, -2, -1 }
		};
	return convolve(input, mat, NULL, False);
}

Image *ImageEdge(Image *input)
{
	static ConvMatrix	mat = {
			{ -1, -2,  0 },
			{ -2,  0,  2 },
			{  0,  2,  1 }
		};
	Image	*image = convolve(input, mat, NULL, True);

	normalize(image);

	return image;
}

Image *ImageEmbose(Image *input)
{
	static ConvMatrix	mat = {
			{ -1, -2,  0 },
			{ -2,  0,  2 },
			{  0,  2,  1 }
		};
	static unsigned char	base[3] = { 128, 128, 128 };
	Image	*image = convolve(input, mat, base, False);

	monochrome(image);
	normalize(image);

	return image;
}

Image *ImageInvert(Image *input)
{
	Image		*output = ImageNewCmap(input->width, input->height, input->cmapSize);
	int		i;
	unsigned char	*ip, *op;
	int		count;
	
	/*
	**   If the input has a colormap, just invert that.
	*/
	if (input->cmapSize != 0) {
		ip = input->cmapData;
		op = output->cmapData;
		count = input->cmapSize;

		memcpy(output->data, input->data, 
			sizeof(char) * input->scale * input->width * input->height);
	} else {
		ip = input->data;
		op = output->data;
		count = input->width * input->height;
	}

	for (i = 0; i < count; i++) {
		*op++ = 255 - *ip++;
		*op++ = 255 - *ip++;
		*op++ = 255 - *ip++;
	}

	return output;
}

#define AREA 	7
#define AREA_2	(AREA/2)

Image *ImageOilPaint(Image *input)
{
	Image		*output = ImageNew(input->width, input->height);
	unsigned char	*op = output->data;
	int		x, y, xx, yy, i;
	int		rVal, gVal, bVal;
	int		rCnt, gCnt, bCnt;
	int		rHist[256], gHist[256], bHist[256];

	for (y = 0; y < input->height; y++) {
		for (x = 0; x < input->width; x++) {
			/*
			**  compute histogram of (on-screen hunk of) n*n 
			**    region centered plane
			*/

			rCnt = gCnt = bCnt = 0;
			rVal = gVal = bVal = 0;
			for (i = 0; i < XtNumber(rHist); i++)
				rHist[i] = gHist[i] = bHist[i] = 0;
			
			for (yy = y - AREA_2; yy < y + AREA_2; yy++) {
				if (yy < 0 || yy >= input->height)
					continue;
				for (xx = x - AREA_2; xx < x + AREA_2; xx++) {
					int		c, p;
					unsigned char	*rgb;

					if (xx < 0 || xx >= input->width)
						continue;

					rgb = ImagePixel(input, xx, yy);

					if ((c= ++rHist[(p=rgb[0])/4]) > rCnt) {
						rVal = p;
						rCnt = c;
					}
					if ((c= ++gHist[(p=rgb[1])/4]) > gCnt) {
						gVal = p;
						gCnt = c;
					}
					if ((c= ++bHist[(p=rgb[2])/4]) > bCnt) {
						bVal = p;
						bCnt = c;
					}
				}
			}

			*op++ = rVal;
			*op++ = gVal;
			*op++ = bVal;
		}

		if (y % 16 == 0)
			StateTimeStep();
	}

	return output;
}
