I have been developing a real time application in C and I need to calculate sine and cosine of a given angle in radians. In respect to the fact that the application is a real time application the functions from the standard library are not suitable for my purposes due to the fact that there is not fixed calculation time. Based on that I have decided to implement those functions in my own.
From the execution time point of view it seems to me that the best approach for real time application is to use the look-up table based method of calculation with improvement of precion based on linear interpolation. I have exploited the symmetry of the sine and cosine functions and I have defined the look-up table covering the quarter of the period.
Look-up table
#define TABLE_SIZE 256
#define PI 3.14
#define STEP_SIZE (PI/(2*(TABLE_SIZE-1))) // N points divide PI/2 into N-1 segments
// look-up table containing values of sine function over one quarter of period
// angel step size is pi/(2*256) i.e. pi/512 approximately 0.35°
double sinLUT[TABLE_SIZE] =
{
0.000000, 0.006160, 0.012320, 0.018479, 0.024637, 0.030795, 0.036951, 0.043107,
0.049260, 0.055411, 0.061561, 0.067708, 0.073853, 0.079994, 0.086133, 0.092268,
0.098400, 0.104528, 0.110653, 0.116773, 0.122888, 0.128999, 0.135105, 0.141206,
0.147302, 0.153392, 0.159476, 0.165554, 0.171626, 0.177691, 0.183750, 0.189801,
0.195845, 0.201882, 0.207912, 0.213933, 0.219946, 0.225951, 0.231948, 0.237935,
0.243914, 0.249883, 0.255843, 0.261793, 0.267733, 0.273663, 0.279583, 0.285492,
0.291390, 0.297277, 0.303153, 0.309017, 0.314870, 0.320710, 0.326539, 0.332355,
0.338158, 0.343949, 0.349727, 0.355491, 0.361242, 0.366979, 0.372702, 0.378411,
0.384106, 0.389786, 0.395451, 0.401102, 0.406737, 0.412356, 0.417960, 0.423549,
0.429121, 0.434676, 0.440216, 0.445738, 0.451244, 0.456733, 0.462204, 0.467658,
0.473094, 0.478512, 0.483911, 0.489293, 0.494656, 0.500000, 0.505325, 0.510631,
0.515918, 0.521185, 0.526432, 0.531659, 0.536867, 0.542053, 0.547220, 0.552365,
0.557489, 0.562593, 0.567675, 0.572735, 0.577774, 0.582791, 0.587785, 0.592758,
0.597707, 0.602635, 0.607539, 0.612420, 0.617278, 0.622113, 0.626924, 0.631711,
0.636474, 0.641213, 0.645928, 0.650618, 0.655284, 0.659925, 0.664540, 0.669131,
0.673696, 0.678235, 0.682749, 0.687237, 0.691698, 0.696134, 0.700543, 0.704926,
0.709281, 0.713610, 0.717912, 0.722186, 0.726434, 0.730653, 0.734845, 0.739009,
0.743145, 0.747253, 0.751332, 0.755383, 0.759405, 0.763398, 0.767363, 0.771298,
0.775204, 0.779081, 0.782928, 0.786745, 0.790532, 0.794290, 0.798017, 0.801714,
0.805381, 0.809017, 0.812622, 0.816197, 0.819740, 0.823253, 0.826734, 0.830184,
0.833602, 0.836989, 0.840344, 0.843667, 0.846958, 0.850217, 0.853444, 0.856638,
0.859800, 0.862929, 0.866025, 0.869089, 0.872120, 0.875117, 0.878081, 0.881012,
0.883910, 0.886774, 0.889604, 0.892401, 0.895163, 0.897892, 0.900587, 0.903247,
0.905873, 0.908465, 0.911023, 0.913545, 0.916034, 0.918487, 0.920906, 0.923289,
0.925638, 0.927951, 0.930229, 0.932472, 0.934680, 0.936852, 0.938988, 0.941089,
0.943154, 0.945184, 0.947177, 0.949135, 0.951057, 0.952942, 0.954791, 0.956604,
0.958381, 0.960122, 0.961826, 0.963493, 0.965124, 0.966718, 0.968276, 0.969797,
0.971281, 0.972728, 0.974139, 0.975512, 0.976848, 0.978148, 0.979410, 0.980635,
0.981823, 0.982973, 0.984086, 0.985162, 0.986201, 0.987202, 0.988165, 0.989092,
0.989980, 0.990831, 0.991645, 0.992421, 0.993159, 0.993859, 0.994522, 0.995147,
0.995734, 0.996284, 0.996795, 0.997269, 0.997705, 0.998103, 0.998464, 0.998786,
0.999070, 0.999317, 0.999526, 0.999696, 0.999829, 0.999924, 0.999981, 1.000000
};
Sine function calculation
double sine(double arg, double lut[TABLE_SIZE]){
double retval;
double rem;
uint8_t index;
if(arg >= 0 && arg <= PI/2){
// first quadrant
index = arg/STEP_SIZE;
rem = arg - index*STEP_SIZE;
if(rem > 0){
// sine value for given argument isn't directly in the lut
if(index == (TABLE_SIZE-1)){
// last point in the lut so the interval for the interpolation
// is the last interval bounded with the index-1 and index
index -= 1;
}
retval = (lut[index+1] - lut[index])/STEP_SIZE*rem + lut[index];
}else{
// sine value for given argument is directly in the lut
retval = lut[index];
}
}else if(arg > PI/2 && arg <= PI){
// second quadrant
index = (PI - arg)/STEP_SIZE;
rem = (PI - arg) - index*STEP_SIZE;
if(rem > 0){
if(index == (TABLE_SIZE-1)){
index -= 1;
}
retval = (lut[index+1] - lut[index])/STEP_SIZE*rem + lut[index];
}else{
retval = lut[index];
}
}else if(arg > PI && arg <= 3*PI/2){
// third quadrant
index = (arg - PI)/STEP_SIZE;
rem = (arg - PI) - index*STEP_SIZE;
if(rem > 0){
if(index == (TABLE_SIZE-1)){
index -= 1;
}
retval = (-lut[index+1] + lut[index])/STEP_SIZE*rem - lut[index];
}else{
retval = -lut[index];
}
}else{
// fourth quadrant
index = (2*PI - arg)/STEP_SIZE;
rem = (2*PI - arg) - index*STEP_SIZE;
if(rem > 0){
if(index == (TABLE_SIZE-1)){
index -= 1;
}
retval = (-lut[index+1] + lut[index])/STEP_SIZE*rem - lut[index];
}else{
retval = -lut[index];
}
}
return retval;
}
Cosine function
double cosine(double arg, double lut[TABLE_SIZE]){
double temp;
temp = (arg + PI/2);
if(temp > 2*PI){
temp -= 2*PI;
}
return sine(temp, lut);
}
#define PI 3.14
I LOLed. \$\endgroup\$