path: root/generic/Geo.c

/*
 * Geo.c -- Implementation of common geometric routines.
 *
 * Authors		: Patrick Lecoanet.
 * Creation date	:
 *
 * $Id$
 */

/*
 *  Copyright (c) 1993 - 1999 CENA, Patrick Lecoanet --
 *
 * This code is free software; you can redistribute it and/or
 * modify it under the terms of the GNU Library General Public
 * License as published by the Free Software Foundation; either
 * version 2 of the License, or (at your option) any later version.
 *
 * This code 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
 * Library General Public License for more details.
 *
 * You should have received a copy of the GNU Library General Public
 * License along with this code; if not, write to the Free
 * Software Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
 *
 */

/*
 * Much of the code here is inspired by (or copied from) the Tk code.
 */


#include "Geo.h"
#include "WidgetInfo.h"

#include <memory.h>


static const char rcsid[] = "$Id$";
static const char compile_id[]="$Compile: " __FILE__ " " __DATE__ " " __TIME__ " $";


/*
 * Compute the origin of the rectangle given
 * by position, anchor, width and height.
 */
void
Anchor2Origin(RadarPoint	*position,
	      RadarDim		width,
	      RadarDim		height,
	      RadarAnchor	anchor,
	      RadarPoint	*origin)
{
  switch (anchor) {
  case RadarAnchorCenter:
    origin->x = position->x - width/2;
    origin->y = position->y - height/2;
    break;
  case RadarAnchorNW:
    *origin = *position;
    break;
  case RadarAnchorN:
    origin->x = position->x - width/2;
    origin->y = position->y;
    break;
  case RadarAnchorNE:
    origin->x = position->x - width;
    origin->y = position->y;
    break;
  case RadarAnchorE:
    origin->x = position->x - width;
    origin->y = position->y - height/2;
    break;
  case RadarAnchorSE:
    origin->x = position->x - width;
    origin->y = position->y - height;
    break;
  case RadarAnchorS:
    origin->x = position->x - width/2;
    origin->y = position->y - height;
    break;
  case RadarAnchorSW:
    origin->x = position->x;
    origin->y = position->y - height;
    break;
  case RadarAnchorW:
    origin->x = position->x;
    origin->y = position->y - height/2;
    break;
  }
}


/*
 * Compute the anchor position given the bbox origin, width,
 * height and the anchor.
 */
void
Origin2Anchor(RadarPoint	*origin,
	      RadarDim		width,
	      RadarDim		height,
	      RadarAnchor	anchor,
	      RadarPoint	*position)
{
  switch (anchor) {
  case RadarAnchorCenter:
    position->x = origin->x + width/2;
    position->y = origin->y + height/2;
    break;
  case RadarAnchorNW:
    *position = *origin;
    break;
  case RadarAnchorN:
    position->x = origin->x + width/2;
    position->y = origin->y;
    break;
  case RadarAnchorNE:
    position->x = origin->x + width;
    position->y = origin->y;
    break;
  case RadarAnchorE:
    position->x = origin->x + width;
    position->y = origin->y + height/2;
    break;
  case RadarAnchorSE:
    position->x = origin->x + width;
    position->y = origin->y + height;
    break;
  case RadarAnchorS:
    position->x = origin->x + width/2;
    position->y = origin->y + height;
    break;
  case RadarAnchorSW:
    position->x = origin->x;
    position->y = origin->y + height;
    break;
  case RadarAnchorW:
    position->x = origin->x;
    position->y = origin->y + height/2;
    break;
  }
}


void
BBox2XRect(RadarBBox	*bbox,
	   XRectangle	*r)
{
  r->x = REAL_TO_INT(bbox->orig.x);
  r->y = REAL_TO_INT(bbox->orig.y);
  r->width = REAL_TO_INT(bbox->corner.x) - r->x;
  r->height = REAL_TO_INT(bbox->corner.y) - r->y;
}


void
GetStringBBox(char	*str,
	      RadarFont	font,
	      RadarPos	x,
	      RadarPos	y,
	      RadarBBox	*str_bbox)
{
  Tk_FontMetrics fm;

  str_bbox->orig.x = x;
  str_bbox->corner.x  = x + RadarTextWidth(font, str, strlen(str));
  Tk_GetFontMetrics(font, &fm);
  str_bbox->orig.y = y - fm.ascent;
  str_bbox->corner.y = str_bbox->orig.y + fm.ascent + fm.descent;
}


void
ResetBBox(RadarBBox *bbox)
{
  bbox->orig.x = bbox->orig.y = 0;
  bbox->corner = bbox->orig;
}


void
CopyBBox(RadarBBox *bbox_from,
	 RadarBBox *bbox_to)
{
  bbox_to->orig = bbox_from->orig;
  bbox_to->corner = bbox_from->corner;
}


void
IntersectBBox(RadarBBox *bbox1,
	      RadarBBox *bbox2,
	      RadarBBox *bbox_inter)
{
  if ((bbox1->corner.x < bbox2->orig.x) ||
      (bbox1->corner.y < bbox2->orig.y) ||
      (bbox2->corner.x < bbox1->orig.x) ||
      (bbox2->corner.y < bbox1->orig.y)) {
    ResetBBox(bbox_inter);
  }
  else {
    bbox_inter->orig.x = MAX(bbox1->orig.x, bbox2->orig.x);
    bbox_inter->orig.y = MAX(bbox1->orig.y, bbox2->orig.y);
    bbox_inter->corner.x = MIN(bbox1->corner.x, bbox2->corner.x);
    bbox_inter->corner.y = MIN(bbox1->corner.y, bbox2->corner.y);
  }
}


RadarBool
IsEmptyBBox(RadarBBox *bbox)
{
  return (bbox->orig.x >= bbox->corner.x) || (bbox->orig.y >= bbox->corner.y);
}


void
AddBBoxToBBox(RadarBBox *bbox,
	      RadarBBox *bbox2)
{
  if (IsEmptyBBox(bbox2)) {
    return;
  }
  if (IsEmptyBBox(bbox)) {
    CopyBBox(bbox2, bbox);
  }
  else {
    bbox->orig.x = MIN(bbox->orig.x, bbox2->orig.x);
    bbox->orig.y = MIN(bbox->orig.y, bbox2->orig.y);
    bbox->corner.x = MAX(bbox->corner.x, bbox2->corner.x);
    bbox->corner.y = MAX(bbox->corner.y, bbox2->corner.y);
  }
}


void
AddPointToBBox(RadarBBox	*bbox,
	       RadarPos		px,
	       RadarPos		py)
{
  if (IsEmptyBBox(bbox)) {
    bbox->orig.x = px;
    bbox->orig.y = py;
    bbox->corner.x = bbox->orig.x + 1;
    bbox->corner.y = bbox->orig.y + 1;
  }
  else {
    bbox->orig.x = MIN(bbox->orig.x, px);
    bbox->orig.y = MIN(bbox->orig.y, py);
    bbox->corner.x = MAX(bbox->corner.x, px + 1);
    bbox->corner.y = MAX(bbox->corner.y, py + 1);
  }
}


