src/px/pmapgen.c

#ifndef lint
static const char       RCSid[] = "$Id: pmapgen.c,v 2.3 2003/02/22 02:07:27 greg Exp $";
#endif
/*
 * pmapgen.c: general routines for 2-D perspective mappings.
 * These routines are independent of the poly structure,
 * so we do not think in terms of texture and screen space.
 *
 * Paul Heckbert        5 Nov 85, 12 Dec 85
 */

#include <stdio.h>
#include "pmap.h"
#include "mx3.h"

#define TOLERANCE 1e-13
#define ZERO(x) ((x)<TOLERANCE && (x)>-TOLERANCE)

#define X(i) qdrl[i][0]         /* quadrilateral x and y */
#define Y(i) qdrl[i][1]

/*
 * pmap_quad_rect: find mapping between quadrilateral and rectangle.
 * The correspondence is:
 *
 *      qdrl[0] --> (u0,v0)
 *      qdrl[1] --> (u1,v0)
 *      qdrl[2] --> (u1,v1)
 *      qdrl[3] --> (u0,v1)
 *
 * This method of computing the adjoint numerically is cheaper than
 * computing it symbolically.
 */
        
extern int
pmap_quad_rect(
        double u0,              /* bounds of rectangle */
        double v0,
        double u1,
        double v1,
        double qdrl[4][2],              /* vertices of quadrilateral */
        double QR[3][3]         /* qdrl->rect transform (returned) */
)
{
    int ret;
    double du, dv;
    double RQ[3][3];            /* rect->qdrl transform */

    du = u1-u0;
    dv = v1-v0;
    if (du==0. || dv==0.) {
        fprintf(stderr, "pmap_quad_rect: null rectangle\n");
        return PMAP_BAD;
    }

    /* first find mapping from unit uv square to xy quadrilateral */
    ret = pmap_square_quad(qdrl, RQ);
    if (ret==PMAP_BAD) return PMAP_BAD;

    /* concatenate transform from uv rectangle (u0,v0,u1,v1) to unit square */
    RQ[0][0] /= du;
    RQ[1][0] /= dv;
    RQ[2][0] -= RQ[0][0]*u0 + RQ[1][0]*v0;
    RQ[0][1] /= du;
    RQ[1][1] /= dv;
    RQ[2][1] -= RQ[0][1]*u0 + RQ[1][1]*v0;
    RQ[0][2] /= du;
    RQ[1][2] /= dv;
    RQ[2][2] -= RQ[0][2]*u0 + RQ[1][2]*v0;

    /* now RQ is transform from uv rectangle to xy quadrilateral */
    /* QR = inverse transform, which maps xy to uv */
    if (mx3d_adjoint(RQ, QR)==0.)
        fprintf(stderr, "pmap_quad_rect: warning: determinant=0\n");
    return ret;
}

/*
 * pmap_square_quad: find mapping between unit square and quadrilateral.
 * The correspondence is:
 *
 *      (0,0) --> qdrl[0]
 *      (1,0) --> qdrl[1]
 *      (1,1) --> qdrl[2]
 *      (0,1) --> qdrl[3]
 */

