--ContourImpl.mesa
-- Written by: Maureen Stone December 1, 1983 12:04 pm
DIRECTORY
Graphics,
Curve USING [StartSamples, defaultHandle, AddSample, ResetContours, AddContour, NextContour, MoveTo, DrawTo],
MessageWindow USING[Append],
GraphicsOps USING[NewAisImage],
TJaMGraphics USING[Painter],
JaMFnsDefs,
Vector USING [Vec],
Rope USING [ROPE, Cat],
Convert USING[ValueToRope],
ConvertUnsafe USING[AppendRope, ToRope],
AIS USING [OpenFile, FRef, WRef, OpenWindow,CloseFile, ReadSample, CloseWindow,Error];

ContourImpl: PROGRAM
IMPORTS Graphics, GraphicsOps, TJaMGraphics, ConvertUnsafe, AIS, JaMFnsDefs, Convert, MessageWindow, Rope, Curve
EXPORTS = BEGIN

Vec: TYPE = Vector.Vec;
VecSequence: TYPE = REF VecSequenceRec;
VecSequenceRec: TYPE = RECORD[element: SEQUENCE length: NAT OF Vec];
ListSamples: TYPE = REF ListSamplesRec;
ListSamplesRec: TYPE = RECORD[next: ListSamples, first, last: Sample];
Sample: TYPE = REF SampleRec;
SampleRec: TYPE = RECORD[next: Sample, value: Vec];
ListHead: TYPE = REF ListHeadRec;
ListHeadRec: TYPE = RECORD[first, last: ListSamples];
Where: TYPE = {head, tail};

ais: AIS.FRef _ NIL;
aisName: Rope.ROPE;
aisWindow: AIS.WRef;
defaultTransitionCount: NAT _ 20;
defaultSampleCount: NAT _ 100;
currentEdges: ListHead;    --samples in edge order
currentContours: ListHead;    --samples in contour order
tValue: REAL _ 128.5;    --global threshhold values
doPaint: BOOLEAN _ TRUE;    --display control boolean for NewTransition
BadCell: PUBLIC SIGNAL=CODE;
transitions: ListHead;    --samples for one scan line
freeList: ListHead;    --freed empty lists from transitions

SetUpAIS: PUBLIC PROCEDURE [name: Rope.ROPE] = {
s: LONG STRING_[64];
ConvertUnsafe.AppendRope[to: s, from: name];    --put name in string
IF ais#NIL THEN FreeAIS[];
ais_AIS.OpenFile[s, FALSE];
aisWindow _ AIS.OpenWindow[ais];
aisName _ name;
};
JSetUpAIS: PROCEDURE = {
ENABLE AIS.Error =>{JaMFnsDefs.JaMExec["(AIS file error) .print"]; CONTINUE};
s: LONG STRING_[64];
JaMFnsDefs.PopString[s];
SetUpAIS[ConvertUnsafe.ToRope[s]];
};
FreeAIS: PROCEDURE = {
AIS.CloseWindow[aisWindow]; aisWindow _ NIL;
AIS.CloseFile[ais];  ais _ NIL;
};
JFreeAIS: PROCEDURE = {
ENABLE AIS.Error =>{JaMFnsDefs.JaMExec["(AIS file error) .print"]; CONTINUE};
FreeAIS[];
};

SetTValue: PUBLIC PROCEDURE[value: REAL] = {tValue _ value};
JSetTValue: PROCEDURE = {SetTValue[JaMFnsDefs.PopReal[]]};

DrawAIS: PUBLIC PROCEDURE= {
Paint: PROC[dc: Graphics.Context] = {
imageRef: Graphics.ImageRef _ GraphicsOps.NewAisImage[aisName];
lx,rx,ty,by: REAL;
mark: Graphics.Mark _ Graphics.Save[dc];
[lx,by] _ AISToScreen[aisWindow.firstPixel,aisWindow.lastScan];
[rx,ty] _ AISToScreen[aisWindow.lastPixel,aisWindow.firstScan];
Graphics.Translate[dc,-lx,-by];
Graphics.ClipBox[dc,[lx,by,rx,ty]];
Graphics.DrawImage[dc,imageRef];
Graphics.Restore[dc,mark];
};
TJaMGraphics.Painter[Paint];
};
AISToScreen: PROC [x,y: REAL] RETURNS[nx,ny: REAL] ={
nx _ x; ny _ aisWindow.lastScan-y;
};

