[Matrix-SIG] Wiener filter using sigtools and Python

[Travis Oliphant]

While this routine will go into the next addition of sigtools, I thought I
would show interested people, how easy it is to do image (volume) 
processing with Python.

Below find a function to do wiener filtering on an array of data.

from sigtools import *
from MLab import *

def wiener(im,mysize=None,noise=None):

    if mysize==None:
        mysize = 3*ones(len(a.shape))
    mysize = array(mysize);

    # Estimate the local mean
    lMean = convolveND(im,ones(mysize),1) / prod(mysize)

    # Estimate the local variance
    lVar = convolveND(im**2,ones(mysize),1) / prod(mysize) - lMean**2

    # Estimate the noise power if needed.
    if noise==None:
        noise = mean(ravel(lVar))

    # Compute result
    # out = lMean + (maximum(0, lVar - noise) ./
    #               maximum(lVar, noise)) * (im - lMean) 
    #
    out = im - lMean
    im = lVar - noise
    im = maximum(im,0)
    lVar = maximum(lVar,noise)
    out = out / lVar
    out = out * im
    out = out + lMean

    return out




Enjoy,

--Travis