<> <> <> <> DIRECTORY Imager USING [Error, ErrorDesc], ImagerTransformation USING [FactoredTransformation, Rectangle, Transformation, TransformationRep], Real USING [FScale, RoundI, RoundLI], RealFns USING [ArcTanDeg, CosDeg, SinDeg, SqRt], Vector2 USING [VEC]; ImagerTransformationImpl: CEDAR PROGRAM IMPORTS Imager, Real, RealFns EXPORTS ImagerTransformation ~ BEGIN OPEN ImagerTransformation; VEC: TYPE ~ Vector2.VEC; < interpress -> Interpress -> imager -> ImagerInterpress works.>> <<>> formErrorDesc: Imager.ErrorDesc ~ [code: $invalidTransformation, explanation: "Transformation has illegal form"]; Validate: PROC [m: Transformation] ~ { form: NAT _ 0; tx, ty: INTEGER _ 0; IF m.form > 10 THEN form _ m.form ELSE IF m.b=0 AND m.d=0 THEN { form _ 1; IF m.a=+1 THEN { IF m.e=+1 THEN form _ 3 ELSE IF m.e=-1 THEN form _ 4 } ELSE IF m.a=-1 THEN { IF m.e=+1 THEN form _ 5 ELSE IF m.e=-1 THEN form _ 6 }; } ELSE IF m.a=0 AND m.e=0 THEN { form _ 2; IF m.b=+1 THEN { IF m.d=+1 THEN form _ 7 ELSE IF m.d=-1 THEN form _ 8 } ELSE IF m.b=-1 THEN { IF m.d=+1 THEN form _ 9 ELSE IF m.d=-1 THEN form _ 10 }; }; IF ABS[m.c]<=INTEGER.LAST/2 THEN tx _ Real.RoundI[m.c]; IF ABS[m.f]<=INTEGER.LAST/2 THEN ty _ Real.RoundI[m.f]; m.form _ form; m.tx _ tx; m.ty _ ty; m.integerTrans _ (tx=m.c AND ty=m.f); }; Create: PUBLIC PROC [a, b, c, d, e, f: REAL] RETURNS [Transformation] ~ { new: Transformation ~ NEW[TransformationRep _ [a: a, b: b, c: c, d: d, e: e, f: f]]; Validate[new]; RETURN[new]; }; Scale: PUBLIC PROC [s: REAL] RETURNS [Transformation] ~ { new: Transformation ~ NEW[TransformationRep _ [a: s, b: 0, c: 0, d: 0, e: s, f: 0]]; Validate[new]; RETURN[new]; }; Scale2: PUBLIC PROC [s: VEC] RETURNS [Transformation] ~ { new: Transformation ~ NEW[TransformationRep _ [a: s.x, b: 0, c: 0, d: 0, e: s.y, f: 0]]; Validate[new]; RETURN[new]; }; Rotate: PUBLIC PROC [r: REAL] RETURNS [Transformation] ~ { cos: REAL ~ RealFns.CosDeg[r]; sin: REAL ~ RealFns.SinDeg[r]; new: Transformation ~ NEW[TransformationRep _ [a: cos, b: -sin, c: 0, d: sin, e: cos, f: 0]]; Validate[new]; RETURN[new]; }; Translate: PUBLIC PROC [t: VEC] RETURNS [Transformation] ~ { new: Transformation ~ NEW[TransformationRep _ [a: 1, b: 0, c: t.x, d: 0, e: 1, f: t.y]]; Validate[new]; RETURN[new]; }; Copy: PUBLIC PROC [m: Transformation] RETURNS [Transformation] ~ { IF m=NIL THEN RETURN[Scale[1]] ELSE RETURN[NEW[TransformationRep _ m^]]; }; Concat: PUBLIC PROC [m, n: Transformation] RETURNS [Transformation] ~ { new: Transformation ~ Copy[m]; ApplyPostConcat[new, n]; RETURN[new]; }; PreScale: PUBLIC PROC [m: Transformation, s: REAL] RETURNS [Transformation] ~ { new: Transformation ~ Copy[m]; ApplyPreScale[new, s]; RETURN[new]; }; PreScale2: PUBLIC PROC [m: Transformation, s: VEC] RETURNS [Transformation] ~ { new: Transformation ~ Copy[m]; ApplyPreScale2[new, s]; RETURN[new]; }; PreRotate: PUBLIC PROC [m: Transformation, r: REAL] RETURNS [Transformation] ~ { new: Transformation ~ Copy[m]; ApplyPreRotate[new, r]; RETURN[new]; }; PreTranslate: PUBLIC PROC [m: Transformation, t: VEC] RETURNS [Transformation] ~ { new: Transformation ~ Copy[m]; ApplyPreTranslate[new, t]; RETURN[new]; }; PostScale: PUBLIC PROC [m: Transformation, s: REAL] RETURNS [Transformation] ~ { new: Transformation ~ Copy[m]; ApplyPostScale[new, s]; RETURN[new]; }; PostScale2: PUBLIC PROC [m: Transformation, s: VEC] RETURNS [Transformation] ~ { new: Transformation ~ Copy[m]; ApplyPostScale2[new, s]; RETURN[new]; }; PostRotate: PUBLIC PROC [m: Transformation, r: REAL] RETURNS [Transformation] ~ { new: Transformation ~ Copy[m]; ApplyPostRotate[new, r]; RETURN[new]; }; PostTranslate: PUBLIC PROC [m: Transformation, t: VEC] RETURNS [Transformation] ~ { new: Transformation ~ Copy[m]; ApplyPostTranslate[new, t]; RETURN[new]; }; TranslateTo: PUBLIC PROC [m: Transformation, t: VEC] RETURNS [Transformation] ~ { new: Transformation ~ Copy[m]; ApplyTranslateTo[new, t]; RETURN[new]; }; Cat: PUBLIC PROC [m1, m2, m3, m4: Transformation _ NIL] RETURNS [Transformation] ~ { new: Transformation ~ NEW[TransformationRep]; ApplyCat[new, m1, m2, m3, m4]; RETURN[new]; }; ApplyPreConcat: PUBLIC PROC [m, p: Transformation] ~ { old: TransformationRep ~ m^; IF p=NIL THEN RETURN; IF p.form > 10 THEN ERROR Imager.Error[formErrorDesc]; m.a _ p.a*old.a + p.d*old.b; m.d _ p.a*old.d + p.d*old.e; m.b _ p.b*old.a + p.e*old.b; m.e _ p.b*old.d + p.e*old.e; m.c _ p.c*old.a + p.f*old.b + old.c; m.f _ p.c*old.d + p.f*old.e + old.f; Validate[m]; }; ApplyPreScale: PUBLIC PROC [m: Transformation, s: REAL] ~ { old: TransformationRep ~ m^; m.a _ s*old.a; m.d _ s*old.d; m.b _ s*old.b; m.e _ s*old.e; Validate[m]; }; ApplyPreScale2: PUBLIC PROC [m: Transformation, s: VEC] ~ { old: TransformationRep ~ m^; m.a _ s.x*old.a; m.d _ s.x*old.d; m.b _ s.y*old.b; m.e _ s.y*old.e; Validate[m]; }; ApplyPreRotate: PUBLIC PROC [m: Transformation, r: REAL] ~ { old: TransformationRep ~ m^; cos: REAL ~ RealFns.CosDeg[r]; sin: REAL ~ RealFns.SinDeg[r]; m.a _ cos*old.a + sin*old.b; m.d _ cos*old.d + sin*old.e; m.b _ cos*old.b - sin*old.a; m.e _ cos*old.e - sin*old.d; Validate[m]; }; ApplyPreTranslate: PUBLIC PROC [m: Transformation, t: VEC] ~ { old: TransformationRep ~ m^; m.c _ t.x*old.a + t.y*old.b + old.c; m.f _ t.x*old.d + t.y*old.e + old.f; Validate[m]; }; ApplyPostConcat: PUBLIC PROC [m, p: Transformation] ~ { old: TransformationRep ~ m^; IF p=NIL THEN RETURN; IF m.form > 10 THEN ERROR Imager.Error[formErrorDesc]; m.form _ p.form; m.a _ old.a*p.a + old.d*p.b; m.d _ old.a*p.d + old.d*p.e; m.b _ old.b*p.a + old.e*p.b; m.e _ old.b*p.d + old.e*p.e; m.c _ old.c*p.a + old.f*p.b + p.c; m.f _ old.c*p.d + old.f*p.e + p.f; Validate[m]; }; ApplyPostScale: PUBLIC PROC [m: Transformation, s: REAL] ~ { old: TransformationRep ~ m^; IF m.form > 10 THEN ERROR Imager.Error[formErrorDesc]; m.a _ old.a*s; m.d _ old.d*s; m.b _ old.b*s; m.e _ old.e*s; m.c _ old.c*s; m.f _ old.f*s; Validate[m]; }; ApplyPostScale2: PUBLIC PROC [m: Transformation, s: VEC] ~ { old: TransformationRep ~ m^; IF m.form > 10 THEN ERROR Imager.Error[formErrorDesc]; m.a _ old.a*s.x; m.d _ old.d*s.y; m.b _ old.b*s.x; m.e _ old.e*s.y; m.c _ old.c*s.x; m.f _ old.f*s.y; Validate[m]; }; ApplyPostRotate: PUBLIC PROC [m: Transformation, r: REAL] ~ { old: TransformationRep ~ m^; cos: REAL ~ RealFns.CosDeg[r]; sin: REAL ~ RealFns.SinDeg[r]; IF m.form > 10 THEN ERROR Imager.Error[formErrorDesc]; m.a _ old.a*cos - old.d*sin; m.d _ old.a*sin + old.d*cos; m.b _ old.b*cos - old.e*sin; m.e _ old.b*sin + old.e*cos; m.c _ old.c*cos - old.f*sin; m.f _ old.c*sin + old.f*cos; Validate[m]; }; ApplyPostTranslate: PUBLIC PROC [m: Transformation, t: VEC] ~ { old: TransformationRep ~ m^; IF m.form > 10 THEN ERROR Imager.Error[formErrorDesc]; m.c _ old.c+t.x; m.f _ old.f+t.y; Validate[m]; }; ApplyTranslateTo: PUBLIC PROC [m: Transformation, t: VEC] ~ { tx, ty: INTEGER _ 0; IF m.form > 10 THEN ERROR Imager.Error[formErrorDesc]; IF ABS[t.x]<=INTEGER.LAST/2 THEN tx _ Real.RoundI[t.x]; IF ABS[t.y]<=INTEGER.LAST/2 THEN ty _ Real.RoundI[t.y]; m.c _ t.x; m.f _ t.y; m.tx _ tx; m.ty _ ty; m.integerTrans _ (tx=t.x AND ty=t.y); }; ApplyCat: PUBLIC PROC [m: Transformation, m1, m2, m3, m4: Transformation _ NIL] ~ { count: NAT _ 0; Apply: PROC [p: Transformation] ~ { IF count=0 THEN m^ _ p^ ELSE ApplyPostConcat[m, p]; count _ count+1; }; IF m1#NIL THEN Apply[m1]; IF m2#NIL THEN Apply[m2]; IF m3#NIL THEN Apply[m3]; IF m4#NIL THEN Apply[m4]; IF count=0 THEN { m^ _ [a: 1, b: 0, c: 0, d: 0, e: 1, f: 0]; Validate[m] }; }; Invert: PUBLIC PROC [m: Transformation] RETURNS [Transformation] ~ { det: REAL ~ m.a*m.e - m.d*m.b; a: REAL ~ m.e/det; d: REAL ~ -m.d/det; b: REAL ~ -m.b/det; e: REAL ~ m.a/det; c: REAL ~ (m.b*m.f - m.e*m.c)/det; f: REAL ~ (m.d*m.c - m.a*m.f)/det; RETURN[Create[a: a, b: b, c: c, d: d, e: e, f: f]]; }; Transform: PUBLIC PROC [m: Transformation, v: VEC] RETURNS [VEC] ~ { IF m=NIL THEN RETURN[v]; SELECT m.form FROM 0 => RETURN[[v.x*m.a + v.y*m.b + m.c, v.x*m.d + v.y*m.e + m.f]]; 1 => RETURN[[v.x*m.a + m.c, v.y*m.e + m.f]]; 2 => RETURN[[v.y*m.b + m.c, v.x*m.d + m.f]]; 3 => RETURN[[m.c + v.x, m.f + v.y]]; 4 => RETURN[[m.c + v.x, m.f - v.y]]; 5 => RETURN[[m.c - v.x, m.f + v.y]]; 6 => RETURN[[m.c - v.x, m.f - v.y]]; 7 => RETURN[[m.c + v.y, m.f + v.x]]; 8 => RETURN[[m.c + v.y, m.f - v.x]]; 9 => RETURN[[m.c - v.y, m.f + v.x]]; 10 => RETURN[[m.c - v.y, m.f - v.x]]; ENDCASE => ERROR Imager.Error[formErrorDesc]; }; TransformVec: PUBLIC PROC [m: Transformation, v: VEC] RETURNS [VEC] ~ { IF m=NIL THEN RETURN[v]; SELECT m.form FROM 0 => RETURN[[v.x*m.a + v.y*m.b, v.x*m.d + v.y*m.e]]; 1 => RETURN[[v.x*m.a, v.y*m.e]]; -- b=0, d=0 2 => RETURN[[v.y*m.b, v.x*m.d]]; -- a=0, e=0 3 => RETURN[[+ v.x, + v.y]]; -- a=+1, b=0, d=0, e=+1 4 => RETURN[[+ v.x, - v.y]]; -- a=+1, b=0, d=0, e=-1 5 => RETURN[[- v.x, + v.y]]; -- a=-1, b=0, d=0, e=+1 6 => RETURN[[- v.x, - v.y]]; -- a=-1, b=0, d=0, e=-1 7 => RETURN[[+ v.y, + v.x]]; -- a=0, b=+1, d=+1, e=0 8 => RETURN[[+ v.y, - v.x]]; -- a=0, b=+1, d=-1, e=0 9 => RETURN[[- v.y, + v.x]]; -- a=0, b=-1, d=+1, e=0 10 => RETURN[[- v.y, - v.x]]; -- a=0, b=-1, d=-1, e=0 ENDCASE => ERROR Imager.Error[formErrorDesc]; }; InverseTransform: PUBLIC PROC [m: Transformation, v: VEC] RETURNS [VEC] ~ { IF m=NIL THEN RETURN[v]; RETURN[InverseTransformVec[m, [v.x - m.c, v.y - m.f]]]; }; InverseTransformVec: PUBLIC PROC [m: Transformation, v: VEC] RETURNS [VEC] ~ { IF m=NIL THEN RETURN[v]; SELECT m.form FROM 0 => { D: REAL ~ m.a*m.e-m.b*m.d; RETURN[[(v.x*m.e-v.y*m.b)/D, (v.y*m.a-v.x*m.d)/D]] }; 1 => RETURN[[v.x/m.a, v.y/m.e]]; -- b=0, d=0 2 => RETURN[[v.y/m.d, v.x/m.b]]; -- a=0, e=0 3 => RETURN[[+ v.x, + v.y]]; -- a=+1, b=0, d=0, e=+1 4 => RETURN[[+ v.x, - v.y]]; -- a=+1, b=0, d=0, e=-1 5 => RETURN[[- v.x, + v.y]]; -- a=-1, b=0, d=0, e=+1 6 => RETURN[[- v.x, - v.y]]; -- a=-1, b=0, d=0, e=-1 7 => RETURN[[+ v.y, + v.x]]; -- a=0, b=+1, d=+1, e=0 8 => RETURN[[- v.y, + v.x]]; -- a=0, b=+1, d=-1, e=0 9 => RETURN[[+ v.y, - v.x]]; -- a=0, b=-1, d=+1, e=0 10 => RETURN[[- v.y, - v.x]]; -- a=0, b=-1, d=-1, e=0 ENDCASE => ERROR Imager.Error[formErrorDesc]; }; InlineRound: PROC [x: REAL] RETURNS [REAL] ~ INLINE { RETURN[IF ABS[x]<100000000B THEN REAL[Real.RoundLI[x]] ELSE x]; }; DRound: PUBLIC PROC [v: VEC] RETURNS [VEC] ~ { RETURN[[InlineRound[v.x], InlineRound[v.y]]] }; RoundXY: PUBLIC PROC [m: Transformation, v: VEC] RETURNS [VEC] ~ { RETURN[InverseTransform[m, DRound[Transform[m, v]]]]; }; RoundXYVec: PUBLIC PROC [m: Transformation, v: VEC] RETURNS [VEC] ~ { RETURN[InverseTransformVec[m, DRound[TransformVec[m, v]]]]; }; TransformRectangle: PUBLIC PROC [m: Transformation, r: Rectangle] RETURNS [Rectangle] ~ { IF m=NIL THEN RETURN[r]; IF m.form=0 THEN { p0: VEC ~ Transform[m, [r.x, r.y]]; p1: VEC ~ Transform[m, [r.x+r.w, r.y]]; p2: VEC ~ Transform[m, [r.x+r.w, r.y+r.h]]; p3: VEC ~ Transform[m, [r.x, r.y+r.h]]; xmin, xmax: REAL _ p0.x; ymin, ymax: REAL _ p0.y; IF p1.xxmax THEN xmax _ p1.x; IF p1.yymax THEN ymax _ p1.y; IF p2.xxmax THEN xmax _ p2.x; IF p2.yymax THEN ymax _ p2.y; IF p3.xxmax THEN xmax _ p3.x; IF p3.yymax THEN ymax _ p3.y; RETURN[[x: xmin, y: ymin, w: xmax-xmin, h: ymax-ymin]]; } ELSE { p0: VEC ~ Transform[m, [r.x, r.y]]; p2: VEC ~ Transform[m, [r.x+r.w, r.y+r.h]]; xmin, xmax: REAL _ p0.x; ymin, ymax: REAL _ p0.y; IF p2.xxmax THEN xmax _ p2.x; IF p2.yymax THEN ymax _ p2.y; RETURN[[x: xmin, y: ymin, w: xmax-xmin, h: ymax-ymin]]; }; }; InverseTransformRectangle: PUBLIC PROC [m: Transformation, r: Rectangle] RETURNS [Rectangle] ~ { IF m=NIL THEN RETURN[r]; IF m.form=0 THEN { p0: VEC ~ InverseTransform[m, [r.x, r.y]]; p1: VEC ~ InverseTransform[m, [r.x+r.w, r.y]]; p2: VEC ~ InverseTransform[m, [r.x+r.w, r.y+r.h]]; p3: VEC ~ InverseTransform[m, [r.x, r.y+r.h]]; xmin, xmax: REAL _ p0.x; ymin, ymax: REAL _ p0.y; IF p1.xxmax THEN xmax _ p1.x; IF p1.yymax THEN ymax _ p1.y; IF p2.xxmax THEN xmax _ p2.x; IF p2.yymax THEN ymax _ p2.y; IF p3.xxmax THEN xmax _ p3.x; IF p3.yymax THEN ymax _ p3.y; RETURN[[x: xmin, y: ymin, w: xmax-xmin, h: ymax-ymin]]; } ELSE IF m.form <= 10 THEN { p0: VEC ~ InverseTransform[m, [r.x, r.y]]; p2: VEC ~ InverseTransform[m, [r.x+r.w, r.y+r.h]]; xmin, xmax: REAL _ p0.x; ymin, ymax: REAL _ p0.y; IF p2.xxmax THEN xmax _ p2.x; IF p2.yymax THEN ymax _ p2.y; RETURN[[x: xmin, y: ymin, w: xmax-xmin, h: ymax-ymin]]; } ELSE ERROR Imager.Error[formErrorDesc]; }; RestrictAbsTo: PROC [limit: REAL, real: REAL] RETURNS [REAL] ~ { <> limit _ ABS[limit]; WHILE real < -limit DO real _ real + 2*limit ENDLOOP; WHILE real >= limit DO real _ real - 2*limit ENDLOOP; RETURN[real] }; NumericalInstability: SIGNAL[real: REAL] ~ CODE; maxRelError: REAL _ 1.0e-4; Factor: PUBLIC PROC [m: Transformation] RETURNS [FactoredTransformation] ~ { IF m=NIL THEN RETURN[[r1: 0, s: [1, 1], r2: 0, t: [0, 0]]]; SELECT m.form FROM 0 => { t: TransformationRep ~ m^; theta: REAL ~ RestrictAbsTo[45, 0.5*RealFns.ArcTanDeg[2*(t.a*t.b+t.d*t.e), (t.a*t.a-t.b*t.b+t.d*t.d-t.e*t.e)]]; cosTheta: REAL ~ RealFns.CosDeg[theta]; sinTheta: REAL ~ RealFns.SinDeg[theta]; aa: REAL ~ cosTheta*t.a + sinTheta*t.b; bb: REAL ~ -sinTheta*t.a + cosTheta*t.b; dd: REAL ~ cosTheta*t.d + sinTheta*t.e; ee: REAL ~ -sinTheta*t.d + cosTheta*t.e; phi: REAL ~ RestrictAbsTo[90, IF ABS[aa]>ABS[ee] THEN RealFns.ArcTanDeg[-dd, aa] ELSE RealFns.ArcTanDeg[bb, ee]]; cosPhi: REAL ~ RealFns.CosDeg[phi]; sinPhi: REAL ~ RealFns.SinDeg[phi]; a: REAL ~ aa*cosPhi - dd*sinPhi; b: REAL ~ bb*cosPhi - ee*sinPhi; d: REAL ~ aa*sinPhi + dd*cosPhi; e: REAL ~ bb*sinPhi + ee*cosPhi; mag: REAL ~ ABS[a]+ABS[d]; IF ABS[b]+ABS[d] > maxRelError*mag THEN SIGNAL NumericalInstability[ABS[b]+ABS[d]]; RETURN[[r1: -theta, s: [a, e], r2: -phi, t: [t.c, t.f]]] }; 1 => RETURN[[r1: 0, s: [m.a, m.e], r2: 0, t: [m.c, m.f]]]; 2 => { IF m.b > 0 THEN RETURN[[r1: 0, s: [-m.d, m.b], r2: -90, t: [m.c, m.f]]] ELSE RETURN[[r1: 0, s: [m.d, -m.b], r2: 90, t: [m.c, m.f]]]; }; 3 => RETURN[[r1: 0, s: [1, 1], r2: 0, t: [m.c, m.f]]]; 4 => RETURN[[r1: 0, s: [1, -1], r2: 0, t: [m.c, m.f]]]; 5 => RETURN[[r1: 0, s: [-1, 1], r2: 0, t: [m.c, m.f]]]; 6 => RETURN[[r1: 0, s: [-1, -1], r2: 0, t: [m.c, m.f]]]; 7 => RETURN[[r1: 0, s: [1, -1], r2: 90, t: [m.c, m.f]]]; 8 => RETURN[[r1: 0, s: [1, 1], r2: -90, t: [m.c, m.f]]]; 9 => RETURN[[r1: 0, s: [1, 1], r2: 90, t: [m.c, m.f]]]; 10 => RETURN[[r1: 0, s: [-1, 1], r2: 90, t: [m.c, m.f]]]; ENDCASE => ERROR Imager.Error[formErrorDesc]; }; Equal: PUBLIC PROC [s, t: Transformation] RETURNS [BOOL] ~ { RETURN [s.a=t.a AND s.b=t.b AND s.c=t.c AND s.d=t.d AND s.e=t.e AND s.f=t.f]; }; WithinQuarterPixel: PROC [p, q: VEC] RETURNS [BOOL] ~ { RETURN [ABS[p.x-q.x] <= 0.25 AND ABS[p.y-q.y] <= 0.25] }; CloseEnough: PUBLIC PROC [s, t: Transformation, range: REAL] RETURNS [BOOL] ~ { p00: VEC _ Transform[s, InverseTransform[t, [0, 0]]]; p10: VEC _ Transform[s, InverseTransform[t, [range, 0]]]; p01: VEC _ Transform[s, InverseTransform[t, [0, range]]]; RETURN [ WithinQuarterPixel[p00, [0, 0]] AND WithinQuarterPixel[p10, [range, 0]] AND WithinQuarterPixel[p10, [0, range]] ] }; CloseToTranslation: PUBLIC PROC [s, t: Transformation, range: REAL] RETURNS [BOOL] ~ { p10: VEC _ TransformVec[s, InverseTransformVec[t, [range, 0]]]; p01: VEC _ TransformVec[s, InverseTransformVec[t, [0, range]]]; RETURN [ WithinQuarterPixel[p10, [range, 0]] AND WithinQuarterPixel[p10, [0, range]] ] }; SingularValues: PUBLIC PROC [m: Transformation] RETURNS [VEC] ~ { t1: REAL _ m.a*m.a + m.d*m.d; t2: REAL _ m.b*m.b + m.e*m.e; t3: REAL _ t1 - t2; t4: REAL _ Real.FScale[m.a*m.b + m.d*m.e, 1]; s: REAL _ RealFns.SqRt[t3*t3+t4*t4]; x: REAL _ RealFns.SqRt[Real.FScale[t1+t2+s, -1]]; y: REAL _ RealFns.SqRt[Real.FScale[t1+t2-s, -1]]; RETURN[[x, y]] }; fuzz: REAL _ 2.0e-6; Different: PROC [a, b: REAL] RETURNS [BOOL] ~ { IF ABS[a-b] < fuzz THEN RETURN [FALSE]; RETURN [ABS[a - b]/MAX[ABS[a], ABS[b]] > fuzz]; }; ShouldBeSameT: PROC [s, t: Transformation] ~ { IF Different[s.a, t.a] THEN ERROR; IF Different[s.b, t.b] THEN ERROR; IF Different[s.c, t.c] THEN ERROR; IF Different[s.d, t.d] THEN ERROR; IF Different[s.e, t.e] THEN ERROR; IF Different[s.f, t.f] THEN ERROR; }; ShouldBeSamePt: PROC [p, q: VEC] ~ { IF Different[p.x, q.x] THEN ERROR; IF Different[p.y, q.y] THEN ERROR; }; TestTransformation: PROC [r: Transformation] ~ { p: VEC _ [3.14159, -265.35]; ShouldBeSameT[PreTranslate[r, p], Concat[Translate[p], r]]; ShouldBeSameT[Rotate[0], Concat[r, Invert[r]]]; ShouldBeSamePt[InverseTransform[r, p], Transform[Invert[r], p]]; ShouldBeSamePt[InverseTransformVec[r, p], TransformVec[Invert[r], p]]; ShouldBeSamePt[TransformVec[r, TransformVec[r, p]], TransformVec[Concat[r, r], p]]; ShouldBeSamePt[SingularValues[r], SingularValues[PreRotate[r, 33.5]]]; }; SelfTest: PROC ~ { ShouldBeSameT[Rotate[90], Concat[Rotate[45], Rotate[45]]]; TestTransformation[PreTranslate[PreScale2[Rotate[61], [1, 5]], [12, 6]]]; TestTransformation[Create[1, 2, 4, 9, 3, 1]]; TestTransformation[Create[3, 1, 4, 1, 5, 9]]; TestTransformation[Create[2, 0, 0, 0, 2, 0]]; TestTransformation[Create[0, 3, 0, 4, 0, 0]]; TestTransformation[Create[1, 0, 0, 0, 1, 0]]; TestTransformation[Create[1, 0, 0, 0, -1, 808]]; TestTransformation[Create[-1, 0, 0, 0, 1, 0]]; TestTransformation[Create[-1, 0, 0, 0, -1, 0]]; TestTransformation[Create[0, 1, 0, 1, 0, 0]]; TestTransformation[Create[0, 1, 0, -1, 0, 0]]; TestTransformation[Create[0, -1, 0, 1, 0, 0]]; TestTransformation[Create[0, -1, 0, -1, 0, 0]]; }; <> <> <> <> <> <> <> <> <> <> <> <> <<};>> <<>> <> <> <> <> <> <> <> <> <new.x DO x0 _ x0-1 ENDLOOP;>> <new.y DO y0 _ y0-1 ENDLOOP;>> <> <> <> <<};>> <<>> END.