InitContours: PROC = {
currentEdges _ NEW[ListHeadRec _ [first: NIL, last: NIL]];
currentContours _ NEW[ListHeadRec _ [first: NIL, last: NIL]];
transitions _ NEW[ListHeadRec _ [first: NIL, last: NIL]];
freeList _ NEW[ListHeadRec _ [first: NIL, last: NIL]];    --don't let it grow without bound
};

lineCount: INTEGER _ 0;
contourCount: INTEGER _ 0;
Outline: PUBLIC PROC = { 
v1,v2,v3,v4,x,y: CARDINAL;
lineCount _ 0;
contourCount _ 0;
InitContours[];
FOR y IN [aisWindow.firstScan..aisWindow.lastScan) DO
IF JaMFnsDefs.GetJaMBreak[] THEN EXIT;
<<this assumes the rd scan direction for ais files>>
FOR x IN [aisWindow.firstPixel..aisWindow.lastPixel) DO
v1 _ AIS.ReadSample[aisWindow,y,x];
v2 _ AIS.ReadSample[aisWindow,y,x+1];
v3 _ AIS.ReadSample[aisWindow,y+1,x+1];
v4 _ AIS.ReadSample[aisWindow,y+1,x];
ContourCell[v1,v2,v3,v4,x,y ! BadCell => EXIT];
ENDLOOP;
IF y MOD 20 = 0 THEN MessageWindow.Append[Rope.Cat["Scan line = ",Convert.ValueToRope[[signed[y]]]], TRUE];
ProcessTransitionBuffer[];
ENDLOOP;
};
Trans: TYPE = RECORD[x,y,d: REAL, dir:{nt,dl,ld}];
<<aisX and aisY relative to current AIS aisWindow>>
ContourCell: PROC[v1,v2,v3,v4:CARDINAL, aisX,aisY:REAL] = {
edgePts: ARRAY [0..3] OF Trans;
i,ldCount,dlCount,index,move: CARDINAL;
FindLD: PROC [start: CARDINAL] RETURNS[CARDINAL] = {
FOR k: CARDINAL IN [start..3] DO
IF edgePts[k].dir=ld THEN RETURN[k];
ENDLOOP;
ERROR;
};
FindDL: PROC [start: CARDINAL] RETURNS[CARDINAL] = {
FOR k: CARDINAL IN [start..3] DO
IF edgePts[k].dir=dl THEN RETURN[k]; ENDLOOP;
FOR k: CARDINAL IN [0..start) DO
IF edgePts[k].dir=dl THEN RETURN[k]; ENDLOOP;
ERROR;
};
<<v1  v2    clockwise>>
<<v4  v3>>
edgePts[0] _ ContourEdge[v1,v2];
edgePts[0].x _ aisX+edgePts[0].d; edgePts[0].y _ aisY;
edgePts[1] _ ContourEdge[v2,v3];
edgePts[1].x _ aisX+1; edgePts[1].y _ aisY+edgePts[1].d;
<<switch the order here for numerical reasons.  Must then invert the dir>>
edgePts[2] _ ContourEdge[v4,v3];
edgePts[2].x _ aisX+edgePts[2].d; edgePts[2].y _ aisY+1;
edgePts[2].dir _
(SELECT edgePts[2].dir FROM ld => dl, dl => ld, ENDCASE => nt);
edgePts[3] _ ContourEdge[v1,v4];
edgePts[3].x _ aisX; edgePts[3].y _ aisY+edgePts[3].d;
edgePts[3].dir _
(SELECT edgePts[3].dir FROM ld => dl, dl => ld, ENDCASE => nt);
ldCount _ dlCount _ 0;
FOR i IN [0..3] DO
SELECT edgePts[i].dir FROM
dl => dlCount _ dlCount+1;
ld => ldCount _ ldCount+1;
ENDCASE;
ENDLOOP;
IF dlCount#ldCount THEN {
MessageWindow.Append["odd number of transitions", TRUE];
SIGNAL BadCell;
};
IF dlCount>2 THEN {
MessageWindow.Append["too many transitions", TRUE];
SIGNAL BadCell;
};
IF dlCount=0 THEN RETURN;--no transitions
<<ok, now have transitions.  Need to connect ld to dl>>
<<this prefers solid dark areas>>
move _ index _ FindLD[0];
index _ FindDL[index];
NewTransition[edgePts[move].x,edgePts[move].y,edgePts[index].x,edgePts[index].y];
lineCount _ lineCount+1;
IF dlCount=2 THEN BEGIN
move _ index _ FindLD[index];
index _ FindDL[index];
NewTransition[edgePts[move].x,edgePts[move].y,
edgePts[index].x,edgePts[index].y];
lineCount _ lineCount+1;
END;
};

