1

I have a function called funcPower3 with body presented below. I would like to plot this function with 3D plot functionnalities in MatplotLib. I saw an example on scipy docs with meshgrid. However, in this example, function is not a defined function, but a simple mathematical operation:

fig = plt.figure()
ax = fig.gca(projection='3d')
X = np.arange(-5, 5, 0.25)
Y = np.arange(-5, 5, 0.25)
X, Y = np.meshgrid(X, Y)
R = np.sqrt(X**2 + Y**2)
Z = np.sin(R)
surf = ax.plot_surface(X, Y, Z, rstride=1, cstride=1, cmap=cm.jet,
    linewidth=0, antialiased=False)

I take inspiration in this sample, but my own function signature is not compatible with data out of meshgrid. I have this error:

 ValueError: zero-size array to minimum.reduce without identity

My function's code is here:

def funcPower3(PARAM):
    inputX1 = open("power3X.txt","r")
    inputY = open("power3Y.txt","r")
    X1=[]
    Y=[]
    for line in inputX1:
        X1.append(float(line))
    for line in inputY:
        Y.append(float(line))
    resTmp_ = 0
    res = 0
    for i in range(len(X1)):
        resTmp_ = Y[i] - (PARAM[0]*(X1[i])**float(PARAM[1]))
        res += resTmp_**2
    return res

And code to plot 3d this function is here:

xmin = 0 ymin = 0 xmax = 10 ymax = 1

fig = plt.figure('Power 3')
ax = fig.gca(projection='3d')
Z = []
Xpl = numpy.arange(xmin, xmax, 0.1).tolist()
Ypl = numpy.arange(ymin, ymax, 0.01).tolist()
Xpl, Ypl = numpy.meshgrid(Xpl, Ypl)
Z=[]
for i in range(len(Xpl[0])):
   for j in range(len(Xpl)):
        Z.append(funcPower3([Xpl[j][i],Ypl[j][i]]))
surf=ax.plot_surface(Xpl, Ypl, Z, rstride=8, cstride=8, alpha=0.3,cmap=cm.jet)
ax.set_xlabel('X Label')
ax.set_ylabel('Y Label')
plt.show()

Thanks and regards for any advice; ))

kiriloff
  • 25,609
  • 37
  • 148
  • 229

1 Answers1

6

Quoting from plot_surface documentation:

X, Y, Z : Data values as 2D arrays

But your Z is 1 dimensional. You need to reshape it to match with the X and Y values. This should work:

Xpl = numpy.arange(xmin, xmax, 0.1).tolist()
Ypl = numpy.arange(ymin, ymax, 0.01).tolist()
Xpl, Ypl = numpy.meshgrid(Xpl, Ypl)

Z=[]
for j in range(len(Xpl)):
   for i in range(len(Xpl[0])):
        # your loop order was backwards
        Z.append(funcPower3([Xpl[j][i],Ypl[j][i]]))

# reshape Z
Z = numpy.array(Z).reshape(Xpl.shape)

fig = plt.figure('Power 3')
ax = fig.gca(projection='3d')
surf=ax.plot_surface(Xpl, Ypl, Z, rstride=8, cstride=8, alpha=0.3,cmap=cm.jet)
ax.set_xlabel('X Label')
ax.set_ylabel('Y Label')
plt.show()

But your code has a couple of flaws. First of all, looking at your funcPower3 function, you are reading two files over and over for each function call. This is way too much wasteful. Instead read those parameters once and give them to your function as parameters. This function could also be simplified a bit (and try to follow the PEP8 naming conventions):

def func_power_3(param, x1, y):
    p1, p2 = param
    res = sum(y_i - (p1*x1_i)**p2 for x1_i, y_i in zip(x1, y))
    return res

and the rest will be

with open("power3X.txt","r") as infile:
    x1 = [float(line) for line in infile]

with open("power3Y.txt","r") as infile:
    y = [float(line) for line in infile]

xpl = numpy.arange(xmin, xmax, 0.1)   # no need for .tolist()
ypl = numpy.arange(ymin, ymax, 0.01)  # meshgrid can work with numpy.array's
xpl, ypl = numpy.meshgrid(xpl, ypl)

# we can form z with list comprehension
z = [[func_power_3([p1,p2], x1, y) for p1, p2 in zip(p1row, p2row)] 
                                   for p1row, p2row in zip(xpl, ypl)]

fig = plt.figure('Power 3')
ax = fig.gca(projection='3d')
surf=ax.plot_surface(xpl, ypl, z, rstride=8, cstride=8, alpha=0.3,cmap=cm.jet)
ax.set_xlabel('X Label')
ax.set_ylabel('Y Label')
plt.show()
Avaris
  • 35,883
  • 7
  • 81
  • 72
  • This looks very promising and code is ultra-clear. I have a little issue with this error: "with open(open("power3X.txt","r")) as infile: TypeError: coercing to Unicode: need string or buffer, file found" – kiriloff Mar 27 '12 at 20:37
  • @dlib: oops too much `open` :). fixed now. – Avaris Mar 27 '12 at 20:38
  • Correcting to with open("power3Y.txt","r") as infile: y = [float(line) for line in infile], works great! – kiriloff Mar 27 '12 at 20:47