[SciPy-user] What is missing in scipy to be a top notch environment for signal processing

Chris Bartley chris at edc.com.au
Mon Nov 20 20:02:38 CST 2006


> From: David Cournapeau <david at ar.media.kyoto-u.ac.jp>
>
>    I was wondering how many people here are using numpy/scipy for signal
processing, 
> and what are their impression compared to eg matlab ? 
> The point is not to do a matlab vs scipy, but more to spot weak points in
scipy, 
> and to change the situation; neither do I want to criticize anything. The
whole 
> point is really to improve scipy.
David,

I use scipy for signal analysis of engineering vibration measurements.
Although I understand it is not strictly speaking part of scipy, the single,
main, thing that stops the other guys here using it is the speed of plotting
large time series datasets (100k-400k data points). Matlab plotting is
incredibly quick! Matplotlib is great for many things, but seems to suffer
for large datasets, especially where the signal covers a lot of the plot
canvas (ie a thick, noisy line, not a thin line). 

The other feature we use a lot is the matlab 'decimate', which I understand
is not currently built into scipy, although I do have code for a port of
this function from Octave.

Otherwise, for our everyday purposes we use: fft, filtering, detrending,
specgram, matrix arithemetic. All freely available in scipy! :)

Chris Bartley

********************************************************
Mechanical Engineer

Engineering Dynamics Consultants
21 Riseley St Applecross WA 6153
+61 8 9316 3577

********************************************************


-----Original Message-----
From: scipy-user-bounces at scipy.org [mailto:scipy-user-bounces at scipy.org] On
Behalf Of scipy-user-request at scipy.org
Sent: Monday, 20 November 2006 2:00 AM
To: scipy-user at scipy.org
Subject: SciPy-user Digest, Vol 39, Issue 34

Send SciPy-user mailing list submissions to
	scipy-user at scipy.org

To subscribe or unsubscribe via the World Wide Web, visit
	http://projects.scipy.org/mailman/listinfo/scipy-user
or, via email, send a message with subject or body 'help' to
	scipy-user-request at scipy.org

You can reach the person managing the list at
	scipy-user-owner at scipy.org

When replying, please edit your Subject line so it is more specific than
"Re: Contents of SciPy-user digest..."


Today's Topics:

   1. What is missing in scipy to be a top notch environment for
      signal processing (lpc and co) ? (David Cournapeau)
   2. Re: Array documentation (Herman Berendsen)
   3. physical constants module (Herman Berendsen)
   4. Re: Unexpected behaviour of convolve2d and correlate2d
      (Gael Varoquaux)


----------------------------------------------------------------------



    I've just finished my second big conversion of matlab -> python code
(~5000 lines of matlab code), and I think there are some "holes" in scipy,
which would be really useful to fill in. I believe they are general enough
so I am not the only one missing them. Here are some functions I missed

    1: linear prediction coefficients computation (the matlab lpc function).
    2: more flexible autocorrelation method (ala xcorr in matlab).
    3: good resampling function in time domain
    4: functions capable of running the same algorithm on some strides of an
array.

More detailed:

    1 requires a method to invert a Toeplitz matrix (durbin levinson) and a
method for autocorrelation.
    2 right now, I believe there is only time domain autocorrelation
implementation, which is expensive for any size exceeding a few
tens/hundreds samples in numpy. Also, it is not possible to select the lag
required. In LPC coding for speech, we often need only a few lags of signals
around a few hundreds samples; just computing what is needed would already
give a 10 times fold speed increase at least for this kind of problems. For
problems where the size of the signal and the number of coefficients is the
same scale, a fft  based autocorrelation would also be beneficial.
    3 basically, we want an equivalent to upfirdn, which is using a
polyphase implementation according to matlab doc + good filter design
methods (which already exist in scipy AFAIK).
    4 I am not sure about this one. A concrete example: if I want to compute
the autocorrelation of some frames of a signal, the obvious thing is to use
a loop. This is expensive. If I had a matrix which each column is a frame,
and an autocorrelation function capable of running an algorithm on each
column, this would be much faster. Incidentally, Matlab offers a function
buffer, which builds a matrix which each column is a frame, with options for
overlapping and border cases. The doc says "Y = BUFFER(X,N) partitions
signal vector X into nonoverlapping data segments (frames) of length N.
Each data frame occupies one column in the output matrix, resulting in a
matrix with N rows." I don't know how, if at all possible, to generalize
that to the numpy array.
 
    Now, because I don't want to look like a whiner, I have some code to
solve partially or totally some of the points:
    - I have C code to compute Levinson Durbin, checked against matlab. 
It expects to have the autocorrelation as an contiguous array; as it is mono
dimensional, adapt it to multiple stride would be trivial
    - I have C code to compute only one side, and a few lags of
autocorrelation in time domain. This also expects contiguous array, and
would be a bit more tricky to adapt.
    - I have code based on fftw to compute the autocorrelation using fft
(can be adapted for cross correlation too). As I think this would be useful
in scipy, I understand that it cannot use fftw. Which fft should I use in
scipy C code ?
    - I would be interested in solving 3, and eventually 4, but I would need