ContourEdge: PROC[v0,v1: REAL] RETURNS[trans: Trans] = { 
trans _ [x:0, y:0, d:0, dir:nt];
IF v0<tValue AND v1<tValue THEN RETURN;
IF v0>tValue AND v1>tValue THEN RETURN;
trans.d _ (v0-tValue)/(v0-v1);
IF v0<tValue THEN trans.dir _ ld ELSE trans.dir _ dl;
};

OutlineBlackCenter: PUBLIC PROC = {OutlineBlack[ContourBlackCenter]};
OutlineBlackEdge: PUBLIC PROC = {OutlineBlack[ContourBlackEdge]};
OutlineBlack: PUBLIC PROC [contourProc: PROC [pts: ARRAY [0..4) OF CARDINAL, aisX,aisY:REAL]] = { 
x,y: CARDINAL;
pts: ARRAY [0..4) OF CARDINAL;
lineCount _ 0;
contourCount _ 0;
InitContours[];
FOR y IN [aisWindow.firstScan..aisWindow.lastScan) DO
IF JaMFnsDefs.GetJaMBreak[] THEN EXIT;
<<this assumes the rd scan direction for ais files>>
FOR x IN [aisWindow.firstPixel..aisWindow.lastPixel) DO
pts[0] _ AIS.ReadSample[aisWindow,y,x];
pts[1]_ AIS.ReadSample[aisWindow,y,x+1];
pts[2] _ AIS.ReadSample[aisWindow,y+1,x+1];
pts[3] _ AIS.ReadSample[aisWindow,y+1,x];
contourProc[pts,x,y ! BadCell => EXIT];
ENDLOOP;
IF y MOD 20 = 0 THEN MessageWindow.Append[Rope.Cat["Scan line = ",Convert.ValueToRope[[signed[y]]]], TRUE];
ProcessTransitionBuffer[];
ENDLOOP;
};
ContourBlackCenter: PROC[pts: ARRAY [0..4) OF CARDINAL, aisX,aisY:REAL] = {
<<0  1    clockwise  goes to 0123 in the code>>
<<3  2>>
trans: PROC [dir: CARDINAL] = {
SELECT dir FROM
01 => NewTransition[aisX,aisY,aisX+1, aisY];
12 => NewTransition[aisX+1,aisY,aisX+1,aisY+1];
23 => NewTransition[aisX+1,aisY+1,aisX,aisY+1];
30 => NewTransition[aisX, aisY+1, aisX,aisY];
ENDCASE => ERROR;
};
i,code, mask, count: CARDINAL _ 0;
mask _ 10B;
FOR i IN [0..4) DO
IF pts[i]=tValue THEN {code _ code+mask; count _ count+1};
mask _ mask/2;
ENDLOOP;
IF count<2 OR count >3 THEN RETURN;    --6 cases
SELECT code FROM
12B, 5B => RETURN;    --the 2 diagonals do not contribute
14B => trans[01];    --4 valid single transitions
11B => trans[30];
6B => trans[12];
3B => trans[23];
16B => {trans[01]; trans[12]};    --4 valid double transitions
13B => {trans[30]; trans[23]};    --give them in left to right order
15B => {trans[30]; trans[01]};
7B => {trans[23]; trans[12]};
ENDCASE => ERROR;
};
ContourBlackEdge: PROC[pts: ARRAY [0..4) OF CARDINAL, aisX,aisY:REAL] = {
<<0  1    clockwise  goes to 0123 in the code>>
<<3  2>>
cx: REAL _ aisX+1;    --shorthand for center point of cell
cy: REAL _ aisY+1;
path: PROC [dx0,dy0,dx1,dy1: REAL] = {  --all paths go thru center
IF dx1 < dx0 THEN {    --make transitions in raster order
NewTransition[cx,cy,cx+dx1,cy+dy1];
NewTransition[cx+dx0,cy+dy0,cx,cy];
}
ELSE {
NewTransition[cx+dx0,cy+dy0,cx,cy];
NewTransition[cx,cy,cx+dx1,cy+dy1];
};
};
i,code, mask: CARDINAL _ 0;
mask _ 10B;
FOR i IN [0..4) DO
IF pts[i]=tValue THEN code _ code+mask;
mask _ mask/2;
ENDLOOP;
<<each case generates 2 or 4 transitions.  However, each transition will
--be seen exactly twice.  So, only generate half the transitions each time>>
SELECT code FROM
10b => path[0,-1,-1,0];    --single black pixel
2 => path[0,1,1,0];
13b => path[0,-1,1,0];    --single white pixel
16b => path[0,1,-1,0];
12b => {path[0,1,-1,0]; path[0,-1,1,0]};  --diagonal
3 => NewTransition[cx,cy,cx+1,cy];    --right
14b => NewTransition[cx,cy,cx-1,cy];    --left
11b => NewTransition[cx,cy-1, cx,cy];    --down
6b => NewTransition[cx,cy+1,cx,cy];    --up
ENDCASE => RETURN;
};
<<we keep all the transitions for one scan line in a list of samples>>
<<at the end of each scan line we add the collected samples to the edges>>
<<this buffer solves the problem of "horizontal" edges>>
NewTransition: PROC[x0,y0,x1,y1: REAL] = {
paint: PROC[dc: Graphics.Context] = {
  Graphics.SetCP[dc,p0.x,p0.y];
  Graphics.DrawTo[dc,p1.x,p1.y];
  };
p0, p1: Vec;
list: ListSamples;
found: BOOLEAN _ FALSE;
[p0.x,p0.y] _ AISToScreen[x0,y0];
[p1.x,p1.y] _ AISToScreen[x1,y1];
IF doPaint THEN TJaMGraphics.Painter[paint];
FOR list _ transitions.first, list.next UNTIL list=NIL DO
IF list.last.value=p0 THEN {
found _ TRUE;
AddSample[p: p1, list: list, where: tail];
EXIT;
};
IF list.first.value=p1 THEN {
found _ TRUE;
AddSample[p: p0, list: list, where: head];
EXIT;
};
ENDLOOP;
IF ~found THEN {
list _ NewList[p: p0];
AddSample[p: p1, list: list, where: tail];
AppendList[head: transitions, list: list];
};
};
ProcessTransitionBuffer: PROC = {
inTailOfEdge: PROC[p: Vec] RETURNS [found: BOOLEAN] = {
FOR edge _ currentEdges.first, edge _ edge.next UNTIL edge=NIL DO
IF edge.last.value=p THEN RETURN[TRUE]; 
ENDLOOP;
RETURN[FALSE];
};
inHeadOfEdge: PROC[p: Vec] RETURNS [found: BOOLEAN] = {
FOR edge _ currentEdges.first, edge _ edge.next UNTIL edge=NIL DO
IF edge.first.value=p THEN RETURN[TRUE]; 
ENDLOOP;
RETURN[FALSE];
};

edge, list: ListSamples;    --set by in*OfEdge
IF transitions.first=NIL THEN RETURN;    --no transitions this scan line
IF currentEdges.first#NIL THEN    --first set of transitions are all disjoint edges
FOR list _ transitions.first, list.next UNTIL list=NIL DO
IF inTailOfEdge[list.first.value] THEN
AddListSamples[from: list, list: edge, where: tail]
ELSE IF inHeadOfEdge[list.last.value]
THEN AddListSamples[from: list, list: edge, where: head];
ENDLOOP;
DO
list _ RemoveFirst[head: transitions ! Empty => EXIT];
IF list.first#NIL THEN AppendList[head: currentEdges, list: list]
ELSE AppendList[head: freeList, list: list];
ENDLOOP;
};
MakeContours: PROC = {
contour: ListSamples;
allFound: SIGNAL = CODE;
findRest: PROC = {
edge: ListSamples _ currentEdges.first;
foundAny: BOOLEAN _ FALSE;
UNTIL edge=NIL OR contour.first.value=contour.last.value DO
IF edge.first.value=contour.last.value THEN {
list: ListSamples _ edge;
edge _ RemoveList[head: currentEdges, list: list];
AddListSamples[from: list, list: contour, where: tail];
foundAny _ TRUE;
}
ELSE IF edge.last.value=contour.first.value THEN {
list: ListSamples _ edge;
edge _ RemoveList[head: currentEdges, list: list];
AddListSamples[from: list, list: contour, where: head];
foundAny _ TRUE;
}
ELSE edge_edge.next;
ENDLOOP;
IF ~foundAny OR contour.first.value=contour.last.value THEN SIGNAL allFound;
};
DO
contour _ RemoveFirst[currentEdges ! Empty=>EXIT];
DO
findRest[! allFound => EXIT];    --adds pieces to contour
ENDLOOP;
IF contour.first.value=contour.last.value
THEN contour.first _ contour.first.next;    --remove duplicates
AppendList[head: currentContours, list: contour];
ENDLOOP;
};

