pipeline_workerd

#!/usr/bin/env python3
#
# This file is part of the Robotic Observatory Control Kit (rockit)
#
# rockit is free software: you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# rockit is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with rockit.  If not, see <http://www.gnu.org/licenses/>.

"""Frame pipeline daemon"""

# pylint: disable=too-many-arguments
# pylint: disable=too-many-branches
# pylint: disable=too-many-locals
# pylint: disable=too-many-statements
# pylint: disable=too-many-lines
# pylint: disable=broad-exception-raised


import argparse
import array
import base64
import datetime
import io
import json
import multiprocessing
import os
import platform
import queue
import re
import shutil
import subprocess
import sys
import tempfile
import threading
import traceback
from astropy.coordinates import Angle
from astropy.io import fits
from astropy.table import Table
from astropy.time import Time
import astropy.units as u
from astropy.wcs import WCS
import numpy as np
from paramiko import SSHClient, RSAKey
from PIL import Image, ImageOps
import Pyro4
from scp import SCPClient
import sep
from rockit.common import daemons, log, IP
from rockit.common.helpers import pyro_client_matches
from rockit.pipeline import CommandStatus, Config

if platform.system() != 'Windows':
    import grp
    import pwd


def window_header_region(header, data, keyword):
    """Returns a cropped view of frame data"""
    if keyword not in header:
        return data, 0, 0

    region = header[keyword]
    r = re.search(r'^\[(\d+):(\d+),(\d+):(\d+)\]$', region).groups()
    x1 = int(r[0]) - 1
    x2 = int(r[1])
    y1 = int(r[2]) - 1
    y2 = int(r[3])

    if x1 == 0 and y1 == 0 and x2 == data.shape[1] and y2 == data.shape[0]:
        return data, 0, 0

    return data[y1:y2, x1:x2], x1, y1


def rescale_image_data(data, clip_low, clip_high):
    """ Returns a normalised array where clip_low percent of the pixels are 0 and
        clip_high percent of the pixels are 255
    """
    high = np.percentile(data, clip_high)
    low = np.percentile(data, clip_low)

    # Avoid a divide-by-zero
    if high == low:
        high = low + 1

    scale = 255. / (high - low)
    data = np.clip(data, low, high)
    return scale * (data - low)


def create_formatted_card(key, value, comment, format_type):
    float_places = {
        'float1dp': 1,
        'float2dp': 2,
        'float3dp': 3,
        'float5dp': 5
    }

    try:
        if format_type in float_places:
            return create_float_card(key, value, comment, float_places[format_type])

        if format_type == 'int':
            value = int(value)
        elif format_type == 'bool':
            value = bool(value)
        elif format_type == 'sexagesimalha':
            value = Angle(value * u.deg).to(u.hourangle).to_string(sep=':', precision=2)
        elif format_type == 'sexagesimaldeg':
            value = Angle(value * u.deg).to_string(sep=':', precision=2)

        return fits.Card(key, value, comment)
    except Exception:
        print(f'failed to create card {key} with error: ')
        traceback.print_exc(file=sys.stdout)
        sys.stdout.flush()
        return fits.Card('COMMENT', f' {key} not available', '')


def create_float_card(key, value, comment, places):
    """Create a fits.Card instance with a float value
       rounded to a specified number of significant figures
    """
    card = fits.Card.fromstring(f'{key:8s}= {value:.0{places}f}')
    card.comment = comment
    return card


def header_to_dict(header):
    ret = {}
    for card in header.cards:
        if card.is_blank or not card.keyword or card.keyword in ['COMMENT', 'HISTORY']:
            continue

        ret[card.keyword] = card.value
    return ret


def move_bscale_bzero(header):
    """Moves the BSCALE and BZERO header keys to the start of the header
       to be grouped with the other image keywords"""
    header.set('BSCALE', header['BSCALE'], after='NAXIS2')
    header.set('BZERO', header['BZERO'], after='BSCALE')


def add_telescope_header(header, info, log_name, telescope_cards, telescope_query_timeout):
    """Queries the telescope status and adds the appropriate header keys to the frame"""
    desc = fits.Card('COMMENT', ' ---          TELESCOPE INFORMATION          --- ')
    header.append(fits.Card(None, None), end=True)
    header.append(desc, end=True)

    data = {}
    # Query data from daemons
    for card in telescope_cards:
        key = card['daemon'].uri + '_' + card['method']
        if key in data:
            continue

        if 'camera' in card and card['camera'] != info['camera_id']:
            continue

        try:
            with card['daemon'].connect(telescope_query_timeout) as d:
                data[key] = getattr(d, card['method'])() or {}
        except Exception as e:
            log.error(log_name, f'Failed to query telescope metadata ({e})')
            sys.stdout.flush()
            data[card['daemon']] = {}

    # Build cards
    for card in telescope_cards:
        if 'camera' in card and card['camera'] != info['camera_id']:
            continue

        key = card['daemon'].uri + '_' + card['method']
        if isinstance(card['parameter'], list):
            value = data.get(key, None)
            for parameter in card['parameter']:
                if value is None:
                    break
                value = value.get(parameter, None)
        else:
            value = data.get(key, {}).get(card['parameter'], None)

        if value is not None:
            header.append(create_formatted_card(card['key'], value, card['comment'], card['type']), end=True)
        else:
            header.append(fits.Card('COMMENT', f' {card["key"]} not available', ''), end=True)


