Merge branch 'main' into mefmess

mirar/catalog/base_catalog.py

@@ -10,7 +10,7 @@
 from mirar.data import Image
 from mirar.paths import BASE_NAME_KEY
 from mirar.processors.candidates.utils import get_image_center_wcs_coords
-from mirar.utils.ldac_tools import save_table_as_ldac
+from mirar.utils.ldac_tools import get_table_from_ldac, save_table_as_ldac
 
 logger = logging.getLogger(__name__)
 
@@ -65,7 +65,8 @@ def write_catalog(self, image: Image, output_dir: str | Path) -> Path:
         :param output_dir: output directory for catalog
         :return: path of catalog
         """
-
+        if isinstance(output_dir, str):
+            output_dir = Path(output_dir)
         ra_deg, dec_deg = get_image_center_wcs_coords(image, origin=1)
 
         base_name = Path(image[BASE_NAME_KEY]).with_suffix(".ldac").name
@@ -107,6 +108,9 @@ def catalog_name(self):
 
     @property
     def projection(self):
+        """
+        projection for kowalski xmatch
+        """
         raise NotImplementedError
 
     @property
@@ -136,3 +140,26 @@ def query(self, coords: dict) -> dict:
         :return: crossmatch
         """
         raise NotImplementedError
+
+
+class CatalogFromFile(BaseCatalog):
+    """
+    Local catalog from file
+    """
+
+    abbreviation = "local"
+
+    def __init__(self, catalog_path: str = None, *args, **kwargs):
+        super().__init__(
+            min_mag=0,
+            max_mag=99,
+            filter_name="None",
+            search_radius_arcmin=0,
+            *args,
+            **kwargs,
+        )
+        self.catalog_path = catalog_path
+
+    def get_catalog(self, ra_deg: float, dec_deg: float) -> astropy.table.Table:
+        catalog = get_table_from_ldac(self.catalog_path)
+        return catalog

mirar/pipelines/winter/blocks.py

@@ -26,6 +26,7 @@
 )
 from mirar.processors.astromatic.swarp.swarp import Swarp
 from mirar.processors.astrometry.anet.anet_processor import AstrometryNet
+from mirar.processors.astrometry.validate import AstrometryStatsWriter
 from mirar.processors.csvlog import CSVLog
 from mirar.processors.dark import DarkCalibrator
 from mirar.processors.mask import MaskPixelsFromPath, WriteMaskedCoordsToFile
@@ -50,7 +51,7 @@
 ]
 
 BOARD_ID = 4
-TARGET_NAME = "EMGW_gal1"
+TARGET_NAME = "m39"
 split = [
     MultiExtParser(
         input_sub_dir="raw/",
@@ -240,9 +241,8 @@
         use_sextractor=True,
         parity="neg",
         search_radius_deg=1.0,
-        # sextractor_config_path=sextractor_autoastrometry_config[
-        #     'config_path'],
-        use_weight=False,
+        sextractor_config_path=sextractor_autoastrometry_config["config_path"],
+        use_weight=True,
     ),
     ImageSaver(output_dir_name=f"anet_{BOARD_ID}", use_existing_weight=False),
     Sextractor(
@@ -362,6 +362,11 @@
         write_regions=True,
         cache=True,
     ),
+    AstrometryStatsWriter(
+        ref_catalog_generator=winter_photometric_catalog_generator,
+        write_regions=True,
+        cache=True,
+    ),
     # ImageSaver(output_dir_name=f"phot_{board_id}_{target_name}")
     ImageSaver(output_dir_name=f"phot_{TARGET_NAME}"),
 ]

mirar/pipelines/winter/config/__init__.py

@@ -1,3 +1,6 @@
+"""
+This file contains the configuration for the winter pipeline.
+"""
 from pathlib import Path
 
 PIPELINE_NAME = "winter"

mirar/pipelines/winter/generator.py

@@ -20,12 +20,14 @@
 from mirar.references.local import RefFromPath
 from mirar.references.ukirt import UKIRTRef
 
+# from mirar.catalog.base_catalog import CatalogFromFile
+
 logger = logging.getLogger(__name__)
 winter_dir = os.path.dirname(__file__)
 astromatic_config_dir = os.path.join(winter_dir, "config/")
 swarp_config_path = os.path.join(astromatic_config_dir, "config.swarp")
 winter_mask_path = os.path.join(winter_dir, "winter_mask.fits")
-scamp_config_path = os.path.join(winter_dir, "scamp.conf")
+scamp_config_path = os.path.join(winter_dir, "config/files/scamp.conf")
 
 
 def winter_reference_generator(image: Image, db_table: Type[BaseDB] = RefStacks):
@@ -101,6 +103,13 @@ def winter_photometric_catalog_generator(image: Image) -> Gaia2Mass:
     search_radius_arcmin = (
         np.max([image["NAXIS1"], image["NAXIS2"]]) * np.abs(image["CD1_1"]) * 60
     ) / 2.0
