--ContourMain.mesa
--Written by M.Stone
--Last changed by M.Stone July 14, 1982 12:14 pm
--Last changed by M.Plass October 14, 1982 2:59 pm


DIRECTORY
        Rope,
        Runtime,
        Graphics,
        GraphicsOps USING[NewAisImage],
        TJaMGraphics USING[Painter],
        AIS USING [OpenFile, FRef, WRef, OpenWindow,CloseFile, ReadSample, CloseWindow,Error],
        CurveDefs,
        Vector USING [Vec],
        ConvertUnsafe,
        Real,
        RealFns,
        Curve,
        JaMFnsDefs;

ContourMain: PROGRAM IMPORTS  CurveDefs, AIS, JaMFnsDefs, Graphics, TJaMGraphics, Real, 
Runtime, ConvertUnsafe, GraphicsOps, Curve
EXPORTS CurveDefs =
BEGIN OPEN JaMFnsDefs, CurveDefs;

Vec: TYPE = Vector.Vec;

ais: AIS.FRef _ NIL;
window: AIS.WRef _ NIL;
aisname: Rope.ROPE;
allSLists,lastSList: SListHandle _ NIL;
allContours,lastContour: ContourHandle _ NIL;

OpenAIS: PROCEDURE = {
        ENABLE AIS.Error =>{JaMFnsDefs.JaMExec["(AIS file error) .print"]; CONTINUE};
        s: LONG STRING_[64];
         PopString[s];
         ais_AIS.OpenFile[s, FALSE];
        window _ AIS.OpenWindow[ais];
        aisname _ ConvertUnsafe.ToRope[s];
        };

SetTValue: PROCEDURE = {tValue _ PopReal[]};

DrawAIS: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[window.firstPixel,window.lastScan];
          [rx,ty] _ AISToScreen[window.lastPixel,window.firstScan];
            Graphics.Translate[dc,-lx,-by];
          Graphics.ClipBox[dc,[lx,by,rx,ty]];
          Graphics.DrawImage[dc,imageRef];
          Graphics.Restore[dc,mark];
          };
    TJaMGraphics.Painter[Paint];
    };

SetAISWindow:PROCEDURE= {
        by,rx,ty,lx,temp: REAL;
        crx,clx,cty,cby: CARDINAL;
        by _ GetReal[]; 
        rx _ GetReal[];
        ty _ GetReal[];
        lx _ GetReal[];
        IF ty < by THEN {temp _ by; by _ ty; ty _ temp};
        IF rx < lx THEN {temp _ rx; rx _ lx; lx _ temp};
        [lx,by] _ ScreenToAIS[lx,by];
        [rx,ty] _ ScreenToAIS[rx,ty];
        cby _ Real.FixC[by];
        cty _ Real.FixC[ty];
        crx _ Real.FixC[rx];
        clx _ Real.FixC[lx];
        AIS.CloseWindow[window];
        window _ AIS.OpenWindow[f: ais, firstScan: cty,lastScan: cby, firstPixel: clx, lastPixel: crx];
        };

AISToScreen: PROC [x,y: REAL] RETURNS[nx,ny: REAL] ={
        nx _ x; ny _ window.lastScan-y;
        };

ScreenToAIS: PROC [x,y: REAL] RETURNS[nx,ny: REAL] ={
        nx _ x; ny _ window.lastScan-y;
        };

AnalyzeContour: PROCEDURE = {
        count: INTEGER _ 0;
        contour:ContourHandle;
        cn: INTEGER _ PopInteger[];
        countProc: PROC[c: ContourHandle] = {IF cn=count THEN contour _ c; count _ count+1};
        ForAllContours[countProc];
        IF contour#NIL THEN {
        CurveDefs.AnalyzeContour[contour];
        PushInteger[cn];
        SetSamples[];
        };
        };

SetSamples: PROCEDURE = {
        count: INTEGER _ 0;
        contour:ContourHandle;
        cn: INTEGER _ PopInteger[];
        countProc: PROC[c: ContourHandle] = {IF cn=count THEN contour _ c; count _ count+1};
        ForAllContours[countProc];
        IF contour#NIL THEN {
            sa: SampleHandle _ contour.sLists.samples;
            Curve.StartSamples[Curve.defaultHandle, sa.xy.x,sa.xy.y];
            FOR sa _ contour.sLists.samples.next, sa.next UNTIL sa=NIL DO
                Curve.AddSample[Curve.defaultHandle, sa.xy.x,sa.xy.y];
                ENDLOOP;
            };
        };

