summaryrefslogtreecommitdiffstats
path: root/plot_stats.py
diff options
context:
space:
mode:
Diffstat (limited to 'plot_stats.py')
-rw-r--r--plot_stats.py36
1 files changed, 32 insertions, 4 deletions
diff --git a/plot_stats.py b/plot_stats.py
index 1802d8d..45eddcd 100644
--- a/plot_stats.py
+++ b/plot_stats.py
@@ -62,33 +62,61 @@ for iout in iouts:
mass_total[ii] = disk.mass_halo
ii += 1
+# Calculate Time Derivatives
+# @todo - What about negative derivatives?
+# @todo - Get them in the plot somehow; color code \pm?
+delta_mass_total_lost_to_star = np.diff(mass_total_lost_to_star) / \
+ np.diff(time)
+delta_mass_total_lost_to_halo = np.diff(mass_total_lost_to_halo) / \
+ np.diff(time)
+delta_mass_disk = np.diff(mass_disk) / np.diff(time)
+delta_mass_halo = np.diff(mass_halo) / np.diff(time)
+delta_mass_total = np.diff(mass_total) / np.diff(time)
+
# Make Plots
-fig = plt.figure()
-ax = fig.add_subplot(1,1,1)
+fig = plt.figure(figsize=(16.0, 6.0))
+ax = fig.add_subplot(1,2,1)
ax2 = ax.twiny()
+ax3 = fig.add_subplot(1,2,2)
+ax4 = ax3.twiny()
if args.disk:
ax.plot(time/2./pi, mass_disk, 'bs-', label='Disk Mass')
+ ax3.plot(time[1::]/2./pi, delta_mass_disk, 'bs-', label='d(Disk Mass)/dt')
if args.halo:
ax.plot(time/2./pi, mass_halo, 'rs-', label='Halo Mass')
+ ax3.plot(time[1::]/2./pi, delta_mass_halo, 'rs-', label='d(Halo Mass)/dt')
if args.box:
ax.plot(time/2./pi, mass_total, 'gs-', label='Box Mass')
+ ax3.plot(time[1::]/2./pi, delta_mass_total, 'gs-', label='d(Box Mass)/dt')
if args.lost_halo:
ax.plot(time/2./pi, mass_total_lost_to_halo, 'ms-', \
label='Total Mass Lost to Halo')
+ ax3.plot(time[1::]/2./pi, delta_mass_total_lost_to_halo, 'ms-', \
+ label='d(Total Mass Lost to Halo)/dt')
if args.lost_star:
ax.plot(time/2./pi, mass_total_lost_to_star, 'ks-', \
label='Total Mass Lost to Star')
+ ax3.plot(time[1::]/2./pi, 2.*pi*delta_mass_total_lost_to_star, 'ks-', \
+ label='d(Total Mass Lost to Star)/dt')
fig.suptitle('Mass Budget')
ax2.plot(time, np.ones_like(time))
ax2.cla()
ax2.set_xlabel('Time [code]')
-ax.grid()
-ax.set_yscale('log')
+ax4.plot(time, np.ones_like(time))
+ax4.cla()
+ax4.set_xlabel('Time [code]')
+ax.grid(); ax.set_yscale('log')
+ax3.grid(); ax3.set_yscale('log')
ax2.set_ylim(ax.get_ylim())
ax2.set_yscale(ax.get_yscale())
+ax4.set_ylim(ax3.get_ylim())
+ax4.set_yscale(ax3.get_yscale())
ax.set_xlabel('Time [yr]')
ax.set_ylabel('Mass [Mstar]')
ax.legend(loc='best')
+ax3.set_xlabel('Time [yr]')
+ax3.set_ylabel('d(Mass)/dt [Mstar/yr]')
+ax3.legend(loc='best')
if args.save:
plt.savefig('Stats.pdf')
if args.show: