Digital Signal Processing Reference
In-Depth Information
result
*
=x;
return result;
}
The above algorithm is based on the Taylor expansion of the sine function evaluated
using a Horner scheme.
x
3
3!
+
x
5
5!
−
x
7
7!
+
x
9
9!
−
x
11
11!
...
sin
(
x
)=
x
−
It can be shown that
r
will always stay in range of the
Accum
if
x
is in the range of
[
−
π
,
π
]
. The precision of this approximation however deteriorates dramatically as
x
approaches
−
π
and
π
. With a 31 bits precision fixed point and 11 iterations, the error
10
−
4
. To improve the results one can use the symmetry of the sine
function and only use the approximation for
at
π
is around 3
.
5
·
[
−
π
/
2
,
π
/
2
]
. Still the approximation
becomes worse as
x
approaches
[
−
π
/
2
]
and
[
π
/
2
]
. For a 31 bits precision fixed
10
−
8
. The results can be
point and 11 iterations, the error at
π
/
2 is around 4
.
5
·
further improved by using the formula
sin
(
2
x
)=
2sin
(
x
)
cos
(
x
)
The approximation for sin
.
For this range only a few iterations are required to obtain high precision. With the
31 bits precision fixed point and 11 iterations, the error is erratic over the whole
range and is maximally 9
(
x
)
and cos
(
x
)
will now only be used for
x
∈
[
−
π
/
4
,
π
/
4
]
10
−
10
, which is around the precision of the fixed point.
·
Note that if
x
∈
[
−
π
/
4
,
π
/
4
]
all values involved in the computation will stay
smaller than 1
0, except for the constant
1.0k
itself. This suggests that using a
Fixed
type instead of
Accum
would suffice. Experimental results show that
replacing the
1.0k
constant by the slightly smaller
1.0r
hurts the precision only
within the precision of the
Fixed
itself.
An optimizing compiler would take
x
*
x
out of the loop and compute it only
once. It could also avoid the integer multiplication and complicated mixed type
division (integer and fixed point) by putting the values for
1.0k/(f
*
(f-1))
in some table. Alternatively the compiler could simply fully unroll the whole
loop. In that case
f
*
(f-1)
will be computed at compile time. Standard compiler
optimizations next can turn the division by a constant into a multiplication by a
constant. The resulting loop then only contains one multiplication and a multiply
subtract operation.
If one takes a careful look at the accuracy of the above algorithm, then one will
see that there is a trend in the error that is the result of the used Taylor expansion.
Variations of the above Horner scheme exist, where the constants 1
.
are
slightly changed to obtain a uniform distribution of the error over the whole range
of
x
. As a final note it must be said that the CORDIC (coordinate rotation digital
computer) algorithm for computing the trigonometric functions is very often used
as it requires only addition and shift operations.
/
(
f
∗
(
f
−
1
))
