# include "psastroInternal.h"
# define ELIXIR_MODE 1

psArray *psastroLoadRefstars (pmConfig *config) {

    int fd;
    bool status;

    // select the current recipe
    psMetadata *recipe  = psMetadataLookupPtr (NULL, config->recipes, PSASTRO_RECIPE);

    // DVO APIs expect decimal degrees
    float RAmin  = DEG_RAD*psMetadataLookupF32(NULL, recipe, "RA_MIN");
    float RAmax  = DEG_RAD*psMetadataLookupF32(NULL, recipe, "RA_MAX");
    float DECmin = DEG_RAD*psMetadataLookupF32(NULL, recipe, "DEC_MIN");
    float DECmax = DEG_RAD*psMetadataLookupF32(NULL, recipe, "DEC_MAX");

    // extra field fraction to add
    double fieldPadding = psMetadataLookupF32 (&status, recipe, "PSASTRO.FIELD.PADDING");
    PS_ASSERT (status, NULL);

    float dRA = RAmax - RAmin;
    RAmin -= dRA * fieldPadding;
    RAmax += dRA * fieldPadding;

    float dDEC = DECmax - DECmin;
    DECmin -= dDEC * fieldPadding;
    DECmax += dDEC * fieldPadding;

    // grab the PSASTRO.CATDIR name from the PSASTRO recipe
    char *catdir_recipe = psMetadataLookupStr(&status, recipe, "PSASTRO.CATDIR");
    psAssert (catdir_recipe, "Need a recipe for the catdir!");

    // substitute abstract name with concrete name, if present in PSASTRO.CATDIRS
    psMetadata *catdirs = psMetadataLookupMetadata(&status, config->complete, "PSASTRO.CATDIRS"); // List of cameras
    if (!catdirs) {
        psError(PS_ERR_UNEXPECTED_NULL, true, "Unable to find PSASTRO.CATDIRS in the system configuration.\n");
        return false;
    }
    char *catdir_virtual = psMetadataLookupStr(&status, catdirs, catdir_recipe);
    char *catdir = (catdir_virtual == NULL) ? catdir_recipe : catdir_virtual;

    // convert the uri to a real filename (ie, path://foobar)
    psString CATDIR = pmConfigConvertFilename(catdir, config, false, false); // Resolved filename
    PS_ASSERT (CATDIR, NULL);

    char *getstarCommand = psStringCopy(psMetadataLookupStr(NULL, recipe, "DVO.GETSTAR"));
    PS_ASSERT (getstarCommand, NULL);

    char *outformat = psMetadataLookupStr(NULL, recipe, "DVO.GETSTAR.OUTFORMAT");
    PS_ASSERT (outformat, NULL);

    char *photcode = psMetadataLookupStr(NULL, recipe, "DVO.GETSTAR.PHOTCODE");
    PS_ASSERT (photcode, NULL);

    float MAGmax = psMetadataLookupF32(NULL, recipe, "DVO.GETSTAR.MAG.MAX");

    // issue the following command:
    // getstar -region RAmin RAmax DECmin DECmax
    char tempFile[64];
    sprintf (tempFile, "/tmp/psastro.XXXXXX");
    if ((fd = mkstemp (tempFile)) == -1) {
        psError(PSASTRO_ERR_REFSTARS, true, "error creating temp output file\n");
        psFree(CATDIR);
        return NULL;
    }
    close (fd);

    psTimerStart ("psastro");

    // supply a known output format (for CATALOG output) so the code below knows what to read
    if (ELIXIR_MODE) {
        psStringAppend (&getstarCommand, " -D CATMODE mef -D CATFORMAT elixir");
    } else {
        psStringAppend (&getstarCommand, " -D CATMODE mef -D CATFORMAT panstarrs");
    }

    // check for default name (use .ptolemyrc), or use specified CATDIR
    if (strcasecmp(CATDIR, "NONE")) {
        psStringAppend (&getstarCommand, " -D CATDIR %s", CATDIR);
    }
    psFree(CATDIR);

    // supply the max magnitude, the output format, and the photcode
    if (strcasecmp (photcode, "NONE")) {
        psStringAppend (&getstarCommand, " -maglim %f", MAGmax);
    }
    if (strcasecmp (photcode, "NONE")) {
        psStringAppend (&getstarCommand, " -photcode %s", photcode);
    }
    psStringAppend (&getstarCommand, " -format %s", outformat);

    // add region and output filename
    psStringAppend (&getstarCommand, " -region %f %f %f %f -o %s", RAmin, DECmin, RAmax, DECmax, tempFile);
    psTrace ("psastro", 3, "%s\n", getstarCommand);

    // XXX use psPipe: catch stderr, stdout, allow for Nsec timeout...
    // use fork to add timeout capability
    status = system (getstarCommand);
    if (status) {
        psError(PSASTRO_ERR_REFSTARS, true, "error loading reference data\n");
        return NULL;
    }
    psFree (getstarCommand);

    psLogMsg ("psastro", 3, "ran getstar : %f sec\n", psTimerMark ("psastro"));

    // the output from getstar is a file with the Average table
    psFits *fits = psFitsOpen (tempFile, "r");

    psTimerStart ("psastro");

    psArray *refstars = NULL;
    if (!strcmp (outformat, "CATALOG")) {
      refstars = psastroReadGetstarCatalog (fits);
    }
    if (!strcmp (outformat, "PS1_DEV_0")) {
      refstars = psastroReadGetstar_PS1_DEV_0 (fits);
    }
    if (refstars == NULL) {
        psError(PSASTRO_ERR_REFSTARS, true, "error reading reference data\n");
        psFitsClose (fits);
        return NULL;
    }

    psLogMsg ("psastro", 3, "loaded %ld reference stars : %f sec\n", refstars->n, psTimerMark ("psastro"));

    psTrace ("psastro", 3, "loaded %ld reference stars from (%10.6f,%10.6f) - (%10.6f,%10.6f)\n",
             refstars->n, RAmin, DECmin, RAmax, DECmax);

    psFitsClose (fits);
    unlink (tempFile);

    // dump or plot the available refstars
    if (psTraceGetLevel("psastro.dump") > 0) {
        psastroDumpRefstars (refstars, "refstars.dat");
    }

    if (psTraceGetLevel("psastro.plot") > 0) {
        psastroPlotRefstars (refstars, recipe);
    }

    return refstars;
}