def add_environment_header(header, info, log_name, environment_cards, environment_daemon, environment_query_timeout):
    """Queries the environment status and adds the appropriate header keys to the frame"""
    desc = fits.Card('COMMENT', ' ---         ENVIRONMENT INFORMATION         --- ')
    header.append(fits.Card(None, None), end=True)
    header.append(desc, end=True)

    try:
        with environment_daemon.connect(environment_query_timeout) as environment:
            data = environment.status()
    except Exception as e:
        data = {}
        log.error(log_name, f'Failed to query environment metadata ({e})')
        sys.stdout.flush()

    for card in environment_cards:
        if 'camera' in card and card['camera'] != info['camera_id']:
            continue

        sensor_data = data.get(card['sensor'], {}).get('parameters', {}).get(card['parameter'], {})
        if sensor_data and sensor_data.get('current', False):
            header.append(create_formatted_card(card['key'], sensor_data['latest'],
                                                card['comment'], card['type']), end=True)
        else:
            header.append(fits.Card('COMMENT', f' {card["key"]} not available', ''), end=True)


def add_pipeline_header(header, info, filename):
    """Adds pipeline configuration to the frame header"""
    desc = fits.Card('COMMENT', ' ---              DATA PIPELINE              --- ')
    header.append(fits.Card(None, None), end=True)
    header.append(desc, end=True)

    header.append(fits.Card('FILESAVE', filename is not None,
                            'image has been archived to disk'), end=True)

    if filename is not None:
        header.append(fits.Card('FILENAME', filename, 'archived image name'), end=True)
    else:
        header.append(fits.Card('COMMENT', ' FILENAME not available', ''))

    header.append(fits.Card('IMAGETYP', info['frame_type'], 'frame type'), end=True)

    if info['frame_type'] == 'SCIENCE':
        header.append(fits.Card('OBJECT', info['frame_object'], 'science target name'), end=True)
    else:
        header.append(fits.Card('COMMENT', ' OBJECT not available', ''))


def add_operations_header(header, cards):
    """Adds custom action information to the frame header"""
    desc = fits.Card('COMMENT', ' ---         OPERATIONS INFORMATION          --- ')
    header.append(fits.Card(None, None), end=True)
    header.append(desc, end=True)

    for card in cards:
        keyword = card.get('keyword', None)
        value = card.get('value', None)
        comment = card.get('comment', None)
        if keyword and value is not None:
            header.append(fits.Card(keyword[0:8], value, comment), end=True)


