I'm studying about comparison between images to determine quality. I've learned about MSE, SNR and PSNR, and now I'm trying to understand WSNR, which I assume is similar to SNR but with weights.
I have a Python code to calculate the WSNR between two images, but I don't fully understand it and I don't know what is the mathematical definition of this measure. I can't find any paper or book that explain how it is defined and calculated.
Question:
How is WSNR (weighted signal to noise ratio) defined?
Python Code:
import numpy as np
from scipy.ndimage.filters import gaussian_filter as __gaussian_filter
from scipy.ndimage.filters import convolve as __convolve
from scipy.ndimage.filters import correlate as __correlate
from scipy.fftpack import fftshift as __fftshift
def wsnr(reference, query):
"""Computes the Weighted Signal to Noise Ratio (WSNR) metric.
value = wsnr(reference, query)
inputs
----------
reference: original image data.
query : modified image data to be compared.
output
----------
value : wsnr value
"""
def __genetate_meshgrid(x, y):
f = lambda u: u / 2 + 0.5 - 1
(H, W) = map(f, (x, y))
return (H, W)
def __create_complex_planes(x, y):
(H, W) = __genetate_meshgrid(x, y)
(xplane, yplane) = np.mgrid[-H:H + 1, -W:W + 1]
return (xplane, yplane)
def __get_evaluated_contrast_sensivity(plane):
w = 0.7
angle = np.angle(plane)
return ((1.0 - w) / 2.0) * np.cos(4.0 * angle) + (1.0 + w) / 2.0
def __get_radial_frequency(x, y):
(xplane, yplane) = __create_complex_planes(x, y)
nfreq = 60
plane = (xplane + 1.0j * yplane) / x * 2.0 * nfreq
s = __get_evaluated_contrast_sensivity(plane)
radfreq = abs(plane) / s
return radfreq
def __generate_CSF(radfreq):
a = -((0.114 * radfreq) ** 1.1)
csf = 2.6 * (0.0192 + 0.114 * radfreq) * np.exp(a)
f = radfreq < 7.8909
csf[f] = 0.9809
return csf
def __weighted_fft_domain(ref, quer, csf):
err = ref.astype('double') - quer.astype('double')
err_wt = __fftshift(np.fft.fft2(err)) * csf
im = np.fft.fft2(ref)
return (err, err_wt, im)
def __get_weighted_error_power(err_wt):
return (err_wt * np.conj(err_wt)).sum()
def __get_signal_power(im):
return (im * np.conj(im)).sum()
def __get_ratio(mss, mse):
if mse != 0:
ratio = 10.0 * np.log10(mss / mse)
else:
ratio = float("inf")
return np.real(ratio)
if not len(reference.shape) < 3:
reference = __convert_to_luminance(reference)
query = __convert_to_luminance(query)
size = reference.shape
(x, y) = (size[0], size[1])
radfreq = __get_radial_frequency(x, y)
csf = __generate_CSF(radfreq)
(err, err_wt, im) = __weighted_fft_domain(reference, query, csf)
mse = __get_weighted_error_power(err_wt)
mss = __get_signal_power(im)
ratio = __get_ratio(mss, mse)
return ratio