圖像矯正-基於opencv實現

一、引言

        上篇文章中四種方法對圖像進行傾角矯正都非常有效。Hough變換和Radon相似，其抗干擾能力比較強，但是運算量大，程序執行慢，其改進方法為：我們可以不對整幅圖像進行操作，可以在圖像中選取一塊（必須含有一條與傾角有關的直線）進行操作，從而減小運算量。這里Hough變換法和Radon變換法進行傾角檢測的最大精度為1度。它們的優點是可以計算有斷點的直線的傾角。最小二乘法的優點就是運算量小，但是其抗干擾能力比較差，容易受到噪聲的影響。兩點法雖然理論簡單，但由於采樣點比較多而且這些點服從隨機分布，計算均值后能有效抑制干擾，實驗表明其矯正效果很好，最大精度可以明顯小於1度，而且計算量也很小。最小二乘法和兩點法不能計算有斷點的直線傾角，這是這兩種方法的缺點。

二、基於opencv的圖像矯正實現

         對圖像進行旋轉矯正，關鍵是獲取旋轉角度是多少，在獲取旋轉角度后，可以用仿射變換對圖像進行矯正。本文是基於opencv的houghline變換實現的圖像旋轉角度獲取，具體代碼為：

        _grayimage = cv2.cvtColor(self._srcimage,cv2.COLOR_RGB2GRAY)
        _cannyimage = cv2.Canny(_grayimage,CANNY_LOW_THRESHOLD, CANNY_HIGH_THRESHOLD, apertureSize=3)
        lines = cv2.HoughLinesP(_cannyimage,1,np.pi/180,160,minLineLength=200, maxLineGap=180)
        
#        尋找長度最長的線
        distance = []
        for line in lines:
            x1,y1,x2,y2 = line[0]
            dis = np.sqrt(pow((x2-x1),2)+pow((y2-y1),2))
            distance.append(dis)
        max_dis_index = distance.index(max(distance))        
        max_line = lines[max_dis_index]
        x1,y1,x2,y2 = max_line[0]

#       獲取旋轉角度
        angle = cv2.fastAtan2((y2-y1),(x2-x1))

根據hough變化獲取旋轉角度后，根據仿射矩陣進行變換，進而得到矯正后的圖像，具體代碼為：

        centerpoint = (self._srcimage.shape[1]/2,self._srcimage.shape[0]/2)
        rotate_mat = cv2.getRotationMatrix2D(centerpoint,angle,1.0)         #獲取旋轉矩陣
        correct_image = cv2.warpAffine(self._srcimage,rotate_mat,(self._srcimage.shape[1],self._srcimage.shape[0]),borderValue =(255,255,255) )

三、霍夫線變換

         上一篇文章已經針對霍夫變換有了一個介紹，本章將結合opencv中的函數進行詳細講解。

         在使用霍夫線變換之前, 首先要對圖像進行邊緣檢測的處理，也即霍夫線變換的直接輸入只能是邊緣二值圖像.

        OpenCV支持三種不同的霍夫線變換，它們分別是：標准霍夫變換(Standard Hough Transform，SHT)和多尺度霍夫變換（Multi-Scale Hough Transform，MSHT）累計概率霍夫變換(Progressive Probabilistic Hough Transform ，PPHT)。

          其中，多尺度霍夫變換（MSHT）為經典霍夫變換（SHT）在多尺度下的一個變種。累計概率霍夫變換(PPHT）算法是標准霍夫變換（SHT）算法的一個改進，它在一定的范圍內進行霍夫變換，計算單獨線段的方向以及范圍，從而減少計算量，縮短計算時間。之所以稱PPHT為“概率”的，是因為並不將累加器平面內的所有可能的點累加，而只是累加其中的一部分，該想法是如果峰值如果足夠高，只用一小部分時間去尋找它就夠了。這樣猜想的話，可以實質性地減少計算時間。

           在OpenCV中，我們可以用HoughLines函數來調用標准霍夫變換SHT和多尺度霍夫變換MSHT。

           而HoughLinesP函數用於調用累計概率霍夫變換PPHT。累計概率霍夫變換執行效率很高，所有相比於HoughLines函數，我們更傾向於使用HoughLinesP函數。

3.1 霍夫線變換原理

        一條直線在圖像二維空間可由兩個變量表示. 如:

<1>在笛卡爾坐標系: 可由參數: 斜率和截距(m,b) 表示。

<2>在極坐標系: 可由參數: 極徑和極角表示。

         對於霍夫變換, 我們將采用第二種方式極坐標系來表示直線. 因此, 直線的表達式可為:

化簡便可得到:

這就意味着每一對代表一條通過點的直線。

         如果對於一個給定點我們在極坐標對極徑極角平面繪出所有通過它的直線, 將得到一條正弦曲線. 例如, 對於給定點X_0= 8 和Y_0= 6 我們可以繪出下圖 (在平面):

         只繪出滿足下列條件的點 和  .

         我們可以對圖像中所有的點進行上述操作. 如果兩個不同點進行上述操作后得到的曲線在平面相交, 這就意味着它們通過同一條直線. 例如,接上面的例子我們繼續對點 和點 繪圖, 得到下圖:

            這三條曲線在平面相交於點 (0.925, 9.6), 坐標表示的是參數對 或者是說點, 點和點組成的平面內的的直線。

         以上的說明表明，一般來說, 一條直線能夠通過在平面 尋找交於一點的曲線數量來檢測。而越多曲線交於一點也就意味着這個交點表示的直線由更多的點組成. 一般來說我們可以通過設置直線上點的閾值來定義多少條曲線交於一點我們才認為檢測到了一條直線。

         這就是霍夫線變換要做的. 它追蹤圖像中每個點對應曲線間的交點. 如果交於一點的曲線的數量超過了閾值, 那么可以認為這個交點所代表的參數對在原圖像中為一條直線。

關於霍夫變換的詳細解釋，可以看此英文頁面：http://homepages.inf.ed.ac.uk/rbf/HIPR2/hough.htm

3.2 HoughLinesP（）函數詳解

因為opencv內部基於c++寫的，在函數詳解中以c++源碼為主。

C++: void HoughLinesP(InputArray image, OutputArray lines, double rho, double theta, int threshold, double minLineLength=0, double maxLineGap=0 )

• 第一個參數，InputArray類型的image，輸入圖像，即源圖像，需為8位的單通道二進制圖像，可以將任意的源圖載入進來后由函數修改成此格式后，再填在這里。
• 第二個參數，InputArray類型的lines，經過調用HoughLinesP函數后后存儲了檢測到的線條的輸出矢量，每一條線由具有四個元素的矢量(x_1,y_1, x_2, y_2）  表示，其中，(x_1, y_1)和(x_2, y_2) 是是每個檢測到的線段的結束點。
• 第三個參數，double類型的rho，以像素為單位的距離精度。另一種形容方式是直線搜索時的進步尺寸的單位半徑。
• 第四個參數，double類型的theta，以弧度為單位的角度精度。另一種形容方式是直線搜索時的進步尺寸的單位角度。
• 第五個參數，int類型的threshold，累加平面的閾值參數，即識別某部分為圖中的一條直線時它在累加平面中必須達到的值。大於閾值threshold的線段才可以被檢測通過並返回到結果中。
• 第六個參數，double類型的minLineLength，有默認值0，表示最低線段的長度，比這個設定參數短的線段就不能被顯現出來。即當檢測出的直線長度大於minLinLength，才認為這是一條直線。
• 第七個參數，double類型的maxLineGap，有默認值0，允許將同一行點與點之間連接起來的最大的距離。即假若直線從中間某處斷開，那么所允許的缺口的最大長度。

//-----------------------------------【頭文件包含部分】---------------------------------------  
//      描述：包含程序所依賴的頭文件  
//----------------------------------------------------------------------------------------------   
#include <opencv2/opencv.hpp>  
#include <opencv2/imgproc/imgproc.hpp>  
  
//-----------------------------------【命名空間聲明部分】---------------------------------------  
//      描述：包含程序所使用的命名空間  
//-----------------------------------------------------------------------------------------------   
using namespace cv;  
//-----------------------------------【main( )函數】--------------------------------------------  
//      描述：控制台應用程序的入口函數，我們的程序從這里開始  
//-----------------------------------------------------------------------------------------------  
int main( )  
{  
    //【1】載入原始圖和Mat變量定義     
    Mat srcImage = imread("1.jpg");  //工程目錄下應該有一張名為1.jpg的素材圖  
    Mat midImage,dstImage;//臨時變量和目標圖的定義  
  
    //【2】進行邊緣檢測和轉化為灰度圖  
    Canny(srcImage, midImage, 50, 200, 3);//進行一此canny邊緣檢測  
    cvtColor(midImage,dstImage, CV_GRAY2BGR);//轉化邊緣檢測后的圖為灰度圖  
  
    //【3】進行霍夫線變換  
    vector<Vec4i> lines;//定義一個矢量結構lines用於存放得到的線段矢量集合  
    HoughLinesP(midImage, lines, 1, CV_PI/180, 80, 50, 10 );  
  
    //【4】依次在圖中繪制出每條線段  
    for( size_t i = 0; i < lines.size(); i++ )  
    {  
        Vec4i l = lines[i];  
        line( dstImage, Point(l[0], l[1]), Point(l[2], l[3]), Scalar(186,88,255), 1, CV_AA);  
    }  
  
    //【5】顯示原始圖    
    imshow("【原始圖】", srcImage);    
  
    //【6】邊緣檢測后的圖   
    imshow("【邊緣檢測后的圖】", midImage);    
  
    //【7】顯示效果圖    
    imshow("【效果圖】", dstImage);    
  
    waitKey(0);    
  
    return 0;    
}

運行截圖：

來一張大圖：

3.3 源碼解釋

void cv::HoughLinesP( InputArray _image,OutputArray _lines,  
                      double rho, double theta,int threshold,  
                      double minLineLength,double maxGap )  
{  
   Ptr<CvMemStorage> storage = cvCreateMemStorage(STORAGE_SIZE);  
   Mat image = _image.getMat();  
   CvMat c_image = image;  
    CvSeq*seq = cvHoughLines2( &c_image, storage, CV_HOUGH_PROBABILISTIC,  
                    rho, theta, threshold,minLineLength, maxGap );  
   seqToMat(seq, _lines);  
}

CV_IMPL CvSeq*  
cvHoughLines2( CvArr* src_image, void*lineStorage, int method,  
               double rho, double theta, intthreshold,  
               double param1, double param2 )  
{  
   CvSeq* result = 0;  
   
   CvMat stub, *img = (CvMat*)src_image;  
   CvMat* mat = 0;  
   CvSeq* lines = 0;  
   CvSeq lines_header;  
   CvSeqBlock lines_block;  
   int lineType, elemSize;  
   int linesMax = INT_MAX;  
   int iparam1, iparam2;  
   
   img = cvGetMat( img, &stub );  
   
   if( !CV_IS_MASK_ARR(img))  
       CV_Error( CV_StsBadArg, "The source image must be 8-bit,single-channel" );  
   
   if( !lineStorage )  
       CV_Error( CV_StsNullPtr, "NULL destination" );  
   
   if( rho <= 0 || theta <= 0 || threshold <= 0 )  
       CV_Error( CV_StsOutOfRange, "rho, theta and threshold must bepositive" );  
   
   if( method != CV_HOUGH_PROBABILISTIC )  
    {  
       lineType = CV_32FC2;  
       elemSize = sizeof(float)*2;  
    }  
   else  
    {  
       lineType = CV_32SC4;  
       elemSize = sizeof(int)*4;  
    }  
   
   if( CV_IS_STORAGE( lineStorage ))  
    {  
       lines = cvCreateSeq( lineType, sizeof(CvSeq), elemSize,(CvMemStorage*)lineStorage );  
    }  
   else if( CV_IS_MAT( lineStorage ))  
    {  
       mat = (CvMat*)lineStorage;  
   
       if( !CV_IS_MAT_CONT( mat->type ) || (mat->rows != 1 &&mat->cols != 1) )  
           CV_Error( CV_StsBadArg,  
           "The destination matrix should be continuous and have a single rowor a single column" );  
   
       if( CV_MAT_TYPE( mat->type ) != lineType )  
           CV_Error( CV_StsBadArg,  
           "The destination matrix data type is inappropriate, see themanual" );  
   
       lines = cvMakeSeqHeaderForArray( lineType, sizeof(CvSeq), elemSize,mat->data.ptr,  
                                        mat->rows + mat->cols - 1, &lines_header, &lines_block );  
       linesMax = lines->total;  
       cvClearSeq( lines );  
    }  
   else  
       CV_Error( CV_StsBadArg, "Destination is not CvMemStorage* norCvMat*" );  
   
   iparam1 = cvRound(param1);  
   iparam2 = cvRound(param2);  
   
   switch( method )  
    {  
   case CV_HOUGH_STANDARD:  
         icvHoughLinesStandard( img, (float)rho,  
               (float)theta, threshold,lines, linesMax );  
         break;  
   case CV_HOUGH_MULTI_SCALE:  
         icvHoughLinesSDiv( img, (float)rho, (float)theta,  
                threshold, iparam1, iparam2,lines, linesMax );  
         break;  
   case CV_HOUGH_PROBABILISTIC:  
         icvHoughLinesProbabilistic( img, (float)rho, (float)theta,  
                threshold, iparam1, iparam2,lines, linesMax );  
         break;  
   default:  
       CV_Error( CV_StsBadArg, "Unrecognized method id" );  
    }  
   
   if( mat )  
    {  
       if( mat->cols > mat->rows )  
           mat->cols = lines->total;  
       else  
           mat->rows = lines->total;  
    }  
   else  
       result = lines;  
   
   return result;  
}

先看Hough檢測直線的代碼，cvHoughLines2也只不過是個對不同Hough方法的封裝，下面是該函數中的部分代碼，選擇不同的Hough變換方法。

/* 這段注釋解釋了函數各個參數的作用
Here image is an input raster;
step is it's step; size characterizes it's ROI;
rho and theta are discretization steps (in pixels and radians correspondingly).
threshold is the minimum number of pixels in the feature for it
to be a candidate for line. lines is the output
array of (rho, theta) pairs. linesMax is the buffer size (number of pairs).
Functions return the actual number of found lines.
*/
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 );  // 確定的直線入隊列輸出
    }
}

