/*
 * polar 2.00 -- image filter plug-in for The GIMP
 * Copyright (C) 1996 Marc Bless
 *
 * E-mail: bless@ai-lab.fh-furtwangen.de
 * WWW:    www.ai-lab.fh-furtwangen.de/~bless
 *
 * This program is free software; you can redistribute it and/or modify
 * it under the terms of the GNU General Public License as published by
 * the Free Software Foundation; either version 2 of the License, or
 * (at your option) any later version.
 *
 * This program is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 * GNU General Public License for more details.
 *
 * You should have received a copy of the GNU General Public License
 * along with this program; if not, write to the Free Software
 * Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
 */

/* 
 * This filter converts an image from rectangular to polar.
 */

/*
 * History:
 * 
 * v 1.00 - 1996 June 01 : standard function works.
 *               June 02 : top in middle, backwards, 
 *                         circle depth, angle
 * v 1.01 - 1996 June 04 : changed PI to M_PI, eliminated rint()
 * v 1.02 - 1996 June 14 : use last calculation
 *
 * v 2.00 - 1996 July 14 : added 'Polar to Rectangular',
 *												 major and minor bugs fixed
 *
 */

#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include "gimp.h"

#define DIALOG_TITLE "Polar/Rectangular (v 2.00)"

#define sqr(x)	((x) * (x))
#define WITHIN(a, b, c) ((((a) <= (b)) && ((b) <= (c))) ? 1 : 0)

typedef unsigned char uchar;

static int 	 dlg_init (void);
static void  polar (Image, Image);
static void  rectangular (Image, Image);
static uchar bilinear (double, double, uchar *);
static void  cb_scale (int, void *, void *);
static void  cb_ok (int, void *, void *);
static void  cb_cancel (int, void *, void *);
static void  cb_toggle (int, void *, void *);
static void  cb_radio (int, void *, void *);
static void  saveimage (void);
static void  freshen (void);
static void  calc_all (void);

static char *prog_name;

static int   dialogID;

static long  pol       = 1;			/* start with Rec2Pol */
static long  rec			 = 0;
static long  inverse	 = 0;
static long  backward  = 0;
static long  aapply    = 0;
static long  uselast   = 0;
static long  circle    = 100;
static long  angle     = 0;

static Image input, output;

static unsigned char *saved;

/*
 * The main function.
 */

int main (int argc, char **argv)
{
  prog_name = argv[0];

  if (!gimp_init (argc, argv))
		return 1;

  input = 0;
  output = 0;
      
  if (!(input = gimp_get_input_image (0)))
		return 1;

	if (!(output = gimp_get_output_image (0)))
	{
		gimp_free_image (input);
		gimp_quit ();
		return 1;
	}

	if ((gimp_image_type (input) != RGB_IMAGE) &&
			(gimp_image_type (input) != GRAY_IMAGE))
	{
  	gimp_message ("polar: can operate on rgb and gray image types only");
		gimp_free_image (input);
		gimp_free_image (output);
		gimp_quit ();
		return 1;
	}

	saveimage ();		/* save our input image before anything happens */

	dialogID = dlg_init ();
						
	if (!gimp_show_dialog (dialogID))
	{
		if (aapply)
		{
			freshen ();
			gimp_update_image (output);
		}
	}
	else
	{
		if (!aapply)
			calc_all ();
	}

	free (saved);

	gimp_free_image (input);
	gimp_free_image (output);

  gimp_quit ();

  return 0;
}

/* 
 * The dialog init function.
 */