extern int
pmap_square_quad(
        register double qdrl[4][2],     /* vertices of quadrilateral */
        register double SQ[3][3]        /* square->qdrl transform */
)
{
    double px, py;

    px = X(0)-X(1)+X(2)-X(3);
    py = Y(0)-Y(1)+Y(2)-Y(3);

    if (ZERO(px) && ZERO(py)) {         /* affine */
        SQ[0][0] = X(1)-X(0);
        SQ[1][0] = X(2)-X(1);
        SQ[2][0] = X(0);
        SQ[0][1] = Y(1)-Y(0);
        SQ[1][1] = Y(2)-Y(1);
        SQ[2][1] = Y(0);
        SQ[0][2] = 0.;
        SQ[1][2] = 0.;
        SQ[2][2] = 1.;
        return PMAP_LINEAR;
    }
    else {                              /* perspective */
        double dx1, dx2, dy1, dy2, del;

        dx1 = X(1)-X(2);
        dx2 = X(3)-X(2);
        dy1 = Y(1)-Y(2);
        dy2 = Y(3)-Y(2);
        del = DET2(dx1,dx2, dy1,dy2);
        if (del==0.) {
            fprintf(stderr, "pmap_square_quad: bad mapping\n");
            return PMAP_BAD;
        }
        SQ[0][2] = DET2(px,dx2, py,dy2)/del;
        SQ[1][2] = DET2(dx1,px, dy1,py)/del;
        SQ[2][2] = 1.;
        SQ[0][0] = X(1)-X(0)+SQ[0][2]*X(1);
        SQ[1][0] = X(3)-X(0)+SQ[1][2]*X(3);
        SQ[2][0] = X(0);
        SQ[0][1] = Y(1)-Y(0)+SQ[0][2]*Y(1);
        SQ[1][1] = Y(3)-Y(0)+SQ[1][2]*Y(3);
        SQ[2][1] = Y(0);
        return PMAP_PERSP;
    }
}
Revision:	2.4
Committed:	Sun Mar 28 20:33:14 2004 UTC (20 years, 1 month ago) by schorsch
Content type:	text/plain
Branch:	MAIN
CVS Tags:	rad5R4, rad5R2, rad4R2P2, rad5R0, rad5R1, rad3R7P2, rad3R7P1, rad4R2, rad4R1, rad4R0, rad3R6, rad3R6P1, rad3R8, rad3R9, rad4R2P1, rad5R3, HEAD
Changes since 2.3:	+15 -9 lines
Log Message:	Continued ANSIfication, and other fixes and clarifications.
#	User	Rev	Content
1	greg	2.1	#ifndef lint
2	schorsch	2.4	static const char RCSid[] = "$Id: pmapgen.c,v 2.3 2003/02/22 02:07:27 greg Exp $";
3	greg	2.1	#endif
4			/*
5			* pmapgen.c: general routines for 2-D perspective mappings.
6			* These routines are independent of the poly structure,
7			* so we do not think in terms of texture and screen space.
8			*
9			* Paul Heckbert 5 Nov 85, 12 Dec 85
10			*/
11
12			#include <stdio.h>
13			#include "pmap.h"
14			#include "mx3.h"
15
16			#define TOLERANCE 1e-13
17			#define ZERO(x) ((x)<TOLERANCE && (x)>-TOLERANCE)
18
19	greg	2.2	#define X(i) qdrl[i][0] /* quadrilateral x and y */
20			#define Y(i) qdrl[i][1]
21	greg	2.1
22			/*
23			* pmap_quad_rect: find mapping between quadrilateral and rectangle.
24			* The correspondence is:
25			*
26	greg	2.2	* qdrl[0] --> (u0,v0)
27			* qdrl[1] --> (u1,v0)
28			* qdrl[2] --> (u1,v1)
29			* qdrl[3] --> (u0,v1)
30	greg	2.1	*
31			* This method of computing the adjoint numerically is cheaper than
32			* computing it symbolically.
33			*/
34
35	schorsch	2.4	extern int
36			pmap_quad_rect(
37			double u0, /* bounds of rectangle */
38			double v0,
39			double u1,
40			double v1,
41			double qdrl[4][2], /* vertices of quadrilateral */
42			double QR[3][3] /* qdrl->rect transform (returned) */
43			)
44	greg	2.1	{
45			int ret;
46			double du, dv;
47	greg	2.2	double RQ[3][3]; /* rect->qdrl transform */
48	greg	2.1
49			du = u1-u0;
50			dv = v1-v0;
51			if (du==0. \|\| dv==0.) {
52			fprintf(stderr, "pmap_quad_rect: null rectangle\n");
53			return PMAP_BAD;
54			}
55
56			/* first find mapping from unit uv square to xy quadrilateral */
57	greg	2.2	ret = pmap_square_quad(qdrl, RQ);
58	greg	2.1	if (ret==PMAP_BAD) return PMAP_BAD;
59
60			/* concatenate transform from uv rectangle (u0,v0,u1,v1) to unit square */
61			RQ[0][0] /= du;
62			RQ[1][0] /= dv;
63			RQ[2][0] -= RQ[0][0]u0 + RQ[1][0]v0;
64			RQ[0][1] /= du;
65			RQ[1][1] /= dv;
66			RQ[2][1] -= RQ[0][1]u0 + RQ[1][1]v0;
67			RQ[0][2] /= du;
68			RQ[1][2] /= dv;
69			RQ[2][2] -= RQ[0][2]u0 + RQ[1][2]v0;
70
71			/* now RQ is transform from uv rectangle to xy quadrilateral */
72			/* QR = inverse transform, which maps xy to uv */
73			if (mx3d_adjoint(RQ, QR)==0.)
74			fprintf(stderr, "pmap_quad_rect: warning: determinant=0\n");
75			return ret;
76			}
77
78			/*
79			* pmap_square_quad: find mapping between unit square and quadrilateral.
80			* The correspondence is:
81			*
82	greg	2.2	* (0,0) --> qdrl[0]
83			* (1,0) --> qdrl[1]
84			* (1,1) --> qdrl[2]
85			* (0,1) --> qdrl[3]
86	greg	2.1	*/
87
88	schorsch	2.4	extern int
89			pmap_square_quad(
90			register double qdrl[4][2], /* vertices of quadrilateral */
91			register double SQ[3][3] /* square->qdrl transform */
92			)
93	greg	2.1	{
94			double px, py;
95
96			px = X(0)-X(1)+X(2)-X(3);
97			py = Y(0)-Y(1)+Y(2)-Y(3);
98
99			if (ZERO(px) && ZERO(py)) { /* affine */
100			SQ[0][0] = X(1)-X(0);
101			SQ[1][0] = X(2)-X(1);
102			SQ[2][0] = X(0);
103			SQ[0][1] = Y(1)-Y(0);
104			SQ[1][1] = Y(2)-Y(1);
105			SQ[2][1] = Y(0);
106			SQ[0][2] = 0.;
107			SQ[1][2] = 0.;
108			SQ[2][2] = 1.;
109			return PMAP_LINEAR;
110			}
111			else { /* perspective */
112			double dx1, dx2, dy1, dy2, del;
113
114			dx1 = X(1)-X(2);
115			dx2 = X(3)-X(2);
116			dy1 = Y(1)-Y(2);
117			dy2 = Y(3)-Y(2);
118			del = DET2(dx1,dx2, dy1,dy2);
119			if (del==0.) {
120			fprintf(stderr, "pmap_square_quad: bad mapping\n");
121			return PMAP_BAD;
122			}
123			SQ[0][2] = DET2(px,dx2, py,dy2)/del;
124			SQ[1][2] = DET2(dx1,px, dy1,py)/del;
125			SQ[2][2] = 1.;
126			SQ[0][0] = X(1)-X(0)+SQ[0][2]*X(1);
127			SQ[1][0] = X(3)-X(0)+SQ[1][2]*X(3);
128			SQ[2][0] = X(0);
129			SQ[0][1] = Y(1)-Y(0)+SQ[0][2]*Y(1);
130			SQ[1][1] = Y(3)-Y(0)+SQ[1][2]*Y(3);
131			SQ[2][1] = Y(0);
132			return PMAP_PERSP;
133			}
134			}