(C语言)最小二乘法的曲线拟合.doc
文本预览下载声明
/*最小二乘法的曲线拟合*/
#includestdio.h
#includemath.h
#includestdlib.h
#define max 100
void main()
{
int i,j,k,m,N,mi;
float mx,temp;
float X[max][max],Y[max],x[max],y[max],a[max];
FILE *fp1;
if((fp1=fopen(in1.txt,r))==NULL) /*输入拟合曲线的次数m以及已知的数据组数N*/
{
printf(Cant open this file!\n);
exit(0);
}
for(i=0;i2;i++)
fscanf(fp1,%d%d,m,N);
fclose(fp1);
FILE *fp2;
if((fp2=fopen(in2.txt,r))==NULL) /*输入已知的N组数据坐标*/
{
printf(Cant open this file!\n);
exit(0);
}
for(i=0;iN;i++)
fscanf(fp2,%f%f,x[i],y[i]);
fclose(fp2);
for(i=0;i=m;i++) /*由给定的点得系数矩阵*/
for(j=i;j=m;j++)
{
temp=0;
for(k=0;kN;k++)
temp=temp+pow(x[k],(i+j));
X[i][j]=temp;
X[j][i]=X[i][j];
}
for(i=0;i=m;i++) /*由给定的点得右端矩阵列向量*/
{
temp=0;
for(k=0;kN;k++)
temp=temp+y[k]*pow(x[k],i);
Y[i]=temp;
}
for(j=0;jm;j++) /*列主元素消去法解该线性方程组,得拟合曲线多项式的系数*/
{
for(i=j+1,mi=j,mx=fabs(X[j][j]);i=m;i++)
if(fabs(X[i][j])mx)
{
mi=i;
mx=fabs(X[i][j]);
}
if(jmi)
{
temp=Y[j];
Y[j]=Y[mi];
Y[mi]=temp;
for(k=j;k=m;k++)
{
temp=X[j][k];
X[j][k]=X[mi][k];
X[mi][k]=temp;
}
}
for(i=j+1;i=m;i++)
{
temp=-X[i][j]/X[j][j];
Y[i]+=Y[j]*temp;
for(k=j;k=m;k++)
X[i][k]+=X[j][k]*temp;
}
}
a[m]=Y[m]/X[m][m];
for(i=m-1;i=0;i--)
{
a[i]=Y[i];
for(j=i+1;j=m;j++)
a[i]-=X[i][j]*a[j];
a[i]/=X[i][i];
}
FILE *fp3; /*输出拟合所得的曲线方程*/
fp3=fopen(out.txt,w);
fprintf(fp3,P(x)=(%f),a[0]);
for(i=1;i=m;i++)
fprintf(fp3,+(%f)*x^%d,a[i],i);
fclose(fp3);
}
显示全部