#ifdef HAVE_CONFIG_H
#include <config.h>
#endif

#include <stdio.h>
#include <pslib.h>

#include "ppMops.h"

static void mopsDetectionsFree(ppMopsDetections *det)
{
    psFree(det->raBoresight);
    psFree(det->decBoresight);
    psFree(det->filter);
    psFree(det->x);
    psFree(det->y);
    psFree(det->ra);
    psFree(det->dec);
    psFree(det->raErr);
    psFree(det->decErr);
    psFree(det->mag);
    psFree(det->magErr);
    psFree(det->extended);
    psFree(det->angle);
    psFree(det->angleErr);
    psFree(det->length);
    psFree(det->lengthErr);
    psFree(det->flags);
    psFree(det->diffSkyfileId);
    psFree(det->naxis1);
    psFree(det->naxis2);
    psFree(det->mask);
    return;
}

ppMopsDetections *ppMopsDetectionsAlloc(long num)
{
    ppMopsDetections *det = psAlloc(sizeof(ppMopsDetections)); // Detections, to return
    psMemSetDeallocator(det, (psFreeFunc)mopsDetectionsFree);

    det->raBoresight = NULL;
    det->decBoresight = NULL;
    det->filter = NULL;
    det->airmass = NAN;
    det->exptime = NAN;
    det->posangle = NAN;
    det->alt = NAN;
    det->az = NAN;
    det->mjd = NAN;
    det->seeing = NAN;
    det->num = 0;
    det->x = psVectorAllocEmpty(num, PS_TYPE_F32);
    det->y = psVectorAllocEmpty(num, PS_TYPE_F32);
    det->ra = psVectorAllocEmpty(num, PS_TYPE_F64);
    det->dec = psVectorAllocEmpty(num, PS_TYPE_F64);
    det->raErr = psVectorAllocEmpty(num, PS_TYPE_F64);
    det->decErr = psVectorAllocEmpty(num, PS_TYPE_F64);
    det->mag = psVectorAllocEmpty(num, PS_TYPE_F32);
    det->magErr = psVectorAllocEmpty(num, PS_TYPE_F32);
    det->extended = psVectorAllocEmpty(num, PS_TYPE_F32);
    det->angle = psVectorAllocEmpty(num, PS_TYPE_F32);
    det->angleErr = psVectorAllocEmpty(num, PS_TYPE_F32);
    det->length = psVectorAllocEmpty(num, PS_TYPE_F32);
    det->lengthErr = psVectorAllocEmpty(num, PS_TYPE_F32);
    det->flags = psVectorAllocEmpty(num, PS_TYPE_U32);
    det->diffSkyfileId = psVectorAllocEmpty(num, PS_TYPE_S64);
    det->naxis1 = psVectorAllocEmpty(num, PS_TYPE_S32);
    det->naxis2 = psVectorAllocEmpty(num, PS_TYPE_S32);
    det->mask = psVectorAllocEmpty(num, PS_TYPE_U8);

    return det;
}


ppMopsDetections *ppMopsDetectionsRealloc(ppMopsDetections *det, long num)
{
    det->x = psVectorRealloc(det->x, num);
    det->y = psVectorRealloc(det->y, num);
    det->ra = psVectorRealloc(det->ra, num);
    det->dec = psVectorRealloc(det->dec, num);
    det->raErr = psVectorRealloc(det->raErr, num);
    det->decErr = psVectorRealloc(det->decErr, num);
    det->mag = psVectorRealloc(det->mag, num);
    det->magErr = psVectorRealloc(det->magErr, num);
    det->extended = psVectorRealloc(det->extended, num);
    det->angle = psVectorRealloc(det->angle, num);
    det->angleErr = psVectorRealloc(det->angleErr, num);
    det->length = psVectorRealloc(det->length, num);
    det->lengthErr = psVectorRealloc(det->lengthErr, num);
    det->flags = psVectorRealloc(det->flags, num);
    det->diffSkyfileId = psVectorRealloc(det->diffSkyfileId, num);
    det->naxis1 = psVectorRealloc(det->naxis1, num);
    det->naxis2 = psVectorRealloc(det->naxis2, num);
    det->mask = psVectorRealloc(det->mask, num);

    return det;
}


