docs: def get_bn(n, mc, dl, F, e):
      """
      Compute b_n from Eq. 22 of Taylor et al. (2015).
      
      :param n: Harmonic number
      :param mc: Chirp mass of binary [Solar Mass]
      :param dl: Luminosity distance [Mpc]
      :param F: Orbital frequency of binary [Hz]
      :param e: Orbital Eccentricity
      
      :returns: b_n
      
      """
      
      # convert to seconds
      mc *= SOLAR2S
      dl *= MPC2S 
      
      omega = 2 * np.pi * F
      
      amp = n * mc**(5/3) * omega**(2/3) / dl
          
      ret = -amp * np.sqrt(1-e**2) *(ss.jn(n-2,n*e) - 2*ss.jn(n,n*e) +
                    ss.jn(n+2,n*e)) 
  
      return ret