// March 31, 1978 2:02 PM *** RESIDENT *** get "zpDefs.bcpl" // outgoing procedures: external [ newSplineID splineType makeSpline remakeSpline computeSpline ] // outgoing statics: external [ @computeRange @numRate @denRate @pseudoCyclic @splineTension ] static [ // spline computation parameters @computeRange=64 @numRate=3 @denRate=2 // flag for pseudo cyclic hack vs periodic splines @pseudoCyclic=false @splineTension=false ] // incoming procedures: external [ MoveBlock // SYSTEM SetBlock Zero Gets typeForm // ZPUTIL giveUp ParametricSpline // PSPLINE FLD; FST; FTR; FLDI // FLOAT FNEG; FAD; FSB; FML FDV; FCM; FSN FSTDP; FLDDP DPAD; DPSB obtainBlock // ZPBLOCK putBlock checkPoint // ZPDRAW startCode endCode checkSplineID curve splineColorSymbol giveMeXY ] // incoming statics: external [ keys // SYSTEM @splineTable // ZPINIT @maxSplineID // globals for encoding: @newX // ZPDRAW @newY ] // local statics: static [ @debugSpline=0 ] // local definitions //***************************************************************** // First: make spline block //***************************************************************** let makeSpline(n, xTable, yTable, brush, color, cyclic; numargs m) = valof [makeSpline // xTable & yTable are FLOATING POINT //returns id of new spline, 0 if error unless n gr 0 resultis 0 if m ls 6 then cyclic=false // check for false curves that actually are lines or dots let minVal=vec 2 let maxVal=vec 2 if n gr 2 & TableConstant(n, yTable) then [ TableMax(n, xTable, maxVal) TableMin(n, xTable, minVal) n=2 FST(FLD(0, minVal), xTable) FST(FLD(0, maxVal), xTable+2) ] if n gr 2 & TableConstant(n, xTable) then [ TableMax(n, yTable, maxVal) TableMin(n, yTable, minVal) n=2 FST(FLD(0, minVal), yTable) FST(FLD(0, maxVal), yTable+2) ] if n eq 2 & TableConstant(2, xTable) & TableConstant(2, yTable) then n=1 // now make the spline let splinePointer=obtainBlock(SPLINEknotBase+4*n) unless splinePointer resultis giveUp("[makeSpline]") Zero(splinePointer, SPLINEknotBase+4*n) //make spline block let knotTable=splinePointer+SPLINEknotBase let dn=2*n MoveBlock(knotTable, xTable, dn) MoveBlock(knotTable+dn, yTable, dn) splinePointer>>SPLINE.nKnots=n splinePointer>>SPLINE.brush=brush splinePointer>>SPLINE.color=color splinePointer>>SPLINE.cyclic=cyclic splineType(splinePointer) giveMeXY(splinePointer, lv splinePointer>>SPLINE.xColor, lv splinePointer>>SPLINE.yColor) let id=remakeSpline(splinePointer) unless id then putBlock(splinePointer) resultis id ]makeSpline and splineType(splinePointer) be [splineType let n=splinePointer>>SPLINE.nKnots let xTable=splinePointer+SPLINEknotBase let yTable=xTable+2*n let x1=FTR(FLD(0, xTable)) let y1=FTR(FLD(0, yTable)) let x2, y2=nil, nil test n eq 1 ifso [ x2=x1; y2=y1 ] ifnot [ x2=FTR(FLD(0, xTable+2)) y2=FTR(FLD(0, yTable+2)) ] splinePointer>>SPLINE.left= (x1 ls x2) ? x1, x2 splinePointer>>SPLINE.right= (x1 gr x2) ? x1, x2 splinePointer>>SPLINE.top= (y1 gr y2) ? y1, y2 splinePointer>>SPLINE.bottom= (y1 ls y2) ? y1, y2 splinePointer>>SPLINE.type = selecton n into [ case 1: dotSpline; case 2: ((y1 eq y2) ? horSpline, ((x1 eq x2) ? verSpline, regSpline)); default: regSpline ] ]splineType //***************************************************************** // Second: get spline id //***************************************************************** and remakeSpline(splinePointer) = valof [remakeSpline let id=newSplineID() unless id resultis 0 test computeSpline(splinePointer) ifso [ splineTable!0=splineTable!0+1 splineTable!id=splinePointer resultis id ] ifnot resultis 0 ]remakeSpline and newSplineID() = valof [newSplineID for id=1 to maxSplineID do unless splineTable!id resultis id typeForm(0, "Sorry, no room for more than ", 10, maxSplineID, 0, "curves*N", 0, "To get a larger work space for curves, start DRAW with switch /S (e.g.: DRAW ", 10, 2*maxSplineID, 0, "/S )*N") resultis 0 ]newSplineID //***************************************************************** // Third: compute [ & generate chain encoding ] //***************************************************************** and computeSpline(splinePointer) = valof [ let r=((splinePointer>>SPLINE.type eq regSpline) ? computeRegularSpline(splinePointer), computeFalseSpline(splinePointer)) splineColorSymbol(splinePointer) resultis r ] and computeFalseSpline(splinePointer) = valof [computeFalseSpline // i.e. single point, or horizontal or vertical line splinePointer>>SPLINE.chain=0 splinePointer>>SPLINE.nBeads=0 splinePointer>>SPLINE.xStart=splinePointer>>SPLINE.left splinePointer>>SPLINE.yStart=splinePointer>>SPLINE.top curve(splinePointer, drawMode) resultis true ]computeFalseSpline and computeRegularSpline(splinePointer) = valof [computeRegularSpline // returns 0 if no storage let procedureName="[computeRegularSpline]" let n=splinePointer>>SPLINE.nKnots let dn=2*n // 4 extra knots for pseudo-cyclic curve let xn=(splinePointer>>SPLINE.cyclic) & pseudoCyclic ? 4, 0 let nt=n+xn let dnt=2*nt // table for derivatives let dTable=obtainBlock(nt*12) unless dTable resultis giveUp(procedureName) let knotTable=splinePointer+SPLINEknotBase let kTable=obtainBlock(2*dnt) unless kTable resultis giveUp(procedureName, dTable) let xTable=kTable+xn let yTable=xTable+dnt MoveBlock(xTable, knotTable, dn) MoveBlock(yTable, knotTable+dn, dn) if xn ne 0 then [ // extra knots MoveBlock(kTable, xTable+dn-6, 4) MoveBlock(xTable+dn, xTable+2, 4) MoveBlock(yTable-4, yTable+dn-6, 4) MoveBlock(yTable+dn, yTable+2, 4) ] manifest [ //local floating Point registers t0=0; t1=1; t2=2; t3=3 //caution: global floating Point registers X=8; Y=9 d1X=10; d2X=11; d3X=12 d1Y=13; d2Y=14; d3Y=15 ] // get storage for chain encoding // (try to estimate block size...) let s=1 let kx=FTR(FLD(t0,xTable)) let ky=FTR(FLD(t1,yTable)) for i=1 to n-1 do [ let kx1=FTR(FLD(t0,xTable+2*i)) let ky1=FTR(FLD(t1,yTable+2*i)) s=s+1+max(abs(kx-kx1),abs(ky-ky1))/16 kx,ky=kx1,ky1 ] let nBeads=(((s*4)/3)/BEADsize)+1 let maxByteCount=16*nBeads let chainSize=nBeads*(BEADsize+2) let chainBlock=obtainBlock(chainSize+maxByteCount/2) unless chainBlock resultis giveUp(procedureName, dTable, (xn ? kTable, 0)) let wTable=0 if splineTension then [ wTable=obtainBlock(nt) if wTable ne 0 then [ SetBlock(wTable, 0, nt) typeForm(0,"Input spline weigths now @", 8, wTable, 1, $*N) Gets(keys) ] ] // spline computation: // caution: ParametricSpline(n,X,Y,X',X'',X''',Y',Y'',Y''') let psDone=ParametricSpline(nt, kTable, kTable+dnt, dTable, dTable+dnt, dTable+2*dnt, dTable+3*dnt, dTable+4*dnt, dTable+5*dnt, ((splinePointer>>SPLINE.cyclic & not pseudoCyclic) ? periodicSpline, naturalSpline), wTable) putBlock(wTable) unless psDone resultis giveUp(procedureName, dTable, chainBlock, (xn ? kTable, 0)) newX=FTR(FLD(t0,xTable)) newY=FTR(FLD(t1,yTable)) splinePointer>>SPLINE.chain=chainBlock chainBlock!(chainSize-1)=maxByteCount splinePointer>>SPLINE.nBeads=nBeads splinePointer>>SPLINE.xStart=newX splinePointer>>SPLINE.yStart=newY startCode(splinePointer,newX,newY) let two=vec 2; FST(FLDI(t0,2),two) let three=vec 2; FST(FLDI(t0,3),three) let six=vec 2; FST(FLDI(t0,6),six) for k=0 to n-2 do [ //knot k let kX,kY=xTable+2*k,yTable+2*k //stepping parameters let m=max(abs(FTR(FLD(t0,kX))-FTR(FLD(t1,kX+2))), abs(FTR(FLD(t2,kY))-FTR(FLD(t3,kY+2)))) m=(m*numRate)/denRate let ni,r=1,m if m gr computeRange then [ r=computeRange; ni=m/r; m=ni*r ] //constants let fni=vec 2; FST(FLDI(t0,ni), fni) test m eq 0 ifso FLDI(t1,0) ifnot FDV(FLDI(t1,1), FLDI(t0,m)) let delta=vec 2; FST(t1,delta) let delta2=vec 2; FST(FML(t1,delta), delta2) let delta3=vec 2; FST(FML(t1,delta), delta3) //derivatives - floating point - let kX1=dTable+xn+2*k let kX2=kX1+dnt let kX3=kX2+dnt let kY1=kX3+dnt let kY2=kY1+dnt let kY3=kY2+dnt // show values if debugSpline then [ typeForm(0, "m=", 10, m, 0, ", delta=", -2, delta, 0, ", X, X', X'', X''' =*N", -2, kX) for i=0 to 2 do typeForm(1, $|, -2, kX1+dnt*i) typeForm(1, $*N, -2, kY) for i=0 to 2 do typeForm(1, $|, -2, kY1+dnt*i) Gets(keys) typeForm(1, $*N) ] //starting derivatives of subintervals - floating point - let X1=vec 2; let X2=vec 2; let X3=vec 2 let Y1=vec 2; let Y2=vec 2; let Y3=vec 2 //start at knot k FLD(X,kX); FLD(Y,kY) X1!0=kX1!0; X1!1=kX1!1 Y1!0=kY1!0; Y1!1=kY1!1 X2!0=kX2!0; X2!1=kX2!1 Y2!0=kY2!0; Y2!1=kY2!1 X3!0=kX3!0; X3!1=kX3!1 Y3!0=kY3!0; Y3!1=kY3!1 //3rd derivatives & differences are constant FML(FLD(d3X,delta3), X3) FML(FLD(d3Y,delta3), Y3) //forward difference computation in subinterval for i=1 to ni do [ //floating point computation of initial values (blahh!) FML(FLD(d2X,delta2),X2) FML(FLD(d1X,delta),X1) FAD(d1X,FDV(FLD(t0,d2X),two)) FAD(d1X,FDV(FLD(t0,d3X),six)) FML(FLD(d2Y,delta2),Y2) FML(FLD(d1Y,delta),Y1) FAD(d1Y,FDV(FLD(t0,d2Y),two)) FAD(d1Y,FDV(FLD(t0,d3Y),six)) FAD(d2X,d3X); FAD(d2Y,d3Y) //double precision fixed point variables let dpX=vec 2; let dpd1X=vec 2; let dpd2X=vec 2; let dpd3X=vec 2 let dpY=vec 2; let dpd1Y=vec 2; let dpd2Y=vec 2; let dpd3Y=vec 2 //initial values FSTDP(X,dpX); FSTDP(d1X,dpd1X); FSTDP(d2X,dpd2X); FSTDP(d3X,dpd3X) FSTDP(Y,dpY); FSTDP(d1Y,dpd1Y); FSTDP(d2Y,dpd2Y); FSTDP(d3Y,dpd3Y) //compute for j=1 to r do [ newX=DPAD(dpX,dpd1X) newY=DPAD(dpY,dpd1Y) DPAD(dpd1X,dpd2X); DPAD(dpd1Y,dpd2Y) DPAD(dpd2X,dpd3X); DPAD(dpd2Y,dpd3Y) checkPoint() ] //now, next interval or next knot test i eq ni ifso [ //next knot newX=FTR(FLD(t0,kX+2)) newY=FTR(FLD(t0,kY+2)) ] ifnot [ //starting point of next interval FDV(FLDI(t1,i),fni) FDV(FLD(t2,t1),two) FDV(FLD(t3,t1),three) FST(FAD(FML(FLD(t0,kX3),t1),kX2),X2) FST(FAD(FML(FLD(t0,kY3),t1),kY2),Y2) FST(FAD(FML(FAD(FML(FLD(t0,kX3),t2),kX2),t1),kX1),X1) FST(FAD(FML(FAD(FML(FLD(t0,kY3),t2),kY2),t1),kY1),Y1) newX=FTR(FAD(FML(FAD(FML(FAD(FML(FLD(X,kX3),t3),kX2),t2),kX1),t1),kX)) newY=FTR(FAD(FML(FAD(FML(FAD(FML(FLD(Y,kY3),t3),kY2),t2),kY1),t1),kY)) ] checkPoint() ] ] //end of chain encoding endCode(splinePointer) //return storage [ chain block is trimmed by endCode() ] putBlock(dTable) putBlock(kTable) resultis true ]computeRegularSpline //***************************************************************** // miscellaneous useful stuff //***************************************************************** and abs(i) = ((i gr 0) ? i, -i) and max(i1,i2) = ((i1 ge i2) ? i1, i2) and TableMax(n, fpTable, fpMax) be [ manifest [ t0=0 ] FST(FLD(t0, fpTable), fpMax) for k=1 to n-1 do if FCM(t0, fpTable+2*k) eq 1 then FST(FLD(t0, fpTable+2*k), fpMax) ] and TableMin(n, fpTable, fpMin) be [ manifest [ t0=0 ] FST(FLD(t0, fpTable), fpMin) for k=1 to n-1 do if FCM(t0, fpTable+2*k) eq -1 then FST(FLD(t0, fpTable+2*k), fpMin) ] and TableConstant(n, fpTable) = valof [ manifest [ t0=0 ] FLD(t0, fpTable) for k=1 to n-1 do if FCM(t0, fpTable+2*k) ne 0 then resultis false resultis true ]