void
AddPointsToBBox(RadarBBox	*bbox,
		RadarPoint	*points,
		int		num_points)
{
  int		x1, y1, x2, y2, cur;

  if (points == NULL) {
    return;
  }
  
  if (num_points == 0) {
    return;
  }
  
  if (IsEmptyBBox(bbox)) {
    x1 = points->x;
    y1 = points->y;
    x2 = x1 + 1;
    y2 = y1 + 1;
    num_points--;
    points++;
  }
  else {
    x1 = bbox->orig.x;
    y1 = bbox->orig.y;
    x2 = bbox->corner.x;
    y2 = bbox->corner.y;
  }

  for ( ; num_points > 0; num_points--, points++) {
    cur = points->x;
    if (cur < x1) {
      x1 = cur;
    }
    if (cur > x2) {
      x2 = cur;
    }
    cur = points->y;
    if (cur < y1) {
      y1 = cur;
    }
    if (cur > y2) {
      y2 = cur;
    }
  }
  bbox->orig.x = x1;
  bbox->orig.y = y1;
  bbox->corner.x = x2 + 1;
  bbox->corner.y = y2 + 1;
}


void
AddStringToBBox(RadarBBox	*bbox,
		char		*str,
		RadarFont	font,
		RadarPos	cx,
		RadarPos	cy)
{
  RadarBBox	str_bbox;
  
  GetStringBBox(str, font, cx, cy, &str_bbox);
  AddBBoxToBBox(bbox, &str_bbox);
}

RadarBool
PointInBBox(RadarBBox	*bbox,
	    RadarPos	x,
	    RadarPos	y)
{
  return ((x >= bbox->orig.x) && (x < bbox->corner.x) &&
	  (y >= bbox->orig.y) && (y < bbox->corner.y));
}


/*
 * Tell where aa area is with respect to another area.
 * Return -1 if the first is entirely outside the second,
 * 1 if it is entirely inside and 0 otherwise.
 */
int
BBoxInBBox(RadarBBox	*bbox1,
	   RadarBBox	*bbox2)
{
  if ((bbox1->corner.x <= bbox2->orig.x) ||
      (bbox1->orig.x >= bbox2->corner.x) ||
      (bbox1->corner.y <= bbox2->orig.y) ||
      (bbox1->orig.y >= bbox2->corner.y)) {
    return -1;
  }
  if ((bbox2->orig.x <= bbox1->orig.x) &&
      (bbox1->corner.x <= bbox2->corner.x) &&
      (bbox2->orig.y <= bbox1->orig.y) &&
      (bbox1->corner.y <= bbox2->corner.y)) {
    return 1;
  }
  return 0;
}

/*
 * Tell where a line is with respect to an area.
 * Return -1 if the line is entirely outside the bbox, 1
 * if it is entirely inside and 0 otherwise.
 */
int
LineInBBox(RadarPoint	*p1,
	   RadarPoint	*p2,
	   RadarBBox	*bbox)
{
  RadarBool	p1_inside = PointInBBox(bbox, p1->x, p1->y);
  RadarBool	p2_inside = PointInBBox(bbox, p2->x, p2->y);

  if (p1_inside != p2_inside) {
    return 0;
  }

  if (p1_inside && p2_inside) {
    return 1;
  }
  
  /*
   * Segment may intersect area, check it more thoroughly.
   */
  /* Vertical line */
  if (p1->x == p2->x) {
    if (((p1->y >= bbox->orig.y) ^ (p2->y >= bbox->orig.y)) &&
	(p1->x >= bbox->orig.x) &&
	(p1->x <= bbox->corner.x)) {
      return 0;
    }
  }
  /* Horizontal line */
  else if (p1->y == p2->y) {
    if (((p1->x >= bbox->orig.x) ^ (p2->x >= bbox->orig.x)) &&
	(p1->y >= bbox->orig.y) &&
	(p1->y <= bbox->corner.y)) {
      return 0;
    }
  }
  /* Diagonal, do it the hard way. */
  else {
    double	slope = ((double) p2->y - p1->y) / ((double) p2->x - p1->x);
    int		low, high, x, y;
    int		bbox_width = bbox->corner.x - bbox->orig.x;
    int		bbox_height = bbox->corner.y - bbox->orig.y;

    /* Check against left edge */
    if (p1->x < p2->x) {
      low = p1->x;
      high = p2->x;
    }
    else {
      low = p2->x;
      high = p1->x;
    }

    y = p1->y + (bbox->orig.x - p1->x) * slope;
    if ((bbox->orig.x >= low) && (bbox->orig.x <= high) &&
	(y >= bbox->orig.y) && (y <= bbox->corner.y))
      return 0;

    /* Check against right edge */
    y += bbox_width * slope;
    if ((y >= bbox->orig.y) && (y <= bbox->corner.y) &&
	(bbox->corner.x >= low) && (bbox->corner.x <= high))
      return 0;
	
    /* Check against bottom edge */
    if (p1->y < p2->y) {
      low = p1->y;
      high = p2->y;
    }
    else {
      low = p2->y;
      high = p1->y;
    }

    x = p1->x + (bbox->orig.y - p1->y) / slope;
    if ((x >= bbox->orig.x) && (x <= bbox->corner.x) &&
	(bbox->orig.y >= low) && (bbox->orig.y <= high))
      return 0;

    /* Check against top edge */
    x += bbox_height / slope;
    if ((x >= bbox->orig.x) && (x <= bbox->corner.x) &&
	(bbox->corner.y >= low) && (bbox->corner.y <= high))
      return 0;
  }

  return -1;
}


/*
 * ShiftLine --
 *	Given two points describing a line and a distance, return
 *	to points describing a line parallel to it at the given distance.
 *	When looking the line from p1 to p2 the new line will be dist away
 *	on its left. Negative values are allowed for dist, resulting in a line
 *	on the right.
 */
void
ShiftLine(RadarPoint	*p1,
	  RadarPoint	*p2,
	  RadarReal	dist,
	  RadarPoint	*p3,
	  RadarPoint	*p4)
{
  static int	shift_table[129];
  RadarBool	dx_neg, dy_neg;
  int		dx, dy;

  /*
   * Initialize the conversion table.
   */
  if (shift_table[0] == 0) {
    int		i;
    double	tangent, cosine;
    
    for (i = 0; i <= 128; i++) {
      tangent = i/128.0;
      cosine = 128/cos(atan(tangent)) + 0.5;
      shift_table[i] = (int) cosine;
    }
  }

  *p3 = *p1;
  dx = p2->x - p1->x;
  dy = p2->y - p1->y;
  if (dx < 0) {
    dx = -dx;
    dx_neg = True;
  }
  else
    dx_neg = False;
  if (dy < 0) {
    dy = -dy;
    dy_neg = True;
  }
  else {
    dy_neg = False;
  }
  if (dy <= dx) {
    dy = ((dist * shift_table[(dy*128)/dx]) + 64) / 128;
    if (!dx_neg) {
      dy = -dy;
    }
    p3->y += dy;
  }
  else {
    dx = ((dist * shift_table[(dx*128)/dy]) + 64) / 128;
    if (dy_neg) {
      dx = -dx;
    }
    p3->x += dx;
  }

  p4->x = p3->x + (p2->x - p1->x);
  p4->y = p3->y + (p2->y - p1->y);
}


