# include "psastroInternal.h"
// leak free 2006.04.27

bool psastroConvertFPA (pmFPA *fpa, psMetadata *recipe) {

    pmChip *chip;
    pmCell *cell;
    pmReadout *readout;
    pmFPAview *view = pmFPAviewAlloc (0);

    while ((chip = pmFPAviewNextChip (view, fpa, 1)) != NULL) {
        psTrace ("psastro", 4, "Chip %d: %x %x\n", view->chip, chip->file_exists, chip->process);
        if (!chip->process || !chip->file_exists) { continue; }

        while ((cell = pmFPAviewNextCell (view, fpa, 1)) != NULL) {
            psTrace ("psastro", 4, "Cell %d: %x %x\n", view->cell, cell->file_exists, cell->process);
            if (!cell->process || !cell->file_exists) { continue; }

            // process each of the readouts
            while ((readout = pmFPAviewNextReadout (view, fpa, 1)) != NULL) {
                if (! readout->data_exists) { continue; }

                psastroConvertReadout (readout, recipe);
            }
        }
    }
    psFree (view);
    return true;
}

// pass/apply the WCS information?
bool psastroConvertReadout (pmReadout *readout, psMetadata *recipe) {

    bool status;

    // PSPHOT.SOURCES carries the pmSource objects (from psphot analysis or loaded externally)
    psArray *sources = psMetadataLookupPtr (NULL, readout->analysis, "PSPHOT.SOURCES");
    if (sources == NULL)
        return false;

    // convert the pmSource objects into pmAstromObj objects (drop !STAR and SATSTAR?)
    psArray *inStars = pmSourceToAstromObj (sources);

    // sort in ascending magnitude order
    // psArraySort (inStars, psastroSortByMag);
    // psVector *index = psArraySortIndex (sources, pmSourceSortBySN);
    psVector *index = psArraySortIndex (NULL, inStars, psastroSortByMag);

    // XXX need to exit gracefully is inStars->n is 0 (or 1?)

    // we are going to select the brighter Nmax subset for astrometry 
    int nMax = psMetadataLookupS32 (&status, recipe, "PSASTRO.MAX.NRAW");
    if (!status || !nMax) nMax = inStars->n;

    // we are going to select the brighter Nmax subset for astrometry 
    float iMagMax = psMetadataLookupF32 (&status, recipe, "PSASTRO.MAX.INST.MAG.RAW");
    float iMagMin = psMetadataLookupF32 (&status, recipe, "PSASTRO.MIN.INST.MAG.RAW");

    // we are going to select the brighter Nmax subset for astrometry 
    pmSourceMode skip = PM_SOURCE_MODE_DEFAULT;
    char *ignoreList = psMetadataLookupStr (&status, recipe, "PSASTRO.IGNORE");
    if (ignoreList != NULL) {
      psArray *list = psStringSplitArray (ignoreList, ",", false);
      for (int i = 0; i < list->n; i++) {
	pmSourceMode mode = pmSourceModeFromString (list->data[i]);
	if (mode == PM_SOURCE_MODE_DEFAULT) {
	  psWarning ("unknown source mode in PSASTRO.IGNORE, skipping");
	  continue;
	}
	skip |= mode;
      }
      psFree (list);
    }

    // choose the first nMax sources
    int j = 0;
    psArray *rawStars = psArrayAlloc (PS_MIN (nMax, inStars->n));

    float mMin = +100.0;
    float mMax = -100.0;
    int nModeSkip = 0;
    int nFaintSkip = 0;
    int nBrightSkip = 0;
    int nInfSkip = 0;

    for (int i = 0; (i < inStars->n) && (j < rawStars->n); i++) {
	int n = index->data.S32[i];
	pmSource *source = sources->data[n];
	// XXX this needs to be flexible
	if (source->mode & skip) { 
	  nModeSkip ++;
	  continue;
	}
	if ((iMagMax != 0.0) && (source->psfMag > iMagMax)) {
	  nFaintSkip ++;
	  continue;
	}
	if ((iMagMin != 0.0) && (source->psfMag < iMagMin)) {
	  nBrightSkip ++;
	  continue;
	}
	if (!isfinite(source->psfMag)) {
	  nInfSkip ++;
	  continue;
	}
	mMin = PS_MIN (mMin, source->psfMag);
	mMax = PS_MAX (mMax, source->psfMag);
	rawStars->data[j] = psMemIncrRefCounter (inStars->data[n]);
	j++;
    }
    rawStars->n = j;

    psMetadataAdd (readout->analysis, PS_LIST_TAIL, "PSASTRO.RAWSTARS", PS_DATA_ARRAY, "astrometry objects", rawStars);
    psLogMsg ("psastro", 4, "loaded %ld sources, using %ld of %ld good sources (inst mag: %f to %f)\n", sources->n, rawStars->n, inStars->n, mMin, mMax);
    psLogMsg ("psastro", 4, "skip reasons: mode: %d, faint: %d, bright: %d, inf: %d\n", nModeSkip, nFaintSkip, nBrightSkip, nInfSkip);

    psFree (index);
    psFree (inStars);
    psFree (rawStars);

    return true;
}

// select a magnitude range?
psArray *pmSourceToAstromObj (psArray *sources) {

    psArray *objects = psArrayAllocEmpty (sources->n);

    for (int i = 0; i < sources->n; i++) {
        pmSource *source = sources->data[i];

        // only accept the PSF sources?
	// XXX drop SATSTAR?
        // if (source->type != PM_SOURCE_TYPE_STAR) continue;

        pmModel *model = source->modelPSF;
        if (model == NULL) continue;

        psF32 *PAR = model->params->data.F32;

        pmAstromObj *obj = pmAstromObjAlloc ();

        // is the source magnitude calibrated in any sense?
        obj->pix->x = PAR[PM_PAR_XPOS];
        obj->pix->y = PAR[PM_PAR_YPOS];
        obj->Mag = source->psfMag;

        // XXX do we have the information giving the readout and cell offset?
        // for the moment, assume chip == cell == readout
        *obj->cell = *obj->pix;
        *obj->chip = *obj->cell;

        psArrayAdd (objects, 100, obj);
        psFree (obj);
    }
    return objects;
}

// sort by Mag (ascending)
int psastroSortByMag (const void *a, const void *b)
{
    const pmAstromObj *A = a;
    const pmAstromObj *B = b;

    psF32 mA = (isfinite(A->Mag)) ? A->Mag : FLT_MAX;
    psF32 mB = (isfinite(B->Mag)) ? B->Mag : FLT_MAX;

    psF32 diff = mA - mB;
    if (diff > +FLT_EPSILON) return (+1);
    if (diff < -FLT_EPSILON) return (-1);
    return (0);
}