+
+    # #TODO Remove this make more generic if needed
+    # catalog_path = '/Users/viraj/winter_data/winter/20230629/phot/' \
+    #                'WINTERcamera_20230630-060524-235_mef_4.tmass.cat'
+    # if os.path.exists(catalog_path):
+    #     return CatalogFromFile(catalog_path=catalog_path)
+
     return Gaia2Mass(
         min_mag=10,
         max_mag=20,

mirar/processors/astromatic/scamp/scamp.py

@@ -59,7 +59,7 @@ def run_scamp(
         f"-VERBOSE_TYPE QUIET -SOLVE_PHOTOM N"
     )
 
-    execute(scamp_cmd, output_dir=output_dir, timeout=60.0)
+    execute(scamp_cmd, output_dir=output_dir, timeout=300.0)
 
 
 class Scamp(BaseImageProcessor):

mirar/processors/astromatic/swarp/__init__.py

@@ -1 +1,4 @@
+"""
+Module for interfacing with the SWarp software.
+"""
 from mirar.processors.astromatic.swarp.swarp import run_swarp

mirar/processors/astrometry/validate.py

@@ -0,0 +1,116 @@
+"""
+Module for validating astrometric solutions
+"""
+import logging
+from collections.abc import Callable
+
+import numpy as np
+from astropy.stats import sigma_clipped_stats
+from astropy.table import Table
+
+from mirar.catalog.base_catalog import BaseCatalog
+from mirar.data import Image, ImageBatch
+from mirar.processors.base_catalog_xmatch_processor import BaseProcessorWithCrossMatch
+
+logger = logging.getLogger(__name__)
+
+# All the Sextractor parameters required for this script to run
+REQUIRED_PARAMETERS = [
+    "X_IMAGE",
+    "Y_IMAGE",
+    "FWHM_WORLD",
+    "FWHM_IMAGE",
+    "FLAGS",
+    "ALPHAWIN_J2000",
+    "DELTAWIN_J2000",
+]
+
+
+def get_fwhm(cleaned_img_cat: Table):
+    """
+    Calculate median FWHM from a ldac path
+    Args:
+
+
+    Returns:
+    """
+    mean_fwhm, med_fwhm, std_fwhm = sigma_clipped_stats(cleaned_img_cat["FWHM_WORLD"])
+
+    mean_fwhm_pix, med_fwhm_pix, std_fwhm_pix = sigma_clipped_stats(
+        cleaned_img_cat["FWHM_IMAGE"]
+    )
+    return med_fwhm, mean_fwhm, std_fwhm, med_fwhm_pix, mean_fwhm_pix, std_fwhm_pix
+
+
+def default_sextractor_catalog_purifier(catalog: Table, image: Image) -> Table:
+    """
+    Default function to purify the photometric image catalog
+    """
+    edge_width_pixels = 100
+    fwhm_threshold_arcsec = 4.0
+    x_lower_limit = edge_width_pixels
+    x_upper_limit = image.get_data().shape[1] - edge_width_pixels
+    y_lower_limit = edge_width_pixels
+    y_upper_limit = image.get_data().shape[0] - edge_width_pixels
+
+    clean_mask = (
+        (catalog["FLAGS"] == 0)
+        & (catalog["FWHM_WORLD"] < fwhm_threshold_arcsec / 3600.0)
+        & (catalog["X_IMAGE"] > x_lower_limit)
+        & (catalog["X_IMAGE"] < x_upper_limit)
+        & (catalog["Y_IMAGE"] > y_lower_limit)
+        & (catalog["Y_IMAGE"] < y_upper_limit)
+    )
+
+    return catalog[clean_mask]
+
+
+class AstrometryStatsWriter(BaseProcessorWithCrossMatch):
+    """
+    Processor to calculate astrometry statistics
+    """
+
+    def __init__(
+        self,
+        ref_catalog_generator: Callable[[Image], BaseCatalog],
+        temp_output_sub_dir: str = "phot",
+        image_photometric_catalog_purifier: Callable[
+            [Table, Image], Table
+        ] = default_sextractor_catalog_purifier,
+        crossmatch_radius_arcsec: float = 3.0,
+        write_regions: bool = False,
+        cache: bool = False,
+    ):
+        super().__init__(
+            ref_catalog_generator=ref_catalog_generator,
+            temp_output_sub_dir=temp_output_sub_dir,
+            crossmatch_radius_arcsec=crossmatch_radius_arcsec,
+            sextractor_catalog_purifier=image_photometric_catalog_purifier,
+            write_regions=write_regions,
+            cache=cache,
+            required_parameters=REQUIRED_PARAMETERS,
+        )
+
+    def _apply_to_images(
+        self,
+        batch: ImageBatch,
+    ) -> ImageBatch:
+        for image in batch:
+            ref_cat, _, cleaned_img_cat = self.setup_catalogs(image)
+            fwhm_med, _, fwhm_std, med_fwhm_pix, _, _ = get_fwhm(cleaned_img_cat)
+            image["FWHM_MED"] = fwhm_med
+            image["FWHM_STD"] = fwhm_std
+            image["FWHM_PIX"] = med_fwhm_pix
+
+            _, _, d2d = self.xmatch_catalogs(
+                ref_cat=ref_cat,
+                image_cat=cleaned_img_cat,
+                crossmatch_radius_arcsec=self.crossmatch_radius_arcsec,
+            )
+
+            image.header["ASTUNC"] = np.nanmedian(d2d.value)
+            image.header["ASTFIELD"] = np.arctan(
+                image.header["CD1_2"] / image.header["CD1_1"]
+            ) * (180 / np.pi)
+
+        return batch