sine.asm

可移植的MicroChip的PIC全系列单片机软件仿真器

        list    p=16C84,t=ON,c=132,n=80,st=off
        radix   dec


        include "p16c84.inc"

  cblock  0x20

        f_hi    ;Low byte of frequency variable
        f_lo
        x1
        x2
        interp
        index
        product

        fstep_lo
        fstep_hi
  endc

ODD_QUADRANT    equ     4       ;In f_hi
NEG_QUADRANT    equ     5       ;In f_hi
INDEX_MASK      equ     0x0f    ;In f_hi
INTERP_MASK     equ     0xf0    ;In f_lo


        ORG     0               ;Reset Vector

        GOTO    Main

        ORG     4               ;Interrupt Vector


Main

        BCF     STATUS,RP0      ;Point to BANK 0

    ;Initialize the step in frequency
        MOVLW   0x10
        MOVWF   fstep_lo        ;One unit of fstep corresponds to 360/16384 degrees
        CLRF    fstep_hi        ;i.e. 0.02197 degrees. Set fstep = 455 (1c7) to step
                                ;10 degrees for example

        CLRF    f_hi            ;Start off at 0 degrees
        CLRF    f_lo

xxx
        CALL    sine
        NOP

    ;Now that the sine has been calculated, we could do something with it other than
    ;throwing it away.

    ;Make a frequency step. Note, the frequency step doesn't necessarily have to be
    ;constant. E.g. it could be changed as part of a PLL algorithm.

        MOVF    fstep_lo,W
        ADDWF   f_lo,F
        SKPNC
          INCF  f_hi,F
        MOVF    fstep_hi,W
        ADDWF   f_hi,F

        GOTO    xxx


;--------------------------------------------------------
;sine
;
;  The purpose of this routine is to take the sine of the
;16-bit frequency variable f_hi:f_lo to produce a signed
;8-bit result that is within +/- 1 count of the true sine
;value.
;  Only the lower 14 bits of the frequency variable are
;actually used. The frequency variable maps into degrees
;by the following transfer function:
;  degrees = (f & 0x3fff) * 360 / 0x4000
;  The sine output is approximately
;sine = int(127*sin( (f & 0x3fff) * 360 / 0x4000) )
; where
; sin() is the true sine function
; int() is the nearest integer function
;
;The technique used to obtain the sine value is a combination
;of table look-up and first order linear interpolation. Sixteen
;equidistant frequency samples of the first quadrant of sin(x) 
;are stored in sine_table.
;
; The frequency variable is broken down as follows:
; xxQQTTTT IIIIPPPP
; where 
;  xx - don't care
;  QQ - Quadrant: 00 = quadrant 1, 01 = quadrant 2, etc.
;  TTTT - Index into the sine_table.
;  IIII - Interpolation between successive entries in the table
;  PPPP - Phase accumulation (not needed in this function, it's
;         only used to increase the dynamic range in frequency
;         steps).
;Once the frequency has been broken down in to these parts, the
;sine function for the first quadrant can be calculated as follows:
;  x1 = sine_table[index]
;  x2 = sine_table[index+1]
;  sine = x1 + ((x2-x1) * interp) / 16
;The first term, x1, is seen to be a first order approximation to
;sine function. The second term improves on that approximation by
;interpolating between x1 and x2. The interpolation variable interp
;is 0 <= interp <= 15 and it is divided by 16. Consequently, the
;interpolation factor ranges between 0 and 15/16.
;
;The sine function in the other three quadrants can be obtained
;from calculations based on the first quadrant by using the following
;trig identities:
; first, let 0 <= f <= 90, i.e. f covers the first quadrant.
;  quadrant 2:  u = 90 + f,   90 < u < 180
;     sin(u-90) = sin(f)
;     x1 = sine_table(16-index), x2 = sine_table(15-index)
;  quadrant 3:  u = 180 + f,  180 < u < 270
;     sin(u) = sin(f+180) = -sin(f)
;     x1 = -sine_table(index), x2 = -sine_table(index+1)
;  quadrant 4:  u = 270 + f,  270 < u < 360
;     sin(u-90) = sin(f+180) = -sin(f)
;     x1 = -sine_table(16-index), x2 = -sine_table(15-index)
;
;Thus, for quadrants 2 and 4, the sine table is indexed in reverse
;order and for quadrants 3 and 4 the values from the sine table 
;are negated. A slight change is made on this indexing and negation
;scheme so that the operation (x2-x1) * interp / 16 only deals with
;positive numbers. This significantly simplifies the multiplication.
;The modification changes the formula for each quadrant as follows:
; quadrant 1:   (no change)
;   x1 = sine_table[index],  x2 = sine_table[index+1]
;   sine = x1 + ((x2-x1) * interp) / 16
; quadrant 2:
;   x1 = sine_table[15-index],  x2 = sine_table[16-index]
;   sine = x2 - ((x2-x1) * interp) / 16
; quadrant 3:
;   x1 = sine_table[index],  x2 = sine_table[index+1]
;   sine = -(x1 + ((x2-x1) * interp) / 16)
; quadrant 4:
;   x1 = sine_table[15-index],  x2 = sine_table[16-index]
;   sine = -(x2 - ((x2-x1) * interp) / 16)
;
;Input
; f_hi:f_lo  -  16-bit frequency variable
;Output
; W = int(127*sin( (f & 0x3fff) * 360 / 0x4000) )
;
;Execution time: 48 Cycles (for all cases)
;

