diff --git a/chandra_aca/aca_image.py b/chandra_aca/aca_image.py
index 6111710..e93b49c 100644
--- a/chandra_aca/aca_image.py
+++ b/chandra_aca/aca_image.py
@@ -279,12 +279,18 @@ def centroid_fm(self, bgd=None, pix_zero_loc="center", norm_clip=None):
         Note that the returned ``norm`` is the sum of the background-subtracted 6x6
         mousebitten image, not the entire image.
 
-        :param bgd: background to subtract, scalar or NxN ndarray (float)
-        :param pix_zero_loc: row/col coords are integral at 'edge' or 'center'
-        :param norm_clip: clip image norm at this min value (default is None and
-                          implies Exception for non-positive norm)
-
-        :returns: row, col, norm float
+        Parameters
+        ----------
+        bgd
+            background to subtract, scalar or NxN ndarray (float)
+        pix_zero_loc
+            row/col coords are integral at 'edge' or 'center'
+        norm_clip : clip image norm at this min value (default is None and
+            implies Exception for non-positive norm)
+
+        Returns
+        -------
+        row, col, norm float
         """
         row, col, norm = centroid_fm(
             self, bgd=bgd, pix_zero_loc=pix_zero_loc, norm_clip=norm_clip
@@ -364,10 +370,14 @@ def flicker_init(self, flicker_mean_time=10000, flicker_scale=1.0, seed=None):
 
         The model was reviewed and approved at SS&AWG on 2019-05-22.
 
-        :param flicker_mean_time: mean flickering time (sec, default=10000)
-        :param flicker_scale: multiplicative factor beyond model default for
-               flickering amplitude (default=1.0)
-        :param seed: random seed for reproducibility (default=None => no seed)
+        Parameters
+        ----------
+        flicker_mean_time
+            mean flickering time (sec, default=10000)
+        flicker_scale : multiplicative factor beyond model default for
+            flickering amplitude (default=1.0)
+        seed
+            random seed for reproducibility (default=None => no seed)
         """
         if not hasattr(self, "flicker_cdf_bins"):
             self._read_flicker_cdfs()
@@ -424,8 +434,12 @@ def flicker_update(self, dt, use_numba=True):
         to use the numba-based version which is about 6 times faster.  The
         vectorized version is left in for reference.
 
-        :param dt: time (secs) to propagate image
-        :param use_numba: use the numba version of updating (default=True)
+        Parameters
+        ----------
+        dt
+            time (secs) to propagate image
+        use_numba
+            use the numba version of updating (default=True)
         """
         if not hasattr(self, "flicker_times"):
             self.flicker_init()
@@ -603,13 +617,20 @@ def centroid_fm(img, bgd=None, pix_zero_loc="center", norm_clip=None):
     (0.0, 0.0) at the center (for pix_zero_loc='center') or the lower-left edge
     (for pix_zero_loc='edge').
 
-    :param img: NxN ndarray
-    :param bgd: background to subtract, float of NXN ndarray
-    :param pix_zero_loc: row/col coords are integral at 'edge' or 'center'
-    :param norm_clip: clip image norm at this min value (default is None and
-                      implies Exception for non-positive norm)
-
-    :returns: row, col, norm float
+    Parameters
+    ----------
+    img
+        NxN ndarray
+    bgd
+        background to subtract, float of NXN ndarray
+    pix_zero_loc
+        row/col coords are integral at 'edge' or 'center'
+    norm_clip : clip image norm at this min value (default is None and
+        implies Exception for non-positive norm)
+
+    Returns
+    -------
+    row, col, norm float
     """
     # Cast to an ndarray (without copying)
     img = img.view(np.ndarray)
