Anti-podal Pairs, Set Width, Set Diameter

In order to use the following procedures with an arbitrary set of points in the plane, the convex hull of the points should first be computed.

/*
 * Anti-podal pairs
 *   Given a convex polygon P, a line of support L is a line intersecting P and
 *   such that the interior of P lies to one side of L
 *   If two points belonging to P admit parallel lines of support, then they 
 *   form an anti-podal pair.
 */

/* anti_podal_pairs(Points, Pairs) is true if Points are the consecutive     */
/*   vertices going clockwise of a convex polygon with the first vertex      */
/*   repeated as the last vertex, and Pairs are the anti-podal pairs of the  */
/*   convex polygon.                                                         */
/*   This uses a "rotating caliper" algorithm.                               */ 
/* e.g.                                                                      */
/* anti_podal_pairs([p(0,4),p(2,6),p(11,4),p(11,2),p(10,0),p(5,0),p(0,4)],   */
/*                  [r(p(2,6),p(10,0)),r(p(11,4),p(5,0)),r(p(11,2),p(0,4)),  */
/*                   r(p(10,0),p(0,4)),r(p(5,0),p(2,6)),r(p(0,4),p(11,4))]). */
anti_podal_pairs([P,P], [r(P,P)]):-!.        % Convex polygon has 1 vertex
anti_podal_pairs([P1,P2,P1], [r(P1,P2)]):-!. % Convex polygon has 2 vertices
anti_podal_pairs(Ps, Pairs):-        % Convex polygon has 3 or more vertices
  Ps=[_,P|Qs],
  anti_podal_pairs_1(Ps, [P|Qs], Qs, Pairs).

anti_podal_pairs_1([_], _, _, []):-!.
anti_podal_pairs_1([P1,P2|Ps], PsCopy, Qs, [r(P2,Q)|Pairs]):-
  furthest_vertex_from_line(Qs, void, PsCopy, P1, P2, 0, [Q|Qs1]),
  anti_podal_pairs_1([P2|Ps], PsCopy, [Q|Qs1], Pairs).

/* furthest_vertex_from_line(Qs, void, PsCopy, P1, P2, 0, [Q|Qs1]) is true   */
/*   if Q is the vertex in Qs which is furthest from the line through P1 and */
/*   P2, and Qs1 are the remaining vertices in Qs. PsCopy is a copy of the   */
/*   original list of vertices of the convex polygon, and allows wraparound  */
/*   from the last vertex to the first vertex.                               */
furthest_vertex_from_line([], Qprev, PsCopy, P1, P2, Area, Qs1):-!,
  furthest_vertex_from_line(PsCopy, Qprev, [], P1, P2, Area, Qs1).
furthest_vertex_from_line([Q|Qs], _, PsCopy, P1, P2, Area, Qs1):-
  area(P2, P1, Q, Area1), /* Always anticlockwise, so always positive */
  Area1 >= Area, !,
  furthest_vertex_from_line(Qs, Q, PsCopy, P1, P2, Area1, Qs1).
furthest_vertex_from_line(Qs, Qprev, _, _, _, _, [Qprev|Qs]).

/* area(Pa, Pb, Pc, Area) is true if Area is twice the signed area of the    */
/*   triangle determined by Pa, Pb and Pc. The area is positive if Pa, Pb    */
/*   and Pc are oriented anti-clockwise, negative if clockwise, and zero if  */
/*   the points are collinear.                                               */
area(p(Xa,Ya), p(Xb,Yb), p(Xc,Yc), Area):-
  Area is (Xb-Xa) * (Yc-Ya) - (Xc-Xa) * (Yb-Ya).

/*
 * Set Width
 *   The width of a set of points is the minimum distance between parallel
 *   lines of support of the points
 */

/* set_width(Points, Point1, Point2, Point3, Width) is true if Points are    */
/*   the consecutive vertices in a clockwise sense of a convex polygon with  */
/*   the first vertex repeated as the last vertex, and the Width is the set  */
/*   width, which is the distance of Point3 from the line through the Point1 */
/*   and Point2.                                                             */
/*   This uses a "rotating caliper" algorithm.                               */ 
/* e.g. set_width([p(0,4),p(2,6),p(11,4),p(11,2),p(10,0),p(5,0),p(0,4)],     */
/*                P1, P2, P3, Width)                                         */
/*      gives P1=p(10,0), P2=p(5,0), P3=p(2,6), Width=6                      */
set_width([P,P], P, P, P, 0.0):-!.         % Convex polygon has 1 vertex
set_width([P1,P2,P1], P1, P2, P1, 0.0):-!. % Convex polygon has 2 vertices
set_width(Ps, P1, P2, P3, Width):- % Convex polygon has 3 or more vertices
  Ps=[_,P|Qs],
  set_width_1(Ps, [P|Qs], Qs, 1.0E9, void, t(P1,P2,P3)),
  area(P2, P1, P3, Area),
  distance_squared(P1, P2, DistSq),
  Width is Area / sqrt(DistSq).