def add_wcs_header(header, objects, binning, add_rotation, wcs_config):
    """Solves frame WCS and adds the appropriate header keys to the frame.
       Returns a list of xpa messages that can be given to update_previews"""

    width = header['NAXIS1']
    height = header['NAXIS2']
    wcs_failed_card = fits.Card('WCSSOLVE', False, 'attempted wcs solution was successful?')

    if objects is None or len(objects) <= 0:
        header.append(wcs_failed_card, end=True)
        return

    if wcs_config is None:
        header.append(wcs_failed_card, end=True)
        return

    try:
        wcs_objects = objects[:50]['X', 'Y', 'FLUX']
        ra = wcs_config.get('tel_ra_card', None)
        dec = wcs_config.get('tel_dec_card', None)
        radius = wcs_config.get('search_radius', None)

        if 'astrometry_daemon' in wcs_config:
            astrometry_daemon = getattr(daemons, wcs_config['astrometry_daemon'])
            table_bytes = io.BytesIO()
            wcs_objects.write(table_bytes, format='fits')
            table_bytes.seek(0)
            wcs_args = {
                'width': width,
                'height': height,
                'scale_high': binning * wcs_config['scale_high'],
                'scale_low': binning * wcs_config['scale_low'],
                'objects': base64.b64encode(table_bytes.getvalue()).decode('ascii'),
                'timeout': wcs_config['timeout'],
                'camera_id': header['CAMID'],
                'exposure_count': header['EXPCNT']
            }

            if ra and ra in header and dec and dec in header and radius:
                wcs_args.update([
                    'ra', header[ra],
                    'dec', header[dec],
                    'radius', radius
                ])

            timeout = astrometry_daemon.default_timeout + wcs_config['timeout']
            with astrometry_daemon.connect(timeout=timeout) as astrometry:
                solution = astrometry.solve(wcs_args)
                if not solution:
                    raise Exception('astrometry failed')

        else:
            with tempfile.TemporaryDirectory() as tmp_dir:
                wcs_path = os.path.join(tmp_dir, 'scratch.wcs')
                table_path = os.path.join(tmp_dir, 'scratch.xyls')
                wcs_objects.write(table_path, format='fits')

                wcs_args = [
                    '/usr/bin/solve-field',
                    '--overwrite', '--no-plots',
                    '--axy', 'none', '--rdls', 'none', '--match', 'none',
                    '--corr', 'none', '--solved', 'none',
                    '--scale-units', 'arcsecperpix', '--scale-high', str(binning * wcs_config['scale_high']),
                    '--scale-low', str(binning * wcs_config['scale_low']), '--width', str(width), '--height', str(height)]

                if ra and ra in header and dec and dec in header and radius:
                    wcs_args.extend([
                        '--ra', header[ra],
                        '--dec', header[dec],
                        '--radius', str(radius)
                    ])

                wcs_args.append(table_path)

                subprocess.check_call(wcs_args, timeout=wcs_config['timeout'],
                                      stdout=subprocess.DEVNULL,
                                      stderr=subprocess.DEVNULL)

                solution = header_to_dict(fits.Header.fromfile(wcs_path))

        header.append(fits.Card('WCSSOLVE', True,
                                'attempted wcs solution was successful?'), end=True)

        wcs_ignore_cards = ['SIMPLE', 'BITPIX', 'NAXIS', 'EXTEND', 'DATE', 'IMAGEW', 'IMAGEH']
        comment = 'astrometry.net wcs solution'
        for key, value in solution.items():
            if key not in wcs_ignore_cards:
                header.append(fits.Card(key, value, comment), end=True)

        if add_rotation:
            wcs = WCS(header)
            coords = wcs.all_pix2world([width // 2, width // 2], [height // 2, height // 2 + 1], 0)
            north_angle = np.degrees(np.arctan2(np.diff(coords[1]), np.diff(coords[0]))).item()
            header.append(create_float_card('FIELDROT', north_angle,
                                            '[deg] direction of north measured ccw from +x', 2), end=True)
    except Exception:
        print('failed to update wcs with error: ')
        traceback.print_exc(file=sys.stdout)
        header.append(wcs_failed_card, end=True)


def detect_objects(data, background_rms, config, binning, crop_x, crop_y):
    """returns an Astropy table containing (X, Y, FLUX, HFD_PX, HFD_ARCSEC) for each
       object detected by sep"""
    try:
        platescale = binning * config['platescale']

        # Note: data is expected to be background subtracted
        if isinstance(data, np.ma.MaskedArray):
            mask = data.mask
            data = data.data
        else:
            mask = None

        objects = sep.extract(data, 5 * background_rms, mask=mask)

        # Work around sep issue #110
        objects['theta'][objects['theta'] > np.pi / 2] -= np.pi / 2

        kronrad, kron_flag = sep.kron_radius(data, objects['x'], objects['y'],
                                             objects['a'], objects['b'], objects['theta'], 6.0, mask=mask)

        flux, _, flux_flag = sep.sum_ellipse(data, objects['x'], objects['y'], objects['a'], objects['b'],
                                             objects['theta'], 2.5 * kronrad, subpix=0, mask=mask)

        r, radius_flag = sep.flux_radius(data, objects['x'], objects['y'],
                                         6.0 * objects['a'], 0.5, normflux=flux, subpix=5, mask=mask)

        # pylint: disable=no-member
        filt = np.logical_and.reduce([
            objects['npix'] >= config['object_minpix'],
            kron_flag == 0,
            flux_flag == 0,
            radius_flag == 0
        ])
        # pylint: enable=no-member

        objects = Table(
            [objects['x'][filt] + crop_x, objects['y'][filt] + crop_y, flux[filt],
             2 * r[filt], 2 * r[filt] * platescale],
            names=('X', 'Y', 'FLUX', 'HFD_PX', 'HFD_ARCSEC'),
            dtype=(float, float, float, float, float)
        )

        objects.sort('FLUX', reverse=True)
        return objects
    except Exception:
        print('failed to extract objects with error: ')
        traceback.print_exc(file=sys.stdout)
        sys.stdout.flush()
        return Table(names=('X', 'Y', 'FLUX', 'HFD_PX', 'HFD_ARCSEC'), dtype=(float, float, float, float, float))


def add_hfd_header(header, objects):
    """Adds the median hfd to the frame header"""
    try:
        if objects is None or len(objects) <= 0:
            return

        header.append(create_float_card('MEDHFD', np.median(objects['HFD_ARCSEC']).item(),
                                        '[arcsec] median estimated half-flux diameter', 1), end=True)
        header.append(fits.Card('HFDCNT', len(objects),
                                'number of sources used to estimate the HFD'), end=True)
    except Exception:
        print('failed to calculate HFD with error:')
        traceback.print_exc(file=sys.stdout)
        sys.stdout.flush()


def add_intensity_stats_header(header, data_raw, data_cropped, config, master_bias_stats):
    """Adds the frame intensity to the frame header"""
    try:
        header.append(create_float_card('MEANCNTS', np.ma.mean(data_cropped),
                                        '[adu] mean frame counts', 2), end=True)
        header.append(create_float_card('MEDCNTS', np.ma.median(data_cropped),
                                        '[adu] median frame counts', 2), end=True)

        if 'overscan_region_card' in config and config['overscan_region_card'] in header:
            bias, *_ = window_header_region(header, data_raw, config['overscan_region_card'])
            header.append(create_float_card('MEANBIAS', np.ma.mean(bias),
                                            '[adu] mean overscan bias counts', 2), end=True)
            header.append(create_float_card('MEDBIAS', np.ma.median(bias),
                                            '[adu] median overscan bias counts', 2), end=True)
        elif master_bias_stats:
            coadds = 1
            if 'coadditions_card' in config and config['coadditions_card'] in header:
                coadds = header[config['coadditions_card']]

            header.append(create_float_card('MEANBIAS', coadds * master_bias_stats['mean'],
                                            '[adu] nominal mean bias counts', 2), end=True)
            header.append(create_float_card('MEDBIAS', coadds * master_bias_stats['median'],
                                            '[adu] nominal median bias counts', 2), end=True)

    except Exception:
        print('failed to calculate intensity stats with error: ')
        traceback.print_exc(file=sys.stdout)
        sys.stdout.flush()


def add_hfd_grid_header(header, objects, config):
    """Adds the hfd grid samples to the frame header"""
    header.append(fits.Card('HFD_ROWS', config['hfd_grid_tiles_y'], 'number of rows in hfd grid'))
    header.append(fits.Card('HFD_COLS', config['hfd_grid_tiles_x'], 'number of columns in hfd grid'))
    x_step = header['NAXIS1'] / config['hfd_grid_tiles_x']
    y_step = header['NAXIS2'] / config['hfd_grid_tiles_y']
    for j in range(config['hfd_grid_tiles_y']):
        for i in range(config['hfd_grid_tiles_x']):
            # pylint: disable=no-member
            filt = np.logical_and.reduce([
                x_step * i < objects['X'],
                objects['X'] < x_step * (i + 1),
                y_step * j < objects['Y'],
                objects['Y'] < y_step * (j + 1)
            ])
            # pylint: enable=no-member

            median = np.median(objects['HFD_ARCSEC'][filt])
            if np.isnan(median):
                header.append(fits.Card('COMMENT', f' HFD_{j:02d}{i:02d} not available', ''), end=True)
            else:
                header.append(create_float_card(f'HFD_{j:02d}{i:02d}', median,
                                                '[arcsec] hfd measurement', 2), end=True)


def update_dashboard(header, data, ssh_root_path, config, process_steps, process_time):
    try:
        size = np.shape(data)
        metadata = {
            'date': header['DATE-OBS'],
            'src': header.get('TIME-SRC', ''),
            'exptime': header['EXPTIME'],
            'filter': header.get('FILTER', ''),
            'type': header.get('IMAGETYP', ''),
            'object': header.get('OBJECT', ''),
            'saved': header['FILESAVE'],
            'filename': header.get('FILENAME', ''),
            'size': [size[1], size[0]],
            'clipsize': config['clip_size'],
            'process_steps': process_steps,
            'process_time': process_time
        }

        scaled = rescale_image_data(data, config['min_threshold'], config['max_threshold'])
        preview = Image.fromarray(scaled).convert('RGB')
        if config['flip_vertical']:
            preview = ImageOps.flip(preview)
        if config['flip_horizontal']:
            preview = ImageOps.mirror(preview)

        preview.thumbnail((config['thumb_size'], config['thumb_size']))
        thumb_bytes = io.BytesIO()
        preview.save(thumb_bytes, 'JPEG', quality=80, optimize=True, progressive=True, clobber=True)
        thumb_bytes.seek(0)

        # Clip the central region from the image
        cx1 = (data.shape[0] - config['clip_size']) // 2
        cx2 = cx1 + config['clip_size']
        cy1 = (data.shape[1] - config['clip_size']) // 2
        cy2 = cy1 + config['clip_size']

        clipped = rescale_image_data(data[cx1:cx2, cy1:cy2],
                                     config['min_threshold'],
                                     config['max_threshold'])
        preview = Image.fromarray(clipped).convert('RGB')
        preview = ImageOps.flip(preview)
        clip_bytes = io.BytesIO()
        preview.save(clip_bytes, 'JPEG', quality=40, optimize=True, progressive=True, clobber=True)
        clip_bytes.seek(0)

        json_bytes = io.BytesIO()
        json_bytes.write(json.dumps(metadata).encode())
        json_bytes.seek(0)

        if 'local_path' in config and platform.system() != 'Windows':
            uid = pwd.getpwnam(config['user']).pw_uid
            gid = grp.getgrnam(config['group']).gr_gid

            path = os.path.join(config['local_path'], config['prefix'] + '-thumb.jpg')
            with open(path, 'wb') as f:
                f.write(thumb_bytes.getbuffer().tobytes())
            os.chown(path, uid, gid)

            path = os.path.join(config['local_path'], config['prefix'] + '-clip.jpg')
            with open(path, 'wb') as f:
                f.write(clip_bytes.getbuffer().tobytes())
            os.chown(path, uid, gid)

            path = os.path.join(config['local_path'], config['prefix'] + '.json')
            with open(path, 'wb') as f:
                f.write(json_bytes.getbuffer().tobytes())
            os.chown(path, uid, gid)

        if 'remote_path' in config:
            with open(os.path.join(ssh_root_path, config['key']), encoding='utf-8') as keyfile:
                key = RSAKey.from_private_key(keyfile)

            with SSHClient() as ssh:
                ssh.load_system_host_keys(os.path.join(ssh_root_path, 'known_hosts'))
                ssh.connect(getattr(IP, config['remote_machine']), username=config['remote_user'], pkey=key)

                with SCPClient(ssh.get_transport()) as scp:
                    scp.putfo(thumb_bytes, remote_path=f'{config["remote_path"]}/{config["prefix"]}-thumb.jpg')
                    scp.putfo(clip_bytes, remote_path=f'{config["remote_path"]}/{config["prefix"]}-clip.jpg')
                    scp.putfo(json_bytes, remote_path=f'{config["remote_path"]}/{config["prefix"]}.json')
    except Exception:
        print('failed to generate dashboard preview with error:')
        traceback.print_exc(file=sys.stdout)
        sys.stdout.flush()


def notify_opsd(header, ops_daemon, log_name):
    try:
        with ops_daemon.connect() as ops:
            return ops.notify_processed_frame(header_to_dict(header))

    except Exception as e:
        log.error(log_name, f'Failed to notify operations daemon of completed frame ({e})')
        sys.stdout.flush()
        return None


def notify_guide_profiles(header, data, log_name, ops_daemon):
    # Collapse into x/y profiles
    # Note: data is expected to already be background-subtracted
    # Note: Pyro doesn't support numpy arrays, so convert to a built-in type
    profile_x = array.array('f', data.mean(axis=0))
    profile_y = array.array('f', data.mean(axis=1))

    try:
        with ops_daemon.connect() as ops:
            return ops.notify_guide_profiles(header_to_dict(header), profile_x, profile_y)

    except Exception as e:
        log.error(log_name, f'Failed to notify operations daemon of guide profiles ({e})')
        sys.stdout.flush()
        return None


def process_previews(preview_queue, ssh_root_path, preview_registrations, preview_lock):
    """A thread running in each processing process that allows ds9 preview updates
       to run concurrently with the main processing work.
     """

    while True:
        # Block until a preview is available
        data, header = preview_queue.get()
        data_bytes = io.BytesIO()

        with preview_lock:
            for r in preview_registrations:
                try:
                    if data is not None:
                        hdu = fits.PrimaryHDU(data, header=header)
                        move_bscale_bzero(hdu.header)
                        fits.HDUList(hdu).writeto(data_bytes)

                    with SSHClient() as ssh:
                        ssh.load_system_host_keys(os.path.join(ssh_root_path, 'known_hosts'))
                        with open(os.path.join(ssh_root_path, r['rockit.key']), encoding='utf-8') as keyfile:
                            key = RSAKey.from_private_key(keyfile)

                        ssh.connect(r['rockit.remote_ip'], username=r['rockit.remote_user'], pkey=key)
                        with SCPClient(ssh.get_transport()) as scp:
                            data_bytes.seek(0)
                            scp.putfo(data_bytes, remote_path=r['rockit.remote_path'])

                        _, stdout, _ = ssh.exec_command(f'{r["rockit.remote_command"]} "{r["samp.hub.xmlrpc.url"]}" "{r["samp.secret"]}" "{r["rockit.client"]}" "{r["rockit.remote_path"]}" "{r["rockit.window_scale"]}"')
                        ret = stdout.channel.recv_exit_status()
                        ssh.exec_command(f'rm "{r["rockit.remote_path"]}"')
                        if ret != 0:
                            raise Exception(f'display command returned code {ret}')

                except Exception as e:
                    preview_registrations.remove(r)
                    print(f'unregistering preview {r["rockit.remote_user"]}@{r["rockit.remote_ip"]}:{r["rockit.remote_path"]}: {e}', flush=True)

def process_frames(processing_queue, camera_config,
                   preview_timestamp, preview_registrations, preview_lock,
                   guide_timestamp, guide_min_interval,
                   dashboard_timestamp, log_name,
                   environment_cards, environment_daemon_name, environment_query_timeout,
                   telescope_cards, telescope_query_timeout,
                   ops_daemon_name):
    """

    Helper process to add frame metadata and calculate image details for the operations daemon.
    This uses a process (rather than a thread) to avoid the GIL bottlenecking throughput,
    and multiple worker processes allow frames to be processed in parallel.
    """
    environment_daemon = ops_daemon = None
    if environment_daemon_name is not None:
        environment_daemon = getattr(daemons, environment_daemon_name)
    if ops_daemon_name is not None:
        ops_daemon = getattr(daemons, ops_daemon_name)

    # Update previews on a separate thread to avoid blocking
    preview_queue = queue.Queue()
    threading.Thread(target=process_previews, daemon=True,
                     args=(preview_queue, camera_config['ssh_root_path'], preview_registrations, preview_lock)
                     ).start()

    # Apply corrections if required
    if 'calibrations_path' in camera_config:
        with fits.open(camera_config['calibrations_path']) as f:
            master_bad_pixel_mask = f['MASK'].data.astype(bool)
            master_bias = np.ma.MaskedArray(f['BIAS'].data, master_bad_pixel_mask)

            master_bias_stats = {
                'median': np.ma.median(master_bias),
                'mean': np.ma.mean(master_bias)
            }

            # Normalize to just the pixel non-uniformity
            master_bias_pixel_nonuniformity = master_bias - master_bias_stats['median']

            master_dark_adu_per_sec = f['DARK'].data

    else:
        master_bias_stats = None
        master_dark_adu_per_sec = None
        master_bias_pixel_nonuniformity = None
        master_bad_pixel_mask = None

    while True:
        # Blocks until a frame is available for processing
        info = processing_queue.get()
        loadpath = os.path.join(camera_config['input_data_path'], info['filename'])

        try:
            start = datetime.datetime.utcnow()

            # We only update the header, not the data.
            # Reduce the overheads from astropy fits by reading what we need at the start,
            # manipulating just the headers in isolation, and then doing a direct byte
            # replacement of the header in the original file.
            # This relies on the camera daemon reserving sufficient blank header keywords
            # to be replaced below. This avoids having to move all the pixel data
            # in the file to create additional header blocks.
            header = fits.getheader(loadpath)
            steps = ['headers']

            # Strip empty rows from the end of the header
            # fits blocks are 2880 bytes, which can hold 36 80 byte header cards.
            # The loaded header does not include the final END card, so in the event
            # that len(header) is a multiple of 36 we know that there must be one
            # extra block
            header_block_capacity = len(header) // 36 + 1
            while header.cards[-1][0] == '' and header.cards[-1][1] == '':
                header.pop()

            # We can skip loading the frame data if we know for sure that
            # we aren't going to need it.
            data_bgsubtracted = data_cropped = data_raw = data_background = None
            crop_x = crop_y = 0
            header_time = Time(header['DATE-OBS'])

            guide_profiles_stale = False
            if info['guide_profiles']:
                with guide_timestamp.get_lock():
                    if header_time > Time(guide_timestamp.value.decode('ascii')) + guide_min_interval * u.s:
                        guide_timestamp.value = header['DATE-OBS'].encode('ascii')
                        guide_profiles_stale = True

            preview_stale = False
            if len(preview_registrations):
                with preview_timestamp.get_lock():
                    next_preview_time = Time(preview_timestamp.value.decode('ascii')) + \
                        camera_config['preview_min_interval'] * u.s
                    if header_time > next_preview_time:
                        preview_timestamp.value = header['DATE-OBS'].encode('ascii')
                        preview_stale = True

            dashboard_stale = False
            if info['dashboard_enabled'] and 'dashboard' in camera_config:
                with dashboard_timestamp.get_lock():
                    next_dashboard_time = Time(dashboard_timestamp.value.decode('ascii')) + \
                        camera_config['dashboard']['min_interval'] * u.s
                    if header_time > next_dashboard_time:
                        dashboard_timestamp.value = header['DATE-OBS'].encode('ascii')
                        dashboard_stale = True

            calculate_objects = info['hfd_enabled'] or info['hfd_grid_enabled'] or \
                info['wcs_enabled'] or info['rotation_enabled']

            calculate_bgsubtracted = calculate_objects or guide_profiles_stale
            calculate_cropped = calculate_bgsubtracted or info['intensity_stats_enabled'] or dashboard_stale
            data_raw_float = None
            if calculate_cropped or preview_stale:
                data_raw = fits.getdata(loadpath)
                steps.append('data')

                window_card = camera_config.get('window_region_card', None)

                if master_bias_pixel_nonuniformity is not None or master_dark_adu_per_sec is not None:
                    data_raw_dtype = data_raw.dtype
                    data_raw_float = data_raw.astype(float)

                    # Correct for bias pattern
                    if master_bias_pixel_nonuniformity is not None:
                        coadds = 1
                        if 'coadditions_card' in camera_config and camera_config['coadditions_card'] in header:
                            coadds = header[camera_config['coadditions_card']]

                        correction, *_ = window_header_region(header, master_bias_pixel_nonuniformity, window_card)
                        data_raw_float -= coadds * correction
                        steps.append('bias')

                    # Correct for hot pixels
                    if master_dark_adu_per_sec is not None:
                        correction, *_ = window_header_region(header, master_dark_adu_per_sec, window_card)
                        data_raw_float -= header['EXPTIME'] * correction
                        steps.append('dark')

                    data_raw = data_raw_float.astype(data_raw_dtype)

                    if master_bad_pixel_mask is not None:
                        mask, *_ = window_header_region(header, master_bad_pixel_mask, window_card)
                        data_raw = np.ma.array(data_raw, mask=mask)
                        data_raw_float = np.ma.array(data_raw_float, mask=mask)

                    # Avoid negative values wrapping to large positive values
                    limits = np.iinfo(data_raw_dtype)
                    data_raw[data_raw_float <= limits.min] = limits.min
                    data_raw[data_raw_float >= limits.max] = limits.max

            if calculate_cropped:
                region = camera_config.get('image_region_card', None)
                data_cropped, crop_x, crop_y = window_header_region(header, data_raw, region)

                if data_cropped is not data_raw:
                    steps.append('crop')

            if calculate_bgsubtracted:
                if data_raw_float is not None:
                    data_bgsubtracted, *_ = window_header_region(header, data_raw_float, region)
                else:
                    data_bgsubtracted = data_cropped.astype(float)

                if isinstance(data_bgsubtracted, np.ma.MaskedArray):
                    data_background = sep.Background(data_bgsubtracted.data, mask=data_bgsubtracted.mask)
                    data_background.subfrom(data_bgsubtracted.data)

                else:
                    data_background = sep.Background(data_bgsubtracted)
                    data_background.subfrom(data_bgsubtracted)

                steps.append('background')

            guide_cards = None
            if guide_profiles_stale and ops_daemon is not None:
                guide_cards = notify_guide_profiles(header, data_bgsubtracted, log_name, ops_daemon)
                steps.append('guiding')

            move_bscale_bzero(header)
            if telescope_cards:
                add_telescope_header(header, info, log_name, telescope_cards, telescope_query_timeout)

            if environment_cards:
                add_environment_header(header, info, log_name, environment_cards,
                                       environment_daemon, environment_query_timeout)

            output_filename = None
            if info['output_archive']:
                date = datetime.datetime.strptime(header['DATE-OBS'], '%Y-%m-%dT%H:%M:%S.%f')
                output_filename = info['output_prefix'] + '-'

                if camera_config['output_filename_add_camera_id']:
                    output_filename += info['camera_id'] + '-'

                if camera_config['output_filename_add_filter']:
                    output_filename += header.get('FILTER', 'UNKNOWN') + '-'

                output_filename += date.strftime('%Y%m%d%H%M%S%f')[:-3] + '.fits'

            add_pipeline_header(header, info, output_filename)

            if info['intensity_stats_enabled']:
                add_intensity_stats_header(header, data_raw, data_cropped, camera_config, master_bias_stats)
                steps.append('intstats')

            binning = 1
            binning_card = camera_config.get('binning_card', None)
            if binning_card and binning_card in header:
                binning = header[binning_card]

            objects = None
            if calculate_objects:
                objects = detect_objects(data_bgsubtracted, data_background.globalrms, camera_config,
                                         binning, crop_x, crop_y)
                steps.append('objdetect')

            if info['hfd_enabled']:
                add_hfd_header(header, objects)
                steps.append('hfd')

            if info['wcs_enabled'] or info['rotation_enabled']:
                add_wcs_header(header, objects, binning, info['rotation_enabled'], camera_config.get('wcs', None))
                steps.append('wcs')
                if info['rotation_enabled']:
                    steps.append('rotation')

            if info['hfd_grid_enabled']:
                add_hfd_grid_header(header, objects, camera_config)
                steps.append('hfd_grid')

            ops_cards = None
            if ops_daemon is not None:
                ops_cards = notify_opsd(header, ops_daemon, log_name)

            if ops_cards or guide_cards:
                cards = []
                if ops_cards:
                    cards.extend(ops_cards)
                if guide_cards:
                    cards.extend(guide_cards)
                add_operations_header(header, cards)

            if output_filename:
                # Replace header directly in the input file before moving to the final output location
                write_fits_header(loadpath, header, header_block_capacity)

                savepath = os.path.join(camera_config['output_data_path'], info['output_subdirectory'], output_filename)
                shutil.move(loadpath, savepath)

                log.info(log_name, 'Saved frame ' + output_filename)
                sys.stdout.flush()

            process_time = round((datetime.datetime.utcnow() - start).total_seconds(), 1)
            print(f'processed {loadpath} in {process_time}s ({", ".join(steps)})', flush=True)

            if preview_stale:
                preview_queue.put((data_raw, header))

            if dashboard_stale:
                update_dashboard(header, data_cropped, camera_config['ssh_root_path'], camera_config['dashboard'],
                                 steps, process_time)

        except Exception:
            print(f'Unexpected exception while processing {loadpath}: ')
            traceback.print_exc(file=sys.stdout)
            sys.stdout.flush()
        finally:
            try:
                if os.path.exists(loadpath):
                    os.remove(loadpath)
            except Exception:
                print(f'Failed to delete temporary frame: {loadpath}', flush=True)


def write_fits_header(loadpath, header, header_block_capacity):
    # Is there enough space to fit all the cards?
    # We need at least one line free for the END card
    truncated = 0
    while len(header) > 36 * header_block_capacity - 1:
        truncated += 1
        header.pop()

    if truncated:
        # Make space for the truncation warning
        header.pop()
        header.append(fits.Card('COMMENT', f' WARNING: {truncated + 1} header cards truncated', ''))

    # It is very likely that we will have more than 36 cards of unused capacity if WCS solutions
    # are disabled. The fits specification requires that the END card is within the block immediately
    # before the data starts, so pad the header with blank cards to do this.
    padding = max(0, 36 * (header_block_capacity - 1) - len(header))
    for _ in range(padding):
        header.append(fits.Card(None, None))

    with open(loadpath, 'r+b') as f:
        f.write(header.tostring().encode('ascii'))


class PipelineWorkerDaemon:
    """Camera-specific pipeline worker"""
    def __init__(self, config, camera_id):
        self._config = config
        self._camera_config = config.cameras[camera_id]
        self._camera_daemon = getattr(daemons, self._camera_config['worker_daemon'])

        # Block if all workers are saturated so the pyro notification
        # can timeout and propagate back to stop the camera exposures
        self._processing_queue = multiprocessing.Queue(maxsize=self._camera_config['worker_processes'])

        # Avoid race conditions between worker processes
        start_timestamp = Time.now().strftime('%Y-%m-%dT%H:%M:%S.%f').encode('ascii')
        self._preview_timestamp = multiprocessing.Array('c', start_timestamp)
        self._guide_timestamp = multiprocessing.Array('c', start_timestamp)
        self._dashboard_timestamp = multiprocessing.Array('c', start_timestamp)

        self._preview_config = {getattr(IP, c['remote_machine']): c for c in config.preview_client_config}
        self._preview_registrations = multiprocessing.Manager().list()
        self._preview_registration_lock = threading.Lock()
        self._preview_lock = multiprocessing.Lock()

        for _ in range(self._camera_config['worker_processes']):
            multiprocessing.Process(target=process_frames, daemon=True, args=(
                self._processing_queue, self._camera_config, self._preview_timestamp,
                self._preview_registrations, self._preview_lock,
                self._guide_timestamp, config.guiding_min_interval,
                self._dashboard_timestamp, config.log_name,
                config.environment_cards, config.environment_daemon_name, config.environment_query_timeout,
                config.telescope_cards, config.telescope_query_timeout,
                config.ops_daemon_name)).start()

    @Pyro4.expose
    def create_output_subdirectory(self, name):
        """Attempts to set the output data subdirectory inside the output_data_path camera config variable
           The given directory name will be created and/or checked to ensure it is writable.
           Returns a CommandStatus indicating the result
        """
        if not pyro_client_matches([self._config.daemon.host]):
            return CommandStatus.InvalidControlIP

        path = os.path.join(self._camera_config['output_data_path'], name)
        try:
            os.makedirs(path, exist_ok=True)
        except FileExistsError:
            pass
        except Exception:
            print('failed to create night dir with exception:')
            traceback.print_exc(file=sys.stdout)
            return CommandStatus.NewDirectoryFailed

        # Throws on error
        try:
            with tempfile.TemporaryFile(dir=path):
                return CommandStatus.Succeeded
        except Exception:
            return CommandStatus.DirectoryNotWritable

    @Pyro4.expose
    def notify_frame(self, info):
        """Called by the master pipeline daemon to notify that a new frame is ady for processing.
           filename is specified relative to the incoming_data_path config value."""
        if not pyro_client_matches([self._config.daemon.host]):
            return

        self._processing_queue.put(info)

    @Pyro4.expose
    def register_preview(self, host, registration):
        """Register a ds9 window for previews"""
        if not pyro_client_matches([self._config.daemon.host]):
            return CommandStatus.InvalidControlIP

        preview_config = self._preview_config.get(host, None)
        if preview_config is None:
            return CommandStatus.InvalidControlIP

        with self._preview_registration_lock:
            if len(self._preview_registrations) >= self._camera_config['preview_max_instances']:
                return CommandStatus.TooManyPreviews

            try:
                # Check whether we are able to connect back to the client, otherwise fail early
                # This requires that key-based ssh access has been configured from the server to the calling client!
                with SSHClient() as ssh:
                    ssh.load_system_host_keys(os.path.join(self._camera_config['ssh_root_path'], 'known_hosts'))
                    key_path = os.path.join(self._camera_config['ssh_root_path'], preview_config['key'])
                    with open(key_path, encoding='utf-8') as keyfile:
                        key = RSAKey.from_private_key(keyfile)

                    remote_ip = getattr(IP, preview_config['remote_machine'])
                    ssh.connect(remote_ip, username=preview_config['remote_user'], pkey=key)

                    registration['rockit.key'] = preview_config['key']
                    registration['rockit.remote_ip'] = remote_ip
                    registration['rockit.remote_user'] = preview_config['remote_user']
                    registration['rockit.remote_command'] = preview_config['remote_command']
                    registration['rockit.remote_path'] = os.path.join(preview_config["remote_path"],
                                                                      f'preview-{registration["rockit.client"]}.fits')
            except Exception:
                return CommandStatus.SSHConnectionFailed

            self._preview_registrations.append(registration)
        return CommandStatus.Succeeded


if __name__ == '__main__':
    parser = argparse.ArgumentParser(description='Data pipeline worker daemon')
    parser.add_argument('config', help='Path to configuration json file')
    parser.add_argument('camera_id', help='Camera ID')
    args = parser.parse_args()
    c = Config(args.config, validate_camera=args.camera_id)
    daemon = getattr(daemons, c.cameras[args.camera_id]['worker_daemon'])
    daemon.launch(PipelineWorkerDaemon(c, args.camera_id))