4. 原文作者有matlab程序,也可以结合大佬的python程序来理解,python程序来源:https://github.com/zsyOAOA/noise_est_ICCV2015
#!/usr/bin/env python
# -*- coding:utf-8 -*-
# Power by Zongsheng Yue 2019-01-07 14:36:55
import sys
import cv2
import numpy as np
from matplotlib import pyplot as plt
from skimage import img_as_float
import time
def im2patch(im, pch_size, stride=8):
'''
Transform image to patches.
Input:
im: 3 x H x W or 1 X H x W image, numpy format
pch_size: (int, int) tuple or integer
stride: (int, int) tuple or integer
'''
if isinstance(pch_size, tuple):
pch_H, pch_W = pch_size
elif isinstance(pch_size, int):
pch_H = pch_W = pch_size
else:
sys.exit('The input of pch_size must be a integer or a int tuple!')
if isinstance(stride, tuple):
stride_H, stride_W = stride
elif isinstance(stride, int):
stride_H = stride_W = stride
else:
sys.exit('The input of stride must be a integer or a int tuple!')
C, H, W = im.shape
num_H = len(range(0, H-pch_H+1, stride_H))
num_W = len(range(0, W-pch_W+1, stride_W))
num_pch = num_H * num_W
pch = np.zeros((C, pch_H*pch_W, num_pch), dtype=im.dtype)
kk = 0
for ii in range(pch_H):
for jj in range(pch_W):
temp = im[:, ii:H-pch_H+ii+1:stride_H, jj:W-pch_W+jj+1:stride_W]
pch[:, kk, :] = temp.reshape((C, num_pch))
kk += 1
return pch.reshape((C, pch_H, pch_W, num_pch))
def noise_estimate(im, pch_size=8):
'''
Implement of noise level estimation of the following paper:
Chen G , Zhu F , Heng P A . An Efficient Statistical Method for Image Noise Level Estimation[C]// 2015 IEEE International Conference
on Computer Vision (ICCV). IEEE Computer Society, 2015.
Input:
im: the noise image, H x W x 3 or H x W numpy tensor, range [0,1]
pch_size: patch_size
Output:
noise_level: the estimated noise level
'''
if im.ndim == 3:
im = im.transpose((2, 0, 1))
else:
im = np.expand_dims(im, axis=0)
# image to patch
pch = im2patch(im, pch_size, 3) # C x pch_size x pch_size x num_pch tensor
num_pch = pch.shape[3]
pch = pch.reshape((-1, num_pch)) # d x num_pch matrix
d = pch.shape[0]
mu = pch.mean(axis=1, keepdims=True) # d x 1
X = pch - mu
sigma_X = np.matmul(X, X.transpose()) / num_pch
sig_value, _ = np.linalg.eigh(sigma_X)
sig_value.sort()
t = 0
for ii in range(-1, -d-1, -1):
tau = np.mean(sig_value[:ii])
if np.sum(sig_value[:ii]>tau) == np.sum(sig_value[:ii] < tau):
t = np.sqrt(tau)
break
return sig_value, t
if __name__ == '__main__':
# file = r'D:\code\data\kodak_IDR\kodim24_gaussian_NOISY_25.0.png'
#file = r'D:\code\data\kodak2\kodim12_gaussian_0.1.png'
# file = r'D:\dataset\noise_img\cap_frame_2016_D65_3000lx_raw.jpg'
# file = r'D:\dataset\noise_img\cap_frame_0001_raw.png'
file = r'a.jpg'
im = cv2.imread(file)
#im = cv2.imread('./lena.png')
im = img_as_float(im)
noise_level = [0, 0.5, 1, 2, 5, 15, 20, 30, 45, 50, 70]
for level in noise_level:
sigma = level / 255
im_noise = im + np.random.randn(*im.shape) * sigma
start = time.time()
eigval, est_level = noise_estimate(im_noise, 8)
end = time.time()
time_elapsed = end -start
str_p = "Time: {0:.4f}, Ture Level: {1:6.4f}, Estimated Level: {2:6.4f}"
print(str_p.format(time_elapsed, level, est_level*255))
plt.plot(255* np.sqrt(eigval)[::-1][:])
plt.grid()
plt.show()
- 作者提到 关于patch size的选择
patch_size比较大时,beta值越大,意味着可以允许比较多的主要特征比例。
但是另一方面, patch的数量越大的时候,冗余特征值的分布才符合
因此size 和 number 需要一个折衷。