// ［このプログラムの目的］カージオイドからハート形曲線への変換3，2009年04月06日（月）

// file name: heart_curve_c.c

#include< stdio.h>
#include< math.h>

void main(void)
{

	double a,pi;// pi は円周率
	double r,f,df;// それぞれ、曲座標の動径、変換後の位相角、および、変換後の位相角の変分
	double t;// 変換前の位相角
	double fmin,fmax;// 位相角の最小値と最大値
	double b,c;// カージオイドからハート形曲線への変換係数
	double p;// 修正関数の定数
	double d;// 縦方向に伸縮する伸縮係数
	int i,imax;
	double x,y;// 修正前の直交座標
	double xx[10001],yy[10001];// メモリ容量の上限に注意

	FILE *fp;

//  定数設定
	pi=3.14159265;
	a=1;

	b=1.8;// 数値は 1.3 から 2.5 の間が良い
	c=0.7;// 0 < beta <1 で、 beta=0.5 付近が良い
	d=1.2;// 0.5 から 1.5 の間が良い
	p=0.1;

//  他のパラメータ設定
	fmin=-pi/2;
	fmax=3*pi/2;
	df=(fmax-fmin)/4000;// 位相角 p のプロット間隔

//  計算実行
	i=0;
	for(f=fmin;f<=fmax;f=f+df)
	{
		i++;

		t=b*sqrt(f+pi/2)-b*sqrt(3*pi/2-f)+(1-b*sqrt(2/pi))*f+b*sqrt(pi/2);
		r=a*(1-sin(t));

		x=r*(1+c*sin(f))*cos(f);
		xx[i]=x;
		y=d*r*(1+c*sin(f))*sin(f);

		if(x>=0)
		{
			yy[i]=y+p*x;
		}
		else
		{
			yy[i]=y-p*x;
		}

			printf("x=%f,y=%f\n",xx[i],yy[i]);
	}

	imax=i;

// 計算データ[ｎ角の星形曲線の(x,y)座標]をテキストファイル（heart_curve_c.txt）に書き込む
	fp=fopen("heart_curve_c.txt","w");
	if(fp==NULL)
		{
		printf("FILE OPEN ERROR\n");
		}
	else
		{
		for(i=1;i<=imax;i++)
			{
			fprintf(fp,"%f,%f\n",xx[i],yy[i]);
			}
		fflush(fp);
		fclose(fp);
		}
	printf("end\n");
}// the end of the program