#include<stdio.h>
#include<stdlib.h>
//M_PI_4使用のためのマクロ定義
#define _USE_MATH_DEFINES
#include<math.h>

//MAXBITSの変更で精度が変動
#define MAXBITS 16

int flag;
static double K;
static double Gain;
static double atanTable[MAXBITS];

static double  gain(void);
void initCordic(void);
void cordic(double* x0, double* y0, double* z0);
double sincosCordic(double a, double* cos);
void dec2bin_for_double(double dn, int bitwidth, int fractional_number, int *bn);
void B_D_output(int bitwidth, int fractional_number,int *c_out);


//データテーブル作成
//初期化関数
void initCordic()
{
   double t = 1;
   int 	i;
   int bitwidth, fractional_number;
   bitwidth = 16;
   fractional_number = 14;

   //MAXBITS分のテーブルを作成
   //printf(" ---------------atanTable----------------\n");
   for (i = 0; i < MAXBITS; ++i){
	    //テーブル作成には標準ライブラリarnctanを使用
        atanTable[i] = atan(t);		
        t /= 2;												
		//テーブルの表示
		//printf(" atanTable[%d] = %.18e\n",i ,atanTable[i]);
   }
   //printf(" ----------------------------------------\n");

	//初期値決定
	Gain = 1/ gain();

}

//cordic関数
void cordic(double* x0, double* y0, double* theta)
{

    double t;
    double x, y, z;
    int i;
	int bitwidth, fractional_number;

    t = 1;
    x = *x0; y = *y0; z = 0;
	bitwidth = 16;
    fractional_number = 14;


	//MAXBITSの精度で繰り返し
    for (i = 0; i < MAXBITS; ++i) {	
        double x1;


		//z >  theta  の時にif文を実行
        if ( z > *theta) {
			//変換式
            x1 = x + y*t;	
            y = y - x*t;
            z = z - atanTable[i];
        }
		//z < theta  の時にelse文を実行
        else {
			//この場合は変換式で正負が逆
            x1 = x - y*t;
            y = y + x*t;
            z = z + atanTable[i];
        }

        x = x1;
        t /= 2;
    }

    *x0 = x;
    *y0 = y;
    *theta = z;
}

//sin cos 計算
double sincosCordic(double a, double* cos)	
{
   //sin変数を初期化
   double sin = 0;
   //cos(0)のx座標値を代入
   *cos =  Gain;
   cordic(cos, &sin, &a);
   return sin;
}

//ゲイン算出用関数
static double  gain()
{
	double x, y, z;
	x = 1;
	y = 0;
	z = 0;
	cordic(&x, &y, &z);
	return x;
}