四、圖像仿射變換

        幾何變換可以看成圖像中物體（或像素）空間位置改變，或者說是像素的移動。

        幾何運算需要空間變換和灰度級差值兩個步驟的算法，像素通過變換映射到新的坐標位置，新的位置可能是在幾個像素之間，即不一定為整數坐標。這時就需要灰度級差值將映射的新坐標匹配到輸出像素之間。最簡單的插值方法是最近鄰插值，就是令輸出像素的灰度值等於映射最近的位置像素，該方法可能會產生鋸齒。這種方法也叫零階插值，相應比較復雜的還有一階和高階插值。

4.1 仿射變換原理

仿射變換的功能是從二維坐標到二維坐標之間的線性變換，且保持二維圖形的“平直性”和“平行性”。仿射變換可以通過一系列的原子變換的復合來實現，包括平移，縮放，翻轉，旋轉和剪切。

空間變換對應矩陣的仿射變換。一個坐標通過函數變換的新的坐標位置：

所以在程序中我們可以使用一個2*3的數組結構來存儲變換矩陣：

4.1、opencv的圖像變換函數

void cvWarpAffine( const CvArr* src,//輸入圖像  CvArr* dst, //輸出圖像  const CvMat* map_matrix, //2*3的變換矩陣  int flags=CV_INTER_LINEAR+CV_WARP_FILL_OUTLIERS, //插值方法的組合  CvScalar fillval=cvScalarAll(0) //用來填充邊界外的值  );

cv2.warpAffine(src, M, dsize[, dst[, flags[, borderMode[, borderValue]]]])

非關鍵字參數有src, M, dsize，分別表示源圖像，變換矩陣，變換后的圖像的長寬。

下面介紹一些典型的仿射變換：

（1）平移，將每一點移到到（x+t , y+t)，變換矩陣為

