A 依据颜色的k均值聚类
def _1rgb_kmeans(picname='', cutnum=50, clus=5):
'''
彩色图像按颜色k-means聚类.
依赖:from scipy.cluster import vq
from scipy.misc import imresize #This function is only available if Python Imaging Library (PIL) is installed.
from PIL import Image
:param picname: 图像名,字符串
:param cutnum: 切成cut*cut个块
:param clus: 分几类
:return: 无,显示图像
'''
from scipy.misc import imresize
import scipy.cluster.vq # 包导入只会在第一次导入时执行
gridsize = cutnum # 多少格
im_a = np.array(Image.open(picname))
x_num = int(im_a.shape[0] / gridsize) # 行方向上的像素数
y_num = int(im_a.shape[1] / gridsize)
grids = []
for x in range(gridsize):
for y in range(gridsize):
R = np.mean(im_a[x*x_num: (x+1)*x_num, y*y_num: (y+1)*y_num, 0]) # 区域内的像素均值,红色
G = np.mean(im_a[x * x_num:(x + 1) * x_num, y * y_num:(y + 1) * y_num, 1])
B = np.mean(im_a[x * x_num:(x + 1) * x_num, y * y_num:(y + 1) * y_num, 2])
grids.append([R,G,B])
grids = np.array(grids, 'f') # 变为一溜数组
centroids, _ = vq.kmeans(grids, clus) # (数据,分几类)->中心点
lable, dist = vq.vq(grids, centroids) # 为原数据根据中心点标签
lable_im_a = lable.reshape(gridsize, gridsize)
lable_im_a = imresize(lable_im_a, im_a.shape[0:2], interp='nearest')#(矩阵,大小)
plt.close('all') #一个习惯
plt.figure(1)
plt.imshow(lable_im_a) # 从图片矩阵中显示
plt.show()
B 图片求导
def _2deriv_im():
sigma = 2
im_a = np.array(Image.open('big.jpg').convert('L'))
im_x = np.zeros(im_a.shape)
filters.gaussian_filter(im_a, (sigma,sigma), (1,1), im_x)
# (1,1): order : int or sequence of ints, optional
# The order of the filter along each axis is given as a sequence
# of integers, or as a single number. An order of 0 corresponds
# to convolution with a Gaussian kernel. A positive order
# corresponds to convolution with that derivative of a Gaussian.
# 可以试试(0,1)或(1,0),效果会不同
# filters.gaussian_filter有两种等价写法:output可以放在参数里也可放在
# 等号左边函数的返回值。注意要求input与output不重名
plt.figure()
plt.imshow(im_x)
plt.show()
C ROF去噪算法
def _3rof_denoise(image, tolerance=0.1, step=0.125, tv_weight=100, U_init=None):
'''
rof去噪算法,目标为min{ ||I-U||**2 + 2*weight* J(U)}, J(U)为全变差
:param image: 灰度图像矩阵
:param tolerance: 停止循环阈值
:param step:
:param tv_weight:正则项J的权值
:param U_init:如果提供初始值。否则初始值就是原矩阵
:return:去噪后新图像和残余
'''
m,n = image.shape
U = U_init if U_init is not None else image
Px = image #x分量
Py = image
error = 1
while(error > tolerance): # 梯度下降,直到小于阈值,实现min{}
Uold = U
gradUx = np.roll(U, -1, axis=1) - U # 原始变量的x梯度
# ~这个式子计算了U的宽度向差分差分
gradUy = np.roll(U, -1, axis=0) - U # roll:在某一或所有坐标轴上“滚动”,使得整体移位。
# >>> x2
# array([[0, 1, 2, 3, 4],
# [5, 6, 7, 8, 9]])
# >>> np.roll(x2, 1)
# array([[9, 0, 1, 2, 3],
# [4, 5, 6, 7, 8]])
# >>> np.roll(x2, 1, axis=0)
# array([[5, 6, 7, 8, 9],
# [0, 1, 2, 3, 4]])
# >>> np.roll(x2, 1, axis=1)
# array([[4, 0, 1, 2, 3],
# [9, 5, 6, 7, 8]])
PxNew = Px + (step / tv_weight) * gradUx
PyNew = Py + (step / tv_weight) * gradUy
normNew = np.maximum(1,np.sqrt(PxNew**2 + PyNew**2))
Px = PxNew / normNew # 更新x分量
Py = PyNew / normNew
rPx = np.roll(Px, 1, axis=1) # 更新原始x分量,坐标轴平移
rPy = np.roll(Py, 1, axis=0)
divP = (Px - rPx) + (Py - rPy) #对偶域散度
U = image + tv_weight * divP # 更新输出变量
error = np.linalg.norm(U - Uold) / np.sqrt(m*n) # norm()衡量矩阵间差异
return U, image-U
D 图像角点值计算与角点数统计
def _4compute_harris_response(ima, sig=3):
'''
计算一张图像(矩阵形式)的harris响应值
:param ima:
:param sig: 越大最后得到值越小
:return:
'''
imx = np.zeros(im_a.shape)
imy = np.zeros(im_a.shape)
filters.gaussian_filter(im_a, (sig,sig), (0,1), imx)
filters.gaussian_filter(im_a, (sig,sig), (1,0), imy) # 计算y方向高斯导数
Wxx = filters.gaussian_filter(imx*imx, sig) # harris矩阵的各分量
Wyy = filters.gaussian_filter(imy*imy, sig)
Wxy = filters.gaussian_filter(imx*imy, sig)
Wdet = Wxx * Wyy - Wxy**2 # 矩阵乘法
Wtr = Wxx + Wyy # 特征值和迹
Wtr = Wtr+ 0.01 # 发现对大块相同的区域,Wxx Wyy会出现导数为0,Wtr中也有0,后面除法做不下去了。
# 万般无奈,试试+0.1吧……验证了一下,影响不大
ans = Wdet / Wtr # 矩阵除法
return ans
def _5get_harris_points(harris_ima, min_dist=10, threshold=0.1):
'''
输入已经计算过harris响应值的图像矩阵,输出相应的角点
:param harris_ima: 以harris响应值为值的图像矩阵
:param min_dist: 分割角点和图像边界的最小像素数,角点的间隔范围的控制
:param threshold: 候选点阈值,越小生成的角点越多
:return: 角点(x,y)坐标列表
'''
corner_thre = harris_ima.max() * threshold
harris_ima_ = (harris_ima > corner_thre) * 1 # 寻找高于阈值的候选点并0-1化
coords = np.array(harris_ima_.nonzero()).T # 得到这些点的坐标,coords第一维为点集序号,每个点的第0维存x、1为y。 @nonzero()
candidate_val = [harris_ima[c[0], c[1]] for c in coords] # 得到这些点的角点值和坐标
indexs = np.argsort(candidate_val) # 按harris响应值排序,得到序列 @argsort
allowed_loca = np.zeros(harris_ima.shape)
allowed_loca[ min_dist:-min_dist, min_dist:-min_dist] = 1 # 标记可行点
selected = []
for i in indexs:
if allowed_loca[ coords[i,0], coords[i,1] ] == 1:
selected.append( (coords[i,1], coords[i,0]) )
allowed_loca[ (coords[i,0] - min_dist):(coords[i,0] + min_dist),
(coords[i, 1] - min_dist):(coords[i, 1] + min_dist) ] = 0 #已定点周围标记为0
return selected