数字图像处理大作业

# -*- coding: utf-8 -*-
import cv2
import numpy as np
import matplotlib.pyplot as plt

# 分别提取图像的红、绿、蓝三个位面并显示
def Q1(img):
    # cv2读进来是bgr
    # 调用通道分离的包实现
    b,g,r = cv2.split(img)

    h, w = r.shape
    # 不调用通道分离的包实现
    b_1 = img.copy()
    g_1 = img.copy()
    r_1 = img.copy()
    for i in range(h):
        for j in range(w):
            # 想要用rgb输出
            r_1[i][j] = (img[i][j][2],0,0)
            g_1[i][j] = (0,img[i][j][1],0)
            b_1[i][j] = (0,0,img[i][j][0])

    # 显示为一张图片
    # 显示中文
    plt.rcParams['font.sans-serif'] = ['SimHei']
    # 显示原图
    plt.subplot(3, 3, 1)
    plt.title('原图')
    plt.imshow(img[:,:,::-1])
    plt.xticks([]), plt.yticks([])
    # 直接按通道输出
    plt.subplot(3, 3, 4)
    plt.title('红色_调包')
    plt.imshow(r)
    plt.xticks([]), plt.yticks([])
    plt.subplot(3, 3, 5)
    plt.title('绿色_调包')
    plt.imshow(g)
    plt.xticks([]), plt.yticks([])
    plt.subplot(3, 3, 6)
    plt.title('蓝色_调包')
    plt.imshow(b)
    plt.xticks([]), plt.yticks([])
    # 不调用
    plt.subplot(3, 3, 7)
    plt.title('红色_循环实现')
    plt.imshow(r_1[:,:,:])
    plt.xticks([]), plt.yticks([])
    plt.subplot(3, 3, 8)
    plt.title('绿色_循环实现')
    plt.imshow(g_1[:,:,:])
    plt.xticks([]), plt.yticks([])
    plt.subplot(3, 3, 9)
    plt.title('蓝色_循环实现')
    plt.imshow(b_1[:,:,:])
    plt.xticks([]), plt.yticks([])
    # 展示
    plt.show()


# 将彩色图像从RGB空间转换到HSI空间，在同一个图像窗口显示两个彩色空间的图像
def Q2(img):
    # 复制一下图片
    img_new = img.copy()
    # 得到长、宽
    h = np.shape(img)[0]
    w = np.shape(img)[1]
    # 1/255 = 0.0039 所以为了防止除以0，给每个点增加了0.001
    B, G, R = cv2.split(img)
    [B, G, R] = [i/255.0+0.001 for i in ([B, G, R])]
    # 公式:https://ask.qcloudimg.com/http-save/7151457/rn6syd7oad.png?imageView2/2/w/1620
    # I可以直接求
    I = (R + G + B) / 3.0  # 计算I通道
    # H
    H = np.zeros((h, w))  # 定义H通道
    for i in range(h):
        den = np.sqrt((R[i] - G[i]) ** 2 + (R[i] - B[i]) * (G[i] - B[i]))+0.001
        # 同理，den把上述的0.001消除了，所以这里重新增加一个不影响的值
        thetha = np.arccos(0.5 * (R[i] - B[i] + R[i] - G[i]) / den)  # 计算夹角
        temp = np.zeros(w)  # 定义临时数组
        # den>0且G>=B的元素h赋值为thetha
        temp[B[i] <= G[i]] = thetha[B[i] <= G[i]]
        # den>0且G<=B的元素h赋值为thetha
        temp[G[i] < B[i]] = 2 * np.pi - thetha[G[i] < B[i]]
        # den<0的元素h赋值为0
        temp[den == 0] = 0
        H[i] = temp / (2 * np.pi)  # 弧度化后赋值给H通道

    # S
    S = np.zeros((h, w))  # 定义H通道

    for i in range(h):
        min = []
        # 找出每组RGB值的最小值
        for j in range(w):
            arr = [B[i][j], G[i][j], R[i][j]]
            min.append(np.min(arr))
        # 计算S通道
        min = np.array(min)
        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]之间
    img_new[:, :, 0] = H * 255
    img_new[:, :, 1] = S * 255
    img_new[:, :, 2] = I * 255
    # 显示中文
    plt.rcParams['font.sans-serif'] = ['SimHei']
    # 显示原图
    plt.subplot(1, 2, 1)
    plt.title('原图')
    plt.imshow(img[:,:,::-1])
    plt.subplot(1, 2, 2)
    plt.title('HSI图像')
    plt.imshow(img_new[:,:,:])
    # 展示
    plt.show()

