/*----------------------------------------------------------------------------*/
/* pamrubber.c - transform images using Rubber Sheeting algorithm
** see: http://www.schaik.com/netpbm/rubber/
**
** Copyright (C) 2011 by Willem van Schaik (willem@schaik.com)
**
** Permission to use, copy, modify, and distribute this software and its
** documentation for any purpose and without fee is hereby granted, provided
** that the above copyright notice appear in all copies and that both that
** copyright notice and this permission notice appear in supporting
** documentation. This software is provided "as is" without express or
** implied warranty.
*/
/*----------------------------------------------------------------------------*/
#include <stdio.h>
#include <stdlib.h>
#include <unistd.h>
#include <assert.h>
#include <limits.h>
#include <math.h>
#include <time.h>
#include "pm_c_util.h"
#include "mallocvar.h"
#include "shhopt.h"
#include "pam.h"
#include "pamdraw.h"
typedef struct {
double x;
double y;
} point;
typedef struct {
point p1;
point p2;
} line;
typedef struct {
point p1;
point p2;
point p3;
} triangle;
typedef struct {
point tl; /* top left */
point tr; /* top right */
point bl; /* bottom left */
point br; /* bottom right */
} quadrilateral;
struct cmdlineInfo {
unsigned int nCP;
point oldCP[4];
point newCP[4];
const char * fileName;
unsigned int quad;
unsigned int tri;
unsigned int frame;
unsigned int linear;
unsigned int verbose;
unsigned int randseedSpec;
unsigned int randseed;
};
static void
parseCmdline(int argc, const char ** argv,
struct cmdlineInfo * const cmdlineP) {
/* parse all parameters from the command line */
unsigned int option_def_index;
char * endptr;
unsigned int i;
unsigned int nCP;
/* instructions to optParseOptions3 on how to parse our options. */
optEntry * option_def;
optStruct3 opt;
MALLOCARRAY_NOFAIL(option_def, 100);
option_def_index = 0; /* incremented by OPTENT3 */
OPTENT3(0, "quad", OPT_FLAG, NULL, &cmdlineP->quad, 0);
OPTENT3(0, "tri", OPT_FLAG, NULL, &cmdlineP->tri, 0);
OPTENT3(0, "frame", OPT_FLAG, NULL, &cmdlineP->frame, 0);
OPTENT3(0, "linear", OPT_FLAG, NULL, &cmdlineP->linear, 0);
OPTENT3(0, "verbose", OPT_FLAG, NULL, &cmdlineP->verbose, 0);
OPTENT3(0, "randseed", OPT_UINT, &cmdlineP->randseed,
&cmdlineP->randseedSpec, 0);
OPTENT3(0, "randomseed", OPT_UINT, &cmdlineP->randseed,
&cmdlineP->randseedSpec, 0);
opt.opt_table = option_def;
opt.short_allowed = FALSE; /* we have no short (old-fashioned) options */
opt.allowNegNum = FALSE; /* we have no parms that are negative numbers */
/* uses and sets argc, argv, and some of *cmdlineP and others. */
pm_optParseOptions3(&argc, (char **)argv, opt, sizeof(opt), 0);
if (!cmdlineP->tri && !cmdlineP->quad)
pm_error("You must specify either -tri or -quad");
if (cmdlineP->tri && cmdlineP->quad)
pm_error("You may not specify both -tri and -quad");
/* Parameters are the control points (in qty of 4) and possibly a file name
*/
nCP = (argc-1) / 4;
if (nCP > 4)
pm_error("Too many arguments: %u. Arguments are "
"control point coordinates and an optional file name, "
"with a maximum of 4 control points", argc-1);
cmdlineP->nCP = nCP;
assert(nCP <= ARRAY_SIZE(cmdlineP->oldCP));
assert(nCP <= ARRAY_SIZE(cmdlineP->newCP));
for (i = 0; i < nCP; ++i) {
cmdlineP->oldCP[i].x = strtol(argv[i * 2 + 1], &endptr, 10);
cmdlineP->oldCP[i].y = strtol(argv[i * 2 + 2], &endptr, 10);
cmdlineP->newCP[i].x = strtol(argv[4 * nCP / 2 + i * 2 + 1],
&endptr, 10);
cmdlineP->newCP[i].y = strtol(argv[4 * nCP / 2 + i * 2 + 2],
&endptr, 10);
}
if (argc - 1 == 4 * nCP)
cmdlineP->fileName = "-";
else if (argc - 2 == 4 * nCP)
cmdlineP->fileName = argv[nCP * 4 + 1];
else
pm_error("Invalid number of arguments. Arguments are "
"control point coordinates and an optional file name, "
"so there must be a multiple of 4 or a multiple of 4 "
"plus 1.");
}
/* global variables */
static int nCP;
static point oldCP[4];
static point newCP[4];
static int nTri;
static triangle tri1s[10];
static triangle tri2s[10];
static quadrilateral quad1;
static quadrilateral quad2;
static tuple black;
static point
makepoint(double const x,
double const y) {
point retval;
retval.x = x;
retval.y = y;
return retval;
}
static double
distance(point const p1,
point const p2) {
return sqrt(SQR(p1.x - p2.x) + SQR(p1.y - p2.y));
}
static line
makeline(point const p1,
point const p2) {
line retval;
retval.p1 = p1;
retval.p2 = p2;
return retval;
}
static bool
intersect(line * const l1P,
line * const l2P,
point * const pP) {
bool cross;
if (((l2P->p2.y - l2P->p1.y) * (l1P->p2.x - l1P->p1.x) -
(l2P->p2.x - l2P->p1.x) * (l1P->p2.y - l1P->p1.y)) == 0) {
/* parallel lines */
cross = false;
if ((l1P->p1.x == l1P->p2.x) && (l2P->p1.x == l2P->p2.x)) {
/* two vertical lines */
pP->x = (l1P->p1.x + l2P->p1.x) / 2.0;
pP->y = 1e10;
} else if ((l1P->p1.y == l1P->p2.y) && (l2P->p1.y == l2P->p2.y)) {
/* two horizontal lines */
pP->x = 1e10;
pP->y = (l1P->p1.y + l2P->p1.y) / 2.0;
} else {
if (fabs(l1P->p2.y - l1P->p1.y) > fabs(l1P->p2.x - l1P->p1.x)) {
/* steep slope */
pP->y = 1e10;
pP->x = (l1P->p2.x - l1P->p1.x) / (l1P->p2.y - l1P->p1.y)
* 1e10;
} else {
/* even slope */
pP->x = 1e10;
pP->y = (l1P->p2.y - l1P->p1.y) / (l1P->p2.x - l1P->p1.x)
* 1e10;
}
}
} else {
/* intersecting lines */
double const ua =
((l2P->p2.x - l2P->p1.x) * (l1P->p1.y - l2P->p1.y)
- (l2P->p2.y - l2P->p1.y) * (l1P->p1.x - l2P->p1.x))
/ ((l2P->p2.y - l2P->p1.y)
* (l1P->p2.x - l1P->p1.x) - (l2P->p2.x - l2P->p1.x)
* (l1P->p2.y - l1P->p1.y));
double const ub =
((l1P->p2.x - l1P->p1.x) * (l1P->p1.y - l2P->p1.y)
- (l1P->p2.y - l1P->p1.y) * (l1P->p1.x - l2P->p1.x))
/ ((l2P->p2.y - l2P->p1.y)
* (l1P->p2.x - l1P->p1.x) - (l2P->p2.x - l2P->p1.x)
* (l1P->p2.y - l1P->p1.y));
pP->x = l1P->p1.x + ua * (l1P->p2.x - l1P->p1.x);
pP->y = l1P->p1.y + ua * (l1P->p2.y - l1P->p1.y);
if ((ua >= 0.0) && (ua <= 1.0) && (ub >= 0.0) && (ub <= 1.0))
cross = true;
else
cross = false;
}
return cross;
}
static triangle
maketriangle(point const p1,
point const p2,
point const p3) {
triangle retval;
retval.p1 = p1;
retval.p2 = p2;
retval.p3 = p3;
return retval;
}
static int
insidetri(triangle * const triP,
point const p) {
int cnt;
cnt = 0; /* initial value */
if ((((triP->p1.y <= p.y) && (p.y < triP->p3.y))
|| ((triP->p3.y <= p.y) && (p.y < triP->p1.y)))
&&
(p.x < (triP->p3.x - triP->p1.x) * (p.y - triP->p1.y)
/ (triP->p3.y - triP->p1.y) + triP->p1.x))
cnt = !cnt;
if ((((triP->p2.y <= p.y) && (p.y < triP->p1.y))
|| ((triP->p1.y <= p.y) && (p.y < triP->p2.y)))
&&
(p.x < (triP->p1.x - triP->p2.x) * (p.y - triP->p2.y)
/ (triP->p1.y - triP->p2.y) + triP->p2.x))
cnt = !cnt;
if ((((triP->p3.y <= p.y) && (p.y < triP->p2.y))
|| ((triP->p2.y <= p.y) && (p.y < triP->p3.y)))
&&
(p.x < (triP->p2.x - triP->p3.x) * (p.y - triP->p3.y)
/ (triP->p2.y - triP->p3.y) + triP->p3.x))
cnt = !cnt;
return cnt;
}
static bool
windtriangle(triangle * const tP,
point const p1,
point const p2,
point const p3) {
point f, c;
line le, lv;
bool cw;
/* find cross of vertical through p3 and the edge p1-p2 */
f.x = p3.x;
f.y = -1.0;
lv = makeline(p3, f);
le = makeline(p1, p2);
intersect(&le, &lv, &c);
/* check for clockwise winding triangle */
if (((p1.x > p2.x) && (p3.y < c.y)) ||
((p1.x < p2.x) && (p3.y > c.y))) {
*tP = maketriangle(p1, p2, p3);
cw = true;
} else { /* p1/2/3 were counter clockwise */
*tP = maketriangle(p1, p3, p2);
cw = false;
}
return cw;
}
static double
tiny(void) {
if (rand() % 2)
return +1E-6 * (double) ((rand() % 90) + 9);
else
return -1E-6 * (double) ((rand() % 90) + 9);
}
static void
angle(point * const p1P,
point * const p2P) {
/*----------------------------------------------------------------------------
Move *p2P slightly if necessary to make sure the line (*p1P, *p2P)
is not horizontal or vertical.
-----------------------------------------------------------------------------*/
if (p1P->x == p2P->x) { /* vertical line */
p2P->x += tiny();
}
if (p1P->y == p2P->y) { /* horizontal line */
p2P->y += tiny();
}
}
static void
sideTriangleVerticalEdge(unsigned int const n,
triangle * const trig1P,
point const p11,
point const p12,
point const p13,
point const p14,
point const r11,
point const r12,
triangle * const trig2P,
point const p21,
point const p22,
point const p23,
point const p24,
point const r21,
point const r22) {
if (((n >= 4) && (r11.x < p11.x)
&& (p14.x < p13.x) && (p14.x < p12.x)
&& (p14.x < p11.x)) /* left edge */
||
((n >= 4) && (r11.x > p11.x)
&& (p14.x > p13.x) && (p14.x > p12.x)
&& (p14.x > p11.x))) /* right edge */ {
*trig1P = maketriangle(r11, r12, p14);
*trig2P = maketriangle(r21, r22, p24);
} else if (((n >= 3) && (r11.x < p11.x) && (p13.x < p12.x)
&& (p13.x < p11.x)) /* left edge */
||
((n >= 3) && (r11.x > p11.x) && (p13.x > p12.x)
&& (p13.x > p11.x))) /* right edge */ {
*trig1P = maketriangle(r11, r12, p13);
*trig2P = maketriangle(r21, r22, p23);
} else if (((n >= 2) && (r11.x < p11.x)
&& (p12.x < p11.x)) /* left edge */
||
((n >= 2) && (r11.x > p11.x)
&& (p12.x > p11.x))) /* right edge */ {
*trig1P = maketriangle(r11, r12, p12);
*trig2P = maketriangle(r21, r22, p22);
} else if (n >= 1) {
*trig1P = maketriangle(r11, r12, p11);
*trig2P = maketriangle(r21, r22, p21);
}
}
static void
sideTriangleHorizontalEdge(unsigned int const n,
triangle * const trig1P,
point const p11,
point const p12,
point const p13,
point const p14,
point const r11,
point const r12,
triangle * const trig2P,
point const p21,
point const p22,
point const p23,
point const p24,
point const r21,
point const r22) {
if (((n >= 4) && (r11.y < p11.y) && (p14.y < p13.y)
&& (p14.y < p12.y) && (p14.y < p11.y)) /* top edge */
||
((n >= 4) && (r11.y > p11.y) && (p14.y > p13.y)
&& (p14.y > p12.y) && (p14.y > p11.y))) /* bottom edge */ {
*trig1P = maketriangle(r11, r12, p14);
*trig2P = maketriangle(r21, r22, p24);
} else if (((n >= 3) && (r11.y < p11.y) && (p13.y < p12.y)
&& (p13.y < p11.y)) /* top edge */
||
((n >= 3) && (r11.y > p11.y) && (p13.y > p12.y)
&& (p13.y > p11.y))) /* bottom edge */ {
*trig1P = maketriangle(r11, r12, p13);
*trig2P = maketriangle(r21, r22, p23);
} else if (((n >= 2) && (r11.y < p11.y)
&& (p12.y < p11.y)) /* top edge */
||
((n >= 2) && (r11.y > p11.y)
&& (p12.y > p11.y))) /* bottom edge */ {
*trig1P = maketriangle(r11, r12, p12);
*trig2P = maketriangle(r21, r22, p22);
} else if (n >= 1) {
*trig1P = maketriangle(r11, r12, p11);
*trig2P = maketriangle(r21, r22, p21);
}
}
static void
sideTriangle(unsigned int const n,
triangle * const trig1P,
point const p11,
point const p12,
point const p13,
point const p14,
point const r11,
point const r12,
triangle * const trig2P,
point const p21,
point const p22,
point const p23,
point const p24,
point const r21,
point const r22) {
if (fabs (r11.x - r12.x) < 1.0)
sideTriangleVerticalEdge(n,
trig1P, p11, p12, p13, p14, r11, r12,
trig2P, p21, p22, p23, p24, r21, r22);
else if (fabs (r11.y - r12.y) < 1.0)
sideTriangleHorizontalEdge(n,
trig1P, p11, p12, p13, p14, r11, r12,
trig2P, p21, p22, p23, p24, r21, r22);
}
static void
edgeTriangle(triangle * const trig1P,
point const p11,
point const p12,
point const tl1,
point const tr1,
point const bl1,
point const br1,
triangle * const trig2P,
point const p21,
point const p22,
point const tl2,
point const tr2,
point const bl2,
point const br2) {
if ((p11.x < p12.x) && (p11.y < p12.y)) {
/* up/left to down/right */
*trig1P = maketriangle(tr1, p12, p11);
*trig2P = maketriangle(tr2, p22, p21);
} else if ((p11.x > p12.x) && (p11.y > p12.y)) {
/* down/right to up/left */
*trig1P = maketriangle(bl1, p12, p11);
*trig2P = maketriangle(bl2, p22, p21);
} else if ((p11.x < p12.x) && (p11.y > p12.y)) {
/* down/left to up/right */
*trig1P = maketriangle(tl1, p12, p11);
*trig2P = maketriangle(tl2, p22, p21);
} else if ((p11.x > p12.x) && (p11.y < p12.y)) {
/* up/right to down/left */
*trig1P = maketriangle(br1, p12, p11);
*trig2P = maketriangle(br2, p22, p21);
}
}
static quadrilateral
quadRect(double const lft,
double const rgt,
double const top,
double const bot) {
quadrilateral retval;
retval.tl = makepoint(lft, top);
retval.tr = makepoint(rgt, top);
retval.bl = makepoint(lft, bot);
retval.br = makepoint(rgt, bot);
return retval;
}
static void
quadCornerSized(point const p0,
point const p1,
point const p2,
point const p3,
point const mid,
quadrilateral * const quadP,
triangle * const triP) {
/* P0-P1 and P2-P3 are the diagonals */
/* P0-P1 are further apart than P2-P3 */
if ((p0.x < p1.x) && (p0.y < p1.y)) {
/* p0 is top-left */
quadP->tl = p0; quadP->br = p1;
if (windtriangle(triP, p0, p2, p1)) {
/* p2 is top-right */
quadP->tr = p2; quadP->bl = p3;
} else {
/* p3 is top-right */
quadP->tr = p3; quadP->bl = p2;
}
} else if ((p0.x > p1.x) && (p0.y < p1.y)) {
/* p0 is top-right */
quadP->tr = p0; quadP->bl = p1;
if (windtriangle(triP, p0, p2, p1)) {
/* p2 is bottom-right */
quadP->br = p2; quadP->tl = p3;
} else {
/* p3 is bottom-right */
quadP->br = p3; quadP->tl = p2;
}
} else if ((p0.x < p1.x) && (p0.y > p1.y)) {
/* p0 is bottom-left */
quadP->bl = p0; quadP->tr = p1;
if (windtriangle(triP, p0, p2, p1)) {
/* p2 is top-left */
quadP->tl = p2; quadP->br = p3;
} else {
/* p3 is top-left */
quadP->tl = p3; quadP->br = p2;
}
} else if ((p0.x > p1.x) && (p0.y > p1.y)) {
/* p0 is bottom-right */
quadP->br = p0; quadP->tl = p1;
if (windtriangle(triP, p0, p2, p1)) {
/* p2 is bottom-left */
quadP->bl = p2; quadP->tr = p3;
} else {
/* p3 is bottom-left */
quadP->bl = p3; quadP->tr = p2;
}
}
}
static void
quadCorner(point const p0,
point const p1,
point const p2,
point const p3,
point const mid,
quadrilateral * const quadP,
triangle * const triP) {
/* p0-p1 and p2-p3 are the diagonals */
if (fabs(p0.x - p1.x) + fabs(p0.y - p1.y) >=
fabs(p2.x - p3.x) + fabs(p2.y - p3.y)) {
quadCornerSized(p0, p1, p2, p3, mid, quadP, triP);
} else {
quadCornerSized(p2, p3, p0, p1, mid, quadP, triP);
}
}
static pamd_drawproc frameDrawproc;
static void
frameDrawproc (tuple ** const tuples,
unsigned int const cols,
unsigned int const rows,
unsigned int const depth,
sample const maxval,
pamd_point const p,
const void * const clientdata) {
int yy;
for (yy = p.y - 1; yy <= p.y + 1; ++yy) {
int xx;
for (xx = p.x - 1; xx <= p.x + 1; ++xx) {
if (xx >= 0 && xx < cols && yy >= 0 && yy < rows) {
unsigned int i;
for (i = 0; i < depth; ++i)
tuples[yy][xx][i] = (sample) *((tuple *) clientdata + i);
}
}
}
}
static void
drawExtendedLine(const struct pam * const pamP,
tuple ** const outTuples,
point const p1,
point const p2) {
pamd_point const p1ext =
pamd_makePoint(p1.x - 10 * (p2.x - p1.x),
p1.y - 10 * (p2.y - p1.y));
pamd_point const p2ext =
pamd_makePoint(p2.x + 10 * (p2.x - p1.x),
p2.y + 10 * (p2.y - p1.y));
pamd_line(outTuples, pamP->width, pamP->height, pamP->depth, pamP->maxval,
p1ext, p2ext, frameDrawproc, black);
}
static pamd_point
clippedPoint(const struct pam * const pamP,
point const p) {
int const roundedX = ROUND(p.x);
int const roundedY = ROUND(p.x);
int clippedX, clippedY;
assert(pamP->width >= 2);
assert(pamP->height >= 2);
if (roundedX <= 0)
clippedX = 1;
else if (roundedX > pamP->width - 1)
clippedX = pamP->width - 2;
else
clippedX = roundedX;
if (roundedY <= 0)
clippedY = 1;
else if (roundedY > pamP->width - 1)
clippedY = pamP->width - 2;
else
clippedY = roundedY;
return pamd_makePoint(clippedX, clippedY);
}
static void drawClippedTriangle(const struct pam * const pamP,
tuple ** const tuples,
triangle const tri) {
pamd_point const p1 = clippedPoint(pamP, tri.p1);
pamd_point const p2 = clippedPoint(pamP, tri.p2);
pamd_point const p3 = clippedPoint(pamP, tri.p3);
pamd_line(tuples, pamP->width, pamP->height, pamP->depth, pamP->maxval,
p1, p2, frameDrawproc, black);
pamd_line(tuples, pamP->width, pamP->height, pamP->depth, pamP->maxval,
p2, p3, frameDrawproc, black);
pamd_line(tuples, pamP->width, pamP->height, pamP->depth, pamP->maxval,
p3, p1, frameDrawproc, black);
}
static void
prepTrig(int const wd,
int const ht) {
/* create triangles using control points */
point rtl1, rtr1, rbl1, rbr1;
point rtl2, rtr2, rbl2, rbr2;
point c1p1, c1p2, c1p3, c1p4;
point c2p1, c2p2, c2p3, c2p4;
line l1, l2;
point p0;
rtl1 = makepoint(0.0 + tiny(), 0.0 + tiny());
rtr1 = makepoint((double) wd - 1.0 + tiny(), 0.0 + tiny());
rbl1 = makepoint(0.0 + tiny(), (double) ht - 1.0 + tiny());
rbr1 = makepoint((double) wd - 1.0 + tiny(), (double) ht - 1.0 + tiny());
rtl2 = makepoint(0.0 + tiny(), 0.0 + tiny());
rtr2 = makepoint((double) wd - 1.0 + tiny(), 0.0 + tiny());
rbl2 = makepoint(0.0 + tiny(), (double) ht - 1.0 + tiny());
rbr2 = makepoint((double) wd - 1.0 + tiny(), (double) ht - 1.0 + tiny());
if (nCP == 1) {
c1p1 = oldCP[0];
c2p1 = newCP[0];
/* connect control point to all corners to get 4 triangles */
/* left side triangle */
sideTriangle(nCP,
&tri1s[0], c1p1, p0, p0, p0, rbl1, rtl1,
&tri2s[0], c2p1, p0, p0, p0, rbl2, rtl2);
/* top side triangle */
sideTriangle(nCP,
&tri1s[1], c1p1, p0, p0, p0, rtl1, rtr1,
&tri2s[1], c2p1, p0, p0, p0, rtl2, rtr2);
/* right side triangle */
sideTriangle(nCP,
&tri1s[2], c1p1, p0, p0, p0, rtr1, rbr1,
&tri2s[2], c2p1, p0, p0, p0, rtr2, rbr2);
/* bottom side triangle */
sideTriangle(nCP,
&tri1s[3], c1p1, p0, p0, p0, rbr1, rbl1,
&tri2s[3], c2p1, p0, p0, p0, rbr2, rbl2);
nTri = 4;
} else if (nCP == 2) {
c1p1 = oldCP[0];
c1p2 = oldCP[1];
c2p1 = newCP[0];
c2p2 = newCP[1];
/* check for hor/ver edges */
angle (&c1p1, &c1p2);
angle (&c2p1, &c2p2);
/* connect two control points to corners to get 6 triangles */
/* left side */
sideTriangle(nCP,
&tri1s[0], c1p1, c1p2, p0, p0, rbl1, rtl1,
&tri2s[0], c2p1, c2p2, p0, p0, rbl2, rtl2);
/* top side */
sideTriangle(nCP,
&tri1s[1], c1p1, c1p2, p0, p0, rtl1, rtr1,
&tri2s[1], c2p1, c2p2, p0, p0, rtl2, rtr2);
/* right side */
sideTriangle(nCP,
&tri1s[2], c1p1, c1p2, p0, p0, rtr1, rbr1,
&tri2s[2], c2p1, c2p2, p0, p0, rtr2, rbr2);
/* bottom side */
sideTriangle(nCP,
&tri1s[3], c1p1, c1p2, p0, p0, rbr1, rbl1,
&tri2s[3], c2p1, c2p2, p0, p0, rbr2, rbl2);
/* edge to corner triangles */
edgeTriangle(&tri1s[4], c1p1, c1p2, rtl1, rtr1, rbl1, rbr1,
&tri2s[4], c2p1, c2p2, rtl2, rtr2, rbl2, rbr2);
edgeTriangle(&tri1s[5], c1p2, c1p1, rtl1, rtr1, rbl1, rbr1,
&tri2s[5], c2p2, c2p1, rtl2, rtr2, rbl2, rbr2);
nTri = 6;
} else if (nCP == 3) {
c1p1 = oldCP[0];
c1p2 = oldCP[1];
c1p3 = oldCP[2];
c2p1 = newCP[0];
c2p2 = newCP[1];
c2p3 = newCP[2];
/* Move vertices slightly if necessary to make sure no edge is
horizontal or vertical.
*/
angle(&c1p1, &c1p2);
angle(&c1p2, &c1p3);
angle(&c1p3, &c1p1);
angle(&c2p1, &c2p2);
angle(&c2p2, &c2p3);
angle(&c2p3, &c2p1);
if (windtriangle(&tri1s[0], c1p1, c1p2, c1p3)) {
tri2s[0] = maketriangle(c2p1, c2p2, c2p3);
} else {
tri2s[0] = maketriangle(c2p1, c2p3, c2p2);
}
c1p1 = tri1s[0].p1;
c1p2 = tri1s[0].p2;
c1p3 = tri1s[0].p3;
c2p1 = tri2s[0].p1;
c2p2 = tri2s[0].p2;
c2p3 = tri2s[0].p3;
/* point to side triangles */
/* left side */
sideTriangle(nCP,
&tri1s[1], c1p1, c1p2, c1p3, p0, rbl1, rtl1,
&tri2s[1], c2p1, c2p2, c2p3, p0, rbl2, rtl2);
/* top side */
sideTriangle(nCP,
&tri1s[2], c1p1, c1p2, c1p3, p0, rtl1, rtr1,
&tri2s[2], c2p1, c2p2, c2p3, p0, rtl2, rtr2);
/* right side */
sideTriangle(nCP,
&tri1s[3], c1p1, c1p2, c1p3, p0, rtr1, rbr1,
&tri2s[3], c2p1, c2p2, c2p3, p0, rtr2, rbr2);
/* bottom side */
sideTriangle(nCP,
&tri1s[4], c1p1, c1p2, c1p3, p0, rbr1, rbl1,
&tri2s[4], c2p1, c2p2, c2p3, p0, rbr2, rbl2);
/* edge to corner triangles */
edgeTriangle(&tri1s[5], c1p1, c1p2, rtl1, rtr1, rbl1, rbr1,
&tri2s[5], c2p1, c2p2, rtl2, rtr2, rbl2, rbr2);
edgeTriangle(&tri1s[6], c1p2, c1p3, rtl1, rtr1, rbl1, rbr1,
&tri2s[6], c2p2, c2p3, rtl2, rtr2, rbl2, rbr2);
edgeTriangle(&tri1s[7], c1p3, c1p1, rtl1, rtr1, rbl1, rbr1,
&tri2s[7], c2p3, c2p1, rtl2, rtr2, rbl2, rbr2);
nTri = 8;
} else if (nCP == 4) {
c1p1 = oldCP[0];
c1p2 = oldCP[1];
c1p3 = oldCP[2];
c1p4 = oldCP[3];
c2p1 = newCP[0];
c2p2 = newCP[1];
c2p3 = newCP[2];
c2p4 = newCP[3];
/* check for hor/ver edges */
angle (&c1p1, &c1p2);
angle (&c1p2, &c1p3);
angle (&c1p3, &c1p4);
angle (&c1p4, &c1p1);
angle (&c1p1, &c1p3);
angle (&c1p2, &c1p4);
angle (&c2p1, &c2p2);
angle (&c2p2, &c2p3);
angle (&c2p3, &c2p4);
angle (&c2p4, &c2p1);
angle (&c2p1, &c2p3);
angle (&c2p2, &c2p4);
/*-------------------------------------------------------------------*/
/* -1- -2- -3- -4- -5- -6- */
/* 1 2 1 3 1 2 1 4 1 3 1 4 */
/* X X X X X X */
/* 3 4 2 4 4 3 2 3 4 2 3 2 */
/*-------------------------------------------------------------------*/
/* center two triangles */
l1 = makeline(c1p1, c1p4);
l2 = makeline(c1p2, c1p3);
if (intersect(&l1, &l2, &p0)) {
if (windtriangle(&tri1s[0], c1p1, c1p2, c1p3)) {
tri1s[1] = maketriangle(c1p4, c1p3, c1p2);
tri2s[0] = maketriangle(c2p1, c2p2, c2p3);
tri2s[1] = maketriangle(c2p4, c2p3, c2p2);
} else {
tri1s[1] = maketriangle(c1p4, c1p2, c1p3);
tri2s[0] = maketriangle(c2p1, c2p3, c2p2);
tri2s[1] = maketriangle(c2p4, c2p2, c2p3);
}
}
l1 = makeline(c1p1, c1p3);
l2 = makeline(c1p2, c1p4);
if (intersect(&l1, &l2, &p0)) {
if (windtriangle(&tri1s[0], c1p1, c1p2, c1p4)) {
tri1s[1] = maketriangle(c1p3, c1p4, c1p2);
tri2s[0] = maketriangle(c2p1, c2p2, c2p4);
tri2s[1] = maketriangle(c2p3, c2p4, c2p2);
} else {
tri1s[1] = maketriangle(c1p3, c1p2, c1p4);
tri2s[0] = maketriangle(c2p1, c2p4, c2p2);
tri2s[1] = maketriangle(c2p3, c2p2, c2p4);
}
}
l1 = makeline(c1p1, c1p2);
l2 = makeline(c1p3, c1p4);
if (intersect(&l1, &l2, &p0)) {
if (windtriangle(&tri1s[0], c1p1, c1p3, c1p4)) {
tri1s[1] = maketriangle(c1p2, c1p4, c1p3);
tri2s[0] = maketriangle(c2p1, c2p3, c2p4);
tri2s[1] = maketriangle(c2p2, c2p4, c2p3);
} else {
tri1s[1] = maketriangle(c1p2, c1p3, c1p4);
tri2s[0] = maketriangle(c2p1, c2p4, c2p3);
tri2s[1] = maketriangle(c2p2, c2p3, c2p4);
}
}
/* control points in correct orientation */
c1p1 = tri1s[0].p1;
c1p2 = tri1s[0].p2;
c1p3 = tri1s[0].p3;
c1p4 = tri1s[1].p1;
c2p1 = tri2s[0].p1;
c2p2 = tri2s[0].p2;
c2p3 = tri2s[0].p3;
c2p4 = tri2s[1].p1;
/* triangle from triangle point to side of image */
/* left side triangle */
sideTriangle(nCP,
&tri1s[2], c1p1, c1p2, c1p3, c1p4, rbl1, rtl1,
&tri2s[2], c2p1, c2p2, c2p3, c2p4, rbl2, rtl2);
/* top side triangle */
sideTriangle(nCP,
&tri1s[3], c1p1, c1p2, c1p3, c1p4, rtl1, rtr1,
&tri2s[3], c2p1, c2p2, c2p3, c2p4, rtl2, rtr2);
/* right side triangle */
sideTriangle(nCP,
&tri1s[4], c1p1, c1p2, c1p3, c1p4, rtr1, rbr1,
&tri2s[4], c2p1, c2p2, c2p3, c2p4, rtr2, rbr2);
/* bottom side triangle */
sideTriangle(nCP,
&tri1s[5], c1p1, c1p2, c1p3, c1p4, rbr1, rbl1,
&tri2s[5], c2p1, c2p2, c2p3, c2p4, rbr2, rbl2);
/*-------------------------------------------------------------------*/
/* -1- -2- -3- -4- -5- -6- */
/* 1 2 1 3 1 2 1 4 1 3 1 4 */
/* X X X X X X */
/* 3 4 2 4 4 3 2 3 4 2 3 2 */
/*-------------------------------------------------------------------*/
/* edge-corner triangles */
edgeTriangle(&tri1s[6], c1p1, c1p2, rtl1, rtr1, rbl1, rbr1,
&tri2s[6], c2p1, c2p2, rtl2, rtr2, rbl2, rbr2);
edgeTriangle(&tri1s[7], c1p2, c1p4, rtl1, rtr1, rbl1, rbr1,
&tri2s[7], c2p2, c2p4, rtl2, rtr2, rbl2, rbr2);
edgeTriangle(&tri1s[8], c1p4, c1p3, rtl1, rtr1, rbl1, rbr1,
&tri2s[8], c2p4, c2p3, rtl2, rtr2, rbl2, rbr2);
edgeTriangle(&tri1s[9], c1p3, c1p1, rtl1, rtr1, rbl1, rbr1,
&tri2s[9], c2p3, c2p1, rtl2, rtr2, rbl2, rbr2);
nTri = 10;
}
}
static void
prepQuad(void) {
/* order quad control points */
double d01, d12, d20;
line l1, l2;
point mid;
triangle tri;
if (nCP == 1) {
/* create a rectangle from top-left corner of image and control
point
*/
quad1 = quadRect(0.0, oldCP[0].x, 0.0, oldCP[0].y);
quad2 = quadRect(0.0, newCP[0].x, 0.0, newCP[0].y);
} else if (nCP == 2) {
/* create a rectangle with the two points as opposite corners */
if ((oldCP[0].x < oldCP[1].x) && (oldCP[0].y < oldCP[1].y)) {
/* top-left and bottom-right */
quad1 = quadRect(oldCP[0].x,oldCP[1].x, oldCP[0].y, oldCP[1].y);
} else if ((oldCP[0].x > oldCP[1].x) && (oldCP[0].y < oldCP[1].y)) {
/* top-right and bottom-left */
quad1 = quadRect(oldCP[1].x, oldCP[0].x, oldCP[0].y, oldCP[1].y);
} else if ((oldCP[0].x < oldCP[1].x) && (oldCP[0].y < oldCP[1].y)) {
/* bottom-left and top-right */
quad1 = quadRect(oldCP[0].x, oldCP[1].x, oldCP[1].y, oldCP[0].y);
} else if ((oldCP[0].x > oldCP[1].x) && (oldCP[0].y < oldCP[1].y)) {
/* bottom-right and top-left */
quad1 = quadRect(oldCP[1].x, oldCP[0].x, oldCP[1].y, oldCP[0].y);
}
if ((newCP[0].x < newCP[1].x) && (newCP[0].y < newCP[1].y)) {
/* top-left and bottom-right */
quad2 = quadRect(newCP[0].x, newCP[1].x, newCP[0].y, newCP[1].y);
} else if ((newCP[0].x > newCP[1].x) && (newCP[0].y < newCP[1].y)) {
/* top-right and bottom-left */
quad2 = quadRect(newCP[1].x, newCP[0].x, newCP[0].y, newCP[1].y);
} else if ((newCP[0].x < newCP[1].x) && (newCP[0].y < newCP[1].y)) {
/* bottom-left and top-right */
quad2 = quadRect(newCP[0].x, newCP[1].x, newCP[1].y, newCP[0].y);
} else if ((newCP[0].x > newCP[1].x) && (newCP[0].y < newCP[1].y)) {
/* bottom-right and top-left */
quad2 = quadRect(newCP[1].x, newCP[0].x, newCP[1].y, newCP[0].y);
}
} else {
if (nCP == 3) {
/* add the fourth corner of a parallelogram */
/* diagonal of the parallelogram is the two control points
furthest apart
*/
d01 = distance(newCP[0], newCP[1]);
d12 = distance(newCP[1], newCP[2]);
d20 = distance(newCP[2], newCP[0]);
if ((d01 > d12) && (d01 > d20)) {
oldCP[3].x = oldCP[0].x + oldCP[1].x - oldCP[2].x;
oldCP[3].y = oldCP[0].y + oldCP[1].y - oldCP[2].y;
} else
if (d12 > d20) {
oldCP[3].x = oldCP[1].x + oldCP[2].x - oldCP[0].x;
oldCP[3].y = oldCP[1].y + oldCP[2].y - oldCP[0].y;
} else {
oldCP[3].x = oldCP[2].x + oldCP[0].x - oldCP[1].x;
oldCP[3].y = oldCP[2].y + oldCP[0].y - oldCP[1].y;
}
if ((d01 > d12) && (d01 > d20)) {
newCP[3].x = newCP[0].x + newCP[1].x - newCP[2].x;
newCP[3].y = newCP[0].y + newCP[1].y - newCP[2].y;
} else
if (d12 > d20) {
newCP[3].x = newCP[1].x + newCP[2].x - newCP[0].x;
newCP[3].y = newCP[1].y + newCP[2].y - newCP[0].y;
} else {
newCP[3].x = newCP[2].x + newCP[0].x - newCP[1].x;
newCP[3].y = newCP[2].y + newCP[0].y - newCP[1].y;
}
} /* end nCP = 3 */
/* nCP = 3 or 4 */
/* find the intersection of the diagonals */
l1 = makeline(oldCP[0], oldCP[1]);
l2 = makeline(oldCP[2], oldCP[3]);
if (intersect(&l1, &l2, &mid)) {
quadCorner(oldCP[0], oldCP[1], oldCP[2], oldCP[3],
mid, &quad1, &tri);
} else {
l1 = makeline(oldCP[0], oldCP[2]);
l2 = makeline(oldCP[1], oldCP[3]);
if (intersect(&l1, &l2, &mid))
quadCorner(oldCP[0], oldCP[2], oldCP[1], oldCP[3],
mid, &quad1, &tri);
else {
l1 = makeline(oldCP[0], oldCP[3]);
l2 = makeline(oldCP[1], oldCP[2]);
if (intersect(&l1, &l2, &mid))
quadCorner(oldCP[0], oldCP[3],
oldCP[1], oldCP[2], mid, &quad1, &tri);
else
pm_error("The four old control points don't seem "
"to be corners.");
}
}
/* repeat for the "to-be" control points */
l1 = makeline(newCP[0], newCP[1]);
l2 = makeline(newCP[2], newCP[3]);
if (intersect(&l1, &l2, &mid))
quadCorner(newCP[0], newCP[1], newCP[2], newCP[3],
mid, &quad2, &tri);
else {
l1 = makeline(newCP[0], newCP[2]);
l2 = makeline(newCP[1], newCP[3]);
if (intersect(&l1, &l2, &mid))
quadCorner(newCP[0], newCP[2], newCP[1], newCP[3],
mid, &quad2, &tri);
else {
l1 = makeline(newCP[0], newCP[3]);
l2 = makeline(newCP[1], newCP[2]);
if (intersect(&l1, &l2, &mid))
quadCorner(newCP[0], newCP[3],
newCP[1], newCP[2], mid, &quad2, &tri);
else
pm_error("The four new control points don't seem "
"to be corners.");
}
}
}
}
static void
warpTrig(point const p2,
point * const p1P) {
/* map target to source by triangulation */
point e1p1, e1p2, e1p3;
point e2p1, e2p2, e2p3;
line lp, le;
line l1, l2, l3;
double d1, d2, d3;
int i;
/* find in which triangle p2 lies */
for (i = 0; i < nTri; i++) {
if (insidetri (&tri2s[i], p2))
break;
}
if (i == nTri)
*p1P = makepoint(0.0, 0.0);
else {
/* where in triangle is point */
d1 = fabs (p2.x - tri2s[i].p1.x) + fabs (p2.y - tri2s[i].p1.y);
d2 = fabs (p2.x - tri2s[i].p2.x) + fabs (p2.y - tri2s[i].p2.y);
d3 = fabs (p2.x - tri2s[i].p3.x) + fabs (p2.y - tri2s[i].p3.y);
/* line through p1 and p intersecting with edge p2-p3 */
lp = makeline(tri2s[i].p1, p2);
le = makeline(tri2s[i].p2, tri2s[i].p3);
intersect (&lp, &le, &e2p1);
/* line through p2 and p intersecting with edge p3-p1 */
lp = makeline(tri2s[i].p2, p2);
le = makeline(tri2s[i].p3, tri2s[i].p1);
intersect (&lp, &le, &e2p2);
/* line through p3 and p intersecting with edge p1-p2 */
lp = makeline(tri2s[i].p3, p2);
le = makeline(tri2s[i].p1, tri2s[i].p2);
intersect (&lp, &le, &e2p3);
/* map target control points to source control points */
e1p1.x = tri1s[i].p2.x
+ (e2p1.x - tri2s[i].p2.x)/(tri2s[i].p3.x - tri2s[i].p2.x)
* (tri1s[i].p3.x - tri1s[i].p2.x);
e1p1.y = tri1s[i].p2.y
+ (e2p1.y - tri2s[i].p2.y)/(tri2s[i].p3.y - tri2s[i].p2.y)
* (tri1s[i].p3.y - tri1s[i].p2.y);
e1p2.x = tri1s[i].p3.x
+ (e2p2.x - tri2s[i].p3.x)/(tri2s[i].p1.x - tri2s[i].p3.x)
* (tri1s[i].p1.x - tri1s[i].p3.x);
e1p2.y = tri1s[i].p3.y
+ (e2p2.y - tri2s[i].p3.y)/(tri2s[i].p1.y - tri2s[i].p3.y)
* (tri1s[i].p1.y - tri1s[i].p3.y);
e1p3.x = tri1s[i].p1.x
+ (e2p3.x - tri2s[i].p1.x)/(tri2s[i].p2.x - tri2s[i].p1.x)
* (tri1s[i].p2.x - tri1s[i].p1.x);
e1p3.y = tri1s[i].p1.y
+ (e2p3.y - tri2s[i].p1.y)/(tri2s[i].p2.y - tri2s[i].p1.y)
* (tri1s[i].p2.y - tri1s[i].p1.y);
/* intersect grid lines in source */
l1 = makeline(tri1s[i].p1, e1p1);
l2 = makeline(tri1s[i].p2, e1p2);
l3 = makeline(tri1s[i].p3, e1p3);
if ((d1 < d2) && (d1 < d3))
intersect (&l2, &l3, p1P);
else if (d2 < d3)
intersect (&l1, &l3, p1P);
else
intersect (&l1, &l2, p1P);
}
}
static void
warpQuad(point const p2,
point * const p1P) {
/* map target to source for quad control points */
point h2, v2;
point c1tl, c1tr, c1bl, c1br;
point c2tl, c2tr, c2bl, c2br;
point e1t, e1b, e1l, e1r;
point e2t, e2b, e2l, e2r;
line l2t, l2b, l2l, l2r;
line lh, lv;
c1tl = quad1.tl;
c1tr = quad1.tr;
c1bl = quad1.bl;
c1br = quad1.br;
c2tl = quad2.tl;
c2tr = quad2.tr;
c2bl = quad2.bl;
c2br = quad2.br;
l2t = makeline(c2tl, c2tr);
l2b = makeline(c2bl, c2br);
l2l = makeline(c2tl, c2bl);
l2r = makeline(c2tr, c2br);
/* find intersections of lines through control points */
intersect (&l2t, &l2b, &h2);
intersect (&l2l, &l2r, &v2);
/* find intersections of axes through P with control point box */
lv = makeline(p2, v2);
intersect (&l2t, &lv, &e2t);
intersect (&l2b, &lv, &e2b);
lh = makeline(p2, h2);
intersect (&l2l, &lh, &e2l);
intersect (&l2r, &lh, &e2r);
/* map target control points to source control points */
e1t.x = c1tl.x + (e2t.x - c2tl.x)/(c2tr.x - c2tl.x) * (c1tr.x - c1tl.x);
if (c1tl.y == c1tr.y)
e1t.y = c1tl.y;
else
e1t.y =
c1tl.y + (e2t.x - c2tl.x)/(c2tr.x - c2tl.x) * (c1tr.y - c1tl.y);
e1b.x = c1bl.x + (e2b.x - c2bl.x)/(c2br.x - c2bl.x) * (c1br.x - c1bl.x);
if (c1bl.y == c1br.y)
e1b.y = c1bl.y;
else
e1b.y =
c1bl.y + (e2b.x - c2bl.x)/(c2br.x - c2bl.x) * (c1br.y - c1bl.y);
if (c1tl.x == c1bl.x)
e1l.x = c1tl.x;
else
e1l.x =
c1tl.x + (e2l.y - c2tl.y)/(c2bl.y - c2tl.y) * (c1bl.x - c1tl.x);
e1l.y = c1tl.y + (e2l.y - c2tl.y)/(c2bl.y - c2tl.y) * (c1bl.y - c1tl.y);
if (c1tr.x == c1br.x)
e1r.x = c1tr.x;
else
e1r.x
= c1tr.x + (e2r.y - c2tr.y)/(c2br.y - c2tr.y) * (c1br.x - c1tr.x);
e1r.y = c1tr.y + (e2r.y - c2tr.y)/(c2br.y - c2tr.y) * (c1br.y - c1tr.y);
/* intersect grid lines in source */
lv = makeline(e1t, e1b);
lh = makeline(e1l, e1r);
intersect (&lh, &lv, p1P);
}
static void
setGlobalCP(struct cmdlineInfo const cmdline) {
unsigned int i;
for (i = 0; i < cmdline.nCP; ++i) {
oldCP[i] = cmdline.oldCP[i];
newCP[i] = cmdline.newCP[i];
}
nCP = cmdline.nCP;
for (i = nCP; i < ARRAY_SIZE(oldCP); ++i)
oldCP[i] = makepoint(-1.0, -1.0);
for (i = nCP; i < ARRAY_SIZE(newCP); ++i)
newCP[i] = makepoint(-1.0, -1.0);
}
static void
createWhiteTuple(const struct pam * const pamP,
tuple * const tupleP) {
tuple white;
unsigned int plane;
white = pnm_allocpamtuple(pamP);
for (plane = 0; plane < pamP->depth; ++plane)
white[plane] = pamP->maxval;
*tupleP = white;
}
static void
makeAllWhite(struct pam * const pamP,
tuple ** const tuples) {
tuple white;
unsigned int row;
createWhiteTuple(pamP, &white);
for (row = 0; row < pamP->height; ++row) {
unsigned int col;
for (col = 0; col < pamP->width; ++col)
pnm_assigntuple(pamP, tuples[row][col], white);
}
pnm_freepamtuple(white);
}
static sample
pix(tuple ** const tuples,
unsigned int const width,
unsigned int const height,
point const p1,
unsigned int const plane,
bool const linear) {
point p;
double pix;
p.x = p1.x + 1E-3;
p.y = p1.y + 1E-3;
if (p.x < 0.0)
p.x = 0.0;
if (p.x > (double) width - 1.0)
p.x = (double) width - 1.0;
if (p.y < 0.0)
p.y = 0.0;
if (p.y > (double) height - 1.0)
p.y = (double) height - 1.0;
if (!linear) {
pix = tuples
[(int) floor(p.y + 0.5)]
[(int) floor(p.x + 0.5)][plane];
} else {
double const rx = p.x - floor(p.x);
double const ry = p.y - floor(p.y);
pix = 0.0;
if (((int) floor(p.x) >= 0) && ((int) floor(p.x) < width) &&
((int) floor(p.y) >= 0) && ((int) floor(p.y) < height)) {
pix += (1.0 - rx) * (1.0 - ry)
* tuples[(int) floor(p.y)][(int) floor(p.x)][plane];
}
if (((int) floor(p.x) + 1 >= 0) && ((int) floor(p.x) + 1 < width) &&
((int) floor(p.y) >= 0) && ((int) floor(p.y) < height)) {
pix += rx * (1.0 - ry)
* tuples[(int) floor(p.y)][(int) floor(p.x) + 1][plane];
}
if (((int) floor(p.x) >= 0) && ((int) floor(p.x) < width) &&
((int) floor(p.y) + 1 >= 0) && ((int) floor(p.y) + 1 < height)) {
pix += (1.0 - rx) * ry
* tuples[(int) floor(p.y) + 1][(int) floor(p.x)][plane];
}
if (((int) floor(p.x) + 1 >= 0) && ((int) floor(p.x) + 1 < width) &&
((int) floor(p.y) + 1 >= 0) && ((int) floor(p.y) + 1 < height)) {
pix += rx * ry
* tuples[(int) floor(p.y) + 1][(int) floor(p.x) + 1][plane];
}
}
return (int) floor(pix);
}
int
main(int argc, const char ** const argv) {
struct cmdlineInfo cmdline;
FILE * ifP;
struct pam inpam, outpam;
tuple ** inTuples;
tuple ** outTuples;
unsigned int p2y;
pm_proginit(&argc, argv);
parseCmdline(argc, argv, &cmdline);
setGlobalCP(cmdline);
srand(cmdline.randseedSpec ? cmdline.randseed : pm_randseed());
ifP = pm_openr(cmdline.fileName);
inTuples = pnm_readpam(ifP, &inpam, PAM_STRUCT_SIZE(tuple_type));
outpam = inpam; /* initial value */
outpam.file = stdout;
outTuples = pnm_allocpamarray(&outpam);
pnm_createBlackTuple(&outpam, &black);
makeAllWhite(&outpam, outTuples);
if (cmdline.tri)
prepTrig(inpam.width, inpam.height);
if (cmdline.quad)
prepQuad();
for (p2y = 0; p2y < inpam.height; ++p2y) {
unsigned int p2x;
for (p2x = 0; p2x < inpam.width; ++p2x) {
point p1, p2;
unsigned int plane;
p2 = makepoint(p2x, p2y);
if (cmdline.quad)
warpQuad(p2, &p1);
if (cmdline.tri)
warpTrig(p2, &p1);
for (plane = 0; plane < inpam.depth; ++plane) {
outTuples[p2y][p2x][plane] =
pix(inTuples, inpam.width, inpam.height, p1, plane,
cmdline.linear);
}
}
}
if (cmdline.frame) {
if (cmdline.quad) {
drawExtendedLine(&outpam, outTuples, quad2.tl, quad2.tr);
drawExtendedLine(&outpam, outTuples, quad2.bl, quad2.br);
drawExtendedLine(&outpam, outTuples, quad2.tl, quad2.bl);
drawExtendedLine(&outpam, outTuples, quad2.tr, quad2.br);
}
if (cmdline.tri) {
unsigned int i;
for (i = 0; i < nTri; ++i)
drawClippedTriangle(&outpam, outTuples, tri2s[i]);
}
}
pnm_writepam(&outpam, outTuples);
pnm_freepamtuple(black);
pnm_freepamarray(outTuples, &outpam);
pnm_freepamarray(inTuples, &inpam);
pm_close(ifP);
pm_close(stdout);
return 0;
}