@@ -707,14 +728,24 @@ def get_psf_image(
         Get interpolated ACA PSF image that corresponds to pixel location
         ``row``, ``col``.
 
-        :param row: (float) row value of PSF centroid
-        :param col: (float) col value of PSF centroid
-        :param norm: (float) summed intensity of PSF image
-        :param pix_zero_loc: row/col coords are integral at 'edge' or 'center'
-        :param interpolation: 'nearest' | 'bilinear' (default)
-        :param aca_image: return ACAImage if True, else return ndarray
-
-        :returns: ACAImage if (aca_image is True) else (ndarray image, row0, col0)
+        Parameters
+        ----------
+        row
+            (float) row value of PSF centroid
+        col
+            (float) col value of PSF centroid
+        norm
+            (float) summed intensity of PSF image
+        pix_zero_loc
+            row/col coords are integral at 'edge' or 'center'
+        interpolation
+            'nearest' | 'bilinear' (default)
+        aca_image
+            return ACAImage if True, else return ndarray
+
+        Returns
+        -------
+        ACAImage if (aca_image is True) else (ndarray image, row0, col0)
         """
         drc = self.drc
 
diff --git a/chandra_aca/attitude.py b/chandra_aca/attitude.py
index fb665e8..22f5427 100644
--- a/chandra_aca/attitude.py
+++ b/chandra_aca/attitude.py
@@ -32,13 +32,20 @@ def calc_roll(yag, zag, yag_obs, zag_obs, sigma=None):
     star inputs then one can supply an array of sigma values corresponding
     to each star.
 
-    :param yag: reference yag (list or array)
-    :param zag: reference zag (list or array)
-    :param yag_obs: observed yag (list or array)
-    :param zag_obs: observed zag (list or array)
-
-    :returns: roll (deg)
-
+    Parameters
+    ----------
+    yag
+        reference yag (list or array)
+    zag
+        reference zag (list or array)
+    yag_obs
+        observed yag (list or array)
+    zag_obs
+        observed zag (list or array)
+
+    Returns
+    -------
+    roll (deg)
     """
     yag = np.asarray(yag)
     zag = np.asarray(zag)
@@ -80,14 +87,22 @@ def calc_roll_pitch_yaw(yag, zag, yag_obs, zag_obs, sigma=None):
     [[cos(th) -sin(th)], [sin(th), cos(th)]] to [[1, -th], [th, 1]].
     In practice anything below 1.0 degree is fine.
 
-    :param yag: reference yag (list or array, arcsec)
-    :param zag: reference zag (list or array, arcsec)
-    :param yag_obs: observed yag (list or array, arcsec)
-    :param zag_obs: observed zag (list or array, arcsec)
-    :param sigma: centroid uncertainties (None or list or array, arcsec)
-
-    :returns: roll, pitch, yaw (degrees)
-
+    Parameters
+    ----------
+    yag
+        reference yag (list or array, arcsec)
+    zag
+        reference zag (list or array, arcsec)
+    yag_obs
+        observed yag (list or array, arcsec)
+    zag_obs
+        observed zag (list or array, arcsec)
+    sigma
+        centroid uncertainties (None or list or array, arcsec)
+
+    Returns
+    -------
+    roll, pitch, yaw (degrees)
     """
     yag = np.array(yag)
     zag = np.array(zag)
@@ -208,15 +223,24 @@ def calc_att(att, yag, zag, yag_obs, zag_obs, sigma=None):
     [[cos(th) -sin(th)], [sin(th), cos(th)]] to [[1, -th], [th, 1]].
     In practice anything below 1.0 degree is fine.
 
-    :param att: reference attitude (Quat-compatible)
-    :param yag: reference yag (list or array, arcsec)
-    :param zag: reference zag (list or array, arcsec)
-    :param yag_obs: observed yag (list or array, arcsec)
-    :param zag_obs: observed zag (list or array, arcsec)
-    :param sigma: centroid uncertainties (None or list or array, arcsec)
-
-    :returns: Quat or list of Quat
-
+    Parameters
+    ----------
+    att
+        reference attitude (Quat-compatible)
+    yag
+        reference yag (list or array, arcsec)
+    zag
+        reference zag (list or array, arcsec)
+    yag_obs
+        observed yag (list or array, arcsec)
+    zag_obs
+        observed zag (list or array, arcsec)
+    sigma
+        centroid uncertainties (None or list or array, arcsec)
+
+    Returns
+    -------
+    Quat or list of Quat
     """
     from Quaternion import Quat
 
diff --git a/chandra_aca/centroid_resid.py b/chandra_aca/centroid_resid.py
index 558b0a3..d19a5fb 100644
--- a/chandra_aca/centroid_resid.py
+++ b/chandra_aca/centroid_resid.py
@@ -94,10 +94,16 @@ def set_centroids(self, source, slot, alg=8, apply_dt=True):
 
         yag, zag, yag_times an zag_times can also be set directly without use of this method.
 
-        :param source: 'ground' | 'obc'
-        :param slot: ACA slot
-        :param alg: for ground processing, use centroids from this algorithm.
-        :param apply_dt: apply centroid time offsets via 'set_offsets'
+        Parameters
+        ----------
+        source
+            'ground' | 'obc'
+        slot
+            ACA slot
+        alg
+            for ground processing, use centroids from this algorithm.
+        apply_dt
+            apply centroid time offsets via 'set_offsets'
         """
         self.centroid_source = source
         self.centroid_dt = None
diff --git a/chandra_aca/dark_model.py b/chandra_aca/dark_model.py
index d67d3d4..0f30833 100644
--- a/chandra_aca/dark_model.py
+++ b/chandra_aca/dark_model.py
@@ -74,11 +74,18 @@ def dark_temp_scale(t_ccd, t_ccd_ref=-19.0, scale_4c=None):
       >>> dark_temp_scale(t_ccd=-18, t_ccd_ref=-10, scale_4c=2.0)
       4.0
 
-    :param t_ccd: actual temperature (degC)
-    :param t_ccd_ref: reference temperature (degC, default=-19.0)
-    :param scale_4c: increase in dark current per 4 degC increase (default=1.0 / 0.7)
-
-    :returns: scale factor
+    Parameters
+    ----------
+    t_ccd
+        actual temperature (degC)
+    t_ccd_ref
+        reference temperature (degC, default=-19.0)
+    scale_4c
+        increase in dark current per 4 degC increase (default=1.0 / 0.7)
+
+    Returns
+    -------
+    scale factor
     """
     if scale_4c is None:
         scale_4c = DARK_SCALE_4C
@@ -90,10 +97,16 @@ def get_dark_hist(date, t_ccd):
     """
     Return the dark current histogram corresponding to ``date`` and ``t_ccd``.
 
-    :param date: date in any DateTime format
-    :param t_ccd: CCD temperature (deg C)
+    Parameters
+    ----------
+    date
+        date in any DateTime format
+    t_ccd
+        CCD temperature (deg C)
 
-    :returns: bin_centers, bins, darkhist
+    Returns
+    -------
+    bin_centers, bins, darkhist
     """
     pars = get_sbp_pars(date)
     x = darkbins.bin_centers
@@ -193,8 +206,14 @@ def get_sbp_pars(dates):
 
     This is only accurate for dates > 2007.0.
 
-    :param dates: one or a list of date(s) in DateTime compatible format
-    :returns: one or a list of parameter lists [g1, g2, g3, x_b, ampl]
+    Parameters
+    ----------
+    dates
+        one or a list of date(s) in DateTime compatible format
+
+    Returns
+    -------
+    one or a list of parameter lists [g1, g2, g3, x_b, ampl]
     """
     dates, is_scalar = as_array(dates)
 
@@ -231,10 +250,18 @@ def get_warm_fracs(warm_threshold, date="2013:001:12:00:00", T_ccd=-19.0):
     Calculate fraction of pixels in modeled dark current distribution
     above warm threshold(s).
 
-    :param warm_threshold: scalar or list of threshold(s) in e-/sec
-    :param date: date to use for modeled dark current distribution/histogram
-    :param T_ccd: temperature (C) of modeled dark current distribution
-    :returns: list or scalar of warm fractions (depends on warm_threshold type)
+    Parameters
+    ----------
+    warm_threshold
+        scalar or list of threshold(s) in e-/sec
+    date
+        date to use for modeled dark current distribution/histogram
+    T_ccd
+        temperature (C) of modeled dark current distribution
+
+    Returns
+    -------
+    list or scalar of warm fractions (depends on warm_threshold type)
     """
 
     x, xbins, y = get_dark_hist(date, T_ccd)
@@ -268,8 +295,12 @@ def synthetic_dark_image(date, t_ccd_ref=None):
     Generate a synthetic dark current image corresponding to the specified
     ``date`` and ``t_ccd``.
 
-    :param date: (DateTime compatible)
-    :param t_ccd_ref: ACA CCD temperature
+    Parameters
+    ----------
+    date
+        (DateTime compatible)
+    t_ccd_ref
+        ACA CCD temperature
     """
 
     from mica.archive.aca_dark import get_dark_cal_image
diff --git a/chandra_aca/drift.py b/chandra_aca/drift.py
index 3940503..914c439 100644
--- a/chandra_aca/drift.py
+++ b/chandra_aca/drift.py
@@ -122,9 +122,16 @@ def calc(self, times, t_ccd):
         SIM DY or DZ values in mm.  This can be converted to a drift in arcsec via
         the scale factor 20.493 arcsec/mm.
 
-        :param times: array of times (CXC secs)
-        :param t_ccd: CCD temperatures (degC)
-        :returns: array of ASOL SIM DY/DZ (mm)
+        Parameters
+        ----------
+        times
+            array of times (CXC secs)
+        t_ccd
+            CCD temperatures (degC)
+
+        Returns
+        -------
+        array of ASOL SIM DY/DZ (mm)
         """
         # The drift model is calibrated assuming t_ccd is in degF, but we want inputs
         # in degC, so convert at this point.
@@ -164,14 +171,24 @@ def get_aca_offsets(detector, chip_id, chipx, chipy, time, t_ccd):
 
     The ``time`` and ``t_ccd`` inputs can be either scalars or arrays.
 
-    :param detector: one of ACIS-I, ACIS-S, HRC-I, HRC-S
-    :param chipx: zero-offset aimpoint CHIPX
-    :param chipy: zero-offset aimpoint CHIPY
-    :param chip_id: zero-offset aimpoint CHIP ID
-    :param time: time(s) of observation (any Chandra.Time compatible format)
-    :param t_ccd: ACA CCD temperature(s) (degC)
-
-    :returns: aca_offset_y, aca_offset_z (arcsec)
+    Parameters
+    ----------
+    detector
+        one of ACIS-I, ACIS-S, HRC-I, HRC-S
+    chipx
+        zero-offset aimpoint CHIPX
+    chipy
+        zero-offset aimpoint CHIPY
+    chip_id
+        zero-offset aimpoint CHIP ID
+    time
+        time(s) of observation (any Chandra.Time compatible format)
+    t_ccd
+        ACA CCD temperature(s) (degC)
+
+    Returns
+    -------
+    aca_offset_y, aca_offset_z (arcsec)
     """
     # Define model instances using calibrated parameters
     drift_y = AcaDriftModel(**DRIFT_PARS["dy"])
@@ -224,7 +241,9 @@ def get_default_zero_offset_table():
     Note the definitive source of this file is:
     https://icxc.harvard.edu/mp/html/aimpoint_table/zero_offset_aimpoints.txt.
 
-    :returns: zero offset aimpoint table as astropy.Table
+    Returns
+    -------
+    zero offset aimpoint table as astropy.Table
     """
     try:
         path = (
@@ -247,14 +266,22 @@ def get_target_aimpoint(date, cycle, detector, too=False, zero_offset_table=None
     """
     Given date, proposal cycle, and detector, return aimpoint chipx, chipy, chip_id
 
-    :param date: observation date
-    :param cycle: proposal cycle of observation
-    :param detector: target detector
-    :param too: boolean. If target is TOO use current cycle not proposal cycle.
-    :param zero_offset_able: table (astropy or numpy) of zero offset aimpoint table
-                             defaults to official SOT MP version if not supplied.
-
-    :returns: tuple of chipx, chipy, chip_id
+    Parameters
+    ----------
+    date
+        observation date
+    cycle
+        proposal cycle of observation
+    detector
+        target detector
+    too
+        boolean. If target is TOO use current cycle not proposal cycle.
+    zero_offset_able : table (astropy or numpy) of zero offset aimpoint table
+        defaults to official SOT MP version if not supplied.
+
+    Returns
+    -------
+    tuple of chipx, chipy, chip_id
     """
     if zero_offset_table is None:
         zero_offset_table = get_default_zero_offset_table()
diff --git a/chandra_aca/maude_decom.py b/chandra_aca/maude_decom.py
index 339818c..e04dc55 100644
--- a/chandra_aca/maude_decom.py
+++ b/chandra_aca/maude_decom.py
@@ -387,8 +387,14 @@ def _aca_header_1(bits):
         """
         Unpack ACA header 1 (ACA User Manual 5.3.2.2.1).
 
-        :param bits: bytes-like object of length 7
-        :return: dict
+        Parameters
+        ----------
+        bits
+            bytes-like object of length 7
+
+        Returns
+        -------
+        dict
         """
         bits = np.unpackbits(np.array(_unpack("BBBbbBB", bits), dtype=np.uint8))
         return {
@@ -413,8 +419,14 @@ def _aca_header_2(bits):
         """
         Unpack ACA header 2 (ACA User Manual 5.3.2.2.2).
 
-        :param bits: bytes-like object of length 7
-        :return: dict
+        Parameters
+        ----------
+        bits
+            bytes-like object of length 7
+
+        Returns
+        -------
+        dict
         """
         bits = _unpack("BbbbbBB", bits)
         c = np.unpackbits(np.array(bits[:2], dtype=np.uint8))
@@ -437,8 +449,14 @@ def _aca_header_3(bits):
         """
         Unpack ACA header 3 (ACA User Manual 5.3.2.2.3).
 
-        :param bits: bytes-like object of length 7
-        :return: dict
+        Parameters
+        ----------
+        bits
+            bytes-like object of length 7
+
+        Returns
+        -------
+        dict
         """
         return {"DIAGNOSTIC": _unpack("BBBBBB", bits[1:])}
 
@@ -447,9 +465,14 @@ def unpack_aca_telemetry(packet):
     """
     Unpack ACA telemetry encoded in 225-byte packets.
 
-    :param packet: bytes
-    :return: list of dict
+    Parameters
+    ----------
+    packet
+        bytes
 
+    Returns
+    -------
+    list of dict
     A list of length 8, one entry per slot, where each entry is a dictionary.
     """
     s1, s2, s3 = _unpack("BBB", packet[5:8])
@@ -512,8 +535,14 @@ def _combine_aca_packets(aca_packets):
 
     This is intended to combine the two 6X6 packets or the four 8X8 packets.
 
-    :param aca_packets: list of dict
-    :return: dict
+    Parameters
+    ----------
+    aca_packets
+        list of dict
+
+    Returns
+    -------
+    dict
     """
     # note that they are reverse-sorted so the first frame overwrites the others if they collide
     aca_packets = sorted(aca_packets, key=lambda p: p["TIME"], reverse=True)
@@ -540,9 +569,16 @@ def _group_packets(packets, discard=True):
     Before decommuting an ACA package we group the ACA-related portion of VCDU frames to form the
     one 225-byte ACA packet.
 
-    :param packets: list of ACA sub-packets
-    :param discard: bool to discard incomplete ACA packets
-    :return: list of ACA packets
+    Parameters
+    ----------
+    packets
+        list of ACA sub-packets
+    discard
+        bool to discard incomplete ACA packets
+
+    Returns
+    -------
+    list of ACA packets
     """
     res = []
     n = None
@@ -579,14 +615,22 @@ def get_raw_aca_packets(start, stop, maude_result=None, **maude_kwargs):
     This function raises an exception if the VCDU frames in maude_result are not contiguous, which
     can also happen if the frames are corrupted in some way.
 
-    :param start: timestamp interpreted as a Chandra.Time.DateTime
-    :param stop: timestamp interpreted as a Chandra.Time.DateTime
-    :param maude_result: the result of calling maude.get_frames. Optional.
-    :param maude_kwargs: keyword args passed to maude.get_frames()
-    :return: dict
-
-        {'flags': int, 'packets': [],
-         'TIME': np.array([]), 'MNF': np.array([]), 'MJF': np.array([])}
+    Parameters
+    ----------
+    start
+        timestamp interpreted as a Chandra.Time.DateTime
+    stop
+        timestamp interpreted as a Chandra.Time.DateTime
+    maude_result
+        the result of calling maude.get_frames. Optional.
+    maude_kwargs
+        keyword args passed to maude.get_frames()
+
+    Returns
+    -------
+    dict
+    {'flags': int, 'packets': [],
+    'TIME': np.array([]), 'MNF': np.array([]), 'MJF': np.array([])}
     """
     date_start, date_stop = DateTime(start), DateTime(
         stop
@@ -721,9 +765,16 @@ def _aca_packets_to_table(aca_packets, dtype=None):
     """
     Store ACA packets in a table.
 
-    :param aca_packets: list of dict
-    :param dtype: dtype to use in the resulting table. Optional.
-    :return: astropy.table.Table
+    Parameters
+    ----------
+    aca_packets
+        list of dict
+    dtype
+        dtype to use in the resulting table. Optional.
+
+    Returns
+    -------
+    astropy.table.Table
     """
     import copy
 
@@ -894,32 +945,39 @@ def get_aca_packets(
          [False False False False False False False False]],
                     fill_value = 1e+20)
 
-
-    :param start: timestamp interpreted as a Chandra.Time.DateTime
-    :param stop: timestamp interpreted as a Chandra.Time.DateTime
-    :param level0: bool.
+    Parameters
+    ----------
+    start
+        timestamp interpreted as a Chandra.Time.DateTime
+    stop
+        timestamp interpreted as a Chandra.Time.DateTime
+    level0 : bool.
         Implies combine=True, adjust_time=True, calibrate=True
-    :param combine: bool.
+    combine : bool.
         If True, multiple ACA packets are combined to form an image (depending on size),
         If False, ACA packets are not combined, resulting in multiple lines for 6x6 and 8x8 images.
-    :param adjust_time: bool
+    adjust_time : bool
         If True, half the integration time is subtracted
-    :param calibrate: bool
+    calibrate : bool
         If True, pixel values will be 'value * imgscale / 32 - 50' and temperature values will
         be: 0.4 * value + 273.15
-    :param blobs: bool or dict
+    blobs : bool or dict
         If set, data is assembled from MAUDE blobs. If it is a dictionary, it must be the
         output of maude.get_blobs ({'blobs': ... }).
-    :param frames: bool or dict
+    frames : bool or dict
         If set, data is assembled from MAUDE frames. If it is a dictionary, it must be the
         output of maude.get_frames ({'data': ... }).
-    :param dtype: np.dtype. Optional.
+    dtype : np.dtype. Optional.
         the dtype to use when creating the resulting table. This is useful to add columns
         including MSIDs that are present in blobs. If used with frames, most probably you will get
         and empty column. This option is intended to augment the default dtype. If a more
         restrictive dtype is used, a KeyError can be raised.
-    :param maude_kwargs: keyword args passed to maude
-    :return: astropy.table.Table
+    maude_kwargs
+        keyword args passed to maude
+
+    Returns
+    -------
+    astropy.table.Table
     """
     if not blobs and not frames:
         frames = True
@@ -1085,10 +1143,18 @@ def get_aca_images(start, stop, **maude_kwargs):
     """
     Fetch ACA image telemetry
 
-    :param start: timestamp interpreted as a Chandra.Time.DateTime
-    :param stop: timestamp interpreted as a Chandra.Time.DateTime
-    :param maude_kwargs: keyword args passed to maude
-    :return: astropy.table.Table
+    Parameters
+    ----------
+    start
+        timestamp interpreted as a Chandra.Time.DateTime
+    stop
+        timestamp interpreted as a Chandra.Time.DateTime
+    maude_kwargs
+        keyword args passed to maude
+
+    Returns
+    -------
+    astropy.table.Table
     """
     return get_aca_packets(start, stop, level0=True, **maude_kwargs)
 
@@ -1111,12 +1177,21 @@ def get_raw_aca_blobs(start, stop, maude_result=None, **maude_kwargs):
 
     This is to blobs what `get_raw_aca_packets` is to frames.
 
-    :param start: timestamp interpreted as a Chandra.Time.DateTime
-    :param stop: timestamp interpreted as a Chandra.Time.DateTime
-    :param maude_result: the result of calling maude.get_blobs. Optional.
-    :param maude_kwargs: keyword args passed to maude.get_frames()
-    :return: dict
-        {'blobs': [], 'names': np.array([]), 'types': np.array([])}
+    Parameters
+    ----------
+    start
+        timestamp interpreted as a Chandra.Time.DateTime
+    stop
+        timestamp interpreted as a Chandra.Time.DateTime
+    maude_result
+        the result of calling maude.get_blobs. Optional.
+    maude_kwargs
+        keyword args passed to maude.get_frames()
+
+    Returns
+    -------
+    dict
+    {'blobs': [], 'names': np.array([]), 'types': np.array([])}
     """
     date_start, date_stop = DateTime(start), DateTime(
         stop
@@ -1172,10 +1247,18 @@ def blob_to_aca_image_dict(blob, imgnum, pea=1):
 
     This does to blobs what unpack_aca_telemetry does to frames, but for a single image.
 
-    :param blob:
-    :param imgnum:
-    :param pea:
-    :return:
+    Parameters
+    ----------
+    blob
+
+    imgnum
+
+    pea
+
+
+    Returns
+    -------
+
     """
     global_msids = ACA_MSID_LIST[pea]
     slot_msids = ACA_SLOT_MSID_LIST[pea][imgnum]
diff --git a/chandra_aca/planets.py b/chandra_aca/planets.py
index 359fe95..10ada92 100644
--- a/chandra_aca/planets.py
+++ b/chandra_aca/planets.py
@@ -115,18 +115,22 @@ def get_planet_angular_sep(
     - 'chandra' (reasonably accurate fast, requires fetching ephemeris)
     - 'chandra-horizons' (most accurate, slow, requires internet access)
 
-    :param body: str
+    Parameters
+    ----------
+    body : str
         Body name (lower case planet name)
-    :param ra: float
+    ra : float
         RA in degrees
-    :param dec: float
+    dec : float
         Dec in degrees
-    :param time: CxoTime-compatible object
+    time : CxoTime-compatible object
         Time or times of observation
-    :param observer_position: str
+    observer_position : str
         Observer position ('earth', 'chandra', or 'chandra-horizons')
 
-    :returns: angular separation (deg)
+    Returns
+    -------
+    angular separation (deg)
     """
     from agasc import sphere_dist
 
@@ -169,9 +173,16 @@ def get_planet_barycentric(body: str, time: CxoTimeLike = None):
     ``body`` must be one of "sun", "mercury", "venus", "earth-moon-barycenter", "earth",
     "moon", "mars", "jupiter", "saturn", "uranus", "neptune", or "pluto".
 
-    :param body: Body name (lower case planet name)
-    :param time: Time or times for returned position (default=NOW)
-    :returns: barycentric position (km) as (x, y, z) or N x (x, y, z)
+    Parameters
+    ----------
+    body
+        Body name (lower case planet name)
+    time
+        Time or times for returned position (default=NOW)
+
+    Returns
+    -------
+    barycentric position (km) as (x, y, z) or N x (x, y, z)
     """
     kernel = KERNEL.val
     if body not in BODY_NAME_TO_KERNEL_SPEC:
@@ -211,10 +222,19 @@ def get_planet_eci(
     - Jupiter: < 1 arcmin with peak around 0.5 arcmin
     - Saturn: < 0.5 arcmin with peak around 0.3 arcmin
 
-    :param body: Body name (lower case planet name)
-    :param time: Time or times for returned position (default=NOW)
-    :param pos_observer: Observer position (default=Earth)
-    :returns: Earth-Centered Inertial (ECI) position (km) as (x, y, z)
+    Parameters
+    ----------
+    body
+        Body name (lower case planet name)
+    time
+        Time or times for returned position (default=NOW)
+    pos_observer
+        Observer position (default=Earth)
+
+    Returns
+    -------
+    ndarray
+        Earth-Centered Inertial (ECI) position (km) as (x, y, z)
         or N x (x, y, z)
     """
     time = CxoTime(time)
@@ -248,9 +268,16 @@ def get_planet_chandra(body: str, time: CxoTimeLike = None):
     - Jupiter: < 0.8 arcsec
     - Saturn: < 0.5 arcsec
 
-    :param body: Body name
-    :param time: Time or times for returned position (default=NOW)
-    :returns: position relative to Chandra (km) as (x, y, z) or N x (x, y, z)
+    Parameters
+    ----------
+    body
+        Body name
+    time
+        Time or times for returned position (default=NOW)
+
+    Returns
+    -------
+    position relative to Chandra (km) as (x, y, z) or N x (x, y, z)
     """
     from cheta import fetch
 
@@ -322,15 +349,22 @@ def get_planet_chandra_horizons(
       2020:001:18:00:00.000 277.15181 -23.19347      34.51       2.03  -1.839         5.408    31.76
       2020:002:00:00:00.000 277.21454 -23.18970      34.69       2.51  -1.839         5.408    31.76
 
-    :param body: one of 'mercury', 'venus', 'mars', 'jupiter', 'saturn',
+    Parameters
+    ----------
+    body : one of 'mercury', 'venus', 'mars', 'jupiter', 'saturn',
         'uranus', 'neptune', or any other body that Horizons supports.
-    :param timestart: start time (any CxoTime-compatible time)
-    :param timestop: stop time (any CxoTime-compatible time)
-    :param n_times: number of time samples (inclusive, default=10)
-    :param timeout: timeout for query to Horizons (secs)
-
-    :returns: Table of information
-
+    timestart
+        start time (any CxoTime-compatible time)
+    timestop
+        stop time (any CxoTime-compatible time)
+    n_times
+        number of time samples (inclusive, default=10)
+    timeout
+        timeout for query to Horizons (secs)
+
+    Returns
+    -------
+    Table of information
     """
     import requests
 
diff --git a/chandra_aca/plot.py b/chandra_aca/plot.py
index 7aa65ad..ba21a8c 100644
--- a/chandra_aca/plot.py
+++ b/chandra_aca/plot.py
@@ -87,10 +87,13 @@ def _plot_catalog_items(ax, catalog):
     Plot catalog items (guide, acq, bot, mon, fid) in yang and zang on the supplied
     axes object in place.
 
-    :param ax: matplotlib axes
-    :param catalog: data structure containing starcheck-style columns/attributes
-                    catalog records.  This can be anything that will work with
-                    astropy.table.Table(catalog).  A list of dicts is the convention.
+    Parameters
+    ----------
+    ax
+        matplotlib axes
+    catalog : data structure containing starcheck-style columns/attributes
+        catalog records.  This can be anything that will work with
+        astropy.table.Table(catalog).  A list of dicts is the convention.
     """
     cat = Table(catalog)
     cat["row"], cat["col"] = yagzag_to_pixels(cat["yang"], cat["zang"], allow_bad=True)
@@ -156,11 +159,18 @@ def _plot_field_stars(ax, stars, attitude, red_mag_lim=None, bad_stars=None):
     Plot plot field stars in yang and zang on the supplied
     axes object in place.
 
-    :param ax: matplotlib axes
-    :param stars: astropy.table compatible set of records of agasc entries of stars
-    :param attitude: Quaternion-compatible attitude
-    :param red_mag_lim: faint limit
-    :param bad_stars: boolean mask of stars to be plotted in red
+    Parameters
+    ----------
+    ax
+        matplotlib axes
+    stars
+        astropy.table compatible set of records of agasc entries of stars
+    attitude
+        Quaternion-compatible attitude
+    red_mag_lim
+        faint limit
+    bad_stars
+        boolean mask of stars to be plotted in red
     """
     stars = Table(stars)
     quat = Quaternion.Quat(attitude)
@@ -248,29 +258,41 @@ def plot_stars(
     Plot a catalog, a star field, or both in a matplotlib figure.
     If supplying a star field, an attitude must also be supplied.
 
-    :param attitude: A Quaternion compatible attitude for the pointing
-    :param catalog: Records describing catalog.  Must be astropy table compatible.
+    Parameters
+    ----------
+    attitude
+        A Quaternion compatible attitude for the pointing
+    catalog : Records describing catalog.  Must be astropy table compatible.
         Required fields are ['idx', 'type', 'yang', 'zang', 'halfw']
-    :param stars: astropy table compatible set of agasc records of stars
+    stars : astropy table compatible set of agasc records of stars
         Required fields are ['RA_PMCORR', 'DEC_PMCORR', 'MAG_ACA', 'MAG_ACA_ERR'].
         If bad_acq_stars will be called (bad_stars is None), additional required fields
         ['CLASS', 'ASPQ1', 'ASPQ2', 'ASPQ3', 'VAR', 'POS_ERR']
         If stars is None, stars will be fetched from the AGASC for the
         supplied attitude.
-    :param title: string to be used as suptitle for the figure
-    :param starcat_time: DateTime-compatible time.  Used in ACASC fetch for proper
+    title
+        string to be used as suptitle for the figure
+    starcat_time : DateTime-compatible time.  Used in ACASC fetch for proper
         motion correction.  Not used if stars is not None.
-    :param red_mag_lim: faint limit for field star plotting.
-    :param quad_bound: boolean, plot inner quadrant boundaries
-    :param grid: boolean, plot axis grid
-    :param bad_stars: boolean mask on 'stars' of those that don't meet minimum requirements
+    red_mag_lim
+        faint limit for field star plotting.
+    quad_bound
+        boolean, plot inner quadrant boundaries
+    grid
+        boolean, plot axis grid
+    bad_stars : boolean mask on 'stars' of those that don't meet minimum requirements
         to be selected as acq stars.  If None, bad_stars will be set by a call
         to bad_acq_stars().
-    :param plot_keepout: plot CCD area to be avoided in star selection (default=False)
-    :param ax: matplotlib axes object to use (optional)
-    :param duration: duration (starting at ``starcat_time``) for plotting planets
+    plot_keepout
+        plot CCD area to be avoided in star selection (default=False)
+    ax
+        matplotlib axes object to use (optional)
+    duration : duration (starting at ``starcat_time``) for plotting planets
         (secs, default=0)
-    :returns: matplotlib figure
+
+    Returns
+    -------
+    matplotlib figure
     """
     if stars is None:
         quat = Quaternion.Quat(attitude)
@@ -373,21 +395,26 @@ def plot_stars(
 
 def _plot_planets(ax, att, date0, duration, lim0, lim1):
     """
-    :param ax: plt.Axes
+
+    Parameters
+    ----------
+    ax : plt.Axes
         Matplotlib axes object to use
-    :param att: Quat
+    att : Quat
         Attitude quaternion
-    :param date0: CxoTime-compatible
+    date0 : CxoTime-compatible
         Date of obs start
-    :param duration: float
+    duration : float
         Duration of plot (secs)
-    :param lim0: float
+    lim0 : float
         Lower limit on x, y axis (row)
-    :param lim1: float
+    lim1 : float
         Upper limit on x, y axis (col)
 
-    :returns: boolean
-        True if planets were plotted
+    Returns
+    -------
+    boolean
+    True if planets were plotted
     """
     if not isinstance(att, Quat):
         att = Quat(att)
@@ -452,9 +479,14 @@ def bad_acq_stars(stars):
     """
     Return mask of 'bad' stars, by evaluating AGASC star parameters.
 
-    :param stars: astropy table-compatible set of agasc records of stars. Required fields
-          are ['CLASS', 'ASPQ1', 'ASPQ2', 'ASPQ3', 'VAR', 'POS_ERR']
-    :returns: boolean mask true for 'bad' stars
+    Parameters
+    ----------
+    stars : astropy table-compatible set of agasc records of stars. Required fields
+        are ['CLASS', 'ASPQ1', 'ASPQ2', 'ASPQ3', 'VAR', 'POS_ERR']
+
+    Returns
+    -------
+    boolean mask true for 'bad' stars
     """
     return (
         (stars["CLASS"] != 0)
@@ -472,8 +504,14 @@ def plot_compass(roll):
     """
     Make a compass plot.
 
-    :param roll: Attitude roll for compass plot.
-    :returns: matplotlib figure
+    Parameters
+    ----------
+    roll
+        Attitude roll for compass plot.
+
+    Returns
+    -------
+    matplotlib figure
     """
     fig = plt.figure(figsize=(3, 3))
     ax = plt.subplot(polar=True)
diff --git a/chandra_aca/star_probs.py b/chandra_aca/star_probs.py
index 5054d4c..8e41c15 100644
--- a/chandra_aca/star_probs.py
+++ b/chandra_aca/star_probs.py
@@ -76,8 +76,14 @@ def get_box_delta(halfwidth):
     Transform from halfwidth (arcsec) to the box_delta value which gets added
     to failure probability (in probit space).
 
-    :param halfwidth: scalar or ndarray of box sizes (halfwidth, arcsec)
-    :returns: box deltas
+    Parameters
+    ----------
+    halfwidth
+        scalar or ndarray of box sizes (halfwidth, arcsec)
+
+    Returns
+    -------
+    box deltas
     """
     # Coefficents for dependence of probability on search box size (halfwidth).  From:
     # https://github.com/sot/skanb/blob/master/pea-test-set/fit_box_size_acq_prob.ipynb
@@ -132,18 +138,30 @@ def t_ccd_warm_limit(
      - Tuple (n, prob): computed probability of acquiring ``n`` or fewer stars
          must not exceed ``prob``.
 
-    :param mags: list of star ACA mags
-    :param date: observation date (any Chandra.Time valid format)
-    :param colors: list of star B-V colors (optional, default=0.0)
-    :param halfwidths: list of acq box halfwidths(optional, default=120)
-    :param min_n_acq: float or tuple (see above)
-    :param cold_t_ccd: coldest CCD temperature to consider (default=-16 C)
-    :param warm_t_ccd: warmest CCD temperature to consider (default=-5 C)
-    :param model: probability model (see acq_success_prob() for allowed values, default)
-
-    :returns: (t_ccd, n_acq | prob_n_or_fewer) tuple with CCD temperature upper limit:
-              - number of expected ACQ stars at that temperature (scalar min_n_acq)
-              - probability of acquiring ``n`` or fewer stars (tuple min_n_acq)
+    Parameters
+    ----------
+    mags
+        list of star ACA mags
+    date
+        observation date (any Chandra.Time valid format)
+    colors
+        list of star B-V colors (optional, default=0.0)
+    halfwidths
+        list of acq box halfwidths(optional, default=120)
+    min_n_acq
+        float or tuple (see above)
+    cold_t_ccd
+        coldest CCD temperature to consider (default=-16 C)
+    warm_t_ccd
+        warmest CCD temperature to consider (default=-5 C)
+    model
+        probability model (see acq_success_prob() for allowed values, default)
+
+    Returns
+    -------
+    (t_ccd, n_acq | prob_n_or_fewer) tuple with CCD temperature upper limit:
+    - number of expected ACQ stars at that temperature (scalar min_n_acq)
+    - probability of acquiring ``n`` or fewer stars (tuple min_n_acq)
     """
 
     if isinstance(min_n_acq, tuple):
@@ -234,9 +252,10 @@ def prob_n_acq(star_probs):
     of length n_stars.  In addition the cumulative sum, which represents
     the probability of acquiring n_acq or fewer stars, is returned.
 
-    :param star_probs: array of star acq probabilities (list or ndarray)
-
-    :returns n_acq_probs, cum_n_acq_probs: tuple of ndarray, ndarray
+    Parameters
+    ----------
+    star_probs : array of star acq probabilities (list or ndarray)
+        :returns n_acq_probs, cum_n_acq_probs: tuple of ndarray, ndarray
     """
     star_probs = np.array(star_probs, dtype=np.float64)
     n_stars = len(star_probs)
@@ -376,14 +395,22 @@ def clip_and_warn(name, val, val_lo, val_hi, model):
     warning if clipping occurs.  The ``name`` and ``model`` are just used in
     the warning.
 
-    :param name: Value name
-    :param val: Value
-    :param val_lo: Minimum
-    :param val_hi: Maximum
-    :param model: Model name
-
-    :returns: Clipped value
-
+    Parameters
+    ----------
+    name
+        Value name
+    val
+        Value
+    val_lo
+        Minimum
+    val_hi
+        Maximum
+    model
+        Model name
+
+    Returns
+    -------
+    Clipped value
     """
     val = np.asarray(val)
     if np.any((val > val_hi) | (val < val_lo)):
@@ -404,8 +431,12 @@ def get_grid_axis_values(hdr, axis):
 
         linspace(hdr[f"{axis}_lo"], hdr[f"{axis}_hi"], n_vals)
 
-    :param hdr: FITS header (dict-like)
-    :param axis: Axis name (e.g. "mag")
+    Parameters
+    ----------
+    hdr
+        FITS header (dict-like)
+    axis
+        Axis name (e.g. "mag")
     """
     n_vals = hdr[f"{axis}_n"]
     if all(f"{axis}_{ii}" in hdr for ii in range(n_vals)):
@@ -442,10 +473,18 @@ def get_grid_func_model(
         "halfw_hi": upper bound of halfw axis
         "info": dict of provenance info for model file
 
-    :param model: Model name (optional), defaults to ``conf.default_model``
-    :param version: Version / tag / branch of ``chandra_models`` repository (optional)
-    :param repo_path: Path to ``chandra_models`` repository (optional)
-    :returns: dict of model data
+    Parameters
+    ----------
+    model
+        Model name (optional), defaults to ``conf.default_model``
+    version
+        Version / tag / branch of ``chandra_models`` repository (optional)
+    repo_path
+        Path to ``chandra_models`` repository (optional)
+
+    Returns
+    -------
+    dict of model data
     """
     if model is None:
         model = conf.default_model
@@ -547,15 +586,24 @@ def grid_model_acq_prob(
     This does a 3-d linear interpolation on mag, t_ccd, and halfwidth using a
     pre-computed gridded model that is stored in a FITS file.
 
-    :param mag: ACA magnitude (float or np.ndarray)
-    :param t_ccd: CCD temperature (degC, float or ndarray)
-    :param color: B-V color to check for B-V=1.5 => red star (float or np.ndarray)
-    :param halfwidth: search box halfwidth (arcsec, default=120, float or ndarray)
-    :param probit: if True then return Probit(p_success). Default=False
-    :param model: Model name, e.g. 'grid-floor-2018-11'
-
-    :returns: Acquisition success probability(s)
-
+    Parameters
+    ----------
+    mag
+        ACA magnitude (float or np.ndarray)
+    t_ccd
+        CCD temperature (degC, float or ndarray)
+    color
+        B-V color to check for B-V=1.5 => red star (float or np.ndarray)
+    halfwidth
+        search box halfwidth (arcsec, default=120, float or ndarray)
+    probit
+        if True then return Probit(p_success). Default=False
+    model
+        Model name, e.g. 'grid-floor-2018-11'
+
+    Returns
+    -------
+    Acquisition success probability(s)
     """
     # Get the grid model function and model parameters from a FITS file. This function
     # call is cached.
@@ -779,12 +827,20 @@ def sota_model_acq_prob(mag, warm_frac, color=0, halfwidth=120):
     https://github.com/sot/skanb/blob/master/pea-test-set/fit_box_size_acq_prob.ipynb
     https://github.com/sot/aca_stats/blob/master/fit_acq_prob_model-2017-07-sota.ipynb
 
-    :param mag: ACA magnitude (float or np.ndarray)
-    :param warm_frac: N100 warm fraction (float or np.ndarray)
-    :param color: B-V color to check for B-V=1.5 => red star (float or np.ndarray)
-    :param halfwidth: search box halfwidth (arcsec, default=120)
-
-    :returns: Acquisition success probability(s)
+    Parameters
+    ----------
+    mag
+        ACA magnitude (float or np.ndarray)
+    warm_frac
+        N100 warm fraction (float or np.ndarray)
+    color
+        B-V color to check for B-V=1.5 => red star (float or np.ndarray)
+    halfwidth
+        search box halfwidth (arcsec, default=120)
+
+    Returns
+    -------
+    Acquisition success probability(s)
     """
     #
     # NOTE: the "WITH_MS" are historical and no longer used in flight
@@ -877,12 +933,18 @@ def mag_for_p_acq(p_acq, date=None, t_ccd=-10.0, halfwidth=120, model=None):
     acquisition probability of ``p_acq``.  Star magnitude is defined/limited
     to the range 5.0 - 12.0 mag.
 
-    :param p_acq: acquisition probability (0 to 1.0)
-    :param date: observation date (any Chandra.Time valid format)
-    :param t_ccd: ACA CCD temperature (deg C)
-    :param halfwidth: search box halfwidth (arcsec, default=120)
-    :param model: probability model (see acq_success_prob() for allowed values, default)
-    :returns mag: star magnitude
+    Parameters
+    ----------
+    p_acq
+        acquisition probability (0 to 1.0)
+    date
+        observation date (any Chandra.Time valid format)
+    t_ccd
+        ACA CCD temperature (deg C)
+    halfwidth
+        search box halfwidth (arcsec, default=120)
+    model : probability model (see acq_success_prob() for allowed values, default)
+        :returns mag: star magnitude
     """
 
     def prob_minus_p_acq(mag):
@@ -928,13 +990,18 @@ def guide_count(mags, t_ccd, count_9th=False):
     of guide_counts to two decimal places (``8 * 0.0005 = 0.004``), but helps with
     minimization.
 
-    :param mags: float, array
+    Parameters
+    ----------
+    mags : float, array
         Star magnitude(s)
-    :param t_ccds: float, array
+    t_ccds : float, array
         CCD temperature(s)
-    :param count_9th: bool
+    count_9th : bool
         Return fractional count of 9th mag or brighter stars
-    :returns: float, fractional count
+
+    Returns
+    -------
+    float, fractional count
     """
     mags = np.atleast_1d(mags)
     mags, t_ccds = np.broadcast_arrays(mags, t_ccd)
@@ -972,12 +1039,20 @@ def t_ccd_warm_limit_for_guide(
     cold end the result may be below ``min_n_acq``, in which case the star catalog
     may be rejected.
 
-    :param mags: list of star ACA mags
-    :param min_guide_count: float minimum fractional guide count
-    :param warm_t_ccd: warmest CCD temperature to consider (default=-5 C)
-    :param cold_t_ccd: coldest CCD temperature to consider (default=-16 C)
-
-    :returns: t_ccd
+    Parameters
+    ----------
+    mags
+        list of star ACA mags
+    min_guide_count
+        float minimum fractional guide count
+    warm_t_ccd
+        warmest CCD temperature to consider (default=-5 C)
+    cold_t_ccd
+        coldest CCD temperature to consider (default=-16 C)
+
+    Returns
+    -------
+    t_ccd
     """
     if guide_count(mags, warm_t_ccd) >= min_guide_count:
         return warm_t_ccd
@@ -1012,12 +1087,20 @@ def binom_ppf(k, n, conf, n_sample=1000):
       >>> binom_ppf(4, 5, [0.17, 0.84])
       array([ 0.55463945,  0.87748177])
 
-    :param k: int, number of successes (0 < k <= n)
-    :param n: int, number of trials
-    :param conf: float or array of floats, percent point values
-    :param n_sample: number of PMF samples for interpolation
-
-    :return: percent point function values corresponding to ``conf``
+    Parameters
+    ----------
+    k
+        int, number of successes (0 < k <= n)
+    n
+        int, number of trials
+    conf
+        float or array of floats, percent point values
+    n_sample
+        number of PMF samples for interpolation
+
+    Returns
+    -------
+    percent point function values corresponding to ``conf``
     """
     ps = np.linspace(0, 1, n_sample)  # Probability values
     pmfs = scipy.stats.binom.pmf(k=k, n=n, p=ps)
diff --git a/chandra_aca/transform.py b/chandra_aca/transform.py
index ec72a34..e05277b 100644
--- a/chandra_aca/transform.py
+++ b/chandra_aca/transform.py
@@ -141,14 +141,24 @@ def pixels_to_yagzag(
     that is the default.  However, it is generally more convenient when doing
     centroids and other manipulations to use the center.
 
-    :param row: ACA pixel row (single value, list, or 1-d numpy array)
-    :param col: ACA pixel column (single value, list, or 1-d numpy array)
-    :param allow_bad: boolean switch.  If true, method will not throw errors
-                         if the row/col values are nominally off the ACA CCD.
-    :param flight: Use flight EEPROM coefficients instead of default ground values.
-    :param t_aca: ACA temperature (degC) for use with flight (default=20C)
-    :param pix_zero_loc: row/col coords are integral at 'edge' or 'center'
-    :rtype: (yang, zang) each vector of the same length as row/col
+    Parameters
+    ----------
+    row
+        ACA pixel row (single value, list, or 1-d numpy array)
+    col
+        ACA pixel column (single value, list, or 1-d numpy array)
+    allow_bad : boolean switch.  If true, method will not throw errors
+        if the row/col values are nominally off the ACA CCD.
+    flight
+        Use flight EEPROM coefficients instead of default ground values.
+    t_aca
+        ACA temperature (degC) for use with flight (default=20C)
+    pix_zero_loc
+        row/col coords are integral at 'edge' or 'center'
+
+    Returns
+    -------
+    (yang, zang) each vector of the same length as row/col
     """
     row = np.asarray(row)
     col = np.asarray(col)
@@ -184,12 +194,20 @@ def yagzag_to_pixels(yang, zang, allow_bad=False, pix_zero_loc="edge"):
     that is the default.  However, it is generally more convenient when doing
     centroids and other manipulations to use ``pix_zero_loc='center'``.
 
-    :param yang: ACA y-angle (single value, list, or 1-d numpy array)
-    :param zang: ACA z-angle (single value, list, or 1-d numpy array)
-    :param allow_bad: boolean switch.  If true, method will not throw errors
-                         if the resulting row/col values are nominally off the ACA CCD.
-    :param pix_zero_loc: row/col coords are integral at 'edge' or 'center'
-    :rtype: (row, col) each vector of the same length as row/col
+    Parameters
+    ----------
+    yang
+        ACA y-angle (single value, list, or 1-d numpy array)
+    zang
+        ACA z-angle (single value, list, or 1-d numpy array)
+    allow_bad : boolean switch.  If true, method will not throw errors
+        if the resulting row/col values are nominally off the ACA CCD.
+    pix_zero_loc
+        row/col coords are integral at 'edge' or 'center'
+
+    Returns
+    -------
+    (row, col) each vector of the same length as row/col
     """
     yang = np.asarray(yang)
     zang = np.asarray(zang)
@@ -276,9 +294,16 @@ def radec_to_eci(ra, dec):
     arrays of length N in which case the output ``ECI`` will be an array with
     shape (N, 3). The N dimension can actually be any multidimensional shape.
 
-    :param ra: Right Ascension (degrees)
-    :param dec: Declination (degrees)
-    :returns: numpy array ECI (3-vector or N x 3 array)
+    Parameters
+    ----------
+    ra
+        Right Ascension (degrees)
+    dec
+        Declination (degrees)
+
+    Returns
+    -------
+    numpy array ECI (3-vector or N x 3 array)
     """
     r = np.deg2rad(ra)
     d = np.deg2rad(dec)
@@ -298,8 +323,14 @@ def eci_to_radec(eci):
     the output RA, Dec will be arrays of shape N. The N dimension can
     actually be any multidimensional shape.
 
-    :param eci: ECI as 3-vector or (N, 3) array
-    :rtype: scalar or array ra, dec (degrees)
+    Parameters
+    ----------
+    eci
+        ECI as 3-vector or (N, 3) array
+
+    Returns
+    -------
+    scalar or array ra, dec (degrees)
     """
     eci = np.asarray(eci)
     if eci.shape[-1] != 3:
@@ -335,8 +366,14 @@ def mag_to_count_rate(mag):
     columns mag0=10.32 and cnt_rate_mag0=5263.0.  To convert to raw telemetered
     counts (ADU), divide e- by the gain of 5.0 e-/ADU.
 
-    :param mag: star magnitude in ACA mag
-    :returns: count rate (e-/sec)
+    Parameters
+    ----------
+    mag
+        star magnitude in ACA mag
+
+    Returns
+    -------
+    count rate (e-/sec)
     """
     count_rate = ACA_CNT_RATE_MAG0 * 10.0 ** ((ACA_MAG0 - mag) / 2.5)
     return count_rate
@@ -349,8 +386,14 @@ def count_rate_to_mag(count_rate):
     Based on $CALDB/data/chandra/pcad/ccd_char/acapD1998-12-14ccdN0002.fits
     columns mag0=10.32 and cnt_rate_mag0=5263.0.
 
-    :param count_rate: count rate (e-/sec)
-    :returns: magnitude (ACA mag)
+    Parameters
+    ----------
+    count_rate
+        count rate (e-/sec)
+
+    Returns
+    -------
+    magnitude (ACA mag)
     """
     mag = ACA_MAG0 - 2.5 * np.log10(count_rate / ACA_CNT_RATE_MAG0)
     return mag
@@ -376,16 +419,21 @@ def snr_mag_for_t_ccd(t_ccd, ref_mag, ref_t_ccd, scale_4c=None):
       counts(mag) = ACA_CNT_RATE_MAG0 * 10.0 ** ((ACA_MAG0 - mag) / 2.5)
       ref_mag = (t_ccd - ref_t_ccd) * np.log10(scale_4c) / 1.6
 
-    :param t_ccd: float, array
+    Parameters
+    ----------
+    t_ccd : float, array
         CCD temperature (degC)
-    :param ref_mag: float, array
+    ref_mag : float, array
         Reference magnitude (mag)
-    :param ref_t_ccd: float, array
+    ref_t_ccd : float, array
         Reference CCD temperature (degC)
-    :param scale_4c: float
+    scale_4c : float
         Scale factor for a 4 degC change in CCD temperature (defaults to DARK_SCALE_4C)
-    :returns: float, array
-        Magnitude(s) with the same expected signal to noise as ref_mag at ref_t_ccd
+
+    Returns
+    -------
+    float, array
+    Magnitude(s) with the same expected signal to noise as ref_mag at ref_t_ccd
     """
     # Allow array inputs
     t_ccd, ref_mag, ref_t_ccd = np.broadcast_arrays(t_ccd, ref_mag, ref_t_ccd)
@@ -400,12 +448,20 @@ def calc_aca_from_targ(targ, y_off, z_off, si_align=None):
     Calculate the PCAD (ACA) pointing attitude from target attitude and
     Y,Z offset (per OR-list).
 
-    :param targ: science target from OR/Obscat (Quat-compliant object)
-    :param y_off: Y offset (deg, sign per OR-list convention)
-    :param z_off: Z offset (deg, sign per OR-list convention)
-    :param si_align: SI ALIGN matrix (default=ODB_SI_ALIGN)
-
-    :rtype: q_aca (Quat)
+    Parameters
+    ----------
+    targ
+        science target from OR/Obscat (Quat-compliant object)
+    y_off
+        Y offset (deg, sign per OR-list convention)
+    z_off
+        Z offset (deg, sign per OR-list convention)
+    si_align
+        SI ALIGN matrix (default=ODB_SI_ALIGN)
+
+    Returns
+    -------
+    q_aca (Quat)
     """
     if si_align is None:
         si_align = ODB_SI_ALIGN
@@ -423,12 +479,20 @@ def calc_targ_from_aca(aca, y_off, z_off, si_align=None):
     Calculate the target attitude from ACA (PCAD) pointing attitude and
     Y,Z offset (per OR-list).
 
-    :param aca: ACA (PCAD) pointing attitude (any Quat-compliant object)
-    :param y_off: Y offset (deg, sign per OR-list convention)
-    :param z_off: Z offset (deg, sign per OR-list convention)
-    :param si_align: SI ALIGN matrix
-
-    :rtype: q_targ (Quat)
+    Parameters
+    ----------
+    aca
+        ACA (PCAD) pointing attitude (any Quat-compliant object)
+    y_off
+        Y offset (deg, sign per OR-list convention)
+    z_off
+        Z offset (deg, sign per OR-list convention)
+    si_align
+        SI ALIGN matrix
+
+    Returns
+    -------
+    q_targ (Quat)
     """
     if si_align is None:
         si_align = ODB_SI_ALIGN
@@ -447,12 +511,20 @@ def calc_target_offsets(aca, ra_targ, dec_targ, si_align=None):
     Calculates required Y and Z offsets (deg) required from a target RA, Dec to
     arrive at the desired PCAD pointing ``aca``.
 
-    :param aca: PCAD attitude (any Quat-compliant initializer)
-    :param ra_targ: RA of science target from OR/Obscat
-    :param dec_targ: Dec of science target from OR/Obscat
-    :param si_align: SI ALIGN matrix (default=ODB)
-
-    :rtype: tuple (y_off, z_off) degrees
+    Parameters
+    ----------
+    aca
+        PCAD attitude (any Quat-compliant initializer)
+    ra_targ
+        RA of science target from OR/Obscat
+    dec_targ
+        Dec of science target from OR/Obscat
+    si_align
+        SI ALIGN matrix (default=ODB)
+
+    Returns
+    -------
+    tuple (y_off, z_off) degrees
     """
     if si_align is None:
         si_align = ODB_SI_ALIGN
@@ -480,11 +552,18 @@ def radec_to_yagzag(ra, dec, q_att):
     input ``ra`` and ``dec`` values can be 1-d arrays in which case the output
     ``yag`` and ``zag`` will be corresponding arrays of the same length.
 
-    :param ra: Right Ascension (degrees)
-    :param dec: Declination (degrees)
-    :param q_att: ACA pointing quaternion (Quat or Quat-compatible input)
-
-    :returns:  yag, zag (arcsec)
+    Parameters
+    ----------
+    ra
+        Right Ascension (degrees)
+    dec
+        Declination (degrees)
+    q_att
+        ACA pointing quaternion (Quat or Quat-compatible input)
+
+    Returns
+    -------
+    yag, zag (arcsec)
     """
     if not isinstance(q_att, Quat):
         q_att = Quat(q_att)
@@ -502,11 +581,18 @@ def yagzag_to_radec(yag, zag, q_att):
     input ``yag`` and ``zag`` values can be 1-d arrays in which case the output
     ``ra`` and ``dec`` will be corresponding arrays of the same length.
 
-    :param yag: ACA Y angle (arcsec)
-    :param zag: ACA Z angle (arcsec)
-    :param q_att: ACA pointing quaternion (Quat or Quat-compatible input)
-
-    :returns: ra, dec (degrees)
+    Parameters
+    ----------
+    yag
+        ACA Y angle (arcsec)
+    zag
+        ACA Z angle (arcsec)
+    q_att
+        ACA pointing quaternion (Quat or Quat-compatible input)
+
+    Returns
+    -------
+    ra, dec (degrees)
     """
     if not isinstance(q_att, Quat):
         q_att = Quat(q_att)
@@ -528,9 +614,16 @@ def normalize_vector(vec, ord=None):
     For an L x M x N input array, this normalizes over the N dimension
     using ``np.linalg.norm``.
 
-    :param vec: input vector or array of vectors
-    :param ord: ord parameter for np.linalg.norm (default=None => 2-norm)
-    :returns: normed array of vectors
+    Parameters
+    ----------
+    vec
+        input vector or array of vectors
+    ord
+        ord parameter for np.linalg.norm (default=None => 2-norm)
+
+    Returns
+    -------
+    normed array of vectors
     """
     norms = np.linalg.norm(vec, axis=-1, keepdims=True, ord=ord)
     return vec / norms
diff --git a/docs/conf.py b/docs/conf.py
index f6a84df..e42c87e 100644
--- a/docs/conf.py
+++ b/docs/conf.py
@@ -27,6 +27,13 @@
 # If your documentation needs a minimal Sphinx version, state it here.
 # needs_sphinx = '1.0'
 
+# Don't show summaries of the members in each class along with the
+# class' docstring
+numpydoc_show_class_members = False
+
+# Don't show type hints
+autodoc_typehints = "none"
+
 # Add any Sphinx extension module names here, as strings. They can be
 # extensions coming with Sphinx (named 'sphinx.ext.*') or your custom
 # ones.
@@ -34,6 +41,7 @@
     "sphinx.ext.autodoc",
     "sphinx.ext.doctest",
     "sphinx.ext.viewcode",
+    "numpydoc",
 ]
 
 # Add any paths that contain templates here, relative to this directory.
diff --git a/docs/dark_model.rst b/docs/dark_model.rst
new file mode 100644
index 0000000..965b71b
--- /dev/null
+++ b/docs/dark_model.rst
@@ -0,0 +1,5 @@
+chandra_aca.dark_model
+==========================
+
+.. automodule:: chandra_aca.dark_model
+   :members:
diff --git a/docs/index.rst b/docs/index.rst
index f5cf2ba..e1d653a 100644
--- a/docs/index.rst
+++ b/docs/index.rst
@@ -37,6 +37,14 @@ Centroid residuals
 
     centroid_resid.rst
 
+Dark model
+----------
+
+.. toctree::
+    :maxdepth: 2
+
+    dark_model.rst
+
 MAUDE Decom
 -----------