/*
 * IntersectLines --
 *	Given two lines described by two points, compute their intersection.
 *	The function returns True if the lines are not parallel and False
 *	otherwise.
 */
RadarBool
IntersectLines(RadarPoint	*a1,
	       RadarPoint	*a2,
	       RadarPoint	*b1,
	       RadarPoint	*b2,
	       RadarPoint	*pi)
{
  int	dxadyb, dxbdya, dxadxb, dyadyb, p, q;

  dxadyb = (a2->x - a1->x)*(b2->y - b1->y);
  dxbdya = (b2->x - b1->x)*(a2->y - a1->y);
  dxadxb = (a2->x - a1->x)*(b2->x - b1->x);
  dyadyb = (a2->y - a1->y)*(b2->y - b1->y);
  
  if (dxadyb == dxbdya) {
    return False;
  }
  
  p = a1->x*dxbdya - b1->x*dxadyb + (b1->y - a1->y)*dxadxb;
  q = dxbdya - dxadyb;
  if (q < 0) {
    p = -p;
    q = -q;
  }
  if (p < 0) {
    pi->x = -(-2*p + q)/(2*q); /*-(-p + q/2)/q;*/
  }
  else {
    pi->x = (2*p + q)/(2*q); /*(p + q/2)/q;*/
  }
  
  p = a1->y*dxadyb - b1->y*dxbdya + (b1->x - a1->x)*dyadyb;
  q = dxadyb - dxbdya;
  if (q < 0) {
    p = -p;
    q = -q;
  }
  if (p < 0) {
    pi->y = -(-2*p + q)/(2*q); /*-(-p + q/2)/q;*/
  }
  else {
    pi->y = (2*p + q)/(2*q); /*(p + q/2)/q;*/
  }
  
  return True;
}


/*
 * InsetPolygon --
 *	Inset the given polygon by the given amount. The
 *     value can be negative, in this case the polygon will
 *     be outset.
 */
/**** A FINIR ****/
void
InsetPolygon(RadarPoint	*p,
	     int	num_points,
	     RadarDim	inset)
{
  RadarPoint	*p1, *p2;
  RadarPoint	new_p1, new_p2;
  RadarPoint	shift1, shift2;
  int		i, processed_points;

  processed_points = 0;
  
  if ((p->x == p[num_points-1].x) && (p->y == p[num_points-1].y)) {
    num_points--;
  }
  for (p1 = p, p2 = p1+1, i = 0; i < num_points; i++, p1 = p2, p2++) {
    /*
     * Wrap to the first point.
     */
    if (i == num_points-1) {
      p2 = p;
    }
    /*
     * Skip duplicate vertices.
     */
    if ((p2->x == p1->x) && (p2->y == p1->y)) {
      continue;
    }
    
    ShiftLine(p1, p2, inset, &new_p1, &new_p2);

    if (processed_points >= 1) {
    }
  }
}


/*
 * Compute the two corner points of a thick line end.
 * Two points describing the line segment and the width
 * are given as input. If projecting is true this function
 * mimics the X11 line projecting behaviour. The computed
 * end is located around p2.
 */
void
GetButtPoints(RadarPoint	*p1,
	      RadarPoint	*p2,
	      int		width,
	      RadarBool		projecting,
	      RadarPoint	*c1,
	      RadarPoint	*c2)
{
  double	w_2 = width/2.0;
  double	length = hypot(p2->x - p1->x, p2->y - p1->y);
  double	delta_x, delta_y;
  
  if (length == 0.0) {
    c1->x = c2->x = p2->x;
    c1->y = c2->y = p2->y;
  }
  else {
    delta_x = -w_2 * (p2->y - p1->y) / length;
    delta_y = w_2 * (p2->x - p1->x) / length;
    c1->x = p2->x + delta_x;
    c2->x = p2->x - delta_x;
    c1->y = p2->y + delta_y;
    c2->y = p2->y - delta_y;
    if (projecting) {
      c1->x += delta_y;
      c2->x += delta_y;
      c1->y -= delta_x;
      c2->y -= delta_x;
    }
  }
}

/*
 * Compute the inside and outside points of the mitered
 * corner formed by a thick line going through 3 points.
 * If the angle formed by the three points is less than
 * 11 degrees, False is returned an no points are computed.
 * Else True is returned and the points are in c1, c2.
 *
 * If someday the code is switched to REAL coordinates, we
 * must round each coordinate to the nearer integer to mimic
 * the way pixels are drawn. Sample code: floor(p->x+0.5);
 *
 * Hmmm, the switch has been done but not the rounding ;-)
 */
RadarBool
GetMiterPoints(RadarPoint	*p1,
	       RadarPoint	*p2,
	       RadarPoint	*p3,
	       int		width,
	       RadarPoint	*c1,
	       RadarPoint	*c2)
{
  static double	deg11 = (11.0*2.0*M_PI)/360.0;
  double	theta1;		/* angle of p2-p1 segment. */
  double	theta2;		/* angle of p2-p3 segment. */
  double	theta;		/* angle of the joint */
  double	theta3;		/* angle of bisector of the joint toward
				 * the external point of the joint. */
  double	dist;		/* distance of the external points
				 * of the corner from the mid point
				 * p2. */
  double	delta_x, delta_y; /* projection of (dist,theta3) on x
				   * and y. */
  
  if (p2->y == p1->y) {
    theta1 = (p2->x < p1->x) ? 0.0 : M_PI;
  }
  else if (p2->x == p1->x) {
    theta1 = (p2->y < p1->y) ? M_PI/2.0 : -M_PI/2.0;
  }
  else {
    theta1 = atan2(p1->y - p2->y, p1->x - p2->x);
  }
  if (p3->y == p2->y) {
    theta2 = (p3->x > p2->x) ? 0.0 : M_PI;
  }
  else if (p3->x == p2->x) {
    theta2 = (p3->y > p2->y) ? M_PI/2.0 : -M_PI/2.0;
  }
  else {
    theta2 = atan2(p3->y - p2->y, p3->x - p2->x);
  }
  theta = theta1 - theta2;
  if (theta > M_PI) {
    theta -= 2.0*M_PI;
  }
  else if (theta < -M_PI) {
    theta += 2*M_PI;
  }
  if ((theta < deg11) && (theta > -deg11)) {
    return False;
  }
  /*
   * Compute the distance of the internal and external
   * corner points from the intersection p2 (considered
   * at 0,0).
   */
  dist = 0.5*width / sin(0.5*theta);
  dist = ABS(dist);

  /*
   * Compute the angle of the bisector of the joint that
   * goes toward the outside of the joint (the left hand
   * when looking from p1-p2).
   */
  theta3 = (theta1 + theta2)/2.0;
  if (sin(theta3 - (theta1 + M_PI)) < 0.0) {
    theta3 += M_PI;
  }
  
  delta_x = dist * cos(theta3);
  c1->x = p2->x + delta_x;
  c2->x = p2->x - delta_x;
  delta_y = dist * sin(theta3);
  c1->y = p2->y + delta_y;
  c2->y = p2->y - delta_y;
  
  return True;
}

