最小二乘法代码(matlab 中polyfit函数原型)

最新推荐文章于 2024-08-09 09:03:13 发布

my_LoveLetter

最新推荐文章于 2024-08-09 09:03:13 发布

阅读量1w

点赞数 2

分类专栏：编程文章标签： matlab ini qt 360

本文链接：https://blog.csdn.net/zhaobwsir/article/details/7078700

版权

编程专栏收录该内容

2 篇文章 0 订阅

订阅专栏


//函数功能:进行最小二乘曲线拟合(拟合y=a0+a1*x+a2*x^2+……+a[poly_n]*x^poly_n),计算出对应的系数a?
//参数说明:
//      n:      给定数据点的个数
//      x[]:    存放给定n个数据点的X坐标
//      y[]:    存放给定n个数据点的Y坐标
//      poly_n: 拟合多项式的系数,表示多项式的最高次数
void polyfit(int n,double x[],double y[],int poly_n,double a[])
{
    void gauss_solve(int n,double A[],double x[],double b[]);
    int i,j;                                            double *tempx,*tempy,*sumxx,*sumxy,*ata;
    tempx=(double *)calloc(double(n),sizeof(double));   sumxx=(double *)calloc(poly_n*2+1,sizeof(double));
    tempy=(double *)calloc(n,sizeof(double));           sumxy=(double *)calloc(poly_n+1,sizeof(double));
    ata=(double *)calloc((poly_n+1)*(poly_n+1),sizeof(double));
    for (i=0;i<n;i++)
    {
        tempx[i]=1;
        tempy[i]=y[i];
    }
    for (i=0;i<2*poly_n+1;i++)
        for (sumxx[i]=0,j=0;j<n;j++)
        {
            sumxx[i]+=tempx[j];
            tempx[j]*=x[j];
        }
    for (i=0;i<poly_n+1;i++)
        for (sumxy[i]=0,j=0;j<n;j++)
        {
            sumxy[i]+=tempy[j];
            tempy[j]*=x[j];
        }
    for (i=0;i<poly_n+1;i++)
        for (j=0;j<poly_n+1;j++)
            ata[i*(poly_n+1)+j]=sumxx[i+j];
    gauss_solve(poly_n+1,ata,a,sumxy);
    free(tempx);    free(sumxx);
    free(tempy);    free(sumxy);
    free(ata);
}
void gauss_solve(int n,double A[],double x[],double b[]){
    int i,j,k,r;
    double max;
    for (k=0;k<n-1;k++)
    {
        max=fabs(A[k*n+k]); /*find maxmum*/
        r=k;
        for (i=k+1;i<n-1;i++)
        if (max<fabs(A[i*n+i]))
        {
            max=fabs(A[i*n+i]);
            r=i;
        }
        if (r!=k)
            for (i=0;i<n;i++)         /*change array:A[k]&A[r]  */
            {
                max=A[k*n+i];
                A[k*n+i]=A[r*n+i];
                A[r*n+i]=max;
            }
            max=b[k];                    /*change array:b[k]&b[r]     */
            b[k]=b[r];
            b[r]=max;
            for (i=k+1;i<n;i++)
            {
                for (j=k+1;j<n;j++)
                    A[i*n+j]-=A[i*n+k]*A[k*n+j]/A[k*n+k];
                b[i]-=A[i*n+k]*b[k]/A[k*n+k];
            }
    }
    for (i=n-1;i>=0;x[i]/=A[i*n+i],i--)
        for (j=i+1,x[i]=b[i];j<n;j++)
            x[i]-=A[i*n+j]*x[j];
}
#define yy y3
void main(void)
{
    QString y1="1935 1877 1826 1804 1740 1665 1648 1611 1556 1497 1477 1443 1371 1320 1281 1246 1185 1118 1082 1034 988 926 870 837 832 769 723 681 667 621 545 517 472 436 370";
    QString y2="-360 -400 -460 -528 -558 -635 -684 -752 -794 -854 -879 -948 -1010 -1020 -1048 -1139 -1179 -1210 -1242 -1328 -1354 -1378 -1452 -1496 -1523 -1559 -1640 -1696 -1707 -1762 -1812 -1880 -1880 -1952 -1958 -2001 -1967 -1969 -1938 -1889 -1874 -1823 -1802 -1805 -1736 -1722 -1713 -1642 -1633 -1634 -1564 -1549 -1512 -1490 -1443 -1404 -1374 -1343 -1305 -1246 -1235 -1235 -1172 -1142 -1123 -1102 -1065 -1029 -999 -971 -929 -888 -900 -871 -815 -770 -779 -747 -680 -672 -656 -604 -551 -542 -518 -466 -404 -385 -380 -344";
    QString y3="364 384 396 440 499 500 542 610 648 675 693 751 791 802 818 897 953 972 962 952 899 851 832 786 743 731 700 681 649 610 573 561 535 497 488 444 412 379 380 365";
    QString y4="-348 -376 -380 -417 -460 -460 -484 -514 -544 -582 -602 -622 -666 -683 -726 -756 -768 -823 -848 -856 -877 -904 -920 -933 -948 -968 -983 -1010 -1066 -1060 -1070 -1108 -1151 -1175 -1166 -1183 -1229 -1273 -1281 -1288 -1311 -1369 -1380 -1373 -1418 -1469 -1521 -1517 -1528 -1565 -1607 -1612 -1607 -1658 -1699 -1711 -1723 -1761 -1789 -1833 -1842 -1846 -1844 -1841 -1807 -1756 -1719 -1684 -1677 -1599 -1559 -1526 -1513 -1466 -1412 -1372 -1359 -1303 -1261 -1213 -1210 -1153 -1072 -1072 -1038 -999 -941 -916 -916 -853 -796 -767 -722 -666 -605 -581 -520 -472 -432 -427 -356 -330";
 void polyfit(int n,double x[],double y[],int poly_n,double a[]);
 int i,j,m,n=7,poly_n=3;
 double x[200],y[200];
 double a[5];
 void polyfit(int n,double *x,double *y,int poly_n,double a[]);
 
 for(i=0;i<yy.split(' ').size();i++)
 {
    x[i]=i+1;
    y[i]=yy.split(' ').at(i).toDouble();
 }
 n=yy.split(' ').size();
 polyfit(n,x,y,poly_n,a);
 for (i=0;i<poly_n+1;i++)/*这里是升序排列，Matlab是降序排列*/
     qDebug()<<a[i];
}
在qt软件中经过验证,输出代码和在matlab相同数据进行仿真时的计算结果一致,相吻合.
转载自:
http://blog.163.com/xing_mu_1/blog/static/661429020097101037114/对里面的部分代码进行了部分修改