Ink
Ink

Reputation: 11

Gaussian filter bug in scipy_filters python

I want to apply a Gaussian filter on an float numpy array using python2. In fact I'm trying to rewrite the code Retrospective Correction using Homomorphic Filtering in python, g(x,y) = exp(LPF(log(f(x,y)))) . C

where LPF(f(x,y)) is the low-pass filter of f(x,y) and C is the normalization coefficient.

Here is my code, which would crash at line:LPFF = scipy_filters.gaussian_filter1d(logF, sigma=10)

fxy = greyscale_matrix
logF = np.log1p(fxy)
print logF
LPFF = scipy_filters.gaussian_filter1d(logF, sigma=10)
print LPFF
expF = np.expm1(LPFF)
tmp_div = fxy / expF
meanF_down = np.mean(tmp_div)
meanF_top = np.mean(fxy)
C = meanF_top / meanF_down
gx_matrix = expF * C
rtn_matrix = greyscale_matrix - gx_matrix

Image.fromarray(gx_matrix).show()
Image.fromarray(rtn_matrix).show()

I try to debug it, when going to

def correlate1d(*args, **kwargs): # real signature unknown
    pass

from _nd_image.py the code seems crash with no idea. I also find a github doing this job as Implementation of homomorphic filter made in Python, but something is wrong with it.

Source img: sourceimg

Using this may give me a filtered img like:filtered img Of course by changing the attributes I can make it lighter, but the excessive shading trend is obviously still exist.

Could anybody shed some light on this problem? Thx a lot!

Upvotes: 0

Views: 1490

Answers (1)

Warren Weckesser
Warren Weckesser

Reputation: 114946

You don't say which version of SciPy that you are using, but I suspect you are hitting a bug that should have been fixed in SciPy version 1.0.1 (and 1.1.0). The bug causes the interpreter to crash when the input to the Gaussian filter function has data type numpy.float16. If your input matrix fxy has type numpy.uint8 (which is typical for an image), then numpy.log1p(fxy) returns an array with data type numpy.float16, e.g.:

In [29]: fxy = np.array([[0, 10], [128, 255]], dtype=np.uint8)

In [30]: logF = np.log1p(fxy)

In [31]: logF
Out[31]: 
array([[0.   , 2.398],
       [4.86 , 5.547]], dtype=float16)

logF has type numpy.float16, and the crash occurs when you pass this to gaussian_filter1d.

A quick fix is to convert fxy to, say, 64 bit floating point before passing it to np.log1p:

In [34]: logF = np.log1p(fxy.astype(np.float64))

In [35]: logF.dtype
Out[35]: dtype('float64')

Upvotes: 1

Related Questions