SetContours: PROCEDURE = {
        fullContour: BOOLEAN _ FALSE;
        SetContour: PROC[contour: ContourHandle] = {
                sa: SampleHandle _ contour.sLists.samples;
            IF fullContour THEN Curve.AddContour[Curve.defaultHandle];
            Curve.StartSamples[Curve.defaultHandle, sa.xy.x,sa.xy.y];
            FOR sa _ contour.sLists.samples.next, sa.next UNTIL sa=NIL DO
                Curve.AddSample[Curve.defaultHandle, sa.xy.x,sa.xy.y];
                        ENDLOOP;
                fullContour _ TRUE;
                };
        Curve.ResetContours[Curve.defaultHandle];
        ForAllContours[SetContour];
        Curve.NextContour[Curve.defaultHandle];
        };

FreeContours: PROC = { allContours _ NIL};
        
Free: PROCEDURE = {
        ENABLE AIS.Error =>{JaMFnsDefs.JaMExec["(AIS file error) .print"]; CONTINUE};
        AIS.CloseWindow[window]; window _ NIL;
        AIS.CloseFile[ais];  ais _ NIL;
        };

LineCount: PROCEDURE = {
        PushInteger[lineCount];
        };

ContourCount: PROCEDURE = {
        PushInteger[contourCount];
        };

RMakeContours: PROCEDURE = {
        MakeContours[];
        DrawContours[];
        };

DrawContours: PROCEDURE = {
    Paint: PROC[dc: Graphics.Context] = {
        drawContour: PROCEDURE[c: ContourHandle] = {
          Graphics.MoveTo[path, 0, 0]; -- flush the old path.
          EnterContour[path, c];
          Graphics.DrawStroke[dc, path];
          };
        ForAllContours[drawContour];
        };
    path: Graphics.Path _ Graphics.NewPath[];
    TJaMGraphics.Painter[Paint];
    };

FillContours: PROCEDURE = {
    Paint: PROC[dc: Graphics.Context] = {
        fillContour: PROCEDURE[c: ContourHandle] = {
          EnterContour[path, c];
          };
        ForAllContours[fillContour];
        Graphics.DrawArea[dc, path];
        };
    path: Graphics.Path _ Graphics.NewPath[];
    TJaMGraphics.Painter[Paint];
    };

DrawSLists: PROCEDURE = {
        ForAllSLists[DrawSList];
        };

DrawSList: PROCEDURE[sl: SListHandle] = {
    path: Graphics.Path _ Graphics.NewPath[];
    Paint: PROC[dc: Graphics.Context] = {Graphics.DrawStroke[dc,path]};
    Graphics.MoveTo[path, sl.samples.xy.x, sl.samples.xy.y];
    LineToSamples[path,sl];
    TJaMGraphics.Painter[Paint];
    };

EnterContour: PROCEDURE[path: Graphics.Path, c: ContourHandle] = {
        sl: SListHandle _ c.sLists;
        Graphics.MoveTo[path, sl.samples.xy.x, sl.samples.xy.y, FALSE];
        DO LineToSamples[path,sl]; IF sl=c.lastSList THEN EXIT ELSE sl _ sl.next; ENDLOOP;
        };


LineToSamples: PROCEDURE[path:Graphics.Path, sl: SListHandle] = {
        s: SampleHandle _ sl.samples.next;
        DO
                Graphics.LineTo[path, s.xy.x, s.xy.y];
                IF s=sl.lastSample THEN EXIT ELSE s _ s.next;
                ENDLOOP;
        };

--the trajectory handling stuff
trajList,lastTraj: TrajHandle _ NIL;

AddTraj: PROC = {
        traj: TrajHandle _ NIL ;-- CurrentTraj went away (MFP) Curve.CurrentTraj[];
        IF traj=NIL THEN RETURN;
        IF lastTraj=NIL THEN trajList _ lastTraj _ traj
        ELSE {
            lastTraj.next _ traj;
            traj.prev _ lastTraj;
            lastTraj _ lastTraj.next;
            };
        };

FreeTrajList: PROC = {trajList _ lastTraj _ NIL};

lineCount: INTEGER _ 0;
contourCount: INTEGER _ 0;
Outline: PROC = { 
        v1,v2,v3,v4,x,y: CARDINAL;
        lineCount _ 0;
        contourCount _ 0;
        InitContours[];
        FOR y IN [window.firstScan..window.lastScan) DO
                IF JaMFnsDefs.GetJaMBreak[] THEN EXIT;
--this assumes the rd scan direction for ais files
                FOR x IN [window.firstPixel..window.lastPixel) DO
                        v1 _ AIS.ReadSample[window,y,x];
                        v2 _ AIS.ReadSample[window,y,x+1];
                        v3 _ AIS.ReadSample[window,y+1,x+1];
                        v4 _ AIS.ReadSample[window,y+1,x];
                        ContourCell[v1,v2,v3,v4,x,y];
                        ENDLOOP;
                ENDLOOP;
        };

Trans: TYPE = RECORD[x,y,d: REAL, dir:{nt,dl,ld}];