bool ppMopsDetectionsAdd(ppMopsDetections *det, float x, float y, double ra, double dec,
                         double raErr, double decErr, float mag, float magErr, float extended,
                         float angle, float angleErr, float length, float lengthErr,
                         psU32 flags, psS64 diffSkyfileId, int naxis1, int naxis2)
{
    psVectorAppend(det->x, x);
    psVectorAppend(det->y, y);
    psVectorAppend(det->ra, ra);
    psVectorAppend(det->dec, dec);
    psVectorAppend(det->raErr, raErr);
    psVectorAppend(det->decErr, decErr);
    psVectorAppend(det->mag, mag);
    psVectorAppend(det->magErr, magErr);
    psVectorAppend(det->extended, extended);
    psVectorAppend(det->angle, angle);
    psVectorAppend(det->angleErr, angleErr);
    psVectorAppend(det->length, length);
    psVectorAppend(det->lengthErr, lengthErr);
    psVectorAppend(det->flags, flags);
    psVectorAppend(det->diffSkyfileId, diffSkyfileId);
    psVectorAppend(det->naxis1, naxis1);
    psVectorAppend(det->naxis2, naxis2);
    psVectorAppend(det->mask, 0);
    return true;
}


bool ppMopsDetectionsCopySingle(ppMopsDetections *target, const ppMopsDetections *source, long index)
{
    psVectorAppend(target->x, source->x->data.F32[index]);
    psVectorAppend(target->y, source->y->data.F32[index]);
    psVectorAppend(target->ra, source->ra->data.F64[index]);
    psVectorAppend(target->dec, source->dec->data.F64[index]);
    psVectorAppend(target->raErr, source->raErr->data.F64[index]);
    psVectorAppend(target->decErr, source->decErr->data.F64[index]);
    psVectorAppend(target->mag, source->mag->data.F32[index]);
    psVectorAppend(target->magErr, source->magErr->data.F32[index]);
    psVectorAppend(target->extended, source->extended->data.F32[index]);
    psVectorAppend(target->angle, source->angle->data.F32[index]);
    psVectorAppend(target->angleErr, source->angleErr->data.F32[index]);
    psVectorAppend(target->length, source->length->data.F32[index]);
    psVectorAppend(target->lengthErr, source->lengthErr->data.F32[index]);
    psVectorAppend(target->flags, source->flags->data.U32[index]);
    psVectorAppend(target->diffSkyfileId, source->diffSkyfileId->data.S64[index]);
    psVectorAppend(target->naxis1, source->naxis1->data.S32[index]);
    psVectorAppend(target->naxis2, source->naxis2->data.S32[index]);
    psVectorAppend(target->mask, 0);
    target->num++;
    return true;
}


bool ppMopsDetectionsPurge(ppMopsDetections *det)
{
    long num = 0;
    for (long i = 0; i < det->num; i++) {
        if (!det->mask->data.U8[i]) {
            if (i == num) {
                // No need to copy
                num++;
                continue;
            }
            det->x->data.F32[num] = det->x->data.F32[i];
            det->y->data.F32[num] = det->y->data.F32[i];
            det->ra->data.F64[num] = det->ra->data.F64[i];
            det->dec->data.F64[num] = det->dec->data.F64[i];
            det->raErr->data.F64[num] = det->raErr->data.F64[i];
            det->decErr->data.F64[num] = det->decErr->data.F64[i];
            det->mag->data.F32[num] = det->mag->data.F32[i];
            det->magErr->data.F32[num] = det->magErr->data.F32[i];
            det->extended->data.F32[num] = det->extended->data.F32[i];
            det->angle->data.F32[num] = det->angle->data.F32[i];
            det->angleErr->data.F32[num] = det->angleErr->data.F32[i];
            det->length->data.F32[num] = det->length->data.F32[i];
            det->lengthErr->data.F32[num] = det->lengthErr->data.F32[i];
            det->flags->data.U32[num] = det->flags->data.U32[i];
            det->diffSkyfileId->data.S64[num] = det->diffSkyfileId->data.S64[i];
            det->naxis1->data.S32[num] = det->naxis1->data.S32[i];
            det->naxis2->data.S32[num] = det->naxis2->data.S32[i];
            det->mask->data.U8[num] = 0;
            num++;
        }
    }
    det->x->n = num;
    det->y->n = num;
    det->ra->n = num;
    det->dec->n = num;
    det->raErr->n = num;
    det->decErr->n = num;
    det->mag->n = num;
    det->magErr->n = num;
    det->extended->n = num;
    det->angle->n = num;
    det->angleErr->n = num;
    det->length->n = num;
    det->lengthErr->n = num;
    det->flags->n = num;
    det->diffSkyfileId->n = num;
    det->naxis1->n = num;
    det->naxis2->n = num;
    det->mask->n = num;
    det->num = num;
    return true;
}

