pygedm.pygedm¶
Python API to YMW16, NE2001, and YT2020 Galactic electron density models
References
[1] Cordes, J. M., & Lazio, T. J. W. (2002), NE2001.I. A New Model for the Galactic Distribution of Free Electrons and its Fluctuations, arXiv e-prints, astro-ph/0207156.
[2] Cordes, J. M., & Lazio, T. J. W. (2003), NE2001. II. Using Radio Propagation Data to Construct a Model for the Galactic Distribution of Free Electrons, arXiv e-prints, astro-ph/0301598.
[3] Yao, J. M., Manchester, R. N., & Wang, N. (2017), A New Electron-density Model for Estimation of Pulsar and FRB Distances, ApJ, 835, 29.
[4] Yamasaki S, Totani T (2020), The Galactic Halo Contribution to the Dispersion Measure of Extragalactic Fast Radio Bursts The Astrophysical Journal, Volume 888, Issue 2, id.105
-
pygedm.pygedm.
calculate_electron_density_lbr
(gl, gb, dist, method='ymw16')¶ Calculate electron density at a point with Galactic coords (ga, gl) at given distance
Parameters: - gl (float, Angle, or Quantity) – Galatic longitude in degrees (or astropy Angle)
- gb (float, Angle, or Quantity) – Galactic latitude in degrees (or astropy Angle)
- dist (float or Quantity) – Distance in pc
Returns: electron density in cm^-3
Return type: N_e (astropy.Quantity)
-
pygedm.pygedm.
calculate_electron_density_xyz
(x, y, z, method='ymw16')¶ Calculate electron density at a point with galactocentric coords (x, y, z)
Parameters: - x (float or Quantity) – galactocentric X coordinate in pc
- y (float or Quantity) – galactocentric Y coordinate in pc
- z (float or Quantity) – galactocentric Z coordinate in pc
Returns: electron density in cm^-3
Return type: N_e (astropy.quantity)
-
pygedm.pygedm.
calculate_halo_dm
(gl, gb, method='yt2020', component='both')¶ Compute halo DM
Parameters: - gl (float, Angle, or Quantity) – Galatic latitude in degrees (or astropy Angle)
- gb (float, Angle, or Quantity) – Galactic latitude in degrees (or astropy Angle)
- method (str) – one of ‘yt2020’ (only YT2020 supported currently)
- component (str) – Compute ‘spherical’ component of halo, ‘disk’, or ‘both’ components.
Returns: Dispersion measure in (pc/cm3)
Return type: DM (float)
-
pygedm.pygedm.
convert_lbr_to_xyz
(gl, gb, dist, method='ymw16')¶ Convert Galactic (l,b,r) coords to Galactocentric (x,y,z) coords
Parameters: - gl (float, Angle, or Quantity) – Galatic longitude in degrees (or astropy Angle)
- gb (float, Angle, or Quantity) – Galactic latitude in degrees (or astropy Angle)
- dist (float or Quantity) – Distance in pc
- method (str) – one of ‘ymw16’, ‘ne2001’, or ‘astropy’
Returns: Galactocentric X, Y, Z coordinates
Return type: xyz (tuple)
Notes
For transform, the Sun is located at (x=0, y=R_sun, z=z_sun) YMW16 uses R_sun of 8300 pc and z_sun of 6.0 pc NE2001 uses R_sun of 8500 pc and z_sun of 0.0 pc Both of these do a basic spherical to cartesian conversion.
astropy does a much more complicated conversion, see https://astropy.readthedocs.io/en/latest/coordinates/galactocentric.html This is the ‘proper’ coordinate system, but note that it is NOT COMPATIBLE WITH NE2001 OR YMW16! (!SEE EXAMPLE OUTPUT BELOW!)
Example output:
pygedm.convert_lbr_to_xyz(0, 0, 0, method='ymw16') (<Quantity 0. pc>, <Quantity 8300. pc>, <Quantity 6. pc>) pygedm.convert_lbr_to_xyz(0, 0, 0, method='ne2001') (<Quantity 0. pc>, <Quantity 8500. pc>, <Quantity 0. pc>) pygedm.convert_lbr_to_xyz(0, 0, 0, method='astropy') (<Quantity -8499.95711754 pc>, <Quantity 0. pc>, <Quantity 27. pc>)
-
pygedm.pygedm.
dist_to_dm
(gl, gb, dist, mode='gal', method='ymw16', nu=1.0)¶ Convert a distance to a DM
Parameters: - gl (float in deg or astropy.Angle) – galactic longitude
- gb (float in deg or astropy.Angle) – galactic latitude
- dist (float or astropy.Quantity) – distance to source (pc) or if in mode IGM use (Mpc)
- method (str) – choose electron density model, either ‘ymw16’ or ‘ne2001’
- mode (str) – Gal, MC, or IGM (for YMW16 only)
- nu (float in GHz or astropy.Quantity) – observing frequency (GHz)
Returns: Dispersion measure (pc / cm3), scattering timescale at 1 GHz (s)
Return type: dm (astropy.Quantity), tau_sc (astropy.Quantity)
-
pygedm.pygedm.
dm_to_dist
(gl, gb, dm, dm_host=0, mode='gal', method='ymw16', nu=1.0)¶ Convert a DM to a distance
Parameters: - gl (float in deg or astropy.Angle) – galactic longitude
- gb (float in deg or astropy.Angle) – galactic latitude
- dm (float in pc/cm3 or astropy.Quantity) – dispersion measure (pc cm^-3)
- method (str) – choose electron density model, either ‘ymw16’ or ‘ne2001’
- mode (str) – Gal, MC, or IGM (for YMW16 only)
- nu (float in GHz or astropy.Quantity) – observing frequency (GHz)
Returns: Distance (pc), scattering timescale at 1 GHz (s)
Return type: dist (astropy.Quantity), tau_sc (astropy.Quantity)
-
pygedm.pygedm.
generate_healpix_dm_map
(dist=1, nside=64, method='ymw16')¶ Generate an all-sky healpix map for a given distance.
Parameters: - dist (float or Quantity) – Distance to integrate EDM out to. 30 kpc will go to edge
- nside (int) – The NSIDE parameter for the healpix map (power of 2, larger => higher resolution)
- method (str) – one of ‘ymw16’, ‘ne2001’, ‘yt2020’ or ‘yt2020_analytic’
Notes
This method takes a fair amount of time to run – tens of seconds for NSIDE=32. YT2020 method is even slower, consider using yt2020_analytic
Returns: Healpix map as a numpy array (1D), which can be viewed using the healpy.mollview() method Return type: hmap (np.array)