(2)縮放變換  將每一點的橫坐標放大或縮小s_x倍，縱坐標放大（縮小）到s_y倍，變換矩陣為

（3）旋轉變換原點：目標圖形圍繞原點順時針旋轉Θ 弧度，變換矩陣為

(4) 旋轉變換  ：目標圖形以（x , y ）為軸心順時針旋轉θ弧度，變換矩陣為

相當於兩次平移與一次原點旋轉變換的復合，即先將軸心（x,y)移到到原點，然后做旋轉變換，最后將圖片的左上角置為圖片的原點,即

有的人可能會說為什么這么復雜呢，那是因為在opencv的圖像處理中，所有對圖像的處理都是從原點進行的，而圖像的原點默認為圖像的左上角，而我們對圖像作旋轉處理時一般以圖像的中點為軸心，因此就需要做如下處理。

其中的變換矩陣由：

M=cv2.getRotationMatrix2D(rotate_center, degree, scale)

rotate_center為一個2元的元組，表示旋轉中心坐標，degree表示逆時針旋轉的角度，scale表示縮放的比例，獲取。

它得到的矩陣是：

其中α = scale * cos( angle ) , β = scale  * sing( angle )  , ( center.x , center.y ) 表示旋轉軸心

五、測試代碼

# -*- coding: utf-8 -*-
"""
Created on Thu Jun 29 14:23:56 2017

@author: Administrator
"""

