// December 19, 1977  11:15 AM


// outgoing procedures:

external [
	InitSpline
	DefineArea
	GetBrush
	DrawSpline
	]
  

// outgoing statics:



// incoming procedures:

external [
	MoveBlock		// SYSTEM
	SetBlock
	Zero
	Allocate
	Free

	ParametricSpline	// SPLINE1 or SPLINE2 or SPLINE3

	FLD; FST; FTR; FLDI	// microFLOAT
	FNEG; FAD; FSB; FML
	FDV; FCM; FSN
	FSTDP; FLDDP
	DPAD; DPSB
	FPSetup

	LoadPackedRAM		// READPRAM
	
	MicroFloatRamImage	// microFLOATMC
	]


// incoming statics:

external [
	sysZone			// SYSTEM

	PSzone			// SPLINE1 or SPLINE2 or SPLINE3

	FPwork			// microFLOAT
	]


// local statics:

static [
	DSnoSpace=0
	// current & next computed points
	DSvec
	DScurX
	DScurY
	DSnewX
	DSnewY
	// clipping limits
	DSminX
	DSmaxX
	DSminY
	DSmaxY
	DSbmh
	// brush height (brush width DSbw is a constant)
	DSbh
	// spline computation parameters
	DScomputeRange=64
	DSnumRate=3
	DSdenRate=2
	DSbbc=0
	]


// local definitions

manifest DSbw=16

structure BBC [
	function word
	spare word
	DBCA word
	DBMR word
	DLX word
	DTY word
	DW word
	DH word
	SBCA word
	SBMR word
	SLX word
	STY word
	gray↑0,3 word
	]

manifest [
	lBBC= size BBC/16

	BBCreplace=0
	BBCpaint=1
	BBCinvert=2
	BBCerase=3

	BBSbitmap=0
	BBSnegBitmap=4
	BBSbrush=8
	BBSgray= 12
	]

structure AREA [
	@BBC
	minX word		// all coordinates relative to bitmap
	maxX word
	minY word
	maxY word
	X0 word			// bitmap coordinate offset
	Y0 word
	bmHeight word
	]

manifest lAREA= size AREA/16

structure DSVEC [
	xTable word
	yTable word
	dTable word
	FPworkSave word
	FPworkNew word
	]

manifest lDSVEC= size DSVEC/16


//	 floating point registers:
manifest [
	// shared by DrawSpline & DrawPoint
	t0=0; t1=1; t2=2; t3=3
	// used by DrawSpline only
	X=8; Y=9
	d1X=10; d2X=11; d3X=12
	d1Y=13; d2Y=14; d3Y=15

	numFPacs=16
	]

let InitSpline(zone; numargs nargs) be [
	PSzone=((nargs ne 1) % (zone eq 0)) ? sysZone, zone
	LoadPackedRAM(MicroFloatRamImage)
	let FPacs=Allocate(PSzone, 32*4+1, lv DSnoSpace)
	if FPacs eq 0 return
	FPacs!0=32
	FPSetup(FPacs)
	]



and DefineArea(bitmap, wordWidth, scanCount, Xw, Yw,
	Xleft, Ybottom, width, height; numargs nargs) = valof [DefineArea
	switchon nargs into [
		case 1:
		case 2:
			resultis 0
		case 3:
			Xw=0
		case 4:
			Yw=0
		case 5:
			Xleft=0
		case 6:
			Ybottom=0
		case 7:
			width=0
		case 8:
			height=0
		]
	if PSzone eq 0 resultis 0
	// even address for BBC block
	let area=Allocate(PSzone, lAREA, lv DSnoSpace, -1)
	if area eq 0 resultis 0
	Zero(area, lAREA)
	
	if width le 0 then width=16*wordWidth
	if height le 0 then height=scanCount
	if Xw ls 0 then Xw=0
	if Yw ls 0 then Yw=0
	let Xright=Xw+width-1
	let Ytop=Yw+height-1
	if Xright ge 16*wordWidth then Xright=16*wordWidth-1
	if Ytop ge height then Ytop=height-1
	area>>AREA.SBMR=1
	area>>AREA.DBCA=bitmap
	area>>AREA.DBMR=wordWidth
	area>>AREA.X0=Xleft-Xw
	area>>AREA.Y0=Ybottom-Yw
	area>>AREA.minX=Xw
	area>>AREA.maxX=Xright
	area>>AREA.minY=Yw
	area>>AREA.maxY=Ytop
	area>>AREA.bmHeight=scanCount-1
	resultis area
	]DefineArea



