[SciPy-User] Matlab imfilter/fspecial equivalent

[Brickle Macho]

I am getting closer.  Here are some observations:

* If I don't imresize() the original image, then visually I can't tell 
the difference in the output image from Matlab or Scipy.
* Empirically the results using different kernels are the same. That is 
the difference between two images is essentially zero. That is doesn't 
seem to matter if I use the imported the "fspecial()" kernel from 
matlab;  generate my own kernel using a fgaussion(), or use 
scipy.ndimage.gaussian_filter().

So look like the problem is with ndimage.resize().

Michael.
--





On 16/05/13 1:25 PM, Travis Oliphant wrote:
> I accidentally hit send before I finished the original email (Gmail 
> interface quirk strikes again...)  I'll including the entire original 
> email....
>
> Hi Michael,
>
> Thanks for sharing your code and example.     The imfilter command is 
> equivalent to scipy.signal.correlate and scipy.ndimage.correlate (the 
> one in scipy.ndimage is faster I believe).   These implement simple 
> correlation-based filtering given a finite kernel.
>
> To get the same output you would need to generate the same kind of 
> kernel in Python as the Matlab fspecial command is producing.  one 
> approach to get you started and to help separate the problem into 2 
> stages (reproducing the imfilter and reproducing the fspecial filter) 
> is to  export the results of the Matlab fspecial command and use that 
> kernel in Python code (save it as a .mat file and read that .mat file 
> with Python).
>
> SciPy does not have the equivalent to the fspecial command but you can 
> generate all kinds of 1-d special filters with 
> scipy.signal.get_window.   You can also "generate" the filters you 
> need directly from code.
>
> Here is some untested code for generating something close to what 
> fspecial('gaussian", [10,10], 2.5) would be producing
>
> import numpy as np
>
> def fgaussian(size, sigma):
>    m,n = size
>      h, k = m//2, n//2
>    x, y = np.mgrid[-h:h, -k:k]
>    return np.exp(-(x**2 + y**2)/(2*sigma**2))
>
> Up to a scaling constant (and possibly shifted a bit) this should 
> return a similar filter to fspecial('gaussian', size, sigma).
>
> I'm not completely sure if the ndimage.gaussian_filter does a 
> finite-window convolution like imfilter (and correlate) or if it does 
> something else.
>
> The uniform filter is the same as an averaging filter (up to a scaling 
> constant).  To other approach is to just use scipy.ndimage.correlate 
> with np.ones(size) / np.product(size)  where size is the size of the 
> kernel.
>
> Perhaps you would be willing to post your code when you get it to work.
>
> Best,
>
> -Travis
>
>
>
>
> On Wed, May 15, 2013 at 9:19 PM, Brickle Macho <bricklemacho at gmail.com 
> <mailto:bricklemacho at gmail.com>> wrote:
>
>     I am porting some Matlab code to python.  When I run the ported
>     version the output differs.   It appears the difference may be due
>     to how I have translated Matlab's imfilter() and/or fspecial().
>     Below are my translation, are there any Scipy/Matlab gurus that
>     can let me know if there are correct or if I have made some
>     errors, or could suggest better translations.
>
>     Matlab:  imfilter(saliencyMap, fspecial('gaussian', [10, 10], 2.5))
>     Scipy:  ndimage.gaussian_filter(saliencyMap, sigma=2.5)
>
>     Also the following.  At the rsik of highlighting my lack of
>     understanding of the temrinology, Is uniform_filter the same as
>     'average'?
>
>     Matlab:  imfilter(myLogAmplitude, fspecial('average', 3),
>     'replicate');
>     Scipy: scipy.ndimage.uniform_filter(mylogAmplitude, size=3,
>     mode="nearest")
>
>
>     Thanks,
>
>     Michael.
>     --
>
>
>     If curious here are the 8 lines of matlab code I am trying to port:
>
>     %% Read image from file
>     inImg = im2double(rgb2gray(imread('yourImage.jpg')));
>     inImg = imresize(inImg, 64/size(inImg, 2));
>
>     %% Spectral Residual
>     myFFT = fft2(inImg);
>     myLogAmplitude = log(abs(myFFT));
>     myPhase = angle(myFFT);
>     mySpectralResidual = myLogAmplitude - imfilter(myLogAmplitude,
>     fspecial('average', 3), 'replicate');
>     saliencyMap = abs(ifft2(exp(mySpectralResidual + i*myPhase))).^2;
>
>     %% After Effect
>     saliencyMap = mat2gray(imfilter(saliencyMap, fspecial('gaussian',
>     [10, 10], 2.5)));
>     imshow(saliencyMap);
>
>
>
>
>
>
>     _______________________________________________
>     SciPy-User mailing list
>     SciPy-User at scipy.org <mailto:SciPy-User at scipy.org>
>     http://mail.scipy.org/mailman/listinfo/scipy-user
>
>
>
>
> -- 
> ---
> Travis Oliphant
> Continuum Analytics, Inc.
> http://www.continuum.io

-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mail.scipy.org/pipermail/scipy-user/attachments/20130517/2ffb4b6e/attachment.html>