some advices from others, as I am not sure how to solve them API-wise.

    I don't know if this is of any interest for other, but I believe some of
those functionalities to be a basic requirement for scipy to be used by
people in signal processing.

    Cheers,

    David
 


------------------------------

Message: 2
Date: Sun, 19 Nov 2006 10:36:43 +0000 (UTC)
From: Herman Berendsen <H.J.C.Berendsen at rug.nl>
Subject: Re: [SciPy-user] Array documentation
To: scipy-user at scipy.org
Message-ID: <loom.20061119T112156-538 at post.gmane.org>
Content-Type: text/plain; charset=us-ascii


Gael Varoquaux <gael.varoquaux <at> normalesup.org> writes:
> 
> Very nice and very complete work. Did you do this be hand or did your 
> generate some part from docstrings. In the later case, would you 
> consider sharing the scripts used ?

Thanks both respondents and those who sent me an email for the nice
comments.
Actually, the array documentation is not yet complete, as I do not have
added examples for all methods and functions. More will follow, so please
check-out in the future. 
I used the docstrings, but modified and extended by hand (programs never add
anything original!). I use the public domain web editor Arachnophilia which
keeps the layout of the html file simple and readable.
Please also check my posted message (Nov 19) on a module physcon.
Herman




------------------------------

Message: 3
Date: Sun, 19 Nov 2006 10:52:02 +0000 (UTC)
From: Herman Berendsen <H.J.C.Berendsen at rug.nl>
Subject: [SciPy-user] physical constants module
To: scipy-user at scipy.org
Message-ID: <loom.20061119T113916-601 at post.gmane.org>
Content-Type: text/plain; charset=us-ascii

If you use physical constants (like Planck's constant h, gas constant R,
etc.) with python, you may be interested to copy my module physcon.py from
http://www.hjcb.nl/python/ It contains a library of physical constants
(CODATA values in SI units, with accuracy) and it defines variables. For
example, after 'import physcon as pc', pc.h will contain Planck's constant
as a float. To see the contents, type help().



------------------------------

Message: 4
Date: Sun, 19 Nov 2006 16:38:10 +0100
From: Gael Varoquaux <gael.varoquaux at normalesup.org>
Subject: Re: [SciPy-user] Unexpected behaviour of convolve2d and
	correlate2d
To: SciPy Users List <scipy-user at scipy.org>
Message-ID: <20061119153810.GD6313 at clipper.ens.fr>
Content-Type: text/plain; charset=iso-8859-1

Hi,

I think you have uncovered a bug, this is taken from a ipython session (I
didn't cheat):

In [42]: signal.convolve2d(ones((4,4)), ones((1,15)), 'same', 'wrap')
Out[42]: 
array([[                           nan,                    13.        ,
                           14.        ,                    15.        ],
       [                   15.        ,                    15.        ,
                           15.        ,                    15.        ],
       [                   15.        ,                    15.        ,
                           15.        ,                    15.        ],
       [                   15.        ,     15892918649396384.        ,
            15892918649396384.        ,     15892918649396384.        ]])

In [43]: signal.convolve2d(ones((4,4)), ones((1,15)), 'same', 'wrap')
Out[43]: 
array([[                           nan,                    13.        ,
                           14.        ,                    15.        ],
       [                   15.        ,                    15.        ,
                           15.        ,                    15.        ],
       [                   15.        ,                    15.        ,
                           15.        ,                    15.        ],
       [                   15.        ,                    14.        ,
                           13.        ,                    12.        ]])

In [44]: signal.convolve2d(ones((4,4)), ones((1,15)), 'same', 'wrap')
Out[44]: 
array([[                           nan,                    13.        ,
                           14.        ,                    15.        ],
       [                   15.        ,                    15.        ,
                           15.        ,                    15.        ],
       [                   15.        ,                    15.        ,
                           15.        ,                    15.        ],
       [                   15.        ,                    14.        ,
                           13.        ,                    12.        ]])

In [45]: signal.convolve2d(ones((4,4)), ones((1,15)), 'same', 'wrap')
Out[45]: 
array([[                           nan,                    13.        ,
                           14.        ,                    15.        ],
       [                   15.        ,                    15.        ,
                           15.        ,                    15.        ],
       [                   15.        ,                    15.        ,
                           15.        ,                    15.        ],
       [                   15.        ,                    14.        ,
                           13.        ,                    12.        ]])

In [46]: signal.convolve2d(ones((4,4)), ones((1,15)), 'same', 'wrap')
Out[46]: 
array([[ 12.,  13.,  14.,  15.],
       [ 15.,  15.,  15.,  15.],
       [ 15.,  15.,  15.,  15.],
       [ 15.,  14.,  13.,  12.]])


Could someone more qualified than me look into this ?

Scipy version is 0.5.1 (ubuntu breezy), and numpy 1.0rc1

 Ga??l


------------------------------

_______________________________________________
SciPy-user mailing list
SciPy-user at scipy.org
http://projects.scipy.org/mailman/listinfo/scipy-user


End of SciPy-user Digest, Vol 39, Issue 34
******************************************



More information about the SciPy-user mailing list