图象滤波器 Coherence-enhance shock filter 实现

一：源码

https://github.com/Itseez/opencv/blob/master/samples/python2/coherence.py

def coherence_filter(img, sigma = 11, str_sigma = 11, blend = 0.5, iter_n = 4):    h, w = img.shape[:2]    for i in xrange(iter_n):        print i,        gray = cv2.cvtColor(img, cv2.COLOR_BGR2GRAY)        eigen = cv2.cornerEigenValsAndVecs(gray, str_sigma, 3)        eigen = eigen.reshape(h, w, 3, 2)  # [[e1, e2], v1, v2]        x, y = eigen[:,:,1,0], eigen[:,:,1,1]        gxx = cv2.Sobel(gray, cv2.CV_32F, 2, 0, ksize=sigma)        gxy = cv2.Sobel(gray, cv2.CV_32F, 1, 1, ksize=sigma)        gyy = cv2.Sobel(gray, cv2.CV_32F, 0, 2, ksize=sigma)        gvv = x*x*gxx + 2*x*y*gxy + y*y*gyy        m = gvv < 0        ero = cv2.erode(img, None)        dil = cv2.dilate(img, None)        img1 = ero        img1[m] = dil[m]        img = np.uint8(img*(1.0 - blend) + img1*blend)    print 'done'    return img

二：

opencv  mat.reshape　c是行优先，　matlab 列优先

三：

http://www.rosoo.net/a/201004/9157.html#CornerEigenValsAndVecs

#!/usr/bin/env python'''Texture flow direction estimation.Sample shows how cv2.cornerEigenValsAndVecs function can be usedto estimate image texture flow direction.Usage:    texture_flow.py [<image>]'''import numpy as npimport cv2if __name__ == '__main__':    import sys    try:        fn = sys.argv[1]    except:        fn = 'data/starry_night.jpg'    img = cv2.imread(fn)    if img is None:        print 'Failed to load image file:', fn        sys.exit(1)    gray = cv2.cvtColor(img, cv2.COLOR_BGR2GRAY)    h, w = img.shape[:2]    eigen = cv2.cornerEigenValsAndVecs(gray, 15, 3)    eigen = eigen.reshape(h, w, 3, 2)  # [[e1, e2], v1, v2]    flow = eigen[:,:,2]    vis = img.copy()    vis[:] = (192 + np.uint32(vis)) / 2    d = 12    points =  np.dstack( np.mgrid[d/2:w:d, d/2:h:d] ).reshape(-1, 2)    for x, y in points:       vx, vy = np.int32(flow[y, x]*d)       cv2.line(vis, (x-vx, y-vy), (x+vx, y+vy), (0, 0, 0), 1, cv2.LINE_AA)    cv2.imshow('input', img)    cv2.imshow('flow', vis)    cv2.waitKey()

int main (int argc, char** argv){    cv::TickMeter tm;    tm.start();    cv::Mat img = cv::imread(argv[1]);    cv::Mat gray = cv::Mat();    cv::cvtColor(img, gray, CV_BGR2GRAY);    // to preserve the original image    cv::Mat flow = gray.clone();    int width = img.cols;    int height = img.rows;    int graySize = width * height;    // "brighten" the flow image     // C++ version of:    // vis[:] = (192 + np.uint32(vis)) / 2    for (unsigned int i=0; i<graySize; ++i)    {         flow.data[i] = (uchar)((192 + (int)flow.data[i]) / 2);    }    cv::Mat eigen = cv::Mat(height, width, CV_32FC(6));    cv::cornerEigenValsAndVecs(gray, eigen, 15, 3);    // this is the equivalent to all the numpy's reshaping etc. to     // generate the flow arrays    // simply use channel 4 and 5 as the actual flow array in C++    std::vector<cv::Mat> channels;    cv::split(eigen, channels);    int d = 12;    cv::Scalar col(0, 0, 0);    // C++ version of:    // points =  np.dstack( np.mgrid[d/2:w:d, d/2:h:d] ).reshape(-1, 2)    // including the actual line drawing part    for (unsigned int y=(d/2); y<flow.rows; y+=d)    {         for (unsigned int x=(d/2); x<flow.cols; x+=d)         {             if (x < flow.cols && y < flow.rows)             {                 cv::Point p(x, y);                 float dx = channels[4].at<float>(p) * (d/2);                 float dy = channels[5].at<float>(p) * (d/2);                 cv::Point p0(p.x - dx, p.y - dy);                 cv::Point p1(p.x + dx, p.y + dy);                 cv::line(flow, p0, p1, col, 1);              }         }    }    tm.stop();    std::cout<<"Flow image generated in "<<tm.getTimeMilli()<<" ms."<<std::endl;    cv::imshow("FLOW", flow);    cv::waitKey();    return 0;}