static int dlg_init (void)
{
  int mainID,   labelID, typeID,    typegrpID,    paramID,   paramgrpID, 
      polID,    recID,   inverseID, backwardID,   aapplyID,  uselastID, 
      circleID, angleID, dlgID,     previewgrpID, previewID;

	dlgID = gimp_new_dialog (DIALOG_TITLE);
  mainID   = gimp_new_row_group (dlgID, DEFAULT, NORMAL, "");
	typeID   = gimp_new_frame (dlgID, mainID, "Type");
		typegrpID = gimp_new_row_group (dlgID, typeID, RADIO, "");
			polID    = gimp_new_radio_button (dlgID, typegrpID, "Rectangular to Polar");
			gimp_change_item (dlgID, polID, sizeof (pol), &pol);
			recID    = gimp_new_radio_button (dlgID, typegrpID, "Polar to Rectangular");
			gimp_change_item (dlgID, recID, sizeof (rec), &rec);
	paramID  = gimp_new_frame (dlgID, mainID, "Parameters");
		paramgrpID = gimp_new_row_group (dlgID, paramID, NORMAL, "");
			labelID  = gimp_new_label (dlgID, paramgrpID, "Circle depth in %");
			circleID = gimp_new_scale (dlgID, paramgrpID, 0, 100, 100, 0);
			labelID  = gimp_new_label (dlgID, paramgrpID, "Angle");
			angleID  = gimp_new_scale (dlgID, paramgrpID, 0, 359, 0, 0);
			inverseID = gimp_new_check_button (dlgID, paramgrpID, "Top in middle");
			gimp_change_item (dlgID, inverseID, sizeof (inverse), &inverse);
			backwardID = gimp_new_check_button (dlgID, paramgrpID, "Backwards");
			gimp_change_item (dlgID, backwardID, sizeof (backward), &backward);
	previewID = gimp_new_frame (dlgID, mainID, "Preview");
		previewgrpID = gimp_new_row_group (dlgID, previewID, NORMAL, "");
			aapplyID = gimp_new_check_button (dlgID, previewgrpID, "Auto Apply");
			gimp_change_item (dlgID, aapplyID, sizeof (aapply), &aapply);
			uselastID = gimp_new_check_button (dlgID, previewgrpID, "Use last calculation");
			gimp_change_item (dlgID, uselastID, sizeof (uselast), &uselast);
			
	gimp_add_callback (dlgID, polID, cb_radio, &pol);
	gimp_add_callback (dlgID, recID, cb_radio, &rec);
	gimp_add_callback (dlgID, circleID, cb_scale, &circle);
	gimp_add_callback (dlgID, angleID, cb_scale, &angle);
	gimp_add_callback (dlgID, inverseID, cb_toggle, &inverse);
	gimp_add_callback (dlgID, backwardID, cb_toggle, &backward);
	gimp_add_callback (dlgID, aapplyID, cb_toggle, &aapply);
	gimp_add_callback (dlgID, uselastID, cb_toggle, &uselast);
	gimp_add_callback (dlgID, gimp_ok_item_id (dlgID), cb_ok, 0);
	gimp_add_callback (dlgID, gimp_cancel_item_id (dlgID), cb_cancel, 0);

	return dlgID;
}

/* 
 * The radio button callback function.
 */

static void cb_radio (int itemID, void *client_data, void *call_data)
{
	if (aapply && (*((long*) client_data) != *((long*) call_data)))
	{
		*((long*) client_data) = *((long*) call_data);
		calc_all ();
	}
	else
		*((long*) client_data) = *((long*) call_data);
}

/* 
 * The check button callback function.
 */

static void cb_toggle (int itemID, void *client_data, void *call_data)
{
	*((long*) client_data) = *((long*) call_data);
	if (aapply)
		calc_all ();
	else
	{
		freshen ();
		gimp_update_image (output);
	}
}

/*
 * The scale callback function.
 */

static void cb_scale (int itemID, void *client_data, void *call_data)
{
	if (aapply && (*((long*) client_data) != *((long*) call_data)))
	{
		*((long*) client_data) = *((long*) call_data);
		calc_all ();
	}
	else
		*((long*) client_data) = *((long*) call_data);
}

/*
 * The saveimage function.
 */