/*
 * Tell where a thick polyline is with respect to an area.
 * Return -1 if the polyline is entirely outside the bbox, 1
 * if it is entirely inside and 0 otherwise. The joints can
 * be specified as JoinMiter, JoinRound, JoinBevel. The ends
 * can be: CapRound, CapButt, CapProjecting.
 */
int
PolylineInBBox(RadarPoint	*points,
	       int		num_points,
	       int		width,
	       int		cap_style,
	       int		join_style,
	       RadarBBox	*bbox)
{
  int		count, inside = -1;
  RadarBool	do_miter_as_bevel;
  RadarPoint	poly[4];
  
  /*
   * If the first point is inside the area, change inside
   * accordingly.
   */
  if ((points[0].x >= bbox->orig.x) && (points[0].x <= bbox->corner.x) &&
      (points[0].y >= bbox->orig.y) && (points[0].y <= bbox->corner.y)) {
    inside = 1;
  }

  /*
   * Now iterate through all the edges. Compute a polygon for
   * each and test it against the area. At each vertex an oval
   * of radius width/2 is also tested to account for round ends
   * and joints.
   */
  do_miter_as_bevel = False;
  for (count = num_points; count >= 2; count--, points++) {

    /*
     * Test a circle around the first point if CapRound or
     * around every joint if JoinRound.
     */
    if (((cap_style == CapRound) && (count == num_points)) ||
	((join_style == JoinRound) && (count != num_points))) {
      if (OvalInBBox(points, width, width, bbox) != inside) {
	return 0;
      }
    }
    /*
     * Build a polygon to represent an edge from the current
     * point to the next. Special cases for the first and
     * last edges to implement the line ends.
     */
    /*
     * First vertex of the edge
     */
    if (count == num_points) {
      GetButtPoints(&points[1], points, width,
		    cap_style == CapProjecting, poly, &poly[1]);
    }
    /*
     * Here we are at a joint starting a new edge. If the
     * joint is mitered, start by carrying over the points
     * from the previous edge. Otherwise compute new points
     * for a butt end.
     */
    else if ((join_style == JoinMiter) && !do_miter_as_bevel) {
      poly[0] = poly[3];
      poly[1] = poly[2];
    }
    else {
      GetButtPoints(&points[1], points, width, 0, poly, &poly[1]);
      /*
       * If the previous joint was beveled (or considered so),
       * check the polygon that fill the bevel. It has more or
       * less an X shape, i.e, it's self intersecting. If this
       * is not ok, it may be necessary to permutte poly[1] &
       * poly[2].
       */
      if ((join_style == JoinBevel) || do_miter_as_bevel) {
	if (PolygonInBBox(poly, 4, bbox) != inside) {
	  return 0;
	}
	do_miter_as_bevel = False;
      }
    }

    /*
     * Opposite vertex of the edge.
     */
    if (count == 2) {
      GetButtPoints(points, &points[1], width, cap_style == CapProjecting,
		    &poly[2], &poly[3]);
    }
    else if (join_style == JoinMiter) {
      if (GetMiterPoints(points, &points[1], &points[2], width,
			 &poly[2], &poly[3]) == False) {
	do_miter_as_bevel = True;
	GetButtPoints(points, &points[1], width, 0, &poly[2], &poly[3]);
      }
    }
    else {
      GetButtPoints(points, &points[1], width, 0, &poly[2], &poly[3]);
    }

    if (PolygonInBBox(poly, 4, bbox) != inside) {
      return 0;
    }
  }

  /*
   * Test a circle around the last point if CapRound.
   */
  if (cap_style == CapRound) {
    if (OvalInBBox(points, width, width, bbox) != inside) {
      return 0;
    }
  }

  return inside;
}


/*
 * Tell where a polygon is with respect to an area.
 * Return -1 if the polygon is entirely outside the bbox, 1
 * if it is entirely inside and 0 otherwise.
 */
int
PolygonInBBox(RadarPoint	*points,
	      int		num_points,
	      RadarBBox		*bbox)
{
  int		inside, count;
  RadarPoint	*p, *head, *first, *second;
  RadarBool	closed;
  
  p = head = points;
  closed = True;
  /*
   * Check to see if closed. If not adjust num_points and
   * record this.
   */
  if ((points[0].x != points[num_points-1].x) ||
      (points[0].y != points[num_points-1].y)) {
    closed = False;
  }

  /*
   * Get the status of the first edge.
   */
  inside = LineInBBox(&points[0], &points[1], bbox);
  if (inside == 0) {
    return 0;
  }
  for (count = num_points; count > 0; p++, count--) {
    first = &p[0];
    /*
     * Pretend the polygon is closed if this is not the case.
     */
    if (!closed && (count == 1)) {
      second = head;
    }
    else {
      second = &p[1];
    }
    
    if (LineInBBox(first, second, bbox) != inside) {
      return 0;
    }
  }

  /*
   * If all the edges are inside the area, the polygon is
   * inside the area. If all the edges are outside, the polygon
   * may completely enclose the area. Test if the origin of
   * the area is inside the polygon to decide this.
   */
  if (inside == 1) {
    return 1;
  }

  printf("PolygonInBBox, np = %d, x = %g, y = %g, dist = %g\n",
	 num_points, bbox->orig.x, bbox->orig.y,
	 PolygonToPointDist(points, num_points, &bbox->orig));
  if (PolygonToPointDist(points, num_points, &bbox->orig) <= 0.0) {
    return 0;
  }
  
  return -1;
}


/*
 * Tell where an oval is with respect to an area.
 * Return -1 if the oval is entirely outside the bbox, 1
 * if it is entirely inside and 0 otherwise.
 */