and DrawSpline(area, n, kXTable, kYTable, brush, drawMode, cyclic; numargs nargs) = valof [DrawSpline
	//returns 0 if it does not work, else |n|

	switchon nargs into [
		case 4:
			brush=0
		case 5:
			drawMode=1		// paint
		case 6:
			cyclic=false
		]

	if (nargs ls 3) % (area eq 0) % (PSzone eq 0) resultis 0
	if brush eq 0 then brush=GetBrush(0, 0)

	let v=vec lDSVEC
	DSvec=v
	Zero(DSvec, lDSVEC)

	let FPworkNew=Allocate(PSzone, 4*numFPacs+1, lv DSnoSpace)
	if FPworkNew eq 0 resultis QuitDrawSpline(0)
	FPworkNew!0=numFPacs
	DSvec>>DSVEC.FPworkNew=FPworkNew
	DSvec>>DSVEC.FPworkSave=FPwork
	FPSetup(FPworkNew)

	let xTable, yTable= 0,0
	let fpX0=vec 2
	let fpY0=vec 2
	let X0=area>>AREA.X0
	let Y0=area>>AREA.Y0
	FST(FLDI(t0, X0), fpX0)
	FST(FLDI(t1, Y0), fpY0)
	test n ls 0
	ifso [
		// convert to floating point & translate to (X0, Y0) origin
		n=-n
		xTable=Allocate(PSzone, 2*n, lv DSnoSpace)
		if xTable eq 0 resultis QuitDrawSpline(0)
		DSvec>>DSVEC.xTable=xTable
		yTable=Allocate(PSzone, 2*n, lv DSnoSpace)
		if yTable eq 0 resultis QuitDrawSpline(0)
		DSvec>>DSVEC.yTable=yTable
		for k=0 to n-1 do [
			FST(FLDI(t0, kXTable!k-X0), xTable+2*k)
			FST(FLDI(t1, kYTable!k-Y0), yTable+2*k)
			]
		]
	ifnot [
		// translate to (X0, Y0) origin
		xTable=kXTable
		yTable=kYTable
		if (X0 ne 0) % (Y0 ne 0) then for k=0 to n-1 do [
			FST(FSB(FLD(t0, xTable+2*k), fpX0), xTable+2*k)
			FST(FSB(FLD(t1, yTable+2*k), fpY0), yTable+2*k)
			]
		]

	DSbh=brush!0
	DSminX=area>>AREA.minX - DSbw
	DSmaxX=area>>AREA.maxX
	DSminY=area>>AREA.minY
	DSmaxY=area>>AREA.maxY + DSbh
	DSbmh=area>>AREA.bmHeight
	DSbbc=area
	DSbbc>>BBC.SBCA=brush+1
	DSbbc>>BBC.function=drawMode & 3

	// is it simply a dot ??
	DSnewX=FTR(FLD(t0, xTable))
	DSnewY=FTR(FLD(t0, yTable))
	if n eq 1 then [
		DrawPoint()
		resultis QuitDrawSpline(1)
		]

	let x2=FTR(FLD(t0, xTable+2))
	let y2=FTR(FLD(t0, yTable+2))
	if DSnewY gr y2 then [ let t=DSnewY; DSnewY=y2; y2=t ]
	if DSnewX gr x2 then [ let t=DSnewX; DSnewX=x2; x2=t ]

	// is it a vertical line ?
	if (n eq 2) & (DSnewX eq x2) then [
		if (DSnewX ls DSminX) % (DSnewX gr DSmaxX) %
			(DSnewY gr DSmaxY) % (y2 ls DSminY) resultis 2
		if DSnewY ls DSminY then DSnewY=DSminY
		if y2 gr DSmaxY then y2=DSmaxY
		for y=DSnewY to y2 do [
			DSnewY=y
			DrawPoint()
			]
		resultis QuitDrawSpline(2)
		]

	// is it a horizontal line ?
	if (n eq 2) & (DSnewY eq y2) then [
		if (DSnewX gr DSmaxX) % (x2 ls DSminX) %
			(DSnewY gr DSmaxY) % (DSnewY ls DSminY) resultis 2
		if DSnewX ls DSminX then DSnewX=DSminX
		if x2 gr DSmaxX then x2=DSmaxX
		for x=DSnewX to x2 do [
			DSnewX=x
			DrawPoint()
			]
		resultis QuitDrawSpline(2)
		]

	// is it a diagonal line ?
	if n eq 2 then [
		DrawPoint()
		DScurX=DSnewX
		DScurY=DSnewY
		DSnewX=x2
		DSnewY=y2
		CheckAndDrawPoint()
		resultis QuitDrawSpline(2)
		]

	// then, it is a curve
	let dn=2*n

	// table for derivatives
	let dTable=Allocate(PSzone, n*12, lv DSnoSpace)
	if dTable eq 0 resultis QuitDrawSpline(0)
	DSvec>>DSVEC.dTable=dTable
	
	// spline computation:
	// ParametricSpline(n,X,Y,X',X'',X''',Y',Y'',Y''')
	let PSdone=ParametricSpline(n, xTable, yTable, dTable, dTable+dn, dTable+2*dn, 
		dTable+3*dn, dTable+4*dn, dTable+5*dn,
		(cyclic ? 1, 0))

	if PSdone eq 0 resultis QuitDrawSpline(0)

	DSnewX=FTR(FLD(t0,xTable))
	DSnewY=FTR(FLD(t1,yTable))
	DScurX=DSnewX
	DScurY=DSnewY
	DrawPoint()

	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
	  // estimate stepping parameters
	  let m1=FTR(FLD(t0,kX))-FTR(FLD(t1,kX+2))
	  if m1 ls 0 then m1=-m1
	  let m2=FTR(FLD(t2,kY))-FTR(FLD(t3,kY+2))
	  if m2 ls 0 then m2=-m2
	  let m=(((m1 gr m2) ? m1, m2)*DSnumRate)/DSdenRate
	  let ni,r=1,m
	  if m gr DScomputeRange then [ r=DScomputeRange; 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+2*k
	  let kX2=kX1+dn
	  let kX3=kX2+dn
	  let kY1=kX3+dn
	  let kY2=kY1+dn
	  let kY3=kY2+dn
	  // 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 [
		DSnewX=DPAD(dpX,dpd1X)
		DSnewY=DPAD(dpY,dpd1Y)
		DPAD(dpd1X,dpd2X); DPAD(dpd1Y,dpd2Y)
		DPAD(dpd2X,dpd3X); DPAD(dpd2Y,dpd3Y)
		CheckAndDrawPoint()
		]
	    // now, next interval or next knot
	    test i eq ni
	    ifso [
		// next knot
		DSnewX=FTR(FLD(t0,kX+2))
		DSnewY=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)
		DSnewX=FTR(FAD(FML(FAD(FML(FAD(FML(FLD(X,kX3),t3),kX2),t2),kX1),t1),kX))
		DSnewY=FTR(FAD(FML(FAD(FML(FAD(FML(FLD(Y,kY3),t3),kY2),t2),kY1),t1),kY))
 		] 
	    CheckAndDrawPoint()
	    ]
	  ]

	resultis QuitDrawSpline(n)
	]DrawSpline