import cv2
import numpy as np
import os

CANNY_LOW_THRESHOLD = 50                # canny算法低閾值
CANNY_HIGH_THRESHOLD = 150              #canny算法高閾值
HOUGH_DELTARHO = 1                      #hough檢測步長
HOUGH_DELTA_THETA = 100
TESTIMAGESDIR = 'C:\\Users\\Administrator\\Desktop\\OCR\\retate'

class RotateByHough(object):
    def __init__(self,image):
        self._srcimage = image

#獲取圖像背景像素值        
    def get_background_pix(self,):
        width = self._srcimage.shape[1]
        height = self._srcimage.shape[0]
        left_coners = self._srcimage[0:9,0:9]
        right_coners = self._srcimage[height-10:height-1,width-10:width-1]
        left_array = np.array(left_coners).reshape(81,1,3)
        left_point = np.mean(left_array,0)        
        right_point = np.mean(np.array(right_coners).reshape(81,1,3),0)
#        print(left_point)
        final = np.array((left_point[0],right_point[0]),np.float32).reshape(2,1,3)        
        final_point = np.mean(final,0)
        
        self._point_scale = (int(final_point[0][0]),int(final_point[0][1]),int(final_point[0][2]))
#        print(self._point_scale)

#       輸入矯正圖像的寬和高        
    def detect_hough_line(self):
        _grayimage = cv2.cvtColor(self._srcimage,cv2.COLOR_RGB2GRAY)
        _cannyimage = cv2.Canny(_grayimage,CANNY_LOW_THRESHOLD, CANNY_HIGH_THRESHOLD, apertureSize=3)
        lines = cv2.HoughLinesP(_cannyimage,1,np.pi/180,160,minLineLength=200, maxLineGap=180)
        
#        尋找長度最長的線
        distance = []
        for line in lines:
            x1,y1,x2,y2 = line[0]
            dis = np.sqrt(pow((x2-x1),2)+pow((y2-y1),2))
            distance.append(dis)
        max_dis_index = distance.index(max(distance))        
        max_line = lines[max_dis_index]
        x1,y1,x2,y2 = max_line[0]

#       獲取旋轉角度
        angle = cv2.fastAtan2((y2-y1),(x2-x1))
        centerpoint = (self._srcimage.shape[1]/2,self._srcimage.shape[0]/2)
        rotate_mat = cv2.getRotationMatrix2D(centerpoint,angle,1.0)         #獲取旋轉矩陣
        correct_image = cv2.warpAffine(self._srcimage,rotate_mat,(self._srcimage.shape[1],self._srcimage.shape[0]),borderValue =(255,255,255) )        
        return correct_image

def get_files(images_dir):
    images = []
    for _file in os.listdir(images_dir):
        filename = os.path.join(images_dir, _file)
        images.append(filename)
        
    return images

#==============================================================================
# def key_board_event(event):
# #    監聽鍵盤輸入值
# 
#     print("Key :", event.Key)
#     if event.Key == 'Down':
#         print('next file')        
#         i = i + 1
#         if i > len(image_files):
#             i = len(image_files)        
#         print(i)
#         print(image_files[i])
#     if event.Key == 'Up':
#         print('pre file')
#         i = i -1
#         if i < 0:
#             i = 0     
#     
#==============================================================================    
def main():
    image_files = get_files(TESTIMAGESDIR)
    files_num = len(image_files)
    while True:
        num = int(input("Enter your input:"))        
        if num > 50:
            break;
        if num < 0 or num > files_num:
            print('resume load')
        else:
            file_name = image_files[num]
            print(file_name)
            image = cv2.imread(file_name)
            roteimage_ins = RotateByHough(image)
            roteimage = roteimage_ins.detect_hough_line()
            width = roteimage.shape[1]
            height = roteimage.shape[0]

            for i in range(0, height,  50):
                for j in range(0, width, 50):
                    cv2.line(roteimage, (0,i),(width, i),(255, 0, 255), 1)                    
            cv2.imshow('test', roteimage)
            cv2.waitKey(0)
#==============================================================================
#     hm = pyHook.HookManager()             #設置一個鈎子管理對象
#     hm.KeyDown = key_board_event          #監聽所有鍵盤事件
#     hm.HookKeyboard()                     #設置鍵盤鈎子
#==============================================================================    
    
if __name__ == "__main__":
    main()