int
OvalInBBox(RadarPoint	*center,
	   int		width,
	   int		height,
	   RadarBBox	*bbox)
{
  RadarPoint	origin, corner;
  int		w_2, h_2;
  double	delta_x, delta_y;
  
  w_2 = (width+1)/2;
  h_2 = (height+1)/2;
  
  origin.x = center->x - w_2;
  origin.y = center->y - h_2;
  corner.x = center->x + w_2;
  corner.y = center->y + h_2;

  /*
   * This if the oval bbox is completely inside or outside
   * of the area. Trivial case first.
   */
  if ((bbox->orig.x <= origin.x) && (bbox->corner.x >= corner.x) &&
      (bbox->orig.y <= origin.y) && (bbox->corner.y >= corner.y)) {
    return 1;
  }
  if ((bbox->corner.x < origin.x) || (bbox->orig.x > corner.x) ||
      (bbox->corner.y < origin.y) || (bbox->orig.y > corner.y)) {
    return -1;
  }

  /*
   * Then test all sides of the area against the oval center.
   * If the point of a side closest to the center is inside
   * the oval, then the oval intersects the area.
   */

  /*
   * Compute the square of the Y axis distance, reducing
   * the oval to a unit circle to take into account the
   * shape factor.
   */
  delta_y = bbox->orig.y - center->y;
  if (delta_y < 0.0) {
    delta_y = center->y - bbox->corner.y;
    if (delta_y < 0.0) {
      delta_y = 0.0;
    }
  }
  delta_y /= h_2;
  delta_y *= delta_y;

  /*
   * Test left and then right edges.
   */
  delta_x = (bbox->orig.x - center->x) / w_2;
  delta_x *= delta_x;
  if ((delta_x + delta_y) <= 1.0) {
    return 0;
  }
  delta_x = (bbox->corner.x - center->x) / w_2;
  delta_x *= delta_x;
  if ((delta_x + delta_y) <= 1.0) {
    return 0;
  }
  
  /*
   * Compute the square of the X axis distance, reducing
   * the oval to a unit circle to take into account the
   * shape factor.
   */
  delta_x = bbox->orig.x - center->x;
  if (delta_x < 0.0) {
    delta_x = center->x - bbox->corner.x;
    if (delta_x < 0.0) {
      delta_x = 0.0;
    }
  }
  delta_x /= w_2;
  delta_x *= delta_x;
  
  /*
   * Test top and then bottom edges.
   */
  delta_y = (bbox->orig.y - center->y) / h_2;
  delta_y *= delta_y;
  if ((delta_x + delta_y) <= 1.0) {
    return 0;
  }
  delta_y = (bbox->corner.y - center->y) / h_2;
  delta_y *= delta_y;
  if ((delta_x + delta_y) <= 1.0) {
    return 0;
  }

  return -1;
}


/*
 * Tell if a point is in an angular range whose center is 0,0.
 * The range is specified by a starting angle and an angle extent.
 * The use of a double here is important, don't change it. In some
 * case we need to normalize a figure to take care of its shape and
 * the result needs precision.
 */
RadarBool
PointInAngle(int	start_angle,
	     int	angle_extent,
	     RadarPoint	*p)
{
  double	point_angle;
  int		angle_diff;

  if ((p->x == 0) && (p->y == 0)) {
    point_angle = 0.0;
  }
  else {
    point_angle = -atan2(p->y, p->x) * 180.0 / M_PI;
  }
  angle_diff = (REAL_TO_INT(point_angle) - start_angle) % 360;
  if (angle_diff < 0) {
    angle_diff += 360;
  }
  return ((angle_diff <= angle_extent) ||
	  ((angle_extent < 0) && ((angle_diff - 360) >= angle_extent)));
}


/*
 * Return the distance of the given point to the rectangle
 * described by rect. Return negative values for points in
 * the rectangle.
 */
double
RectangleToPointDist(RadarBBox	*bbox,
		     RadarPoint	*p)
{
  double	new_dist, dist;
  RadarPoint	p1, p2;

  p1.x = bbox->orig.x;
  p1.y = p2.y = bbox->orig.y;
  p2.x = bbox->corner.x;
  dist = LineToPointDist(&p1, &p2, p);
  if (dist == 0.0) {
    return 0.0;
  }

  p1 = p2;
  p2.y = bbox->corner.y;
  new_dist = LineToPointDist(&p1, &p2, p);
  dist = MIN(dist, new_dist);
  if (dist == 0.0) {
    return 0.0;
  }

  p1 = p2;
  p2.x = bbox->orig.x;
  new_dist = LineToPointDist(&p1, &p2, p);
  dist = MIN(dist, new_dist);
  if (dist == 0.0) {
    return 0.0;
  }
  
  p1 = p2;
  p2.y = bbox->orig.y;
  new_dist = LineToPointDist(&p1, &p2, p);
  dist = MIN(dist, new_dist);

  if (PointInBBox(bbox, p->x, p->y)) {
    return -dist;
  }
  else {
    return dist;
  }
}


/*
 * Return the distance of the given point to the line
 * described by <xl1,yl1>, <xl2,yl2>..
 */
double
LineToPointDist(RadarPoint	*p1,
		RadarPoint	*p2,
		RadarPoint	*p)
{
  double	x, y;
  int		x_int, y_int;

  /*
   * First compute the closest point on the line. This is done
   * separately for vertical, horizontal, other lines.
   */

  /* Vertical */
  if (p1->x == p2->x) {
    x = p1->x;
    if (p1->y >= p2->y) {
      y_int = MIN(p1->y, p->y);
      y_int = MAX(y_int, p2->y);
    }
    else {
      y_int = MIN(p2->y, p->y);
      y_int = MAX(y_int, p1->y);
    }
    y = y_int;
  }

  /* Horizontal */
  else if (p1->y == p2->y) {
    y = p1->y;
    if (p1->x >= p2->x) {
      x_int = MIN(p1->x, p->x);
      x_int = MAX(x_int, p2->x);
    }
    else {
      x_int = MIN(p2->x, p->x);
      x_int = MAX(x_int, p1->x);
    }
    x = x_int;
  }

  /*
   * Other. Compute its parameters of form y = a1*x + b1 and
   * then compute the parameters of the perpendicular passing
   * through the point y = a2*x + b2, last find the closest point
   * on the segment.
   */
  else {
    double	a1, a2, b1, b2;

    a1 = ((double) (p2->y - p1->y)) / ((double) (p2->x - p1->x));
    b1 = p1->y - a1*p1->x;

    a2 = -1.0/a1;
    b2 = p->y - a2*p->x;

    x = (b2 - b1) / (a1 - a2);
    y = a1*x + b1;

    if (p1->x > p2->x) {
      if (x > p1->x) {
	x = p1->x;
	y = p1->y;
      }
      else if (x < p2->x) {
	x = p2->x;
	y = p2->y;
      }
    }
    else {
      if (x > p2->x) {
	x = p2->x;
	y = p2->y;
      }
      else if (x < p1->x) {
	x = p1->x;
	y = p1->y;
      }
    }
  }

  /* Return the distance */
  return hypot(p->x - x, p->y - y);
}


/*
 * Return the distance of the polygon described by
 * points, to the given point. If the point is
 * inside return values are negative.
 */
