>>> import astropy.units as u
>>> from sbpy.calib import Sun
>>> from sbpy.photometry import bandpass
>>> from sbpy.spectroscopy import SpectralGradient
>>> import sbpy.units as sbu

>>> S = SpectralGradient(18 * u.percent / sbu.hundred_nm, wave0=550 * u.nm)
>>> sun = Sun.from_builtin("calspec")
>>> comet = sun.redden(S)

>>> bp_g = bandpass("LSST g")
>>> bp_r = bandpass("LSST r")

>>> _, r = comet.observe_bandpass(bp_r, unit=u.ABmag)
>>> _, g = comet.observe_bandpass(bp_g, unit=u.ABmag)
>>> print("g-r =", g - r)
g-r = 0.7007308548908533 mag