main()
{
   //値演算用変数
   double sin_v;
   double theta;
   double cos_v;
   double true_sin_v, true_cos_v;
   double difference;
   double i=1,l=0;
   int sign=0;

   //二進数変換用変数
   int  c_out[16];
   int bitwidth, fractional_number;

   //arctan,gain()を実行
   initCordic();
   flag = 0;
   
   //入力thetaの範囲　[-π/2 : π/2]
	//#define M_PI_2     1.57079632679489661923
   //#define M_PI_4     0.785398163397448309616
   
   
   //theta = M_PI_2;
   //theta = M_PI_4 * 3 / 2;
     theta = M_PI_4;
   //theta = M_PI_4 / 2;
   //theta = 0;
   //theta = -M_PI_4 / 2;
   //theta = -M_PI_4;
   //theta = -M_PI_4 * 3 / 2;
   //theta = -M_PI_2;


   //theta = M_PI_2 - M_PI_2 / 128;
   //theta = M_PI_2 - M_PI_2 / 64;
   //theta = M_PI_2 / 128;
   //theta = M_PI_2 / 64;
   //theta = -M_PI_2 / 64;
   //theta = -M_PI_2 / 128;
   //theta = -M_PI_2 + M_PI_2 / 64;
   //theta = -M_PI_2 + M_PI_2 / 128;
   

   //入力値正規化用
   //theta = theta/M_PI;

   flag = 1;
   sin_v = sincosCordic(theta, &cos_v);

   
   // fixed point singed 16bit (fractional number is 14bit)
   bitwidth = 16;
   fractional_number = 14;

   //-------------------------------------------------------------------
   //演算値(sin)

   // decimal (double) to binary
   dec2bin_for_double(sin_v, bitwidth, fractional_number, c_out);

   //decimal output
   printf("sin (%f) = %.10e\n",theta ,sin_v);

   // binary output
   B_D_output(bitwidth, fractional_number, c_out);

   //-------------------------------------------------------------------
   //真値(sin)
   true_sin_v = sin(theta);

   // decimal (double) to binary
   dec2bin_for_double(true_sin_v, bitwidth, fractional_number, c_out);

   // decimal output
   printf("true sin (%f) = %.10e\n",theta ,true_sin_v);

   // binary output
   B_D_output(bitwidth, fractional_number, c_out);

   //-------------------------------------------------------------------
   //真値と演算値の差(sin)

   //演算値が負にならないように調整
   difference = true_sin_v - sin_v;
   if(difference < 0)
	   difference = sin_v - true_sin_v;

   // decimal (double) to binary
   dec2bin_for_double(difference, bitwidth, fractional_number, c_out);

   // decimal output
   printf("difference sin (%f) = %.10e\n",theta ,difference);

   // binary output
   B_D_output(bitwidth, fractional_number, c_out);

   //-------------------------------------------------------------------
   //演算値(cos)

   // decimal (double) to binary
   dec2bin_for_double(cos_v, bitwidth, fractional_number, c_out);

   // decimal output
   printf("cos (%f) = %.10e\n",theta ,cos_v);

   // binary output
   B_D_output(bitwidth, fractional_number, c_out);

   //-------------------------------------------------------------------
   //真値(cos)
   true_cos_v = cos(theta);

   // decimal (double) to binary
   dec2bin_for_double(true_cos_v, bitwidth, fractional_number, c_out);

   // decimal output
   printf("true cos (%f) = %.10e\n",theta ,true_cos_v);

   // binary output
   B_D_output(bitwidth, fractional_number, c_out);

   //-------------------------------------------------------------------
   //真値と演算値の差(cos)

   //演算値が負にならないように調整
   difference = true_cos_v - cos_v;
   if(difference < 0)
	   difference = cos_v - true_cos_v;

   // decimal (double) to binary
   dec2bin_for_double(difference, bitwidth, fractional_number, c_out);

   // decimal output
   printf("difference cos (%f) = %.10e\n",theta ,difference);

   // binary output
   B_D_output(bitwidth, fractional_number, c_out);
   //-------------------------------------------------------------------

   //他パラメータ確認用
   //-------------------------------------------------------------------
   //αの固定小数点
   //decimal (double) to binary
   //dec2bin_for_double(atanTable[14], bitwidth, fractional_number, c_out);
   //binary output
   //B_D_output(bitwidth, fractional_number, c_out);
   //-------------------------------------------------------------------
   //1/Kの固定小数点
   //decimal (double) to binary
   //l = 1/K;
   //dec2bin_for_double(l, bitwidth, fractional_number, c_out);
   //binary output
   //B_D_output(bitwidth, fractional_number, c_out);
   //-------------------------------------------------------------------
   //decimal (double) to binary
   //dec2bin_for_double(theta, bitwidth, fractional_number, c_out);
   //binary output
   //B_D_output(bitwidth, fractional_number, c_out);
   //-------------------------------------------------------------------
   //deltaの固定小数点
   //decimal (double) to binary
   //i /= 8;
   //dec2bin_for_double(i, bitwidth, fractional_number, c_out);
   //binary output
   //B_D_output(bitwidth, fractional_number, c_out);
   //-------------------------------------------------------------------

   getchar();										//visualStudio用コンソール待機

}

void dec2bin_for_double(double dn, int bitwidth, int fractional_number, int *bn)
{
	int dn_int;
	double shift_v;
	int i;
	int temp_bn;

	//fractional_numberの値だけ左シフトを行い，
	//shift_vに整数値をdouble型として代入する
	shift_v = (double)(1 << fractional_number);

	//演算結果dnとshift_vを乗算(fractional_number精度)
	//またint型より値は1,0で表現されている
	dn_int = (int)(dn * shift_v);

	//bitwidth分繰り返し二進数変換を行う
	for(i = 0; i < bitwidth; i++)
	{
		//int型で表現された1,0とANDをとる
		temp_bn = dn_int & 0x1;

		bn[i] = temp_bn;
		//i=0が最大値が入っていくため表示は逆から格納
		dn_int = (dn_int >> 1);
	}
}


void B_D_output(int bitwidth, int fractional_number, int *c_out){
	int j;
	

	// binary output
  	 for( j=bitwidth-1; j>=0; j-- )
  	 {
   	    printf("%1d",c_out[j]);
	
		  if (j == fractional_number)
		   {
			   printf(".");
		   }
   	}	
  	 printf("\n");


}