# 将彩色图像转换成灰度图像，并进行直方图规定化操作，显示规定化操作前后的图像；
# 累积直方图的制作
def leiji(img):
    # 初始化
    x = [0] * 256
    y = [0] * 256
    prob = [0] * 256
    for i in range(256):
        x[i] = i
        y[i] = 0

    for rv in img:
        # 算出现次数
        for cv in rv:
            y[cv] += 1

    # 长和宽
    h, w = img.shape
    for i in range(256):
        # 算出现概率
        prob[i] = y[i] / (h * w)

    # 累积直方图准备工作
    prob_sum = [0] * 256
    # 第0个值
    prob_sum[0] = prob[0]
    for i in range(1, 256):
        prob_sum[i] = prob_sum[i - 1] + prob[i]
    return prob_sum

def Q3(img1,img2):
    # 将彩色图像转化为灰度图像
    gray_img1 = cv2.cvtColor(img1, cv2.COLOR_BGR2GRAY) # 待转换图像
    gray_img2 = cv2.cvtColor(img2, cv2.COLOR_BGR2GRAY) # 目标转换成这样的图像

    # 颜色累计直方图
    x = range(256)
    img1_leiji = leiji(gray_img1)
    img2_leiji = leiji(gray_img2)

    # 构造一个256*256的空表
    abs_chazhi = [[0 for i in range(256)] for j in range(256)]
    for i in range(256):
        for j in range(256):
            abs_chazhi[i][j] = abs(img1_leiji[i] - img2_leiji[j])

    # 构造灰度级映射
    img_map = [0] * 256
    for i in range(256):
        zuixiao_n = abs_chazhi[i][0]
        index = 0
        for j in range(256):
            if zuixiao_n > abs_chazhi[i][j]:
                zuixiao_n = abs_chazhi[i][j]
                index = j
        img_map[i] = ([i, index])
    # 进行映射，实现灰度规定化
    h, w = gray_img1.shape
    img_new = gray_img1.copy()
    for i in range(h):
        for j in range(w):
            img_new[i, j] = img_map[gray_img1[i, j]][1]
    # 显示中文
    plt.rcParams['font.sans-serif'] = ['SimHei']
    # 显示原图
    plt.subplot(4, 2, 1)
    plt.title('原图1')
    plt.imshow(img1[:,:,::-1])
    plt.subplot(4, 2, 2)
    plt.title('原图2')
    plt.imshow(img2[:,:,::-1])
    plt.subplot(4, 2, 3)
    plt.title('灰度图1')
    plt.imshow(gray_img1,'gray')
    plt.subplot(4, 2, 4)
    plt.title('灰度图2')
    plt.imshow(gray_img2,'gray')
    plt.subplot(4, 2, 5)
    plt.title('累积直方图1')
    plt.bar(x, img1_leiji, color='orange')
    plt.subplot(4, 2, 6)
    plt.title('累积直方图2')
    plt.bar(x, img2_leiji, color='green')
    plt.subplot(4, 2, 7)
    plt.title('灰度图1')
    plt.imshow(gray_img1,'gray')
    plt.subplot(4, 2, 8)
    plt.title('规定化后')
    plt.imshow(img_new,'gray')
    # 展示
    plt.show()

# 给图像外层加一圈 这里采用+最外层一圈而不是补0
def jiayiquan(h,w,img):
    a = []
    for i in range(h + 2):
        a.append([0] * (w + 2))
    for i in range(h):
        for j in range(w):
            a[i + 1][j + 1] = img[i][j]

    # 完成在原图img1 周围补一圈 变成a
    # 四个角
    a[0][0] = a[1][1]  # 左上角
    a[0][w + 1] = a[1][w]  # 右上角
    a[h + 1][0] = a[h][1]  # 左下角
    a[h + 1][w + 1] = a[h][w]  # 右下角

    # 左右
    for i in range(1, h + 1):
        a[i][0] = a[i][1]
        a[i][w + 1] = a[i][w]
    # 上下
    for j in range(1, w + 1):
        a[0][j] = a[1][j]
        a[h + 1][j] = a[h][j]
    return a

