最小二乘法拟合直线 C++/OpenCV

问题:

　　我们在拥有一系列散列的点(x1,y1),(x2,y2)... (xm,ym),这些点在一条直线附近，通过点拟合直线。

　　我在工程中是要拟合一系列线段，其实一条线段就对应着两个要拟合的点，算法上稍有区别，原理完全一致。

 

思路一：

　　用OpenCV自带的最小二乘法拟合

 

　　函数：cvFitLine()

　　　　void cvFitLine( const CvArr* points,    int dist_type,    double param,  double reps,   double aeps,   float* line );

　　

　　points 

　　2D 或 3D 点集，32-比特整数或浮点数坐标 ，存放输入点。

　　　　Char* cvSeqPush(CvSeq* seq,void* element=NULL)

　　　　功能：存放元素到序列尾部。            

　　　　 　　

     dist_type 
　　      拟合的距离类型.

　　　　dist_type=CV_DIST_L2 (L2): ρ(r)=r2/2 (最小二乘法，我们要用的)

　　　　dist_type=CV_DIST_L1 (L1): ρ(r)=r 
　　　　dist_type=CV_DIST_L12 (L1-L2): ρ(r)=2?[sqrt(1+r2/2) - 1] 
　　　　dist_type=CV_DIST_FAIR (Fair): ρ(r)=C2?[r/C - log(1 + r/C)], C=1.3998 
　　　　dist_type=CV_DIST_WELSCH (Welsch): ρ(r)=C2/2?[1 - exp(-(r/C)2)], C=2.9846 
　　　　dist_type=CV_DIST_HUBER (Huber): ρ(r)= r2/2, if r < C; C?(r-C/2), otherwise; C=1.345  

　　param 
　　对某些距离的数字参数，如果是 0， 则选择某些最优值 

　　reps, aeps 
　　半径 (坐标原点到直线的距离) 和角度的精度，一般设为0.01。 

　　line 
　　输出的直线参数。2D 拟合情况下，它是包含 4 个浮点数的数组 (vx, vy, x0, y0)，其中 (vx, vy) 是线的单位向量而 (x0, y0) 是线上的某个点. 

对 3D 拟合，它是包含 6 个浮点数的数组 (vx, vy, vz, x0, y0, z0)， 其中 (vx, vy, vz) 是线的单位向量，而 (x0, y0, z0) 是线上某点。 

        float *line = new float[4];
        CvMemStorage* storage = cvCreateMemStorage(0);
        
        //往点序列中存放需要参与拟合直线的点
        CvSeq* point_seq = cvCreateSeq( CV_32FC2, sizeof(CvSeq), sizeof(CvPoint2D32f), storage );
        for(int j=0;j<temp.size();j++)
        {
            cvSeqPush(point_seq, &cvPoint2D32f(temp[j].stPot().x,temp[j].stPot().y));
            cvSeqPush(point_seq, &cvPoint2D32f(temp[j].enPot().x,temp[j].enPot().y));
        }

        cvFitLine(point_seq,CV_DIST_L2,0,0.01,0.01,line);   //CV_DIST_L2表示最小二乘法

        //line[0],line[1]为x，y的单位方向向量
        //line[2],line[3]为直线经过某点的X,Y值
        //正好组成参数方程
        //X=line[0]t+line[2]
        //Y=line[1]t+line[3]
        
        //根据line[4]计算自己需要的直线或线段
        Point A(0,line[3]-(line[1]*line[2]/line[0]));
        Point B(src->width,(src->width-line[2])*line[1]/line[0]+line[3]);
        cvLine(src,A,B,CV_RGB(255,0,0),3,8,0);

        cvClearSeq(point_seq);
        cvReleaseMemStorage(&storage);        

思路二：

　　自己写算法，但是需要注意，openCV中原点在屏幕的左上角，与我们数学中的直角坐标系不太一样。     

//最小二乘法拟合直线
void line_fit(vector<MyLine>& h,IplImage* src)
{
    int n = 2*h.size();
    int k;                       //目标直线斜率
    int b;                       //目标直线截距
    vector<MyLine>::iterator it = h.begin();
    int sumx=0,sumy=0,sumxy=0,sumxsq=0;
    while(it != h.end())
    {
        sumx += it->stPot().x;
        sumx += it->enPot().x;
        sumy += it->stPot().y;
        sumy += it->enPot().y;
        sumxy += it->stPot().x*it->stPot().y;
        sumxy += it->enPot().x*it->enPot().y;
        sumxsq += it->stPot().x*it->stPot().x;
        sumxsq += it->enPot().x*it->enPot().x;
        it++;
    }
    
    if(sumxsq == (sumx*sumx/n))
        k = atan2(0,1.0);
    else
        k = (sumxy-((sumx*sumy)/n))/(sumxsq-(sumx*sumx/n));
    b = (sumy-k*sumx)/n;
    cvLine(src,Point(0,b),Point(src->width,k*src->width+b),CV_RGB(255,0,0),1,8,0);
}