Skip to content

Move spex prism loader to specutils #104

@kelle

Description

@kelle

Move Spex loader to specutils. Related issue astropy/specutils#1259

def _identify_spex(filename):
    """
    Check whether the given file is a SpeX data product.
    """
    try:
        with fits.open(filename, memmap=False) as hdulist:
            return "spex" in hdulist[0].header["INSTRUME"].lower() and "irtf" in hdulist[0].header["TELESCOP"].lower()
    except Exception:  # pylint: disable=broad-except,
        return False


def identify_spex_prism(origin, *args, **kwargs):
    """
    Confirm this is a SpeX Prism FITS file.
    See FITS keyword reference at http://irtfweb.ifa.hawaii.edu/~spex/observer/
    Notes: GRAT has values of: ShortXD, Prism, LXD_long, LXD_short, SO_long, SO_short
    """
    is_spex = _identify_spex(args[0])
    if is_spex:
        with fits.open(args[0], memmap=False) as hdulist:
            return (
                isinstance(args[0], str)
                and os.path.splitext(args[0].lower())[1] == ".fits"
                and is_spex
                and ("lowres" in hdulist[0].header["GRAT"].lower() or "prism" in hdulist[0].header["GRAT"].lower())
            )
    else:
        return is_spex


@data_loader("Spex Prism", identifier=identify_spex_prism, extensions=["fits"], dtype=Spectrum)
def spex_prism_loader(filename, **kwargs):
    """Open a SpeX Prism file and convert it to a Spectrum1D object"""

    with fits.open(filename, **kwargs) as hdulist:
        header = hdulist[0].header

        tab = hdulist[0].data

        # Handle missing/incorrect units
        try:
            flux_unit = header["YUNITS"].replace("ergs", "erg ").strip()
            wave_unit = header["XUNITS"].replace("Microns", "um")
        except (KeyError, ValueError):
            # For now, assume some default units
            flux_unit = "erg"
            wave_unit = "um"

        wave, data = tab[0] * Unit(wave_unit), tab[1] * Unit(flux_unit)

        if tab.shape[0] == 3:
            uncertainty = StdDevUncertainty(tab[2])
        else:
            uncertainty = None

        meta = {"header": header}

    return Spectrum(flux=data, spectral_axis=wave, uncertainty=uncertainty, meta=meta)

Metadata

Metadata

Assignees

No one assigned

    Labels

    No labels
    No labels

    Type

    No type
    No fields configured for issues without a type.

    Projects

    No projects

    Milestone

    No milestone

    Relationships

    None yet

    Development

    No branches or pull requests

    Issue actions