<commit_before>"""
legacypipe.galex
================

Code to generate GALEX custom coadds / mosaics.

"""
import os
import numpy as np

import fitsio

from astrometry.util.util import Tan
from astrometry.util.fits import fits_table
from tractor.psf import PixelizedPSF

import logging
logger = logging.getLogger('legacypipe.galex')
def info(*args):
    from legacypipe.utils import log_info
    log_info(logger, args)
def debug(*args):
    from legacypipe.utils import log_debug
    log_debug(logger, args)

def galex_psf(band, galex_dir):
    band2file = {
        'f': os.path.join(galex_dir, 'PSFfuv.fits'),
        'n': os.path.join(galex_dir, 'PSFnuv_faint.fits')
        }
    debug('Reading {}'.format(band2file[band]))
    psfimg = fitsio.read(band2file[band])
    psfimg = psfimg[:psfimg.shape[0] -1, :psfimg.shape[0] -1] # make odd
    psfimg /= psfimg.sum()
    return psfimg

def stage_galex_forced(
    survey=None,
    cat=None,
    T=None,
    targetwcs=None,
    targetrd=None,
    W=None, H=None,
    pixscale=None,
    brickname=None,
    galex_dir=None,
    brick=None,
    galex_ceres=True,
    version_header=None,
    maskbits=None,
    mp=None,
    record_event=None,
    ps=None,
    plots=False,
    **kwargs):
    '''
    After the model fits are finished, we can perform forced
    photometry of the GALEX coadds.
    '''
    from legacypipe.runbrick import _add_stage_version
    from legacypipe.bits import MASKBITS
    from legacypipe.survey import wise_apertures_arcsec
    from tractor import NanoMaggies

    record_event and record_event('stage_galex_forced: starting')
    _add_stage_version(version_header, 'GALEX', 'galex_forced')

    galex_apertures_arcsec = wise_apertures_arcsec

    if not plots:
        ps = None

    # Skip if $GALEX_DIR or --galex-dir not set.
    if galex_dir is None:
        info('GALEX_DIR not set -- skipping GALEX')
        return None

    tiles = galex_tiles_touching_wcs(targetwcs, galex_dir)
    info('Cut to', len(tiles), 'GALEX tiles')

    # the way the roiradec box is used, the min/max order doesn't matter
    roiradec = [targetrd[0,0], targetrd[2,0], targetrd[0,1], targetrd[2,1]]

    # Sources to photometer
    do_phot = np.ones(len(cat), bool)

    # Drop sources within the CLUSTER mask from forced photometry.
    Icluster = None
    if maskbits is not None:
        incluster = (maskbits & MASKBITS['CLUSTER'] > 0)
        if np.any(incluster):
            print('Checking for sources inside CLUSTER mask')
            ra  = np.array([src.getPosition().ra  for src in cat])
            dec = np.array([src.getPosition().dec for src in cat])
            ok,xx,yy = targetwcs.radec2pixelxy(ra, dec)
            xx = np.round(xx - 1).astype(int)
            yy = np.round(yy - 1).astype(int)
            I = np.flatnonzero(ok * (xx >= 0)*(xx < W) * (yy >= 0)*(yy < H))
            if len(I):
                Icluster = I[incluster[yy[I], xx[I]]]
                print('Found', len(Icluster), 'of', len(cat), 'sources inside CLUSTER mask')
                do_phot[Icluster] = False
    Nskipped = len(do_phot) - np.sum(do_phot)

    gcat = []
    for i in np.flatnonzero(do_phot):
        src = cat[i]
        src = src.copy()
        src.setBrightness(NanoMaggies(galex=1.))
        gcat.append(src)

    # use pixelized PSF model
    pixpsf = True

    # Photometer the two bands in parallel
    args = []
    for band in ['n','f']:
        btiles = tiles[tiles.get('has_%s' % band)]
        if len(btiles) == 0:
            continue
        args.append((galex_dir, gcat, btiles, band, roiradec, pixpsf, ps, galex_ceres))
    # Run the forced photometry!
    record_event and record_event('stage_galex_forced: photometry')
    phots = mp.map(galex_phot, args)
    record_event and record_event('stage_galex_forced: results')

    # Initialize coadds even if we don't have any overlapping tiles -- we'll write zero images.
    # Create the WCS into which we'll resample the tiles.
    # Same center as "targetwcs" but bigger pixel scale.
    gpixscale = 1.5
    gra  = np.array([src.getPosition().ra  for src in cat])
    gdec = np.array([src.getPosition().dec for src in cat])
    rc,dc = targetwcs.radec_center()
    ww = int(W * pixscale / gpixscale)
    hh = int(H * pixscale / gpixscale)
    gcoadds = GalexCoadd(rc, dc, ww, hh, gpixscale)

    # Unpack results...
    GALEX = None
    if len(phots):
        # The "phot" results for the full-depth coadds are one table per
        # band.  Merge all those columns.
        galex_models = []
        for i,p in enumerate(phots[:len(args)]):
            if p is None:
                (_,_,tiles,band) = args[i][:4]
                print('"None" result from GALEX forced phot:', tiles, band)
                continue
            galex_models.extend(p.models)
            if GALEX is None:
                GALEX = p.phot
            else:
                GALEX.add_columns_from(p.phot)
        gcoadds.add(galex_models)

        if Nskipped > 0:
            GALEX = _fill_skipped_values(GALEX, Nskipped, do_phot)
            assert(len(GALEX) == len(cat))
            assert(len(GALEX) == len(T))
<commit_msg>galex crash; import _fill_skipped_values<commit_after>"""
legacypipe.galex
================

Code to generate GALEX custom coadds / mosaics.

"""
import os
import numpy as np

import fitsio

from astrometry.util.util import Tan
from astrometry.util.fits import fits_table
from tractor.psf import PixelizedPSF

