基于《数字图像处理 第四版》冈萨雷斯 这本书进行算法操作,有些实现的可能不是很完美,请辩证看待。
图像增强:
一、附件中的20幅噪声图像中均加入了均值为0,方差为0.03的高斯噪声,请对图像进行去噪增强,使增强效果尽可能接近图1gray.jpg的效果,请上传增强后的效果图,并写出原理。
由于20幅噪声图像中均加入了均值为0,方差为0.03的高斯噪声。可以认为当这些噪声像素全部进行线性叠加求平均后,得到的图像即为去噪图像。
算法原理:1、像素的线性叠加。测试图片读入是三通道的ndarray格式数据,输入图片后,由于数据格式是unit8,范围为0-255,如果进行简单的叠加会使数据溢出,造成误差,首先对像素点进行归一化,归一化后再进行相加求和求平均,再将求得的均值乘回原来的0-255范围,输出图像即为去噪后的图像。
2、图像的平均:由于老师提供了原图,所以选择采用处理后的图像与原图像进行像素对比。选择均方根误差 RMSE和归一化均方根误差 NRMSE进行判断,指标越小,说明图像去噪的效果越好。
关键代码:
def cv_imread(file_path):
"""
批量读取图片,由于路径中有中文名称于是将其转化后读取
:param file_path:图片路径
:return:
"""
cv_img = cv2.imdecode(np.fromfile(file_path, dtype=np.uint8), cv2.IMREAD_COLOR)
return cv_img
def read_path(file_pathname):
"""
遍历目标目录下的所有图片文件
:param file_pathname: 目录的地址
:return: 返回一个图像列表
"""
images = []
for filename in os.listdir(file_pathname):
img = cv_imread(r"D:\desk\图像增强1\noise" + "\\" + filename)
images.append(img)
return images
def cal_mean(img_list):
"""
线性图像相加去噪
:param img_list: 含有多个图像
:return:
"""
sum = img_list[0] - img_list[0]
for i in range(len(img_list)):
sum = sum + img_list[i] / 255.0
sum_final = sum / 20 * 255.0
aver = np.around(sum_final)
pic_out = aver.astype('uint8')
return pic_out
def multi_pic():
path = 'D:\desk\图像增强1\\noise'
image_pp = read_path(path)
image_after = cal_mean(image_pp)
img_raw = cv2.imread('pic_gray_1.jpg')
# print(image_after)
# cv2.imshow('raw', image_pp[0])
# cv2.imshow('after', image_after)
# cv2.imshow('raw_orign',img_raw)
cv2.waitKey(0)
cv2.imwrite('denoise_after2.jpg', image_after)
print(np.sqrt(skimage.metrics.mean_squared_error(img_raw,image_after)))
print(np.sqrt(skimage.metrics.mean_squared_error(img_raw,image_pp[0])))
print(skimage.metrics.normalized_root_mse(img_raw,image_after))
print(skimage.metrics.normalized_root_mse(img_raw,image_pp[0]))
结果图:
4、结果分析:
将结果图和数据集中加噪声的图片与原始图片进行对比,得到的结果如下:
二、请对下图做简单的空域均值滤波,模板为15*15的矩阵,用三种不同的方法考虑边界滤波后的效果,并提交所采用的三种方法和对应的效果。(内容请包含 1. 算法原理; 2.关键代码;3.结果图;4.结果分析)
def MeanFilter(img, K_size):#扩展边缘版
# 均值滤波
h, w, c = img.shape
# 计算边缘点的零填充
pad = K_size // 2
out1 = np.zeros((h + 2 * pad, w + 2 * pad, c), dtype=np.float)
out1[pad:pad + h, pad:pad + w] = img.copy().astype(np.float)
# 卷积的过程
tmp = out1.copy()
for y in range(h):
for x in range(w):
for ci in range(c):
out1[pad + y, pad + x, ci] = np.mean(tmp[y:y + K_size, x:x + K_size, ci])
out_1 = out1[pad:pad + h, pad:pad + w].astype(np.uint8)
return out_1
第二种方法(不扩充,保留边界)
def MeanFilter_2(img, K_size):#不对边缘点进行0填充
# 均值滤波
h, w, c = img.shape
# 计算边缘点的零填充
pad = K_size // 2
out1 = np.zeros((h + 2 * pad, w + 2 * pad, c), dtype=np.float)
out1[pad:pad + h, pad:pad + w] = img.copy().astype(np.float)
# 卷积的过程
tmp = out1.copy()
for y in range(h-2*pad):
for x in range(w-2*pad):
for ci in range(c):
out1[pad + y+pad, pad + x+pad, ci] = np.mean(tmp[y:y + K_size, x:x + K_size, ci])
out_1 = out1[pad:pad + h, pad:pad + w].astype(np.uint8)
return out_1
第三种:扩充,但填的是边界像素版
Python
def MeanFilter_3(img, K_size):#对边缘点外的边界点也进行0填充 | |
# 均值滤波 | |
h, w, c = img.shape | |
# 计算边缘点的零填充 | |
pad = K_size // 2 | |
out1 = np.zeros((h + 4 * pad, w + 4 * pad, c), dtype=np.float) | |
out1[pad:pad + h, pad:pad + w] = img.copy().astype(np.float) | |
out1[pad-1,0:w-1]=img[0,0:w-1].copy().astype(np.float) | |
out1[h+pad +1,0:w-1] = img[h-1,0:w-1].copy().astype(np.float) | |
out1[0:h-1,pad-1] = img[0:h-1,0].copy().astype(np.float) | |
out1[0:h-1,w + pad + 1] = img[0:h-1,w-1].copy().astype(np.float) | |
# 卷积的过程 | |
tmp = out1.copy() | |
for y in range(h+2*pad): | |
for x in range(w+2*pad): | |
for ci in range(c): | |
out1[y+pad, x+pad, ci] = np.mean(tmp[y:y + K_size, x:x + K_size, ci]) | |
out_1 = out1[pad:pad + h, pad:pad + w].astype(np.uint8) | |
return out_1 |
三、请用直方图均衡对下图进行增强
算法原理:1、利用图像的直方图对图像进行对比度修正,由于给的图像是三通道的,首先要将图像通道拆开来,变成三个单通道,然后计算每个通道的直方图。直方图的计算过程:首先统计不同像素值,再统计不同像素值在图像中对应的像素点个数,在计算不同像素值在所有像素中的比例,最终画出直方图。
2、得到直方图后,进行直方图均衡化,首先建立映射表,计算直方图累计概率,将映射表中不同像素值对应的累计概率存入映射表中,计算完后将结果存入映射表中,再将映射表*255,返回真实像素值,最终得到直方图均衡后的图像。
代码:
1、画直方图与计算直方图的代码:
Python
img = cv2.imread('2.jpg') | |
# R-G-B三种通道直方图 | |
countb = np.zeros(256, np.float32) | |
countg = np.zeros(256, np.float32) | |
countr = np.zeros(256, np.float32) | |
for i in range(img.shape[0]): | |
for j in range(img.shape[1]): | |
(b, g, r) = img[i, j] | |
b = int(b) | |
g = int(g) | |
r = int(r) | |
countb[b] += 1 # 统计该像素出现的次数 | |
countg[g] += 1 | |
countr[r] += 1 | |
countb = countb / (img.shape[0] * img.shape[1]) # 得到比例 | |
countg = countg / (img.shape[0] * img.shape[1]) | |
countr = countr / (img.shape[0] * img.shape[1]) | |
x = np.linspace(0, 255, 256) | |
plt.figure() | |
plt.bar(x, countb, color='b') | |
plt.figure() | |
plt.bar(x, countg, color='g') | |
plt.figure() | |
plt.bar(x, countr, color='r') | |
plt.show() |
结果
请对下图(二选一)做伽马变换,以增强其对比度。
算法原理:做伽马变换,增加对比度。由于所给的照片整体的灰度值偏暗,需要扩展像素的灰度大小,使其变亮,于是采用幂次变换进行。设计的算法流程是将图像像素点首先归一化后使输入值在[0,1]间,通过求函数的幂
,其中s为变化后输出,c设置为1,r为输出,y为伽马值,由于图像整体值偏暗,要让其变亮,选择y值小于1进行变换。
关键代码:
def gamma_bianhuan(image, gamma):
image = image / 255.0
New = np.power(image, gamma)
return New
def question_4():
a = cv2.imread('1.jpg')
image1 = cv2.split(a)[0] # 蓝
image2 = cv2.split(a)[1] # 绿
image3 = cv2.split(a)[2] # 红
cv2.imshow("raw_3", image1)
# image_1 = gamma_bianhuan(image1, 1.5)
image_1_2 = gamma_bianhuan(image1, 0.5)
image_2 = gamma_bianhuan(image2, 0.6)
image_3 = gamma_bianhuan(image3, 0.4)
# merged = cv2.merge([image1, image2, image3])
cv2.imshow("chuli_1_2", image_1_2)
# cv2.imshow("chuli_1_1", image_1)
cv2.imshow("chuli2",image_2)
cv2.imshow("chuli3", image_3)
# cv2.imshow("chuli_pp",merged)
cv2.waitKey(0)
五、请对图像分别采用均值滤波器和中值滤波器进行滤波处理,并对处理后的结果进行分析。
算法原理:
均值滤波:对目标像素及周边像素取平均值后再填目标像素来实现滤波目的。
中值滤波:非线性的处理方法,在去噪的同时可以兼顾到边界信息的保留。选一个奇数*奇数的模板,将这个模板与像素点比对,把模板中所含的像素点按灰度级的升或降序排列,取灰度值的中位数来代替该点的灰度值。
设置模板为3*3的均值滤波器和中值滤波器,使各个权重和为1,对像素进行填充卷积,最终得到处理后的图像。
def MeanFilter(img, K_size):
"""
空域均值滤波
:param img:待处理图像
:param K_size: 卷积核大小
:return:
"""
# 均值滤波
h, w, c = img.shape
# 计算边缘点的零填充
pad = K_size // 2
out1 = np.zeros((h + 2 * pad, w + 2 * pad, c), dtype=np.float)
out1[pad:pad + h, pad:pad + w] = img.copy().astype(np.float)
# 卷积的过程
tmp = out1.copy()
for y in range(h):
for x in range(w):
for ci in range(c):
out1[pad + y, pad + x, ci] = np.mean(tmp[y:y + K_size, x:x + K_size, ci])
out_1 = out1[pad:pad + h, pad:pad + w].astype(np.uint8)
return out_1
def MedianFilter(img, K_size):
"""
中值滤波器
:param img: 目标图像
:param K_size: 卷积核大小
:return:
"""
h, w, c = img.shape
# 零填充
pad = K_size // 2
out = np.zeros((h + 2 * pad, w + 2 * pad, c), dtype=np.float)
out[pad:pad + h, pad:pad + w] = img.copy().astype(np.float)
# 卷积的过程
tmp = out.copy()
for y in range(h):
for x in range(w):
for ci in range(c):
out[pad + y, pad + x, ci] = np.median(tmp[y:y + K_size, x:x + K_size, ci])
out = out[pad:pad + h, pad:pad + w].astype(np.uint8)
return out
def question_5():
img_5=cv2.imread('image_salt_pepper2.jpg')
img_mean=MeanFilter(img_5,3)
img_median=MedianFilter(img_5,3)
cv2.imshow('raw_5',img_5)
cv2.imshow('mean',img_mean)
cv2.imshow('median',img_median)
cv2.waitKey(0)
六、请对下图进行锐化处理,并写出所用算子和分析结果
算法原理:锐化空间滤波器主要有基于一阶微分和基于二阶微分的锐化滤波器,而将一阶微分定义为
,二阶微分定义为
在对于图像进行锐化时,我们使用sobel算子和拉普拉斯算子作为滤波模板,根据公式:
得到拉普拉斯矩阵的模板为,
最终的计算公式为
则根据模板对图像进行卷积操作,由于中间系数为负数,需要原图像减去拉普拉斯变换后的图像,即c=-1,若中间系数为正数,则为加上,即c=1。
关键代码:
拉普拉斯锐化过程函数:
def laplace(img,kernel,c):
des_8U = cv2.filter2D(img, -1, kernel=kernel, borderType=cv2.BORDER_DEFAULT)
des_16S = cv2.filter2D(img, ddepth=cv2.CV_16SC1, kernel=kernel, borderType=cv2.BORDER_DEFAULT)
if c==-1:
g = img - des_16S
if c==1:
g = img + des_16S
g[g < 0] = 0
g[g > 255] = 255
# g=cv2.merge([g,g,g])
return g
定义四种拉普拉斯算子模板:
kernel1 = np.asarray([[0, 1, 0],
[1, -4, 1],
[0, 1, 0]])
im1_laplace1=laplace(f, kernel1,-1)
cv2.imwrite("laplace1.jpg", im1_laplace1)
kernel2 = np.asarray([[1, 1, 1],
[1, -8, 1],
[1, 1, 1]])
im1_laplace2=laplace(f, kernel2,-1)
kernel3 = np.asarray([[0, -1, 0],
[-1, 4, -1],
[0, -1, 0]])
im1_laplace3=laplace(f, kernel3,1)
cv2.imwrite("laplace3.jpg", im1_laplace3)
kernel4 = np.asarray([[-1, -1, -1],
[-1, 8, -1],
[-1, -1, -1]])
im1_laplace4=laplace(f, kernel4,1)
cv2.imwrite("laplace4.jpg", im1_laplace4)
#sobel模板
x = cv2.Sobel(binary, cv2.CV_16S, 1, 0)
y = cv2.Sobel(binary, cv2.CV_16S, 0, 1)
absX = cv2.convertScaleAbs(x)
absY = cv2.convertScaleAbs(y)
Sobel = cv2.addWeighted(absX, 0.5, absY, 0.5, 0)
img_sobel=cv2.merge([Sobel,Sobel,Sobel])+img
cv2.imwrite("sobel.jpg",img_sobel)
七、请对下图进行增强,以使其达到较佳的图像显示效果
(内容请包含 1. 算法原理; 2.关键代码;3.结果图;4.结果分析)
算法原理:待处理图像整体偏亮,对于图像的细节突出的不是很好,为使其达到较佳的图像显示效果,采用拉普拉斯法突出图像细节,采用经过sobel梯度算子突出图像边缘,将sobel梯度算子产生的图像进行均值滤波后与拉普拉斯法产生的图像进行叠加,得到最终的图像。
关键代码:
均值滤波:
def laplace(img,kernel,c):
des_8U = cv2.filter2D(img, -1, kernel=kernel, borderType=cv2.BORDER_DEFAULT)
des_16S = cv2.filter2D(img, ddepth=cv2.CV_16SC1, kernel=kernel, borderType=cv2.BORDER_DEFAULT)
if c==-1:
g = img - des_16S
if c==1:
g = img + des_16S
g[g < 0] = 0
g[g > 255] = 255
# g=cv2.merge([g,g,g])
return g
def pic_enhence_1():
img_raw=cv2.imread("body_x_ray.jpg")
img1=cv2.imread("body_x_ray.jpg",cv2.IMREAD_GRAYSCALE)
kernel2 = np.asarray([[1, 1, 1],
[1, -8, 1],
[1, 1, 1]])
im1_laplace=laplace(img1, kernel2,-1)
cv2.imwrite("laplace_q2.jpg", cv2.merge([im1_laplace,im1_laplace,im1_laplace]))
sobel_1(img1,img_raw)
def pic_enhence_2():
img_laplace=cv2.imread("laplace2.jpg")
img_sobel=cv2.imread("sobel_q2.jpg")
img_sobel=img_sobel
img_laplace=img_laplace
g=img_sobel+img_laplace
# g=cv2.multiply(img_laplace,img_sobel)
cv2.imshow("g",g)
cv2.imwrite("q2.jpg",g)
cv2.waitKey(0)
八、请对下面图像进行祛斑或去皱处理,二选一即可
算法原理:对于图像进行去斑处理,首先采用sobel梯度算法提取出边缘像素点,然后通过连通域分析确定出斑点位置,由于五官也算是连通域,所以可以通过面积的阈值来区分五官与斑点位置,连通域面积小的是斑点,连通域面积较大的是五官区域。得到斑点的连通域后再通过图像修复函数inpaint函数对于斑点进行修复,最后在对图像进行幂律变换,使其总体灰度值变亮,最终得到祛斑处理后的图像。
关键代码:
幂律变换:
Python
def gamma_bianhuan(image, gamma): | |
image = image / 255.0 | |
New = np.power(image, gamma) | |
return New | |
a = cv2.imread('1.jpg') | |
image1 = cv2.split(a)[0] # 蓝 | |
image2 = cv2.split(a)[1] # 绿 | |
image3 = cv2.split(a)[2] # 红 | |
# cv2.imshow("raw_3", image1) | |
# image_1 = gamma_bianhuan(image1, 1.5) | |
image_1_2 = cv2.merge([gamma_bianhuan(image1, 0.5),gamma_bianhuan(image1, 0.5),gamma_bianhuan(image1, 0.5)]) |
sobel算子找到可能点:
Plain Text
src=cv2.imread("3 (2).jpg") | |
src1 = cv2.cvtColor(src,cv2.COLOR_BGR2GRAY) | |
x_gray = cv2.Sobel(src1,cv2.CV_32F,1,0) | |
y_gray = cv2.Sobel(src1,cv2.CV_32F,0,1) | |
x_gray = cv2.convertScaleAbs(x_gray) | |
y_gray = cv2.convertScaleAbs(y_gray) | |
dst = cv2.add(x_gray,y_gray,dtype=cv2.CV_16S) | |
dst = cv2.convertScaleAbs(dst) |
经过高斯滤波后,通过连通域分析找到斑点:
Plain Text
gaussianBlur = cv2.GaussianBlur(src1, (3, 3), 0) | |
ret, binary = cv2.threshold(gaussianBlur, 127, 255, cv2.THRESH_BINARY) | |
num_labels,labels,stats,centers = cv2.connectedComponentsWithStats(binary, connectivity=8,ltype=cv2.CV_32S) | |
for t in range(1, num_labels, 1): | |
x, y, w, h, area = stats[t] | |
if area>100: | |
index = np.where(labels==t) | |
labels[index[0], index[1]] = 0 |
去除斑点:
Plain Text
dst_TELEA = cv2.inpaint(src1,labels,3,cv2.INPAINT_TELEA) | |
cv2.imshow("final",dst_TELEA) |
九、请将下图进行增强,以能显示阴影区域内细节
算法原理:
由于需要显示原图像中阴影区域细节,在频率域中,图像边缘及尖锐细节往往是高频信号,于是可以采用高通滤波器对图像进行处理。高通滤波器有理想高通滤波器(IHPF),巴特沃斯高通滤波器(BHPF),高斯高通滤波器(GHPF)三种类型。但是理想高通滤波器有明显振铃现象,巴特沃斯高通滤波器较平滑,边缘失真小,二阶BHPF只有轻微振铃现象,高斯高通滤波器没有振铃现象,完全平滑。结合三种滤波器的特点,理论上应该选择使用高斯滤波器进行图像增强,这里将三种滤波器都实现了,以观察增强后的效果。
巴特沃斯高通滤波器的数学表达式
高斯滤波器的数学表达式:
关键代码:
1、离散傅里叶变换
dft = cv.dft(np.float32(img), flags=cv.DFT_COMPLEX_OUTPUT)
# 中心化
dtf_shift = np.fft.fftshift(dft)
rows, cols = img.shape
crow, ccol = int(rows-1 / 2), int(cols-1 / 2) # 计算频谱中心
mask = np.zeros((rows, cols, 2)) # 生成rows行cols列的二维矩阵
2、 三种滤波器的实现:
for i in range(rows):
for j in range(cols):
D = np.sqrt((i - crow) ** 2 + (j - ccol) ** 2) # 计算D(u,v)
if (filter.lower() == 'butterworth'): # 巴特沃斯滤波器
if (type == 'lp'):
mask[i, j] = 1 / (1 + (D / D0) ** (2 * N))
elif (type == 'hp'):
mask[i, j] = 1 / (1 + (D0 / D) ** (2 * N))
elif (filter.lower() == 'ideal'): # 理想滤波器
if (type == 'lp'):
if (D <= D0):
mask[i, j] = 1
elif (type == 'hp'):
if (D > D0):
mask[i, j] = 1
elif (filter.lower() == 'gaussian'): # 高斯滤波器
if (type == 'lp'):
mask[i, j] = np.exp(-(D * D) / (2 * D0 * D0))
elif (type == 'hp'):
mask[i, j] = (1 - np.exp(-(D * D) / (2 * D0 * D0)))
3、傅里叶反变换与反变换后参数的归一化(由于值很大,需要归一化后*255)
fshift = dtf_shift * mask
f_ishift = np.fft.ifftshift(fshift)
img_back = cv.idft(f_ishift)
img_back = cv.magnitude(img_back[:, :, 0], img_back[:, :, 1]) # 计算像素梯度的绝对值
img_back = np.abs(img_back)
img_back = 255*(img_back - np.amin(img_back)) / (np.amax(img_back) - np.amin(img_back))
4、主代码
img=cv2.imread("HF.jpg",0)
img_back1 = filter(img, 20, type='hp', filter='ideal')
cv2.imwrite("img_pp1.jpg",cv2.merge([img_back1,img_back1,img_back1]))
img_back2 = filter(img,20, type='hp', filter='butterworth')
cv2.imwrite("img_pp2.jpg",cv2.merge([img_back2,img_back2,img_back2]))
img_back3 = filter(img, 20, type='hp', filter='gaussian')
cv2.imwrite("img_pp3.jpg",cv2.merge([img_back3,img_back3,img_back3]))
十、请去除下面图中规律变化的噪声
算法原理:
由于图像中包含莫尔波纹,故应该采用陷波滤波减少莫尔波纹。陷波滤波器时一种选择滤波器,由于莫尔波纹是类似周期性的冲激信号,故采用陷波带阻滤波器对图像进行处理。
形式为:
距离计算公式为:
实现流程:首先对图像进行傅里叶变换,并且查看其傅里叶频谱图,效果如下:
可以看出,由于图像的莫尔波纹主要是竖直方向的,对应于频域里就是水平方向的条纹点,在频域图里面将对应位置置零,即可消除莫尔条纹。
关键代码:
def notch_filter(img, radius,uk,vk):
dft = cv.dft(np.float32(img), flags=cv.DFT_COMPLEX_OUTPUT)
# 中心化
dtf_shift = np.fft.fftshift(dft)
rows, cols = img.shape
crow, ccol = int(rows - 1 / 2), int(cols - 1 / 2) # 计算频谱中心
mask = np.zeros((rows, cols, 2), np.float32)
for i in range(rows):
for j in range(cols):
D = np.sqrt((i - crow) ** 2 + (j - ccol) ** 2) # 计算D(u,v)
D1=np.sqrt((i - rows//2 - vk)**2 + (j - cols//2 - uk)**2)
D2=np.sqrt((i - rows//2 + crow)**2 + (j - cols//2 + ccol)**2)
D0=radius
mask[i, j] = (1 / (1 + (D0 / (D1+1e-5))**2)) * (1 / (1 + (D0 / (D2+1e-5))**2))
for i in range(cols):
mask[331,i]=0
mask[331,317]= (1 / (1 + (D0 / (D1+1e-5))**2)) * (1 / (1 + (D0 / (D2+1e-5))**2))
fshift = dtf_shift * mask
f_ishift = np.fft.ifftshift(fshift)
img_back = cv.idft(f_ishift)
img_back = cv.magnitude(img_back[:, :, 0], img_back[:, :, 1]) # 计算像素梯度的绝对值
img_back = np.abs(img_back)
img_back = 255 * (img_back - np.amin(img_back)) / (np.amax(img_back) - np.amin(img_back))
return img_back
img=cv2.imread("1-2.png",0)
img1=notch_filter(img,10,50,60)
cv2.imwrite("pp2.png",cv2.merge([img1,img1,img1]))
结果图:
结果分析:
处理后的图像莫尔波纹明显减少,呈清晰图像显示,如果仔细调参数,可以把波纹完全去除,后续可以仔细寻找对应的波纹的频域坐标点,完全消除。
几何变换和彩色空间
一请分别用最近邻法和双线性插值法将附件中80*60的小图放大至800*600,并分别截取边缘区域与原图相比较。
(内容请包含 1. 算法原理; 2.关键代码;3.结果图;4.结果分析)
算法原理:
1、最近邻法:放大后图像的像素点值存在空缺的值,为了 对覆盖的每一个点赋灰度值,在原图像中寻找和放大后空缺点最接近的像素点所在位置(即最近邻点),将这个最近邻点的灰度值赋给空缺值,得到放大后的图像。
2、双线性插值:用四个最近邻去估计给定位置的灰度,公式为
,则(x,y)点的灰度值可以由公式计算得到,但是计算量很大,因为要求出abcd参数的值,列方程求解。
关键代码:
1、最近邻插值:输入参数,原图像,目标图像高度与宽度,输出放大后的图像,首先产生一个放大后的图像的空值模板,然后找到与待填充位置最接近的像素点,将其复制给待填充位置。
2、双线性插值:输入输入参数,原图像,目标图像高度与宽度,输出放大后的图像
def NN_interpolation(img, transformH, transformW): # 最近邻插值
scrH, scrW, _ = img.shape
retimg = np.zeros((transformH, transformW, 3), dtype=np.uint8)
for i in range(transformH):
for j in range(transformW):
scrx = round((i + 1) * (scrH / transformH))#寻找最接近像素点
scry = round((j + 1) * (scrW / transformW))
retimg[i, j] = img[scrx - 1, scry - 1]
return retimg
def BiLinear_interpolation(img, transformH, transformW): # 双线性插值法
scrH, scrW, _ = img.shape
img = np.pad(img, ((0, 1), (0, 1), (0, 0)), 'constant') #给边缘像素填充连续值
retimg = np.zeros((transformH, transformW, 3), dtype=np.uint8)
for i in range(transformH):
for j in range(transformW):
scrx = (i + 1) * (scrH / transformH) - 1
scry = (j + 1) * (scrW / transformW) - 1
x = math.floor(scrx)#返回一个向下舍入的数字
y = math.floor(scry)
u = scrx - x
v = scry - y
retimg[i, j] = (1 - u) * (1 - v) * img[x, y] + u * (1 - v) * img[x + 1, y] + (1 - u) * v * img[
x, y + 1] + u * v * img[x + 1, y + 1] # ax+by+cxy+d的反解方程代入
return retimg
二、
考察下面的一幅500x500的彩色图像,其中的方块部分分别为纯的红、绿和蓝色。
1. 如果将此图像转换到HSI空间,对H分量图像用一个25x25的算术平均掩模进行处理,再转换回到RGB空间,得到的结果将是怎么样的?
2.采取上面同样的步骤,只是这次处理的是S分量,结果又会怎样?
算法原理:图像转换到HSI空间,对H分量图像用一个25x25的算术平均掩模进行处理,对S分量转换进行处理同理。
首先是HSI模型的数学原理:H色度, S饱和度,I强度,三者的数学公式为:
其中为了防止出现除零错误,我们在分母里加一个极小的值
来避免这种情况。
实现思路:输入图像进行RGB的分离,分离出三个通道后通过RGB三通道的值根据上述公式计算出HSI模型的三个值,将HSI模型中的H分量与S分量进行25x25的算术平均掩模进行处理后,将HSI模型转化为RGB模型,得到结果。
其中HSI转RGB的数学公式如下所示:
def RGB2HSI(rgb_img):
#保存原始图像的行列数
row = np.shape(rgb_img)[0]
col = np.shape(rgb_img)[1]
#对原始图像进行复制
hsi_img = rgb_img.copy()
#对图像进行通道拆分
B,G,R = cv2.split(rgb_img)
#把通道归一化到[0,1]
[B,G,R] = [ i/ 255.0 for i in ([B,G,R])]
H = np.zeros((row, col)) #定义H通道
I = (R + G + B) / 3.0 #计算I通道
S = np.zeros((row,col)) #定义S通道
for i in range(row):
den = np.sqrt((R[i]-G[i])**2+(R[i]-B[i])*(G[i]-B[i]))
theta = np.arccos(0.5*(R[i]-B[i]+R[i]-G[i])/(den+np.finfo(np.float).eps)) #计算夹角,加个小角度,防止出现除0错误
h = np.zeros(col) #定义临时数组
#den>0且G>=B的元素h赋值为theta
h[B[i]<=G[i]] = theta[B[i]<=G[i]]
#den>0且G<=B的元素h赋值为theta
h[G[i]<B[i]] = 2*np.pi-theta[G[i]<B[i]]
#den<0的元素h赋值为0
h[den == 0] = 0
H[i] = h/(2*np.pi) #弧度化后赋值给H通道
#计算S通道
for i in range(row):
min = []
#找出每组RGB值的最小值
for j in range(col):
arr = [B[i][j],G[i][j],R[i][j]]
min.append(np.min(arr))
min = np.array(min)
#计算S通道
S[i] = 1 - min*3/(R[i]+B[i]+G[i])
#I为0的值直接赋值0
S[i][R[i]+B[i]+G[i] == 0] = 0
#扩充到255以方便显示,一般H分量在[0,2pi]之间,S和I在[0,1]之间
# hsi_img[:,:,0] = H/np.max(H)*255
hsi_img[:, :, 0] = H * 255
hsi_img[:,:,1] = S*255
hsi_img[:,:,2] = I*255
return hsi_img
HSI转RGB代码:
def HSI2RGB(hsi_img):
# 保存原始图像的行列数
row = np.shape(hsi_img)[0]
col = np.shape(hsi_img)[1]
#对原始图像进行复制
rgb_img = hsi_img.copy()
#对图像进行通道拆分
H,S,I = cv2.split(hsi_img)
#把通道归一化到[0,1]
[H,S,I] = [ i/ 255.0 for i in ([H,S,I])]
R,G,B = H,S,I
for i in range(row):
h = H[i]*2*np.pi
#H大于等于0小于120度时
a1 = h >=0
a2 = h < 2*np.pi/3
a = a1 & a2 #采用numpy的花式索引根据索引整型数组的值作为目标数组的某个轴的下标来取值
tmp = np.cos(np.pi / 3 - h)
b = I[i] * (1 - S[i])
r = I[i]*(1+S[i]*np.cos(h)/tmp)
g = 3*I[i]-r-b
B[i][a] = b[a]
R[i][a] = r[a]
G[i][a] = g[a]
#H大于等于120度小于240度
a1 = h >= 2*np.pi/3
a2 = h < 4*np.pi/3
a = a1 & a2
tmp = np.cos(np.pi - h)
r = I[i] * (1 - S[i])
g = I[i]*(1+S[i]*np.cos(h-2*np.pi/3)/tmp)
b = 3 * I[i] - r - g
R[i][a] = r[a]
G[i][a] = g[a]
B[i][a] = b[a]
#H大于等于240度小于360度
a1 = h >= 4 * np.pi / 3
a2 = h < 2 * np.pi
a = a1 & a2
tmp = np.cos(5 * np.pi / 3 - h)
g = I[i] * (1-S[i])
b = I[i]*(1+S[i]*np.cos(h-4*np.pi/3)/tmp)
r = 3 * I[i] - g - b
B[i][a] = b[a]
G[i][a] = g[a]
R[i][a] = r[a]
rgb_img[:,:,0] = B*255
rgb_img[:,:,1] = G*255
rgb_img[:,:,2] = R*255
return rgb_img
#利用opencv读入图片
rgb_img = cv2.imread('image1.bmp',cv2.IMREAD_COLOR)
#进行颜色空间转换
hsi_img = RGB2HSI(rgb_img)
H,S,I=cv2.split(hsi_img)
H_only=cv2.blur(H,(25,25))
S_only=cv2.blur(S,(25,25))
hsi_img_S_only=cv2.merge([H,S_only,I])
hsi_img_H_only=cv2.merge([H_only,S,I])
#进行滤波
hsi_img_mean=cv2.blur(hsi_img,(25,25))
rgb_mean=HSI2RGB(hsi_img_mean)
rgb_S_only=HSI2RGB((hsi_img_S_only))
rgb_H_only=HSI2RGB(hsi_img_H_only)
#opencv库的颜色空间转换结果
hsi_img_cv = cv2.cvtColor(rgb_img,cv2.COLOR_BGR2HSV)
hsi_img_cv=cv2.blur(hsi_img_cv,(25,25))
rgb_img_cv = cv2.cvtColor(hsi_img_cv,cv2.COLOR_HSV2BGR)
cv2.imshow("HSI", hsi_img_mean)
cv2.imshow("RGB",rgb_mean)
cv2.imshow('hsi_2',rgb_S_only)
cv2.imshow('hsi_3',rgb_H_only)
cv2.imwrite("HSI.jpeg",hsi_img)
cv2.imwrite("RGB.jpeg", rgb_mean)
cv2.imwrite("rgb_H_only.bmp",rgb_H_only)
cv2.imwrite("rgb_S_only.bmp",rgb_S_only)
cv2.imwrite("OpenCV_HSI.jpeg", hsi_img_cv)
cv2.imwrite("OpenCV_RGB.jpeg", rgb_img_cv)
cv2.waitKey()
cv2.destroyAllWindows()
结果分析:
opencv自带的库的结果和自编算法的HSI转换做出来的结果不同,经过调试与资料查找,发现自带的库函数是转为HSV模型,与HSI模型的I不同,计算公式为V = Max(R, G, B),与i取rgb的平均值不同。所以最终的结果不一样。
对于改变h分量进行均值掩模处理,可以看到,主要效果在边界处呈现。使用自编函数进行滤波后,中间线条是颜色叠加,绿色相连。
对改变I分量的操作,由于改变的是图片的强度分量,不会改变图片的色度,所以图片颜色保持不变。
彩色分割
1. (简答题)
请完整地将图中的小丑鱼分割出来(包括橙色,白色和黑色区域),分割后请用其它颜色做背景色。
算法原理:
题目要求将图中的小丑鱼分割出来,并且分割后用其他颜色做背景色,思路是首先将小丑鱼的轮廓给分割出来,轮廓作为一个mask,然后再生成一个其他颜色的背景图,将其他颜色背景图中小丑鱼轮廓对应的像素点扣出来,将两幅图像进行融合,最后得到结果图。
图像分割部分:彩色图像分割主要有基于RGB模型和HSV模型的图像分割,在RGB空间中,主要给定一个感兴趣的有代表性色彩的彩色样点集,可得到我们希望分割的颜色的“平均”估计。平均彩色用RGB向量α 表示,分割的目的是将图像中每个RGB像素分类,即在指定的区域内是否有一种颜色。主要通过欧式距离或者曼哈顿距离对图像进行分割,公式如下:
欧式距离:
曼哈顿距离:
由于题目给的这幅图小丑鱼的橙色与红色部分非常容易受光照影响(人眼也比较难分辨橙色与红色黄色),如果采用RGB空间进行图像分割,很暗分割出理想的小丑鱼,于是我决定采用HSI空间的分割(在实际应用中,由于转换成HSI模型需要手写算子,运算速率较低,于是采用opencv内置库的HSV模型,与HSI模型相近,除了强度I变成了亮度V,但是运算速率得到了大幅的提升)。
在进行HSV图像分割时,主要使用要分割区域的颜色阈值进行分割,比如橙色的HSV区域为(1,190,200)-(35,255,255),像素点的取值在这块区域的可以认为是橙色区域,通过对HSV颜色进行区域设置达到分割的效果。
图片叠加模块:生成一个纯色背景图,将纯色背景图中小丑鱼轮廓对应的像素点扣出来,将两幅图像进行融合,最后得到结果图。
关键代码:
图像分割模块:
#读取图片,并对图片进行处理,转成HSV格式
img = cv2.imread('fish.jpg')
img = cv2.cvtColor(img, cv2.COLOR_BGR2RGB)
hsv_img = cv2.cvtColor(img, cv2.COLOR_RGB2HSV)
#hsv颜色的阈值区间
light_red_1=(0, 43, 46)
dark_red_1=(14, 255, 255)
light_red_2=(156,43,46)
dark_red_2=(180,255,255)
light_orange = (1,190,200)
dark_orange = (35,255,255)
light_white = (0,0,200)
dark_white = (180,60,255)
light_black_1 = (0,0,46)
dark_black_1= (150,43,220)
light_black = (0,0,0)
dark_black = (180,255,46)
light_brown = (19,47,38)
dark_brown = (37,51,30)
light_blue = (100,43,56)
dark_blue = (124,255,255)
#设置不同颜色的mask,并将这些mask叠加,得到要分割的小丑鱼区域
red_mask=cv2.inRange(hsv_img,light_red_2,dark_red_2)
shade_mask=cv2.inRange(hsv_img, light_red_1, dark_red_1)
blue_mask=cv2.inRange(hsv_img,light_blue,dark_blue)
orange_mask = cv2.inRange(hsv_img,light_orange,dark_orange)
white_mask = cv2.inRange(hsv_img,light_white,dark_white)
black_mask=cv2.inRange(hsv_img,light_black,dark_black)
brown_mask=cv2.inRange(hsv_img,light_brown,dark_brown)
mask = orange_mask+white_mask+black_mask+blue_mask+shade_mask
#进行抠图,把小丑鱼区域从原图中抠出来
mask_inv=cv2.bitwise_not(mask)
result = cv2.bitwise_and(img, img,mask=mask)
#生成纯色图像的函数
def create_image(h,w,c):
img_create = np.zeros([h, w, c], np.uint8)
img_create[:, :,0] = np.zeros([h, w]) + 255#生成绿色 BGR
return img_create
#调用函数,并将两幅图像合成
green=create_image(470,640,3)
img1=cv2.bitwise_and(green,green,mask=mask_inv)
img2=cv2.bitwise_and(img,img,mask=mask)
final=cv2.add(img1,img2)
final=cv2.cvtColor(final,cv2.COLOR_RGB2BGR)
final1=final
cv2.imshow('final',final)
cv2.waitKey(0)
图像分割
1. (简答题)
请检测出下图中的不同大小的硬币数目与线条数
算法原理:
霍夫变换(Hough Transform)是图像处理中的一种特征提取技术,它通过一种投票算法检测具有特定形状的物体。Hough变换是图像处理中从图像中识别几何形状的基本方法之一。Hough变换的基本原理在于利用点与线的对偶性,将原始图像空间的给定的曲线通过曲线表达形式变为参数空间的一个点。这样就把原始图像中给定曲线的检测问题转化为寻找参数空间中的峰值问题。也即把检测整体特性转化为检测局部特性。比如直线、椭圆、圆、弧线等。
霍夫变换可以检测任何形状,但复杂的形状需要的参数就多,霍夫空间的维数就多,因此在程序实现上所需的内存空间以及运行效率上都不利于把标准霍夫变换应用于实际复杂图形的检测中。霍夫梯度法是霍夫变换的改进(圆检测),它的目的是减小霍夫空间的维度,提高效率。
霍夫圆变换:是将二维图像空间中一个圆转换为该圆半径、圆心横纵坐标所确定的三维参数空间中一个点的过程,因此,圆周上任意三点所确定的圆,经Hough变换后在三维参数空间应对应一点。该过程类似于选举投票过程,圆周上任意三个点为一选举人,而这三个点所确定的圆则为一侯选人(以下称为候选圆)。遍历圆周上所有点,任意三个点所确定的候选圆进行投票。遍历结束后,得票数最高点(理论上圆周上任意三点确定的圆在Hough变换后均对应三维参数空间中的同一点)所确定的圓即为该圆周上,绝大多数点所确定的圆(以下称为当选圆),即绝大多数点均在该当选圆的圆周上,以此确定该圆。
霍夫线变换:任何一条线都可以用(ρ,θ)这两个术语表示。首先创建2D数组或累加器(以保存两个参数的值),并将其初始设置为0。让行表示ρ,列表示θ。阵列的大小取决于所需的精度。假设您希望角度的精度为1度,则需要180列。对于ρ,最大距离可能是图像的对角线长度。因此,以一个像素精度为准,行数可以是图像的对角线长度。
实现思路:采用霍夫圆检测与霍夫线检测的函数进行图像中圆与线条的提取,得到圆与线条并计数。在实际操作中由于图像中的一条线条被硬币挡住后容易被检测成多条不同的线条,故定义了一个检测函数,当两条直线的ρ
和θ值的差值在一定范围内时,可以认为这两条直线为一条直线,此时将其合并成一条。之后再进行计数与画图操作。
关键代码:
读取并转换成灰度图:
img1 = cv2.imread("2-2.jpg")
gray = cv2.cvtColor(img1,cv2.COLOR_BGR2GRAY)#灰度图像
result=cv2.GaussianBlur(gray,(5,5),0)
def check(list1,theta,rou):
flag=0
if list1==[]:
return flag
for line in list1:
rou1,theta1=line
if abs(theta1-theta)<=0.05 and abs(rou1-rou)<=50 :
flag=1
return flag
edges = cv2.Canny(gray,50,150,apertureSize = 3)
line1=cv2.HoughLines(edges,1.3,np.pi/180,200)
line_test=[]
for line in line1[0:,]:
rou,theta=line[0]
line_test.append([rou,theta])
line_22=[]
for line in line_test:
rou,theta = line
if check(line_22,theta,rou)==0:
line_22.append(line)
line_num=0
num_coins=0
for i in circles[:]:
cv2.circle(img1,(i[0],i[1]),i[2],(255,0,0),5)#画圆
cv2.circle(img1,(i[0],i[1]),2,(255,0,255),10)#画圆心
num_coins+=1
for line in line_22:
rho = line[0]
theta=line[1]
a = np.cos(theta)
b = np.sin(theta)
x0 = a*rho
y0 = b*rho
x1 = int(x0 + 1000*(-b))
y1 = int(y0 + 1000*(a))
x2 = int(x0 - 1000*(-b))
y2 = int(y0 - 1000*(a))
cv2.line(img1,(x1,y1),(x2,y2),(0,0,255),2)
line_num+=1
print("硬币数量为",num_coins)
print("线条数量为",line_num)
cv2.imshow('result',img1)
cv2.waitKey(0)
import cv2 as cv
import numpy as np
img = cv.imread(cv.samples.findFile('2-2.jpg'))
gray = cv.cvtColor(img,cv.COLOR_BGR2GRAY)
edges = cv.Canny(gray,50,150,apertureSize = 3)
lines = cv.HoughLinesP(edges,1,np.pi/180,100,minLineLength=5,maxLineGap=20)
for line in lines:
x1,y1,x2,y2 = line[0]
cv.line(img,(x1,y1),(x2,y2),(0,255,0),2)
cv.imshow('result',img)
cv.waitKey(0)
形态学
请将图像中的小星系去除,只保留最大的星云
算法原理:
开运算就是先腐蚀后膨胀的过程,用来消除背景中的小点前景色噪声、平滑形状边界、断开物体之间的细小连接,当然不同的核矩阵效果会不同,有效元素为圆的核矩阵可以平滑边界、去除突刺。
闭运算有助于关闭前景物体上的小孔,或者小黑点。闭运算用来填充物体内的小空洞,连接邻近的物体,连接断开的轮廓线,平滑其边界的同时不改变面积。闭运算是先膨胀后腐蚀的过程。合理选择卷积核大小,太小了无法去除前景图的黑点。
题目要求将图像中的小星系去除,只保留最大的星云,星云可看作白色小点前景色噪声,背景可看作主要图像,所以将原图二值化后,对图像采用开运算去除背景中的小点前景色噪声,最终得到结果图。
关键代码
import cv2
import numpy as np
original_img = cv2.imread('star.png')
gray=cv2.cvtColor(original_img,cv2.COLOR_BGR2GRAY)
dst=255-gray
gray_res = gray
kernel = cv2.getStructuringElement(cv2.MORPH_ELLIPSE,(5,5)) #定义椭圆形结构元素
opened1 = cv2.morphologyEx(gray, cv2.MORPH_OPEN, kernel,iterations=1) #开运算,经过一轮迭代的
opened2 = cv2.morphologyEx(gray, cv2.MORPH_OPEN, kernel,iterations=5) #开运算,经过5轮迭代的
cv2.imshow("Open1", opened1)
cv2.imshow("Open2", opened2)
cv2.waitKey(0)
cv2.destroyAllWindows()
2. (简答题)
(灰度图形态学)请对下图分别进行开操作和闭操作,并对处理结果进行分析
算法原理:灰度图像与二值图像的区别在于其记录了灰度信息,所以,形态学处理的定义与二值图像有些不同,因为二值图像可以用一系列的二维坐标来表示图像信息,而灰度图需要一个三维坐标表示,而且二值图像中结构元SE是平坦的,没有灰度信息的,但灰度图中,结构元是可以带有第三维信息的,即结构元也是灰度的,这就带来了一些问题,因为二值图像中,形态学的输出结果完全由输入图像产生,但是结构元一旦引入灰度信息,那么输出结果将不再由输入图像唯一确定。所以,一般情况下,结构元都使用平坦的结构元。
腐蚀与膨胀是形态学的基本操作,在灰度图像中也是如此,在二值图像中腐蚀和膨胀定义为对图像进行 translation 以后的“与”和“或”的逻辑操作结果,在灰度图像中,为了保存灰度信息,“与”和“或”操作被对应的替换成了“最大值”和“最小值”操作这样就给出了灰度图像中腐蚀和膨胀的操作定义。
Erosion(腐蚀),使暗的地方会更暗(整体亮度减小),公式如下:
Dilation(膨胀),公式如下:
与腐蚀不同的就是与 SE(卷积核) 相加,取最大值。直观上,就是原图亮的地方会更亮(整体亮度增加),且范围扩大。
开操作是先进行腐蚀,再进行膨胀的过程。可以去除图像中小的明亮区域,使得整个图像的边界更加清晰。
闭操作是先进行膨胀,再进行腐蚀的过程。可以填补图像中一些小的暗区域,使得整个图像的连通性更加完整。
关键代码:
根据公式定义灰度图下的膨胀算子(取最大值):
def dilation_1(image, kernel):
# 获取结构元素的高度和宽度
kh, kw = kernel.shape
h,w=image.shape
# 创建一个与原始图像相同大小的全零矩阵
result = np.zeros(image.shape)
# 边界填充
padded_image = np.pad(image, (kh // 2, kw // 2), 'edge')
# 迭代遍历每个像素,进行膨胀操作
for i in range(kh // 2, padded_image.shape[0] - kh // 2):
for j in range(kw // 2, padded_image.shape[1] - kw // 2):
a=[]
for k in range(kh):
for l in range(kw):
if -1 < (i - kw + k) < h and -1 < (j - kh + l) < w:
a.append(image[i - kh + k, j - kw + l])
k2=max(a)
result[i-kw // 2, j-kw // 2] = k2
return result
def erosion_1(image, kernel):
# 获取结构元素的高度和宽度
kh, kw = kernel.shape
h,w=image.shape
# 创建一个与原始图像相同大小的全零矩阵
result = np.zeros(image.shape)
# 边界填充
padded_image = np.pad(image, (kh // 2, kw // 2), 'edge')
# 迭代遍历每个像素,进行腐蚀操作
for i in range(kh // 2, padded_image.shape[0] - kh // 2):
for j in range(kw // 2, padded_image.shape[1] - kw // 2):
a=[]
for k in range(kh):
for l in range(kw):
if -1 < (i - kw + k) < h and -1 < (j - kh + l) < w:
a.append(image[i - kh + k, j - kw + l])
k2=min(a)
result[i-kw // 2,j-kw // 2]=k2
return result
# 定义开操作
def opening(image, kernel):
return dilation_1(erosion_1(image, kernel), kernel)
# 定义闭操作
def closing(image, kernel):
return erosion_1(dilation_1(image, kernel), kernel)
#总程序
# 读入一幅灰度图像
img = cv2.imread('blb1.jpg', 0)
# 创建一个方形结构元素
kernel = np.ones((7, 7), np.uint8)
# 对图像进行开操作和闭操作
opening_img = opening(img/255.0 , kernel)*255.0
closing_img = closing(img / 255.0, kernel)*255.0
opening_img=opening_img.astype(np.uint8)
closing_img=closing_img.astype(np.uint8)
print(img)
print(opening_img)
print(closing_img)
cv2.imshow('raw',img)
cv2.imwrite('open.jpg',opening_img)
cv2.imwrite('close.jpg',closing_img)
cv2.waitKey(0)
综合题 珍珠:
算法原理:原图有些珍珠颜色和背景颜色很接近,影响珍珠图像的提取,对珍珠图像进行直方图绘制,结果如下所示:
可以看出珍珠颜色分布不均匀,大多数分布在200左右,背景颜色和珍珠颜色相似,难以检测,采用直方图均衡方法对图像进行处理。处理后的直方图图像和处理后图像如下所示:
图像预处理后采用霍夫圆检测对图像进行圆提取,提取珍珠图像并画出圆图像。
代码:
# 对图像进行灰度图处理
gray = cv2.cvtColor(img,cv2.COLOR_BGR2GRAY) # 将彩色图片转换为灰度图
# 绘制灰度直方图
hist1 = cv2.calcHist([img],[0],None,[256],[0,256]) # 计算灰度直方图
plt.figure() # 新建一个图像
plt.title("Grayscale Histogram") # 设置图像的标题
plt.xlabel("Bins") # 设置X轴标签
plt.ylabel("# of Pixels") # 设置Y轴标签
plt.plot(hist1) # 画图
plt.xlim([0,256]) # 设置x坐标轴范围
# 对直方图进行均衡化
equ = cv2.equalizeHist(gray) # 对灰度图进行直方图均衡化处理
# 绘制灰度均衡化以后的直方图
hist2 = cv2.calcHist([equ],[0],None,[256],[0,256]) # 计算均衡化后的灰度直方图
plt.figure() # 新建一个图像
plt.title("Grayscale Histogram") # 设置图像的标题
plt.xlabel("Bins") # 设置X轴标签
plt.ylabel("# of Pixels") # 设置Y轴标签
plt.plot(hist2) # 画图
plt.xlim([0,256]) # 设置x坐标轴范围
# 对图像进行二值化处理
ret,img1 = cv2.threshold(equ,200,255,cv2.THRESH_BINARY) # 对均衡化后的灰度图进行二值化处理
# 用霍夫变换检测圆形物体
circles1 = cv2.HoughCircles(equ,cv2.HOUGH_GRADIENT,1,55,param1=100,param2=25,minRadius=20,maxRadius=40)
print(circles1) # 打印出检测到的圆形物体的信息
circles = circles1[0,:,:] # 提取出检测到的圆形物体的二维数组
print(circles)
circles = np.uint16(np.around(circles)) # 对所有元素四舍五入并转换为无符号整型数据类型
for i in circles[:]:
# 画出每个圆形物体的外接圆和圆心
cv2.circle(img,(i[0],i[1]),i[2],(255,0,0),2) # 画圆
cv2.circle(img,(i[0],i[1]),2,(255,0,255),5) # 画圆心
目标检测
请用下面鸡蛋模板在原图中检测出对应的鸡蛋图像。
算法原理:基于特征的匹配算法,提取图像的特征生成特征描述子,根据描述子的相似程度对两幅图像的特征之间进行匹配。采用sift算法建立特征匹配算子,然后使用BFMatcher(暴力匹配法)特征匹配算法进行特征匹配,检测出鸡蛋的图像。
sift描述子:
暴力匹配法:先计算待匹配图像中每一对特征点的距离,然后按照距离从小到大的顺序排序,选取排名靠前的N个特征点作为匹配结果。通过设置匹配距离的阈值来筛选出合适的匹配结果,从而得到更加精准的匹配效果。
import cv2
import numpy as np
from matplotlib import pyplot as plt
# 加载原图像和目标图像
img1 = cv2.imread('eggs.png')
img2 = cv2.imread('origin.jpg')
# 将图像转换为灰度图像
gray1 = cv2.cvtColor(img1, cv2.COLOR_BGR2GRAY)
gray2 = cv2.cvtColor(img2, cv2.COLOR_BGR2GRAY)
# 创建SIFT特征检测器
sift = cv2.SIFT_create()
# 检测原图像和目标图像的关键点和特征描述符
kp1, des1 = sift.detectAndCompute(gray1, None)
kp2, des2 = sift.detectAndCompute(gray2, None)
# 创建BFMatcher匹配器
bf = cv2.BFMatcher(cv2.NORM_L2, crossCheck=True)
# 对原图像和目标图像的特征描述符进行匹配
matches = bf.match(des1, des2)
# 将匹配点按照距离排序
matches = sorted(matches, key=lambda x: x.distance)
# 仅保留良好的匹配点
good_matches = matches[:40]
# 将匹配点绘制在原图像和目标图像中
img3 = cv2.drawMatches(img1, kp1, img2, kp2, good_matches, None, flags=cv2.DrawMatchesFlags_NOT_DRAW_SINGLE_POINTS)
# 显示最终的结果
plt.imshow(cv2.cvtColor(img3, cv2.COLOR_BGR2RGB))
plt.show()
结果分析:图中和原图完全相似的鸡蛋被检测出来了,但是倾斜的鸡蛋还有缩小的鸡蛋没有被匹配出来,探究原因是因为采用的SIFT特征匹配算法在图像有旋转、缩放时匹配效果较差,此时可以采用相似性测度模板匹配的算法进行匹配。
请计算得到下图的HOG特征,并对其可视化。
算法原理:HOG(方向梯度直方图),通过计算和统计图像局部区域的梯度方向直方图来构成特征,在一幅图像中,局部目标的表象和形状能够被梯度或边缘的方向密度分布很好地描述。
实现方法:首先将图像分成小的连通区域,我们把它叫细胞单元(cell)。然后采集 细胞单元中各像素点的梯度的或边缘的方向直方图。最后把这些直方图组合起来就可 以构成特征描述器。
流程:1. 图像预处理:将图像进行灰度化,并对图像进行归一化处理。
2. 计算梯度:对图像进行求导,得到图像中每个像素点的梯度向量和梯度大小。
3. 划分cell:将图像分成多个cell,对每个cell内的像素点的梯度进行累加,统计每个cell内梯度角度的频率。这个过程可以采用直方图的方式实现。
4. 划分block:将相邻若干个cell组成block,对每个block内的cell的直方图进行归一化处理(L2范数),然后将不同block内的特征进行串联。
5. 特征向量归一化:对整个特征向量进行L2范数归一化。
可视化HOG特征的流程:
1. 将图像划分为多个cell,对每个cell内的梯度角度进行可视化,用箭头表示梯度方向。
2. 将相邻的cell组成block,对每个block的梯度频率直方图进行可视化,常用的方式是将不同方向的梯度频率反映在颜色上,例如灰度值。
3. 将整个图像的HOG特征向量可视化为一个二维图,例如使用PCA得到的二维坐标系表示特征向量的主要方向。
关键代码:
import cv2
import numpy as np
import math
import matplotlib.pyplot as plt
class Hog_descriptor():
# 初始化
# cell_size每个细胞单元的像素数
# bin_size表示把360分为多少边
def __init__(self, img, cell_size=16, bin_size=8):
self.img = img
self.img = np.sqrt(img / np.max(img))
self.img = img * 255
self.cell_size = cell_size
self.bin_size = bin_size
self.angle_unit = 360 / self.bin_size
# 获取hog向量和图片
def extract(self):
# 获得原图的shape
height, width = self.img.shape
# 计算原图的梯度大小
gradient_magnitude, gradient_angle = self.global_gradient()
gradient_magnitude = abs(gradient_magnitude)
# cell_gradient_vector用来保存每个细胞的梯度向量
cell_gradient_vector = np.zeros((int(height / self.cell_size), int(width / self.cell_size), self.bin_size))
height_cell,width_cell,_ = np.shape(cell_gradient_vector)
# 计算每个细胞的梯度直方图
for i in range(height_cell):
for j in range(width_cell):
# 获取这个细胞内的梯度大小
cell_magnitude = gradient_magnitude[i * self.cell_size:(i + 1) * self.cell_size,
j * self.cell_size:(j + 1) * self.cell_size]
# 获得这个细胞内的角度大小
cell_angle = gradient_angle[i * self.cell_size:(i + 1) * self.cell_size,
j * self.cell_size:(j + 1) * self.cell_size]
# 转化为梯度直方图格式
cell_gradient_vector[i][j] = self.cell_gradient(cell_magnitude, cell_angle)
# hog图像
hog_image = self.render_gradient(np.zeros([height, width]), cell_gradient_vector)
hog_vector = []
# block为2x2
for i in range(height_cell - 1):
for j in range(width_cell - 1):
block_vector = []
block_vector.extend(cell_gradient_vector[i][j])
block_vector.extend(cell_gradient_vector[i][j + 1])
block_vector.extend(cell_gradient_vector[i + 1][j])
block_vector.extend(cell_gradient_vector[i + 1][j + 1])
mag = lambda vector: math.sqrt(sum(i ** 2 for i in vector))
magnitude = mag(block_vector)
if magnitude != 0:
normalize = lambda block_vector, magnitude: [element / magnitude for element in block_vector]
block_vector = normalize(block_vector, magnitude)
hog_vector.append(block_vector)
return hog_vector, hog_image
# 计算原图的梯度大小
# 角度大小
def global_gradient(self):
gradient_values_x = cv2.Sobel(self.img, cv2.CV_64F, 1, 0, ksize=5)
gradient_values_y = cv2.Sobel(self.img, cv2.CV_64F, 0, 1, ksize=5)
gradient_magnitude = cv2.addWeighted(gradient_values_x, 0.5, gradient_values_y, 0.5, 0)
gradient_angle = cv2.phase(gradient_values_x, gradient_values_y, angleInDegrees=True)
return gradient_magnitude, gradient_angle
# 分解角度信息到不同角度的直方图上
def cell_gradient(self, cell_magnitude, cell_angle):
orientation_centers = [0] * self.bin_size
for i in range(cell_magnitude.shape[0]):
for j in range(cell_magnitude.shape[1]):
gradient_strength = cell_magnitude[i][j]
gradient_angle = cell_angle[i][j]
min_angle, max_angle, mod = self.get_closest_bins(gradient_angle)
orientation_centers[min_angle] += (gradient_strength * (1 - (mod / self.angle_unit)))
orientation_centers[max_angle] += (gradient_strength * (mod / self.angle_unit))
return orientation_centers
# 计算每个像素点所属的角度
def get_closest_bins(self, gradient_angle):
idx = int(gradient_angle / self.angle_unit)
mod = gradient_angle % self.angle_unit
return idx, (idx + 1) % self.bin_size, mod
# 将梯度直方图进行绘图
def render_gradient(self, image, cell_gradient):
cell_width = self.cell_size / 2
max_mag = np.array(cell_gradient).max()
for x in range(cell_gradient.shape[0]):
for y in range(cell_gradient.shape[1]):
cell_grad = cell_gradient[x][y]
cell_grad /= max_mag
angle = 0
angle_gap = self.angle_unit
for magnitude in cell_grad:
angle_radian = math.radians(angle)
x1 = int(x * self.cell_size + magnitude * cell_width * math.cos(angle_radian))
y1 = int(y * self.cell_size + magnitude * cell_width * math.sin(angle_radian))
x2 = int(x * self.cell_size - magnitude * cell_width * math.cos(angle_radian))
y2 = int(y * self.cell_size - magnitude * cell_width * math.sin(angle_radian))
cv2.line(image, (y1, x1), (y2, x2), int(255 * math.sqrt(magnitude)))
angle += angle_gap
return image
img = cv2.imread('cats.jpg', cv2.IMREAD_GRAYSCALE)
hog = Hog_descriptor(img, cell_size=30, bin_size=12)
vector, image = hog.extract()
plt.imshow(image, cmap=plt.cm.gray)
plt.show()
import cv2
from skimage import feature as ft
from matplotlib import pyplot as plt
if __name__=='__main__':
img = cv2.imread('cats.jpg', cv2.IMREAD_GRAYSCALE)
features = ft.hog(img,orientations=6,pixels_per_cell=[20,20],cells_per_block=[2,2],visualize=True)
cv2.imshow("original",img)
cv2.waitKey(0)
plt.imshow(features[1],cmap=plt.cm.gray)
plt.show()