第5章 Python 数字图像处理(DIP) - 图像复原与重建7 - 周期噪声 余弦噪声生成方法

周期噪声

周期噪声通常是在获取图像期间由电气或机电干扰产生的

def add_sin_noise(img, scale=1, angle=0):
    """
    add sin noise for image
    param: img: input image, 1 channel, dtype=uint8
    param: scale: sin scaler, smaller than 1, will enlarge, bigger than 1 will shrink
    param: angle: angle of the rotation
    return: output_img: output image is [0, 1] image which you could use as mask or any you want to
    """
    
    height, width = img.shape[:2]  # original image shape
    
    # convert all the angle
    if int(angle / 90) % 2 == 0:
        rotate_angle = angle % 90
    else:
        rotate_angle = 90 - (angle % 90)
    
    rotate_radian = np.radians(rotate_angle)    # convert angle to radian
    
    # get new image height and width
    new_height = int(np.ceil(height * np.cos(rotate_radian) + width * np.sin(rotate_radian)))
    new_width = int(np.ceil(width * np.cos(rotate_radian) + height * np.sin(rotate_radian))) 
    
    # if new height or new width less than orginal height or width, the output image will be not the same shape as input, here set it right
    if new_height < height:
        new_height = height
    if new_width < width:
        new_width = width
    
    # meshgrid
    u = np.arange(new_width)
    v = np.arange(new_height)
    u, v = np.meshgrid(u, v)
    
    # get sin noise image, you could use scale to make some difference, better you could add some shift
#     noise = abs(np.sin(u * scale))
    noise = 1 - np.sin(u * scale)
    
    # here use opencv to get rotation, better write yourself rotation function
    C1 = cv2.getRotationMatrix2D((new_width/2.0, new_height/2.0), angle, 1)
    new_img = cv2.warpAffine(noise, C1, (int(new_width), int(new_height)), borderValue=0)
    
    # ouput image should be the same shape as input, so caculate the offset the output image and the new image
    # I make new image bigger so that it will cover all output image
    
    offset_height = abs(new_height - height) // 2
    offset_width = abs(new_width - width) // 2
    img_dst = new_img[offset_height:offset_height + height, offset_width:offset_width+width]
    output_img = normalize(img_dst)
    
    return output_img 

def spectrum_fft(fft):
    """
    return FFT spectrum
    """
    return np.sqrt(np.power(fft.real, 2) + np.power(fft.imag, 2))

# 周期噪声
img_ori = cv2.imread('DIP_Figures/DIP3E_Original_Images_CH05/Fig0507(a)(ckt-board-orig).tif', 0) #直接读为灰度图像

# 正弦噪声
noise = add_sin_noise(img_ori, scale=0.35, angle=-20)
img = np.array(img_ori / 255, np.float32)
img_noise = img + noise
img_noise = np.uint8(normalize(img_noise)*255)

# 频率域中的其他特性
# FFT
img_fft = np.fft.fft2(img_noise.astype(np.float32))
# 中心化
fshift = np.fft.fftshift(img_fft)            # 将变换的频率图像四角移动到中心
# 中心化后的频谱
spectrum_fshift = spectrum_fft(fshift)
spectrum_fshift_n = np.uint8(normalize(spectrum_fshift) * 255)

# 对频谱做对数变换
spectrum_log = np.log(1 + spectrum_fshift)

plt.figure(figsize=(15, 10))
plt.subplot(121), plt.imshow(img_noise, 'gray'), plt.title('With Sine noise'), plt.xticks([]),plt.yticks([])
plt.subplot(122), plt.imshow(spectrum_log, 'gray'), plt.title('Spectrum'), plt.xticks([]),plt.yticks([])
# 在图像上加上箭头
plt.arrow(180, 180, 25, 30, width=5,length_includes_head=True, shape='full')
plt.arrow(285, 265, -25, -30, width=5,length_includes_head=True, shape='full')
plt.tight_layout()
plt.show()