import logging
logger = logging.getLogger('legacypipe.galex')
def info(*args):
    from legacypipe.utils import log_info
    log_info(logger, args)
def debug(*args):
    from legacypipe.utils import log_debug
    log_debug(logger, args)

def galex_psf(band, galex_dir):
    band2file = {
        'f': os.path.join(galex_dir, 'PSFfuv.fits'),
        'n': os.path.join(galex_dir, 'PSFnuv_faint.fits')
        }
    debug('Reading {}'.format(band2file[band]))
    psfimg = fitsio.read(band2file[band])
    psfimg = psfimg[:psfimg.shape[0] -1, :psfimg.shape[0] -1] # make odd
    psfimg /= psfimg.sum()
    return psfimg

def stage_galex_forced(
    survey=None,
    cat=None,
    T=None,
    targetwcs=None,
    targetrd=None,
    W=None, H=None,
    pixscale=None,
    brickname=None,
    galex_dir=None,
    brick=None,
    galex_ceres=True,
    version_header=None,
    maskbits=None,
    mp=None,
    record_event=None,
    ps=None,
    plots=False,
    **kwargs):
    '''
    After the model fits are finished, we can perform forced
    photometry of the GALEX coadds.
    '''
    from legacypipe.runbrick import _add_stage_version
    from legacypipe.bits import MASKBITS
    from legacypipe.survey import wise_apertures_arcsec
    from tractor import NanoMaggies

    record_event and record_event('stage_galex_forced: starting')
    _add_stage_version(version_header, 'GALEX', 'galex_forced')

    galex_apertures_arcsec = wise_apertures_arcsec

    if not plots:
        ps = None

    # Skip if $GALEX_DIR or --galex-dir not set.
    if galex_dir is None:
        info('GALEX_DIR not set -- skipping GALEX')
        return None

    tiles = galex_tiles_touching_wcs(targetwcs, galex_dir)
    info('Cut to', len(tiles), 'GALEX tiles')

    # the way the roiradec box is used, the min/max order doesn't matter
    roiradec = [targetrd[0,0], targetrd[2,0], targetrd[0,1], targetrd[2,1]]

    # Sources to photometer
    do_phot = np.ones(len(cat), bool)

    # Drop sources within the CLUSTER mask from forced photometry.
    Icluster = None
    if maskbits is not None:
        incluster = (maskbits & MASKBITS['CLUSTER'] > 0)
        if np.any(incluster):
            print('Checking for sources inside CLUSTER mask')
            ra  = np.array([src.getPosition().ra  for src in cat])
            dec = np.array([src.getPosition().dec for src in cat])
            ok,xx,yy = targetwcs.radec2pixelxy(ra, dec)
            xx = np.round(xx - 1).astype(int)
            yy = np.round(yy - 1).astype(int)
            I = np.flatnonzero(ok * (xx >= 0)*(xx < W) * (yy >= 0)*(yy < H))
            if len(I):
                Icluster = I[incluster[yy[I], xx[I]]]
                print('Found', len(Icluster), 'of', len(cat), 'sources inside CLUSTER mask')
                do_phot[Icluster] = False
    Nskipped = len(do_phot) - np.sum(do_phot)

    gcat = []
    for i in np.flatnonzero(do_phot):
        src = cat[i]
        src = src.copy()
        src.setBrightness(NanoMaggies(galex=1.))
        gcat.append(src)

    # use pixelized PSF model
    pixpsf = True

    # Photometer the two bands in parallel
    args = []
    for band in ['n','f']:
        btiles = tiles[tiles.get('has_%s' % band)]
        if len(btiles) == 0:
            continue
        args.append((galex_dir, gcat, btiles, band, roiradec, pixpsf, ps, galex_ceres))
    # Run the forced photometry!
    record_event and record_event('stage_galex_forced: photometry')
    phots = mp.map(galex_phot, args)
    record_event and record_event('stage_galex_forced: results')

    # Initialize coadds even if we don't have any overlapping tiles -- we'll write zero images.
    # Create the WCS into which we'll resample the tiles.
    # Same center as "targetwcs" but bigger pixel scale.
    gpixscale = 1.5
    gra  = np.array([src.getPosition().ra  for src in cat])
    gdec = np.array([src.getPosition().dec for src in cat])
    rc,dc = targetwcs.radec_center()
    ww = int(W * pixscale / gpixscale)
    hh = int(H * pixscale / gpixscale)
    gcoadds = GalexCoadd(rc, dc, ww, hh, gpixscale)

    # Unpack results...
    GALEX = None
    if len(phots):
        # The "phot" results for the full-depth coadds are one table per
        # band.  Merge all those columns.
        galex_models = []
        for i,p in enumerate(phots[:len(args)]):
            if p is None:
                (_,_,tiles,band) = args[i][:4]
                print('"None" result from GALEX forced phot:', tiles, band)
                continue
            galex_models.extend(p.models)
            if GALEX is None:
                GALEX = p.phot
            else:
                GALEX.add_columns_from(p.phot)
        gcoadds.add(galex_models)

        if Nskipped > 0:
            from legacypipe.runbrick import _fill_skipped_values
            GALEX = _fill_skipped_values(GALEX, Nskipped, do_phot)
            assert(len(GALEX) == len(cat))
            assert(len(GALEX) == len(T))

{'fraction_non_alphanumeric': 0.11536979969183359, 'fraction_numerical': 0.005200308166409861, 'mean_word_length': 2.5773337926283157, 'pattern_counts': {'":': 0, '<': 7, '<?xml version=': 0, '>': 11, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'dir': 'git-commits-cleaned', 'id': '4474972', 'n_tokens_mistral': 3717, 'n_tokens_neox': 3530, 'n_words': 899}