Hough变换原理

博客里面编辑数学公式太麻烦了，还是截图吧，自己写的，有不严谨的地方还望多多指教。

标准hough变换源码：

static void icvHoughLinesStandard(const CvMat* img, float rho, float theta,int threshold, CvSeq *lines, int linesMax)
{
    cv::AutoBuffer<int> _accum, _sort_buf;    // _accum：计数用数组，_sort_buf，排序用数组
    cv::AutoBuffer<float> _tabSin, _tabCos;   // 提前计算sin与cos值，避免重复计算带来的计算性能下降

    const uchar* image;
    int step, width, height;
    int numangle, numrho;
    int total = 0;
    float ang;
    int r, n;
    int i, j;
    float irho = 1 / rho;   // rho指像素精度，常取1，因此irho常为1
    double scale;

    CV_Assert(CV_IS_MAT(img) && CV_MAT_TYPE(img->type) == CV_8UC1);

    image = img->data.ptr;
    step = img->step;
    width = img->cols;
    height = img->rows;

    numangle = cvRound(CV_PI / theta);  // 根据th精度计算th维度的长度
    numrho = cvRound(((width + height) * 2 + 1) / rho);  // 根据r精度计算r维度的长度

    _accum.allocate((numangle + 2) * (numrho + 2));
    _sort_buf.allocate(numangle * numrho);
    _tabSin.allocate(numangle);
    _tabCos.allocate(numangle);
    int *accum = _accum, *sort_buf = _sort_buf;
    float *tabSin = _tabSin, *tabCos = _tabCos;

    memset(accum, 0, sizeof(accum[0]) * (numangle + 2) * (numrho + 2));

    for (ang = 0, n = 0; n < numangle; ang += theta, n++)   // 计算三角函数表，避免重复计算
    {
        tabSin[n] = (float)(sin(ang) * irho);
        tabCos[n] = (float)(cos(ang) * irho);
    }

    // stage 1. fill accumulator 
    for (i = 0; i < height; i++)
    for (j = 0; j < width; j++)
    {
        if (image[i * step + j] != 0)
        for (n = 0; n < numangle; n++)
        {
            r = cvRound(j * tabCos[n] + i * tabSin[n]);  // Hough极坐标变换式
            r += (numrho - 1) / 2;
            accum[(n + 1) * (numrho + 2) + r + 1]++;  // 计数器统计
        }
    }

    // stage 2. find local maximums
    for (r = 0; r < numrho; r++)
    for (n = 0; n < numangle; n++)
    {
        int base = (n + 1) * (numrho + 2) + r + 1;
        if (accum[base] > threshold       // 大于阈值，且是局部极大值
            && accum[base] > accum[base - 1]
            && accum[base] >= accum[base + 1]
            && accum[base] > accum[base - numrho - 2]
            && accum[base] >= accum[base + numrho + 2])
            sort_buf[total++] = base;
    }

    // stage 3. sort the detected lines by accumulator value
    icvHoughSortDescent32s(sort_buf, total, accum);

    // stage 4. store the first min(total,linesMax) lines to the output buffer
    linesMax = MIN(linesMax, total);  // linesMax是输入参数，表示最多输出多少个直线参数
    scale = 1. / (numrho + 2);
    for (i = 0; i < linesMax; i++)
    {
        CvLinePolar line;           // 输出结构，就是(r,theta)
        int idx = sort_buf[i];
        int n = cvFloor(idx*scale) - 1;
        int r = idx - (n + 1)*(numrho + 2) - 1;
        line.rho = (r - (numrho - 1)*0.5f) * rho;
        line.angle = n * theta;
        cvSeqPush(lines, &line);  // 确定的直线入队列输出
    }
}