bits 16

org 100h

; shorter but has artifacts due to paging / banking (small black squares)
;les bp,[bx]
; get VRAM addr (compat + remove artifacts)
push 0xa000
pop es

fninit ; required on bochs for clean output

mov di,fp_data + 2040
;mov cx,0xff ; compat
fldz ; compat (don't works in bochs / freedos otherwise)
init:
    fld1
    faddp
    fst dword [di]
    ; compat (bochs / freedos)
    fldz
    fstp dword [di + 4]
    sub di,8
    loop init

mov ax,$4f02
mov bx,$4213 ; 640x400 32 bpp
;mov bx,$4141 ; BOCHS
;mov bx,$4115 ; 800x600x4
;mov bx,$4143
int 10h

; upper ECX init to prevent frame to be 0
dec cx
inc ecx
lp:
    mov dword [di - 8],ecx
    mov al,8
    mul cl
    xchg bx,ax

    fld dword [di + bx] ; x[i]
    fild word [di - 6] ; frame
    fld dword [di + bx + 4] ; y[i]
    fld st0 ; dup
    fmul dword [di - 30] ; -0.995683 *
    fld st2 ; frame
    fidiv word [di - 18] ; 128.0 /
    fld st2 ; y[i]
    fld st5 ; x[i]
    fidiv word [di - 16] ; 32.0 /
    fld1 ; 1.0
    faddp
    fdivp
    faddp
    fcos
    fsubr st4 ; x[i] - TOS
    fcos
    fdivr dword [di - 28] ; 0.009429633 /r
    fsub dword [di - 26] ; 5.1594524
    faddp
    fst dword [di + bx] ; x[i]

    fimul word [di - 16] ; 32.0 *
    fmulp
    fptan
    fxch st2 ; swap with frame (due to fptan pushing 1)
    fdivp
    faddp ; 1 at TOS-1 now here to add
    faddp
    fst dword [di + bx + 4]

    ; sy
    fldpi
    faddp
    ftst
    fnstsw ax
    sahf
    fdivr dword [di - 36] ; 0.8561283
    fadd dword [di - 24] ; 0.36063576
    ; sx
    fld dword [di + bx]
    fiadd word [di - 14] ; 14
    fidiv word [di - 16] ; 32.0
    jb sky
ground:
    fimul word [di - 12] ; 4
    fxch st1
    fmul dword [di - 22] ; 0.98328006
    jmp ssky 
sky:
    fld st1
    fidiv word [di - 69] ; 2
    fsubp
    fxch st1
ssky:

    ; compute coords
    fimul word [di - 10] ; SCREEN_WIDTH
    frndint
    faddp
    fimul word [di - 10] ; SCREEN_WIDTH
    frndint
    fimul word [di - 12] ; 4
    fistp dword [di - 4] ; index
    mov si,[di - 4] ; offset
    mov dx,[di - 2] ; bank

    cmp dx,15 ; check boundaries -> ceil(640x400x4 / 65536)
    ja short skp
    xor bx,bx
    mov ax,$4f05
    int 10h

    px_sum:
        add al,[es:si+bx]
        adc ah,bh
        inc bx
        jpo short px_sum

    ; compute color increments index
    div byte [di - 100] ; 192 (from code data)
    cbw
    xchg ax,bp
    mov ah,[bp + di - 33] ; pal value

    dec bx
    px_sat:
        mov al,ah
        ;and al,3
        aam 4

        add [es:si+bx],al
        salc
        or [es:si+bx],al
        ;shr ah,2
        dec bx
        jns short px_sat
skp:
    inc ecx
    jmp lp

db 0x5b ; ~0.8561283
pal:
    db 0x3f,0x39,0x2b,0x15

; truncated using https://www.h-schmidt.net/FloatConverter/IEEE754.html
fp_constants:
    db 0xe5   ; 0xe6
    dw 0xbf7e ; ~-0.995683
    dw 0x3c1a ; ~0.009429633 noise spread, low value smooth things out
    dw 0x40a5 ; ~5.1594524
    dw 0x3eb8 ; ~0.36063576
    dw 0x3f7b ; ~0.98328006
    dw 0x0080 ; 128
    dw 0x0020 ; 32
    dw 0x000e ; 14
    dw 0x0004 ; 4
    dw 0x0280 ; SCREEN_WIDTH
    ;dw 0x0320 ; 800x600

fp_data equ $+8