static void saveimage (void)
{
	saved = (unsigned char *) malloc (gimp_image_width (input) *
	                          gimp_image_height (input) *
														gimp_image_channels (input));
	memcpy (saved, gimp_image_data (input),
	               gimp_image_width (input) *
								 gimp_image_height (input) *
								 gimp_image_channels (input));
}

/*
 * The freshen function.
 */

static void freshen (void)
{
	memcpy (gimp_image_data (output), saved,
		      gimp_image_width (input) *
					gimp_image_height (input) *
					gimp_image_channels (input));
}

/*
 * The ok button callback function.
 */

static void cb_ok (int itemID, void *client_data, void *call_data)
{
	gimp_close_dialog (dialogID, 1);
}

/*
 * The cancel button callback function.
 */

static void cb_cancel (int itemID, void *client_data, void *call_data)
{
	gimp_close_dialog (dialogID, 0);
}

/*
 * The polar filter function.
 */

static void polar (Image linput, Image loutput)
{
  long width, height;
  long channels, rowstride;
  unsigned char *src_row, *dest_row;
  unsigned char *src, *dest;
  unsigned char *src_org;
  uchar *p;
  uchar values[4];
  int k;
  uchar val;
  short row, col;
  int x1, y1, x2, y2;
  int xi, yi;
	int xdiff, ydiff;
  double x, y, xm, ym, r, rmax, m;
  double xmax, ymax;
	double t;
  double phi;
	double angl;
  int progress, max_progress;

  gimp_image_area (linput, &x1, &y1, &x2, &y2);

  width     = gimp_image_width (linput);
  height    = gimp_image_height (linput);
  channels  = gimp_image_channels (linput);
  rowstride = width * channels;

	xdiff = x2 - x1;
	ydiff = y2 - y1;

  progress     = 0;
  max_progress = ydiff;
	gimp_init_progress ("Rectangular to Polar");

	if (uselast)
		src_org = src_row = gimp_image_data (linput);
	else
		src_org = src_row = saved;
  dest_row = gimp_image_data (loutput);

  src_row  += rowstride * y1 + x1 * channels;
  dest_row += rowstride * y1 + x1 * channels;

  xm = xdiff / 2.0;
  ym = ydiff / 2.0;

	angl = (double)angle / 180.0 * M_PI;

  for (row = y1; row < y2; row++)
  {
    src     = src_row;
    dest    = dest_row;

    for (col = x1; col < x2; col++)
    {
			if (col >= xm)
			{
				if (row > ym)
				{
					phi = M_PI - atan (((double)(col - xm))/((double)(row - ym)));
					r   = sqrt (sqr (col - xm) + sqr (row - ym));
				}
				else
				if (row < ym)
				{
					phi = atan (((double)(col - xm))/((double)(ym - row)));
					r   = sqrt (sqr (col - xm) + sqr (ym - row));
				}
				else
				{
					phi = M_PI / 2;
					r   = col - xm; /* xm - x1; */
				}
			}
			else
			if (col < xm)
			{
				if (row < ym)
				{
					phi = 2 * M_PI - atan (((double)(xm - col))/((double)(ym - row)));
					r   = sqrt (sqr (xm - col) + sqr (ym - row));
				}
				else
				if (row > ym)
				{
					phi = M_PI + atan (((double)(xm - col))/((double)(row - ym)));
					r   = sqrt (sqr (xm - col) + sqr (row - ym));
				}
				else
				{
					phi = 1.5 * M_PI;
					r   = xm - col; /* xm - x1; */
				}
			}

			if (col != xm)
				m = fabs (((double)(row - ym)) / ((double)(col - xm)));
			else
				m = 0;

			if (m <= ((double)(ydiff) / (double)(xdiff)))
			{
				if (col == xm)
				{
					xmax = 0;
					ymax = ym - y1;
				}
				else
				{
	  			xmax = xm - x1;
	  			ymax = m * xmax;
				}
			}
			else
			{
	  		ymax = ym - y1;
	  		xmax = ymax / m;
			}
	
			rmax = sqrt ( (double)(sqr (xmax) + sqr (ymax)) );

			t = ((ym - y1) < (xm - x1)) ? (ym - y1) : (xm - x1);
			rmax = (rmax - t) / 100 * (100 - circle) + t;

			phi = fmod (phi + angl, 2*M_PI);

			if (backward)
				x = x2 - 1 - (xdiff - 1)/(2*M_PI) * phi;
			else
				x = (xdiff - 1)/(2*M_PI) * phi + x1;
			if (inverse)
				y = (ydiff)/rmax   * r   + y1;
			else
				y = y2 - (ydiff)/rmax * r;

			xi = (int) (x+0.5);
			yi = (int) (y+0.5);

			p = src_org + rowstride * yi + xi * channels;

			for (k = 0; k < channels; k++)
			{
/*				if (WITHIN(0, xi, width - 1) && WITHIN(0, yi, height - 1))
*/				if (WITHIN(x1, xi, x2 - 1) && WITHIN(y1, yi, y2 - 1))
					values[0] = *(p + k);
				else
					values[0] = 0;

/*				if (WITHIN(0, xi+1, width - 1) && WITHIN(0, yi, height - 1))
*/				if (WITHIN(x1, xi+1, x2 - 1) && WITHIN(y1, yi, y2 - 1))
					values[1] = *(p + channels + k);
				else
					values[1] = 0;

/*				if (WITHIN(0, xi, width - 1) && WITHIN(0, yi + 1, height - 1))
*/				if (WITHIN(x1, xi, x2 - 1) && WITHIN(y1, yi + 1, y2 - 1))
					values[2] = *(p + rowstride + k);
				else
					values[2] = 0;

/*				if (WITHIN(0, xi + 1, width - 1) && WITHIN(0, yi + 1, height - 1))
*/				if (WITHIN(x1, xi + 1, x2 - 1) && WITHIN(y1, yi + 1, y2 - 1))
					values[3] = *(p + rowstride + channels + k);
				else
					values[3] = 0;

				val = bilinear (x, y, values);

				*dest++ = val;
			}

		}

    src_row += rowstride;
    dest_row += rowstride;

    if (++progress % 2 == 0)
    	gimp_do_progress (progress, max_progress);
  }

	gimp_do_progress (1, 1);
}

