add various math.atan() routines

This commit is contained in:
Irmen de Jong 2023-06-17 00:43:23 +02:00
parent 3de10adac2
commit c0b398e0ce
7 changed files with 434 additions and 14 deletions

View File

@ -53,7 +53,8 @@ class VarConstantValueTypeAdjuster(private val program: Program, private val err
}
if(to==null) {
if(!range.to.inferType(program).isInteger)
val toType = range.to.inferType(program)
if(toType.isKnown && !range.to.inferType(program).isInteger)
errors.err("range expression to value must be integer", range.to.position)
} else if(to-to.toInt()>0) {
errors.err("range expression to value must be integer", range.to.position)

View File

@ -95,4 +95,293 @@ _sinecosR8 .char trunc(127.0 * sin(range(180+45) * rad(360.0/180.0)))
rts
}}
}
sub atan_coarse_sgn(byte x1, byte y1, byte x2, byte y2) -> ubyte {
; From a pair of signed coordinates around the origin, calculate discrete direction between 0 and 23 into A.
cx16.r0L = 3 ; quadrant
cx16.r1sL = x2-x1 ; xdelta
if_neg {
cx16.r0L--
cx16.r1sL = -cx16.r1sL
}
cx16.r2sL = y2-y1 ; ydelta
if_neg {
cx16.r0L-=2
cx16.r2sL = -cx16.r2sL
}
return atan_coarse_qd(cx16.r0L, cx16.r1L, cx16.r2L)
}
sub atan_coarse(ubyte x1, ubyte y1, ubyte x2, ubyte y2) -> ubyte {
; From a pair of positive coordinates, calculate discrete direction between 0 and 23 into A.
cx16.r0L = 3 ; quadrant
if x2>=x1 {
cx16.r1L = x2-x1
} else {
cx16.r1L = x1-x2
cx16.r0L--
}
if y2>=y1 {
cx16.r2L = y2-y1
} else {
cx16.r2L = y1-y2
cx16.r0L -= 2
}
return atan_coarse_qd(cx16.r0L, cx16.r1L, cx16.r2L)
}
asmsub atan_coarse_qd(ubyte quadrant @A, ubyte xdelta @X, ubyte ydelta @Y) -> ubyte @A {
;Arctan https://github.com/dustmop/arctan24
; From a pair of X/Y deltas (both >=0), and quadrant 0-3, calculate discrete direction between 0 and 23 into A.
; .reg:a @in quadrant Number 0 to 3.
; .reg:x @in x_delta Delta for x direction.
; .reg:y @in y_delta Delta for y direction.
; Returns A as the direction (0-23).
%asm {{
x_delta = cx16.r0L
y_delta = cx16.r1L
quadrant = cx16.r2L
half_value = cx16.r3L
region_number = cx16.r4L
small = cx16.r5L
large = cx16.r5H
sta quadrant
sty y_delta
stx x_delta
cpx y_delta
bcs _XGreaterOrEqualY
_XLessY:
lda #16
sta region_number
stx small
sty large
bne _DetermineRegion
_XGreaterOrEqualY:
lda #0
sta region_number
stx large
sty small
_DetermineRegion:
; set A = small * 2.5
lda small
lsr a
sta half_value
lda small
asl a
bcs _SmallerQuotient
clc
adc half_value
bcs _SmallerQuotient
cmp large
bcc _LargerQuotient
; S * 2.5 > L
_SmallerQuotient:
; set A = S * 1.25
lsr half_value
lda small
clc
adc half_value
cmp large
bcc _Region1 ; if S * 1.25 < L then goto Region1 (L / S > 1.25)
bcs _Region0 ; (L / S < 1.25)
; S * 2.5 < L
_LargerQuotient:
; set A = S * 7.5
lda small
asl a
asl a
asl a
bcs _Region2
sec
sbc half_value
cmp large
bcc _Region3 ; if S * 7.5 < L then goto Region3 (L / S > 7.5)
jmp _Region2 ; (L / S < 7.5)
_Region0:
; L / S < 1.25. d=3,9,15,21
jmp _LookupResult
_Region1:
; 1.25 < L / S < 2.5. d=2,4,8,10,14,16,20,22
lda region_number
clc
adc #4
sta region_number
bpl _LookupResult
_Region2:
; 2.5 < L / S < 7.5. d=1,5,7,11,13,17,19,23
lda region_number
clc
adc #8
sta region_number
bpl _LookupResult
_Region3:
; 7.5 < L / S. d=0,6,12,18
lda region_number
clc
adc #12
sta region_number
_LookupResult:
lda quadrant
clc
adc region_number
tax
lda _quadrant_region_to_direction,x
rts
_quadrant_region_to_direction:
.byte 9, 3,15,21
.byte 10, 2,14,22
.byte 11, 1,13,23
.byte 12, 0,12, 0
.byte 9, 3,15,21
.byte 8, 4,16,20
.byte 7, 5,17,19
.byte 6, 6,18,18
}}
}
asmsub atan(ubyte x1 @R0, ubyte y1 @R1, ubyte x2 @R2, ubyte y2 @R3) -> ubyte @A {
;; Calculate the angle, in a 256-degree circle, between two points into A.
;; The points (x1, y1) and (x2, y2) have to use *unsigned coordinates only* from the positive quadrant in the carthesian plane!
;; https://www.codebase64.org/doku.php?id=base:8bit_atan2_8-bit_angle
;; This uses 2 large lookup tables so uses a lot of memory but is super fast.
%asm {{
x1 = cx16.r0L
y1 = cx16.r1L
x2 = cx16.r2L
y2 = cx16.r3L
octant = cx16.r4L ;; temporary zeropage variable
lda x1
sbc x2
bcs *+4
eor #$ff
tax
rol octant
lda y1
sbc y2
bcs *+4
eor #$ff
tay
rol octant
lda log2_tab,x
sbc log2_tab,y
bcc *+4
eor #$ff
tax
lda octant
rol a
and #%111
tay
lda atan_tab,x
eor octant_adjust,y
rts
octant_adjust
.byte %00111111 ;; x+,y+,|x|>|y|
.byte %00000000 ;; x+,y+,|x|<|y|
.byte %11000000 ;; x+,y-,|x|>|y|
.byte %11111111 ;; x+,y-,|x|<|y|
.byte %01000000 ;; x-,y+,|x|>|y|
.byte %01111111 ;; x-,y+,|x|<|y|
.byte %10111111 ;; x-,y-,|x|>|y|
.byte %10000000 ;; x-,y-,|x|<|y|
;;;;;;;; atan(2^(x/32))*128/pi ;;;;;;;;
atan_tab
.byte $00,$00,$00,$00,$00,$00,$00,$00
.byte $00,$00,$00,$00,$00,$00,$00,$00
.byte $00,$00,$00,$00,$00,$00,$00,$00
.byte $00,$00,$00,$00,$00,$00,$00,$00
.byte $00,$00,$00,$00,$00,$00,$00,$00
.byte $00,$00,$00,$00,$00,$00,$00,$00
.byte $00,$00,$00,$00,$00,$00,$00,$00
.byte $00,$00,$00,$00,$00,$00,$00,$00
.byte $00,$00,$00,$00,$00,$00,$00,$00
.byte $00,$00,$00,$00,$00,$00,$00,$00
.byte $00,$00,$00,$00,$00,$01,$01,$01
.byte $01,$01,$01,$01,$01,$01,$01,$01
.byte $01,$01,$01,$01,$01,$01,$01,$01
.byte $01,$01,$01,$01,$01,$01,$01,$01
.byte $01,$01,$01,$01,$01,$02,$02,$02
.byte $02,$02,$02,$02,$02,$02,$02,$02
.byte $02,$02,$02,$02,$02,$02,$02,$02
.byte $03,$03,$03,$03,$03,$03,$03,$03
.byte $03,$03,$03,$03,$03,$04,$04,$04
.byte $04,$04,$04,$04,$04,$04,$04,$04
.byte $05,$05,$05,$05,$05,$05,$05,$05
.byte $06,$06,$06,$06,$06,$06,$06,$06
.byte $07,$07,$07,$07,$07,$07,$08,$08
.byte $08,$08,$08,$08,$09,$09,$09,$09
.byte $09,$0a,$0a,$0a,$0a,$0b,$0b,$0b
.byte $0b,$0c,$0c,$0c,$0c,$0d,$0d,$0d
.byte $0d,$0e,$0e,$0e,$0e,$0f,$0f,$0f
.byte $10,$10,$10,$11,$11,$11,$12,$12
.byte $12,$13,$13,$13,$14,$14,$15,$15
.byte $15,$16,$16,$17,$17,$17,$18,$18
.byte $19,$19,$19,$1a,$1a,$1b,$1b,$1c
.byte $1c,$1c,$1d,$1d,$1e,$1e,$1f,$1f
;;;;;;;; log2(x)*32 ;;;;;;;;
log2_tab
.byte $00,$00,$20,$32,$40,$4a,$52,$59
.byte $60,$65,$6a,$6e,$72,$76,$79,$7d
.byte $80,$82,$85,$87,$8a,$8c,$8e,$90
.byte $92,$94,$96,$98,$99,$9b,$9d,$9e
.byte $a0,$a1,$a2,$a4,$a5,$a6,$a7,$a9
.byte $aa,$ab,$ac,$ad,$ae,$af,$b0,$b1
.byte $b2,$b3,$b4,$b5,$b6,$b7,$b8,$b9
.byte $b9,$ba,$bb,$bc,$bd,$bd,$be,$bf
.byte $c0,$c0,$c1,$c2,$c2,$c3,$c4,$c4
.byte $c5,$c6,$c6,$c7,$c7,$c8,$c9,$c9
.byte $ca,$ca,$cb,$cc,$cc,$cd,$cd,$ce
.byte $ce,$cf,$cf,$d0,$d0,$d1,$d1,$d2
.byte $d2,$d3,$d3,$d4,$d4,$d5,$d5,$d5
.byte $d6,$d6,$d7,$d7,$d8,$d8,$d9,$d9
.byte $d9,$da,$da,$db,$db,$db,$dc,$dc
.byte $dd,$dd,$dd,$de,$de,$de,$df,$df
.byte $df,$e0,$e0,$e1,$e1,$e1,$e2,$e2
.byte $e2,$e3,$e3,$e3,$e4,$e4,$e4,$e5
.byte $e5,$e5,$e6,$e6,$e6,$e7,$e7,$e7
.byte $e7,$e8,$e8,$e8,$e9,$e9,$e9,$ea
.byte $ea,$ea,$ea,$eb,$eb,$eb,$ec,$ec
.byte $ec,$ec,$ed,$ed,$ed,$ed,$ee,$ee
.byte $ee,$ee,$ef,$ef,$ef,$ef,$f0,$f0
.byte $f0,$f1,$f1,$f1,$f1,$f1,$f2,$f2
.byte $f2,$f2,$f3,$f3,$f3,$f3,$f4,$f4
.byte $f4,$f4,$f5,$f5,$f5,$f5,$f5,$f6
.byte $f6,$f6,$f6,$f7,$f7,$f7,$f7,$f7
.byte $f8,$f8,$f8,$f8,$f9,$f9,$f9,$f9
.byte $f9,$fa,$fa,$fa,$fa,$fa,$fb,$fb
.byte $fb,$fb,$fb,$fc,$fc,$fc,$fc,$fc
.byte $fd,$fd,$fd,$fd,$fd,$fd,$fe,$fe
.byte $fe,$fe,$fe,$ff,$ff,$ff,$ff,$ff
}}
}
}

