Hi folks,
I'm trying to estimate smooth contours with scipy's parametric
splines, given a set of (x, y) points known to lie along the contour.
Now, because contours are closed, special care must be taken to
ensure that the interpolated value at the last point be the same as
the interpolated value at the first point.
I assume that the 'per' option to scipy.interpolate.splprep is
designed to allow for this, as per the documentation:
> per -- If non-zero, data points are considered periodic with period
> x[m-1] - x[0] and a smooth periodic spline
> approximation is returned.
> Values of y[m-1] and w[m-1] are not used.
However, I cannot get this to work on my computer. Setting 'per =
true' always results in memory errors or other problems.
Here is a simple example to reproduce the problem:
# first make x and y points along the unit circle, from zero to just
below two pi
In [1]: import numpy, scipy, scipy.interpolate
In [2]: twopi = numpy.arange(0, 2 * numpy.pi, 0.1)
In [3]: xs = numpy.cos(twopi)
In [4]: ys = numpy.sin(twopi)
In [5]: tck, uout = scipy.interpolate.splprep([xs, ys], u = twopi,
per = True)
Warning: Setting x[0][63]=x[0][0]
Warning: Setting x[1][63]=x[1][0]
...[here my machine grinds for 2-3 minutes]...
Warning: The required storage space exceeds the available
strorage space.
Probably causes: nest to small or s is too small. (fp>s)
At this point, the returned tck arrays are just all zeros. Sometimes
I get other malloc errors printed to stdout and memory error
exceptions, e.g.:
Python(6820,0xa000ed68) malloc: *** vm_allocate(size=3716243456)
failed (error code=3)
Python(6820,0xa000ed68) malloc: *** error: can't allocate region
Python(6820,0xa000ed68) malloc: *** set a breakpoint in szone_error
to debug
During the time that my machine is grinding, python is using very
little CPU -- the grind is all because python is allocating huge
amounts of memory, causing the pager to go nuts.
If I explicitly make the last value of the x and y input arrays equal
to the first value (as the warnings say that the function is doing),
I get the same problem:
In [6]: xs[-1] = xs[0]
In [7]: ys[-1] = ys[0]
In [8]: tck, uout = scipy.interpolate.splprep([xs, ys], u = twopi,
per = True)
#same thing
Any thoughts?
Zach Pincus
Program in Biomedical Informatics and Department of Biochemistry
Stanford University School of Medicine