--aisX and aisY relative to current AIS window
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 Runtime.CallDebugger["odd transitions"];
        IF dlCount>2 THEN Runtime.CallDebugger["too many transitions"];
        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];
        NewSample[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];
                NewSample[edgePts[move].x,edgePts[move].y,
                        edgePts[index].x,edgePts[index].y];
                lineCount _ lineCount+1;
                END;

        };

tValue: REAL _ 128.5;
nPoints: INTEGER _ 4;
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;
        };

InitContours: PROC = {
        allSLists _ lastSList _ NIL;
        allContours _ lastContour _ NIL;
        };

--There exist a set of SLists  Connect these edges into closed contours
--holes must have opposite wrap than areas.
MakeContours: PROC = {
        contour: ContourHandle _ NewContour[FirstSList[]];
        found: BOOLEAN _ FALSE;
        close: PROC[c: ContourHandle] = { 
                IF c.sLists.samples.xy=c.sLists.lastSample.xy THEN {
                    c.sLists.lastSample _ c.sLists.lastSample.prev;
                    c.sLists.closed _ TRUE;
                    }
                };
        find: PROC[sl: SListHandle] = {
                IF sl.samples.xy=contour.lastSList.lastSample.xy THEN {
                  AddSamplesToContour[contour,sl];
                  found _ TRUE};
                  };
        
        UNTIL FirstSList[]=NIL DO
                found _ FALSE;
                ForAllSLists[find];        --finds pieces for contour
                IF ~found THEN contour _ NewContour[FirstSList[]];
                ENDLOOP;
        ForAllContours[close];
        };

--list manipulating routines

--these are for the stand alone set of SLists which are created as the image is outlined
NewSample: PROC[x0,y0,x1,y1: REAL] = {OPEN Real;
    Paint: PROC[dc: Graphics.Context] = {
      [p0.x,p0.y] _ AISToScreen[x0,y0];
      [p1.x,p1.y] _ AISToScreen[x1,y1];
      Graphics.SetCP[dc,p0.x,p0.y];
      Graphics.DrawTo[dc,p1.x,p1.y];
      };
        p0,p1: Vec;
        slist: SListHandle _ NIL;
        sample: SampleHandle _ NEW[Sample];
        found: BOOLEAN _ FALSE;
        find: PROC[sl: SListHandle] = {
          IF found THEN RETURN;
          IF sl.lastSample.xy=p0 THEN
              {sample.xy _ p1; AddSample[sl,sample,TRUE]; found _ TRUE}  --add to tail of sl
          ELSE IF sl.samples.xy=p1 THEN
              {AddSample[sl,sample,FALSE]; found _ TRUE};  --add to head of sl
          };
        TJaMGraphics.Painter[Paint];
        sample.xy _ p0;
        ForAllSLists[find];
        IF ~found THEN {
          slist _ NewSList[sample];
          sample _ NEW[Sample];
          sample.xy _ p1;
          AddSample[slist,sample]}; 
        };

NewSList: PROC[sample: SampleHandle] RETURNS[SListHandle] = {
        slist: SListHandle _ NEW[SList];
        slist.samples _ slist.lastSample _ sample;
        IF allSLists=NIL THEN allSLists _ lastSList _ slist
        ELSE {        slist.prev _ lastSList; slist.next _ NIL;
                lastSList.next _ slist; lastSList _ slist};
        RETURN[slist];
        };

--atTail determines which end of the slist the sample will be added
AddSample: PROC [sl: SListHandle, sample: SampleHandle,atTail: BOOLEAN _ TRUE]= {
        IF atTail THEN {
          sample.prev _ sl.lastSample; sample.next _ NIL;
          sl.lastSample.next _ sample; sl.lastSample _ sample}
        ELSE {
          sample.next _ sl.samples; sample.prev _ NIL;
          sl.samples.prev _ sample; sl.samples _ sample};
        };

FirstSList: PROC RETURNS[SListHandle] = {
        RETURN[allSLists];
        };

ForAllSLists: PROC[proc:PROC[sl: SListHandle]] = {
        ptr: SListHandle _ allSLists;
        next: SListHandle _ NIL;
        UNTIL ptr=NIL DO next _ ptr.next; proc[ptr]; ptr _ next; ENDLOOP;
        };

RemoveSList: PROC [sl: SListHandle]= {
        IF sl.next#NIL THEN sl.next.prev _ sl.prev;
        IF sl.prev#NIL THEN sl.prev.next _ sl.next;
        IF sl=allSLists THEN allSLists _ sl.next;
        IF sl=lastSList THEN lastSList _ sl.prev;
        sl.next _ sl.prev _ NIL;
        };