double
PolygonToPointDist(RadarPoint	*points,
		   int		num_points,
		   RadarPoint	*p)
{
  double	best_distance;
  int		intersections;
  int		x_int, y_int;
  RadarPoint	*first_point;
  double	x, y, dist;
  RadarPoint	p1, p2;

  /*
   * The algorithm iterates through all the edges of the polygon
   * computing for each the distance to the point and whether a vertical
   * ray starting at the point, intersects the edge. The smallest
   * distance of all edges is stored in best_distance while intersections
   * hold the count of edges to rays intersections. For more informations
   * on how the distance is computed see LineToPointDist.
   */
  best_distance = 1.0e40;
  intersections = 0;

  first_point = points;

  /*
   * Check to see if closed. Adjust num_points to open it (the
   * algorithm always consider a set of points as a closed polygon).
   */
  if ((points[0].x == points[num_points-1].x) &&
      (points[0].y == points[num_points-1].y)) {
    num_points--;
  }

  for ( ; num_points >= 1; num_points--, points++) {
    p1 = points[0];
    /*
     * Wrap over to the first point.
     */
    if (num_points == 1) {
      p2 = *first_point;
    }
    else {
      p2 = points[1];
    }
    
    /*
     * First try to find the closest point on this edge.
     */

    /* Vertical edge */
    if (p1.x == p2.x) {
      x = p1.x;
      if (p1.y >= p2.y) {
	y_int = MIN(p1.y, p->y);
	y_int = MAX(y_int, p2.y);
      }
      else {
	y_int = MIN(p2.y, p->y);
	y_int = MAX(y_int, p1.y);
      }
      y = y_int;
    }

    /* Horizontal edge */
    else if (p1.y == p2.y) {
      y = p1.y;
      if (p1.x >= p2.x) {
	x_int = MIN(p1.x, p->x);
	x_int = MAX(x_int, p2.x);
	if ((p->y < y) && (p->x < p1.x) && (p->x >= p2.x)) {
	  intersections++;
	}
      }
      else {
	x_int = MIN(p2.x, p->x);
	x_int = MAX(x_int, p1.x);
	if ((p->y < y) && (p->x < p2.x) && (p->x >= p1.x)) {
	  intersections++;
	}
      }
      x = x_int;
    }

    /* Other */
    else {
      double	a1, b1, a2, b2;

      a1 = ((double) (p2.y - p1.y)) / ((double) (p2.x - p1.x));
      b1 = p1.y - a1 * p1.x;

      a2 = -1.0/a1;
      b2 = p->y - a2 * p->x;

      x = (b2 - b1)/(a1 - a2);
      y = a1 * x + b1; 

      if (p1.x > p2.x) {
	if (x > p1.x) {
	  x = p1.x;
	  y = p1.y;
	}
	else if (x < p2.x) {
	  x = p2.x;
	  y = p2.y;
	}
      }
      else {
	if (x > p2.x) {
	  x = p2.x;
	  y = p2.y;
	}
	else if (x < p1.x) {
	  x = p1.x;
	  y = p1.y;
	}
      }

      if (((a1 * p->x + b1) > p->y) &&	/* True if point is lower */
	  (p->x >= MIN(p1.x, p2.x)) &&
	  (p->x < MAX(p1.x, p2.x))) {
	intersections++;
      }
    }

    /*
     * Now compute the distance to the closest point and
     * keep it if it is the shortest.
     */
    dist = hypot(p->x - x, p->y - y);
    best_distance = MIN(best_distance, dist);
    /*
     * We can safely escape here if distance is zero.
     */
    if (best_distance == 0.0) {
      return 0.0;
    }
  }

  /*
   * Well, all the edges are processed, if the
   * intersection count is odd the point is inside.
   */
  if (intersections & 0x1) {
    return -best_distance;
  }
  else {
    return best_distance;
  }
}


/*
 * Return the distance of a thick polyline to the
 * given point. Cap and Join parameters are considered
 * in the process.
 */
double
PolylineToPointDist(RadarPoint	*points,
		    int		num_points,
		    int		width,
		    int		cap_style,
		    int		join_style,		   
		    RadarPoint	*p)
{
  RadarBool	miter2bevel = False;
  int		count;
  RadarPoint	*ptr;
  RadarPoint	outline[5];
  double	dist, best_dist, h_width;

  best_dist = 1.0e36;
  h_width = width/2.0;

  for (count = num_points, ptr = points; count >= 2; count--, ptr++) {
    if (((cap_style == CapRound) && (count == num_points)) ||
	((join_style == JoinRound) && (count != num_points))) {
      dist = hypot(ptr->x - p->x, ptr->y - p->y) - h_width;
      if (dist <= 0.0) {
	best_dist = 0.0;
	goto done;
      }
      else if (dist < best_dist) {
	best_dist = dist;
      }
    }
    /*
     * Build the polygonal outline of the current edge.
     */
    if (count == num_points) {
      GetButtPoints(&ptr[1], ptr, width, cap_style==CapProjecting, outline, &outline[1]);
    }
    else if ((join_style == JoinMiter) && !miter2bevel) {
      outline[0] = outline[3];
      outline[1] = outline[2];
    }
    else {
      GetButtPoints(&ptr[1], ptr, width, 0, outline, &outline[1]);
      /*
       * If joints are beveled, check the distance to the polygon
       * that fills the joint.
       */
      if ((join_style == JoinBevel) || miter2bevel) {
	outline[4] = outline[0];
	dist = PolygonToPointDist(outline, 5, p);
	if (dist <= 0.0) {
	  best_dist = 0.0;
	  goto done;
	}
	else if (dist < best_dist) {
	  best_dist = dist;
	}
	miter2bevel = False;
      }
    }
    if (count == 2) {
      GetButtPoints(ptr, &ptr[1], width, cap_style==CapProjecting,
		    &outline[2], &outline[3]);
    }
    else if (join_style == JoinMiter) {
      if (GetMiterPoints(ptr, &ptr[1], &ptr[2], width,
			 &outline[2], &outline[3]) == False) {
	miter2bevel = True;
	GetButtPoints(ptr, &ptr[1], width, 0, &outline[2], &outline[3]);
      }
      /*printf("2=%g+%g, 3=%g+%g\n",
	outline[2].x, outline[2].y, outline[3].x, outline[3].y);*/
    }
    else {
      GetButtPoints(ptr, &ptr[1], width, 0, &outline[2], &outline[3]);
    }
    outline[4] = outline[0];
    /*printf("0=%g+%g, 1=%g+%g, 2=%g+%g, 3=%g+%g, 4=%g+%g\n",
	   outline[0].x, outline[0].y, outline[1].x, outline[1].y,
	   outline[2].x, outline[2].y, outline[3].x, outline[3].y,
	   outline[4].x, outline[4].y);*/
    dist = PolygonToPointDist(outline, 5, p);
    if (dist <= 0.0) {
      best_dist = 0.0;
      goto done;
    }
    else if (dist < best_dist) {
      best_dist = dist;
    }
  }

  /*
   * Test the final point if cap style is round. The code so far
   * has only handled the butt and projecting cases.
   */
  if (cap_style == CapRound) {
    dist = hypot(ptr->x - p->x, ptr->y - p->y) - h_width;
    if (dist <= 0.0) {
      best_dist = 0.0;
      goto done;
    }
    else if (dist < best_dist) {
      best_dist = dist;
    }
  }
  
 done:
  return best_dist;
}


