Subversion Repositories Projects

Compare Revisions

Ignore whitespace Rev 1686 → Rev 1687

/NGVideo5_8/tags/v1.33/mymath.c
0,0 → 1,184
/****************************************************************/
/* */
/* NG-Video 5,8GHz */
/* */
/* Copyright (C) 2011 - gebad */
/* */
/* This code is distributed under the GNU Public License */
/* which can be found at http://www.gnu.org/licenses/gpl.txt */
/* */
/****************************************************************/
 
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <avr/pgmspace.h>
 
#include "mymath.h"
 
// http://www.mikrocontroller.net/articles/AVR_Arithmetik#avr-gcc_Implementierung_.2832_Bit.29
unsigned int sqrt32(uint32_t q)
{ uint16_t r, mask;
uint16_t i = 8*sizeof (r) -1;
r = mask = 1 << i;
for (; i > 0; i--) {
mask >>= 1;
if (q < (uint32_t) r*r)
r -= mask;
else
r += mask;
}
if (q < (uint32_t) r*r)
r -= 1;
return r;
}
 
// aus Orginal MK source code
// sinus with argument in degree at an angular resolution of 1 degree and a discretisation of 13 bit.
const uint16_t sinlookup[91] = {0, 143, 286, 429, 571, 714, 856, 998, 1140, 1282, 1423, 1563, 1703, \
1843, 1982, 2120, 2258, 2395, 2531, 2667, 2802, 2936, 3069, 3201, 3332, \
3462, 3591, 3719, 3846, 3972, 4096, 4219, 4341, 4462, 4581, 4699, 4815, \
4930, 5043, 5155, 5266, 5374, 5482, 5587, 5691, 5793, 5893, 5991, 6088, \
6183, 6275, 6366, 6455, 6542, 6627, 6710, 6791, 6870, 6947, 7022, 7094, \
7165, 7233, 7299, 7363, 7424, 7484, 7541, 7595, 7648, 7698, 7746, 7791, \
7834, 7875, 7913, 7949, 7982, 8013, 8041, 8068, 8091, 8112, 8131, 8147, \
8161, 8172, 8181, 8187, 8191, 8192};
 
//double c_cos_8192(int16_t angle)
int16_t c_cos_8192(int16_t angle)
{
int8_t m,n;
int16_t sinus;
 
angle = 90 - angle;
// avoid negative angles
if (angle < 0)
{
m = -1;
angle = -angle;
}
else m = +1;
 
// fold angle to interval 0 to 359
angle %= 360;
 
// check quadrant
if (angle <= 90) n = 1; // first quadrant
else if ((angle > 90) && (angle <= 180)) {angle = 180 - angle; n = 1;} // second quadrant
else if ((angle > 180) && (angle <= 270)) {angle = angle - 180; n = -1;} // third quadrant
else {angle = 360 - angle; n = -1;} //fourth quadrant
// get lookup value
sinus = sinlookup[angle];
// calculate sinus value
return (sinus * m * n);
}
 
/* ----------------------------------------------------------------------------------------------
 
atan2.S
 
Author: Reiner Patommel
atan2.S uses CORDIC, an algorithm which was developed in 1956 by Jack Volder.
CORDIC can be used to calculate trigonometric, hyperbolic and linear
functions and is until today implemented in most pocket calculators.
The algorithm makes only use of simple integer arithmetic.
 
The CORDIC equations in vectoring mode for trigonometric functions are:
Xi+1 = Xi - Si * (Yi * 1 / 2^i)
Yi+1 = Yi + Si * (Xi * 1 / 2^i)
Zi+1 = Zi - Si * atan(1 / 2^i)
where Si = +1 if Yi < 0 else Si = -1
The rotation angles are limited to -PI/2 and PI/2 i.e.
-90 degrees ... +90 degrees
 
For the calculation of atan2(x,y) we need a small pre-calculated table of
angles or radians with the values Tn = atan(1/2^i).
We rotate the vector(Xo,Yo) in steps to the x-axis i.e. we drive Y to 0 and
accumulate the rotated angles or radians in Z. The direction of the rotation
will be positive or negative depending on the sign of Y after the previous
rotation and the rotation angle will decrease from step to step. (The first
step is 45 degrees, the last step is 0.002036 degrees for n = 15).
 
After n rotations the variables will have the following values:
Yn = ideally 0
Xn = c*sqrt(Xo^2+Yo^2) (magnitude of the vector)
Zn = Zo+atan(Yo/Xo) (accumulated rotation angles or radians)
 
c, the cordic gain, is the product Pn of sqrt(1+2^(-2i)) and amounts to
approx. 1.64676 for n = 15.
 
In order to represent X, Y and Z as integers we introduce a scaling factor Q
which is chosen as follows:
1. We normalize Xo and Yo (starting values) to be less than or equal to 1*Q and
set Zo = 0 i.e. Xomax = 1*Q, Yomax = 1*Q, Zo = 0.
2. With Xo = 1*Q and Yo = 1*Q, Xn will be Xnmax = Q*c*sqrt(2) = 2.329*Q
3. In order to boost accuracy we only cover the rotation angles between 0 and PI/2
i.e. X and Z are always > 0 and therefore can be unsigned.
(We do the quadrant correction afterwards based on the initial signs of x and y)
4. If X is an unsigned int, Xnmax = 65536, and Q = 65536/2.329 = 28140.
i.e.
Xnmax = 65536 --> unsigned int
Ynmax = +/- 28140 --> signed int
Znmax = PI/2 * 28140 = 44202 --> unsigned int
The systematic error is 90/44202 = 0.002036 degrees or 0.0000355 rad.
Below is atan2 and atan in C: */
 
#define getAtanVal(x) (uint16_t)pgm_read_word(&atantab[x])
//uint16_t atantab[T] PROGMEM = {22101,13047,6894,3499,1756,879,440,220,110,55,27,14,7,3,2,1}; // T 16
// Genauigkeit etwas eingeschränkt, da Servos nur in 1° Abständen positioniert werden
uint16_t atantab[T] PROGMEM = {22101,13047,6894,3499,1756,879,440,220,110,55,27,14,7,3};
 
// max. Wert für y oder x ist 76314
int16_t my_atan2(int32_t y, int32_t x)
{ unsigned char i;
uint32_t x1;
int32_t y1;
uint32_t angle = 0;
uint32_t tmp;
 
x1 = abs(x);
y1 = abs(y);
 
if (y1 == 0) {
if (x < 0)
return (180);
else
return 0;
}
if (x1 < y1) {
x1 = SCALED(x1) / y1;
y1 = Q;
}
else {
y1 = SCALED(y1) / x1;
x1 = Q;
}
for (i = 0; i < T-1; i++) {
tmp = x1 >> i;
if (y1 < 0) {
x1 -= y1 >> i;
y1 += tmp;
angle -= getAtanVal(i);
}
else {
x1 += y1 >> i;
y1 -= tmp;
angle += getAtanVal(i);
}
}
 
angle = angle/RQ; // RAD_TO_DEG * angle / Q
if (x <= 0)
angle = 180 - angle;
if (y <= 0)
return(-angle);
return(angle);
}