sine

    ;Get the 4-bit index and add 1 to it.
        MOVF    f_hi,W
        ANDLW   INDEX_MASK
        ADDLW   1
        BTFSC   f_hi,ODD_QUADRANT
         SUBLW  17                      ;Odd quadrants, index = 16 - index
                                        ;Actually: (index + 1) = 17 - (index + 1)
                                        ;                      = 16 - index
        MOVWF   index
        CALL    sine_table              ;Get x2=sin(index+1)
        MOVWF   x2

        DECF    index,W
        CALL    sine_table              ;Get x1=sin(index)
        MOVWF   x1

        SUBWF   x2,W                    ;W=x2-x1, This is always positive. 

    ;Initialize the product of (x2-x1)*interp/16 to 1/2. Note 8/16 == 1/2
    ;(This rounds the product to the nearest integer.)
        CLRF    product
        BSF     product,3       ;(note, product and index could be aliased to
                                ; save one byte of ram).

    ;multiply interp and x2 - x1 and divide by 16. This is actually a 4 by 8
    ;bit multiplication. The division by 16 is implemented with a shift right
    ;one position for each of the four multiplication iterations.

	clrc
	btfsc	f_lo,4
         addwf  product,f       ;Then add (x2-x1) to the product
        rrf     product,f       ;Divide the product by two

	clrc
	btfsc	f_lo,5
         addwf  product,f
        rrf     product,f	;Divide the product by four

	clrc
	btfsc	f_lo,6
         addwf  product,f
        rrf     product,f	;Divide the product by eight

	clrc
	btfsc	f_lo,7
         addwf  product,f
        rrf     product,w	;Divide the product by sixteen


        BTFSS   f_hi,ODD_QUADRANT
          ADDWF   x1,W
        BTFSC   f_hi,ODD_QUADRANT
          SUBWF   x2,W
        BTFSC   f_hi,NEG_QUADRANT
          SUBLW 0
        RETURN

sine_table
        ADDWF   PCL,F
        RETLW   0       ;127*sin(0 * 90/16)
        RETLW   12      ;127*sin(1 * 90/16)
        RETLW   25      ;127*sin(2 * 90/16)
        RETLW   37      ;127*sin(3 * 90/16)
        RETLW   49      ;127*sin(4 * 90/16)
        RETLW   60      ;127*sin(5 * 90/16)
        RETLW   71      ;127*sin(6 * 90/16)
        RETLW   81      ;127*sin(7 * 90/16)
        RETLW   90      ;127*sin(8 * 90/16)
        RETLW   98      ;127*sin(9 * 90/16)
        RETLW   106     ;127*sin(10 * 90/16)
        RETLW   112     ;127*sin(11 * 90/16)
        RETLW   117     ;127*sin(12 * 90/16)
        RETLW   122     ;127*sin(13 * 90/16)
        RETLW   125     ;127*sin(14 * 90/16)
        RETLW   126     ;127*sin(15 * 90/16)
        RETLW   127     ;127*sin(16 * 90/16)

        END