/*
 * The rectangular filter function.
 */

static void
rectangular (linput, loutput)
     Image linput, loutput;
{
  long width, height;
  long channels, rowstride;
  unsigned char *src_row, *dest_row;
  unsigned char *src, *dest;
  unsigned char *src_org;
  uchar *p;
  uchar values[4];
  int k;
  uchar val;
  short row, col;
  int x1, y1, x2, y2;
  int xi, yi;
	int xdiff, ydiff;
  double x, y, xm, ym, r, rmax, m, xx, yy;
  double xmax, ymax;
	double t;
  double phi;
	double phi2;
	double angl;
  int progress, max_progress;
	
  gimp_image_area (linput, &x1, &y1, &x2, &y2);

  width     = gimp_image_width (linput);
  height    = gimp_image_height (linput);
  channels  = gimp_image_channels (linput);
  rowstride = width * channels;

	xdiff = x2 - x1;
	ydiff = y2 - y1;

  progress     = 0;
  max_progress = ydiff;
	gimp_init_progress ("Polar to Rectangular");

	if (uselast)
		src_org = src_row = gimp_image_data (linput);
	else
		src_org = src_row = saved;
  dest_row = gimp_image_data (loutput);

  src_row  += rowstride * y1 + x1 * channels;
  dest_row += rowstride * y1 + x1 * channels;

  xm = xdiff / 2.0;
  ym = ydiff / 2.0;

	angl = (double)angle / 180.0 * M_PI;

  for (row = y1; row < y2; row++)
  {
    src     = src_row;
    dest    = dest_row;

    for (col = x1; col < x2; col++)
    {
			if (backward)
				phi = (2 * M_PI) * (x2 - col) / xdiff;
			else
				phi = (2 * M_PI) * (col - x1) / xdiff;

			phi = fmod (phi + angl, 2 * M_PI);

			if (phi >= 1.5 * M_PI)
				phi2 = 2 * M_PI - phi;
			else
			if (phi >= M_PI)
				phi2 = phi - M_PI;
			else
			if (phi >= 0.5 * M_PI)
				phi2 = M_PI - phi;
			else
				phi2 = phi;

			xx = tan (phi2);
			if (xx != 0)
				m = (double) 1.0 / xx;
			else
				m = 0;

			if (m <= ((double)(ydiff) / (double)(xdiff)))
			{
				if (phi2 == 0)
				{
					xmax = 0;
					ymax = ym - y1;
				}
				else
				{
					xmax = xm - x1;
					ymax = m * xmax;
				}
			}
			else
			{
				ymax = ym - y1;
				xmax = ymax / m;
			}

			rmax = sqrt ((double)(sqr (xmax) + sqr (ymax)));
			
			t = ((ym - y1) < (xm - x1)) ? (ym - y1) : (xm - x1);

			rmax = (rmax - t) / 100.0 * (100 - circle) + t;

			if (inverse)
				r = rmax * (double)((row - y1) / (double)(ydiff));
			else
				r = rmax * (double)((y2 - row) / (double)(ydiff));

			xx = r * sin (phi2);
			yy = r * cos (phi2);
			
			if (phi >= 1.5 * M_PI)
			{
				x = (double)xm - xx;
				y = (double)ym - yy;
			}
			else
			if (phi >= M_PI)
			{
				x = (double)xm - xx;
				y = (double)ym + yy;
			}
			else
			if (phi >= 0.5 * M_PI)
			{
				x = (double)xm + xx;
				y = (double)ym + yy;
			}
			else
			{
				x = (double)xm + xx;
				y = (double)ym - yy;
			}

			xi = (int)(x + 0.5);
			yi = (int)(y + 0.5);

			p = src_org + rowstride * yi + xi * channels;

			for (k = 0; k < channels; k++)
			{
				if (WITHIN(0, xi, width - 1) && WITHIN(0, yi, height - 1))
					values[0] = *(p + k);
				else
					values[0] = 0;

				if (WITHIN(0, xi+1, width - 1) && WITHIN(0, yi, height - 1))
					values[1] = *(p + channels + k);
				else
					values[1] = 0;

				if (WITHIN(0, xi, width - 1) && WITHIN(0, yi + 1, height - 1))
					values[2] = *(p + rowstride + k);
				else
					values[2] = 0;

				if (WITHIN(0, xi + 1, width - 1) && WITHIN(0, yi + 1, height - 1))
					values[3] = *(p + rowstride + channels + k);
				else
					values[3] = 0;

				val = bilinear (x, y, values);

				*dest++ = val;
			}

		}

    src_row  += rowstride;
    dest_row += rowstride;

    if (++progress % 2 == 0)
    	gimp_do_progress (progress, max_progress);
  }

	gimp_do_progress (1, 1);
}

static uchar
bilinear(double x, double y, uchar *v)
{
    double xx, yy, m0, m1;

    xx = fmod(x, 1.0);
    yy = fmod(y, 1.0);

    m0 = (1.0 - xx) * v[0] + xx * v[1];
    m1 = (1.0 - xx) * v[2] + xx * v[3];

    return (uchar) ((1.0 - yy) * m0 + yy * m1);
} 

static void calc_all (void)
{
	if (pol && !rec)
	{
		polar (input, output);
		gimp_update_image (output);
	}
	else
	if (!pol && rec)
	{
		rectangular (input, output);
		gimp_update_image (output);
	}
}