and QuitDrawSpline(r) = valof [
	FPSetup(DSvec>>DSVEC.FPworkSave)
	DSvec>>DSVEC.FPworkSave=0
	for i=0 to lDSVEC-1 do
		if DSvec!i ne 0 then Free(PSzone, DSvec!i)
	resultis r
	]


//*****************************************************************
// point draw/erase procedures
//*****************************************************************


and CheckAndDrawPoint() be [CheckAndDrawPoint
	let deltaX=(DSnewX gr DScurX) ? DSnewX-DScurX, DScurX-DSnewX
	let deltaY=(DSnewY gr DScurY) ? DSnewY-DScurY, DScurY-DSnewY
	// same points ?
	unless deltaX % deltaY return
	// gap ?
	test deltaX gr 1 % deltaY gr 1
	ifso [
		//linear interpolation
		let x=vec 2; let dpdX=vec 2
		let y=vec 2; let dpdY=vec 2
		let m=(deltaX gr deltaY) ? deltaX, deltaY
		FLDI(t0, m)
		FSTDP(FDV(FLDI(t1, DSnewX-DScurX), t0), dpdX)
		FSTDP(FDV(FLDI(t2, DSnewY-DScurY), t0), dpdY)
		x!0=DScurX; x!1=0
		y!0=DScurY; y!1=0
		for i=1 to m do [
			DSnewX=DPAD(x, dpdX)
			DSnewY=DPAD(y, dpdY)
			DrawPoint()
			]
		]
	ifnot DrawPoint()
	DScurX=DSnewX
	DScurY=DSnewY
	]CheckAndDrawPoint




and DrawPoint() be [DrawPoint
	let x=DSnewX-DSbw/2
	let y=DSnewY+DSbh/2
	if (x le DSminX) % (x gr DSmaxX)
	 % (y ls DSminY) % (y ge DSmaxY) return
	test x ls (DSminX+DSbw)
	ifso [
		DSbbc>>BBC.DLX= DSminX+DSbw
		DSbbc>>BBC.DW= x-DSminX
		DSbbc>>BBC.SLX= DSminX+DSbw-x
		]
	ifnot [
		DSbbc>>BBC.DLX= x
		DSbbc>>BBC.DW= ((x gr (DSmaxX-DSbw)) ? (DSmaxX-x+1), DSbw)
		DSbbc>>BBC.SLX= 0
		]
	test y gr (DSmaxY-DSbh)
	ifso [
		DSbbc>>BBC.DTY= DSbmh-DSmaxY+DSbh
		DSbbc>>BBC.STY= y-DSmaxY+DSbh
		DSbbc>>BBC.DH= DSmaxY-y
		]
	ifnot [
		DSbbc>>BBC.DTY= DSbmh-y
		DSbbc>>BBC.STY= 0
		DSbbc>>BBC.DH= (y ls (DSminY+DSbh)) ? y-DSminY+1, DSbh
		]
	DSBitBlt(DSbbc)
	]DrawPoint