--these are for the stand alone set of contours created from the SLists by MakeContours
NewContour: PROC[sl: SListHandle] RETURNS[ContourHandle _ NIL] = {
        IF sl # NIL THEN {
                contour: ContourHandle _ NEW[Contour];
                RemoveSList[sl];
                contour.sLists _ contour.lastSList _ sl;
                AddContour[contour];
                RETURN[contour];
                };
        };

AddContour: PROC[contour: ContourHandle]= {
        IF allContours=NIL THEN allContours _ lastContour _ contour
        ELSE {        contour.prev _ lastContour; contour.next _ NIL;
                lastContour.next _ contour; lastContour _ contour};
        };

DeleteContour: PROC [contour: ContourHandle]= {
        IF contour.next#NIL THEN contour.next.prev _ contour.prev;
        IF contour.prev#NIL THEN contour.prev.next _ contour.next;
        contour.next _ contour.prev _ NIL;
        };

AddSamplesToContour: PROC [contour: ContourHandle, sl: SListHandle]= {
        c: SListHandle _ contour.sLists;
        RemoveSList[sl];
        c.lastSample _ c.lastSample.prev;  --skip duplicate point
        IF c.lastSample.next#NIL
          THEN c.lastSample.next _ sl.samples
          ELSE c.samples _ sl.samples;  --single sample sLists
        sl.samples.prev _ c.lastSample;
        c.lastSample _ sl.lastSample;
        };

ForAllContours: PROC[proc:PROC[c: ContourHandle]] = {
        ptr: ContourHandle _ allContours;
        next: ContourHandle _ NIL;
        UNTIL ptr=NIL DO next _ ptr.next; proc[ptr]; ptr _ next; ENDLOOP;
        };

--these are the calls to the output routines to write out the data
OpenOutputFile: PROC = {
        s: STRING _ [40];
        PopString[s];
        CurveDefs.OpenFile[ConvertUnsafe.ToRope[s]];
        };
        
CloseOutputFile: PROC = {
        CurveDefs.CloseFile[];
        };
        
RWriteContours: PROC = {
        s: REAL;
        write: PROC[c: ContourHandle] = {CurveDefs.WriteContour[c,1/s,[-minY,-minX]]};
        minX _ 10000; minY _ 10000;
        maxX _ 0; maxY _ 0;
        ForAllContours[SetMinMax];
        s _ maxY-minY;
        IF maxX-minX > s THEN s _ maxX-minX;
        ForAllContours[write];
        };

RWriteTrajs: PROC = {
        FOR traj: TrajHandle _ trajList, traj.next UNTIL traj=NIL DO CurveDefs.WriteTraj[traj,TRUE]; ENDLOOP;
        };

minX,maxX,minY,maxY: REAL _ 0;
SetMinMax: PROC [c: ContourHandle] = {
        compare: PROC [sa: SampleHandle] = {
                IF sa.xy.x>maxX THEN maxX _ sa.xy.x;
                IF sa.xy.y>maxY THEN maxY _ sa.xy.y;
                IF sa.xy.x<minX THEN minX _ sa.xy.x;
                IF sa.xy.y<minY THEN minY _ sa.xy.y;
                };
        ForAllSamples[c,compare];
        };

ForAllSamples: PROC [c: ContourHandle,proc: PROC[sa: SampleHandle]] = {
        sl: SListHandle _ c.sLists;
        sa: SampleHandle _ NIL;
        UNTIL sl=NIL DO
                sa _ sl.samples;
                UNTIL sa=NIL DO
                        proc[sa];
                        IF sa=sl.lastSample THEN EXIT ELSE sa _ sa.next;
                        ENDLOOP;
                IF sl=c.lastSList THEN EXIT ELSE sl _ sl.next;
                ENDLOOP
        };

Register[".setupais"L,OpenAIS];
Register[".freeais"L,Free];
Register[".tvalue"L,SetTValue];
Register[".outline"L,Outline];
Register[".drawmyais"L,DrawAIS];
--Register[".window"L,SetAISWindow];
--Register[".linecount"L,LineCount];
--Register[".contourcount"L,ContourCount];
Register[".makecontours"L,RMakeContours];
--Register[".freecontours"L,FreeContours];

--Register[".drawcontours"L,DrawContours];
--Register[".fillcontours"L,FillContours];
--Register[".drawslists"L,DrawSLists];
Register[".analyzecontour"L,AnalyzeContour];
--Register[".writecontours"L,RWriteContours];
--Register[".writetrajs"L,RWriteTrajs];
--Register[".outputfile"L,OpenOutputFile];
--Register[".closeout"L,CloseOutputFile];

Register[".addtraj"L,AddTraj];
Register[".freetrajlist"L,FreeTrajList];
Register[".setsa"L,SetSamples];
Register[".setcontours"L,SetContours];

        
END.

M.Plass August 13, 1982 10:38 am: Added SetContours.