Bilinear Resampling

For my current project I need a bilinear resampling function which is NOT based on Tformfilter and Resizeimage. The goal is to stretch a 64x64 bank to a 256x256 or higher dimensioned bank with bilinear interpolated values. I've setup a quick test with a 4x4 "data image" but there is an error I cant identify. The target size is shown in the green rectangle and you'll notice that the blob is distorted a little bit (it will get even worse if you increase the s# scale variable)

Where is the error? I added some kind of debug logging and the correct results should be:

0, 32, 96,128,128, 96, 32, 0
32, 64,128,160,160,128, 64, 32
96,128,191,223,223,191,128, 96
128,160,223,255,255,223,160,128
128,160,223,255,255,223,160,128
96,128,191,223,223,191,128, 96
32, 64,128,160,160,128, 64, 32
0, 32, 96,128,128, 96, 32, 0

but they are

0, 32, 64, 64, 64, 32, 0, 0
32, 96,160,160,160, 96, 32, 16
64,160,255,255,255,160, 64, 32
64,160,255,255,255,160, 64, 32
64,160,255,255,255,160, 64, 32
32, 96,160,160,160, 96, 32, 16
0, 32, 64, 64, 64, 32, 0, 0
0, 16, 32, 32, 32, 16, 0, 0

Please help, I'm stuck there :-/

Graphics 400,300,32,2

width=4
Height=4
s#=2.0

Restore Image
Dim Heightmap(width,Height)

For Y=0 To width-1
	For X=0 To Height-1
		Read Heightmap(x,y)
	Next
Next

Scale(s,width,Height)
Color 0,128,0 : Rect 16,16,64*s,64*s,0

While Not KeyHit(1)
	
	rgb=(ReadPixel(MouseX(),MouseY(),FrontBuffer()) And $ff0000) / $10000
	
	Color   0,  0,  0 : Rect 0,0,24,16,1
	Color 255,255,255 : Text 0,0,rgb
	
	Flip
	
Wend

End

Function Bilinear#(x#,y#)
	
	Local x0#=Floor(x)
	Local y0#=Floor(y)
	Local dx#=x-x0
	Local dy#=y-y0
	
	Local r0=Heightmap(x0,y0)
	Local r1=Heightmap(x0+1,y0)
	Local r2=Heightmap(x0+1,y0+1)
	Local r3=Heightmap(x0,y0+1)
	
	Return LinearInterpolate(LinearInterpolate(r0,r1,dx),LinearInterpolate(r3,r2,dx),dy)
	
End Function

Function Scale(ratio#,width%,height%)
	
	Local x%,y%,c%
	
	For y=0 To (height*ratio)-1
		
		out$=""
		
		For x=0 To (width*ratio)-1
			
			Local u# = x*1.0/ratio
			Local v# = y*1.0/ratio
			
			c=Bilinear(u,v)
			
			Color c,c,c
			Rect (x+1)*16,(y+1)*16,16,16,1
			
			out$=out$+Fill(c,3)
			If x<(width*ratio)-1 Then out$=out$+","
			
		Next
		
		DebugLog out
		
	Next
	
End Function

Function LinearInterpolate#(x1#,x2#,mu#=0.5)
	
	Return (x1*(1.0-mu)+x2*mu)
	
End Function

Function Fill$(number%,lenght%=2)
	
	Local r$=""
	
	For i=1 To lenght-Len(Str(number))
		
		r$=r$+" "
		
	Next
	
	Return r$+Str(number)
	
End Function

.Image
Data   0, 64, 64,  0
Data  64,255,255, 64
Data  64,255,255, 64
Data   0, 64, 64,  0