[SciPy-user] Re: data smoothing: interpolate.splrep ignores
barrett at stsci.edu
Thu Feb 26 09:15:58 CST 2004
george young wrote:
> [[reposting with "rough" data file appended, sorry]]
> [SciPy-0.2.0_alpha_200.4161, Numeric-23.1, Python 2.3.3, x86 Linux]
> My goal really is to smooth some noisy measurement data without messing
> up it's *shape*. My first attempt was 1d splines. I did:
> from Numeric import *
> import scipy
> f = file('rough')
> data = array([(float(x),float(y)) for (x,y) in [l.split() for l in f]])
> rep = scipy.interpolate.splrep(data[:,0], data[:,1], s=s)
> smooth_y = scipy.interpolate.splev(data[:,0], rep)
> "s" is supposed to vary the amount of smoothing. For s=0, i get
> my original data back, as expected. But for all other values, including
> the recommended range of (m-sqrt(2.0*m)) to (m+sqrt(2.0*m), I get a
> single, too much smoothed, result. It seems like splrep is not using
> the "s" value to adjust the splines, except to sense that it's not zero.
> Splines may not be the right method anyway, since they tend to warp the
> shape of the curve, and I need to get the data's derivatives. Is there
> some way to fashion a low pass filter? It seems like fft should be useful
> here, but I have very little experience with fft's.
Have you considered wavelet smoothing or denoising?
The 'a-trous' wavelet method is a simple, understandable technique that can be
done by convolving a function with the data. The difference between your raw
data and your smoothed (convolved) data is the wavelet at that particular
spatial resolution. You can filter the wavelet values using thresholding and
then add them back into the smoothed data to get your smoothed or denoised
result. This can be done at different spatial resolutions depending on the what
scale you wish to smooth.
If you set the filtering threshold correctly, you'll be able to keep the
significant features of the data without loosing resolution as is common with
averaging techniques, or adding spurious features as is common with FFT
techniques. A B3-spline is good kernel function and usually does a good job in
this situation. Its 1-dimensional representation is [1, 4, 6, 4, 1]/16.
Paul Barrett, PhD Space Telescope Science Institute
Phone: 410-338-4475 ESS/Science Software Branch
FAX: 410-338-4767 Baltimore, MD 21218
More information about the SciPy-user