def power_law(nu, nu_ref_s, beta_s):
    """ Function to compute power law SED.
    Parameters
    ----------
    nu: float, or array_like(float)
        Frequency in GHz.
    beta_s: float
        Power law index in RJ units.
    Returns
    -------
    array_like(float)
        SED relative to reference frequency.
    """
    x = nu / nu_ref_s
    sed = x ** beta_s
    return sed