-- CGMatrixImpl.mesa
-- Last changed by Doug Wyatt, August 23, 1982 4:14 pm

DIRECTORY
CGMatrix USING [ErrorType, Id, Ref, Rep],
CGStorage USING [qZone],
GraphicsBasic USING [Transformation, Vec],
Real USING [RealException];

CGMatrixImpl: CEDAR MONITOR
IMPORTS CGStorage, Real
EXPORTS CGMatrix = { OPEN CGMatrix, GraphicsBasic;

repZone: ZONE = CGStorage.qZone;

lastid: Id ← 1; -- Monitored global!

identityMatrix: Rep = [m: [a: 1, b: 0, c: 0, d: 1, e: 0, f: 0],
id: 0, trivial: TRUE, rectangular: TRUE, determinant: 1];

Error: PUBLIC ERROR[type: ErrorType] = CODE;

NewID: ENTRY PROC RETURNS[Id] = INLINE {
IF lastid<LAST[Id] THEN lastid ← lastid + 1
ELSE lastid ← 1; -- just wrap around
RETURN[lastid];
};

SetHints: PROC[self: Ref] = {
trivial: BOOLEAN ← self.m.a=1 AND self.m.b=0 AND self.m.c=0 AND self.m.d=1;
self.id ← IF trivial THEN 1 ELSE NewID[];
self.trivial ← trivial;
self.rectangular ← trivial OR (self.m.b=0 AND self.m.c=0) OR (self.m.a=0 AND self.m.d=0);
self.determinant ← self.m.a*self.m.d - self.m.b*self.m.c;
};

New: PUBLIC PROC RETURNS[Ref] = {
self: Ref ← repZone.NEW[Rep ← identityMatrix];
RETURN[self];
};

Make: PUBLIC PROC[m: Transformation] RETURNS[Ref] = {
self: Ref ← repZone.NEW[Rep ← [m: m,
id: 0, trivial: FALSE, rectangular: FALSE, determinant: 0]];
SetHints[self];
RETURN[self];
};

MapH: PUBLIC PROC[self: Ref, v: Vec, h: [0..1]] RETURNS[Vec] = {
x,y: REAL;
x ← v.x*self.m.a + v.y*self.m.c;
y ← v.x*self.m.b + v.y*self.m.d;
RETURN[IF h=1 THEN [x + self.m.e, y + self.m.f] ELSE [x,y]];
};

InvH: PUBLIC PROC[self: Ref, v: Vec, h: [0..1]] RETURNS[Vec] = {
ENABLE Real.RealException => GOTO Singular;
u: Vec ← IF h=1 THEN Inv0[self,v.x-self.m.e,v.y-self.m.f] ELSE Inv0[self,v.x,v.y];
RETURN[u]; EXITS Singular => ERROR Error[singularMatrix];
};

Inv0: PROC[self: Ref, x,y: REAL] RETURNS[Vec] = {
det: REAL ← self.determinant;
RETURN[[(x*self.m.d-y*self.m.c)/det,(y*self.m.a-x*self.m.b)/det]];
};

Concat: PUBLIC PROC[self: Ref, a,b,c,d: REAL] = {
IF self.trivial THEN {
self.m.a ← a; self.m.b ← b;
self.m.c ← c; self.m.d ← d;
}
ELSE {
aa: REAL ← self.m.a; bb: REAL ← self.m.b;
cc: REAL ← self.m.c; dd: REAL ← self.m.d;
self.m.a ← a*aa + b*cc; self.m.b ← a*bb + b*dd;
self.m.c ← c*aa + d*cc; self.m.d ← c*bb + d*dd;
};
SetHints[self];
};

ConcatTransformation: PUBLIC PROC[self: Ref, m: Transformation] = {
IF self.trivial THEN {
self.m.a ← m.a; self.m.b ← m.b;
self.m.c ← m.c; self.m.d ← m.d;
self.m.e ← self.m.e + m.e;
self.m.f ← self.m.f + m.f;
}
ELSE {
aa: REAL ← self.m.a; bb: REAL ← self.m.b;
cc: REAL ← self.m.c; dd: REAL ← self.m.d;
self.m.a ← m.a*aa + m.b*cc; self.m.b ← m.a*bb + m.b*dd;
self.m.c ← m.c*aa + m.d*cc; self.m.d ← m.c*bb + m.d*dd;
self.m.e ← self.m.e + m.e*aa + m.f*cc;
self.m.f ← self.m.f + m.e*bb + m.f*dd;
};
SetHints[self];
};

Invert: PUBLIC PROC[self: Ref] = {
aa: REAL ← self.m.a; bb: REAL ← self.m.b;
cc: REAL ← self.m.c; dd: REAL ← self.m.d;
ee: REAL ← self.m.e; ff: REAL ← self.m.f;
det: REAL ← aa*dd - bb*cc;
self.m.a ← dd/det;
self.m.b ← -bb/det;
self.m.c ← -cc/det;
self.m.d ← aa/det;
self.m.e ← (cc*ff-dd*ee)/det;
self.m.f ← (bb*ee-aa*ff)/det;
SetHints[self];
};

Copy: PUBLIC PROC[self: Ref] RETURNS[Ref] = {
RETURN[repZone.NEW[Rep ← self^]];
};

}.