Maths: Sine16Bit

       Name: Sine16Bit                                               [Show more]
       Type: Subroutine
   Category: Maths
    Summary: Calculate the sine of a 16-bit number
  Deep dive: Trigonometry
    Context: See this subroutine in context in the source code
 References: This subroutine is called as follows:
             * ProjectAxisAngle calls Sine16Bit


 Sets (A Y) = sin(X W)

 where (X W) is a 16-bit angle, with 0 to 65535 representing a quarter circle.

 The sine lookup table contains 256 entries that represent a quarter circle, so
 we work out the result by first looking up the sine for the high byte X, and
 then interpolating W/256 of the way between the results for X and X + 1. It
 might help to think of X being an integer (0 to 255) and W being a fraction
 (0 to 255) and we're trying to map X.W onto the sine table by finding the
 entry for X and doing linear interpolation between X and X + 1 for the
 fractional amount.


 Arguments:

   (X W)                The 16-bit angle, representing a quarter-circle in the
                        range 0 to 255


 Returns:

   (A Y)                sin(X W)

   (X W)                Unchanged


.Sine16Bit

 STX H                  \ Store X in H for later

 SEC                    \ Set (A I) = sin(X + 1) - sin(X)
 LDA sinLo+1,X          \
 SBC sinLo,X            \ starting with the low bytes
 STA I

 LDA sinHi+1,X          \ And then the high bytes
 SBC sinHi,X

                        \ Let's call this dX, the difference between sin(X + 1)
                        \ sin(X) and, so (A I) = dX

 LSR A                  \ The maximum differential value we can have is 402, so
 ROR I                  \ we halve (A I) to fit the result into one byte, like
 LDX I                  \ this:
                        \
                        \   (A I) = dX / 2
                        \
                        \ and then setting X to the low byte in I, so:
                        \
                        \   X = dX / 2

 LDY W                  \ Set Y = W

 JSR Multiply8x8        \ Set (A V) = X * Y
                        \           = (dX / 2) * W

 LDX H                  \ Set X to the original argument we stored above

 ASL V                  \ Set (A V) = (A V) * 2
 ROL A                  \           = (dX / 2) * W * 2
                        \           = dX * W
                        \
                        \ and put bit 7 of (A V) into the C flag, so in effect
                        \ we have:
                        \
                        \   (C A V) = dX * W

 PHP                    \ Store the C flag on the stack

 CLC                    \ Set (A Y) = (C A) + sin(X)
 ADC sinLo,X            \           = (dX * W / 256) + sin(X)
 TAY                    \
                        \ starting with the low bytes

 LDA #0                 \ And then adding the carry to the high byte, so now we
 ADC sinHi,X            \ have the interim result for (0 A) + sin(X)

 PLP                    \ And finally adding C (which we retrieve from the
 ADC #0                 \ stack) to give the final result for (C A) + sin(X)

                        \ So at this point we have:
                        \
                        \   (A Y) = sin(X) + (dX * W / 256)
                        \
                        \ which is the result we want, as it takes sin(X) and
                        \ adds on W/256 of the difference between sin(X) and
                        \ sin(X + 1)

 RTS                    \ Return from the subroutine