/*
 * Return the distance of the given oval to the point given.
 * The oval is described by its bounding box <xbb,ybb,wbb,hbb>,
 * the thickness of its outline <width>. Return values are negative
 * if the point is inside.
 */
double
OvalToPointDist(RadarPoint	*center,
		int		width,
		int		height,
		unsigned int	line_width,
		RadarPoint	*p)
{
  double	x_delta, y_delta;
  /*  double	x_diameter, y_diameter;*/
  double	scaled_distance;
  double	distance_to_outline;
  double	distance_to_center;

  /*
   * Compute the distance from the point given to the center
   * of the oval. Then compute the same distance in a coordinate
   * system where the oval is a circle with unit radius.
   */

  x_delta = p->x - center->x;
  y_delta = p->y - center->y;
  distance_to_center = hypot(x_delta, y_delta);
  scaled_distance = hypot(x_delta / ((width + line_width) / 2.0),
			  y_delta / ((height + line_width) / 2.0));

  /*
   * If the scaled distance is greater than 1.0 the point is outside
   * the oval. Compute the distance to the edge and convert it back
   * to the original coordinate system. This distance is not much
   * accurate and can overestimate the real distance if the oval is
   * very eccentric.
   */
  if (scaled_distance > 1.0) {
    distance_to_outline = (distance_to_center / scaled_distance) * (scaled_distance - 1.0);
    return distance_to_outline;
  }

  /*
   * The point is inside the oval. Compute the distance as above and check
   * if the point is within the outline.
   */
  if (scaled_distance > 1.0e-10) {
    distance_to_outline = (distance_to_center / scaled_distance) * (1.0 - scaled_distance) - line_width;
  }
  else {
    /*
     * If the point is very close to the center avoid dividing by a
     * very small number, take another method.
     */
    if (width < height)
      distance_to_outline = ((double) (width - line_width)) / 2;
    else
      distance_to_outline = ((double) (height - line_width)) / 2;
  }

  if (distance_to_outline < 0.0)
    return 0.0;
  else
    return -distance_to_outline;
}


/*
 **********************************************************************************
 *
 * GetBezierPoints --
 *	Use recursive subdivision to approximate the curve. The subdivision stops
 *	when the error is under eps.
 *	This algorithm is adaptive, meaning that it computes the minimum number
 *	of segments needed to render each curve part.
 *
 **********************************************************************************
 */
static void
GetBezierPoints(RadarPoint	*controls,
		RadarList	to_points,
		double		eps)
{
  RadarReal	dist2;
  RadarPoint	mid_segm, mid_cord, delta;

  /*
   * Compute distance between cord center and curve at t = 0.5
   */
  mid_segm.x = (controls[0].x + 3*controls[1].x + 3*controls[2].x + controls[3].x) / 8.0;
  mid_segm.y = (controls[0].y + 3*controls[1].y + 3*controls[2].y + controls[3].y) / 8.0;
  mid_cord.x = (controls[0].x + controls[3].x) / 2.0;
  mid_cord.y = (controls[0].y + controls[3].y) / 2.0;
  delta.x = mid_segm.x - mid_cord.x;
  delta.y = mid_segm.y - mid_cord.y;
  dist2 = delta.x*delta.x + delta.y*delta.y;
  
  if (dist2 > eps) {
    RadarPoint	new_controls[4];
    /*
     * Subdivide the curve at t = 0.5
     * and compute each new curve.
     */
    new_controls[0] = controls[0];
    new_controls[1].x = (controls[0].x + controls[1].x) / 2.0;
    new_controls[1].y = (controls[0].y + controls[1].y) / 2.0;
    new_controls[2].x = (controls[0].x + 2*controls[1].x + controls[2].x) / 4.0;
    new_controls[2].y = (controls[0].y + 2*controls[1].y + controls[2].y) / 4.0;
    new_controls[3] = mid_segm;
    GetBezierPoints(new_controls, to_points, eps);
    
    new_controls[0] = mid_segm;
    new_controls[1].x = (controls[1].x + 2*controls[2].x + controls[3].x) / 4.0;
    new_controls[1].y = (controls[1].y + 2*controls[2].y + controls[3].y) / 4.0;
    new_controls[2].x = (controls[2].x + (controls[3].x)) / 2.0;
    new_controls[2].y = (controls[2].y + (controls[3].y)) / 2.0;
    new_controls[3] = controls[3];
    GetBezierPoints(new_controls, to_points, eps);
  }
  else {
    /*
     * Flat enough add the end to the current path.
     * The start should already be there.
     */
    RadarListAdd(to_points, &controls[3], RadarListTail);
  }
}


/*
 **********************************************************************************
 *
 * GetBezierPath --
 *	Compute in to_points a new set of points describing a Bezier path based
 *	on the control points given in from_points.
 *	If more than four points are given, the algorithm iterate over the
 *	set using the last point of a segment as the first point of the next.
 *	If 3 points are left, they are interpreted as a Bezier segment with
 *	coincident internal control points. If 2 points are left a straight
 *	is emitted.
 *
 **********************************************************************************
 */
void
GetBezierPath(RadarList	from_points,
	      RadarList	to_points)
{
  RadarPoint	*fp;
  int		num_fp, i;
  RadarPoint	s[4];
  
  fp = (RadarPoint *) RadarListArray(from_points);
  num_fp = RadarListSize(from_points);

  /*
   * make sure the output vector is empty, then add the first point.
   */
  RadarListEmpty(to_points);
  RadarListAdd(to_points, &fp[0], RadarListTail);

  for (i = 0; i < num_fp; ) {
    if (i < (num_fp-3)) {
      GetBezierPoints(fp, to_points, 1.0);
      fp += 3;
      i += 3;
    }
    else if (i == (num_fp-3)) {
      s[0] = fp[0];
      s[1] = s[2] = fp[1];
      s[3] = fp[2];
      GetBezierPoints(s, to_points, 1.0);
      break;
    }
    else if (i == (num_fp-2)) {
      RadarListAdd(to_points, &fp[1], RadarListTail);
      break;
    }
  }
}

/*
 **********************************************************************************
 *
 * SmoothPathWithBezier --
 *	Compute in to_points a new set of points describing a smoothed path based
 *	on the path given in from_points. The algorithm use Bezier cubic curves.
 *
 **********************************************************************************
 */
