2

I'm trying to plot elevation cross section and extract elevations along a cross section using scipy's "interpolate" function. I'm converting a point feature class of profile elevation points to 2 numpy arrays representing "distance" and "elevation", respectively.

Here is my code:

import numpy
import arcpy
from scipy import interpolate
import matplotlib.pyplot as plt

arcpy.env.workspace = "D:\\NJBPN\\volume_calc_scripts\\test_njbpn.gdb"
#interpolate profile
x = arcpy.da.FeatureClassToNumPyArray("test_prof",["Distance"])
y = arcpy.da.FeatureClassToNumPyArray("test_prof",["Z"])
f = interpolate.interp1d(x,y,kind='linear')

xnew = numpy.arange(numpy.amin(x, axis = 0),numpy.amax(y, axis = 1))
ynew = f(xnew)
interpolate.interp1d._call_(xnew)
plt.plot(x,y,'o',xnew,f(xnew),'-')
plt.legend(['data', 'linear'], loc='best')
plt.show()

It keeps giving the following error:

Runtime error 
Traceback (most recent call last):
  File "<string>", line 12, in <module>
  File "C:\Python27\ArcGIS10.2\lib\site-packages\numpy\core\fromnumeric.py", line 1987, in amin
    out=out, keepdims=keepdims)
  File "C:\Python27\ArcGIS10.2\lib\site-packages\numpy\core\_methods.py", line 14, in _amin
    out=out, keepdims=keepdims)
TypeError: cannot perform reduce with flexible type

Here is an updated code to extract the xnew value at ynew=0

xnew = numpy.arange(numpy.amin(x),numpy.amax(x)).astype(numpy.int)
ynew = f(xnew).astype(numpy.float)
shorelinedist = numpy.where(ynew==0.00,[xnew]

plt.plot(x,y,'o',xnew,ynew,'-')
plt.legend(['data', 'linear'], loc='best')
plt.show()
  • I don't understand your xnew line did you mean >>> x = np.arange(10) >>> y = np.arange(2,12) >>> xnew = np.arange(np.amin(x, axis = 0),np.amax(y, axis = 0)) >>> xnew array([ 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10]) >>> – user681 May 13 '15 at 15:30
  • I've figured out the error - I needed to specify the values at "float". Sorry I'm new to this. – thecurrentwave May 13 '15 at 17:21
  • The next piece of code I'm having trouble with. I can plot the interpolated values, but now I need to find the xnew value where ynew = 0 thanks for your help!!! – thecurrentwave May 13 '15 at 17:24
  • If your first issue has been resolved through your own efforts then you could either edit your question to revise it so that it is only about the current issue (overhauling a question is fine to do while no answers have been offered) or add an answer to your own first question (self-answering and later self-accepting is fine) and then ask your second question as a new one. – PolyGeo May 14 '15 at 7:28