add more refs

nanograv · Apr 12, 2020 · 25846da · 25846da
1 parent 2c4a81f
commit 25846da
Showing 1 changed file with 13 additions and 6 deletions.
diff --git a/src/pint/models/solar_wind_dispersion.py b/src/pint/models/solar_wind_dispersion.py
@@ -22,6 +22,7 @@ class SolarWindDispersion(Dispersion):
     References
     ----------
     Madison et al. 2019, ApJ, 872, 150; Section 3.1
+    Edwards et al. 2006, MNRAS, 372, 1549; Setion 2.5.4
 
     """
 
@@ -54,6 +55,12 @@ def setup(self):
     def solar_wind_delay(self, toas, acc_delay=None):
         """Return the solar wind dispersion delay for a set of frequencies
         Eventually different solar wind models will be supported
+
+        Implements equations 29, 30 of Edwards et al. 2006,
+        where their rho is given as theta here
+
+        rvec: radial vector from observatory to the center of the Sun
+        pos: pulsar position
         """
         if self.SWM.value == 0:
             tbl = toas.table
@@ -63,10 +70,10 @@ def solar_wind_delay(self, toas, acc_delay=None):
                 warn("Using topocentric frequency for dedispersion!")
                 bfreq = tbl["freq"]
 
-            rsa = tbl["obs_sun_pos"].quantity
+            rvec = tbl["obs_sun_pos"].quantity
             pos = self.ssb_to_psb_xyz_ICRS(epoch=tbl["tdbld"].astype(np.float64))
-            r = np.sqrt(np.sum(rsa * rsa, axis=1))
-            cos_theta = (np.sum(rsa * pos, axis=1) / r).to(u.Unit("")).value
+            r = np.sqrt(np.sum(rvec * rvec, axis=1))
+            cos_theta = (np.sum(rvec * pos, axis=1) / r).to(u.Unit("")).value
             ret = (
                 const.au ** 2.0
                 * np.arccos(cos_theta)
@@ -91,10 +98,10 @@ def d_delay_d_ne_sw(self, toas, param_name, acc_delay=None):
                 warn("Using topocentric frequency for solar wind dedispersion!")
                 bfreq = tbl["freq"]
 
-            rsa = tbl["obs_sun_pos"].quantity
+            rvec = tbl["obs_sun_pos"].quantity
             pos = self.ssb_to_psb_xyz_ICRS(epoch=tbl["tdbld"].astype(np.float64))
-            r = np.sqrt(np.sum(rsa * rsa, axis=1))
-            cos_theta = np.sum(rsa * pos, axis=1) / r
+            r = np.sqrt(np.sum(rvec * rvec, axis=1))
+            cos_theta = np.sum(rvec * pos, axis=1) / r
 
             # ret = AUdist**2.0 / const.c * np.arccos(cos_theta) * DMconst / \
             ret = (