4-python图像处理之Scipy

SciPy（http://scipy.org/） 是建立在NumPy 基础上， 用于数值运算的开源工具包。
SciPy 提供很多高效的操作，可以实现数值积分、优化、统计、信号处理，以及对
我们来说最重要的图像处理功能。

---

一. 图像模糊

1.图像的高斯模糊是非常经典的图像卷积例子本质上，图像模糊就是
将（灰度）图像I 和一个高斯核进行卷积操作：
I = I*G
其中* 表示卷积操作。
高斯模糊通常是其他图像处理操作的一部分，比如图像插值操作、兴趣点计算以及
很多其他应用。

SciPy 有用来做滤波操作的scipy.ndimage.filters 模块。该模块使用快速一维分离
的方式来计算卷积。

测试案例：

#coding=utf-8
#
# 图像模糊

#高斯模糊
'''
高斯模糊通常是其他图像处理操作的一部分，比如图像插值操作、兴趣点计算以及
很多其他应用。
'''
from PIL import Image
from numpy import *
import pylab as plb
from scipy.ndimage import filters
import matplotlib.pyplot as plt
im = array(Image.open('test1.jpg').convert('L'))
#随着第二个参数的增大模糊程度也逐渐增大
plb.gray()
for i in range(4):
    im2 = filters.gaussian_filter(im,i*2)
    plt.subplot(2,2,i)
    plt.title(i)
    plb.imshow(im2)
plb.show()

二.图像导数

在很多应用中图像强度的变化情况是非常重要的信息。强度的
变化可以用灰度图像I（对于彩色图像，通常对每个颜色通道分别计算导数）的x
和y 方向导数Ix 和Iy 进行描述。

测试案例：

#coding=utf-8
#
# 图像倒数

from PIL import Image
from numpy import *
from scipy.ndimage import filters
import pylab as plb
im = array(Image.open('test1.jpg').convert('L'))
# Sobel 导数滤波器
imx = zeros(im.shape)   #对x求导
img2 = Image.fromarray(imx)

#img2.show()
filters.sobel(im,1,imx)
imy = zeros(im.shape)   #对y求导
#img2 = Image.fromarray(imy)

#img2.show()
filters.sobel(im,0,imy)
#和，，将x y求导结果合一
magnitude = 255-sqrt(imx**2+imy**2)
#img2 = Image.fromarray(magnitude)
plb.gray()
plb.imshow(magnitude)
plb.show()

---

三.形态学处理计数

形态学（或数学形态学）是度量和分析基本形状的图像处理方法的基本框架与集合。
形态学通常用于处理二值图像，但是也能够用于灰度图像。二值图像是指图像的每
个像素只能取两个值，通常是0 和1。二值图像通常是，在计算物体的数目，或者
度量其大小时，对一幅图像进行阈值化后的结果。你可以从http://en.wikipedia.org/
wiki/Mathematical_morphology 大体了解形态学及其处理图像的方式。

测试案例：

#coding=utf-8
#
# 形态学：对象计数

from PIL import Image
from numpy import *
from scipy.ndimage import measurements,morphology
import pylab as plb
# 载入图像，然后使用阈值化操作，以保证处理的图像为二值图像
im = array(Image.open('test3.png').convert('L'))
im = 1*(im<225)#暂时还不太清楚其所以然，0.0
labels, nbr_objects = measurements.label(im)
print "Number of objects:", nbr_objects


'''
imgLabels = Image.fromarray(labels)
imgLabels.show()
'''
#形态学的开操作，iterations=5进行五次
im_open = morphology.binary_opening(im,ones((2,2)),iterations=5)
labels_open, nbr_objects_open = measurements.label(im_open)
print "Number of objects:", nbr_objects_open
plb.gray()
plb.imshow(labels_open)
plb.show()

四. 图像去噪简单实例：

#coding=utf-8
#
# 图像去噪小案例


'''
使用A. Chambolle（2005）在公式（11）中的计算步骤实现Rudin-Osher-Fatemi（ROF）去噪模型
输入：含有噪声的输入图像（灰度图像）、U 的初始值、TV 正则项权值、步长、停业条件
输出：去噪和去除纹理后的图像、纹理残留
'''
from numpy import *
def denoise(im,U_init,tolerance=0.1,tau=0.125,tv_weight=100):
    m,n = im.shape # 噪声图像的大小
    # 初始化
    U = U_init
    Px = im # 对偶域的x 分量
    Py = im # 对偶域的y 分量
    error = 1
    while (error > tolerance):
        Uold = U
        # 原始变量的梯度
        GradUx = roll(U,-1,axis=1)-U # 变量U 梯度的x 分量
        GradUy = roll(U,-1,axis=0)-U # 变量U 梯度的y 分量
        # 更新对偶变量
        PxNew = Px + (tau/tv_weight)*GradUx
        PyNew = Py + (tau/tv_weight)*GradUy
        NormNew = maximum(1,sqrt(PxNew**2+PyNew**2))
        Px = PxNew/NormNew # 更新x 分量（对偶）
        Py = PyNew/NormNew # 更新y 分量（对偶）
        # 更新原始变量
        RxPx = roll(Px,1,axis=1) # 对x 分量进行向右x 轴平移
        RyPy = roll(Py,1,axis=0) # 对y 分量进行向右y 轴平移
        DivP = (Px-RxPx)+(Py-RyPy) # 对偶域的散度
        U = im + tv_weight*DivP # 更新原始变量
        # 更新误差
        error = linalg.norm(U-Uold)/sqrt(n*m);
    return U,im-U # 去噪后的图像和纹理残余


from PIL import Image
from pylab import *
import matplotlib.pyplot as plt
im = array(Image.open('test3.jpg').convert('L'))
U,T = denoise(im,im)
gray()
plt.subplot(121)
imshow(im)
plt.subplot(122)
imshow(U)
axis('equal')
axis('off')
show()