Modelling an asteroid with Matplotlib using surface and wireframe

Question

I am trying to model an asteroid using plot_surface and plot_wireframe. I have x y and z values for the points on the surface of the asteroid. The wireframe is accurate to the shape of the asteroid but the surface plot does not fit the wireframe. How can I get the surface plot to fit the wireframe or how could i use the wireframe to get a 3d solid model? Here is my code for the model:

from mpl_toolkits.mplot3d import axes3d
import matplotlib.pyplot as plt
import numpy as np
from matplotlib import cm
from matplotlib.mlab import griddata

fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
data = np.genfromtxt('data.txt')
x = data[:,0]
y = data[:,1]
z = data[:,2]

ax.plot_wireframe(x, y, z, rstride=1, cstride=1, alpha=1)

xi = np.linspace(min(x), max(x))
yi = np.linspace(min(y), max(y))

X, Y = np.meshgrid(xi, yi)
Z = griddata(x, y, z, xi, yi)


surf = ax.plot_surface(X, Y, Z, rstride=1, cstride=1, cmap=cm.coolwarm,
        linewidth=0, antialiased=False)
ax.set_zlim(-1.01, 1.01)

plt.show()

The data is in this format, although there are a lot more lines in the original file:

-1.7738946051191869E-002   4.3461451610545973E-002   1.3393057231408241     
-0.29733561550902488       0.32305812106837900        1.3393057231408241     
-0.29733561550902488       0.16510132228266330        1.3548631099230350     
-0.21872587865015569        2.4170900455101410E-002   1.3610011616437809     
1.4452975249810950E-002 -0.20900795344486520        1.3610011616437809    
1.5732454381265970E-002 -0.20900795344486520        1.3608751439485580    
-0.34501536374240321       0.51320241386595655        1.3158820995876130     
-0.40193014435941982       0.45628763324893978        1.3158820995876130     
-0.42505849480150409       0.28183419537116011        1.3307863198123011     
-0.18994178462386799      -0.19294290416565860        1.3424523041534830     
1.4452975249810939E-002 -0.39733766403933751        1.3424523041534830     
5.8021940902131752E-002 -0.57108837516584876        1.3210481842104100     
9.3746267961881152E-002 -0.61017602710257668        1.3136798474111200     
0.26609469681891229      -0.43782759824554562        1.3136798474111200     
  0.17938460413447810       0.39179924148155021        1.2357401964919650     
   8.9613011902522258E-002  0.42818009222325598        1.2584008460875080     
  0.33671539027096409      -0.47165177581327772        1.2965073126705291     
  0.53703772594296528      -0.47165177581327777        1.2357401964919561     
 -0.19242375014122229       0.71021685426700043        1.2584008460875080     
 -0.34501536374240321       0.66763766324752027        1.2904902860951690 

Hope you can help

Bedmate

I would recomend looking into code.enthought.com/projects/mayavi which leverages vtk and provides much better 3D support than matplotlib. – Almagest

Interesting problem, can you post a link to the complete data? – Cartridge

I agree with @tcaswell: matplotlib's 3D isn't really up to this; see #6030598 for how mayavi was used to render isosurfaces of 3D volume data (not quite the same problem). – Tiannatiara

Answer 1

Could you provide an image of what is happening? I'm guessing that you are getting an image where the surface of the asteroid seems to be jumping all over the place. Is that correct? If so it could be caused by the plotter not knowing the order of the points.

If we have a set of points, which contains all the points of a unit circle, then we would expect a drawing of those points to create a unit circle. If, however, you decided to connect each point to two other points, then it wouldn't necessarily look like a circle. If (for some reason) you connected one point to another point on the other side of the circle and continued to do that until each point was connected to two other points, it may or may not look like a circle because each point isn't necessarily connected to the adjacent points.

The same is true for the your asteroid. You need to come up with some scheme so that the plotter knows how to connect the points, otherwise you will continue to have the same problem.

The following example of a circle should illustrate my point:

import math
import matplotlib.pylab as plt
import random

thetaList = range(360)
random.shuffle(thetaList)
degToRad = lambda x: float(x) * math.pi / float(180)
x = [math.cos(degToRad(theta)) for theta in thetaList]
y = [math.sin(degToRad(theta)) for theta in thetaList]       

#plot the cirlce
plt.plot(x,y)
plt.show()

Precambrian