<<>> <> <> <> <<>> DIRECTORY G2dBasic, G2dMatrix, Imager, RealFns; G2dMatrixImpl: CEDAR MONITOR IMPORTS RealFns EXPORTS G2dMatrix ~ BEGIN Pair: TYPE ~ G2dBasic.Pair; Triple: TYPE ~ G2dBasic.Triple; Matrix: TYPE ~ G2dMatrix.Matrix; Rectangle: TYPE ~ Imager.Rectangle; singular: PUBLIC ERROR ~ CODE; <> Identity: PUBLIC PROC RETURNS [Matrix] ~ { RETURN[[[1.0, 0.0, 0.0], [0.0, 1.0, 0.0], [0.0, 0.0, 1.0]]]; }; Transpose: PUBLIC PROC [m: Matrix] RETURNS [t: Matrix] ~ { t Ź [[m.row1.x, m.row2.x, m.row3.x], [m.row1.y, m.row2.y, m.row3.y], [m.row1.z, m.row2.z, m.row3.z]]; }; Determinant: PUBLIC PROC [m: Matrix] RETURNS [d: REAL] ~ { d Ź m.row1.x*(m.row2.y*m.row3.z-m.row3.y*m.row2.z)- m.row1.y*(m.row2.x*m.row3.z-m.row3.x*m.row2.z)+ m.row1.z*(m.row2.x*m.row3.y-m.row3.x*m.row2.y); }; Adjoint: PUBLIC PROC [m: Matrix] RETURNS [r: Matrix] ~ { r.row1.x Ź m.row2.y*m.row3.z-m.row2.z*m.row3.y; r.row2.x Ź -(m.row2.x*m.row3.z-m.row2.z*m.row3.x); r.row3.x Ź m.row2.x*m.row3.y-m.row2.y*m.row3.x; r.row1.y Ź -(m.row1.y*m.row3.z-m.row1.z*m.row3.y); r.row2.y Ź m.row1.x*m.row3.z-m.row1.z*m.row3.x; r.row3.y Ź -(m.row1.x*m.row3.y-m.row1.y*m.row3.x); r.row1.z Ź m.row1.y*m.row2.z-m.row1.z*m.row2.y; r.row2.z Ź -(m.row1.x*m.row2.z-m.row1.z*m.row2.x); r.row3.z Ź m.row1.x*m.row2.y-m.row1.y*m.row2.x; }; Invert: PUBLIC PROC [m: Matrix] RETURNS [r: Matrix] ~ { det: REAL; r Ź Adjoint[m]; IF (det Ź r.row1.x*m.row1.x+r.row1.y*m.row2.x+r.row1.z*m.row3.x) = 0.0 THEN ERROR singular ELSE { di: REAL Ź 1.0/det; r.row1.x Ź di*r.row1.x; r.row1.y Ź di*r.row1.y; r.row1.z Ź di*r.row1.z; r.row2.x Ź di*r.row2.x; r.row2.y Ź di*r.row2.y; r.row2.z Ź di*r.row2.z; r.row3.x Ź di*r.row3.x; r.row3.y Ź di*r.row3.y; r.row3.z Ź di*r.row3.z; }; }; Mul: PUBLIC PROC [left, rite: Matrix] RETURNS [ret: Matrix] ~ { ret.row1.x Ź left.row1.x*rite.row1.x+left.row1.y*rite.row2.x+left.row1.z*rite.row3.x; ret.row1.y Ź left.row1.x*rite.row1.y+left.row1.y*rite.row2.y+left.row1.z*rite.row3.y; ret.row1.z Ź left.row1.x*rite.row1.z+left.row1.y*rite.row2.z+left.row1.z*rite.row3.z; ret.row2.x Ź left.row2.x*rite.row1.x+left.row2.y*rite.row2.x+left.row2.z*rite.row3.x; ret.row2.y Ź left.row2.x*rite.row1.y+left.row2.y*rite.row2.y+left.row2.z*rite.row3.y; ret.row2.z Ź left.row2.x*rite.row1.z+left.row2.y*rite.row2.z+left.row2.z*rite.row3.z; ret.row3.x Ź left.row3.x*rite.row1.x+left.row3.y*rite.row2.x+left.row3.z*rite.row3.x; ret.row3.y Ź left.row3.x*rite.row1.y+left.row3.y*rite.row2.y+left.row3.z*rite.row3.y; ret.row3.z Ź left.row3.x*rite.row1.z+left.row3.y*rite.row2.z+left.row3.z*rite.row3.z; }; Rotate: PUBLIC PROC [mat: Matrix, radians: REAL] RETURNS [Matrix] ~ { cos: REAL _ RealFns.Cos[radians]; sin: REAL _ RealFns.Sin[radians]; RETURN[Mul[mat, [[cos, sin, 0.0], [-sin, cos, 0.0], [0.0, 0.0, 1.0]]]]; }; Scale: PUBLIC PROC [mat: Matrix, scale: Pair] RETURNS [Matrix] ~ { RETURN[Mul[mat, [[scale.x, 0.0, 0.0], [0.0, scale.y, 0.0], [0.0, 0.0, 1.0]]]]; }; Translate: PUBLIC PROC [mat: Matrix, translate: Pair] RETURNS [Matrix] ~ { RETURN[Mul[mat, [[1.0, 0.0, 0.0], [0.0, 1.0, 0.0], [translate.x, translate.y, 1.0]]]]; }; Transform: PUBLIC PROC [p: Triple, mat: Matrix] RETURNS [Triple] ~ { RETURN[[ p.x*mat.row1.x+p.y*mat.row2.x+p.z*mat.row3.x, p.x*mat.row1.y+p.y*mat.row2.y+p.z*mat.row3.y, p.x*mat.row1.z+p.y*mat.row2.z+p.z*mat.row3.z]]; }; <> UnitSquareToQuadrilateral: PUBLIC PROC [p1, p2, p3, p4: Pair] RETURNS [Matrix] ~ { <> <<{S0, S1, S2, S3}) to an arbitrary quadrilateral, {P0, P1, P2, P3}. Thus,>> << [Artwork node; type 'Artwork on' to command tool] >> <> <
> << [Artwork node; type 'Artwork on' to command tool] >> <> <> <> <<(0, 0, 1) . (N00, N10, N20) = x1w, or N20 = x1w.>> <<(0, 0, 1) . (N01, N11, N21) = y1w, or N21 = y1w.>> <<(0, 0, 1) . (N02, N12, N22) = w, or N22 = w = 1.>> <> <> <> <<>> Det2: PROC [a, b, c, d: REAL] RETURNS [REAL] ~ INLINE {RETURN[a*d-b*c]}; x0: REAL Ź p1.x; y0: REAL Ź p1.y; x1: REAL Ź p2.x; y1: REAL Ź p2.y; x2: REAL Ź p3.x; y2: REAL Ź p3.y; x3: REAL Ź p4.x; y3: REAL Ź p4.y; px: REAL Ź x0-x1+x2-x3; py: REAL Ź y0-y1+y2-y3; IF RealFns.AlmostZero[px, -21] AND RealFns.AlmostZero[py, -21] THEN RETURN[[[x2-x1, y2-y1, 0.0], [x1-x0, y1-y0, 0.0], [x0, y0, 1.0]]] -- affine ELSE { -- perspective dx1: REAL Ź x3-x2; dx2: REAL Ź x1-x2; dy1: REAL Ź y3-y2; dy2: REAL Ź y1-y2; del: REAL Ź Det2[dx1, dx2, dy1, dy2]; a, b, c, d, e, f, g, h: REAL; IF del = 0.0 THEN ERROR; -- should be Error["Degenerate quadrilateral"] <> << |b e h|>> << |c f 1|>> g Ź Det2[px, dx2, py, dy2]/del; <> h Ź Det2[dx1, px, dy1, py]/del; a Ź x3-x0+g*x3; b Ź x1-x0+h*x1; c Ź x0; d Ź y3-y0+g*y3; e Ź y1-y0+h*y1; f Ź y0; RETURN[[[a, d, g], [b, e, h], [c, f, 1.0]]]; }; }; QuadrilateralToRectangle: PUBLIC PROC [p1, p2, p3, p4: Pair, rectangle: Rectangle] RETURNS [transform: Matrix] ~ { transform Ź UnitSquareToQuadrilateral[p1, p2, p3, p4]; transform Ź Invert[transform]; transform Ź Scale[transform, [rectangle.w, rectangle.h]]; transform Ź Translate[transform, [rectangle.x, rectangle.y]]; }; END.