set_width_1([_], _, _, _, Triple, Triple):-!.
set_width_1([P1,P2|Ps], PsCopy, Qs, DistSq0, Triple0, Triple):-
  furthest_vertex_from_line(Qs, void, PsCopy, P1, P2, 0, Qs1),
  set_width_2([P1,P2|Ps], PsCopy, Qs1, DistSq0, Triple0, Triple).
  
set_width_2([P1,P2|Ps], PsCopy, [Q|Qs], DistSq0, _, Triple):-
  distance_squared(P2, Q, DistSq1),
  DistSq1 < DistSq0, !,
  set_width_1([P2|Ps], PsCopy, [Q|Qs], DistSq1, t(P1,P2,Q), Triple).
set_width_2([_|Ps], PsCopy, Qs, DistSq0, Triple0, Triple):-
  set_width_1(Ps, PsCopy, Qs, DistSq0, Triple0, Triple).

/* distance_squared(p(X1,Y1), p(X2,Y2), DistanceSquared) is true if          */
/*   DistanceSquared is the square of the distance between the points        */
/*   (X1,Y1) and (X2,Y2).                                                    */
distance_squared(p(X1,Y1), p(X2,Y2), DistanceSquared):-
  X1_X2 is X1 - X2,
  Y1_Y2 is Y1 - Y2,
  DistanceSquared is X1_X2*X1_X2 + Y1_Y2*Y1_Y2.

/*
 * Set Diameter
 *   The set diameter of a set of points is the maximum distance between any
 *   two points in the set
 *   The set diameter of a set of points is the maximum distance between
 *   parallel lines of support of the points
 */  

/* set_diameter(Points, Point1, Point2, Diameter) is true if Points are the  */
/*   consecutive vertices in a clockwise sense of a convex polygon with the  */
/*   first vertex repeated as the last vertex, and Diameter is the set       */
/*   diameter, which is the distance between Point1 and Point2.              */
/*   This uses a "rotating caliper" algorithm.                               */
/* e.g. set_diameter([p(0,4),p(2,6),p(11,4),p(11,2),p(10,0),p(5,0),p(0,4)],  */
/*                   P1, P2, Diameter)                                       */
/*      gives P1=p(11,2), P2=p(0,4), Diameter=11.180339887                   */
set_diameter([P,P], P, P, 0.0):-!.             % Convex polygon has 1 vertex
set_diameter([P1,P2,P1], P1, P2, Diameter):-!, % Convex polygon has 2 vertices
  distance_squared(P1, P2, DistSq),
  Diameter is sqrt(DistSq).
set_diameter(Ps, P1, P2, Diameter):-   % Convex polygon has 3 or more vertices
  Ps=[_,P|Qs],
  set_diameter_1(Ps, [P|Qs], Qs, diam(void,void,0.0), diam(P1,P2,DistSq)),
  Diameter is sqrt(DistSq).

set_diameter_1([_], _, _, Best, Best):-!.
set_diameter_1([P1,P2|Ps], PsCopy, Qs, BestSoFar, Best):-
  furthest_vertex_from_line(Qs, void, PsCopy, P1, P2, 0, Qs1),
  set_diameter_2([P2|Ps], PsCopy, Qs1, BestSoFar, Best).

set_diameter_2([P2|Ps], PsCopy, [Q|Qs], diam(_,_,DistSq0), Best):-
  distance_squared(P2, Q, DistSq1),
  DistSq1 > DistSq0, !,
  set_diameter_1([P2|Ps], PsCopy, [Q|Qs], diam(P2,Q,DistSq1), Best).
set_diameter_2(Ps, PsCopy, Qs, BestSoFar, Best):-
  set_diameter_1(Ps, PsCopy, Qs, BestSoFar, Best).

LPA Index     Home Page