View File

@ -182,4 +182,58 @@ math {
return
}}
}
sub atan_coarse_sgn(byte x1, byte y1, byte x2, byte y2) -> ubyte {
; From a pair of signed coordinates around the origin, calculate discrete direction between 0 and 23 into A.
cx16.r0L = 3 ; quadrant
cx16.r1sL = x2-x1 ; xdelta
if_neg {
cx16.r0L--
cx16.r1sL = -cx16.r1sL
}
cx16.r2sL = y2-y1 ; ydelta
if_neg {
cx16.r0L-=2
cx16.r2sL = -cx16.r2sL
}
return atan_coarse_qd(cx16.r0L, cx16.r1L, cx16.r2L)
}
sub atan_coarse(ubyte x1, ubyte y1, ubyte x2, ubyte y2) -> ubyte {
; From a pair of positive coordinates, calculate discrete direction between 0 and 23 into A.
cx16.r0L = 3 ; quadrant
if x2>=x1 {
cx16.r1L = x2-x1
} else {
cx16.r1L = x1-x2
cx16.r0L--
}
if y2>=y1 {
cx16.r2L = y2-y1
} else {
cx16.r2L = y1-y2
cx16.r0L -= 2
}
return atan_coarse_qd(cx16.r0L, cx16.r1L, cx16.r2L)
}
sub atan_coarse_qd(ubyte quadrant, ubyte xdelta, ubyte ydelta) -> ubyte {
; From a pair of X/Y deltas (both >=0), and quadrant 0-3, calculate discrete direction between 0 and 23.
return lsb(mkword(atan(0, 0, xdelta, ydelta), 0) / 2730)
}
sub atan(ubyte x1, ubyte y1, ubyte x2, ubyte y2) -> ubyte {
;; Calculate the angle, in a 256-degree circle, between two points into A.
;; The points (x1, y1) and (x2, y2) have to use *unsigned coordinates only* from the positive quadrant in the carthesian plane!
%ir {{
loadm.b r65532,math.atan.x1
loadm.b r65533,math.atan.y1
loadm.b r65534,math.atan.x2
loadm.b r65535,math.atan.y2
syscall 44 (r65532.b, r65533.b, r65534.b, r65535.b): r0.b
returnr.b r0
}}
}
}