CountEdges: PROC = {
count: NAT _ 0;
FOR edge: ListSamples _ currentEdges.first, edge.next UNTIL edge=NIL DO
count _ count+1;
ENDLOOP;
JaMFnsDefs.PushInteger[count];
};
CountContours: PROC = {
count: NAT _ 0;
FOR contour: ListSamples _ currentContours.first, contour.next UNTIL contour=NIL DO
count _ count+1;
ENDLOOP;
JaMFnsDefs.PushInteger[count];
};
JDrawContour: PROCEDURE = {
cn: INTEGER _ JaMFnsDefs.PopInteger[];
DrawContourNumber[cn];
};

DrawContourNumber: PUBLIC PROCEDURE [cn: NAT] = {
count: NAT _ 0;
FOR contour: ListSamples _ currentContours.first, contour.next UNTIL contour=NIL DO
IF count=cn THEN DrawContour[contour];
count _ count+1;
ENDLOOP;
};
DrawContour: PROCEDURE [contour: ListSamples] = {
Paint: PROC[dc: Graphics.Context] = {
Curve.MoveTo[dc, [contour.first.value.x, contour.first.value.y]];
FOR sa: Sample _ contour.first.next, sa.next UNTIL sa=NIL DO
Curve.DrawTo[dc, [sa.value.x,sa.value.y]];
ENDLOOP;
};
TJaMGraphics.Painter[Paint];
};
JSetSamples: PROCEDURE = {
cn: INTEGER _ JaMFnsDefs.PopInteger[];
SetSamples[cn];
};