void
SmoothPathWithBezier(RadarList	from_points,
		     RadarList	to_points)
{
  RadarPoint		*fp;
  int			num_fp;
  RadarBool		closed;
  RadarPoint		s[4];
  int			i;

  fp = (RadarPoint *) RadarListArray(from_points);
  num_fp = RadarListSize(from_points);

  /*
   * make sure the output vector is empty
   */
  RadarListEmpty(to_points);

  /*
   * If the curve is closed, generates a Bezier curve that
   * spans the closure. Else simply add the first point to
   * the path.
   */
  if ((fp[0].x == fp[num_fp-1].x) && (fp[0].y == fp[num_fp-1].y)) {
    closed = True;
    s[0].x = REAL_TO_INT(0.5*fp[num_fp-2].x + 0.5*fp[0].x);
    s[0].y = REAL_TO_INT(0.5*fp[num_fp-2].y + 0.5*fp[0].y);
    s[1].x = REAL_TO_INT(0.167*fp[num_fp-2].x + 0.833*fp[0].x);
    s[1].y = REAL_TO_INT(0.167*fp[num_fp-2].y + 0.833*fp[0].y);
    s[2].x = REAL_TO_INT(0.833*fp[0].x + 0.167*fp[1].x);
    s[2].y = REAL_TO_INT(0.833*fp[0].y + 0.167*fp[1].y);
    s[3].x = REAL_TO_INT(0.5*fp[0].x + 0.5*fp[1].x);
    s[3].y = REAL_TO_INT(0.5*fp[0].y + 0.5*fp[1].y);
    RadarListAdd(to_points, s, RadarListTail);
    GetBezierPoints(s, to_points, 1.0);
  }
  else {
    closed = False;
    RadarListAdd(to_points, &fp[0], RadarListTail);
  }

  for (i = 2; i < num_fp; i++, fp++) {
    /*
     * Setup the first two control points. This differ
     * for first segment of open curves.
     */
    if ((i == 2) && !closed) {
      s[0] = fp[0];
      s[1].x = REAL_TO_INT(0.333*fp[0].x + 0.667*fp[1].x);
      s[1].y = REAL_TO_INT(0.333*fp[0].y + 0.667*fp[1].y);
    }
    else {
      s[0].x = REAL_TO_INT(0.5*fp[0].x + 0.5*fp[1].x);
      s[0].y = REAL_TO_INT(0.5*fp[0].y + 0.5*fp[1].y);
      s[1].x = REAL_TO_INT(0.167*fp[0].x + 0.833*fp[1].x);
      s[1].y = REAL_TO_INT(0.167*fp[0].y + 0.833*fp[1].y);
    }

    /*
     * Setup the last two control points. This also differ
     * for last segment of open curves.
     */
    if ((i == num_fp-1) && !closed) {
      s[2].x = REAL_TO_INT(0.667*fp[1].x + 0.333*fp[2].x);
      s[2].y = REAL_TO_INT(0.667*fp[1].y + 0.333*fp[2].y);
      s[3] = fp[2];
    }
    else {
      s[2].x = REAL_TO_INT(0.833*fp[1].x + 0.167*fp[2].x);
      s[2].y = REAL_TO_INT(0.833*fp[1].y + 0.167*fp[2].y);
      s[3].x = REAL_TO_INT(0.5*fp[1].x + 0.5*fp[2].x);
      s[3].y = REAL_TO_INT(0.5*fp[1].y + 0.5*fp[2].y);
    }

    /*
     * If the first two points or the last two are equal
     * output the last control point. Else generate the
     * Bezier curve.
     */
    if (((fp[0].x == fp[1].x) && (fp[0].y == fp[1].y)) ||
	((fp[1].x == fp[2].x) && (fp[1].y == fp[2].y))) {
      RadarListAdd(to_points, &s[3], RadarListTail);
    }
    else {
      GetBezierPoints(s, to_points, 1.0);
    }
  }
}


/*
 **********************************************************************************
 *
 * GetLineEnd --
 *	Compute the points describing the given line end style at point p1 for
 *	the line p1,p2. Point p1 is adjusted to fit the line end.
 *	If bbox is non null, it is filled with the bounding box of the end.
 *
 * For the time being this procedure handles open/filled arrows.
 *
 * Here are the parameters describing arrows.
 *
 *              *  | ARROW_SHAPE_C
 *            **   |
 *          * ***************************
 *        *  *
 *      *   *  +p1                      +p2
 *      | * |*
 *      |   * ***************************
 *      |   | **
 *      |   |   *
 *      |   |   |
 *      |---|   | ARROW_SHAPE_A
 *      |       |
 *      |-------| ARROW_SHAPE_B
 *
 **********************************************************************************
 */
void
GetLineEnd(RadarPoint	*p1,
	   RadarPoint	*p2,
	   unsigned int	line_width,
	   int		cap_style,
	   LineEnd	end_style,
	   RadarPoint	*points)
{
  RadarReal	dx, dy, length, temp, backup;
  RadarReal	frac_height, sin_theta, cos_theta;
  RadarReal	vert_x, vert_y;
  RadarReal	shape_a, shape_b, shape_c;
  
  if (end_style != NULL) {
    shape_a = end_style->shape_a + 0.001;
    shape_b = end_style->shape_b + 0.001;
    shape_c = end_style->shape_c + line_width/2.0 + 0.001;

    frac_height = (line_width/2.0) / shape_c;
    dx = p1->x - p2->x;
    dy = p1->y - p2->y;
    length = hypot(dx, dy);
    if (length == 0) {
      sin_theta = cos_theta = 0.0;
    }
    else {
      sin_theta = dy/length;
      cos_theta = dx/length;
    }

    if (cap_style != CapProjecting) {
      temp = frac_height;
    }
    else {
      temp = line_width / shape_c;
    }
    backup = temp * shape_b + shape_a * (1.0 - temp) / 2.0;
    points[0].x = points[5].x = p1->x + backup * cos_theta;
    points[0].y = points[5].y = p1->y + backup * sin_theta;
    
    vert_x = points[0].x - shape_a*cos_theta;
    vert_y = points[0].y - shape_a*sin_theta;
    temp = shape_c*sin_theta;
    points[1].x = REAL_TO_INT(points[0].x - shape_b*cos_theta + temp);
    points[4].x = REAL_TO_INT(points[1].x - 2*temp);
    temp = shape_c*cos_theta;
    points[1].y = REAL_TO_INT(points[0].y - shape_b*sin_theta - temp);
    points[4].y = REAL_TO_INT(points[1].y + 2*temp);
    points[2].x = REAL_TO_INT(points[1].x*frac_height + vert_x*(1.0-frac_height));
    points[2].y = REAL_TO_INT(points[1].y*frac_height + vert_y*(1.0-frac_height));
    points[3].x = REAL_TO_INT(points[4].x*frac_height + vert_x*(1.0-frac_height));
    points[3].y = REAL_TO_INT(points[4].y*frac_height + vert_y*(1.0-frac_height));
  }
}