View File

@ -372,6 +372,23 @@ but perhaps the provided ones can be of service too.
Fast 8-bit byte cosine of angle 0..179 (each is a 2 degree step), result is in range -127..127
Angles 180..255 will yield a garbage result!
``atan(ubyte x1, ubyte y1, ubyte x2, ubyte y2)``
Fast arctan routine that uses more memory because of large lookup tables.
Calculate the angle, in a 256-degree circle, between two points in the positive quadrant.
``atan_coarse_sgn(byte x1, byte y1, byte x2, byte y2)``
Small and fast, but imprecise, arctan routine
From a pair of signed coordinates around the origin, calculate discrete direction between 0 and 23.
``atan_coarse(ubyte x1, ubyte y1, ubyte x2, ubyte y2)``
Small and fast, but imprecise, arctan routine
From a pair of positive coordinates, calculate discrete direction between 0 and 23.
``atan_coarse_qd(ubyte quadrant, ubyte xdelta, ubyte ydelta)``
Small and fast, but imprecise, arctan routine
If you already know the quadrant and x/y deltas, calculate discrete direction between 0 and 23.
cx16logo
--------

View File

@ -1,6 +1,12 @@
TODO
====
- this should give a compiler error because word returnvalue:
sub atan_coarse_qd(ubyte quadrant, ubyte xdelta, ubyte ydelta) -> ubyte {
return mkword(math.atan(0, 0, xdelta, ydelta), 0) / 2730
}
- vm: fix syscall.ATAN calculation
- document some library modules better (diskio, etc)
...