and DSBitBlt(bbc) be [
	DSBitBlt= table [
		#055001;		// sta 3,1,2
		#145000;		// mov 2,1
		#111000;		// mov 0,2
		#045001;		// sta 1,1,2
		#126400;		// sub 1,1
		#061024;		// BitBlt
		#031001;		// lda 2,1,2
		#035001;		// lda 3,1,2
		#001401;		//  jmp 1,3
		]
	DSBitBlt(bbc)
	]



and GetBrush(brushShape, brushSize) = valof [
	// returns a pointer to the brush pattern (16 wide, H height)
	// The brush pattern is stored as:
	//			H word
	//			bitPattern↑1,H word
	//
	// Values for brushShape are:
	//		0= round brush
	//		1= rectangular brush
	//		2= horizontal bar brush
	//		3= vertical bar brush
	//		4= diagonal bar brush
	//
	// brushSize is any value between 1 and 16, which will be rounded
	// to the values 1, 2, 4, 8, 16.

	let brushTable= table [
		// dot 1:
		1; #200;
		// dot 2:
		2; #600; #600;
		// dot 4:
		4; #600; #1700; #1700; #600;
		// dot 8:
		8; #1700; #3740; #7760; #7760; #7760; #7760; #3740; #1700;
		// dot 16:
		16; #3740; #17770; #37774; #77776; #77776; #177777; #177777; #177777; 
		#177777; #177777; #177777; #77776; #77776; #37774; #17770; #3740;

		// rect 4:
		4; #1700; #1700; #1700; #1700;
		// rect 8:
		8; #7760; #7760; #7760; #7760; #7760; #7760; #7760; #7760;
		// rect 16:
		16; #177777; #177777; #177777; #177777; #177777; #177777; #177777; #177777; 
		#177777; #177777; #177777; #177777; #177777; #177777; #177777; #177777;

		// hor 2:
		1; #600;
		// hor 4:
		1; #1700;
		// hor 8:
		1; #7760;
		// hor 16:
		1; #17777;

		// ver 2:
		2; #400; #400;
		// ver 4:
		4; #400; #400; #400; #400;
		// ver 8:
		8; #400; #400; #400; #400; #400; #400; #400; #400;
		// ver 16:
		16; #400; #400; #400; #400; #400; #400; #400; #400;
		#400; #400; #400; #400; #400; #400; #400; #400;

		// diag 2:
		2; #200; #400;
		// diag 4;
		4; #100; #200; #400; #1000;
		// diag 8:
		8; #20; #40; #100; #200; #400; #1000; #2000; #4000;
		// diag 16:
		16; 1; 2; 4; #10; #20; #40; #100; #200;
		#400; #1000; #2000; #4000; #10000; #20000; #40000; #100000
		]

	manifest [
		dot1=0; dot2=dot1+2; dot4=dot2+3; dot8=dot4+5; dot16=dot8+9
		rect1=dot1; rect2=dot2; rect4=dot16+17; rect8=rect4+5; rect16=rect8+9
		hor1=dot1; hor2=rect16+17; hor4=hor2+2; hor8=hor4+2; hor16=hor8+2
		ver1=dot1; ver2=hor16+2; ver4=ver2+3; ver8=ver4+5; ver16=ver8+9
		diag1=dot1; diag2=ver16+17; diag4=diag2+3; diag8=diag4+5; diag16=diag8+9
		]

	let brushIndexTable = table [
		dot1; dot2; dot4; dot8; dot16;
		rect1; rect2; rect4; rect8; rect16;
		hor1; hor2; hor4; hor8; hor16;
		ver1; ver2; ver4; ver8; ver16;
		diag1; diag2; diag4; diag8; diag16;
		]

	let widthRoundingTable = table [
		0; 0; 1; 1; 2; 2; 2; 3; 3; 3; 3; 3; 4; 4; 4; 4; 4 ]

	if brushSize gr 16 then brushSize=16
	if (brushShape ls 0) % (brushShape gr 4) then brushShape=0
	resultis (brushTable +
		brushIndexTable!(5*brushShape + widthRoundingTable!brushSize))
	]