[Scipy-tickets] [SciPy] #1260: scipy.signal.cspline1d_eval cannot return complex results
SciPy Trac
scipy-tickets@scipy....
Wed Aug 11 09:59:00 CDT 2010
#1260: scipy.signal.cspline1d_eval cannot return complex results
--------------------------+-------------------------------------------------
Reporter: anielsen | Owner: somebody
Type: defect | Status: new
Priority: normal | Milestone: 0.9.0
Component: scipy.signal | Version: 0.7.0
Keywords: |
--------------------------+-------------------------------------------------
I apologize in advance if I misunderstand what's going on here. I'm new to
python and scipy.
I've found that cspline1d_eval cannot return complex results. The forward
transform (cspline1d) will create complex coefficients when given a
complex series input. I did some simple tests that confirm these changes
produce acceptable results for appropriately smooth complex input.
the code is in scipy/signal/bsplines.py
I made a simple hack which produces the behavior I'd anticipated by
replacing
res = zeros_like(newx)
with
res = zeros(newx.shape,dtype=cj.dtype)
in two locations, that I marked in the code snippet below.
I've attached a very simple script which demonstrates the results of using
these changes to interpolate a simple complex signal.
def cspline1d_eval(cj, newx, dx=1.0, x0=0):
"""Evaluate a spline at the new set of points.
dx is the old sample-spacing while x0 was the old origin.
In other-words the old-sample points (knot-points) for which the cj
represent spline coefficients were at equally-spaced points of
oldx = x0 + j*dx j=0...N-1
N=len(cj)
edges are handled using mirror-symmetric boundary conditions.
"""
newx = (asarray(newx)-x0)/float(dx)
### change made here
#res = zeros_like(newx)
res = zeros(newx.shape,dtype=cj.dtype)
if (res.size == 0):
return res
N = len(cj)
cond1 = newx < 0
cond2 = newx > (N-1)
cond3 = ~(cond1 | cond2)
# handle general mirror-symmetry
res[cond1] = cspline1d_eval(cj, -newx[cond1])
res[cond2] = cspline1d_eval(cj, 2*(N-1)-newx[cond2])
newx = newx[cond3]
if newx.size == 0:
return res
### change made here
#result = zeros_like(newx)
result = zeros(newx.shape,dtype=cj.dtype)
jlower = floor(newx-2).astype(int)+1
for i in range(4):
thisj = jlower + i
indj = thisj.clip(0,N-1) # handle edge cases
result += cj[indj] * cubic(newx - thisj)
res[cond3] = result
return res
--
Ticket URL: <http://projects.scipy.org/scipy/ticket/1260>
SciPy <http://www.scipy.org>
SciPy is open-source software for mathematics, science, and engineering.
More information about the Scipy-tickets
mailing list