View File

@ -1,17 +1,44 @@
%import math
%import textio
%import diskio
%zeropage basicsafe
main {
sub start() {
txt.print("pwd: ")
txt.print(diskio.curdir())
txt.print("\ndisk name: ")
uword name = diskio.diskname()
if name
txt.print(name)
else
txt.print("!error!")
const ubyte HEIGHT = 30 ; txt.DEFAULT_HEIGHT-1
const ubyte WIDTH = 80 ; txt.DEFAULT_WIDTH-1
const ubyte HALFWIDTH = 40 ; txt.DEFAULT_WIDTH/2
const ubyte HALFHEIGHT = 15 ; txt.DEFAULT_HEIGHT/2
txt.print_ub(math.atan(0, 0, 10, 20))
; ubyte @zp value
; ubyte xx
; ubyte yy
; for yy in 0 to HEIGHT {
; for xx in 0 to WIDTH {
; value = math.atan(HALFWIDTH, HALFHEIGHT, xx, yy)
; txt.setchr(xx,yy,value)
; }
; }
;
; byte sx
; byte sy
; for sy in 0 to HEIGHT as byte {
; for sx in 0 to WIDTH as byte {
; value = math.atan_coarse_sgn(0, 0, sx-HALFWIDTH, sy-HALFHEIGHT)
; txt.setchr(sx as ubyte,sy as ubyte,value)
; }
; }
;
; for yy in 0 to HEIGHT {
; for xx in 0 to WIDTH {
; value = math.atan_coarse(HALFWIDTH, HALFHEIGHT, xx, yy)
; txt.setchr(xx,yy,value)
; }
; }
;
; goto start
}
}

View File

@ -42,6 +42,18 @@ SYSCALLS:
30 = gfx_getpixel ; get byte pixel value at coordinates r0.w/r1.w
31 = rndseed
32 = rndfseed
33 = RND
34 = RNDW
35 = RNDF
36 = STRING_CONTAINS
37 = BYTEARRAY_CONTAINS
38 = WORDARRAY_CONTAINS
39 = CLAMP_BYTE
40 = CLAMP_UBYTE
41 = CLAMP_WORD
42 = CLAMP_UWORD
43 = CLAMP_FLOAT
43 = ATAN
*/
enum class Syscall {
@ -88,7 +100,8 @@ enum class Syscall {
CLAMP_UBYTE,
CLAMP_WORD,
CLAMP_UWORD,
CLAMP_FLOAT
CLAMP_FLOAT,
ATAN
;
companion object {
@ -454,12 +467,25 @@ object SysCalls {
}
Syscall.CLAMP_FLOAT -> {
val (valueU, minimumU, maximumU) = getArgValues(callspec.arguments, vm)
val value = (valueU as Float)
val minimum = (minimumU as Float)
val maximum = (maximumU as Float)
val value = valueU as Float
val minimum = minimumU as Float
val maximum = maximumU as Float
val result = min(max(value, minimum), maximum)
returnValue(callspec.returns!!, result, vm)
}
Syscall.ATAN -> {
val (x1, y1, x2, y2) = getArgValues(callspec.arguments, vm)
val x1f = (x1 as UByte).toDouble()
val y1f = (y1 as UByte).toDouble()
val x2f = (x2 as UByte).toDouble()
val y2f = (y2 as UByte).toDouble()
val xd = x2f-x1f
val yd = y2f-y1f
TODO("calculate atan the same way as the 6502 routine does 0-255")
// val radians = atan2(yd, xd) + PI // 0 to 2*PI
// val result = floor(2*PI/radians*256.0)
// returnValue(callspec.returns!!, result, vm)
}
else -> throw AssemblyError("missing syscall ${call.name}")
}
}