def rho_NFW(r, rs, rhos):
    """
    The density profile [GeV/cm**3] of an NFW halo.

    Parameters
    ----------
    r : the distance from the center [kpc]
    rs : the NFW r_s parameter [kpc]
    rhos : the NFW rho_s parameter [GeV/cm**3]
    """
    res = rhos/(r/rs)/(1.+r/rs)**2
    return res