console-dev.de

Home of VsTortoise, VisualHAM, N3D and HEL Library

Archive for July, 2009

sine approximation with fixed point math part 2

Posted by Peter Schraut under NDS, Programming

In my previous post sine approximation with fixed point math I wrote the output precision of f32sin is getting worse when the incoming radians parameter grows and guessed it is due to using fixed point.

Well, Jasper “cearn” Vijn suggested to use a higher resolution fixed point format for 1.0f/(2*pi) and this solves the problem indeed! But it does not end here!

He also created an awesome article about the ins and outs of fixed point sine approximation and presents a couple of highly precise and optimized routines.

Head over to Jasper Vijn’s document Another fast fixed-point sine approximation now!

I’ll post a corrected version of the f32sin routine from my previous post for completeness, but I highly encourage to vists Jasper’s site and get one of his!

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31

// Gets sin of specified radians.
// Input and output format is fixed point sign.19.12
int32 f32sin(int32 radians)
{
  const int32 PI2_PRECISION_BITS = 20;
  const int32 PI2_RECIPROCAL     = (1.0f / (2*M_PI)) * (1<<PI2_PRECISION_BITS);
 
  // Divide incoming radians by 2*PI, so it is a "normalized"
  // angle between -1 and +1, where 1 equals 2*PI (360 degree).
  int32 value = (radians * PI2_RECIPROCAL) >> PI2_PRECISION_BITS;
 
  // Now that we have our normalized angle,
  // we just discard all numbers which are not in the range -1..+1.
  // Since this is fixed point, a simple AND operation can be used.
  value &= floattof32(1)-1;
 
  // Always wrap normalized angle to -0.5(-PI)..+0.5(+PI)
  if(value > floattof32(0.5f))
    value -= floattof32(1); // subtract 2*PI in normalized form
  else
    if(value < floattof32(-0.5f))
      value += floattof32(1); // add 2*PI in normalized form
 
  // Convert normalized angle (-1..1) back to radians (-2*PI..2*PI)
  value = f32mul(value, floattof32(2*M_PI));
 
  // Approximate sin value
  const int32 B = floattof32(1.2732395447351f);
  const int32 C = floattof32(-0.405284734569f);
  return f32mul(B, value) + f32mul(f32mul(C, value), f32abs(value));
}
1 COMMENT

sine approximation with fixed point math

Posted by Peter Schraut under NDS, Programming

When I was writing the Nintendo DS 4k intro, I obviously needed sin and cos functions to set up rotation matrices. Since using the standard C library math module would have blown the 4kb limit to mars, I was looking for an alternative that does not use much memory.

I was searching for an approximation method using taylor series and came across the following sites:

The only task that remained was porting it from float to fixed point, which is rather trivial, yay!

Most amazing for me is how precise the approximation actually is. Here is a screenshot of the original sinf routine (green) and the approximated one (red) using the source code below:

Sine approximation

Sine approximation


Well, when the incoming radians parameter is getting bigger (1000 and above), precision becomes is not so good anymore, but this is rather a fixed point problem and the way I implemented it, I guess.

The routine from polygon labs also can’t really handle radians greater than 2*PI, which I tried to fix by normalizing it to an angle between -1..+1 and then just use a bit-AND to discard all numbers outside this range.

I wanted to keep the radians format for the incoming parameter, because it’s how the standard sinf routine works. Using another range, like 0..4096 to represent a full rotation would probably solve that precision problem with large radians because this would remove the need for the first two code lines, but I wanted 2*PI. ;-)

Here comes the source code, to be compiled with libnds and devkitARM:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44

// Gets sin of specified radians.
// Input and output format is fixed point sign.19.12
int32 f32sin(int32 radians)
{
  // Offset to slightly restore precision when radians is large.
  const int32 offset= f32mul(radians, floattof32(0.00147f));
 
  // Divide incoming radians by 2*PI, so it is a "normalized"
  // angle between -1 and +1, where 1 equals 2*PI (360 degree).
  int32 value = f32mul(radians + offset, floattof32(1.0f / (2*M_PI)));
 
  // Now that we have our normalized angle,
  // we just discard all numbers which are not in the range -1..+1.
  // Since this is fixed point, a simple AND operation can be used.
  value &= floattof32(1)-1;
 
  // Always wrap normalized angle to -0.5(-PI)..+0.5(+PI)
  if(value > floattof32(0.5f))
    value -= floattof32(1); // subtract 2*PI in normalized form
  else
  if(value < floattof32(-0.5f))
    value += floattof32(1); // add 2*PI in normalized form
 
  // Convert normalized angle (-1..1) back to radians (-2*PI..2*PI)
  value = f32mul(value, floattof32(2*M_PI));
 
  // Approximate sin value
  const int32 B = floattof32(1.2732395447351f);
  const int32 C = floattof32(-0.405284734569f);
  return f32mul(B, value) + f32mul(f32mul(C, value), f32abs(value));
}
 
inline int32 f32abs(int32 a)
{
  return a >= 0 ? a : -a;
}
 
// multiplies two fixed point numbers.
// since the result is shifted rather than the two
// operands, make sure a and b are small.
inline int32 f32mul(int32 a, int32 b)
{
  return (a*b)>>12;
}

Let me know when you fixed the precision problem with large radian values :-)

2 COMMENTS
Subscribe to console-dev.de
This site uses a Hackadelic PlugIn, Hackadelic SEO Table Of Contents 1.6.0.