SetSamples: PUBLIC PROCEDURE [cn: NAT] = {
count: NAT _ 0;
FOR contour: ListSamples _ currentContours.first, contour.next UNTIL contour=NIL DO
IF count=cn THEN SetContourSamples[contour];
count _ count+1;
ENDLOOP;
};
SetContourSamples: PROCEDURE [contour: ListSamples] = {
Curve.StartSamples[Curve.defaultHandle, contour.first.value.x,contour.first.value.y];
FOR sa: Sample _ contour.first.next, sa.next UNTIL sa=NIL DO
Curve.AddSample[Curve.defaultHandle, sa.value.x,sa.value.y];
ENDLOOP;
};
SetContours: PUBLIC PROCEDURE = {
first: BOOLEAN _ TRUE;
FOR contour: ListSamples _ currentContours.first, contour.next UNTIL contour=NIL DO
IF first THEN {Curve.ResetContours[Curve.defaultHandle]; first _ FALSE}
ELSE Curve.AddContour[Curve.defaultHandle];
SetContourSamples[contour];
ENDLOOP;
Curve.NextContour[Curve.defaultHandle];    --select first contour
};

NewList: PROC [p: Vec] RETURNS [list: ListSamples] = {
new: Sample _ NEW[SampleRec _ [NIL, p]];
list _ RemoveFirst[freeList ! Empty => {list _ NEW[ListSamplesRec]; CONTINUE};];
list.first _ list.last _ new;
};
<<next two procedures can assume that the lists are not empty>>
AddSample: PROC [p: Vec, list: ListSamples, where: Where] = {
new: Sample _ NEW[SampleRec _ [NIL, p]];
IF where=head THEN {new.next _ list.first; list.first _ new}
ELSE {list.last.next _ new; list.last _ new};
};
<<when adding samples, may need to eliminate duplicate points>>
AddListSamples: PROC [from: ListSamples, list: ListSamples, where: Where] = {
IF where = head THEN {
IF from.last.value=list.first.value THEN from.last.next _ list.first.next    --remove duplicate
ELSE from.last.next _ list.first;
list.first _ from.first;
}
ELSE {
IF list.last.value=from.first.value THEN list.last.next _ from.first.next    --remove duplicate
ELSE list.last.next _ from.first;
list.last _ from.last;
};
from.first _ from.last _ NIL;    --keep those reference counts down
};
<<head may be empty>>
AppendList: PROC [head: ListHead, list: ListSamples] = {
IF head.last=NIL THEN head.first _ list    --first time
ELSE head.last.next _ list;
head.last _ list;
};
Empty: SIGNAL=CODE;
RemoveFirst: PROC [head: ListHead] RETURNS [first: ListSamples] = {
IF head.first=NIL THEN SIGNAL Empty;
first _ head.first;
head.first _ head.first.next;
IF head.first=NIL THEN head.last _ NIL;
first.next _ NIL;
RETURN[first];
};
RemoveList: PROC [head: ListHead, list: ListSamples] RETURNS[next: ListSamples] = {
found: BOOLEAN _ FALSE;
IF head.first=NIL THEN SIGNAL Empty;
IF list=head.first THEN {
found _ TRUE;
head.first _ list.next;
IF head.first=NIL THEN head.last_NIL;
}
ELSE {
prev: ListSamples;
FOR prev _ head.first, prev.next UNTIL prev.next=head.last DO
IF prev.next=list THEN {found _ TRUE; prev.next _ list.next; EXIT};
ENDLOOP;
IF list=head.last THEN {found _ TRUE; prev.next _ NIL; head.last _ prev}
};
IF ~found THEN ERROR ELSE {next _ list.next; list.next _ NIL};
RETURN[next];
};
JaMFnsDefs.Register[".setupais"L,JSetUpAIS];    --open the ais file (will close current file)
JaMFnsDefs.Register[".freeais"L,JFreeAIS];    --close the ais file
JaMFnsDefs.Register[".drawmyais"L,DrawAIS];    --use Graphics to display the current file
JaMFnsDefs.Register[".countEdges"L,CountEdges];    
JaMFnsDefs.Register[".countContours"L,CountContours];    
JaMFnsDefs.Register[".drawContour"L,JDrawContour];    
JaMFnsDefs.Register[".tvalue"L,JSetTValue];    --threshhold value.  Make sure it is not xx.0
JaMFnsDefs.Register[".outline"L,Outline];    --find the contours
JaMFnsDefs.Register[".outlineBlackCenter"L,OutlineBlackCenter];    --find the edge through the black region.  Makes sense for 
binary files.
JaMFnsDefs.Register[".outlineBlackEdge"L,OutlineBlackEdge];    --find the edge around the black region.  Makes sense for binary 
files.
JaMFnsDefs.Register[".makecontours"L,MakeContours];    --turn the edges into contours
JaMFnsDefs.Register[".setsa"L,JSetSamples];    --send numbered contour's samples to Curve
JaMFnsDefs.Register[".setcontours"L,SetContours];    --send all contours to Curve
END.