psArray *psastroReadGetstarCatalog (psFits *fits) {

    bool status;

    if (ELIXIR_MODE) {
        psFitsMoveExtName (fits, "DVO_AVERAGE_ELIXIR");
    } else {
        psFitsMoveExtName (fits, "DVO_AVERAGE_PANSTARRS");
    }

    long numSources = psFitsTableSize(fits); // Number of sources in table

    // convert the Average table to the pmAstromObj entries
    psArray *refstars = psArrayAllocEmpty (numSources);
    for (int i = 0; i < numSources; i++) {
        pmAstromObj *ref = pmAstromObjAlloc ();

        psMetadata *row = psFitsReadTableRow(fits, i); // Table row

        // DVO tables are stored in degrees
        if (ELIXIR_MODE) {
            ref->sky->r   = RAD_DEG*psMetadataLookupF32 (&status, row, "RA");
            ref->sky->d   = RAD_DEG*psMetadataLookupF32 (&status, row, "DEC");
            ref->Mag      = 0.001*psMetadataLookupS32 (&status, row, "MAG");  // ELIXIR uses millimags
        } else {
            ref->sky->r   = RAD_DEG*psMetadataLookupF64 (&status, row, "RA");
            ref->sky->d   = RAD_DEG*psMetadataLookupF64 (&status, row, "DEC");
            ref->Mag      = psMetadataLookupF32 (&status, row, "MAG"); // PANSTARRS uses mags
        }

        // XXX VERY temporary hack to avoid M31 bulge
        if ((fabs(ref->sky->r - 0.186438) < 0.002) && (fabs(ref->sky->d - 0.720270) < 0.002)) {
          psFree (ref);
          psFree (row);
          continue;
        }

        psArrayAdd (refstars, 100, ref);
        psFree (ref);
        psFree (row);
    }
    return refstars;
}

psArray *psastroReadGetstar_PS1_DEV_0 (psFits *fits) {

    bool status;

    psFitsMoveExtName (fits, "GETSTAR_PS1_DEV_0");

    long numSources = psFitsTableSize(fits); // Number of sources in table

    // convert the Average table to the pmAstromObj entries
    psArray *refstars = psArrayAllocEmpty (numSources);
    for (int i = 0; i < numSources; i++) {
        pmAstromObj *ref = pmAstromObjAlloc ();

        psMetadata *row = psFitsReadTableRow(fits, i); // Table row

        ref->sky->r   = RAD_DEG*psMetadataLookupF32 (&status, row, "RA");
        ref->sky->d   = RAD_DEG*psMetadataLookupF32 (&status, row, "DEC");
        ref->Mag      = psMetadataLookupF32 (&status, row, "MAG");

        // XXX VERY temporary hack to avoid M31 bulge
        if ((fabs(ref->sky->r - 0.186438) < 0.002) && (fabs(ref->sky->d - 0.720270) < 0.002)) {
          psFree (ref);
          psFree (row);
          continue;
        }

        psArrayAdd (refstars, 100, ref);
        psFree (ref);
        psFree (row);
    }
    return refstars;
}
