[Matplotlib-users] Attempting to get animated plot to work in this example
Samuel Dupree
sdupree at computational-astronomer.com
Fri Apr 22 14:06:20 EDT 2016
I can't reduce the example to a point that demonstrates my problem
without the "orbital" module. However, I've done some deeper digging
into the "orbital" module to find the functions that support the
plotting. I've attached that source file to this note, and if you could
point out any issues you may see dealing with the animate function, I
would appreciate it.
Sam Dupree.
On 4/22/16 13:48:45, Paul Hobson wrote:
> Can you reduce this problem down to an example that doesn't involve
> this "orbital" module? For all we know, the issue is in there
>
> On Fri, Apr 22, 2016 at 10:41 AM, Samuel Dupree
> <sdupree at computational-astronomer.com
> <mailto:sdupree at computational-astronomer.com>> wrote:
>
> There are two plots that should output plots showing a satellite
> orbiting about the Earth. Instead the plots are static which they
> sjouldn't be.
>
> Hope this helps.
>
> Sam Dupree.
>
>
>
>
>
> On 4/22/16 13:37:17, Paul Hobson wrote:
>> On Fri, Apr 22, 2016 at 10:07 AM, Samuel Dupree
>> <sdupree at computational-astronomer.com
>> <mailto:sdupree at computational-astronomer.com>> wrote:
>>
>> I'm running matplotlib ver. 1.5.1 in an anaconda distribution
>> on a MacBook Pro running Mac OS X ver. 10.11.4 (El Capitan).
>> I'm attempting to get the animation portion of the example
>> attached to this note running. I've looked through the
>> on-line documentation on plot() and plot3d() but have turned
>> up empty.
>>
>> Any thoughts.
>>
>> Sam Dupree.
>>
>>
>> What specific problems are you encountering? Does the script halt
>> with some error messages? If so, what are they?
>>
>> Does the script run, but fails to produce output?
>>
>> Does the script produce output, but it is not as expected?
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mail.python.org/pipermail/matplotlib-users/attachments/20160422/7cddaaed/attachment-0001.html>
-------------- next part --------------
"""plotting module for orbital
This implementation was inspired by poliastro (c) 2012 Juan Luis Cano (BSD License)
"""
# encoding: utf-8
from __future__ import absolute_import, division, print_function
from copy import copy
import matplotlib as mpl
import matplotlib.pyplot as plt
import numpy as np
from matplotlib import animation
from matplotlib.patches import Circle
from mpl_toolkits.mplot3d.axes3d import Axes3D
from numpy import cos, sin
from scipy.constants import kilo, pi
from .maneuver import TimeOperation
from .utilities import (
lookahead, orbit_radius, saved_state, uvw_from_elements)
__all__ = [
'plot2d',
'plot3d',
'plot',
'Plotter2D',
'Plotter3D',
]
def plot2d(orbit, title='', maneuver=None, animate=False, speedup=5000):
"""Convenience function to 2D plot orbit in a new figure."""
plotter = Plotter2D()
if animate:
return plotter.animate(orbit, title=title, speedup=speedup)
else:
plotter.plot(orbit, title=title, maneuver=maneuver)
def plot3d(orbit, title='', maneuver=None, animate=False, speedup=5000):
"""Convenience function to 3D plot orbit in a new figure."""
plotter = Plotter3D()
if animate:
return plotter.animate(orbit, title=title, speedup=speedup)
else:
plotter.plot(orbit, title=title, maneuver=maneuver)
plot = plot2d
class Plotter2D():
"""2D Plotter
Handles still and animated plots of an orbit.
"""
def __init__(self, axes=None, num_points=100):
if axes:
self.fig = axes.get_figure()
else:
self.fig = plt.figure()
axes = self.fig.add_subplot(111)
self.axes = axes
self.axes.set_aspect(1)
self.axes.set_xlabel("$p$ [km]")
self.axes.set_ylabel("$q$ [km]")
self.points_per_rad = num_points / (2 * pi)
def plot(self, orbit, maneuver=None, title=''):
self._plot_body(orbit)
if maneuver is None:
self._plot_orbit(orbit)
self.pos_dot = self._plot_position(orbit)
else:
self._plot_orbit(orbit, label='Initial orbit')
self.propagate_counter = 1
states = lookahead(
orbit.apply_maneuver(maneuver, iter=True, copy=True),
fillvalue=(None, None))
with saved_state(orbit):
for (orbit, operation), (_, next_operation) in states:
with saved_state(orbit):
operation.plot(orbit, self, next_operation)
self.axes.legend()
self.axes.set_title(title)
def animate(self, orbit, speedup=5000, title=''):
# Copy orbit, because it will be modified in the animation callback.
orbit = copy(orbit)
self.plot(orbit)
p = orbit.a * (1 - orbit.e ** 2)
def fpos(f):
pos = np.array([cos(f), sin(f), 0 * f]) * p / (1 + orbit.e * cos(f))
pos /= kilo
return pos
time_per_orbit = orbit.T / speedup
interval = 1000 / 30
times = np.linspace(orbit.t, orbit.t + orbit.T, time_per_orbit * 30)
def animate(i):
orbit.t = times[i - 1]
pos = fpos(orbit.f)
self.pos_dot.set_data(pos[0], pos[1])
return self.pos_dot
self.axes.set_title(title)
# blit=True causes an error on OS X, disable for now.
ani = animation.FuncAnimation(
self.fig, animate, len(times), interval=interval, blit=False)
return ani
@staticmethod
def _perifocal_coords(orbit, f):
p = orbit.a * (1 - orbit.e ** 2)
pos = np.array([cos(f), sin(f), 0 * f]) * p / (1 + orbit.e * cos(f))
pos /= kilo
return pos
def _plot_orbit(self, orbit, f1=0, f2=2 * pi, label=None):
if f2 < f1:
f2 += 2 * pi
num_points = self.points_per_rad * (f2 - f1)
f = np.linspace(f1, f2, num_points)
pos = self._perifocal_coords(orbit, f)
self.axes.plot(pos[0, :], pos[1, :], '--', linewidth=1, label=label)
def _plot_position(self, orbit, f=None, propagated=False, label=None):
if f is None:
f = orbit.f
pos = self._perifocal_coords(orbit, f)
if propagated:
if label is not None:
raise TypeError('propagated flag sets label automatically')
label = 'Propagated position {}'.format(self.propagate_counter)
self.propagate_counter += 1
pos_dot, = self.axes.plot(
pos[0], pos[1], 'o', label=label)
return pos_dot
def _plot_body(self, orbit):
color = '#EBEBEB'
if orbit.body.plot_color is not None:
color = orbit.body.plot_color
self.axes.add_patch(Circle((0, 0), orbit.body.mean_radius / kilo,
linewidth=0, color=color))
class Plotter3D(object):
"""3D Plotter
Handles still and animated plots of an orbit.
"""
def __init__(self, axes=None, num_points=100):
if axes:
self.fig = axes.get_figure()
else:
self.fig = plt.figure()
axes = self.fig.add_subplot(111, projection='3d')
self.axes = axes
self.axes.set_xlabel("$x$ [km]")
self.axes.set_ylabel("$y$ [km]")
self.axes.set_zlabel("$z$ [km]")
# These are used to fix aspect ratio of final plot.
# See Plotter3D._force_aspect()
self._coords_x = np.array(0)
self._coords_y = np.array(0)
self._coords_z = np.array(0)
self.points_per_rad = num_points / (2 * pi)
def plot(self, orbit, maneuver=None, title=''):
self._plot_body(orbit)
if maneuver is None:
self._plot_orbit(orbit)
self.pos_dot = self._plot_position(orbit)
else:
self._plot_orbit(orbit, label='Initial orbit')
self.propagate_counter = 1
states = lookahead(
orbit.apply_maneuver(maneuver, iter=True, copy=True),
fillvalue=(None, None))
with saved_state(orbit):
for (orbit, operation), (_, next_operation) in states:
with saved_state(orbit):
operation.plot(orbit, self, next_operation)
self.axes.legend()
self.axes.set_title(title)
self._force_aspect()
def animate(self, orbit, speedup=5000, title=''):
# Copy orbit, because it will be modified in the animation callback.
orbit = copy(orbit)
self.plot(orbit)
num_points = self.points_per_rad * 2 * pi
f = np.linspace(0, 2 * pi, num_points)
def fpos(f):
U, _, _ = uvw_from_elements(orbit.i, orbit.raan, orbit.arg_pe, f)
pos = orbit_radius(orbit.a, orbit.e, f) * U
pos /= kilo
return pos[0], pos[1], pos[2]
time_per_orbit = orbit.T / speedup
interval = 1000 / 30
times = np.linspace(orbit.t, orbit.t + orbit.T, time_per_orbit * 30)
def animate(i):
orbit.t = times[i - 1]
x, y, z = fpos(orbit.f)
self.pos_dot.set_data([x], [y])
self.pos_dot.set_3d_properties([z])
return self.pos_dot
self.axes.set_title(title)
# blit=True causes an error on OS X, disable for now.
ani = animation.FuncAnimation(
self.fig, animate, len(times), interval=interval, blit=False)
return ani
@staticmethod
def _xyz_coords(orbit, f):
U, _, _ = uvw_from_elements(orbit.i, orbit.raan, orbit.arg_pe, f)
pos = orbit_radius(orbit.a, orbit.e, f) * U
pos /= kilo
return pos
def _plot_orbit(self, orbit, f1=0, f2=2 * pi, label=None):
if f2 < f1:
f2 += 2 * pi
num_points = self.points_per_rad * (f2 - f1)
f = np.linspace(f1, f2, num_points)
pos = self._xyz_coords(orbit, f)
x, y, z = pos[0, :], pos[1, :], pos[2, :]
self.axes.plot(x, y, z, '--', linewidth=1, label=label)
self._append_coords_for_aspect(x, y, z)
def _plot_position(self, orbit, f=None, propagated=False, label=None):
if f is None:
f = orbit.f
pos = self._xyz_coords(orbit, f)
x, y, z = pos[0], pos[1], pos[2]
if propagated:
if label is not None:
raise TypeError('propagated flag sets label automatically')
label = 'Propagated position {}'.format(self.propagate_counter)
self.propagate_counter += 1
pos_dot, = self.axes.plot(
[x], [y], [z], 'o', label=label)
return pos_dot
def _plot_body(self, orbit):
color = '#EBEBEB'
if orbit.body.plot_color is not None:
color = orbit.body.plot_color
u = np.linspace(0, 2 * np.pi, 100)
v = np.linspace(0, np.pi, 50)
cx = orbit.body.mean_radius * np.outer(np.cos(u), np.sin(v))
cy = orbit.body.mean_radius * np.outer(np.sin(u), np.sin(v))
cz = orbit.body.mean_radius * np.outer(np.ones(np.size(u)), np.cos(v))
cx, cy, cz = cx / kilo, cy / kilo, cz / kilo
self.axes.plot_surface(cx, cy, cz, rstride=5, cstride=5, color=color,
edgecolors='#ADADAD', shade=False)
def _append_coords_for_aspect(self, x, y, z):
self._coords_x = np.append(self._coords_x, x)
self._coords_y = np.append(self._coords_y, y)
self._coords_z = np.append(self._coords_z, z)
def _force_aspect(self):
# Thanks to the following SO answer, we can make sure axes are equal
# http://stackoverflow.com/a/13701747/2093785
# Create cubic bounding box to simulate equal aspect ratio
x = self._coords_x
y = self._coords_y
z = self._coords_z
max_range = np.array([x.max() - x.min(),
y.max() - y.min(),
z.max() - z.min()]).max()
Xb = (0.5 * max_range * np.mgrid[-1:2:2, -1:2:2, -1:2:2][0].flatten() +
0.5 * (x.max() + x.min()))
Yb = (0.5 * max_range * np.mgrid[-1:2:2, -1:2:2, -1:2:2][1].flatten() +
0.5 * (y.max() + y.min()))
Zb = (0.5 * max_range * np.mgrid[-1:2:2, -1:2:2, -1:2:2][2].flatten() +
0.5 * (z.max() + z.min()))
for xb, yb, zb in zip(Xb, Yb, Zb):
self.axes.plot([xb], [yb], [zb], 'w')
More information about the Matplotlib-users
mailing list