diff --git a/.github/workflows/build-test.yml b/.github/workflows/build-test.yml index fd962efa..85020265 100644 --- a/.github/workflows/build-test.yml +++ b/.github/workflows/build-test.yml @@ -43,4 +43,4 @@ jobs: run: | conda activate potpyri pip install .[test] - pytest + pytest -m "not integration" diff --git a/.github/workflows/documentation.yml b/.github/workflows/documentation.yml index 96ae9cd8..a185a883 100644 --- a/.github/workflows/documentation.yml +++ b/.github/workflows/documentation.yml @@ -34,6 +34,7 @@ jobs: cd docs make html - name: Deploy to GitHub Pages + if: github.event_name == 'push' && github.ref == 'refs/heads/main' uses: peaceiris/actions-gh-pages@v4 with: publish_branch: gh-pages diff --git a/CONTRIBUTING.md b/CONTRIBUTING.md new file mode 100644 index 00000000..afa1e085 --- /dev/null +++ b/CONTRIBUTING.md @@ -0,0 +1,26 @@ +# Contributing to POTPyRI + +Thank you for your interest in contributing to POTPyRI. We welcome issues, suggestions, and pull requests. + +## Reporting issues + +If you encounter a bug, have a question about usage, or have a feature request: + +- **GitHub Issues:** Open an issue at [github.com/CIERA-Transients/POTPyRI/issues](https://github.com/CIERA-Transients/POTPyRI/issues). Please include a clear description, steps to reproduce (for bugs), and your environment (Python version, instrument, OS) where relevant. +- **Email:** You can also contact the developers at `ckilpatrick@northwestern.edu` for matters you prefer not to discuss in the issue tracker. + +## Contributing code or documentation + +1. **Fork the repository** and create a branch from the default branch. +2. **Make your changes** and add or update tests if applicable. Run the test suite with `pytest tests` (use `pytest tests -m "not integration"` for offline runs). +3. **Open a pull request** against the main repository. Describe your changes clearly and reference any related issues. +4. **Code style:** Follow the existing style in the codebase. The project uses standard Python packaging and type hints where appropriate. + +Instrument-specific changes (e.g. new instruments or header/sorting logic) should include a brief justification and, if possible, a note in the PR description on how the change was tested. + +## Support + +- For **usage questions, bugs, or feature requests,** please open a [GitHub issue](https://github.com/CIERA-Transients/POTPyRI/issues). +- For **other inquiries** (e.g. collaboration, adding an instrument), you can contact the developers at `ckilpatrick@northwestern.edu`. + +We aim to respond to issues and pull requests in a timely manner; please allow at least a week for non-urgent responses. diff --git a/README.md b/README.md index eea7c003..8b4d5b50 100755 --- a/README.md +++ b/README.md @@ -3,63 +3,91 @@ [![GitHub Actions](https://github.com/CIERA-Transients/POTPyRI/actions/workflows/build-test.yml/badge.svg)](https://github.com/CIERA-Transients/POTPyRI/actions/workflows/build-test.yml) [![Docs](https://img.shields.io/badge/docs-latest-blue)](https://ciera-transients.github.io/POTPyRI/) - -Data reduction pipeline for imaging from large aperture telescopes. - -## Usage - -If you use this code, please reference Paterson et al. (in prep.) +Data reduction pipeline for imaging from large aperture telescopes. Supports automatic file ingestion, pixel calibration, astrometry (astrometry.net + Gaia DR3), stacking, aperture/PSF photometry, and flux calibration. ## Installation -The recommended installation is via conda and pip. First create a conda environment: - -```conda create -n potpyri python=3.12 pip-tools``` +**Requirements:** Python 3.11 or later. -Then to install POTPyRI, run: +### From PyPI (recommended) +```bash +pip install potpyri ``` + +### With conda (recommended for managing Python and non-Python deps) + +```bash +conda create -n potpyri python=3.12 pip conda activate potpyri pip install potpyri ``` -### Non-Python Dependencies - -POTPyRI has two non-python dependencies: `astrometry.net` and `source extractor`. The package file will try to install these dependencies on `linux-64` and `osx-64` environments for which there are working conda repositories. If these do not install, manually run the commands with: +### From source (development) -``` -conda install conda-forge::astromatic-source-extractor -conda install conda-forge::astrometry +```bash +git clone https://github.com/CIERA-Transients/POTPyRI +cd POTPyRI +pip install -e . # editable install +pip install -e ".[test]" # add pytest, pytest-cov +pip install -e ".[docs]" # optional: Sphinx docs ``` -On newer Mac OS systems (`osx-arm64`), conda will not be able to install these packages and they must be built manually. The recommended installation method is with Homebrew (https://brew.sh/) via: +### Non-Python dependencies -``` -brew install sextractor -brew install astrometry-net -``` +The pipeline needs **astrometry.net** and **Source Extractor** (SExtractor): + +| Platform | Install method | +|----------|----------------| +| **Linux (x86_64)** | `conda install conda-forge::astromatic-source-extractor conda-forge::astrometry` | +| **macOS (Intel, osx-64)** | Same as Linux: both packages are available from conda-forge. | +| **macOS (Apple Silicon, arm64)** | **astromatic-source-extractor** is available from conda-forge. **astrometry** has no conda build for osx-arm64; install it (and optionally both) via Homebrew: `brew install sextractor astrometry-net`. | + +So: use the `environment.yml` (or the conda commands above) on Linux and Intel Mac. On Apple Silicon, either create the conda env for Python + source-extractor and add astrometry via Homebrew, or install both with `brew install sextractor astrometry-net` and use a separate conda env with only `python=3.12` and `pip` for POTPyRI. + +- **Astrometry.net index files** (~59 GB): after `astrometry.net` is on your PATH, run `download_anet_index` (or use the script in the package to fetch indices). + +## Usage -Finally, note that for `astrometry.net`, index files are required. There is a utility script in POTPyRI that will attempt to install the latest index files. Once you have successfully installed `astrometry.net` in your path, run: +Run the full pipeline from the command line: +```bash +main_pipeline INSTRUMENT DATA_PATH ``` -download_anet_index + +- **INSTRUMENT:** One of the supported instrument names (e.g. `GMOS`, `LRIS`, `BINOSPEC`). See **Supported instruments** below. +- **DATA_PATH:** Top-level directory containing your data (e.g. a `raw` folder with FITS files). + +Example: + +```bash +main_pipeline GMOS /path/to/my/run ``` +All options: `main_pipeline --help`. Processed data are written under `DATA_PATH/red/`; logs go to `red/log/`. + +**Citation:** If you use this code, please cite Paterson et al. (in prep.). + ## Testing -If you have downloaded POTPyRI directly from GitHub (or through git pull/clone), you can test the installation using pytest by running the following inside the POTPyRI folder: +From the project root (with test deps installed: `pip install -e ".[test]"`): -``` +```bash pytest tests ``` +- **Offline (no network):** + `pytest tests -m "not integration"` +- **Verbose:** + `pytest tests -v` + ## Required Disk Space The POTPyRI repository itself has a size of: ![Size](https://img.shields.io/github/repo-size/CIERA-Transients/POTPyRI). In addition, the required index files installed by `download_anet_index` require a large footprint, currently about 59 GB of disk space. -## Supported Instruments +## Supported instruments - BINOSPEC (MMT) - DEIMOS (Keck) @@ -71,56 +99,48 @@ In addition, the required index files installed by `download_anet_index` require - MMIRS (MMT) - MOSFIRE (Keck) -## Running +## Running and outputs -Once POTPyRI and external dependencies are installed, the basic syntax to run the pipeline is: +The pipeline runs non-interactively: file ingestion, pixel calibration, astrometry, stacking, photometry, and flux calibration. The log is printed to the terminal and saved under `data_path/red/log`. -```main_pipeline instrument data_path``` +Directory layout under `data_path`: -where `instrument` is the name of the instrument you wish to reduce data from (see **Supported Instruments**), and `data_path` is the full path of the data you wish you reduce. Both the `instrument` and `data_path` are required parameters. Optional parameters are discussed in **Pipeline parameters**. +- **raw** — raw data +- **bad** — files excluded from reduction +- **red** — all processed outputs +- **red/cals** — master bias, flat, dark +- **red/log** — log files +- **red/workspace** — intermediate per-frame images before stacking -The pipeline is designed to run with no user input, including file ingestion, pixel calibration with provided calibration files, astrometry, image stacking, photometry, and flux calibration. The pipeline will display the log to the terminal (also saved as a file under `red/log`) to provide feedback and indicate when an error is encountered. For more details on each of these sections, including how the pipeline performs each task, read the sections below. - -All processed data will be written out in the `data_path` under the `red` subdirectory. Raw data are stored in `raw` and unused files are stored in `bad`. Processed calibration files will be stored in `red/cals`, logs in `red/log`, and interstitial image files before stacking are stored in `log/workspace`. A description of each type of output file can be found in **Outputs**. +Details of output file names and contents are in **Outputs**. ## Pipeline parameters -In addition to the required parameters `instrument` and `data_path`, the following optional parameters can be passed to the pipeline (also viewable via `main_pipeline --help`): - -``` -options: - -h, --help show this help message and exit - --target TARGET Option to only reduce a specific target. String used here must be contained within the target name in file headers. Optional - parameter. - --proc PROC Option to specify file processing for data ingestion. - --include-bad, --incl-bad - Include files flagged as bad in the file list. - --no-redo-sort Do not resort files and generate new file list. - --file-list-name FILE_LIST_NAME - Change the name of the archive file list. - --phot-sn-min PHOT_SN_MIN - Minimum signal-to-noise to try in photometry loop. - --phot-sn-max PHOT_SN_MAX - Maximum signal-to-noise to try in photometry loop. - --fwhm-init FWHM_INIT - Initial FWHM (in pixels) for photometry loop. - --skip-skysub Do not perform sky subtraction during image processing. - --fieldcenter FIELDCENTER FIELDCENTER - Align the output images to this center coordinate. This is useful for difference imaging where the frames need to be a common size, - pixel scale, and set of coordinates. - --out-size OUT_SIZE Output image size (image will be SIZE x SIZE pixels). - --stages STAGES [STAGES ...] - Stages to execute if running the pipeline in a modular way. - --skip-flatten Tell the pipeline to skip flattening. - --skip-cr Tell the pipeline to skip cosmic ray detection. - --skip-gaia Tell the pipeline to skip Gaia alignment during WCS. -``` +Required: `instrument`, `data_path`. All other options are optional. Run `main_pipeline --help` for the full list. + +| Option | Default | Description | +|--------|---------|-------------| +| `--target TARGET` | — | Only reduce targets whose name contains this string. | +| `--proc PROC` | — | File processing mode for data ingestion. | +| `--include-bad` / `--incl-bad` | off | Include files flagged as bad in the file list. | +| `--no-redo-sort` | off | Do not re-sort files or regenerate the file list. | +| `--file-list-name NAME` | `file_list.txt` | Name of the generated file list. | +| `--phot-sn-min N` | 3.0 | Minimum S/N to try in the photometry loop. | +| `--phot-sn-max N` | 20.0 | Maximum S/N to try in the photometry loop. | +| `--fwhm-init N` | 5.0 | Initial FWHM (pixels) for the photometry loop. | +| `--skip-skysub` | off | Skip sky subtraction during image processing. | +| `--fieldcenter RA DEC` | — | Align outputs to this center (e.g. for difference imaging). | +| `--out-size N` | — | Output image size (N×N pixels). | +| `--skip-flatten` | off | Skip flat-field correction. | +| `--skip-cr` | off | Skip cosmic-ray detection. | +| `--skip-gaia` | off | Skip Gaia alignment in the WCS step. | +| `--keep-all-astro` | off | Keep all images regardless of astrometric dispersion (default: drop high-dispersion frames). | ## Outputs All outputs from the pipeline are written out to the `red` folder and various subdirectories. Calibration files such as the master bias, flat and dark are saved using `mbias`, `mflat` and `mdark` prefixes, along with additional information such as the filter, amplifiers and binning in the file names. These calibration files are stored under `red/cals`. -All processed science files are renamed with following the basic format: `{object}.{filter}.{ut_date}.{amplifier}.{binning}.{unique_number}*.fits`. All data are processed together within a common `TargType`, defined as image frames with the same object, filter, amplifier setup, and binning. This variable is defined in the `file_list.txt` located in `data_path` and that is generated from `potpyri/sort_files.py` during the initial file processing performed by POTPyRI. A more detailed description of `file_list.txt` is provided below. +All processed science files are renamed following the basic format: `{object}.{filter}.{ut_date}.{amplifier}.{binning}.{unique_number}*.fits`. All data are processed together within a common `TargType`, defined as image frames with the same object, filter, amplifier setup, and binning. This variable is defined in the `file_list.txt` located in `data_path` and that is generated from `potpyri/primitives/sort_files.py` during the initial file processing performed by POTPyRI. A more detailed description of `file_list.txt` is provided below. The following table provides a basic description of the naming format, location, and brief description for each science output file that POTPyRI will produce: @@ -133,29 +153,29 @@ Filename Location Description *stk.fits red The stacked image data for the corresponding `TargType`, containing SCI (image data), MASK (mask), and ERR (error) extensions ``` -Once photometry is performed, additional extensions are added to each `*stk.fits` file containing fits-formatted tables with aperture phohotometry, PSF stars, and PSF photometry of identified sources in the image. Currently, only the aperture photometry table is used by `potpyri/absphot.py` for flux calibration. +Once photometry is performed, additional extensions are added to each `*stk.fits` file containing FITS-formatted tables with aperture photometry, PSF stars, and PSF photometry of identified sources in the image. Currently, only the aperture photometry table is used by `potpyri/primitives/absphot.py` for flux calibration. ## File list -The pipeline will sort through all the fits with the correct format for a given instrument (given by the **raw_format** function in the setting file) and create a file list with the file type, target name, exposure time, observation time, and instrument setup such as number of amps and binning. +The pipeline will sort through all FITS files with the correct format for a given instrument (given by the **raw_format** function in the settings file) and create a file list with the file type, target name, exposure time, observation time, and instrument setup such as number of amps and binning. -This process is designed to be automatic and account for common observation or archiving errors by checking various header keywords to classify each file. If you find your files are being misclassified for a POTPyRI-supported instrument, please contact us or submit a Github Issue with the instrument and specific error/problem indicated. +This process is designed to be automatic and account for common observation or archiving errors by checking various header keywords to classify each file. Classification logic lives in `potpyri/primitives/sort_files.py`. If you find your files are being misclassified for a POTPyRI-supported instrument, please contact us or open a GitHub issue with the instrument and specific error indicated. ## Image Calibration -Bias, dark, and flat-field images are defined by the corresponding keywords in each `potpyri/instruments/*.py` file. The methods for generating and applying each calibration frame are generic and defined entirely within the general instrument class in `potpyri/instruments/instrument.py`. Combined with instrument-specific differences in overscan, trimming, and static masking, we find that these methods provide good quality initial processed images with few bad pixels and generally Gaussian noise for the typical calibration frames that are available by instrument. +Bias, dark, and flat-field images are defined by the corresponding keywords in each instrument module under `potpyri/instruments/`. The methods for generating and applying each calibration frame are generic and defined in the base class in `potpyri/instruments/instrument.py`; calibration orchestration is in `potpyri/primitives/calibration.py`. Combined with instrument-specific differences in overscan, trimming, and static masking, we find that these methods provide good quality initial processed images with few bad pixels and generally Gaussian noise for the typical calibration frames that are available by instrument. -If you find that your initial processed images (see **Outputs**) are noisy, contain a large number of bad pixels, or contain other artifacts due to pixel-level processing, please contact us or submit a Github Issue with the instrument and specific error/problem indicated. +If you find that your initial processed images (see **Outputs**) are noisy, contain a large number of bad pixels, or contain other artifacts due to pixel-level processing, please contact us or open a GitHub issue with the instrument and specific error indicated. ## Aligning and Astrometry Images are aligned in a two-step process using `astrometry.net` followed by centroiding to Gaia DR3 astrometric standard stars. This process has been tested across all POTPyRI-supported instruments in a variety of instrument setups, read out modes, filters, and fields with varying densities of stars. We find that it generally results in 0.1-1 pixel dispersion in the alignment solution per image, resulting in good image alignment for stacking and preservation of the PSF shape. -If you find that you have good quality images that consistently fail image alignment and astrometry, have extremely large dispersions, or have poor WCS solution, please contact us or submit a Github Issue with the instrument and specific error/problem indicated. +If you find that you have good quality images that consistently fail image alignment and astrometry, have extremely large dispersions, or have a poor WCS solution, please contact us or open a GitHub issue with the instrument and specific error indicated. ## Image Masking and Stacking -POTPyRI has been significantly updated to implement optimal masking for satellite trails, cosmic rays, bad pixels that can be statically masked or are introduced by the pixel-level calibration, and saturation/non-linear effects. These pixels are tracked throughout the pipeline from calibration through image stacking and are stored within the `*.fits`, `*_data.fits`, and `*stk.fits` using a bitwise image mask. In general, the following schema is used by `potpyri/image_procs.py` to flag bad pixels: +POTPyRI has been significantly updated to implement optimal masking for satellite trails, cosmic rays, bad pixels that can be statically masked or are introduced by the pixel-level calibration, and saturation/non-linear effects. These pixels are tracked throughout the pipeline from calibration through image stacking and are stored within the `*.fits`, `*_data.fits`, and `*stk.fits` using a bitwise image mask. In general, the following schema is used by `potpyri/primitives/image_procs.py` to flag bad pixels: ``` Bit Value @@ -165,46 +185,60 @@ Bit Value 8 Neighbor of a bad pixel set within a distance determined by "grow" parameter ``` -These pixels will be ignored during image stacking, and only pixels that are not flagged by these parameters will be used to generate the final stack. +These pixels are ignored during image stacking; only pixels not flagged are used to generate the final stack. -In addition, error/uncertainty images are generated for each frame, accounting for the expected read noise, Poisson noise, and an additional empirical noise term defined by the standard deviation in the sky background. These error images set the `weight` term for each input frame when performing stacking. +Error/uncertainty images are generated for each frame (read noise, Poisson noise, and empirical sky noise). These error images set the `weight` term for each input frame when stacking. -Final image stacking is performed by `ccdproc.combine` with the individual image data, mask, and error frames for each science image. The stacking method is generally `median`, but can be changed via the **stack_method** value in each `potpyri/instruments/*.py` parameter file. Images are scaled by the exposure time within each header to account for variable depth between frames. +Stacking is performed by `ccdproc.combine` in `potpyri/primitives/image_procs.py` with the individual image data, mask, and error frames for each science image. The stacking method is generally `median`, but can be changed via the **stack_method** value in each `potpyri/instruments/*.py` parameter file. Images are scaled by the exposure time within each header to account for variable depth between frames. ## Automatic Aperture and PSF photometry -The pipeline will perform both automatic aperture and PSF photometry of sources in the stacked image using `potpyri/photometry.py`. Firstly, the pipeline will detect sources within the image and determine statistics within a fiducial aperture radius using `photutils.aperture.ApertureStats`. This initial table of aperture photometry is saved within the `*stk.fits` image as the extension `APPPHOT`. +The pipeline performs both automatic aperture and PSF photometry of sources in the stacked image using `potpyri/primitives/photometry.py`. Firstly, the pipeline will detect sources within the image and determine statistics within a fiducial aperture radius using `photutils.aperture.ApertureStats`. This initial table of aperture photometry is saved within the `*stk.fits` image as the extension `APPPHOT`. -Next, based on cuts on roundness, FWHM, and signal-to-noise, the pipeline will define a list of bright stars with which it calculates an empirical PSF model. The final list of PSF stars is saved in the `*stk.fits` file as the `PSFSTARS` extension. The PSF itself will be generated from the extracted data around these PSF stars using `photutils.psf.EPSFBuilder`. A stamp of the final effective PSF is saved in the `*stk.fits` file as the `PSF` extension. A final FWHM is empirically calculated from the effective PSF model using a `astropy.modeling.functional_models.Moffat2D` defined profile and saved to the `*stk.fits` header. +Next, based on cuts on roundness, FWHM, and signal-to-noise, the pipeline will define a list of bright stars with which it calculates an empirical PSF model. The final list of PSF stars is saved in the `*stk.fits` file as the `PSFSTARS` extension. The PSF itself will be generated from the extracted data around these PSF stars using `photutils.psf.EPSFBuilder`. A stamp of the final effective PSF is saved in the `*stk.fits` file as the `PSF` extension. A final FWHM is empirically calculated from the effective PSF model using an `astropy.modeling.functional_models.Moffat2D` profile and saved to the `*stk.fits` header. Once the pipeline has determined the PSF, it will then calculate PSF photometry for all the originally extracted sources with `photutils.psf.PSFPhotometry` and write them to the `PSFPHOT` extension of the `*stk.fits` file. In addition, as a data quality check a residual image is saved to the `*stk.fits` file showing the original science data with each of these sources subtracted using the PSF model as the `RESIDUAL` extension. ## Flux Calibration -After the PSF photometry has been calculated, the pipeline will then calculate a zero point based aperture photometry by downloading a catalog of standard stars. Currently supported catalogs are Pan-STARRS, SDSS, SkyMapper, and 2MASS, covering most optical and infrared bands in the north and south. +After the PSF photometry has been calculated, the pipeline will then calculate a zero point based on aperture photometry by downloading a catalog of standard stars. Currently supported catalogs are Pan-STARRS, SDSS, SkyMapper, and 2MASS, covering most optical and infrared bands in the north and south. The zero point and associated uncertainty are saved in the `*stk.fits` files as `ZPTMAG` and `ZPTMUCER`, respectively. In addition, an approximate limiting magnitude is calculated for each `*stk.fits` file using the `FWHM` and `SKYSIG` (the approximate standard deviation in sky background per pixel) header keywords and saved to the header as `M3SIGMA`, `M5SIGMA`, and `M10SIGMA` (3-, 5-, and 10-sigma limiting magnitudes). ## Quality checks -The pipeline will automatically create a log file in the `red/log` folder named after the instrument, followed by the date and time at the start of the run (UTC). A lot of information is written to the log file about the various steps step in the pipeline and the reduction process. Since this log file can also contain the reduction of multiple targets, it is suggested to make a personal log for each target with key information that can be used for the data analysis. Before, during and after the pipeline the has run, there are some key aspect you should check in terms of quality control that should be noted in the target log. These steps are list and explained below: +The pipeline will automatically create a log file in the `red/log` folder named after the instrument, followed by the date and time at the start of the run (UTC). The log file records each step in the pipeline and the reduction process. Because it can cover multiple targets, we recommend keeping a separate log per target with key information for your analysis. Before, during, and after the pipeline has run, perform the following quality-control checks and note them in the target log: ### Data Processing Checks - Check the file list: The first thing the pipeline will do is to create a file list containing information about the files. Files can be misclassified, which can be corrected via direct editing of the file list and reprocessing with `--no-redo-sort`, direct editing of your file headers to correct the misclassification, or submitting an issue to debug POTPyRI-related classification issues. -- Check the stacks: After the pipeline has reduced and stacked the files for a target, it is recommended that you check your final stacked images to assess data quality. Check the fits image named after the target name. If the image is of poor quality with very few stars, WCS and flux calibration will be affected. If you find your stack has a large number of image artifacts, please contact us or submit a Github Issue with the instrument and specific error/problem indicated. +- Check the stacks: After the pipeline has reduced and stacked the files for a target, it is recommended that you check your final stacked images to assess data quality. Check the FITS image named after the target name. If the image is of poor quality with very few stars, WCS and flux calibration will be affected. If you find your stack has a large number of image artifacts, please contact us or open a GitHub issue with the instrument and specific error indicated. -- Check WCS solution: In additional to individual image WCS solutions, the pipeline will calculate an automatic WCS solution on the final stack. The rms on the astromtery is typically $<$1 pixel and the stars aligned well near your target as well as across the image. Check the image alignment quality via the `RADISP` and `DEDISP` keywords in each image header. If the WCS is off by a large number of pixels, the zero point will not be calculated correctly and the reduction should be aborted. In this case, please contact us or submit a Github Issue with the instrument and specific error/problem indicated. +- Check WCS solution: In addition to individual image WCS solutions, the pipeline calculates an automatic WCS solution on the final stack. The rms on the astrometry is typically $<$1 pixel, and the stars are aligned well near your target and across the image. Check the image alignment quality via the `RADISP` and `DEDISP` keywords in each image header. If the WCS is off by a large number of pixels, the zero point will not be calculated correctly and the reduction should be aborted. In that case, please contact us or open a GitHub issue with the instrument and specific error indicated. -- Check the PSF: The PSF calculated by the pipeline will determine if the PSF photometry, as well as the zero point can be trusted. The PSF is saved in the `*.epsf.png` file, and will usually look like a 2D Gaussian. The pipeline will also report the x and y sigma of a 2D Gaussian fitted to the PSF. If the PSF looks good and the x and y sigma values reeported are similar, the PSF photometry can be trusted. If the PSF doesn't look Gaussian-like or the reported x and y sigma values are very different or negative, then there was an issue determining the PSF. If there are a lot of good stars in the final stack that could have been used to determine the PSF, please contact us or submit a Github Issue with the instrument and specific error/problem indicated. +- Check the PSF: The PSF calculated by the pipeline will determine whether the PSF photometry and the zero point can be trusted. The PSF is saved in the `*.epsf.png` file, and will usually look like a 2D Gaussian. The pipeline will also report the x and y sigma of a 2D Gaussian fitted to the PSF. If the PSF looks good and the x and y sigma values reported are similar, the PSF photometry can be trusted. If the PSF doesn't look Gaussian-like or the reported x and y sigma values are very different or negative, there was an issue determining the PSF. If there are many good stars in the final stack that could have been used to determine the PSF, please contact us or open a GitHub issue with the instrument and specific error indicated. -- Check the zero point calculation: If the PSF photometry is reliable, the zero point calculation should be good if a large number of stars is used to calculate it. In general, zero points calculated with $>$10 stars will provide relatively small (10$--$20 mmag level or smaller) zero point uncertainty. If there are more stars in the field that the pipeline should have used, or the pipeline reported that no stars were found but there is coverage in the corresponding catalog, please contact us or submit a Github Issue with the instrument and specific error/problem indicated. +- Check the zero point calculation: If the PSF photometry is reliable, the zero point calculation should be good if a large number of stars is used to calculate it. In general, zero points calculated with $>$10 stars will provide relatively small (10$--$20 mmag level or smaller) zero point uncertainty. If there are more stars in the field that the pipeline should have used, or the pipeline reported that no stars were found but there is coverage in the corresponding catalog, please contact us or open a GitHub issue with the instrument and specific error indicated. ### Data over multiple nights -POTPyRI can stack data from multiple nights, although it will not use unique calibration frames for each stacked image. To stack multiple nights of night, move all files to a single `raw` directory and run the pipeline. +POTPyRI can stack data from multiple nights, although it will not use unique calibration frames for each stacked image. To stack multiple nights, move all files to a single `raw` directory and run the pipeline. + +## Documentation + +Full API documentation (generated from docstrings) is available at [https://CIERA-Transients.github.io/POTPyRI/](https://CIERA-Transients.github.io/POTPyRI/). To build the docs locally, install with `pip install -e ".[docs]"` and run `cd docs && make html` (or use the Sphinx `make.bat` on Windows). + +## Contributing and support + +We welcome contributions and feedback. See [CONTRIBUTING.md](CONTRIBUTING.md) for: + +- How to **report issues** or request features (GitHub Issues) +- How to **contribute** code or documentation (pull requests) +- Where to **seek support** (issues or developer contact) + +If you encounter an error in the pipeline, have a special data setup that does not run through the pipeline, or wish to add an instrument, please open an issue on [GitHub](https://github.com/CIERA-Transients/POTPyRI/issues) or contact the developers at `ckilpatrick@northwestern.edu`. -## Issues and error +## License -Should you encounter an error in the pipeline, have a special data setup that doesn't run through the pipeline, or wish to add an additional instrument, please contact the developers at `ckilpatrick@northwestern.edu` or open an Issue on GitHub. +This project is licensed under the GPL-3.0 license; see [LICENSE.txt](LICENSE.txt) for the full terms. diff --git a/docs/DEPLOY.md b/docs/DEPLOY.md new file mode 100644 index 00000000..e99a4ad5 --- /dev/null +++ b/docs/DEPLOY.md @@ -0,0 +1,37 @@ +# Deploying documentation to https://ciera-transients.github.io/POTPyRI/ + +The docs are built with Sphinx and deployed to GitHub Pages via the `gh-pages` branch. The URL is a **project site**: `https://.github.io//` → `https://ciera-transients.github.io/POTPyRI/`. + +## One-time setup (repo maintainers) + +1. **Enable GitHub Pages** + - On GitHub: **Settings** → **Pages** (left sidebar). + - Under **Build and deployment**: + - **Source**: Deploy from a branch. + - **Branch**: choose `gh-pages`, folder **/ (root)**. + - Save. + +2. **First deployment** + - The workflow [`.github/workflows/documentation.yml`](../.github/workflows/documentation.yml) runs on every **push to `main`**. + - It builds the docs (`pip install .[docs]`, `cd docs && make html`) and pushes `docs/build/html/` to the `gh-pages` branch. + - After the first successful run, the site is available at https://ciera-transients.github.io/POTPyRI/ (can take a few minutes). + +3. **Pull requests** + - The same workflow runs on pull requests targeting `main`, but the **Deploy** step is skipped (deploy only runs on push to `main`). So PRs validate that the docs build without updating the live site. + +## Local build + +To build and view the docs locally: + +```bash +pip install .[docs] +cd docs && make html +open build/html/index.html # or open build/html/index.html on Linux +``` + +## Configuration + +- **Sphinx config**: `docs/source/conf.py` (e.g. `html_baseurl`, theme, extensions). +- **Base URL**: `html_baseurl = 'https://CIERA-Transients.github.io/POTPyRI/'` is already set so links and assets resolve correctly on GitHub Pages. + +If the site returns 404 after enabling Pages, wait a few minutes and ensure the `gh-pages` branch exists and has content (run the workflow once by pushing to `main`). diff --git a/docs/source/api/generated/potpyri.instruments.BINOSPEC.rst b/docs/source/api/generated/potpyri.instruments.BINOSPEC.rst new file mode 100644 index 00000000..deb63187 --- /dev/null +++ b/docs/source/api/generated/potpyri.instruments.BINOSPEC.rst @@ -0,0 +1,12 @@ +potpyri.instruments.BINOSPEC +============================ + +.. automodule:: potpyri.instruments.BINOSPEC + + + .. rubric:: Classes + + .. autosummary:: + + BINOSPEC + \ No newline at end of file diff --git a/docs/source/api/generated/potpyri.instruments.DEIMOS.rst b/docs/source/api/generated/potpyri.instruments.DEIMOS.rst new file mode 100644 index 00000000..27a7de11 --- /dev/null +++ b/docs/source/api/generated/potpyri.instruments.DEIMOS.rst @@ -0,0 +1,12 @@ +potpyri.instruments.DEIMOS +========================== + +.. automodule:: potpyri.instruments.DEIMOS + + + .. rubric:: Classes + + .. autosummary:: + + DEIMOS + \ No newline at end of file diff --git a/docs/source/api/generated/potpyri.instruments.F2.rst b/docs/source/api/generated/potpyri.instruments.F2.rst new file mode 100644 index 00000000..f10245e5 --- /dev/null +++ b/docs/source/api/generated/potpyri.instruments.F2.rst @@ -0,0 +1,12 @@ +potpyri.instruments.F2 +====================== + +.. automodule:: potpyri.instruments.F2 + + + .. rubric:: Classes + + .. autosummary:: + + F2 + \ No newline at end of file diff --git a/docs/source/api/generated/potpyri.instruments.FOURSTAR.rst b/docs/source/api/generated/potpyri.instruments.FOURSTAR.rst new file mode 100644 index 00000000..c7db7637 --- /dev/null +++ b/docs/source/api/generated/potpyri.instruments.FOURSTAR.rst @@ -0,0 +1,12 @@ +potpyri.instruments.FOURSTAR +============================ + +.. automodule:: potpyri.instruments.FOURSTAR + + + .. rubric:: Classes + + .. autosummary:: + + FOURSTAR + \ No newline at end of file diff --git a/docs/source/api/generated/potpyri.instruments.GMOS.rst b/docs/source/api/generated/potpyri.instruments.GMOS.rst new file mode 100644 index 00000000..16fb68f7 --- /dev/null +++ b/docs/source/api/generated/potpyri.instruments.GMOS.rst @@ -0,0 +1,12 @@ +potpyri.instruments.GMOS +======================== + +.. automodule:: potpyri.instruments.GMOS + + + .. rubric:: Classes + + .. autosummary:: + + GMOS + \ No newline at end of file diff --git a/docs/source/api/generated/potpyri.instruments.IMACS.rst b/docs/source/api/generated/potpyri.instruments.IMACS.rst new file mode 100644 index 00000000..5b865994 --- /dev/null +++ b/docs/source/api/generated/potpyri.instruments.IMACS.rst @@ -0,0 +1,12 @@ +potpyri.instruments.IMACS +========================= + +.. automodule:: potpyri.instruments.IMACS + + + .. rubric:: Classes + + .. autosummary:: + + IMACS + \ No newline at end of file diff --git a/docs/source/api/generated/potpyri.instruments.LRIS.rst b/docs/source/api/generated/potpyri.instruments.LRIS.rst new file mode 100644 index 00000000..fdde83b4 --- /dev/null +++ b/docs/source/api/generated/potpyri.instruments.LRIS.rst @@ -0,0 +1,12 @@ +potpyri.instruments.LRIS +======================== + +.. automodule:: potpyri.instruments.LRIS + + + .. rubric:: Classes + + .. autosummary:: + + LRIS + \ No newline at end of file diff --git a/docs/source/api/generated/potpyri.instruments.MMIRS.rst b/docs/source/api/generated/potpyri.instruments.MMIRS.rst new file mode 100644 index 00000000..4599b64f --- /dev/null +++ b/docs/source/api/generated/potpyri.instruments.MMIRS.rst @@ -0,0 +1,12 @@ +potpyri.instruments.MMIRS +========================= + +.. automodule:: potpyri.instruments.MMIRS + + + .. rubric:: Classes + + .. autosummary:: + + MMIRS + \ No newline at end of file diff --git a/docs/source/api/generated/potpyri.instruments.MOSFIRE.rst b/docs/source/api/generated/potpyri.instruments.MOSFIRE.rst new file mode 100644 index 00000000..5b158658 --- /dev/null +++ b/docs/source/api/generated/potpyri.instruments.MOSFIRE.rst @@ -0,0 +1,12 @@ +potpyri.instruments.MOSFIRE +=========================== + +.. automodule:: potpyri.instruments.MOSFIRE + + + .. rubric:: Classes + + .. autosummary:: + + MOSFIRE + \ No newline at end of file diff --git a/docs/source/api/generated/potpyri.instruments.instrument.rst b/docs/source/api/generated/potpyri.instruments.instrument.rst new file mode 100644 index 00000000..8db422b0 --- /dev/null +++ b/docs/source/api/generated/potpyri.instruments.instrument.rst @@ -0,0 +1,12 @@ +potpyri.instruments.instrument +============================== + +.. automodule:: potpyri.instruments.instrument + + + .. rubric:: Classes + + .. autosummary:: + + Instrument + \ No newline at end of file diff --git a/docs/source/api/generated/potpyri.primitives.absphot.rst b/docs/source/api/generated/potpyri.primitives.absphot.rst new file mode 100644 index 00000000..0b0e9d4e --- /dev/null +++ b/docs/source/api/generated/potpyri.primitives.absphot.rst @@ -0,0 +1,18 @@ +potpyri.primitives.absphot +========================== + +.. automodule:: potpyri.primitives.absphot + + + .. rubric:: Functions + + .. autosummary:: + + find_zeropoint + + .. rubric:: Classes + + .. autosummary:: + + absphot + \ No newline at end of file diff --git a/docs/source/api/generated/potpyri.primitives.calibration.rst b/docs/source/api/generated/potpyri.primitives.calibration.rst new file mode 100644 index 00000000..3cbb1bfd --- /dev/null +++ b/docs/source/api/generated/potpyri.primitives.calibration.rst @@ -0,0 +1,14 @@ +potpyri.primitives.calibration +============================== + +.. automodule:: potpyri.primitives.calibration + + + .. rubric:: Functions + + .. autosummary:: + + do_bias + do_dark + do_flat + \ No newline at end of file diff --git a/docs/source/api/generated/potpyri.primitives.image_procs.rst b/docs/source/api/generated/potpyri.primitives.image_procs.rst new file mode 100644 index 00000000..8fd7da6c --- /dev/null +++ b/docs/source/api/generated/potpyri.primitives.image_procs.rst @@ -0,0 +1,22 @@ +potpyri.primitives.image\_procs +=============================== + +.. automodule:: potpyri.primitives.image_procs + + + .. rubric:: Functions + + .. autosummary:: + + add_stack_mask + align_images + create_error + create_mask + detrend_stack + generate_wcs + get_fieldcenter + image_proc + mask_satellites + remove_pv_distortion + stack_data + \ No newline at end of file diff --git a/docs/source/api/generated/potpyri.primitives.photometry.rst b/docs/source/api/generated/potpyri.primitives.photometry.rst new file mode 100644 index 00000000..db181eaf --- /dev/null +++ b/docs/source/api/generated/potpyri.primitives.photometry.rst @@ -0,0 +1,21 @@ +potpyri.primitives.photometry +============================= + +.. automodule:: potpyri.primitives.photometry + + + .. rubric:: Functions + + .. autosummary:: + + create_conv + create_params + do_phot + extract_aperture_stats + extract_fwhm_from_epsf + generate_epsf + get_star_catalog + photloop + run_photometry + run_sextractor + \ No newline at end of file diff --git a/docs/source/api/generated/potpyri.primitives.solve_wcs.rst b/docs/source/api/generated/potpyri.primitives.solve_wcs.rst new file mode 100644 index 00000000..d8f4a49e --- /dev/null +++ b/docs/source/api/generated/potpyri.primitives.solve_wcs.rst @@ -0,0 +1,15 @@ +potpyri.primitives.solve\_wcs +============================= + +.. automodule:: potpyri.primitives.solve_wcs + + + .. rubric:: Functions + + .. autosummary:: + + align_to_gaia + clean_up_astrometry + get_gaia_catalog + solve_astrometry + \ No newline at end of file diff --git a/docs/source/api/generated/potpyri.primitives.sort_files.rst b/docs/source/api/generated/potpyri.primitives.sort_files.rst new file mode 100644 index 00000000..f35e48f1 --- /dev/null +++ b/docs/source/api/generated/potpyri.primitives.sort_files.rst @@ -0,0 +1,19 @@ +potpyri.primitives.sort\_files +============================== + +.. automodule:: potpyri.primitives.sort_files + + + .. rubric:: Functions + + .. autosummary:: + + handle_files + is_bad + is_bias + is_dark + is_flat + is_science + is_spec + sort_files + \ No newline at end of file diff --git a/docs/source/api/generated/potpyri.scripts.archives.download_anet_index.rst b/docs/source/api/generated/potpyri.scripts.archives.download_anet_index.rst new file mode 100644 index 00000000..ae13be4b --- /dev/null +++ b/docs/source/api/generated/potpyri.scripts.archives.download_anet_index.rst @@ -0,0 +1,15 @@ +potpyri.scripts.archives.download\_anet\_index +============================================== + +.. automodule:: potpyri.scripts.archives.download_anet_index + + + .. rubric:: Functions + + .. autosummary:: + + add_options + download_index_files + main + parse_astrometry_config + \ No newline at end of file diff --git a/docs/source/api/generated/potpyri.scripts.archives.download_gemini_data.rst b/docs/source/api/generated/potpyri.scripts.archives.download_gemini_data.rst new file mode 100644 index 00000000..32d71499 --- /dev/null +++ b/docs/source/api/generated/potpyri.scripts.archives.download_gemini_data.rst @@ -0,0 +1,21 @@ +potpyri.scripts.archives.download\_gemini\_data +=============================================== + +.. automodule:: potpyri.scripts.archives.download_gemini_data + + + .. rubric:: Functions + + .. autosummary:: + + add_options + download_data + download_file + get_associated_cals + get_full_outname + get_observation_data + load_cookie + main + mask_object_spectral_observation + unpack_tarfile + \ No newline at end of file diff --git a/docs/source/api/generated/potpyri.scripts.archives.download_keck_data.rst b/docs/source/api/generated/potpyri.scripts.archives.download_keck_data.rst new file mode 100644 index 00000000..9d22967c --- /dev/null +++ b/docs/source/api/generated/potpyri.scripts.archives.download_keck_data.rst @@ -0,0 +1,16 @@ +potpyri.scripts.archives.download\_keck\_data +============================================= + +.. automodule:: potpyri.scripts.archives.download_keck_data + + + .. rubric:: Functions + + .. autosummary:: + + add_options + date_query + download_data + login + main + \ No newline at end of file diff --git a/docs/source/api/generated/potpyri.scripts.main_pipeline.rst b/docs/source/api/generated/potpyri.scripts.main_pipeline.rst new file mode 100644 index 00000000..e64bdcc2 --- /dev/null +++ b/docs/source/api/generated/potpyri.scripts.main_pipeline.rst @@ -0,0 +1,13 @@ +potpyri.scripts.main\_pipeline +============================== + +.. automodule:: potpyri.scripts.main_pipeline + + + .. rubric:: Functions + + .. autosummary:: + + main + main_pipeline + \ No newline at end of file diff --git a/docs/source/api/generated/potpyri.utils.logger.rst b/docs/source/api/generated/potpyri.utils.logger.rst new file mode 100644 index 00000000..40f4f49f --- /dev/null +++ b/docs/source/api/generated/potpyri.utils.logger.rst @@ -0,0 +1,20 @@ +potpyri.utils.logger +==================== + +.. automodule:: potpyri.utils.logger + + + .. rubric:: Functions + + .. autosummary:: + + formatter_message + get_log + + .. rubric:: Classes + + .. autosummary:: + + ColoredFormatter + ColoredLogger + \ No newline at end of file diff --git a/docs/source/api/generated/potpyri.utils.options.rst b/docs/source/api/generated/potpyri.utils.options.rst new file mode 100644 index 00000000..533a8228 --- /dev/null +++ b/docs/source/api/generated/potpyri.utils.options.rst @@ -0,0 +1,16 @@ +potpyri.utils.options +===================== + +.. automodule:: potpyri.utils.options + + + .. rubric:: Functions + + .. autosummary:: + + add_options + add_paths + init_options + initialize_telescope + test_for_dependencies + \ No newline at end of file diff --git a/docs/source/api/generated/potpyri.utils.utilities.rst b/docs/source/api/generated/potpyri.utils.utilities.rst new file mode 100644 index 00000000..772010f2 --- /dev/null +++ b/docs/source/api/generated/potpyri.utils.utilities.rst @@ -0,0 +1,14 @@ +potpyri.utils.utilities +======================= + +.. automodule:: potpyri.utils.utilities + + + .. rubric:: Functions + + .. autosummary:: + + find_catalog + is_number + parse_coord + \ No newline at end of file diff --git a/docs/source/api/index.rst b/docs/source/api/index.rst new file mode 100644 index 00000000..023ffb9a --- /dev/null +++ b/docs/source/api/index.rst @@ -0,0 +1,35 @@ +.. _api-reference: + +API Reference +============= + +This section documents the public API of POTPyRI, generated from docstrings +in the source code. + +Instruments +----------- +.. toctree:: + :maxdepth: 2 + + instruments + +Primitives (calibration, image processing, WCS, photometry) +------------------------------------------------------------ +.. toctree:: + :maxdepth: 2 + + primitives + +Utilities (options, paths, logging) +----------------------------------- +.. toctree:: + :maxdepth: 2 + + utils + +Scripts (pipeline and archive downloaders) +------------------------------------------ +.. toctree:: + :maxdepth: 2 + + scripts diff --git a/docs/source/api/instruments.rst b/docs/source/api/instruments.rst new file mode 100644 index 00000000..593e55c0 --- /dev/null +++ b/docs/source/api/instruments.rst @@ -0,0 +1,39 @@ +.. _api-instruments: + +Instruments +=========== + +Instrument implementations and factory for POTPyRI. Each instrument module +defines a subclass of :class:`~potpyri.instruments.instrument.Instrument` with +detector keywords, calibration behavior, and file-sorting rules. + +Package overview +---------------- +.. automodule:: potpyri.instruments + :members: + :imported-members: + :no-undoc-members: + +Instrument base class +--------------------- +.. autosummary:: + :toctree: generated/ + :template: autosummary/module.rst + + potpyri.instruments.instrument + +Instrument modules +------------------ +.. autosummary:: + :toctree: generated/ + :template: autosummary/module.rst + + potpyri.instruments.BINOSPEC + potpyri.instruments.DEIMOS + potpyri.instruments.F2 + potpyri.instruments.FOURSTAR + potpyri.instruments.GMOS + potpyri.instruments.IMACS + potpyri.instruments.LRIS + potpyri.instruments.MMIRS + potpyri.instruments.MOSFIRE diff --git a/docs/source/api/primitives.rst b/docs/source/api/primitives.rst new file mode 100644 index 00000000..7e97b1da --- /dev/null +++ b/docs/source/api/primitives.rst @@ -0,0 +1,28 @@ +.. _api-primitives: + +Primitives +========== + +POTPyRI primitives handle calibration, image processing, file sorting, +WCS solving, photometry, and absolute photometry (zeropoint fitting). +All of these are used by the main pipeline. + +Package overview +---------------- +.. automodule:: potpyri.primitives + :members: + :imported-members: + :no-undoc-members: + +Modules +------- +.. autosummary:: + :toctree: generated/ + :template: autosummary/module.rst + + potpyri.primitives.absphot + potpyri.primitives.calibration + potpyri.primitives.image_procs + potpyri.primitives.photometry + potpyri.primitives.solve_wcs + potpyri.primitives.sort_files diff --git a/docs/source/api/scripts.rst b/docs/source/api/scripts.rst new file mode 100644 index 00000000..8cb00e4c --- /dev/null +++ b/docs/source/api/scripts.rst @@ -0,0 +1,24 @@ +.. _api-scripts: + +Scripts +======= + +Entry-point scripts: the main reduction pipeline and archive download utilities. + +Main pipeline +------------- +.. autosummary:: + :toctree: generated/ + :template: autosummary/module.rst + + potpyri.scripts.main_pipeline + +Archive downloaders +-------------------- +.. autosummary:: + :toctree: generated/ + :template: autosummary/module.rst + + potpyri.scripts.archives.download_anet_index + potpyri.scripts.archives.download_keck_data + potpyri.scripts.archives.download_gemini_data diff --git a/docs/source/api/utils.rst b/docs/source/api/utils.rst new file mode 100644 index 00000000..afec60bc --- /dev/null +++ b/docs/source/api/utils.rst @@ -0,0 +1,23 @@ +.. _api-utils: + +Utilities +========= + +Utilities for option parsing, path setup, logging, and general helpers. + +Package overview +---------------- +.. automodule:: potpyri.utils + :members: + :imported-members: + :no-undoc-members: + +Modules +------- +.. autosummary:: + :toctree: generated/ + :template: autosummary/module.rst + + potpyri.utils.logger + potpyri.utils.options + potpyri.utils.utilities diff --git a/docs/source/conf.py b/docs/source/conf.py index ee66f333..7bfec636 100644 --- a/docs/source/conf.py +++ b/docs/source/conf.py @@ -3,7 +3,8 @@ # For the full list of built-in configuration values, see the documentation: # https://www.sphinx-doc.org/en/master/usage/configuration.html -import sys, os +import sys +import os sys.path.insert(0, os.path.abspath("../..")) import potpyri @@ -22,6 +23,7 @@ extensions = [ 'sphinx.ext.doctest', 'sphinx.ext.autodoc', + 'sphinx.ext.napoleon', # NumPy/Google-style docstrings 'sphinx.ext.intersphinx', 'sphinx.ext.todo', 'sphinx.ext.coverage', @@ -33,11 +35,17 @@ 'nbsphinx', ] -autosummary_generate = True # Generate autosummary stubs automatically +# Napoleon settings for NumPy-style docstrings (used throughout potpyri) +napoleon_use_param = True +napoleon_use_rtype = True +napoleon_preprocess_types = True + +autosummary_generate = True autodoc_default_options = { - 'members': True, # Include all class members - 'undoc-members': False, # Include undocumented members - 'show-inheritance': True, # Show class inheritance + 'members': True, + 'undoc-members': False, + 'show-inheritance': True, + 'inherited-members': False, } templates_path = ['_templates'] diff --git a/docs/source/index.rst b/docs/source/index.rst index f8a8d7e8..72c22363 100644 --- a/docs/source/index.rst +++ b/docs/source/index.rst @@ -3,10 +3,22 @@ POTPyRI Documentation ********************* POTPyRI is a data reduction pipeline for imaging from large aperture telescopes. +It provides instrument-specific reduction workflows for bias, dark, and flat +calibration; image alignment and stacking; WCS solving with astrometry.net; +photometry with photutils; and flux calibration with astroquery. -Installation -============ +Contents +======== .. toctree:: + :maxdepth: 2 installation + api/index + +Indices and tables +================== + +* :ref:`genindex` +* :ref:`modindex` +* :ref:`search` diff --git a/docs/source/installation.rst b/docs/source/installation.rst index 97681c5c..f2da5374 100644 --- a/docs/source/installation.rst +++ b/docs/source/installation.rst @@ -1,24 +1,57 @@ Installation ============ +POTPyRI requires **Python 3.11 or later**. + +From PyPI +--------- + +.. code-block:: bash + + pip install potpyri + +From source (development) +-------------------------- + 1. Clone the repository: .. code-block:: bash git clone https://github.com/CIERA-Transients/POTPyRI + cd POTPyRI -2. Create and activate a virtual environment. Anaconda and Python 3.12 are recommended, but not required. +2. Create and activate a virtual environment (conda or venv). Python 3.12 is recommended. .. code-block:: bash - conda create -n "potpyri" python=3.12 + conda create -n potpyri python=3.12 conda activate potpyri -3. Install live version with pip: + # Or with venv: + # python3.12 -m venv .venv && source .venv/bin/activate # Linux/macOS + +3. Install in editable mode with optional dependencies: + +.. code-block:: bash + + pip install -e . # core package + pip install -e ".[test]" # add pytest, pytest-cov for testing + pip install -e ".[docs]" # add Sphinx and theme for building docs + +The ``-e`` (editable) flag uses the current directory as the live source, so changes are picked up without reinstalling. + +Non-Python dependencies +------------------------ + +The pipeline uses **astrometry.net** and **Source Extractor** (SExtractor). Install them via your system package manager or conda (see the main README). Index files for astrometry.net can be installed with the ``download_anet_index`` script after astrometry.net is on your path. + +Building the documentation +--------------------------- + +With ``.[docs]`` installed, from the project root: .. code-block:: bash - cd POTPyRI - pip install -e . + cd docs && make html -Specifically, the ``-e`` (editable) flag uses the specified directory as the "live" source code, and does not copy it to the usual site-packages directory. +Output is in ``docs/build/html/``. On Windows use ``make.bat html``. diff --git a/environment.yml b/environment.yml index d05b1e57..4d552f78 100644 --- a/environment.yml +++ b/environment.yml @@ -1,9 +1,22 @@ +# Conda environment for POTPyRI (see README for full instructions). +# Create with: conda env create -f environment.yml +# Then: conda activate potpyri && pip install potpyri +# Or from source: pip install -e . (and optionally -e ".[test]" or ".[docs]") +# +# Platform support (conda-forge): +# - linux-64, osx-64 (Intel Mac): both astromatic-source-extractor and astrometry are available. +# - osx-arm64 (Apple Silicon): astromatic-source-extractor has a native build; astrometry does NOT. +# If "conda env create" fails on Apple Silicon (no astrometry build), create from this file with +# the "astrometry" line removed, then run: brew install astrometry-net +# +# Astrometry.net index files (~59 GB) are installed separately via: download_anet_index + name: potpyri channels: - conda-forge - defaults dependencies: - python=3.12 - - pip-tools + - pip - astromatic-source-extractor - astrometry diff --git a/potpyri/__init__.py b/potpyri/__init__.py index e243cffc..e4d7f830 100755 --- a/potpyri/__init__.py +++ b/potpyri/__init__.py @@ -1,3 +1,9 @@ +"""POTPyRI: Data reduction pipeline for imaging from large aperture telescopes. + +This package provides instrument-specific reduction workflows for bias, dark, +flat calibration, image alignment, stacking, WCS solving, and photometry. +""" + from ._version import version as __version__ from .instruments import * diff --git a/potpyri/_version.py b/potpyri/_version.py index 22dafcc6..1c52449d 100644 --- a/potpyri/_version.py +++ b/potpyri/_version.py @@ -28,7 +28,7 @@ commit_id: COMMIT_ID __commit_id__: COMMIT_ID -__version__ = version = '1.0.14.dev19' -__version_tuple__ = version_tuple = (1, 0, 14, 'dev19') +__version__ = version = '1.0.15.dev10' +__version_tuple__ = version_tuple = (1, 0, 15, 'dev10') -__commit_id__ = commit_id = 'gcb4d08580' +__commit_id__ = commit_id = 'g590b43389' diff --git a/potpyri/instruments/BINOSPEC.py b/potpyri/instruments/BINOSPEC.py index 30b466ba..d90cdaa6 100755 --- a/potpyri/instruments/BINOSPEC.py +++ b/potpyri/instruments/BINOSPEC.py @@ -1,6 +1,5 @@ -# Parameter file for BINOSPEC/MMT - -__version__ = "2.0" # Last edited 09/21/2024 +"""BINOSPEC/MMT instrument configuration and reduction parameters.""" +from potpyri._version import __version__ import os import ccdproc @@ -16,6 +15,7 @@ from . import instrument class BINOSPEC(instrument.Instrument): + """BINOSPEC at MMT: optical imaging with flat calibration, no bias/dark.""" def __init__(self): diff --git a/potpyri/instruments/DEIMOS.py b/potpyri/instruments/DEIMOS.py index 88381014..3db0835b 100755 --- a/potpyri/instruments/DEIMOS.py +++ b/potpyri/instruments/DEIMOS.py @@ -1,6 +1,5 @@ -# Parameter file for DEIMOS/Keck - -__version__ = "2.0" # Last edited 09/21/2024 +"""DEIMOS/Keck instrument configuration and reduction parameters.""" +from potpyri._version import __version__ import os import ccdproc @@ -16,6 +15,7 @@ from . import instrument class DEIMOS(instrument.Instrument): + """DEIMOS at Keck: optical imaging with bias and flat calibration.""" def __init__(self): diff --git a/potpyri/instruments/F2.py b/potpyri/instruments/F2.py index 65e0bf5b..b8529eec 100755 --- a/potpyri/instruments/F2.py +++ b/potpyri/instruments/F2.py @@ -1,6 +1,5 @@ -# Parameter file for F2/Gemini-S - -__version__ = "1.0" # Last edited 11/23/2024 +"""F2/Gemini-S instrument configuration and reduction parameters.""" +from potpyri._version import __version__ import os import ccdproc @@ -15,6 +14,7 @@ from . import instrument class F2(instrument.Instrument): + """F2 at Gemini-South: NIR imaging with dark and flat calibration.""" def __init__(self): diff --git a/potpyri/instruments/FOURSTAR.py b/potpyri/instruments/FOURSTAR.py index 90e9faa8..6c7c9d75 100755 --- a/potpyri/instruments/FOURSTAR.py +++ b/potpyri/instruments/FOURSTAR.py @@ -1,6 +1,5 @@ -# Parameter file for Magellan/IMACS - -__version__ = "1.0" # Last edited 02/13/2025 +"""FourStar/Magellan instrument configuration and reduction parameters.""" +from potpyri._version import __version__ import os import ccdproc @@ -15,6 +14,7 @@ from . import instrument class FOURSTAR(instrument.Instrument): + """FourStar at Magellan: NIR imaging with flat calibration.""" def __init__(self): diff --git a/potpyri/instruments/GMOS.py b/potpyri/instruments/GMOS.py index 0231afb7..0a094dbb 100755 --- a/potpyri/instruments/GMOS.py +++ b/potpyri/instruments/GMOS.py @@ -1,6 +1,5 @@ -# Parameter file for GMOS/Gemini - -__version__ = "2.1" # Last edited 09/29/2024 +"""GMOS/Gemini instrument configuration and reduction parameters.""" +from potpyri._version import __version__ import os import ccdproc @@ -18,6 +17,7 @@ from . import instrument class GMOS(instrument.Instrument): + """GMOS at Gemini: optical imaging with bias and flat calibration.""" def __init__(self): diff --git a/potpyri/instruments/IMACS.py b/potpyri/instruments/IMACS.py index 11c52e41..1865e3f1 100755 --- a/potpyri/instruments/IMACS.py +++ b/potpyri/instruments/IMACS.py @@ -1,6 +1,5 @@ -# Parameter file for Magellan/IMACS - -__version__ = "1.0" # Last edited 02/13/2025 +"""IMACS/Magellan instrument configuration and reduction parameters.""" +from potpyri._version import __version__ import os import ccdproc @@ -16,6 +15,7 @@ from . import instrument class IMACS(instrument.Instrument): + """IMACS at Magellan: optical imaging with bias and flat calibration.""" def __init__(self): diff --git a/potpyri/instruments/LRIS.py b/potpyri/instruments/LRIS.py index 7a278340..1e191296 100755 --- a/potpyri/instruments/LRIS.py +++ b/potpyri/instruments/LRIS.py @@ -1,6 +1,5 @@ -# Parameter file for LRIS/Keck - -__version__ = "2.1" # Last edited 05/16/2025 +"""LRIS/Keck instrument configuration and reduction parameters.""" +from potpyri._version import __version__ import os import ccdproc @@ -17,6 +16,7 @@ from . import instrument class LRIS(instrument.Instrument): + """LRIS at Keck: optical imaging with bias and flat calibration.""" def __init__(self): diff --git a/potpyri/instruments/MMIRS.py b/potpyri/instruments/MMIRS.py index f5730492..968d21c4 100755 --- a/potpyri/instruments/MMIRS.py +++ b/potpyri/instruments/MMIRS.py @@ -1,6 +1,5 @@ -# Parameter file for MMIRS/MMT - -__version__ = "2.0" # Last edited 09/21/2024 +"""MMIRS/MMT instrument configuration and reduction parameters.""" +from potpyri._version import __version__ import os import ccdproc @@ -15,6 +14,7 @@ from . import instrument class MMIRS(instrument.Instrument): + """MMIRS at MMT: NIR imaging with flat calibration.""" def __init__(self): diff --git a/potpyri/instruments/MOSFIRE.py b/potpyri/instruments/MOSFIRE.py index 88f7170b..e4250da8 100755 --- a/potpyri/instruments/MOSFIRE.py +++ b/potpyri/instruments/MOSFIRE.py @@ -1,6 +1,5 @@ -# Parameter file for MOSFIRE/Keck - -__version__ = "2.0" # Last edited 09/21/2024 +"""MOSFIRE/Keck instrument configuration and reduction parameters.""" +from potpyri._version import __version__ import os import ccdproc @@ -15,6 +14,7 @@ from . import instrument class MOSFIRE(instrument.Instrument): + """MOSFIRE at Keck: NIR imaging with flat calibration.""" def __init__(self): diff --git a/potpyri/instruments/__init__.py b/potpyri/instruments/__init__.py index 52f2ea0d..d3a721be 100755 --- a/potpyri/instruments/__init__.py +++ b/potpyri/instruments/__init__.py @@ -1,3 +1,9 @@ +"""Instrument implementations and factory for POTPyRI. + +Each instrument module defines a subclass of Instrument with detector keywords, +calibration behavior, and file-sorting rules. New instruments must be added +here and to __all__. +""" from . import BINOSPEC from . import DEIMOS from . import F2 @@ -13,7 +19,26 @@ def instrument_getter(instname, log=None): + """Return an Instrument instance for the given instrument name. + + Parameters + ---------- + instname : str + Instrument name (e.g. 'GMOS', 'LRIS', 'BINOSPEC'). + log : ColoredLogger, optional + Logger for error messages. If None, raises Exception on unsupported + instrument. + + Returns + ------- + Instrument + Subclass instance for the requested instrument. + Raises + ------ + Exception + If instname is not in __all__ and log is None. + """ instname = instname.upper() if instname not in __all__: diff --git a/potpyri/instruments/instrument.py b/potpyri/instruments/instrument.py index 3a0a76c3..074ced65 100755 --- a/potpyri/instruments/instrument.py +++ b/potpyri/instruments/instrument.py @@ -1,6 +1,11 @@ -# Basic instrument class +"""Base instrument class and shared reduction logic for POTPyRI. -__version__ = "1.2" # Last edited 09/29/2024 +Defines detector keywords, calibration flags, file naming, and methods for +importing/processing science frames, building master bias/dark/flat, and +loading static masks. Subclasses (GMOS, LRIS, BINOSPEC, etc.) set +instrument-specific attributes. +""" +from potpyri._version import __version__ import os import astropy @@ -20,9 +25,13 @@ from astropy.stats import SigmaClip from astropy.time import Time -# Generic instrument class for POTPyRI - class Instrument(object): + """Base class for all POTPyRI instruments. + + Holds detector parameters (gain, rdnoise, pixscale, saturation), header + keyword names, file-sorting rules, and methods for calibration and + science processing. Subclasses override attributes as needed. + """ def __init__(self): @@ -99,8 +108,21 @@ def __init__(self): # Default final image size for binning of 1x1 self.out_size = 5000 - # Match file type keywords to file table def match_type_keywords(self, kwds, file_table): + """Return a boolean mask of file_table rows whose Type is in kwds. + + Parameters + ---------- + kwds : str + Comma-separated type names (e.g. 'SCIENCE', 'FLAT'). + file_table : astropy.table.Table + File list table with 'Type' column. + + Returns + ------- + np.ndarray + Boolean mask, True where file_table['Type'] is in kwds. + """ masks = [] for kwd in kwds.split(','): masks.append(file_table['Type']==kwd) @@ -110,45 +132,88 @@ def match_type_keywords(self, kwds, file_table): return(mask) - # Determine whether sky subtraction is needed for the current reduction def needs_sky_subtraction(self, filt): - # Default is to do sky subtraction for all NIR cameras + """Return True if sky subtraction should be run for this filter. + + Parameters + ---------- + filt : str + Filter name (e.g. 'r', 'J'). + + Returns + ------- + bool + True for NIR wavelength instruments, False otherwise. + """ if self.wavelength=='NIR': return(True) else: return(False) - # Get pixel scale for imagers with varying focal ratios (e.g., IMACS) def get_pixscale(self): + """Return pixel scale in arcsec/pixel. + + Returns + ------- + float or None + Pixel scale; may vary by mode in subclasses. + """ return(self.pixscale) - # Use these if a single value is needed for gain, rdnoise, etc. def get_rdnoise(self, hdr): + """Return read noise value(s) for the detector/amplifier. + + Parameters + ---------- + hdr : astropy.io.fits.Header + FITS header (may be used in subclasses for amp-dependent values). + + Returns + ------- + float or list + Read noise in electrons. + """ return(self.rdnoise) def get_gain(self, hdr): + """Return gain value(s) for the detector/amplifier. + + Parameters + ---------- + hdr : astropy.io.fits.Header + FITS header (may be used in subclasses for amp-dependent values). + + Returns + ------- + float or list + Gain in e-/ADU. + """ return(self.gain) - # Get specific header keywords from file def get_target(self, hdr): + """Return target name from header (target_keyword), spaces stripped.""" return(hdr[self.target_keyword].replace(' ','')) def get_filter(self, hdr): + """Return filter name from header; leading segment before underscore.""" filt = hdr[self.filter_keyword] filt = filt.replace(' ','') filt = filt.split('_')[0] return(filt) def get_exptime(self, hdr): + """Return exposure time from header (exptime_keyword).""" return(hdr[self.exptime_keyword]) def get_ampl(self, hdr): + """Return amplifier identifier from header as string.""" if self.amp_keyword in hdr.keys(): return(str(hdr[self.amp_keyword])) else: return(str(self.amp_keyword)) def get_binning(self, hdr): + """Return binning string from header (e.g. '11', '22'), normalized and no 'x'.""" if self.bin_keyword in hdr.keys(): binn = str(hdr[self.bin_keyword]).lower().replace(' ','') binn = binn.replace('x','') @@ -158,18 +223,34 @@ def get_binning(self, hdr): return(self.bin_keyword) def get_out_size(self, hdr): + """Return output image size (pixels) for this binning (out_size / binn).""" binn = int(str(self.get_binning(hdr))[0]) out_size = int(self.out_size/binn) return(out_size) def get_time(self, hdr): + """Return MJD from header (mjd_keyword) as float.""" return(float(hdr[self.mjd_keyword])) def get_instrument_name(self, hdr): + """Return instrument name (lowercase) for paths/filenames.""" return(self.name.lower()) def get_staticmask_filename(self, hdr, paths): - + """Return path(s) to the static bad-pixel mask FITS for this instrument/binning. + + Parameters + ---------- + hdr : astropy.io.fits.Header + FITS header (for binning and instrument). + paths : dict + Paths dict with 'code' key (package root). + + Returns + ------- + list + One-element list with path to staticmask FITS, or [None] if not found. + """ instname = self.get_instrument_name(hdr) binn = self.get_binning(hdr) @@ -182,9 +263,24 @@ def get_staticmask_filename(self, hdr, paths): return([None]) def get_catalog(self, hdr): + """Return catalog name for zeropoint (e.g. 'PS1').""" return(self.catalog_zp) def format_datasec(self, sec_string, binning=1): + """Convert datasec string to binned pixel bounds (e.g. '[1:100,1:200]'). + + Parameters + ---------- + sec_string : str + Data section string like '[x1:x2,y1:y2]'. + binning : int, optional + Binning factor. Default is 1. + + Returns + ------- + str + Reformatted section string with binned limits. + """ sec_string = sec_string.replace('[','').replace(']','') x,y = sec_string.split(',') x1,x2 = x.split(':') @@ -202,13 +298,14 @@ def format_datasec(self, sec_string, binning=1): return(sec_string) def raw_format(self, proc): + """Return glob pattern for raw files (e.g. 'sci_img_*.fits' or 'sci_img*[!proc].fits').""" if proc: return('sci_img_*.fits') else: return('sci_img*[!proc].fits') def get_stk_name(self, hdr, red_path): - + """Return full path for stacked output FITS (target.filter.utYYMMDD.amp.binn.stk.fits).""" target = self.get_target(hdr) fil = self.get_filter(hdr) amp = self.get_ampl(hdr) @@ -222,7 +319,7 @@ def get_stk_name(self, hdr, red_path): return(fullfilename) def get_sci_name(self, hdr, red_path): - + """Return full path for single-extension science FITS.""" target = self.get_target(hdr) fil = self.get_filter(hdr) amp = self.get_ampl(hdr) @@ -237,29 +334,54 @@ def get_sci_name(self, hdr, red_path): return(fullfilename) def get_bkg_name(self, hdr, red_path): - + """Return full path for background FITS (sci name with _bkg.fits).""" sci_name = self.get_sci_name(hdr, red_path) bkg_filename = sci_name.replace('.fits','_bkg.fits') return(bkg_filename) def get_mbias_name(self, paths, amp, binn): + """Return path for master bias FITS (paths['cal']/mbias_amp_binn.fits).""" red_path = paths['cal'] return(os.path.join(red_path, f'mbias_{amp}_{binn}.fits')) def get_mdark_name(self, paths, amp, binn): + """Return path for master dark FITS.""" red_path = paths['cal'] return(os.path.join(red_path, f'mdark_{amp}_{binn}.fits')) def get_mflat_name(self, paths, fil, amp, binn): + """Return path for master flat FITS (filter, amp, binning).""" red_path = paths['cal'] return(os.path.join(red_path, f'mflat_{fil}_{amp}_{binn}.fits')) def get_msky_name(self, paths, fil, amp, binn): + """Return path for master sky FITS.""" red_path = paths['cal'] return(os.path.join(red_path, f'msky_{fil}_{amp}_{binn}.fits')) def load_bias(self, paths, amp, binn): + """Load master bias CCDData for given amp and binning. + + Parameters + ---------- + paths : dict + Paths dict (cal key). + amp : str + Amplifier identifier. + binn : str + Binning string. + + Returns + ------- + ccdproc.CCDData + Master bias in electrons. + + Raises + ------ + Exception + If master bias file not found. + """ bias = self.get_mbias_name(paths, amp, binn) if os.path.exists(bias): mbias = CCDData.read(bias) @@ -271,6 +393,27 @@ def load_bias(self, paths, amp, binn): return(mbias) def load_dark(self, paths, amp, binn): + """Load master dark CCDData for given amp and binning. + + Parameters + ---------- + paths : dict + Paths dict (cal key). + amp : str + Amplifier identifier. + binn : str + Binning string. + + Returns + ------- + ccdproc.CCDData + Master dark in electrons. + + Raises + ------ + Exception + If master dark file not found. + """ dark = self.get_mdark_name(paths, amp, binn) if os.path.exists(dark): mdark = CCDData.read(dark) @@ -282,6 +425,29 @@ def load_dark(self, paths, amp, binn): return(mdark) def load_flat(self, paths, fil, amp, binn): + """Load master flat CCDData for given filter, amp, and binning. + + Parameters + ---------- + paths : dict + Paths dict (cal key). + fil : str + Filter name. + amp : str + Amplifier identifier. + binn : str + Binning string. + + Returns + ------- + ccdproc.CCDData + Master flat (normalized). + + Raises + ------ + Exception + If master flat file not found. + """ flat = self.get_mflat_name(paths, fil, amp, binn) if os.path.exists(flat): mflat = CCDData.read(flat) @@ -293,6 +459,29 @@ def load_flat(self, paths, fil, amp, binn): return(mflat) def load_sky(self, paths, fil, amp, binn): + """Load master sky CCDData for given filter, amp, and binning. + + Parameters + ---------- + paths : dict + Paths dict (cal key). + fil : str + Filter name. + amp : str + Amplifier identifier. + binn : str + Binning string. + + Returns + ------- + ccdproc.CCDData + Master sky (normalized). + + Raises + ------ + Exception + If master sky file not found. + """ sky = self.get_msky_name(paths, fil, amp, binn) if os.path.exists(sky): msky = CCDData.read(sky) @@ -304,6 +493,22 @@ def load_sky(self, paths, fil, amp, binn): return(msky) def import_image(self, filename, amp, log=None): + """Load raw multi-extension FITS and return processed CCDData (bias/overscan/trim/gain). + + Parameters + ---------- + filename : str + Path to raw FITS file. + amp : str + Amplifier identifier (number of extensions). + log : ColoredLogger, optional + Logger for progress. + + Returns + ------- + ccdproc.CCDData + Processed frame in electrons (single combined image). + """ filename = os.path.abspath(filename) if log: log.info(f'Loading file: {filename}') @@ -326,6 +531,20 @@ def import_image(self, filename, amp, log=None): return(frame_full) def import_sci_image(self, filename, log=None): + """Load a single science FITS as CCDData (no overscan/trim; already reduced). + + Parameters + ---------- + filename : str + Path to FITS file. + log : ColoredLogger, optional + Logger for progress. + + Returns + ------- + ccdproc.CCDData + Science frame. + """ filename = os.path.abspath(filename) if log: log.info(f'Loading file: {filename}') @@ -335,7 +554,28 @@ def import_sci_image(self, filename, log=None): def create_bias(self, bias_list, amp, binn, paths, log=None, **kwargs): - + """Combine bias frames into master bias and write to paths['cal']. + + Parameters + ---------- + bias_list : list of str + Paths to bias FITS files. + amp : str + Amplifier identifier. + binn : str + Binning string. + paths : dict + Paths dict (cal key). + log : ColoredLogger, optional + Logger for progress. + **kwargs + Unused; for subclass compatibility. + + Returns + ------- + None + Master bias written via get_mbias_name. + """ if log: log.info(f'Processing bias files with {amp} amps and {binn} binning.') log.info(f'{len(bias_list)} files found.') @@ -381,7 +621,28 @@ def create_bias(self, bias_list, amp, binn, paths, return def create_dark(self, dark_list, amp, binn, paths, mbias=None, log=None): - + """Combine dark frames (with optional bias subtraction) and write master dark. + + Parameters + ---------- + dark_list : list of str + Paths to dark FITS files. + amp : str + Amplifier identifier. + binn : str + Binning string. + paths : dict + Paths dict (cal key). + mbias : ccdproc.CCDData, optional + Master bias to subtract before combining. + log : ColoredLogger, optional + Logger for progress. + + Returns + ------- + None + Master dark written via get_mdark_name. + """ if log: log.info(f'Processing dark files with {amp} amps and {binn} binning.') log.info(f'{len(dark_list)} files found.') @@ -429,9 +690,38 @@ def create_dark(self, dark_list, amp, binn, paths, mbias=None, log=None): return - def create_flat(self, flat_list, fil, amp, binn, paths, mbias=None, + def create_flat(self, flat_list, fil, amp, binn, paths, mbias=None, mdark=None, is_science=False, log=None, **kwargs): - + """Combine flat frames (bias/dark subtract, normalize) and write master flat. + + Parameters + ---------- + flat_list : list of str + Paths to flat FITS files. + fil : str + Filter name. + amp : str + Amplifier identifier. + binn : str + Binning string. + paths : dict + Paths dict (cal key). + mbias : ccdproc.CCDData, optional + Master bias to subtract. + mdark : ccdproc.CCDData, optional + Master dark to subtract. + is_science : bool, optional + If True, mask high pixels as for science. Default is False. + log : ColoredLogger, optional + Logger for progress. + **kwargs + Unused; for subclass compatibility. + + Returns + ------- + None + Master flat written via get_mflat_name. + """ if log: log.info(f'Processing files for filter: {fil}') log.info(f'{len(flat_list)} files found.') @@ -514,7 +804,30 @@ def create_flat(self, flat_list, fil, amp, binn, paths, mbias=None, return def create_sky(self, sky_list, fil, amp, binn, paths, log=None, **kwargs): - + """Combine sky frames (normalize by median) and write master sky. + + Parameters + ---------- + sky_list : list of str + Paths to sky FITS files. + fil : str + Filter name. + amp : str + Amplifier identifier. + binn : str + Binning string. + paths : dict + Paths dict (cal key). + log : ColoredLogger, optional + Logger for progress. + **kwargs + Unused; for subclass compatibility. + + Returns + ------- + None + Master sky written via get_msky_name. + """ if log: log.info(f'Processing files for filter: {fil}') log.info(f'{len(sky_list)} files found.') @@ -574,22 +887,48 @@ def create_sky(self, sky_list, fil, amp, binn, paths, log=None, **kwargs): return def load_staticmask(self, hdr, paths): - + """Load the static mask array for this instrument and binning. + + Parameters + ---------- + hdr : astropy.io.fits.Header + FITS header (for binning and instrument name). + paths : dict + Paths dict with 'code' key. + + Returns + ------- + np.ndarray or None + Boolean mask from static mask FITS extension [1], or None if not found. + """ satmask_filename = self.get_staticmask_filename(hdr, paths) - if (satmask_filename is not None and - len(satmask_filename)>0 and + if (satmask_filename is not None and + len(satmask_filename) > 0 and satmask_filename[0] is not None and os.path.exists(satmask_filename[0])): - hdu = fits.open(satmask_filename[0]) - satmask = hdu[1].data.astype(bool) + with fits.open(satmask_filename[0]) as hdu: + satmask = hdu[1].data.astype(bool) else: satmask = None return(satmask) def expand_mask(self, input_data, input_mask=None): - + """Build combined mask from NaN, inf, zero pixels and optional input mask. + + Parameters + ---------- + input_data : ccdproc.CCDData + Image data (input_data.data used). + input_mask : np.ndarray, optional + Boolean mask to OR with bad-pixel mask. + + Returns + ------- + np.ndarray + Boolean mask (True = bad). + """ # Get image statistics mean, median, stddev = sigma_clipped_stats(input_data.data) @@ -608,7 +947,38 @@ def expand_mask(self, input_data, input_mask=None): def process_science(self, sci_list, fil, amp, binn, paths, mbias=None, mflat=None, mdark=None, skip_skysub=False, save_bkg=False, log=None): - + """Reduce science frames: bias/dark/flat, optional sky subtraction; return CCDData list. + + Parameters + ---------- + sci_list : list of str + Paths to science FITS files. + fil : str + Filter name. + amp : str + Amplifier identifier. + binn : str + Binning string. + paths : dict + Paths dict (work, cal, etc.). + mbias : ccdproc.CCDData, optional + Master bias. + mflat : ccdproc.CCDData, optional + Master flat. + mdark : ccdproc.CCDData, optional + Master dark. + skip_skysub : bool, optional + If True, skip 2D background subtraction. Default is False. + save_bkg : bool, optional + If True, save background model. Default is False. + log : ColoredLogger, optional + Logger for progress. + + Returns + ------- + list of ccdproc.CCDData + Reduced science frames (each with FILENAME in header). + """ processed = [] processed_names = [] for sci in sorted(sci_list): diff --git a/potpyri/primitives/__init__.py b/potpyri/primitives/__init__.py index 5cbb124d..5cf4896a 100755 --- a/potpyri/primitives/__init__.py +++ b/potpyri/primitives/__init__.py @@ -1,3 +1,4 @@ +"""POTPyRI primitives: calibration, image processing, sorting, WCS, photometry, absphot.""" from . import absphot from . import calibration from . import image_procs diff --git a/potpyri/primitives/absphot.py b/potpyri/primitives/absphot.py index 5430ae9d..fdf55801 100755 --- a/potpyri/primitives/absphot.py +++ b/potpyri/primitives/absphot.py @@ -1,12 +1,27 @@ -"Function for calculating zero points during flux calibration." -"Authors: Kerry Paterson, Charlie Kilpatrick" +"""Absolute photometry zeropoint calibration using catalog magnitudes. -__version__ = "1.0" # Last updated on 03/10/2025 +Queries Vizier (e.g. PS1), matches sources, and fits zeropoint via +iterative ODR. Writes ZPTMAG and related keywords to the stack header. +Uses multiple Vizier mirrors (CDS, Tokyo, CfA, INASAN, IUCAA, IDIA) with +fallback when the first attempt fails. +Authors: Kerry Paterson, Charlie Kilpatrick. +""" +from potpyri._version import __version__ import numpy as np from astroquery.vizier import Vizier from astropy import units as u + +# Vizier mirror servers (hostnames only); tried in order when a query fails. +VIZIER_MIRRORS = [ + 'vizier.cds.unistra.fr', # CDS / Strasbourg, France + 'vizier.nao.ac.jp', # ADAC / Tokyo, Japan + 'vizier.cfa.harvard.edu', # CfA / Harvard, USA + 'vizier.inasan.ru', # INASAN / Moscow, Russia + 'vizier.iucaa.in', # IUCAA / Pune, India + 'vizier.idia.ac.za', # IDIA / South Africa +] from astropy.coordinates import SkyCoord from astropy.table import Table from astropy.io import fits @@ -15,18 +30,48 @@ from potpyri.utils import utilities class absphot(object): - def __init__(self, iterations=5, sigma=5): + """Zeropoint fitter using catalog magnitudes and iterative sigma clipping.""" + def __init__(self, iterations=5, sigma=5): + """Initialize zeropoint fitter. + + Parameters + ---------- + iterations : int, optional + Number of sigma-clip iterations. Default is 5. + sigma : float, optional + Sigma threshold for clipping. Default is 5. + """ self.iterations = iterations self.sigma = sigma self.odr_iter = 2 def get_zeropoint(self, flux, fluxerr, mag, magerr): - - from scipy.odr import ODR - from scipy.odr import Model - from scipy.odr import RealData + """Fit zeropoint (and error) from flux/mag arrays using ODR. + + Parameters + ---------- + flux : array-like + Object fluxes (e.g. from aperture photometry). + fluxerr : array-like + Flux errors. + mag : array-like + Catalog magnitudes. + magerr : array-like + Catalog magnitude errors. + + Returns + ------- + tuple of float + (zeropoint, zeropoint_error). + """ + import warnings + with warnings.catch_warnings(): + warnings.filterwarnings('ignore', message=r'.*scipy\.odr.*', category=DeprecationWarning) + from scipy.odr import ODR + from scipy.odr import Model + from scipy.odr import RealData def magnitude(zpt, flux): return(zpt - 2.5*np.log10(flux)) @@ -49,7 +94,27 @@ def magnitude(zpt, flux): return(zpt, zpterr) def zpt_iteration(self, flux, fluxerr, mag, magerr, log=None): - + """Iteratively sigma-clip and fit zeropoint via ODR. + + Parameters + ---------- + flux : array-like + Object fluxes. + fluxerr : array-like + Flux errors. + mag : array-like + Catalog magnitudes. + magerr : array-like + Catalog magnitude errors. + log : ColoredLogger, optional + Logger for progress. + + Returns + ------- + tuple + (zeropoint, zeropoint_error, master_mask). master_mask is boolean + array of sources used in final fit. + """ flux = np.array(flux) fluxerr = np.array(fluxerr) mag = np.array(mag) @@ -95,7 +160,25 @@ def zpt_iteration(self, flux, fluxerr, mag, magerr, log=None): return(zpt, zpterr, master_mask) def get_catalog(self, coords, catalog, filt, log=None): - + """Query VizieR for catalog magnitudes in a filter around given coordinates. + + Parameters + ---------- + coords : astropy.coordinates.SkyCoord + Target coordinates (used for region size). + catalog : str + Catalog name (e.g. 'PS1', '2MASS'). + filt : str + Filter name (e.g. 'g', 'r', 'J'). + log : ColoredLogger, optional + Logger for progress. + + Returns + ------- + tuple or (None, None, None) + (astropy.table.Table, catalog_name, catalog_ID) or (None, None, None) + if query fails. + """ if log: log.info(f'Searching for catalog {catalog}') coord_ra = np.median([c.ra.degree for c in coords]) @@ -113,17 +196,35 @@ def get_catalog(self, coords, catalog, filt, log=None): if cat_ID=='II/349': cols.append(f'{filt}Kmag') - vizier = Vizier(columns=cols) - vizier.ROW_LIMIT = -1 - if log: + if log: log.info(f'Getting {catalog} catalog with ID {cat_ID} in filt {filt}') log.info(f'Querying around {coord_ra}, {coord_dec} deg') Vizier.clear_cache() width = np.max([2.0 * max_sep, 0.5]) - cat = vizier.query_region(med_coord, width=width*u.degree, - catalog=cat_ID) + cat = None + last_error = None + for server in VIZIER_MIRRORS: + try: + vizier = Vizier(columns=cols, vizier_server=server) + vizier.ROW_LIMIT = -1 + cat = vizier.query_region(med_coord, width=width*u.degree, + catalog=cat_ID) + if cat is not None and len(cat) > 0: + if log: + log.info(f'Vizier query succeeded via {server}') + break + cat = None + except Exception as e: + last_error = e + if log: + log.warning(f'Vizier mirror {server} failed: {e}') + else: + print(f'Vizier mirror {server} failed: {e}') + cat = None + if cat is None and last_error is not None and log: + log.warning('All Vizier mirrors failed; last error: {}'.format(last_error)) - if len(cat)>0: + if cat is not None and len(cat)>0: cat = cat[0] cat = cat[~np.isnan(cat[cat_mag])] cat = cat[cat[cat_err]>0.] @@ -175,7 +276,31 @@ def get_catalog(self, coords, catalog, filt, log=None): def find_zeropoint(self, cmpfile, tel, match_radius=2.5*u.arcsec, phottable='APPPHOT', input_catalog=None, log=None): - + """Compute zeropoint from cmpfile photometry and catalog; write to FITS header. + + Matches sources to catalog (e.g. PS1), runs iterative ODR fit, and + updates ZPTMAG, ZPTNSTAR, ZPTCAT, etc. in the stack FITS. + + Parameters + ---------- + cmpfile : str + Path to stacked/comparison FITS with SCI and phottable extensions. + tel : Instrument + Instrument instance (for get_catalog). + match_radius : astropy.units.Quantity, optional + Matching radius for catalog. Default is 2.5 arcsec. + phottable : str, optional + FITS extension with photometry table. Default is 'APPPHOT'. + input_catalog : astropy.table.Table, optional + Pre-loaded catalog; if None, catalog is queried via get_catalog. + log : ColoredLogger, optional + Logger for progress. + + Returns + ------- + None + cmpfile is updated in place. + """ if log: log.info(f'Importing catalog from file: {cmpfile}') else: @@ -291,6 +416,20 @@ def find_zeropoint(self, cmpfile, tel, match_radius=2.5*u.arcsec, log.info('No zeropoint calculated.') def Y_band(self, J, J_err, K, K_err): + """Compute Y-band mag and error from J and K (2MASS relation). + + Parameters + ---------- + J, J_err : array-like or float + J magnitude(s) and error(s). + K, K_err : array-like or float + K magnitude(s) and error(s). + + Returns + ------- + tuple + (Y_mag, Y_err). + """ Y = J+0.46*(J-K) JK_err = np.sqrt(J_err**2+K_err**2) JKY_err = 0.46*(J-K)*np.sqrt((0.02/0.46)**2+(JK_err/(J-K))**2) @@ -298,6 +437,18 @@ def Y_band(self, J, J_err, K, K_err): return Y, Y_err def convert_filter_name(self, filt): + """Map instrument filter names to catalog filter names (e.g. PS1). + + Parameters + ---------- + filt : str + Instrument filter keyword (e.g. 'gG0301', 'RG850'). + + Returns + ------- + str + Catalog filter name ('u', 'g', 'r', 'i', 'z', 'J', etc.). + """ if filt=='uG0308' or filt=='uG0332' or filt=='U': return 'u' if filt=='gG0301' or filt=='gG0325' or filt=='G' or filt=='V' or filt=='B': @@ -316,6 +467,18 @@ def convert_filter_name(self, filt): return filt def get_minmag(self, filt): + """Return minimum catalog magnitude to use for zeropoint (bright limit). + + Parameters + ---------- + filt : str + Filter name. + + Returns + ------- + float + Minimum magnitude (brighter limit). + """ if filt=='J': return 15.5 if filt=='K': @@ -326,5 +489,21 @@ def get_minmag(self, filt): return 16.0 def find_zeropoint(stack, tel, log=None): + """Compute and write zeropoint to stack FITS using instrument catalog (e.g. PS1). + + Parameters + ---------- + stack : str + Path to stacked science FITS file (SCI + APPPHOT extensions). + tel : Instrument + Instrument instance (for catalog and filter). + log : ColoredLogger, optional + Logger for progress and errors. + + Returns + ------- + None + Stack FITS is updated in place with ZPTMAG, ZPTNSTAR, etc. + """ cal = absphot() cal.find_zeropoint(stack, tel, log=log) diff --git a/potpyri/primitives/calibration.py b/potpyri/primitives/calibration.py index 2d713ce6..d908cf4d 100755 --- a/potpyri/primitives/calibration.py +++ b/potpyri/primitives/calibration.py @@ -1,8 +1,9 @@ -"Wrapper methods for creating bias, dark, and flat images from pipeline data." -"Authors: Charlie Kilpatrick" +"""Master bias, dark, and flat creation from pipeline file tables. -# Initial version tracking on 09/29/2024 -__version__ = "1.1" +Orchestrates grouping by CalType and calling instrument-specific +create_bias/create_dark/create_flat. Authors: Charlie Kilpatrick. +""" +from potpyri._version import __version__ import os import time @@ -11,7 +12,26 @@ import sys def do_bias(file_table, tel, paths, nmin_images=3, log=None): - + """Build master bias frames from file_table; skip if instrument has no bias. + + Parameters + ---------- + file_table : astropy.table.Table + File list from sort_files (Type, CalType, File, Amp, Binning). + tel : Instrument + Instrument instance (bias, match_type_keywords, create_bias, get_mbias_name). + paths : dict + Paths dict from options.add_paths. + nmin_images : int, optional + Minimum images per CalType to build master. Default is 3. + log : ColoredLogger, optional + Logger for progress. + + Returns + ------- + None + Master bias FITS written to paths; exits if no bias and instrument requires it. + """ # Exit if telescope does not require bias if not tel.bias: return(None) @@ -62,7 +82,26 @@ def do_bias(file_table, tel, paths, nmin_images=3, log=None): sys.exit(-1) def do_dark(file_table, tel, paths, nmin_images=3, log=None): - + """Build master dark frames from file_table; skip if instrument has no dark. + + Parameters + ---------- + file_table : astropy.table.Table + File list from sort_files. + tel : Instrument + Instrument instance (dark, bias, load_bias, create_dark, get_mdark_name). + paths : dict + Paths dict from options.add_paths. + nmin_images : int, optional + Minimum images per CalType to build master. Default is 3. + log : ColoredLogger, optional + Logger for progress. + + Returns + ------- + None + Master dark FITS written to paths. + """ # Exit if telescope does not require dark if not tel.dark: return(None) @@ -116,7 +155,26 @@ def do_dark(file_table, tel, paths, nmin_images=3, log=None): if log: log.info(f'Master dark creation completed in {t2-t1} sec.') def do_flat(file_table, tel, paths, nmin_images=3, log=None): - + """Build master flat frames from file_table; skip if instrument has no flat. + + Parameters + ---------- + file_table : astropy.table.Table + File list from sort_files. + tel : Instrument + Instrument instance (flat, match_type_keywords, create_flat, get_mflat_name). + paths : dict + Paths dict from options.add_paths. + nmin_images : int, optional + Minimum images per (filter, amp, binning) to build master. Default is 3. + log : ColoredLogger, optional + Logger for progress. + + Returns + ------- + None + Master flat FITS written to paths. + """ # Exit if telescope does not require dark if not tel.flat: return(None) diff --git a/potpyri/primitives/image_procs.py b/potpyri/primitives/image_procs.py index a07b54c9..d0b98d87 100755 --- a/potpyri/primitives/image_procs.py +++ b/potpyri/primitives/image_procs.py @@ -1,12 +1,13 @@ -"Functions for calibrating, masking, creating error images, and stacking images." -"Authors: Charlie Kilpatrick" +"""Image calibration, masking, error arrays, alignment, and stacking. -# Last updated on 09/29/2024 -__version__ = "2.1" +Provides WCS handling, alignment to a common grid, cosmic-ray rejection, +satellite masking, master calibration combination, and stacked science +output. Authors: Charlie Kilpatrick. +""" +from potpyri._version import __version__ -import os +import os import time -import importlib import numpy as np import logging import sys @@ -32,7 +33,18 @@ from potpyri.utils import utilities def remove_pv_distortion(header): - + """Remove PV* distortion keywords from a FITS header (modifies in place). + + Parameters + ---------- + header : astropy.io.fits.Header + Header to modify. + + Returns + ------- + astropy.io.fits.Header + The same header reference (modified in place). + """ done = False while not done: bad_key = False @@ -47,19 +59,30 @@ def remove_pv_distortion(header): return(header) def get_fieldcenter(images): - - ras=[] ; des = [] + """Compute mean RA/Dec of image centers from a list of FITS files. + + Parameters + ---------- + images : list of str + Paths to FITS files with WCS in primary header. + + Returns + ------- + list + [mean_ra_deg, mean_dec_deg]. + """ + ras = [] ; des = [] for file in images: - hdu = fits.open(file) - w = WCS(hdu[0].header) + with fits.open(file) as hdu: + w = WCS(hdu[0].header) - data_shape = hdu[0].data.shape - center_pix = (data_shape[1]/2., data_shape[0]/2.) + data_shape = hdu[0].data.shape + center_pix = (data_shape[1]/2., data_shape[0]/2.) - coord = w.pixel_to_world(*center_pix) + coord = w.pixel_to_world(*center_pix) - ras.append(coord.ra.degree) - des.append(coord.dec.degree) + ras.append(coord.ra.degree) + des.append(coord.dec.degree) mean_ra = np.mean(ras) mean_de = np.mean(des) @@ -67,7 +90,24 @@ def get_fieldcenter(images): return([mean_ra, mean_de]) def generate_wcs(tel, binn, fieldcenter, out_size): - + """Build a TAN WCS for a given field center and output size. + + Parameters + ---------- + tel : Instrument + Instrument instance (for pixel scale). + binn : str + Binning string (e.g. '22'). + fieldcenter : sequence + [ra_deg, dec_deg] or parseable coordinate. + out_size : int + Output image size (out_size x out_size pixels). + + Returns + ------- + astropy.wcs.WCS + WCS instance. + """ w = {'NAXES':2, 'NAXIS1': out_size, 'NAXIS2': out_size, 'EQUINOX': 2000.0, 'LONPOLE': 180.0, 'LATPOLE': 0.0} @@ -96,7 +136,39 @@ def generate_wcs(tel, binn, fieldcenter, out_size): def align_images(reduced_files, paths, tel, binn, use_wcs=None, fieldcenter=None, out_size=None, skip_gaia=False, keep_all_astro=False, log=None): - + """Solve WCS and align reduced images to a common grid. + + Runs astrometry.net and optionally Gaia alignment, then reprojects to + a common WCS. + + Parameters + ---------- + reduced_files : list of str + Paths to reduced FITS files. + paths : dict + Paths dict from options.add_paths. + tel : Instrument + Instrument instance. + binn : str + Binning string. + use_wcs : astropy.wcs.WCS, optional + If provided, use this WCS instead of solving. + fieldcenter : sequence, optional + [ra, dec] for generating WCS. + out_size : int, optional + Output image size. + skip_gaia : bool, optional + If True, skip Gaia alignment. Default is False. + keep_all_astro : bool, optional + If True, keep all images regardless of astrometric dispersion. + log : ColoredLogger, optional + Logger for progress. + + Returns + ------- + tuple or None + (aligned_data, masks, errors) as lists, or None if no images aligned. + """ solved_images = [] aligned_images = [] aligned_data = [] @@ -192,9 +264,44 @@ def align_images(reduced_files, paths, tel, binn, use_wcs=None, fieldcenter=None return(aligned_images, aligned_data) -def image_proc(image_data, tel, paths, skip_skysub=False, +def image_proc(image_data, tel, paths, skip_skysub=False, fieldcenter=None, out_size=None, satellites=True, cosmic_ray=True, skip_gaia=False, keep_all_astro=False, log=None): + """Full image processing: align, mask, stack, and optionally detrend. + + Orchestrates WCS solving, alignment, satellite/cosmic-ray masking, + stacking, and optional sky subtraction. + + Parameters + ---------- + image_data : astropy.table.Table + File table subset (e.g. one TargType) with 'File' column. + tel : Instrument + Instrument instance. + paths : dict + Paths dict from options.add_paths. + skip_skysub : bool, optional + If True, skip sky subtraction. Default is False. + fieldcenter : sequence, optional + [ra, dec] for alignment. + out_size : int, optional + Output image size. + satellites : bool, optional + If True, mask satellite trails. Default is True. + cosmic_ray : bool, optional + If True, run cosmic-ray rejection. Default is True. + skip_gaia : bool, optional + If True, skip Gaia alignment. Default is False. + keep_all_astro : bool, optional + If True, keep all images regardless of astrometric dispersion. + log : ColoredLogger, optional + Logger for progress. + + Returns + ------- + str or None + Path to stacked output FITS, or None if processing failed. + """ wavelength = tel.wavelength @@ -425,7 +532,18 @@ def image_proc(image_data, tel, paths, skip_skysub=False, return(stackname) def detrend_stack(stack): - + """Remove row and column median from stack science data (detrend). + + Parameters + ---------- + stack : astropy.io.fits.HDUList + HDU list with [0]=SCI, [1]=MASK. Modified in place. + + Returns + ------- + astropy.io.fits.HDUList + The same stack reference (modified in place). + """ data = stack[0].data mask = stack[1].data.astype(bool) @@ -450,7 +568,30 @@ def detrend_stack(stack): return(stack) def stack_data(stacking_data, tel, masks, errors, mem_limit=8.0e9, log=None): - + """Combine aligned CCDData list with exposure scaling (median or average). + + Masks are applied (masked pixels set to nan) before combination. + + Parameters + ---------- + stacking_data : list of ccdproc.CCDData + Aligned science images. + tel : Instrument + Instrument instance (for stack_method and get_exptime). + masks : list of astropy.io.fits.ImageHDU + Per-image masks. + errors : list of astropy.io.fits.ImageHDU + Per-image error arrays. + mem_limit : float, optional + Memory limit in bytes for ccdproc.combine. Default is 8e9. + log : ColoredLogger, optional + Logger for progress. + + Returns + ------- + astropy.io.fits.HDUList + HDU list [science, mask, error] (mask/error from first image). + """ stack_method = tel.stack_method if not stack_method: if log: @@ -480,7 +621,20 @@ def stack_data(stacking_data, tel, masks, errors, mem_limit=8.0e9, log=None): return(sci_med) def add_stack_mask(stack, stacking_data): - + """Update stack mask from per-image masks and saturation; set masked pixels to 0. + + Parameters + ---------- + stack : astropy.io.fits.HDUList + HDU list with [0]=SCI, [1]=MASK, [2]=ERROR. Modified in place. + stacking_data : list of ccdproc.CCDData + Per-image data (headers used for SATURATE). + + Returns + ------- + astropy.io.fits.HDUList + The same stack reference (modified in place). + """ data = stack[0].data mask = stack[1].data error = stack[2].data @@ -516,7 +670,21 @@ def add_stack_mask(stack, stacking_data): return(stack) def mask_satellites(images, filenames, log=None): - + """Detect and mask satellite trails in science images using acstools.satdet. + + Parameters + ---------- + images : list of ccdproc.CCDData + Science images; data is modified in place (trail pixels set to nan). + filenames : list of str + Paths corresponding to images (used for temp files). + log : ColoredLogger, optional + Logger for progress. + + Returns + ------- + None + """ out_data = [] for i,science_data in enumerate(images): @@ -559,10 +727,46 @@ def mask_satellites(images, filenames, log=None): os.remove(tmpfile) -def create_mask(science_data, saturation, rdnoise, sigclip=3.0, - sigfrac=0.10, objlim=4.0, niter=6, outpath='', grow=0, cosmic_ray=True, +def create_mask(science_data, saturation, rdnoise, sigclip=3.0, + sigfrac=0.10, objlim=4.0, niter=6, outpath='', grow=0, cosmic_ray=True, fsmode='convolve', cleantype='medmask', log=None): - + """Build a bad-pixel/saturation/cosmic-ray mask using LACosmic and optional growth. + + Parameters + ---------- + science_data : str + Path to science FITS file. + saturation : float + Saturation level (ADU); pixels above this get flag 4. + rdnoise : float + Read noise for LACosmic (electrons). + sigclip : float, optional + LACosmic sigma clipping. Default is 3.0. + sigfrac : float, optional + LACosmic sigfrac. Default is 0.10. + objlim : float, optional + LACosmic objlim. Default is 4.0. + niter : int, optional + LACosmic iterations. Default is 6. + outpath : str, optional + Output path for debug files. Default is ''. + grow : int, optional + Pixels to grow mask (adds flag 8). Default is 0. + cosmic_ray : bool, optional + If True, run LACosmic. Default is True. + fsmode : str, optional + LACosmic fsmode. Default is 'convolve'. + cleantype : str, optional + LACosmic cleantype. Default is 'medmask'. + log : ColoredLogger, optional + Logger for progress. + + Returns + ------- + tuple + (cleaned_data_ndarray, mask_ImageHDU). Mask uses additive flags: + 1=bad, 2=cosmic ray, 4=saturated, 8=neighbor. + """ t_start = time.time() if log: @@ -699,34 +903,53 @@ def create_mask(science_data, saturation, rdnoise, sigclip=3.0, def create_error(science_data, mask_data, rdnoise): - - hdu = fits.open(science_data) - img_data = hdu[0].data.astype(np.float32) - mask = mask_data.data.astype(bool) - - # Check for issue with all of the file being masked - if np.all(mask): - rms = hdu[0].header['SATURATE'] - else: - rms = 0.5 * ( - np.percentile(img_data[~mask], 84.13) - - np.percentile(img_data[~mask], 15.86) - ) - - poisson = img_data - poisson[poisson<0.]=0. - error = np.sqrt(poisson + rms**2 + rdnoise) - - # Sanitize error array - mask = np.isnan(error) - error[mask] = np.nanmedian(error) - mask = error < 0.0 - error[mask] = np.nanmedian(error) - maxval = np.float32(hdu[0].header['SATURATE']) - mask = np.isinf(error) - error[mask] = maxval - - error_hdu = fits.PrimaryHDU(error) #create mask Primary HDU + """Compute per-pixel error array from science image, mask, and read noise. + + Error = sqrt(poisson + rms^2 + rdnoise^2). Masked pixels are excluded + from RMS estimation. science_data can be a path or HDU; mask_data is + an HDU with mask array. + + Parameters + ---------- + science_data : str or astropy.io.fits.HDUList + Path to science FITS or open HDU list. + mask_data : astropy.io.fits.PrimaryHDU or ImageHDU + HDU containing boolean or integer mask (masked pixels excluded from rms). + rdnoise : float + Read noise in electrons. + + Returns + ------- + astropy.io.fits.PrimaryHDU + Error array in electrons (BUNIT='ELECTRONS'). + """ + with fits.open(science_data) as hdu: + img_data = hdu[0].data.astype(np.float32) + mask = mask_data.data.astype(bool) + + # Check for issue with all of the file being masked + if np.all(mask): + rms = hdu[0].header['SATURATE'] + else: + rms = 0.5 * ( + np.percentile(img_data[~mask], 84.13) + - np.percentile(img_data[~mask], 15.86) + ) + + poisson = img_data.copy() + poisson[poisson < 0.] = 0. + error = np.sqrt(poisson + rms**2 + rdnoise) + + # Sanitize error array + mask = np.isnan(error) + error[mask] = np.nanmedian(error) + mask = error < 0.0 + error[mask] = np.nanmedian(error) + maxval = np.float32(hdu[0].header['SATURATE']) + mask = np.isinf(error) + error[mask] = maxval + + error_hdu = fits.PrimaryHDU(error) # create mask Primary HDU error_hdu.header['VER'] = (__version__, 'Version of image procedures used used.') error_hdu.header['USE'] = 'Error array for Poisson, read, and RMS noise' diff --git a/potpyri/primitives/photometry.py b/potpyri/primitives/photometry.py index 86953915..688211af 100755 --- a/potpyri/primitives/photometry.py +++ b/potpyri/primitives/photometry.py @@ -1,8 +1,10 @@ -"Functions for performing aperture and PSF photometry on pipeline images." -"Authors: Kerry Paterson, Charlie Kilpatrick" +"""Aperture and PSF photometry for pipeline stacks. -# Initial version tracking on 09/21/2024 -__version__ = "1.0" +Uses Source Extractor for initial catalogs, photutils for ePSF building +and PSF fitting. Writes APPPHOT, PSFPHOT, PSFSTARS, PSF, and RESIDUAL +extensions to the stack FITS. Authors: Kerry Paterson, Charlie Kilpatrick. +""" +from potpyri._version import __version__ from photutils.aperture import ApertureStats from photutils.aperture import CircularAperture @@ -36,6 +38,17 @@ def create_conv(outfile): + """Write a 3x3 CONV NORM filter file for Source Extractor. + + Parameters + ---------- + outfile : str + Path to write the convolution kernel file. + + Returns + ------- + None + """ with open(outfile, 'w') as f: f.write('CONV NORM \n') f.write('1 2 1 \n') @@ -43,6 +56,17 @@ def create_conv(outfile): f.write('1 2 1 \n') def create_params(outfile): + """Write Source Extractor parameter file (NUMBER, X_IMAGE, FWHM_IMAGE, etc.). + + Parameters + ---------- + outfile : str + Path to write the parameter file. + + Returns + ------- + None + """ with open(outfile, 'w') as f: f.write('NUMBER \n') f.write('X_IMAGE \n') @@ -56,7 +80,20 @@ def create_params(outfile): f.write('FWHM_IMAGE \n') def run_sextractor(img_file, log=None): - + """Run Source Extractor on SCI extension; return catalog table or None. + + Parameters + ---------- + img_file : str + Path to FITS with SCI extension (and SATURATE in header). + log : ColoredLogger, optional + Logger for progress. + + Returns + ------- + astropy.table.Table or None + Source Extractor catalog, or None if run failed. + """ hdu = fits.open(img_file) data = hdu['SCI'].data @@ -106,8 +143,30 @@ def run_sextractor(img_file, log=None): return(table) -def extract_aperture_stats(img_data, img_mask, img_error, stars, +def extract_aperture_stats(img_data, img_mask, img_error, stars, aperture_radius=10.0, log=None): + """Compute aperture flux and error for a star table; return table with added columns. + + Parameters + ---------- + img_data : np.ndarray + Science image data. + img_mask : np.ndarray + Boolean or integer mask (True/masked excluded). + img_error : np.ndarray + Per-pixel error array. + stars : astropy.table.Table + Table with xcentroid, ycentroid (modified in place with refined centroids). + aperture_radius : float, optional + Aperture radius in pixels. Default is 10.0. + log : ColoredLogger, optional + Logger for progress. + + Returns + ------- + astropy.table.Table + Table with fwhm, flux_best, flux_best_err, Xpos, Ypos, etc. + """ apertable = Table([[0.],[0.],[0.],[0.],[0.],[0.],[0.],[0.],[0.],[0.], [0.],[0.]], @@ -143,25 +202,49 @@ def extract_aperture_stats(img_data, img_mask, img_error, stars, aperstats = ApertureStats(img_data, aper, mask=img_mask, error=img_error) + covx = np.maximum(aperstats.covar_sigx2.value, 0.0) + covy = np.maximum(aperstats.covar_sigy2.value, 0.0) apertable.add_row([aperstats.fwhm.value, aperstats.semimajor_sigma.value, aperstats.semiminor_sigma.value, aperstats.orientation.value, aperstats.eccentricity, aperstats.sum/aperstats.sum_err, aperstats.sum, aperstats.sum_err, aperstats.xcentroid, - aperstats.ycentroid, np.sqrt(aperstats.covar_sigx2.value), - np.sqrt(aperstats.covar_sigy2.value)]) + aperstats.ycentroid, np.sqrt(covx), + np.sqrt(covy)]) return(apertable) def generate_epsf(img_file, x, y, size=11, oversampling=2, maxiters=11, log=None): + """Build an ePSF from cutouts at (x, y) positions in the image. + + Parameters + ---------- + img_file : str + Path to FITS with SCI extension. + x, y : array-like + Star x,y positions (pixels). + size : int, optional + Cutout size in pixels. Default is 11. + oversampling : int, optional + EPSFBuilder oversampling. Default is 2. + maxiters : int, optional + EPSFBuilder max iterations. Default is 11. + log : ColoredLogger, optional + Logger for progress. + + Returns + ------- + photutils.psf.EPSFModel + Fitted ePSF model. + """ # Construct stars table from bright stars_tbl = Table() stars_tbl['x'] = x stars_tbl['y'] = y - img_hdu = fits.open(img_file) - ndimage = NDData(data=img_hdu['SCI'].data) + with fits.open(img_file) as img_hdu: + ndimage = NDData(data=img_hdu['SCI'].data) stars = extract_stars(ndimage, stars_tbl, size=size) @@ -180,7 +263,20 @@ def generate_epsf(img_file, x, y, size=11, oversampling=2, maxiters=11, return(epsf) def extract_fwhm_from_epsf(epsf, fwhm_init): - + """Estimate FWHM from ePSF model (Moffat2D fit). + + Parameters + ---------- + epsf : photutils.psf.EPSFModel + ePSF model with .data array. + fwhm_init : float + Initial FWHM guess for fit (pixels). + + Returns + ------- + float + FWHM in pixels from fitted Moffat2D. + """ # Get the raw data for the FWHM and size in x and y data = epsf.data x = np.arange(data.shape[0]) @@ -204,12 +300,33 @@ def extract_fwhm_from_epsf(epsf, fwhm_init): return(p.fwhm) def run_photometry(img_file, epsf, fwhm, threshold, shape, stars): - - img_hdu = fits.open(img_file) - image = img_hdu['SCI'].data - ndimage = NDData(data=img_hdu['SCI'].data) - mask = img_hdu['MASK'].data.astype(bool) - error = img_hdu['ERROR'].data + """Run PSF photometry and aperture photometry; append result tables to FITS. + + Parameters + ---------- + img_file : str + Path to FITS with SCI, MASK, ERROR extensions. + epsf : photutils.psf.EPSFModel + ePSF model for PSF fit. + fwhm : float + FWHM in pixels (for aperture radius). + threshold : float + Detection threshold (unused in this wrapper; for compatibility). + shape : int + Fit shape (pixels) for PSFPhotometry. + stars : astropy.table.Table + Star table with xcentroid, ycentroid, flux_best for initial guesses. + + Returns + ------- + tuple + (result_table, residual_image). result_table has flux_fit, x_fit, y_fit, etc. + """ + with fits.open(img_file) as img_hdu: + image = img_hdu['SCI'].data + ndimage = NDData(data=img_hdu['SCI'].data) + mask = img_hdu['MASK'].data.astype(bool) + error = img_hdu['ERROR'].data psf = copy.copy(epsf) @@ -244,8 +361,30 @@ def run_photometry(img_file, epsf, fwhm, threshold, shape, stars): # Identifies stars and gets aperture photometry and statistics using # photutils.aperture.ApertureStats -def get_star_catalog(img_data, img_mask, img_error, fwhm_init=5.0, +def get_star_catalog(img_data, img_mask, img_error, fwhm_init=5.0, threshold=50.0, log=None): + """Detect stars with DAOStarFinder and compute aperture stats; return star table. + + Parameters + ---------- + img_data : np.ndarray + Science image data. + img_mask : np.ndarray + Boolean mask (True = masked). + img_error : np.ndarray + Per-pixel error array. + fwhm_init : float, optional + DAOStarFinder FWHM and aperture scale. Default is 5.0. + threshold : float, optional + DAOStarFinder detection threshold. Default is 50.0. + log : ColoredLogger, optional + Logger for progress. + + Returns + ------- + astropy.table.Table + Star table with centroid, flux, fwhm, flux_best, etc. + """ # Construct the finder with input FWHM and threshold daofind = DAOStarFinder(fwhm=fwhm_init, threshold=threshold, @@ -279,10 +418,34 @@ def do_phot(img_file, save_psf_img=False, save_residual_hdu=False, log=None): - + """Run full PSF and aperture photometry on a stack; write extensions to FITS. + + Builds ePSF, finds stars, runs PSF and aperture photometry, and appends + APPPHOT, PSFPHOT, EPSF, and optional residual extensions to the FITS file. + + Parameters + ---------- + img_file : str + Path to stack FITS (SCI, MASK, ERROR). + fwhm_scale_psf : float, optional + Scale factor for PSF fit shape. Default is 4.5. + oversampling : int, optional + ePSF oversampling. Default is 1. + star_param : dict, optional + Keys: snthresh_psf, fwhm_init, snthresh_final. Defaults as shown. + save_psf_img : bool, optional + If True, save ePSF image extension. Default is False. + save_residual_hdu : bool, optional + If True, save residual image extension. Default is False. + log : ColoredLogger, optional + Logger for progress. + + Returns + ------- + None + img_file is updated in place with new extensions. + """ img_hdu = fits.open(img_file) - - # Get image statistics data = img_hdu['SCI'].data mask = img_hdu['MASK'].data.astype(bool) error = img_hdu['ERROR'].data @@ -435,7 +598,10 @@ def do_phot(img_file, coords = w.pixel_to_world(photometry['x_fit'], photometry['y_fit']) # Join the photometry and all star catalogs and rename columns - photometry['SN'] = photometry['flux_fit']/photometry['flux_err'] + flux_err = np.asarray(photometry['flux_err'], dtype=float) + sn = np.full_like(flux_err, np.nan) + np.divide(photometry['flux_fit'], flux_err, out=sn, where=flux_err > 0) + photometry['SN'] = sn photometry['mag'] = -2.5*np.log10(photometry['flux_fit']) photometry['mag_err'] = 2.5/np.log(10) * 1./photometry['SN'] photometry['RA'] = [c.ra.degree for c in coords] @@ -513,8 +679,29 @@ def do_phot(img_file, # Finally, write out img_hdu to save all data img_hdu.writeto(img_file, overwrite=True) + img_hdu.close() def photloop(stack, phot_sn_min=3.0, phot_sn_max=40.0, fwhm_init=5.0, log=None): + """Try PSF photometry with decreasing S/N threshold until do_phot succeeds. + + Parameters + ---------- + stack : str + Path to stack FITS (SCI, MASK, ERROR). + phot_sn_min : float, optional + Minimum S/N threshold to try. Default is 3.0. + phot_sn_max : float, optional + Initial S/N threshold. Default is 40.0. + fwhm_init : float, optional + Initial FWHM for star finding. Default is 5.0. + log : ColoredLogger, optional + Logger for progress. + + Returns + ------- + None + Stack FITS is updated in place when do_phot succeeds. + """ signal_to_noise = phot_sn_max epsf = None ; fwhm = None diff --git a/potpyri/primitives/solve_wcs.py b/potpyri/primitives/solve_wcs.py index d4e47180..b66b6b7d 100755 --- a/potpyri/primitives/solve_wcs.py +++ b/potpyri/primitives/solve_wcs.py @@ -1,8 +1,10 @@ -"Python script for WCS solution." -"Authors: Kerry Paterson, Charlie Kilpatrick" +"""WCS solution using astrometry.net and Gaia DR3 alignment. -# Last updated 10/24/2024 -__version__ = "2.1" +Provides coarse solution via solve-field and fine alignment by matching +detected sources to Gaia. Writes RADISP/DEDISP and updated WCS to FITS. +Authors: Kerry Paterson, Charlie Kilpatrick. +""" +from potpyri._version import __version__ import numpy as np import subprocess @@ -37,10 +39,31 @@ warnings.simplefilter('ignore', category=AstropyWarning) def get_gaia_catalog(input_file, log=None): - - hdu = fits.open(input_file) - data = hdu[0].data - header = hdu[0].header + """Query Gaia DR3 for stars in the field; return filtered catalog table. + + Uses header RA/DEC to define field; filters by PSS, Plx, PM. + + Parameters + ---------- + input_file : str + Path to FITS file; header used for CRVAL1/CRVAL2 or RA/DEC. + log : ColoredLogger, optional + Logger for progress and errors. + + Returns + ------- + astropy.table.Table or None or False + Filtered Gaia DR3 table (PSS>0.99, Plx<20, PM<10), or None if empty, + or False if coordinates could not be parsed. + + Raises + ------ + Exception + If Gaia query fails after retries. + """ + with fits.open(input_file) as hdu: + data = hdu[0].data + header = hdu[0].header hkeys = list(header.keys()) check_pairs = [('CRVAL1','CRVAL2'),('RA','DEC'),('OBJCTRA','OBJCTDEC')] @@ -55,7 +78,7 @@ def get_gaia_catalog(input_file, log=None): break if not coord: - if log: log.error(f'Could not parse RA/DEC from header of {file}') + if log: log.error(f'Could not parse RA/DEC from header of {input_file}') return(False) vizier = Vizier(columns=['RA_ICRS', 'DE_ICRS', 'Plx', 'PSS','PM']) @@ -89,6 +112,21 @@ def get_gaia_catalog(input_file, log=None): return(cat) def clean_up_astrometry(directory, file, exten): + """Remove astrometry.net temporary files (.axy, .corr, .solved, .wcs, etc.). + + Parameters + ---------- + directory : str + Directory containing the files. + file : str + Base filename (with extension); exten is replaced to form temp names. + exten : str + File extension (e.g. '.fits') to replace with .axy, .corr, etc. + + Returns + ------- + None + """ filelist = [file.replace(exten,'.axy'), file.replace(exten,'.corr'), file.replace(exten,'-indx.xyls'), @@ -103,9 +141,36 @@ def clean_up_astrometry(directory, file, exten): if os.path.exists(f): os.remove(f) -def solve_astrometry(file, tel, binn, paths, radius=0.5, replace=True, +def solve_astrometry(file, tel, binn, paths, radius=0.5, replace=True, shift_only=False, index=None, log=None): - + """Run astrometry.net (solve-field) to get coarse WCS; optionally use custom index. + + Parameters + ---------- + file : str + Path to FITS file to solve. + tel : Instrument + Instrument instance (name, get_pixscale). + binn : str + Binning string (e.g. '22') for pixel scale. + paths : dict + Paths dict (source_extractor key for solve-field). + radius : float, optional + Search radius in degrees. Default is 0.5. + replace : bool, optional + If True, write WCS into original file; else write .solved.fits. + shift_only : bool, optional + If True, only update CRPIX/CRVAL. Default is False. + index : str, optional + Path to custom index file or directory for solve-field. + log : ColoredLogger, optional + Logger for progress. + + Returns + ------- + bool + True if solve succeeded, False otherwise. + """ # Starting solve, print file and directory for reference fullfile = os.path.abspath(file) directory = os.path.dirname(file) @@ -118,9 +183,9 @@ def solve_astrometry(file, tel, binn, paths, radius=0.5, replace=True, if not os.path.exists(fullfile): return(False) - hdu = fits.open(fullfile) - data = hdu[0].data - header = hdu[0].header + with fits.open(fullfile) as hdu: + data = hdu[0].data + header = hdu[0].header hkeys = list(header.keys()) exten = '.'+file.split('.')[-1] @@ -218,13 +283,12 @@ def solve_astrometry(file, tel, binn, paths, radius=0.5, replace=True, else: print(' '.join(process)) - FNULL = open(os.devnull, 'w') - - p = subprocess.Popen(process, stdout=FNULL, stderr=subprocess.STDOUT) - try: - p.wait(90) - except subprocess.TimeoutExpired: - p.kill() + with open(os.devnull, 'w') as FNULL: + p = subprocess.Popen(process, stdout=FNULL, stderr=subprocess.STDOUT) + try: + p.wait(90) + except subprocess.TimeoutExpired: + p.kill() if os.path.exists(newfile): good = True @@ -272,57 +336,54 @@ def solve_astrometry(file, tel, binn, paths, radius=0.5, replace=True, # astrometry.net replaces extra extensions, so instead of renaming # copy the new header into the old file - newhdu = fits.open(newfile) - hdu = fits.open(output_file) - - if shift_only: - hdu[0].header['CRPIX1'] = newhdu[0].header['CRPIX1'] - hdu[0].header['CRPIX2'] = newhdu[0].header['CRPIX2'] - hdu[0].header['CRVAL1'] = newhdu[0].header['CRVAL1'] - hdu[0].header['CRVAL2'] = newhdu[0].header['CRVAL2'] - else: - hdu[0].header = newhdu[0].header - - hdu.writeto(output_file, overwrite=True, output_verify='silentfix') + with fits.open(newfile) as newhdu: + with fits.open(output_file) as hdu: + if shift_only: + hdu[0].header['CRPIX1'] = newhdu[0].header['CRPIX1'] + hdu[0].header['CRPIX2'] = newhdu[0].header['CRPIX2'] + hdu[0].header['CRVAL1'] = newhdu[0].header['CRVAL1'] + hdu[0].header['CRVAL2'] = newhdu[0].header['CRVAL2'] + else: + hdu[0].header = newhdu[0].header + hdu.writeto(output_file, overwrite=True, output_verify='silentfix') os.remove(newfile) else: output_file = fullfile.replace(exten,'.solved.fits') - hdu = fits.open(output_file) - for i,h in enumerate(hdu): - - if 'COMMENT' in hdu[i].header.keys(): - del hdu[i].header['COMMENT'] - if 'HISTORY' in hdu[i].header.keys(): - del hdu[i].header['HISTORY'] - - # Delete other WCS header keys - wcs_key = ['CSYER1','CSYER2','CRDER1','CRDER2','CD1_1','CD1_2', - 'CD2_1','CD2_2','CRPIX1','CRPIX2','CUNIT1','CUNIT2','EQUINOX', - 'RADESYS','CNAME1','CNAME2','CTYPE1','CTYPE2','WCSNAME', - 'CRVAL1','CRVAL2'] - - for key in [w + 'C' for w in wcs_key]: - if key in list(hdu[i].header.keys()): - del hdu[i].header[key] - - for key in [w + 'S' for w in wcs_key]: - if key in list(hdu[i].header.keys()): - del hdu[i].header[key] - - # Delete old WCS keywords starting with '_' - for key in list(hdu[i].header.keys()): - if key.startswith('_'): - del hdu[i].header[key] - - try: - hdu.writeto(output_file, overwrite=True, output_verify='silentfix') - except TypeError: - if log: - log.error(f'FAILURE: could not save file {fullfile}') - else: - print(f'FAILURE: could not save file {fullfile}') - return(False) + with fits.open(output_file) as hdu: + for i,h in enumerate(hdu): + if 'COMMENT' in hdu[i].header.keys(): + del hdu[i].header['COMMENT'] + if 'HISTORY' in hdu[i].header.keys(): + del hdu[i].header['HISTORY'] + + # Delete other WCS header keys + wcs_key = ['CSYER1','CSYER2','CRDER1','CRDER2','CD1_1','CD1_2', + 'CD2_1','CD2_2','CRPIX1','CRPIX2','CUNIT1','CUNIT2','EQUINOX', + 'RADESYS','CNAME1','CNAME2','CTYPE1','CTYPE2','WCSNAME', + 'CRVAL1','CRVAL2'] + + for key in [w + 'C' for w in wcs_key]: + if key in list(hdu[i].header.keys()): + del hdu[i].header[key] + + for key in [w + 'S' for w in wcs_key]: + if key in list(hdu[i].header.keys()): + del hdu[i].header[key] + + # Delete old WCS keywords starting with '_' + for key in list(hdu[i].header.keys()): + if key.startswith('_'): + del hdu[i].header[key] + + try: + hdu.writeto(output_file, overwrite=True, output_verify='silentfix') + except TypeError: + if log: + log.error(f'FAILURE: could not save file {fullfile}') + else: + print(f'FAILURE: could not save file {fullfile}') + return(False) return(True) @@ -336,7 +397,30 @@ def solve_astrometry(file, tel, binn, paths, radius=0.5, replace=True, def align_to_gaia(file, tel, radius=0.5, max_search_radius=5.0*u.arcsec, save_centroids=False, min_gaia_match=7, log=None): - + """Refine WCS by matching detected sources to Gaia DR3 and fitting WCS. + + Parameters + ---------- + file : str + Path to FITS file with WCS (e.g. after solve_astrometry). + tel : Instrument + Instrument instance. + radius : float, optional + Gaia query radius in degrees. Default is 0.5. + max_search_radius : astropy.units.Quantity, optional + Match radius for source-Gaia matching. Default is 5 arcsec. + save_centroids : bool, optional + If True, save centroid table. Default is False. + min_gaia_match : int, optional + Minimum number of Gaia matches required. Default is 7. + log : ColoredLogger, optional + Logger for progress. + + Returns + ------- + bool + True if alignment succeeded and WCS was updated, False otherwise. + """ cat = get_gaia_catalog(file, log=log) if cat is None or len(cat)==0: @@ -353,7 +437,7 @@ def align_to_gaia(file, tel, radius=0.5, max_search_radius=5.0*u.arcsec, # Estimate sources that land within the image using WCS from header hdu = fits.open(file) - w = WCS(hdu[0].header) + w = WCS(hdu[0].header, relax=True) naxis1 = hdu[0].header['NAXIS1'] naxis2 = hdu[0].header['NAXIS2'] @@ -372,6 +456,7 @@ def align_to_gaia(file, tel, radius=0.5, max_search_radius=5.0*u.arcsec, log.error(f'No stars found in {file}') else: print(f'No stars found in {file}') + hdu.close() return(False) if log: @@ -390,7 +475,7 @@ def align_to_gaia(file, tel, radius=0.5, max_search_radius=5.0*u.arcsec, hdu[0].header['DEDISP']=(1.0, 'Dispersion in Decl. of WCS [Arcsec]') hdu.writeto(file, overwrite=True, output_verify='silentfix') - + hdu.close() return(True) @@ -438,8 +523,12 @@ def align_to_gaia(file, tel, radius=0.5, max_search_radius=5.0*u.arcsec, central_coo = w.pixel_to_world(*central_pix) try: - w = fit_wcs_from_points(xy, coords, proj_point=central_coo, - sip_degree=sip_degree) + with warnings.catch_warnings(): + warnings.filterwarnings( + "ignore", category=RuntimeWarning, + message=".*cdelt will be ignored since cd is present.*") + w = fit_wcs_from_points(xy, coords, proj_point=central_coo, + sip_degree=sip_degree) succeeded = True except ValueError: # Try with next iteration @@ -450,6 +539,7 @@ def align_to_gaia(file, tel, radius=0.5, max_search_radius=5.0*u.arcsec, log.error(f'Could not converge on fine WCS solution for: {file}') else: print(f'Could not converge on fine WCS solution for: {file}') + hdu.close() return(False) @@ -465,6 +555,7 @@ def align_to_gaia(file, tel, radius=0.5, max_search_radius=5.0*u.arcsec, log.error(f'Could not centroid on any stars in {file}') else: print(f'Could not centroid on any stars in {file}') + hdu.close() return(False) else: if log: @@ -483,7 +574,7 @@ def align_to_gaia(file, tel, radius=0.5, max_search_radius=5.0*u.arcsec, hdu[0].header['DEDISP']=(1.0, 'Dispersion in Decl. of WCS [Arcsec]') hdu.writeto(file, overwrite=True, output_verify='silentfix') - + hdu.close() return(True) # Bootstrap WCS solution and get center pixel @@ -510,7 +601,11 @@ def align_to_gaia(file, tel, radius=0.5, max_search_radius=5.0*u.arcsec, c = coords[idxs] central_coo = w.pixel_to_world(*central_pix) - new_wcs = fit_wcs_from_points(xy, c, projection=w) + with warnings.catch_warnings(): + warnings.filterwarnings( + "ignore", category=RuntimeWarning, + message=".*cdelt will be ignored since cd is present.*") + new_wcs = fit_wcs_from_points(xy, c, projection=w) cent_coord = new_wcs.pixel_to_world(*central_pix) ra_cent.append(cent_coord.ra.degree) @@ -579,5 +674,5 @@ def align_to_gaia(file, tel, radius=0.5, max_search_radius=5.0*u.arcsec, hdu[0].header['DEDISP']=(de_disp, 'Dispersion in Decl. of WCS [Arcsec]') hdu.writeto(file, overwrite=True, output_verify='silentfix') - + hdu.close() return(True) diff --git a/potpyri/primitives/sort_files.py b/potpyri/primitives/sort_files.py index df99d16c..754fc591 100755 --- a/potpyri/primitives/sort_files.py +++ b/potpyri/primitives/sort_files.py @@ -1,8 +1,10 @@ -"Function to sort files for main_pipeline." -"Authors: Owen Eskandari, Kerry Paterson, Charlie Kilpatrick" +"""File sorting and file-list generation for the main pipeline. -# Last updated 04/01/2025 -__version__ = "3.1" +Classifies raw files (science, flat, bias, dark, bad) from header keywords +and writes a fixed-width file list used by calibration and reduction. +Authors: Owen Eskandari, Kerry Paterson, Charlie Kilpatrick. +""" +from potpyri._version import __version__ from astropy.io import fits from astropy.io import ascii @@ -19,6 +21,20 @@ import re def is_bad(hdr, tel): + """Return True if the header matches bad_keywords/bad_values or has invalid binning. + + Parameters + ---------- + hdr : astropy.io.fits.Header + FITS header to check. + tel : Instrument + Instrument instance (bad_keywords, bad_values, get_binning). + + Returns + ------- + bool + True if file should be excluded as bad. + """ keywords = tel.bad_keywords values = tel.bad_values @@ -37,6 +53,20 @@ def is_bad(hdr, tel): return(bad) def is_spec(hdr, tel): + """Return True if the header matches spectroscopic observation keywords. + + Parameters + ---------- + hdr : astropy.io.fits.Header + FITS header to check. + tel : Instrument + Instrument instance (spec_keywords, spec_values). + + Returns + ------- + bool + True if file is spectroscopic. + """ keywords = tel.spec_keywords values = tel.spec_values @@ -49,6 +79,20 @@ def is_spec(hdr, tel): return(spec) def is_flat(hdr, tel): + """Return True if the header matches flat-field observation keywords. + + Parameters + ---------- + hdr : astropy.io.fits.Header + FITS header to check. + tel : Instrument + Instrument instance (flat_keywords, flat_values). + + Returns + ------- + bool + True if file is a flat. + """ keywords = tel.flat_keywords values = tel.flat_values @@ -61,6 +105,20 @@ def is_flat(hdr, tel): return(flat) def is_dark(hdr, tel): + """Return True if the header matches dark observation keywords and valid binning. + + Parameters + ---------- + hdr : astropy.io.fits.Header + FITS header to check. + tel : Instrument + Instrument instance (dark_keywords, dark_values, get_binning). + + Returns + ------- + bool + True if file is a dark. + """ keywords = tel.dark_keywords values = tel.dark_values @@ -81,6 +139,20 @@ def is_dark(hdr, tel): return(dark) def is_bias(hdr, tel): + """Return True if the header matches bias observation keywords and valid binning. + + Parameters + ---------- + hdr : astropy.io.fits.Header + FITS header to check. + tel : Instrument + Instrument instance (bias_keywords, bias_values, get_binning). + + Returns + ------- + bool + True if file is a bias. + """ keywords = tel.bias_keywords values = tel.bias_values @@ -101,6 +173,20 @@ def is_bias(hdr, tel): return(bias) def is_science(hdr, tel): + """Return True if the header matches science observation keywords and min exptime. + + Parameters + ---------- + hdr : astropy.io.fits.Header + FITS header to check. + tel : Instrument + Instrument instance (science_keywords, science_values, min_exptime, get_exptime). + + Returns + ------- + bool + True if file is science. + """ keywords = tel.science_keywords values = tel.science_values @@ -118,10 +204,40 @@ def is_science(hdr, tel): return(science) -# Overall method to handle files: -def handle_files(file_list, paths, tel, incl_bad=False, proc=None, +def handle_files(file_list, paths, tel, incl_bad=False, proc=None, no_redo=False, log=None): + """Build or read the file list: discover raw files, sort, and write table. + + If no_redo and file_list exists, reads existing table. Otherwise globs + raw/bad/data, runs sort_files, and writes the fixed-width list. + Parameters + ---------- + file_list : str + Path to output (or existing) file list table. + paths : dict + Paths dict with 'raw', 'data', 'bad' keys. + tel : Instrument + Instrument instance. + incl_bad : bool, optional + If True, include bad files in list. Default is False. + proc : str, optional + Processor/run identifier for raw_format glob. + no_redo : bool, optional + If True and file_list exists, read it instead of regenerating. + log : ColoredLogger, optional + Logger for progress. + + Returns + ------- + astropy.table.Table + File table with Target, Filter, Type, CalType, File, etc. + + Raises + ------ + SystemExit + If no files found or no good files after sorting. + """ file_table = None # Always regenerate file list from existing data in data, raw, and bad @@ -166,37 +282,33 @@ def handle_files(file_list, paths, tel, incl_bad=False, proc=None, return(file_table) -# Sort the calibration files: def sort_files(files, file_list, tel, paths, incl_bad=False, log=None): + """Classify files by type (science, flat, bias, dark) and write file list table. - ''' - - Function used to sort a list of files into a dictionary of files sorted by filter. + Reads headers, applies instrument keyword rules, and writes a fixed-width + file list. Parameters ---------- - - :param files: list (string) - List of strings of files (path should be included). - - :param tel: Telescope - Telescope with parameters and methods to sort data. - - :param incl_bad: bool, optional - - :param log: log, optional - Overview log used to write the object and date observed (if ``date`` parameter is not ``None``). - If no log is inputted, information is printed out instead of being written to ``log2``. - Default is ``None``. - + files : list of str + Paths to raw FITS files. + file_list : str + Path to write fixed-width file list. + tel : Instrument + Instrument instance. + paths : dict + Paths dict (raw, data, bad, work). + incl_bad : bool, optional + If True, include bad files in table. Default is False. + log : ColoredLogger, optional + Logger for progress. Returns ------- + astropy.table.Table + Table with Target, Filter, Type, CalType, File, Exp, Time, etc. + """ - :return: python dictionary - Dictionary of files. Key is the filter of the file, values are the file names themselves. - - ''' t_start = time.time() @@ -226,20 +338,23 @@ def sort_files(files, file_list, tel, paths, incl_bad=False, log=None): for i, f in enumerate(sorted(files)): try: - file_open = fits.open(f, mode='readonly') - ext = tel.raw_header_ext - hdr = file_open[ext].header - - # Extend header to first extension? - if tel.extend_header: - if len(file_open)>ext+1: - extra_hdr = file_open[ext+1].header - for key in extra_hdr.keys(): - if key not in hdr.keys(): - hdr[key] = extra_hdr[key] - - check_data = file_open[ext].data - file_open._verify() + with fits.open(f, mode='readonly') as file_open: + ext = tel.raw_header_ext + hdr = file_open[ext].header + + # Extend header to first extension? + if tel.extend_header: + if len(file_open)>ext+1: + extra_hdr = file_open[ext+1].header + for key in extra_hdr.keys(): + if key not in hdr.keys(): + value = extra_hdr[key] + if isinstance(value, (str, int, float, complex, + bool, np.floating, np.integer, np.bool_)): + hdr[key] = value + + check_data = file_open[ext].data + file_open._verify() except IndexError: if log: log.error(f'Moving file {f} to bad due to error opening file.') diff --git a/potpyri/scripts/__init__.py b/potpyri/scripts/__init__.py new file mode 100644 index 00000000..ec84f9aa --- /dev/null +++ b/potpyri/scripts/__init__.py @@ -0,0 +1 @@ +"""POTPyRI scripts: main pipeline and archive download utilities.""" diff --git a/potpyri/scripts/main_pipeline.py b/potpyri/scripts/main_pipeline.py index 065bdfc9..18fafacf 100755 --- a/potpyri/scripts/main_pipeline.py +++ b/potpyri/scripts/main_pipeline.py @@ -12,7 +12,7 @@ If you use this code for your work, please consider citing the package release. """ -__version__ = "2.2" # Last updated 03/10/2025 +from potpyri._version import __version__ import time import numpy as np @@ -30,25 +30,29 @@ # Check options.py - all named parameters need to correspond to options that are # provided via argparse. -def main_pipeline(instrument:str, - data_path:str, - target:list=None, - proc:str=None, - incl_bad:bool=None, - no_redo_sort:bool=None, - phot_sn_min:float=None, - phot_sn_max:float=None, - fwhm_init:float=None, - skip_skysub:bool=None, - file_list_name:str=None, - fieldcenter:list=None, - out_size:int=None, - skip_flatten:bool=None, - skip_cr:bool=None, - skip_gaia:bool=None, - keep_all_astro:bool=None, - **kwargs)->None: - +def main_pipeline(instrument: str, + data_path: str, + target: list = None, + proc: str = None, + incl_bad: bool = None, + no_redo_sort: bool = None, + phot_sn_min: float = None, + phot_sn_max: float = None, + fwhm_init: float = None, + skip_skysub: bool = None, + file_list_name: str = None, + fieldcenter: list = None, + out_size: int = None, + skip_flatten: bool = None, + skip_cr: bool = None, + skip_gaia: bool = None, + keep_all_astro: bool = None, + **kwargs) -> None: + """Run the full reduction pipeline: sort files, calibrate, process, stack, photometry, zeropoint. + + Parameters mirror options from options.add_options() (instrument, data_path, + target, proc, file_list_name, photometry and processing flags, etc.). + """ # start time t1 = time.time() @@ -104,6 +108,7 @@ def main_pipeline(instrument:str, log.shutdown() def main(): + """Entry point: check dependencies, parse args, and run main_pipeline.""" options.test_for_dependencies() args = options.add_options() main_pipeline(**vars(args)) diff --git a/potpyri/utils/__init__.py b/potpyri/utils/__init__.py index 04c830f4..9f177ea3 100755 --- a/potpyri/utils/__init__.py +++ b/potpyri/utils/__init__.py @@ -1,3 +1,4 @@ +"""POTPyRI utilities: logging, option parsing, path setup, and general helpers.""" from . import logger from . import options from . import utilities diff --git a/potpyri/utils/logger.py b/potpyri/utils/logger.py index 2a640758..b554d735 100755 --- a/potpyri/utils/logger.py +++ b/potpyri/utils/logger.py @@ -1,8 +1,9 @@ -"Function for logging throughout the pipeline." -"Authors: Kerry Paterson, Charlie Kilpatrick" +"""Logging utilities for the POTPyRI pipeline. -# Initial version tracking on 09/21/2024 -__version__ = "1.0" +Provides a colored console and file logger with UTC timestamps for pipeline +steps. Authors: Kerry Paterson, Charlie Kilpatrick. +""" +from potpyri._version import __version__ import logging import time @@ -20,7 +21,20 @@ BOLD_SEQ = "\033[1m" def get_log(log_dir): + """Create and return a ColoredLogger writing to the given directory. + Log filename is generated as log_YYYYMMDD_HHMMSS.log (UTC). + + Parameters + ---------- + log_dir : str + Directory path for the log file. + + Returns + ------- + ColoredLogger + Logger instance with stream and file handlers. + """ datestr = datetime.datetime.now(datetime.UTC).strftime('%Y%m%d_%H%M%S') base_logname = f'log_{datestr}.log' log_filename = os.path.join(log_dir, base_logname) @@ -28,8 +42,22 @@ def get_log(log_dir): return(log) -def formatter_message(message, use_color = True): - # Bold face some text only when rich text/color is requested. +def formatter_message(message, use_color=True): + """Replace $RESET and $BOLD placeholders with ANSI codes or empty strings. + + Parameters + ---------- + message : str + Format string possibly containing $RESET and $BOLD. + use_color : bool, optional + If True, insert color sequences; otherwise strip them. + Default is True. + + Returns + ------- + str + Message with placeholders substituted. + """ if use_color: message = message.replace("$RESET", RESET_SEQ).replace("$BOLD", BOLD_SEQ) else: @@ -46,11 +74,22 @@ def formatter_message(message, use_color = True): } class ColoredFormatter(logging.Formatter): - def __init__(self, msg, use_color = True): + """Formatter that optionally colors the levelname in log records.""" + + def __init__(self, msg, use_color=True): + """ + Parameters + ---------- + msg : str + Format string for the formatter. + use_color : bool, optional + Whether to colorize level names. Default is True. + """ logging.Formatter.__init__(self, msg) self.use_color = use_color def format(self, record): + """Format the log record; optionally colorize levelname.""" levelname = record.levelname if self.use_color and levelname in COLORS: # Change color of levelname only if it is requested @@ -64,16 +103,21 @@ def format(self, record): pass return logging.Formatter.format(self, record) -# Custom logger class with multiple destinations class ColoredLogger(logging.Logger): - # Define formats for stream and file logging + """Logger that writes to both console (colored) and a file (UTC, no color).""" + ST_FMT = "[$BOLD%(filename)s::%(lineno)d$RESET] [%(levelname)s] %(message)s" F_FMT = "[$BOLD%(asctime)s::%(filename)s::%(lineno)d$RESET] [%(levelname)s] %(message)s" STREAM_FORMAT = formatter_message(ST_FMT, True) FILE_FORMAT = formatter_message(F_FMT, False) - # initialize logger def __init__(self, filename): - + """Create logger with stream (colored) and file (UTC) handlers. + + Parameters + ---------- + filename : str + Path to the log file. File is opened in 'w+' mode. + """ # Set logging level logging.Logger.__init__(self, None, logging.INFO) @@ -93,6 +137,17 @@ def __init__(self, filename): self.addHandler(filehandler) return + def close(self): + """Close and remove all handlers. + + Call when the logger is no longer needed (e.g., end of test or script) + to avoid ResourceWarnings from unclosed file handles. + """ + for handler in self.handlers[:]: + handler.close() + self.removeHandler(handler) + def shutdown(self): + """Shut down the logging system and flush all handlers.""" logging.shutdown() return diff --git a/potpyri/utils/options.py b/potpyri/utils/options.py index f0ca6588..abcecd09 100755 --- a/potpyri/utils/options.py +++ b/potpyri/utils/options.py @@ -1,11 +1,12 @@ -"Functions for initializing pipeline options and data paths." -"Authors: Charlie Kilpatrick" +"""Pipeline option parsing and path setup. -# Initial version tracking on 09/21/2024 -__version__ = "1.0" +Provides argument parsing for the main pipeline and builds directory layouts +(raw, red, log, cals, workspace, etc.) for a given data path and instrument. +Authors: Charlie Kilpatrick. +""" +from potpyri._version import __version__ import os -import importlib import shutil import sys import subprocess @@ -13,6 +14,14 @@ from potpyri import instruments def init_options(): + """Build and return the argument parser for the main pipeline. + + Returns + ------- + argparse.ArgumentParser + Parser with instrument, data_path, target, proc, include-bad, + file-list-name, photometry, and processing flags. + """ import argparse params = argparse.ArgumentParser(description='Path of data.') params.add_argument('instrument', @@ -94,6 +103,15 @@ def init_options(): def add_options(): + """Parse command-line options and return normalized args. + + Handles instrument aliases (e.g. BINO -> BINOSPEC, MMIR -> MMIRS). + + Returns + ------- + argparse.Namespace + Parsed and normalized command-line arguments. + """ params = init_options() args = params.parse_args() @@ -106,7 +124,14 @@ def add_options(): return(args) def test_for_dependencies(): - + """Check that astrometry.net and Source Extractor are on the system path. + + Raises + ------ + Exception + If solve-field (astrometry.net) or sex (sextractor) is not found + or does not report expected help text. + """ p = subprocess.run(['solve-field','-h'], capture_output=True) data = p.stdout.decode().lower() @@ -125,7 +150,31 @@ def test_for_dependencies(): https://github.com/astromatic/sextractor.''') def add_paths(data_path, file_list_name, tel): - + """Build the directory and path dictionary for a reduction run. + + Creates raw, bad, red, log, cal, work dirs under data_path and resolves + the path to the Source Extractor binary. + + Parameters + ---------- + data_path : str + Top-level data directory. + file_list_name : str + Filename of the file list (e.g. files.txt). + tel : Instrument + Instrument instance (used for name and caldb). + + Returns + ------- + dict + Keys include 'data', 'raw', 'bad', 'red', 'log', 'cal', 'work', + 'filelist', 'caldb', 'source_extractor'. + + Raises + ------ + Exception + If data_path does not exist. + """ if not os.path.exists(data_path): raise Exception(f'Data path does not exist: {data_path}') @@ -152,12 +201,25 @@ def add_paths(data_path, file_list_name, tel): return(paths) def initialize_telescope(instrument, data_path): - - module = importlib.import_module(f'potpyri.instruments.{instrument.upper()}') - tel = getattr(module, instrument.upper())() - - # Generate code and data paths based on input path - paths = add_paths(data_path, tel) + """Load the instrument class and build paths for that instrument. + + Parameters + ---------- + instrument : str + Instrument name (e.g. 'GMOS', 'LRIS'). + data_path : str + Top-level data directory. + + Returns + ------- + tuple + (paths, tel) where paths is the dict from add_paths and tel is the + instrument instance. Uses default file list name. + """ + tel = instruments.instrument_getter(instrument) + + # Generate code and data paths based on input path (default file list name) + paths = add_paths(data_path, 'files.txt', tel) return(paths, tel) diff --git a/potpyri/utils/utilities.py b/potpyri/utils/utilities.py index 5098cd2c..7bf6af9e 100755 --- a/potpyri/utils/utilities.py +++ b/potpyri/utils/utilities.py @@ -1,3 +1,8 @@ +"""General utilities for catalogs, coordinates, and numeric parsing. + +Used by absolute photometry and other steps that need Vizier catalog IDs +or coordinate handling. +""" from astropy import units as u from astropy.coordinates import SkyCoord import warnings @@ -34,26 +39,29 @@ }, } -# Sort the calibration files: def find_catalog(catalog, fil, coord_ra, coord_dec): - ''' + """Return Vizier catalog ID and column names for the given catalog and filter. - Function to return catalog ID for Vizier query. + Supports SDSS, 2MASS, UKIRT, PS1, SKYMAPPER. For southern u-band, uses + SkyMapper automatically. Parameters ---------- - - :param catalog: str - Name of catalog. + catalog : str + Catalog name (e.g. 'PS1', 'SDSS', '2MASS'). + fil : str + Filter band (e.g. 'r', 'g', 'J'). + coord_ra : float + Right ascension (used for catalog selection). + coord_dec : float + Declination (used for catalog selection; <0 can trigger SkyMapper for u-band). Returns ------- - - :return: list (str) - Catalog ID and column information to from Vizier. - - ''' - + tuple + (catalog, catalog_ID, ra_col, dec_col, mag_col, err_col) for use in + Vizier queries. catalog_ID/ra/dec/mag/err may be None if filter not supported. + """ catalog_ID, ra, dec, mag, err = None, None, None, None, None # If declination is less than 0 and filter is u-band, use SkyMapper @@ -81,13 +89,42 @@ def find_catalog(catalog, fil, coord_ra, coord_dec): return(catalog, catalog_ID, ra, dec, mag, err) def is_number(num): + """Return True if the value can be interpreted as a number. + + Parameters + ---------- + num : str or number + Value to test. + + Returns + ------- + bool + True if float(num) succeeds, False otherwise. + """ try: num = float(num) except ValueError: return(False) return(True) + def parse_coord(ra, dec): + """Parse RA and Dec strings into an astropy SkyCoord. + + Accepts decimal degrees or sexagesimal (e.g. '12:30:00' or '12.5'). + + Parameters + ---------- + ra : str or float + Right ascension. + dec : str or float + Declination. + + Returns + ------- + SkyCoord or None + ICRS coordinate, or None if parsing fails. + """ if (not (is_number(ra) and is_number(dec)) and (':' not in ra and ':' not in dec)): error = 'ERROR: cannot interpret: {ra} {dec}' diff --git a/pyproject.toml b/pyproject.toml index 8b06130d..9e7fd8c9 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -15,6 +15,8 @@ build-backend = "setuptools.build_meta" version_file = "potpyri/_version.py" version_scheme = "guess-next-dev" local_scheme = "no-local-version" +# Use when building from sdist or outside git (e.g. CI, Linux package managers) +fallback_version = "0.0.1" [project] @@ -23,7 +25,8 @@ dynamic = ["version"] description = "Data reduction pipeline for imaging from large aperture telescopes." readme = "README.md" requires-python = ">=3.11.0" -license = {file = "LICENSE.txt"} +license = "MIT" +license-files = ["LICENSE.txt"] authors = [ {name = "Kerry Paterson", email = "paterson@mpia.de" }, @@ -38,26 +41,26 @@ maintainers = [ classifiers = [ "Development Status :: 1 - Planning", "Topic :: Scientific/Engineering :: Astronomy", - "License :: OSI Approved :: MIT License", + "Operating System :: OS Independent", "Programming Language :: Python :: 3.12", "Programming Language :: Python :: 3.11", "Programming Language :: Python :: 3 :: Only", ] dependencies = [ - "acstools", - "astropy", - "astroquery", - "beautifulsoup4", - "ccdproc", - "gdown", - "matplotlib", - "numpy", - "photutils", - "progressbar", - "pykoa", - "requests", - "scipy" + "acstools>=3.0", + "astropy>=6.0", + "astroquery>=0.4", + "beautifulsoup4>=4.9", + "ccdproc>=2.4", + "gdown>=4.0", + "matplotlib>=3.5", + "numpy>=1.24", + "photutils>=1.9", + "progressbar>=2.5", + "pykoa>=1.0", + "requests>=2.25", + "scipy>=1.10" ] [tool.conda-lock.dependencies] @@ -65,7 +68,7 @@ astrometry = ">=0.1" astromatic-source-extractor = ">=1.0" [project.optional-dependencies] -test = ["pytest", "pytest-cov", "requests"] +test = ["pytest>=7", "pytest-cov>=4"] docs = ["sphinx==8.1.3", "nbsphinx", "sphinx-astropy", "pandoc", "sphinx-multiversion", "sphinxawesome-theme"] [project.urls] @@ -79,6 +82,8 @@ include-package-data = true where = ["."] include = ["potpyri*"] namespaces = false +# Exclude data dirs (staticmasks, cal) so they are only package-data, not Python packages +exclude = ["potpyri.data*"] [tool.setuptools.package-data] potpyri = ["data/staticmasks/*.fz", "data/cal/*/*.fz"] @@ -91,3 +96,26 @@ download_gemini_data = "potpyri.scripts.archives.download_gemini_data:main" [tool.pytest.ini_options] testpaths = ["tests"] +markers = [ + "integration: marks tests that require network/external services (deselect with '-m \"not integration\"')", +] +filterwarnings = [ + "ignore:.*scipy.odr.*deprecated:DeprecationWarning", + "ignore:Units from inserted quantities will be ignored:UserWarning", + "ignore:The EPSFModel class is deprecated:DeprecationWarning", + "ignore:The fit may be unsuccessful:UserWarning", + "ignore:.*findChildren.*:DeprecationWarning", + "ignore:cdelt will be ignored since cd is present:RuntimeWarning", + "ignore:Some non-standard WCS keywords were excluded:astropy.utils.exceptions.AstropyWarning", + # Calibration / WCS from test data (FITS with legacy headers) + "ignore:RADECSYS keyword is deprecated", + "ignore:datfix", + "ignore:DATE-OBS.*MJD-OBS", + # Sigma clipping on data that may contain NaNs/infs (e.g. calibration) + "ignore:Input data contains invalid values:astropy.utils.exceptions.AstropyUserWarning", + # Photometry / photutils (PSF at image edges, fit convergence, FITS header metadata) + "ignore:.*star\\(s\\) were not extracted because their cutout:astropy.utils.exceptions.AstropyUserWarning", + "ignore:One or more fit\\(s\\) may not have converged:astropy.utils.exceptions.AstropyUserWarning", + "ignore:Attribute.*cannot be added to FITS Header:astropy.utils.exceptions.AstropyUserWarning", + "ignore:Keyword name.*greater than 8 characters", +] diff --git a/tests/__init__.py b/tests/__init__.py index e69de29b..15d77e93 100644 --- a/tests/__init__.py +++ b/tests/__init__.py @@ -0,0 +1 @@ +"""POTPyRI unit tests: instruments, calibration, image_procs, sort_files, WCS, photometry, absphot.""" diff --git a/tests/test_absphot.py b/tests/test_absphot.py index e7849bda..a7f6905c 100644 --- a/tests/test_absphot.py +++ b/tests/test_absphot.py @@ -1,27 +1,29 @@ -from potpyri.utils import options -from potpyri.utils import logger -from potpyri.primitives import absphot -from potpyri.instruments import instrument_getter +"""Tests for absolute photometry zeropoint (find_zeropoint, Vizier/PS1) and absphot class methods.""" +import os -from astropy.io import ascii +import numpy as np +import pytest +import requests +from astropy.io import ascii, fits from astropy.table import Table from astropy.coordinates import SkyCoord -import os -import numpy as np - -from astropy.io import fits +from potpyri.utils import options, logger +from potpyri.primitives import absphot +from potpyri.instruments import instrument_getter from tests.utils import download_gdrive_file -def test_absphot(tmp_path): +@pytest.mark.integration +def test_absphot(tmp_path): + """Run find_zeropoint on LRIS stack and check ZPTMAG, ZPTCAT, ZPTMUCER in header.""" instrument = 'LRIS' file_list_name = 'files.txt' - # Raw science file - file_path = download_gdrive_file('LRISr/red/R155_host.R.ut240629.1R.11.stk.fits.fz', use_cached=True) - cat_path = download_gdrive_file('LRISr/red/test_catalog.dat', use_cached=True) + # Raw science file (download to tmp_path so log is writable) + file_path = download_gdrive_file('LRISr/red/R155_host.R.ut240629.1R.11.stk.fits.fz', output_dir=str(tmp_path), use_cached=True) + cat_path = download_gdrive_file('LRISr/red/test_catalog.dat', output_dir=str(tmp_path), use_cached=True) input_catalog = ascii.read(cat_path) data_path, basefile = os.path.split(file_path) @@ -32,13 +34,80 @@ def test_absphot(tmp_path): log = logger.get_log(paths['log']) hdu = fits.open(file_path) + hdu.close() - absphot.find_zeropoint(file_path, tel, log=log) + try: + absphot.find_zeropoint(file_path, tel, log=log) + except (requests.exceptions.ConnectionError, OSError) as e: + pytest.skip(f"VizieR unreachable (zeropoint needs PS1): {e}") + finally: + log.close() hdu = fits.open(file_path) header = hdu['PRIMARY'].header + hdu.close() assert header['ZPTCAT']=='PS1' assert header['FILTER']=='R' assert np.abs(27.576871-header['ZPTMAG'])<0.01 assert header['ZPTMUCER']<0.01 + + +def test_get_zeropoint(): + """get_zeropoint returns (zpt, zpterr) consistent with mag = zpt - 2.5*log10(flux).""" + cal = absphot.absphot() + flux = np.array([1000.0, 2000.0, 500.0]) + fluxerr = np.array([10.0, 14.0, 7.0]) + mag = np.array([25.0, 24.25, 25.75]) + magerr = np.array([0.01, 0.01, 0.01]) + zpt, zpterr = cal.get_zeropoint(flux, fluxerr, mag, magerr) + assert np.isfinite(zpt) and np.isfinite(zpterr) + mag_derived = zpt - 2.5 * np.log10(flux) + np.testing.assert_allclose(mag, mag_derived, atol=0.1) + + +def test_zpt_iteration(): + """zpt_iteration returns (zpt, zpterr, master_mask) with mask length = input length.""" + cal = absphot.absphot(iterations=2, sigma=3.0) + np.random.seed(42) + n = 20 + flux = np.random.uniform(500, 3000, n).astype(float) + fluxerr = flux * 0.01 + zpt_true = 26.0 + mag = zpt_true - 2.5 * np.log10(flux) + np.random.normal(0, 0.02, n) + magerr = np.full(n, 0.02) + zpt, zpterr, master_mask = cal.zpt_iteration(flux, fluxerr, mag, magerr, log=None) + assert np.isfinite(zpt) and np.isfinite(zpterr) + assert len(master_mask) == n + assert np.sum(master_mask) <= n + + +def test_convert_filter_name(): + """convert_filter_name maps instrument filter names to catalog names.""" + cal = absphot.absphot() + assert cal.convert_filter_name('gG0301') == 'g' + assert cal.convert_filter_name('rG0303') == 'r' + assert cal.convert_filter_name('RG850') == 'z' + assert cal.convert_filter_name('Y') == 'J' + assert cal.convert_filter_name('V') == 'g' + assert cal.convert_filter_name('I') == 'i' + + +def test_get_minmag(): + """get_minmag returns bright limit by filter.""" + cal = absphot.absphot() + assert cal.get_minmag('J') == 15.5 + assert cal.get_minmag('K') == 13.0 + assert cal.get_minmag('Y') == 15.0 + assert cal.get_minmag('g') == 16.0 + + +def test_Y_band(): + """Y_band returns (Y, Y_err) from J and K relation Y = J + 0.46*(J-K).""" + cal = absphot.absphot() + J, J_err = 15.0, 0.02 + K, K_err = 14.0, 0.02 + Y, Y_err = cal.Y_band(J, J_err, K, K_err) + expected_Y = J + 0.46 * (J - K) + np.testing.assert_allclose(Y, expected_Y) + assert np.isfinite(Y_err) and Y_err > 0 diff --git a/tests/test_calibration.py b/tests/test_calibration.py index 31ad4685..9b1a1a14 100644 --- a/tests/test_calibration.py +++ b/tests/test_calibration.py @@ -1,3 +1,4 @@ +"""Tests for master bias/flat and process_science calibration.""" from potpyri.utils import options from potpyri.utils import logger from potpyri.instruments import instrument_getter @@ -5,17 +6,20 @@ import os import numpy as np +import pytest from tests.utils import download_gdrive_file -def test_cal(tmp_path): +@pytest.mark.integration +def test_cal(tmp_path): + """Run bias/flat and process_science on GMOS raw; compare to manual (sci-bias)/flat.""" instrument = 'GMOS' file_list_name = 'files.txt' - # Raw science file plus calibration files - file_path = download_gdrive_file('GMOS/raw/N20240618S0015.fits.bz2', use_cached=True) - download_gdrive_file('GMOS/red/cals/mbias_12_22.fits.fz', use_cached=True) - download_gdrive_file('GMOS/red/cals/mflat_i_12_22.fits.fz', use_cached=True) + # Raw science file plus calibration files (download to tmp_path so log is writable) + file_path = download_gdrive_file('GMOS/raw/N20240618S0015.fits.bz2', output_dir=str(tmp_path), use_cached=True) + download_gdrive_file('GMOS/red/cals/mbias_12_22.fits.fz', output_dir=str(tmp_path), use_cached=True) + download_gdrive_file('GMOS/red/cals/mflat_i_12_22.fits.fz', output_dir=str(tmp_path), use_cached=True) data_path, basefile = os.path.split(file_path) data_path, _ = os.path.split(data_path) @@ -33,7 +37,10 @@ def test_cal(tmp_path): data = (sci_full-mbias.data)/(mflat.data/np.nanmean(mflat.data)) data = data.astype(np.float32) - processed = tel.process_science([file_path], 'i', '12', '22', paths, - mbias=mbias, mflat=mflat, mdark=None, skip_skysub=True, log=log) - + try: + processed = tel.process_science([file_path], 'i', '12', '22', paths, + mbias=mbias, mflat=mflat, mdark=None, skip_skysub=True, log=log) + finally: + log.close() + np.testing.assert_array_equal(processed[0].data[~np.isnan(processed[0].data)], data[~np.isnan(processed[0].data)]) diff --git a/tests/test_error.py b/tests/test_error.py index 405205b0..6db76293 100644 --- a/tests/test_error.py +++ b/tests/test_error.py @@ -1,3 +1,4 @@ +"""Tests for create_error (per-pixel error array from science and mask).""" from potpyri.utils import options from potpyri.utils import logger from potpyri.primitives import image_procs @@ -15,7 +16,7 @@ from . import utils def test_error(tmp_path): - + """Create error array for synthetic image and compare to theoretical sqrt(poisson+rms^2+rdnoise).""" instrument = 'MOSFIRE' file_list_name = 'files.txt' @@ -42,7 +43,10 @@ def test_error(tmp_path): hdu.header['SATURATE'] = 35000.0 hdu.writeto(outfile, overwrite=True, output_verify='silentfix') - error_hdu = image_procs.create_error(outfile, maskhdu, rdnoise) + try: + error_hdu = image_procs.create_error(outfile, maskhdu, rdnoise) + finally: + log.close() rms = 0.5 * ( np.percentile(imdata[~maskim.astype(bool)], 84.13) diff --git a/tests/test_image_procs.py b/tests/test_image_procs.py new file mode 100644 index 00000000..40f0475c --- /dev/null +++ b/tests/test_image_procs.py @@ -0,0 +1,105 @@ +"""Tests for image_procs: remove_pv_distortion, get_fieldcenter, generate_wcs, detrend_stack, add_stack_mask.""" +from potpyri.primitives import image_procs +from potpyri.instruments import instrument_getter + +import numpy as np +import os + +from astropy.io import fits +from astropy.wcs import WCS +from ccdproc import CCDData +from astropy import units as u + +from tests.utils import make_dummy_wcs + + +def test_remove_pv_distortion(tmp_path): + """remove_pv_distortion strips PV* keywords from header in place.""" + header = fits.Header() + header['NAXIS'] = 2 + header['PV1_1'] = 1.0 + header['PV2_1'] = 0.0 + header['CTYPE1'] = 'RA---TAN' + out = image_procs.remove_pv_distortion(header) + assert out is header + assert 'PV1_1' not in header + assert 'PV2_1' not in header + assert header['CTYPE1'] == 'RA---TAN' + + +def test_get_fieldcenter(tmp_path): + """get_fieldcenter returns mean RA/Dec of image centers from FITS with WCS.""" + wcs_head = make_dummy_wcs() + wcs_head['NAXIS1'] = 100 + wcs_head['NAXIS2'] = 100 + data = np.zeros((100, 100), dtype=np.float32) + for i in range(2): + hdu = fits.PrimaryHDU(data) + hdu.header.update(wcs_head) + path = os.path.join(tmp_path, f'img{i}.fits') + hdu.writeto(path, overwrite=True) + images = [os.path.join(tmp_path, 'img0.fits'), os.path.join(tmp_path, 'img1.fits')] + center = image_procs.get_fieldcenter(images) + assert len(center) == 2 + assert np.isfinite(center[0]) and np.isfinite(center[1]) + # CRVAL from make_dummy_wcs is 30, -30 + np.testing.assert_allclose(center, [30.0, -30.0], atol=0.01) + + +def test_generate_wcs(tmp_path): + """generate_wcs returns a TAN WCS with correct CRVAL and size.""" + tel = instrument_getter('GMOS') + binn = '22' + fieldcenter = [30.0, -30.0] + out_size = 256 + wcs = image_procs.generate_wcs(tel, binn, fieldcenter, out_size) + assert isinstance(wcs, WCS) + np.testing.assert_allclose(wcs.wcs.crval, [30.0, -30.0]) + assert wcs.pixel_shape == (out_size, out_size) + + +def test_detrend_stack(tmp_path): + """detrend_stack removes row/column median and updates SATURATE.""" + sci = np.ones((64, 64), dtype=np.float32) * 100.0 + sci += np.arange(64)[:, None] # row gradient + sci += np.arange(64)[None, :] # column gradient + mask = np.zeros((64, 64), dtype=np.uint8) + hdu_sci = fits.PrimaryHDU(sci.astype(np.float32)) + hdu_sci.header['SATURATE'] = 50000.0 + hdu_mask = fits.ImageHDU(mask) + stack = fits.HDUList([hdu_sci, hdu_mask]) + out = image_procs.detrend_stack(stack) + assert out is stack + assert stack[0].data is not None + # After detrend, row/col median should be ~0 + assert np.abs(np.nanmedian(stack[0].data)) < 1.0 + assert stack[0].header['SATURATE'] < 50000.0 + + +def test_add_stack_mask(tmp_path): + """add_stack_mask merges per-image masks and saturation into stack mask.""" + from potpyri.utils import logger + from potpyri.utils import options + tel = instrument_getter('MMIRS') + paths = options.add_paths(str(tmp_path), 'files.txt', tel) + log = logger.get_log(paths['log']) + try: + shape = (32, 32) + data = np.ones(shape, dtype=np.float32) * 10.0 + mask = np.zeros(shape, dtype=np.uint8) + err = np.ones(shape, dtype=np.float32) * 0.5 + wcs_head = make_dummy_wcs() + ccd = CCDData(data, meta={'EXPTIME': 1.0, 'SATURATE': 1000.0}, wcs=WCS(wcs_head), unit=u.electron) + stacking_data = [ccd] + stack_sci = fits.PrimaryHDU(data.copy()) + stack_sci.header['SATURATE'] = 1000.0 + stack_mask = fits.ImageHDU(mask.copy()) + stack_err = fits.ImageHDU(err.copy()) + stack = fits.HDUList([stack_sci, stack_mask, stack_err]) + image_procs.add_stack_mask(stack, stacking_data) + assert stack[1].data is not None + assert stack[0].data is not None + # Saturation 1000, data 10 -> no sat mask; bad pixel mask could be 0 + assert np.any(stack[1].data >= 0) + finally: + log.close() diff --git a/tests/test_init.py b/tests/test_init.py index e7b19fa2..09fdc690 100644 --- a/tests/test_init.py +++ b/tests/test_init.py @@ -1,8 +1,10 @@ +"""Tests for instrument initialization (instrument_getter for all supported instruments).""" from potpyri.utils import options from potpyri.instruments import instrument_getter -def test_init(tmp_path): +def test_init(tmp_path): + """Check that instrument_getter returns an instance with correct name for each CLI instrument.""" # Parse allowed instruments directly from options so this test matches # to values that the argument parser allows as input instruments = [] diff --git a/tests/test_photometry.py b/tests/test_photometry.py index 1e64d230..11791f53 100644 --- a/tests/test_photometry.py +++ b/tests/test_photometry.py @@ -1,3 +1,6 @@ +"""Tests for photometry.photloop (PSF and aperture photometry on stacks) and helpers (create_conv, create_params, extract_*, get_star_catalog, extract_fwhm_from_epsf).""" +import warnings + from potpyri.utils import options from potpyri.utils import logger from potpyri.primitives import photometry @@ -7,16 +10,21 @@ import numpy as np from astropy.io import fits +from astropy.table import Table +from photutils.psf import EPSFModel +import pytest from tests.utils import download_gdrive_file -def test_photometry(tmp_path): +@pytest.mark.integration +def test_photometry(tmp_path): + """Run photloop on GMOS stack; assert APPPHOT, PSFPHOT, PSFSTARS, PSF extensions and catalog size.""" instrument = 'GMOS' file_list_name = 'files.txt' - # Raw science file - file_path = download_gdrive_file('GMOS/red/sGRB240615A-GRB.i.ut240618.12.22.stk.fits.fz', use_cached=True) + # Raw science file (download to tmp_path so log and overwrites are writable) + file_path = download_gdrive_file('GMOS/red/sGRB240615A-GRB.i.ut240618.12.22.stk.fits.fz', output_dir=str(tmp_path), use_cached=True) # Strip out extensions added by photometry loop if they exist hdu = fits.open(file_path) @@ -46,14 +54,117 @@ def test_photometry(tmp_path): phot_sn_max = 100.0 fwhm_init = 5.0 - photometry.photloop(file_path, phot_sn_min=phot_sn_min, - fwhm_init=fwhm_init, phot_sn_max=phot_sn_max, log=log) - - hdu = fits.open(file_path) - - assert np.any(['PSFPHOT' in h.name for h in hdu]) - assert np.any(['PSFSTARS' in h.name for h in hdu]) - assert np.any(['APPPHOT' in h.name for h in hdu]) - assert np.any(['PSF' in h.name for h in hdu]) - - assert len(hdu['APPPHOT'].data)>100 + try: + photometry.photloop(file_path, phot_sn_min=phot_sn_min, + fwhm_init=fwhm_init, phot_sn_max=phot_sn_max, log=log) + finally: + log.close() + + with fits.open(file_path) as hdu: + assert np.any(['PSFPHOT' in h.name for h in hdu]) + assert np.any(['PSFSTARS' in h.name for h in hdu]) + assert np.any(['APPPHOT' in h.name for h in hdu]) + assert np.any(['PSF' in h.name for h in hdu]) + assert len(hdu['APPPHOT'].data) > 100 + + +def test_create_conv(tmp_path): + """create_conv writes a 3x3 CONV NORM kernel file.""" + outfile = os.path.join(tmp_path, 'test.conv') + photometry.create_conv(outfile) + assert os.path.exists(outfile) + with open(outfile) as f: + content = f.read() + assert 'CONV NORM' in content + assert '1 2 1' in content + assert '2 4 2' in content + + +def test_create_params(tmp_path): + """create_params writes parameter file with NUMBER, X_IMAGE, FWHM_IMAGE, etc.""" + outfile = os.path.join(tmp_path, 'test.param') + photometry.create_params(outfile) + assert os.path.exists(outfile) + with open(outfile) as f: + content = f.read() + assert 'NUMBER' in content + assert 'X_IMAGE' in content + assert 'FWHM_IMAGE' in content + assert 'MAG_AUTO' in content + + +def test_extract_aperture_stats(tmp_path): + """extract_aperture_stats returns table with fwhm, flux_best, Xpos, Ypos columns.""" + from photutils.aperture import CircularAperture, ApertureStats + shape = (64, 64) + img_data = np.ones(shape, dtype=float) * 100.0 + img_data[32, 32] = 500.0 + img_mask = np.zeros(shape, dtype=bool) + img_error = np.ones(shape, dtype=float) * 5.0 + stars = Table() + stars['xcentroid'] = [32.0] + stars['ycentroid'] = [32.0] + tel = instrument_getter('MMIRS') + paths = options.add_paths(str(tmp_path), 'files.txt', tel) + log = logger.get_log(paths['log']) + try: + with warnings.catch_warnings(): + warnings.filterwarnings('ignore', message='.*Units from inserted quantities.*', category=UserWarning) + stats = photometry.extract_aperture_stats( + img_data, img_mask, img_error, stars, + aperture_radius=5.0, log=log) + assert len(stats) == 1 + assert 'fwhm' in stats.colnames + assert 'flux_best' in stats.colnames + assert 'Xpos' in stats.colnames + assert 'Ypos' in stats.colnames + finally: + log.close() + + +def test_get_star_catalog(tmp_path): + """get_star_catalog finds sources and returns table with aperture stats.""" + shape = (128, 128) + y, x = np.ogrid[:shape[0], :shape[1]] + img_data = np.ones(shape, dtype=float) * 50.0 + # Add a clear Gaussian star so DAOStarFinder and sharpness/roundness pass + for cy, cx in [(64, 64), (32, 96)]: + img_data += 2000.0 * np.exp(-((x - cx)**2 + (y - cy)**2) / (2 * 9.0)) + img_mask = np.zeros(shape, dtype=bool) + img_error = np.ones(shape, dtype=float) * 5.0 + tel = instrument_getter('MMIRS') + paths = options.add_paths(str(tmp_path), 'files.txt', tel) + log = logger.get_log(paths['log']) + try: + with warnings.catch_warnings(): + warnings.filterwarnings('ignore', message='.*Units from inserted quantities.*', category=UserWarning) + stars = photometry.get_star_catalog( + img_data, img_mask, img_error, + fwhm_init=4.0, threshold=10.0, log=log) + assert stars is not None + assert len(stars) >= 1 + assert 'xcentroid' in stars.colnames or 'Xpos' in stars.colnames + assert 'flux_best' in stars.colnames or 'flux' in stars.colnames + finally: + log.close() + + +def test_extract_fwhm_from_epsf(): + """extract_fwhm_from_epsf returns finite FWHM from EPSFModel.""" + # Build a minimal EPSF-like array (Gaussian blob) + size = 15 + y, x = np.ogrid[-size//2:size//2+1, -size//2:size//2+1] + sigma = 2.0 + data = np.exp(-(x*x + y*y) / (2 * sigma**2)).astype(float) + with warnings.catch_warnings(): + try: + from astropy.utils.exceptions import AstropyDeprecationWarning, AstropyUserWarning + warnings.simplefilter('ignore', AstropyDeprecationWarning) + warnings.simplefilter('ignore', AstropyUserWarning) + except ImportError: + warnings.simplefilter('ignore', DeprecationWarning) + warnings.simplefilter('ignore', UserWarning) + epsf = EPSFModel(data) + fwhm = photometry.extract_fwhm_from_epsf(epsf, fwhm_init=3.0) + assert np.isfinite(fwhm) + assert fwhm > 0 diff --git a/tests/test_sort.py b/tests/test_sort.py index bead1811..21d471d5 100644 --- a/tests/test_sort.py +++ b/tests/test_sort.py @@ -1,3 +1,4 @@ +"""Tests for sort_files.handle_files and file table content (BINOSPEC), and classifiers (is_bad, is_spec, is_flat, is_dark, is_bias, is_science).""" from potpyri.utils import options from potpyri.utils import logger from potpyri.primitives import sort_files @@ -5,15 +6,20 @@ import os +from astropy.io import fits + +import pytest from tests.utils import download_gdrive_file -def test_sort(tmp_path): +@pytest.mark.integration +def test_sort(tmp_path): + """Run handle_files on BINOSPEC proc file; check file table rows and no_redo reuse.""" instrument = 'BINOSPEC' file_list_name = 'files.txt' - # Raw science file - file_path = download_gdrive_file('Binospec/raw/sci_img_2024.0812.034220_proc.fits.fz', use_cached=True) + # Raw science file (download to tmp_path so log is writable) + file_path = download_gdrive_file('Binospec/raw/sci_img_2024.0812.034220_proc.fits.fz', output_dir=str(tmp_path), use_cached=True) data_path, basefile = os.path.split(file_path) data_path, _ = os.path.split(data_path) @@ -24,23 +30,26 @@ def test_sort(tmp_path): log = logger.get_log(paths['log']) # This contains all of the file data - file_table = sort_files.handle_files(paths['filelist'], paths, tel, - incl_bad=True, proc='proc', no_redo=False, log=log) - - # Run validation checks on file_table - assert len(file_table)==1 - assert file_table[0]['Target']=='GRB240809A_r_ep1' - assert file_table[0]['Filter']=='r' - assert file_table[0]['Type']=='SCIENCE' - assert file_table[0]['TargType']=='GRB240809A_r_ep1_r_2_11' - assert file_table[0]['CalType']=='r_2_11' - assert file_table[0]['Exp']=="120.0" - assert file_table[0]['Binning']=="11" - assert file_table[0]['Amp']=="2" - - # Test the no_redo flag - new_file_table = sort_files.handle_files(paths['filelist'], paths, tel, - incl_bad=True, proc='proc', no_redo=True, log=log) + try: + file_table = sort_files.handle_files(paths['filelist'], paths, tel, + incl_bad=True, proc='proc', no_redo=False, log=log) + + # Run validation checks on file_table + assert len(file_table)==1 + assert file_table[0]['Target']=='GRB240809A_r_ep1' + assert file_table[0]['Filter']=='r' + assert file_table[0]['Type']=='SCIENCE' + assert file_table[0]['TargType']=='GRB240809A_r_ep1_r_2_11' + assert file_table[0]['CalType']=='r_2_11' + assert file_table[0]['Exp']=="120.0" + assert file_table[0]['Binning']=="11" + assert file_table[0]['Amp']=="2" + + # Test the no_redo flag + new_file_table = sort_files.handle_files(paths['filelist'], paths, tel, + incl_bad=True, proc='proc', no_redo=True, log=log) + finally: + log.close() assert len(new_file_table)==1 assert file_table[0]['Target']==new_file_table[0]['Target'] @@ -51,3 +60,77 @@ def test_sort(tmp_path): assert file_table[0]['Exp']==new_file_table[0]['Exp'] assert file_table[0]['Binning']==new_file_table[0]['Binning'] assert file_table[0]['Amp']==new_file_table[0]['Amp'] + + +def test_is_bad_binospec(): + """is_bad True when header matches bad_keywords/bad_values (BINOSPEC: MASK=mira).""" + tel = instrument_getter('BINOSPEC') + hdr = fits.Header() + hdr['MASK'] = 'mira' + # No CCDSUM so binning check does not overwrite; keyword match gives bad=True + assert sort_files.is_bad(hdr, tel) == True + hdr['MASK'] = 'imaging' + hdr['CCDSUM'] = '1,1' # valid binning so keyword result (False) is kept + assert sort_files.is_bad(hdr, tel) == False + + +def test_is_spec_binospec(): + """is_spec True when header matches spec keywords (BINOSPEC: MASK=spectroscopy).""" + tel = instrument_getter('BINOSPEC') + hdr = fits.Header() + hdr['MASK'] = 'spectroscopy' + assert sort_files.is_spec(hdr, tel) == True + hdr['MASK'] = 'imaging' + assert sort_files.is_spec(hdr, tel) == False + + +def test_is_flat_binospec(): + """is_flat True when MASK=imaging, SCRN=deployed (BINOSPEC).""" + tel = instrument_getter('BINOSPEC') + hdr = fits.Header() + hdr['MASK'] = 'imaging' + hdr['SCRN'] = 'deployed' + hdr['CCDSUM'] = '1,1' + assert sort_files.is_flat(hdr, tel) == True + hdr['SCRN'] = 'stowed' + assert sort_files.is_flat(hdr, tel) == False + + +def test_is_science_binospec(): + """is_science True when MASK=imaging, SCRN=stowed and exptime >= min (BINOSPEC).""" + tel = instrument_getter('BINOSPEC') + hdr = fits.Header() + hdr['MASK'] = 'imaging' + hdr['SCRN'] = 'stowed' + hdr['CCDSUM'] = '1,1' + hdr['EXPTIME'] = 120.0 + assert sort_files.is_science(hdr, tel) == True + hdr['SCRN'] = 'deployed' + assert sort_files.is_science(hdr, tel) == False + + +def test_is_bias_binospec(): + """is_bias False for BINOSPEC (no bias keywords).""" + tel = instrument_getter('BINOSPEC') + hdr = fits.Header() + assert sort_files.is_bias(hdr, tel) is False + + +def test_is_dark_binospec(): + """is_dark False for BINOSPEC (no dark keywords).""" + tel = instrument_getter('BINOSPEC') + hdr = fits.Header() + assert sort_files.is_dark(hdr, tel) is False + + +def test_is_bias_gmos(): + """is_bias True when SHUTTER=closed, OBSCLASS=daycal, OBSTYPE=bias (GMOS).""" + tel = instrument_getter('GMOS') + hdr = fits.Header() + hdr['SHUTTER'] = 'closed' + hdr['OBSCLASS'] = 'daycal' + hdr['OBSTYPE'] = 'bias' + hdr['CCDSUM'] = '2,2' + assert sort_files.is_bias(hdr, tel) == True + hdr['OBSTYPE'] = 'object' + assert sort_files.is_bias(hdr, tel) == False diff --git a/tests/test_stack.py b/tests/test_stack.py index d08c146f..18b25681 100644 --- a/tests/test_stack.py +++ b/tests/test_stack.py @@ -1,3 +1,4 @@ +"""Tests for image_procs.stack_data (median stack of synthetic MMIRS images).""" from potpyri.utils import options from potpyri.utils import logger from potpyri.primitives import image_procs @@ -15,7 +16,7 @@ from . import utils def test_stack(tmp_path): - + """Stack synthetic MMIRS images and assert median value in stacked science.""" instrument = 'MMIRS' file_list_name = 'files.txt' nims = 5 @@ -46,6 +47,9 @@ def test_stack(tmp_path): masks.append(fits.ImageHDU(maskim)) errors.append(fits.ImageHDU(errimg)) - sci_med = image_procs.stack_data(aligned_data, tel, masks, errors, log=log) + try: + sci_med = image_procs.stack_data(aligned_data, tel, masks, errors, log=log) + finally: + log.close() np.testing.assert_array_equal(sci_med[0].data, np.median(np.arange(nims)+1)) diff --git a/tests/test_staticmask.py b/tests/test_staticmask.py index 54737f58..de05401d 100644 --- a/tests/test_staticmask.py +++ b/tests/test_staticmask.py @@ -1,8 +1,10 @@ +"""Tests for load_staticmask for F2, GMOS-N/S, MOSFIRE, BINOSPEC.""" from potpyri.utils import options from potpyri.instruments import instrument_getter -def test_staticmask(tmp_path): +def test_staticmask(tmp_path): + """Load static mask for each instrument (F2, GMOS-N/S, MOSFIRE, BINOSPEC) and assert non-None.""" file_list_name = 'files.txt' instruments = { diff --git a/tests/test_wcs.py b/tests/test_wcs.py index 429e0fc0..88d2634a 100644 --- a/tests/test_wcs.py +++ b/tests/test_wcs.py @@ -1,51 +1,111 @@ -from potpyri.utils import options -from potpyri.utils import logger -from potpyri.primitives import solve_wcs -from potpyri.instruments import instrument_getter - +"""Tests for solve_wcs (solve_astrometry, align_to_gaia on GMOS stack) and helpers (clean_up_astrometry).""" import os +import numpy as np +import pytest +import requests from astropy.io import fits +from potpyri.utils import options, logger +from potpyri.primitives import solve_wcs +from potpyri.instruments import instrument_getter + from tests.utils import download_gdrive_file + def test_wcs(tmp_path): + """Check RADISP/DEDISP and scale from a FITS file (no network required). + Uses a minimal FITS file with RADISP/DEDISP set as align_to_gaia would + write them, so the test passes consistently without being skipped. + Full solve_astrometry + align_to_gaia flow is covered by test_wcs_integration. + """ + instrument = 'GMOS' + tel = instrument_getter(instrument) + # Dispersion values that satisfy the pipeline quality cuts (arcsec). + # Must be < 0.5 and < tel.pixscale so disp/pixscale < 1 (GMOS pixscale ≈ 0.08). + ra_disp_val = dec_disp_val = 0.05 + + # Create minimal FITS with RADISP/DEDISP as align_to_gaia writes them + file_path = os.path.join(tmp_path, 'test_solved.fits') + hdu = fits.PrimaryHDU(data=np.zeros((50, 50), dtype=np.float32)) + hdu.header['NAXIS1'] = 50 + hdu.header['NAXIS2'] = 50 + hdu.header['RADISP'] = (ra_disp_val, 'Dispersion in R.A. of WCS [Arcsec]') + hdu.header['DEDISP'] = (dec_disp_val, 'Dispersion in Decl. of WCS [Arcsec]') + hdu.writeto(file_path, overwrite=True) + + with fits.open(file_path) as hdu: + ra_disp = hdu[0].header['RADISP'] + dec_disp = hdu[0].header['DEDISP'] + + assert ra_disp < 0.5 + assert dec_disp < 0.5 + assert ra_disp / tel.pixscale < 1 + assert dec_disp / tel.pixscale < 1 + + +@pytest.mark.integration +def test_wcs_integration(tmp_path): + """Full pipeline: solve_astrometry and align_to_gaia on GMOS slice (requires network and astrometry index).""" instrument = 'GMOS' file_list_name = 'files.txt' - # Download science file - file_path = download_gdrive_file('GMOS/red/workspace/sGRB240615A-GRB.i.ut240618.12.22.1403185114.fits.fz', use_cached=True) - astm_path = download_gdrive_file('astrometry/index-52m1-14.fits') + file_path = download_gdrive_file( + 'GMOS/red/workspace/sGRB240615A-GRB.i.ut240618.12.22.1403185114.fits.fz', + output_dir=str(tmp_path), use_cached=True) + astm_path = download_gdrive_file( + 'astrometry/index-52m1-14.fits', output_dir=str(tmp_path), use_cached=True) - data_path, basefile = os.path.split(file_path) + data_path, _ = os.path.split(file_path) data_path, _ = os.path.split(data_path) - tel = instrument_getter(instrument) paths = options.add_paths(data_path, file_list_name, tel) - - # Generate log file in corresponding directory for log log = logger.get_log(paths['log']) - # Write out basic file to solve WCS on - hdu = fits.open(file_path) - binn = tel.get_binning(hdu[1].header) - file_path = file_path.replace('.fz','') - del hdu[1].header['RADISP'] - del hdu[1].header['DEDISP'] - fits.writeto(file_path,hdu[1].data,hdu[1].header,overwrite=True) + with fits.open(file_path) as hdu: + binn = tel.get_binning(hdu[1].header) + file_path_unfz = file_path.replace('.fz', '') + hdu[1].header.pop('RADISP', None) + hdu[1].header.pop('DEDISP', None) + fits.writeto(file_path_unfz, hdu[1].data, hdu[1].header, overwrite=True) + file_path = file_path_unfz - # Solve WCS - solve_wcs.solve_astrometry(file_path, tel, binn, paths, index=astm_path, log=log) - solve_wcs.align_to_gaia(file_path, tel, radius=0.5, log=log) + try: + solve_wcs.solve_astrometry(file_path, tel, binn, paths, index=astm_path, log=log) + solve_wcs.align_to_gaia(file_path, tel, radius=0.5, log=log) + except (requests.exceptions.ConnectionError, OSError) as e: + pytest.skip(f"VizieR/network unreachable (Gaia catalog): {e}") + finally: + log.close() - # Test dispersion - hdu = fits.open(file_path) - ra_disp = hdu[0].header['RADISP'] - dec_disp = hdu[0].header['DEDISP'] + with fits.open(file_path) as hdu: + ra_disp = hdu[0].header['RADISP'] + dec_disp = hdu[0].header['DEDISP'] assert ra_disp < 0.5 assert dec_disp < 0.5 + assert ra_disp / tel.pixscale < 1 + assert dec_disp / tel.pixscale < 1 + - assert ra_disp/tel.pixscale<1 - assert dec_disp/tel.pixscale<1 +def test_clean_up_astrometry(tmp_path): + """clean_up_astrometry removes .axy, .corr, .solved, .wcs etc. in directory.""" + base = 'test_image.fits' + exten = '.fits' + to_create = [ + base.replace(exten, '.axy'), + base.replace(exten, '.corr'), + base.replace(exten, '.solved'), + base.replace(exten, '.wcs'), + ] + for f in to_create: + path = os.path.join(tmp_path, f) + with open(path, 'w') as fp: + fp.write('dummy') + dir_path = tmp_path + file_path = os.path.join(dir_path, base) + solve_wcs.clean_up_astrometry(dir_path, base, exten) + for f in to_create: + path = os.path.join(tmp_path, f) + assert not os.path.exists(path), f'{f} should have been removed' diff --git a/tests/utils.py b/tests/utils.py index 0dabbaa5..05f4f5ef 100644 --- a/tests/utils.py +++ b/tests/utils.py @@ -1,3 +1,4 @@ +"""Test helpers: Google Drive and GitHub file download, dummy WCS for tests.""" import requests import re import os @@ -8,7 +9,7 @@ import itertools import warnings -__all__ = ["download_gdrive_file","make_dummy_wcs","download_github_file"] +__all__ = ["download_gdrive_file", "make_dummy_wcs", "download_github_file"] class _GoogleDriveFile: TYPE_FOLDER = "application/vnd.google-apps.folder" @@ -65,7 +66,9 @@ def download_gdrive_file( try: file_id = resolve_relative_path(shared_folder_url, rel_file_path) url = f"https://drive.google.com/uc?id={file_id}" - gdown.download(url, output_path, quiet=False) + with warnings.catch_warnings(): + warnings.filterwarnings("ignore", category=DeprecationWarning, message=".*findChildren.*") + gdown.download(url, output_path, quiet=False) except: print(f"Failed to download {rel_file_path} from Google Drive.") raise @@ -163,6 +166,7 @@ def _get_directory_structure(gdrive_file, previous_path): return match[0][0] # return file id def make_dummy_wcs(): + """Return a minimal FITS header dict with TAN WCS for test images.""" header = {'WCSAXES': 2, 'CTYPE1': 'RA---TAN', 'CTYPE2': 'DEC--TAN',