This is how I chose to solve the problem of making trig calculations with the PICmicro MCU.
Note that this does not use the CORDIC method because I imagine that would hog too much valuable computing time. This solution probably should not be used to try to get someone to Mars or anything — I’ll work on a more complete solution at a later date. Until then, I hope this will do.
First of all, I attacked this problem with a few basic assumptions:
- 256-degree circle – It makes sense, when using a computer, to base everything on powers of two. Here I’m going to consider circle in terms of only 256 degrees. This is because 256 can be represented nicely by one byte and it’s close enough to 360 degrees that results shouldn’t be too far off base. Also, common angles can be easily represented as well, with a half-circle, quarter, eighth, sixteenth, etc, all settling nicely at a power of two. If it’s difficult to think of a circle in 256 degrees, just remember that the number 360 was completely arbitrary in the first place. (You can thank the ancient Babylonians for it, as well as for the number of minutes in an hour and seconds in a minute. They used a base-60 counting system.)
- Unsigned magnitude, zero to 255 – Sines and cosines are usually represented by fractions. Since computers can’t tolerate fractions very easily, I chose to return results as a positive number from zero to 255, where zero is, well, zero, and 255 represents 1. Why only positive numbers? Because when you take the sine of an angle in a 256-degree circle, the most significant bit of the angle you’re looking for will quickly tell you if the result should be positive or negative. This is because quadrants 3 and 4 (angles 128-255) generate negative sine values. (As it turns out, the MSB for all these numbers is 1.) This allows me to use the full byte for sine values, effectively doubling the granularity of the result.
- One-quadrant sine-only look-up table – The sine values for any one quadrant can be used to determine any of the other trig values for any other quadrant. I’ll provide a refresher on trigonometric identities below. Suffice it to say that one quadrant of sine values is all we need.
- Non-CORDIC Solution – CORDIC would take too long to implement with a negligable increase in accuracy for typical applications. This solution is a look-up table that essentially digitizes the first quadrant of a sine wave as shown below.
The code listed below was written for mid-range PIC MPU’s like the 16f84, but can be ported down with a few minor changes. (For example, you may need to use comf and incf to negate values in the lower-end divices.) It finds the unsigned sine of an angle in eight instruction cycles (32 clock cycles), or 8us on a 4MHz microprocessor. With a few extra instructions (about 13 total) you can determine the sign as well as the magnitude. Determining the not-as-accurate signed sine requires 17 instruction cycles and returns a value ranging from -127 to 127. It takes about 15 cycles to determine the 8-bit cosine magnitude and discover its sign, while the signed cosine requires 18 instruction cycles. If you see a more efficient way of doing this, please let me know. Because instruction cycles vary, I don’t recommend using this for bit-banging routines where timing is critical. Maybe later I’ll do a version that takes the same number of cycles regardless of what function its performing.Before you forget, please help spread the word about this site by referring friends and associates!
The sinw subroutine does not call any other routines. This should be helpful if you’re working with a device that has a short call stack. If I were to return a signed value for the sine, I would have needed to make an additional call on the stack. I believe the way I’ve done it works out nicer. As it is, only ten or so instructions are required to determine both the sine of the angle and it’s sign (positive or negative).
The code below is well commented and includes examples of how to compute both sine and cosine (signed and unsigned) using the sine look-up table. Also, in the comments I’ve listed all the other trig equations I can think of — you can use these to figure out tangents, cosecants, etc.
It might be helpful to remember when sines and cosines are positive or negative. In our 256-degree circle, each quadrant is represented by 64 degrees. So Q1 is angles 0-63, Q2 is 64-127, etc.
- Q1 (0-63 deg) – Sine and cosine are both positive. Bits 6 and 7 are both clear.
- Q2 (64-127 deg) – Sine is positive, cosine is negative. Bit 6 is set and 7 is clear.
- Q3 (128-191 deg) – Sine and cosine are both negative. Bit 6 is clear and 7 is set.
- Q4 (192-255 deg) – Sine is negative, cosine is positive. Bits 6 and 7 are both set.
Bit seven (the MSB) of the 8-bit angle indicates if the sine is positive or negative (zero or one, respectively). Bit six indicates whether or not the sine wave is a mirror image of the preceeding quadrant. (Q2 is a mirror image of Q1 and Q4 is a mirror image of Q3.) This bit can therefor be used to indicate whether or not we need to read from the start or the end of the look-up table. The following rules always apply:
- If bit 7 (MSB) is set, the angle is in Q3 or Q4, where sine is always negative.
- If bit 7 (MSB) is clear, the angle is in Q1 or Q2, where sine is always positive.
- If bit 6 is set, the angle is in Q2 or Q4, which are mirror-images of Q1 and Q3, so read the table backwards.
- If bit 6 is clear, the angle is in Q1 or Q3, which can be read directly from the look-up table.
We only use the table to look up a sine value. Rather than trying to use the table to also look up a cosine, figure out what angle’s sine returns the same value and go to the table for that. This will ensure you can determine the correct sign (positive or negative) of the magnitude returned. To really drive this point home, let’s take a look at each case for the two most significant bits in an 8-bit unsigned angle:
- 00 – Angle is in Q1, the sine is positive, and it can be looked up directly from the table
- 01 – Angle is in Q2, the sine is positive, and it can be looked up by indexing the table backwards
- 10 – Angle is in Q3, the sine is negative, and it can be looked up directly from the table
- 11 – Angle is in Q4, the sine is negative, and it can be looked up by indexing the table backwards
Skip this paragraph if you’re not a beginner. Why does rolling (shifting) bits to the left or right magically double or halve whatever number your shifting? It works in a similar way for any counting system — just think about the numbers you’re used to using. If you “shift” the digits in the decimal number 123 to the left you get 1230 — effectively multiplying the original number by ten, which is the base of the decimal counting system. Likewise, shifting the hexadecimal value F3 to the left you get F30. F3 in decimal notation is 243. 243 times 16 (the base of hex) is 3,888, which in hexadecimal is F30.
Feel free to use this code, but be sure to leave the copyright information intact.
#define _version "0.01"
; Author: Alex Franke
; Copyright: (c)2002 by Alex Franke. All rights reserved.
;
; Title: sin
; Description: Computes the sine of an angle given in degrees of a
; 256-degree circle.
;
; Device: Microchip PIC MCU, p16f84
; or other mid-range PIC MCU's
;
; Update History: 6/11/02 Created
;
LIST R=DEC
INCLUDE "p16f84.inc"
; Registers
CBLOCK 0x020
sinTemp, theta, temp, temp2
ENDC
__CONFIG _CP_OFF & _WDT_OFF & _RC_OSC & _PWRTE_ON
; Main
org 0
Start
; Set up
movlw 84 ; specify an angle value
nop ; pause here to change value if desired for testing
movwf theta ; save that angle for use later
; now start the examples...
; Assume the angle you're trying to find is in theta and you want to find...
; ...its sine (0 to 255) and its sign (positive or negative)
movf theta, w ; copy to w so we can work with it
call sinw ; we now have the sine in w
btfsc theta, 7 ; test the original angle to see if it's in Q3 or Q4
nop ; ... if so, we know it's a negative, so handle that here
; ...its sine (-127 to 127)
movf theta, w ; copy to w so we can work with it
call sinw ; we now have the sine in w
movwf temp ; store result in temp register
bcf STATUS, C ; clear the carry bit if it's set -- we're about to roll
rrf temp, w ; roll temp result right to cut in half, store in w
btfsc theta, 7 ; test the original angle to see if it's in Q3 or Q4
sublw 0 ; ...it is, so negate w
nop ; sine (-127 to 127) is now in w
; ...its cosine (0 to 255) and its sign (positive or negative)
; A cosine is basically a sine, but shifted to the left one quadrant. This means
; that a sine is a cosine shifted to the right. So the cosine of any angle is
; equal to the sine at the same position one quadrant to the right, or
; cos(x)=sin(90+x). So to convert to sine we add 90 deg (64 degrees in our 256-deg circle)
; to the original angle.
movf theta, w ; copy to w so we can work with it
addlw 64 ; add pi/2 so we can just take the sine
movwf temp ; move new angle to temp because we test this later
call sinw ; we now have the sine in w
btfsc temp, 7 ; test the original angle to see if it's in Q3 or Q4
nop ; ... if so, we know it's a negative, so handle that here
; ...its cosine (-127 to 127)
movf theta, w ; copy to w so we can work with it
addlw 64 ; convert cos to sin
movwf temp2 ; store in temp2 register so we can compare for negative
call sinw ; we now have the sine in w
movwf temp ; store result in temp register
bcf STATUS, C ; clear the carry bit if it's set -- we're about to roll
rrf temp, w ; roll temp result right to cut in half, store in w
btfsc temp2, 7 ; test the original angle to see if it's in Q3 or Q4
sublw 0 ; ...it is, so negate w
nop ; sine (-127 to 127) is now in w
; ...its tangent, cotangent, secant, cosecant, etc
; Rather than going into a dissertation on how computers divide numbers,
; let me just recommend using a good divide macro and any of the following
; equations...
;
; tan( x ) = sin( x ) / cos( x )
; cot( x ) = 1 / tan( x ) = cos( x ) / sin( x )
; csc( x ) = 1 / sin( x )
; sec( x ) = 1 / cos( x )
;
; These Pythagorean identities may also be useful
;
; sin^2( x ) + cos^2( x ) = 1
; sec^2( x ) = 1 + tan^2( x )
; csc^2( x ) = 1 + cot^2( x )
;
; These will help with adding and subtracting
;
; sin( x + y ) = ( sin( x ) * cos( y ) ) + ( cos( x ) * sin( y ) )
; sin( x - y ) = ( sin( x ) * cos( y ) ) - ( cos( x ) * sin( y ) )
; cos( x + y ) = ( cos( x ) * cos( y ) ) - ( sin( x ) * sin( y ) )
; cos( x - y ) = ( cos( x ) * cos( y ) ) + ( sin( x ) * sin( y ) )
; tan( x + y ) = ( tan( x ) + tan( y ) ) / ( 1 - ( tan( x ) * tan( y ) ) )
; tan( x - y ) = ( tan( x ) - tan( y ) ) / ( 1 + ( tan( x ) * tan( y ) ) )
;
; And for doubled angles...
;
; sin( 2x ) = 2 * sin( x ) * cos( x )
; cos( 2x ) = cos^2( x ) - sin^2( x ) = ( 2 * cos^2( x ) ) - 1 = 1 - ( 2 * sin^2( x ) )
; tan( 2x ) = ( 2 * tan( x ) ) / ( 1 - tan^2( x ) )
;
; And as we've used above... What's a cosine? It's on of several cofunctions:
;
; sin( 90 deg - x ) = cos( x ) cos ( x - 90 deg ) = sin( x )
; sec( 90 deg - x ) = csc( x ) csc ( x - 90 deg ) = sec( x )
; tan( 90 deg - x ) = cot( x ) cot ( x - 90 deg ) = tan( x )
;
; And finally... What's an arcsine, arccosecant, arccosine, arcsecant, arctangent,
; and arccotangent? They're simply the inverses...
;
; arcsin ( sin(x) ) = x arccos ( cos(x) ) = x
; arcsec ( sec(x) ) = x arccsc ( csc(x) ) = x
; arctan ( tan(x) ) = x arccot ( cot(x) ) = x
;
goto Start
Finished
goto $
sinw
; If bit 7 is clear, we're in quadrant 1 or 3
; If bit 8 is clear, we're in quadrant 1 or 2
;
; bit 8 bit 7
; Quadrant 1 0 0 read from begining of table
; Quadrant 2 0 1 read from end of table
; Quadrant 3 1 0 read from begining of table and negate
; Quadrant 4 1 1 read from end of table and negate
;
; Angles in quadrants 3 and 4 result in negative sines. Because returning
; a negative value from this subroutine would require an additional entry
; on the call stack (a call to the sine table to get the value before
; negating it), we'll just leave it up to the user to negate the result
; if necessary. This saves us an additional call.
; This table returns a value sclaed 0 to 255
movwf sinTemp ; Move input to temp variable to test it
andlw 0x03F ; Clear the first two bits of input - number needs to be <= 64
btfsc sinTemp, 6 ; read from the begining or the end of the table?
sublw 64 ; ...read from end, so w = TableSize - Input
addwf PCL, f ; Position program counter to appropriate result
dt 0, 6, 13, 19, 25, 31, 37, 44
dt 50, 56, 62, 68, 74, 80, 86, 92
dt 98, 103, 109, 115, 120, 126, 131, 136
dt 142, 147, 152, 157, 162, 167, 171, 176
dt 180, 185, 189, 193, 197, 201, 205, 208
dt 212, 215, 219, 222, 225, 228, 231, 233
dt 236, 238, 240, 242, 244, 246, 247, 249
dt 250, 251, 252, 253, 254, 254, 255, 255, 255
end