radialDensity.py 1.7 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748
  1. """Plot a profile of the radial density"""
  2. import matplotlib.pyplot as plt
  3. import numpy as np
  4. from analysis import utils
  5. def plotRadialDensity(data):
  6. """Plots a profile of the radial density, divided into components.
  7. Parameters:
  8. data: a data dictionary in the format saved by the simulation
  9. """
  10. f, ax = plt.subplots(1, 1)
  11. # We measure the radial distance from the Center Of Mass
  12. COM = np.sum(data['r_vec'] * data['mass'][:,np.newaxis], axis=0) / np.sum(data['mass'])
  13. r = np.linalg.norm((data['r'] - COM), axis=-1)
  14. # Divide the particles into groups
  15. diskMask = data['types'][:,1] == 'disk'
  16. bulgeMask = data['types'][:,1] == 'bulge'
  17. luminousMask = np.logical_or(diskMask, bulgeMask)
  18. darkMask = data['types'][:,1] == 'dark'
  19. # Plot each group separately
  20. style = {'histtype':'step', 'cumulative':True, 'bins':np.logspace(-2, 2, 1000),
  21. 'color':'black', 'linestyle':'dashed', 'label':'disk'}
  22. plt.hist(r[diskMask], weights=data['mass'][diskMask],
  23. label='disk', **style)
  24. plt.hist(r[bulgeMask], weights=data['mass'][bulgeMask],
  25. label='bulge', **style)
  26. plt.hist(r[luminousMask], weights=data['mass'][luminousMask],
  27. label='luminous', **style)
  28. plt.hist(r[darkMask], weights=data['mass'][darkMask],
  29. label='dark', **style)
  30. # Plotting styling
  31. plt.legend(loc='upper left')
  32. plt.yscale('log'); plt.xscale('log')
  33. utils.setAxes(ax, x=(3E-2, 1E2), y=(1E-2, 2))
  34. utils.setAxes(ax, x=r'$r$', y=r'Cumulative mass',
  35. xcoords=(.9,-0.08), ycoords=(-0.1,.7))
  36. utils.stylizePlot([ax])
  37. plt.show()
  38. data = utils.loadData('hyperbolic_halo', 30000)
  39. plotRadialDensity(data)