def huifu(h,w,img):
    a = np.zeros((h,w))
    for i in range(h-2):
        for j in range(w-2):
            a[i][j]=img[i+1][j+1]
    return a

def juanji(is_ditong,h,w,img_a,mod):
    b = np.zeros((h+2,w+2))
    if is_ditong:
        ditong_sum = sum(mod)
    else:
        ditong_sum = 1

    for i in range(1,h+1):
        for j in range(1,w+1):
            c = int((mod[0] * img_a[i-1][j-1] + mod[1] * img_a[i-1][j] + mod[2] * img_a[i-1][j+1] +
                    mod[3] * img_a[i][j - 1] + mod[4] * img_a[i][j]   + mod[5] * img_a[i][j+1] +
                    mod[6] * img_a[i+1][j-1] + mod[7] * img_a[i+1][j] + mod[8] * img_a[i+1][j+1])/ditong_sum-0.5)
            if c < 0:
                b[i][j] = 0
            elif c > 255:
                b[i][j] = 255
            else:
                b[i][j] = c
        # 恢复为图片原大小
    b = huifu(h,w,b)
    return b

# 分别设计低通、高通滤波模板，并对图像进行滤波操作，显示原图、平滑图像和锐化图像；
def Q4(img):
    # 转化为灰度图
    img1 = cv2.cvtColor(img, cv2.COLOR_BGR2GRAY)
    h, w = img1.shape

    # 创建模板
    # model_d是低通滤波模板
    # 均值卷积核
    model_d_1 = [1, 1, 1,
                1, 1, 1,
                1, 1, 1] # 还要除以1/9
    # 高斯卷积核
    model_d_2 = [1, 2, 1,
                 2, 4, 2,
                 1, 2, 1] # 还需要除1/16

    # model_g是低通滤波模板
    # 中间加1来给
    model_g_1 = [-1, -1, -1,
                 -1,  9, -1,
                 -1, -1, -1] # 锐化卷积核1，中间+1在原图的基础上叠加了卷积核2
    model_g_2 = [-1, -1, -1,
                 -1,  8, -1,
                 -1, -1, -1] # 锐化卷积核2

    # 给图片周围增加一圈,使得 高[0,h+1]（h+2） 宽[0,w+1](w+1)
    a = jiayiquan(h,w,img1)

    # 卷积操作
    # 均值卷积核 111 111 111
    d1 = juanji(True,h,w,a,model_d_1)
    # 高斯卷积核 121 242 121
    d2 = juanji(True,h,w,a,model_d_2)
    # 锐化卷积1 -1-1-1 -1 9-1 -1-1-1
    g1 = juanji(False,h,w,a,model_g_1)
    # 锐化卷积2 -1-1-1 -1 8-1 -1-1-1
    g2 = juanji(False,h,w,a,model_g_2)

    # 显示中文
    plt.rcParams['font.sans-serif'] = ['SimHei']
    # 显示原图
    plt.subplot(3, 2, 1)
    plt.title('原图')
    plt.imshow(img[:,:,::-1])
    plt.subplot(3, 2, 2)
    plt.title('灰度图')
    plt.imshow(img1,'gray')
    plt.subplot(3, 2, 3)
    plt.title('111 111 111滤波')
    plt.imshow(d1,'gray')
    plt.subplot(3, 2, 4)
    plt.title('121 242 121滤波')
    plt.imshow(d2, 'gray')
    plt.subplot(3, 2, 5)
    plt.title('-1绕9滤波')
    plt.imshow(g1, 'gray')
    plt.subplot(3, 2, 6)
    plt.title('-1绕8滤波')
    plt.imshow(g2, 'gray')
    # 展示
    plt.show()

# 设计高斯高通滤波器，对图像在频域进行高通滤波，显示原始图像、频谱图及滤波后图像和频谱图；
def Q5(img):
    img1 = cv2.cvtColor(img, cv2.COLOR_BGR2GRAY)
    h,w = img1.shape

    # 快速傅里叶变换算法得到频率分布
    f = np.fft.fft2(img1)
    # 频谱图1
    fimg1 = np.log(np.abs(f))

    # 调用fftshift()函数转移到中间位置
    fshift = np.fft.fftshift(f)
    # 频谱图2
    fimg2 = np.log(np.abs(fshift))

    # 高通滤波 ---> 低频部分置0
    # 获得中心点
    half_h = int(h / 2)
    half_w = int(w / 2)
    # 中心点附近置零
    fshift1 = fshift.copy()
    fshift1[half_h - 30: half_h + 30, half_w - 30: half_w + 30] = 0
    # 高通滤波频谱图
    fimg_g = np.abs(fshift1)
    fimg_g = np.log(fimg_g)
    # 傅里叶逆变换
    ishift = np.fft.ifftshift(fshift1)
    img_new1 = np.fft.ifft2(ishift)
    img_new1 = np.abs(img_new1)


    # 低通滤波---> 高频部分 置0
    fshift2 = fshift.copy()
    for i in range(h):
        for j in range(w):
            if (i<half_h-30 or i >half_h+30) or (j<half_w-30 or j>half_w+30):
                fshift2[i][j]=0
    # 低通滤波频谱图
    fimg_d = np.abs(fshift2)
    fimg_d = np.log(fimg_d)
    # 傅里叶逆变换
    ishift = np.fft.ifftshift(fshift2)
    img_new2 = np.fft.ifft2(ishift)
    img_new2 = np.abs(img_new2)

    # 显示中文
    plt.rcParams['font.sans-serif'] = ['SimHei']
    # 显示原图
    plt.subplot(2, 4, 1)
    plt.title('原图')
    plt.imshow(img[:,:,::-1])
    plt.subplot(2, 4, 2)
    plt.title('灰度图')
    plt.imshow(img1,'gray')
    plt.subplot(2, 4, 3)
    plt.title('频谱图_1')
    plt.imshow(fimg1,'gray')
    plt.subplot(2, 4, 4)
    plt.title('频谱图_2')
    plt.imshow(fimg2,'gray')
    plt.subplot(2, 4, 5)
    plt.title('高通-频谱图')
    plt.imshow(fimg_g,'gray')
    plt.subplot(2, 4, 6)
    plt.title('高通-处理后')
    plt.imshow(img_new1,'gray')
    plt.subplot(2, 4, 7)
    plt.title('低通-频谱图')
    plt.imshow(fimg_d,'gray')
    plt.subplot(2, 4, 8)
    plt.title('低通-处理后')
    plt.imshow(img_new2,'gray')

    # 展示
    plt.show()

# 对图像采用canny算子进行边缘检测，并对检测后的图像进行开和闭数学形态学运算，显示边缘图像及开闭运算图像。
def Q6(img):
    img1 = cv2.cvtColor(img, cv2.COLOR_BGR2GRAY)
    # 高斯滤波
    gaussian = cv2.GaussianBlur(img1,(3,3),0)
    # Canny算子
    Canny = cv2.Canny(gaussian, 50, 150)
    # 使用闭运算连接中断的图像前景,迭代运算三次
    # 调包实现
    kai = cv2.morphologyEx(Canny, cv2.MORPH_CLOSE, kernel=(3, 3), iterations=3)
    bi = cv2.morphologyEx(Canny, cv2.MORPH_OPEN, kernel=(3, 3), iterations=3)
    # 显示中文
    plt.rcParams['font.sans-serif'] = ['SimHei']
    # 显示原图
    plt.subplot(2, 3, 1)
    plt.title('原图')
    plt.imshow(img[:,:,::-1])
    plt.subplot(2, 3, 2)
    plt.title('灰度图')
    plt.imshow(img1,'gray')
    plt.subplot(2, 3, 4)
    plt.title('Canny边缘')
    plt.imshow(Canny,'gray')
    plt.subplot(2, 3, 5)
    plt.title('开运算')
    plt.imshow(kai,'gray')
    plt.subplot(2, 3, 6)
    plt.title('闭运算')
    plt.imshow(bi,'gray')
    # 展示
    plt.show()


img = cv2.imread('b.jpg',1)
Q1(img) # finished
Q2(img) # finished

img1 = cv2.imread('img1.png',1)
img2 = cv2.imread('img2.png',1)
Q3(img1,img2) # finished

Q4(img) # finished
Q5(img) # finished

img3 = cv2.imread('jiaoyan2.png',1)
Q6(img3) # finished
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452