Changeset 34281

Timestamp:

Aug 2, 2012, 2:47:03 PM (14 years ago)

Author:

Serge CHASTEL

Message:

Rollback to version 34278

Location:

tags/ipp-20120802/ppTranslate

Files:

• . (modified) (1 prop)
• src/ppMops.c (modified) (2 diffs)
• src/ppMops.h (modified) (3 diffs)
• src/ppMopsArguments.c (modified) (1 diff)
• src/ppMopsGetSkyChipPsfVersion.c (modified) (1 diff)
• src/ppMopsMerge.c (modified) (1 diff)
• src/ppMopsRead.c (modified) (1 diff)
• src/ppMopsWrite.c (modified) (1 diff)

tags/ipp-20120802/ppTranslate/src/ppMops.c

@@ -67 +67 @@
     }
 
+
     psArray *detections = ppMopsRead(args); // Detections from each input
     if (!detections) {
@@ -91 +92 @@
     psLibFinalize();
 
-/*     fprintf (stderr, "found %d leaks at %s\n",  */
-/*      psMemCheckLeaks2 (0, */
-/*              NULL, stdout, false, 500), "ppMops"); */
+    fprintf (stderr, "found %d leaks at %s\n", 
+        psMemCheckLeaks2 (0,
+                NULL, stdout, false, 500), "ppMops");
 
     return PS_EXIT_SUCCESS;
 }
 
+
+#if 0
+    ps
+
+
+    psArray *detections = NULL;         // Detections
+    psMetadata *header = NULL;          // Header for detections
+    {
+        psFits *fits = psFitsOpen(data->detections, "r"); // FITS file
+        header = psFitsReadHeader(NULL, fits);
+        if (!header) {
+            psErrorStackPrint(stderr, "Unable to read header");
+            psFitsClose(fits);
+            psFree(data);
+            exit(PS_EXIT_DATA_ERROR);
+        }
+        if (!psFitsMoveExtName(fits, IN_EXTNAME)) {
+            psErrorStackPrint(stderr, "Unable to move to extension %s", IN_EXTNAME);
+            psFitsClose(fits);
+            psFree(header);
+            psFree(data);
+            exit(PS_EXIT_DATA_ERROR);
+        }
+        detections = psFitsReadTable(fits);
+        psFitsClose(fits);
+        if (!detections) {
+            psErrorStackPrint(stderr, "Unable to read detections");
+            psFree(data);
+            psFree(header);
+            exit(PS_EXIT_DATA_ERROR);
+        }
+    }
+
+    // Translate the columns
+    int numIn = detections->n;        // Number of rows in input
+    int numOut = 0;                   // Number of rows in output
+    psArray *output = psArrayAllocEmpty(numIn); // Output table
+    double plateScale = 0.0;                 // Average plate scale
+    for (int i = 0; i < numIn; i++) {
+        psMetadata *inRow = detections->data[i]; // Input row
+
+        double ra = psMetadataLookupF64(NULL, inRow, "RA_PSF");
+        double dec = psMetadataLookupF64(NULL, inRow, "DEC_PSF");
+        double mag = psMetadataLookupF64(NULL, inRow, "PSF_INST_MAG") + data->zp;
+        double magErr = psMetadataLookupF64(NULL, inRow, "PSF_INST_MAG_SIG");
+        float ext = psMetadataLookupF32(NULL, inRow, "EXT_NSIGMA");
+        double xErr = psMetadataLookupF64(NULL, inRow, "X_PSF_SIG");
+        double yErr = psMetadataLookupF64(NULL, inRow, "Y_PSF_SIG");
+        double scale = psMetadataLookupF64(NULL, inRow, "PLTSCALE");
+        double angle = psMetadataLookupF64(NULL, inRow, "POSANGLE");
+        psU32 flags = psMetadataLookupU32(NULL, inRow, "FLAGS");
+
+        if (!isfinite(mag) || !isfinite(magErr) || (flags & SOURCE_MASK) ||
+            ((flags & PM_SOURCE_MODE_DEFECT) && !(flags & PM_SOURCE_MODE_MOMENTS_FAILURE))
+            ) {
+            // DEFECT can be due to bad moments
+            continue;
+        }
+
+        psMetadata *outRow = output->data[numOut] = psMetadataAlloc(); // Output row
+
+        numOut++;
+        plateScale += scale;
+
+        // XXX Not at all sure I've got the angles around the right way here...
+        double cosAngle = cos(angle), sinAngle = sin(angle);
+        double cosAngle2 = PS_SQR(cosAngle), sinAngle2 = PS_SQR(sinAngle);
+        double xErr2 = PS_SQR(xErr), yErr2 = PS_SQR(yErr);
+        double errScale = scale / 3600.0;
+        double raErr = errScale * sqrt(cosAngle2 * xErr2 + sinAngle2 * yErr2);
+        double decErr = errScale * sqrt(sinAngle2 * xErr2 + cosAngle2 * yErr2);
+
+        psMetadataAddF64(outRow, PS_LIST_TAIL, "RA_DEG", 0, "Right ascension (degrees)", ra);
+        psMetadataAddF64(outRow, PS_LIST_TAIL, "DEC_DEG", 0, "Declination (degrees)", dec);
+        psMetadataAddF64(outRow, PS_LIST_TAIL, "RA_SIG", 0, "Right ascension error (degrees)", raErr);
+        psMetadataAddF64(outRow, PS_LIST_TAIL, "DEC_SIG", 0, "Declination error (degrees)", decErr);
+        psMetadataAddF64(outRow, PS_LIST_TAIL, "MAG", 0, "Magnitude", mag);
+        psMetadataAddF64(outRow, PS_LIST_TAIL, "MAG_SIG", 0, "Magnitude error", magErr);
+        psMetadataAddU32(outRow, PS_LIST_TAIL, "FLAGS", 0, "Detection bit flags", flags);
+
+        // The below need fixing
+        psMetadataAddF64(outRow, PS_LIST_TAIL, "STARPSF", 0, "EXT_NSIGMA", ext);
+        psMetadataAddF64(outRow, PS_LIST_TAIL, "ANG", 0, "Position angle of trail (degrees)", 0.0);
+        psMetadataAddF64(outRow, PS_LIST_TAIL, "ANG_SIG", 0, "Position angle error (degrees)", 0.0);
+        psMetadataAddF64(outRow, PS_LIST_TAIL, "LEN", 0, "Length of trail (arcsec)", 0.0);
+        psMetadataAddF64(outRow, PS_LIST_TAIL, "LEN_SIG", 0, "Length error (arcsec)", 0.0);
+    }
+    output->n = numOut;
+    plateScale /= numOut;
+    psFree(detections);
+
+    // Translate the header
+    psMetadata *outHeader = psMetadataAlloc(); // Output header
+    ppMopsVersionHeader(outHeader);
+    {
+        const char *ra = psMetadataLookupStr(NULL, header, "FPA.RA");
+        const char *dec = psMetadataLookupStr(NULL, header, "FPA.DEC");
+        const char *filter = psMetadataLookupStr(NULL, header, "FPA.FILTER");
+        float airmass = psMetadataLookupF32(NULL, header, "FPA.AIRMASS");
+        float exptime = psMetadataLookupF32(NULL, header, "EXPTIME");
+        double angle = psMetadataLookupF64(NULL, header, "FPA.POSANGLE");
+        double alt = psMetadataLookupF64(NULL, header, "FPA.ALT");
+        double az = psMetadataLookupF64(NULL, header, "FPA.AZ");
+        psS64 imageid = psMetadataLookupS64(NULL, header, "IMAGEID");
+        double mjd = psMetadataLookupF64(NULL, header, "MJD-OBS") + exptime / 2.0 / 3600 / 24;
+
+        float psf = plateScale * 0.5 * (psMetadataLookupF32(NULL, header, "FWHM_MAJ") +
+                                        psMetadataLookupF32(NULL, header, "FWHM_MIN"));
+
+        psMetadataAddStr(outHeader, PS_LIST_TAIL, "RA", 0, "Right ascension of boresight", ra);
+        psMetadataAddStr(outHeader, PS_LIST_TAIL, "DEC", 0, "Declination of boresight", dec);
+        psMetadataAddF32(outHeader, PS_LIST_TAIL, "AIRMASS", 0, "Airmass of exposure", airmass);
+        psMetadataAddF64(outHeader, PS_LIST_TAIL, "MJD-OBS", 0, "MJD of exposure midpoint", mjd);
+        psMetadataAddStr(outHeader, PS_LIST_TAIL, "FILTER", 0, "Filter name", filter);
+        psMetadataAddF64(outHeader, PS_LIST_TAIL, "EXPTIME", 0, "Exposure time (sec)", exptime);
+        psMetadataAddF64(outHeader, PS_LIST_TAIL, "ROTANGLE", 0, "Rotator position angle", angle);
+        psMetadataAddF64(outHeader, PS_LIST_TAIL, "TEL_ALT", 0, "Telescope altitude", alt);
+        psMetadataAddF64(outHeader, PS_LIST_TAIL, "TEL_AZ", 0, "Telescope azimuth", az);
+        psMetadataAddS64(outHeader, PS_LIST_TAIL, "DIFFIMID", 0, "Difference image identifier", imageid);
+        psMetadataAddStr(outHeader, PS_LIST_TAIL, "FPA_ID", 0, "Exposure name", data->exp_name);
+        psMetadataAddS64(outHeader, PS_LIST_TAIL, "EXP_ID", 0, "Exposure identifier", data->exp_id);
+        psMetadataAddBool(outHeader, PS_LIST_TAIL, "POSITIVE", 0, "Positive subtraction?", data->direction);
+        psMetadataAddStr(outHeader, PS_LIST_TAIL, "OBSCODE", 0, "IAU Observatory code", OBSERVATORY_CODE);
+        psMetadataAddF32(outHeader, PS_LIST_TAIL, "STARPSF", 0, "Stellar PSF (arcsec)", psf);
+
+        // These are completely fake
+        psMetadataAddF32(outHeader, PS_LIST_TAIL, "LIMITMAG", 0, "Limiting magnitude (FAKE)", 99.0);
+        psMetadataAddF32(outHeader, PS_LIST_TAIL, "DE1", 0, "Detection efficiency (FAKE)", 0.0);
+        psMetadataAddF32(outHeader, PS_LIST_TAIL, "DE2", 0, "Detection efficiency (FAKE)", 0.0);
+        psMetadataAddF32(outHeader, PS_LIST_TAIL, "DE3", 0, "Detection efficiency (FAKE)", 0.0);
+        psMetadataAddF32(outHeader, PS_LIST_TAIL, "DE4", 0, "Detection efficiency (FAKE)", 0.0);
+        psMetadataAddF32(outHeader, PS_LIST_TAIL, "DE5", 0, "Detection efficiency (FAKE)", 0.0);
+        psMetadataAddF32(outHeader, PS_LIST_TAIL, "DE6", 0, "Detection efficiency (FAKE)", 0.0);
+        psMetadataAddF32(outHeader, PS_LIST_TAIL, "DE7", 0, "Detection efficiency (FAKE)", 0.0);
+        psMetadataAddF32(outHeader, PS_LIST_TAIL, "DE8", 0, "Detection efficiency (FAKE)", 0.0);
+        psMetadataAddF32(outHeader, PS_LIST_TAIL, "DE9", 0, "Detection efficiency (FAKE)", 0.0);
+        psMetadataAddF32(outHeader, PS_LIST_TAIL, "DE10", 0, "Detection efficiency (FAKE)", 0.0);
+    }
+    psFree(header);
+
+    // Write the new table
+    {
+        psFits *fits = psFitsOpen(data->output, "w"); // FITS file
+        if (!fits) {
+            psErrorStackPrint(stderr, "Unable to open %s for writing", data->output);
+            psFree(outHeader);
+            psFree(output);
+            psFree(data);
+            exit(PS_EXIT_SYS_ERROR);
+        }
+
+        if (numOut == 0) {
+            // Write dummy table
+            psMetadata *outRow = psMetadataAlloc(); // Dummy output row
+            psMetadataAddF64(outRow, PS_LIST_TAIL, "RA_DEG", 0, "Right ascension (degrees)", NAN);
+            psMetadataAddF64(outRow, PS_LIST_TAIL, "DEC_DEG", 0, "Declination (degrees)", NAN);
+            psMetadataAddF64(outRow, PS_LIST_TAIL, "RA_SIG", 0, "Right ascension error (degrees)", NAN);
+            psMetadataAddF64(outRow, PS_LIST_TAIL, "DEC_SIG", 0, "Declination error (degrees)", NAN);
+            psMetadataAddF64(outRow, PS_LIST_TAIL, "MAG", 0, "Magnitude", NAN);
+            psMetadataAddF64(outRow, PS_LIST_TAIL, "MAG_SIG", 0, "Magnitude error", NAN);
+            psMetadataAddU32(outRow, PS_LIST_TAIL, "FLAGS", 0, "Detection bit flags", 0);
+
+            // The below need fixing
+            psMetadataAddF64(outRow, PS_LIST_TAIL, "STARPSF", 0, "EXT_NSIGMA", NAN);
+            psMetadataAddF64(outRow, PS_LIST_TAIL, "ANG", 0, "Position angle of trail (degrees)", NAN);
+            psMetadataAddF64(outRow, PS_LIST_TAIL, "ANG_SIG", 0, "Position angle error (degrees)", NAN);
+            psMetadataAddF64(outRow, PS_LIST_TAIL, "LEN", 0, "Length of trail (arcsec)", NAN);
+            psMetadataAddF64(outRow, PS_LIST_TAIL, "LEN_SIG", 0, "Length error (arcsec)", NAN);
+            if (!psFitsWriteTableEmpty(fits, outHeader, outRow, OUT_EXTNAME)) {
+                psErrorStackPrint(stderr, "Unable to write table.");
+                psFree(outHeader);
+                psFree(output);
+                psFree(outRow);
+                psFree(data);
+                exit(PS_EXIT_SYS_ERROR);
+            }
+            psFree(outRow);
+        } else if (!psFitsWriteTable(fits, outHeader, output, OUT_EXTNAME)) {
+            psErrorStackPrint(stderr, "Unable to write table.");
+            psFree(outHeader);
+            psFree(output);
+            psFree(data);
+            exit(PS_EXIT_SYS_ERROR);
+        }
+        psFitsClose(fits);
+    }
+
+    psFree(outHeader);
+    psFree(output);
+    psFree(data);
+
+#endif
+

tags/ipp-20120802/ppTranslate/src/ppMops.h

@@ -30 +30 @@
 } ppMopsArguments;
 
+#if 0
+#warning "IS THERE ANYTHING TO BE MODIFIED HERE?"
+TTYPE19 = 'PSF_CHISQ'          / label for field  19
+TFORM19 = '1E      '           / data format of field: 4-byte REAL
+TTYPE20 = 'CR_NSIGMA'          / label for field  20
+TFORM20 = '1E      '           / data format of field: 4-byte REAL
+TTYPE21 = 'EXT_NSIGMA'         / label for field  21
+TFORM21 = '1E      '           / data format of field: 4-byte REAL
+TTYPE22 = 'PSF_MAJOR'          / label for field  22
+TFORM22 = '1E      '           / data format of field: 4-byte REAL
+TTYPE23 = 'PSF_MINOR'          / label for field  23
+TFORM23 = '1E      '           / data format of field: 4-byte REAL
+TTYPE24 = 'PSF_THETA'          / label for field  24
+TFORM24 = '1E      '           / data format of field: 4-byte REAL
+TTYPE25 = 'PSF_QF  '           / label for field  25
+TFORM25 = '1E      '           / data format of field: 4-byte REAL
+TTYPE26 = 'PSF_NDOF'           / label for field  26
+TFORM26 = '1J      '           / data format of field: 4-byte INTEGER
+TTYPE27 = 'PSF_NPIX'           / label for field  27
+TFORM27 = '1J      '           / data format of field: 4-byte INTEGER
+TTYPE28 = 'MOMENTS_XX'         / label for field  28
+TFORM28 = '1E      '           / data format of field: 4-byte REAL
+TTYPE29 = 'MOMENTS_XY'         / label for field  29
+TFORM29 = '1E      '           / data format of field: 4-byte REAL
+TTYPE30 = 'MOMENTS_YY'         / label for field  30
+TFORM30 = '1E      '           / data format of field: 4-byte REAL
+#endif
+
 /// Parse arguments
 ppMopsArguments *ppMopsArgumentsParse(int argc, char *argv[]);
@@ -47 +75 @@
   long numGood;                       // Number of "good" detections
   psS64 diffSkyfileId;                // unique id for input skyfile
-  psMetadata *table;                  // Columns from the input file (SkyChip.psf extension)
-  psMetadata *fittedTrails;           // Columns from the fitted trails (SkyChip.xfit extension)
-  long fittedTrailsSize;       // Size of the fitted trails (SkyChip.xfit extension)
+  psMetadata *table;                  // Columns from the input file
   psVector *x, *y;                    // Image coordinates
   psVector *ra, *dec;                 // Sky coordinates
@@ -77 +103 @@
 /// @returns 1 if EXTTYPE of "SkyChip.psf" is PS1_DV1
 /// @returns 2 if EXTTYPE of "SkyChip.psf" is PS1_DV2
-/// @returns 3 if EXTTYPE of "SkyChip.psf" is PS1_DV3
 /// @returns 0 otherwise
 int ppMopsGetSkyChipPsfVersion(const psFits* fits);

tags/ipp-20120802/ppTranslate/src/ppMopsArguments.c

@@ -34 +34 @@
     fprintf(stderr, "\n");
     psArgumentHelp(arguments);
-    fprintf(stderr, "\t\tCMF file version can be set to either 1 for PS1_DV1, 2 for PS1_DV2, or 3 for PS1_DV3\n");
+    fprintf(stderr, "\t\tCMF file version can be set to either 1 for PS1_DV1 or 2 for PS1_DV1\n");
     fprintf(stderr, "\t\tSee IPP-MOPS ICD for details\n");
     psLibFinalize();

tags/ipp-20120802/ppTranslate/src/ppMopsGetSkyChipPsfVersion.c

@@ -14 +14 @@
     psFree(headerSkyChip);
     return 2;
-  } else if (strcmp(version, "PS1_DV3") == 0) {
-    psFree(headerSkyChip);
-    return 2;
   }
   psWarning("Unsupported EXTTYPE in SkyChip.psf table: [%s]", version);

tags/ipp-20120802/ppTranslate/src/ppMopsMerge.c

@@ -71 +71 @@
 
     for (int i = 0; i < numInputs; i++) {
-      psTrace("ppMops.merge", 1, "Processing batch %d\n", i);
         ppMopsDetections *det = detections->data[i]; // Detections of interest
         if (!det) {

tags/ipp-20120802/ppTranslate/src/ppMopsRead.c

@@ -8 +8 @@
 #include "ppMops.h"
 
-static psMetadata* sort_by_increasing_idet(psMetadata *fittedTrail, long length, long fillLength);
-
 /*
   ppMopsRead possibly modifies the args->version if the user did not
   set it explicitely.
-*/
+ */
 psArray *ppMopsRead(ppMopsArguments *args)
 {
-  psTrace("ppMops.read", 1, "Reading input detections\n");
-
-  psArray *inNames = args->input;          // Input names
-  long num = inNames->n;                   // Number of inputs
-  psArray *detections = psArrayAlloc(num); // Array of detections, to return
-  for (int i = 0; i < num; i++) {
-    const char *name = inNames->data[i];
-
-    psFits *fits = psFitsOpen(name,  "r"); // FITS file
-    if (!fits) {
-      psError(PS_ERR_IO, false, "Unable to open input %d", i);
-      return false;
+    psTrace("ppMops.read", 1, "Reading input detections\n");
+
+    psArray *inNames = args->input;          // Input names
+    long num = inNames->n;                   // Number of inputs
+    psArray *detections = psArrayAlloc(num); // Array of detections, to return
+    for (int i = 0; i < num; i++) {
+        const char *name = inNames->data[i];
+
+        psFits *fits = psFitsOpen(name,  "r"); // FITS file
+        if (!fits) {
+            psError(PS_ERR_IO, false, "Unable to open input %d", i);
+            return false;
+        }
+
+        psMetadata *header = psFitsReadHeader(NULL, fits); // Primary header
+        if (!header) {
+            psError(PS_ERR_IO, false, "Unable to read header %d", i);
+            return false;
+        }
+
+        psS64 diffSkyfileId = psMetadataLookupS64(NULL, header, "IMAGEID"); // Identifier for image
+        if (diffSkyfileId == 0) {
+            psError(PS_ERR_UNEXPECTED_NULL, false, "Unable to find identifier for image %d", i);
+            return false;
+        }
+
+        if (!psFitsMoveExtName(fits, "SkyChip.psf")) {
+            psError(PS_ERR_IO, false, "Unable to move to HDU with detections");
+            return false;
+        }
+        int skyChipPsfVersion = ppMopsGetSkyChipPsfVersion(fits);
+        if (args->version == 0) {
+          psTrace("ppMops.read", 1, "Changing args->version to %d\n", skyChipPsfVersion);
+          args->version = (unsigned short) skyChipPsfVersion;
+        }
+        if (skyChipPsfVersion == 0) {
+          // Try to read with the user specified version?
+          skyChipPsfVersion = args->version;
+        }
+        /* Display a warning message if there are version
+           inconsistencies between the file and the flag (note that
+           those inconsistencies might be wanted) */
+        if (skyChipPsfVersion != args->version) {
+          if (skyChipPsfVersion > args->version) {
+            psWarning("The FITS data will be downgraded from PS1_DV%d to PS1_DV%d\n",
+                      skyChipPsfVersion, args->version);
+          } else { // Necessarily: skyChipPsfVersion > args->version
+            psWarning("The FITS data will be upgraded from PS1_DV%d to PS1_DV%d (new values set to default 0, NaN...)\n",
+                      skyChipPsfVersion, args->version);            
+          }
+        }
+
+        long size = psFitsTableSize(fits); // Size of table
+        if (size <= 0) {
+            psErrorStackPrint(stderr, "Unable to determine size of table %d", i);
+            psErrorClear();
+            psWarning("Ignoring input %d", i);
+            psFree(header);
+            psFitsClose(fits);
+            continue;
+        }
+        ppMopsDetections *det = ppMopsDetectionsAlloc(size);
+        detections->data[i] = det;
+        det->component = psStringNCopy(name, strrchr(name, '.') - name); // Strip off extension
+        det->num = size;
+        det->diffSkyfileId = diffSkyfileId;
+
+        psTrace("ppMops.read", 3, "Reading %ld rows from %s\n", size, (const char*)inNames->data[i]);
+
+        det->raBoresight = psMemIncrRefCounter(psMetadataLookupStr(NULL, header, "FPA.RA"));
+        det->decBoresight = psMemIncrRefCounter(psMetadataLookupStr(NULL, header, "FPA.DEC"));
+        det->filter = psMemIncrRefCounter(psMetadataLookupStr(NULL, header, "FPA.FILTER"));
+        det->airmass = psMetadataLookupF32(NULL, header, "AIRMASS");
+        det->exptime = psMetadataLookupF32(NULL, header, "EXPTIME");
+        det->posangle = psMetadataLookupF64(NULL, header, "FPA.POSANGLE");
+        det->alt = psMetadataLookupF64(NULL, header, "FPA.ALT");
+        det->az = psMetadataLookupF64(NULL, header, "FPA.AZ");
+        det->mjd = psMetadataLookupF64(NULL, header, "MJD-OBS") + det->exptime / 2.0 / 3600 / 24;
+
+        det->seeing = (float) 0.5 * (psMetadataLookupF32(NULL, header, "FWHM_MAJ") +
+                                     psMetadataLookupF32(NULL, header, "FWHM_MIN"));
+
+        det->naxis1 = psMetadataLookupS32(NULL, header, "IMNAXIS1"); // Number of columns
+        det->naxis2 = psMetadataLookupS32(NULL, header, "IMNAXIS2"); // Number of rows
+
+        psFree(header);
+
+        psMetadata *table = psFitsReadTableAllColumns(fits); // Table of interest
+        if (!table) {
+            psError(PS_ERR_IO, false, "Unable to read table %d", i);
+            return NULL;
+        }
+        det->table = table;
+        psFitsClose(fits);
+        if (args->version == 0) {
+          if (skyChipPsfVersion < 2) {
+            // XXX: TODO: Do we need to add dummy vectors for the missing columns?
+           }
+        }
+
+        psVector *ra = psMetadataLookupVector(NULL, table, "RA_PSF");
+        psVector *dec = psMetadataLookupVector(NULL, table, "DEC_PSF");
+
+        det->raErr = psVectorAlloc(size, PS_TYPE_F64);
+        det->decErr = psVectorAlloc(size, PS_TYPE_F64);
+        det->mask = psVectorAlloc(size, PS_TYPE_U8);
+
+        // convert ra and dec to radians for use in the purge duplicates function
+        det->ra = (psVector*)psBinaryOp(NULL, ra, "*", psScalarAlloc(DEG_TO_RAD(1.0), PS_TYPE_F64));
+        det->dec = (psVector*)psBinaryOp(NULL, dec, "*", psScalarAlloc(DEG_TO_RAD(1.0), PS_TYPE_F64));
+        det->x = psMemIncrRefCounter(psMetadataLookupVector(NULL, table, "X_PSF"));
+        det->y = psMemIncrRefCounter(psMetadataLookupVector(NULL, table, "Y_PSF"));
+        if (!det->ra || !det->dec || !det->x || !det->y) {
+            psError(PS_ERR_UNEXPECTED_NULL, true, "Unable to find all of RA, Dec, X and Y columns");
+            return NULL;
+        }
+
+        // Add our new vectors to the table so that duplicates and masked items may be purged
+        psMetadataAddVector(table, PS_LIST_HEAD, "DEC_ERR", 0, NULL, det->decErr);
+        psMetadataAddVector(table, PS_LIST_HEAD, "RA_ERR", 0, NULL, det->raErr);
+
+        psTrace("ppMops.read", 2, "Read %ld rows from %s\n", det->num, det->component);
+
+        psVector *mag    = psMetadataLookupVector(NULL, table, "PSF_INST_MAG");
+        psVector *magErr = psMetadataLookupVector(NULL, table, "PSF_INST_MAG_SIG");
+        psVector *xErrV  = psMetadataLookupVector(NULL, table, "X_PSF_SIG");
+        psVector *yErrV  = psMetadataLookupVector(NULL, table, "Y_PSF_SIG");
+        psVector *scaleV = psMetadataLookupVector(NULL, table, "PLTSCALE");
+        psVector *angleV = psMetadataLookupVector(NULL, table, "POSANGLE");
+        psVector *flagsV = psMetadataLookupVector(NULL, table, "FLAGS");
+
+        double plateScale = 0.0;        // Plate scale
+        long numGood = 0;               // Number of good rows
+        for (long row = 0; row < size; row++) {
+
+            psU32 flags = flagsV->data.U32[row]; // psFitsTableGetU32(NULL, table, row, "FLAGS");
+            if (flags & SOURCE_MASK) {
+                psTrace("ppMops.read", 10, "Discarding row %ld from input %d because of flags: %ud", row, i, flags);
+                det->mask->data.U8[row] = 0xFF;
+                continue;
+            }
+
+            // Calculate error in RA, Dec
+
+            double xErr = xErrV->data.F32[row];
+            double yErr = yErrV->data.F32[row];
+            double scale = scaleV->data.F32[row];
+            double angle = angleV->data.F32[row];
+
+            if (!isfinite(det->x->data.F32[row]) || !isfinite(det->y->data.F32[row]) ||
+                !isfinite(det->ra->data.F64[row]) || !isfinite(det->dec->data.F64[row]) ||
+                !isfinite(mag->data.F32[row]) || !isfinite(magErr->data.F32[row]) ||
+                !isfinite(xErr) || !isfinite(yErr) || !isfinite(scale) || !isfinite(angle)) {
+                psTrace("ppMops.read", 10,
+                        "Discarding row %ld from input %d because of non-finite values: "
+                        "%f %f %lf %lf %f %f %f %f %f %f",
+                        row, i,
+                        det->x->data.F32[row], det->y->data.F32[row],
+                        det->ra->data.F64[row], det->dec->data.F64[row],
+                        mag->data.F32[row], magErr->data.F32[row],
+                        xErr, yErr, scale, angle);
+                det->mask->data.U8[row] = 0xFF;
+                continue;
+            }
+
+            // XXX Not at all sure I've got the angles around the right way here...
+            double cosAngle = cos(angle), sinAngle = sin(angle);
+            double cosAngle2 = PS_SQR(cosAngle), sinAngle2 = PS_SQR(sinAngle);
+            double xErr2 = PS_SQR(xErr), yErr2 = PS_SQR(yErr);
+            double errScale = scale / 3600.0;
+            det->raErr->data.F64[row] = errScale * sqrt(cosAngle2 * xErr2 + sinAngle2 * yErr2);
+            det->decErr->data.F64[row] = errScale * sqrt(sinAngle2 * xErr2 + cosAngle2 * yErr2);
+
+            det->mask->data.U8[row] = 0;
+            plateScale += scale;
+            numGood++;
+        }
+        det->seeing *= ((float) plateScale) / ((float) numGood);
+
+        // Are we using numGood for anything outside of this function?
+        det->numGood = numGood;
+
+        if (isfinite(args->zp) && numGood > 0) {
+            psBinaryOp(mag, mag, "+", psScalarAlloc(args->zp, PS_TYPE_F32));
+        }
+
+        psTrace("ppMops.read", 2, "Read %ld good rows from %s\n", numGood, (const char*)name);
     }
 
-    psMetadata *header = psFitsReadHeader(NULL, fits); // Primary header
-    if (!header) {
-      psError(PS_ERR_IO, false, "Unable to read header %d", i);
-      return false;
-    }
-
-    psS64 diffSkyfileId = psMetadataLookupS64(NULL, header, "IMAGEID"); // Identifier for image
-    if (diffSkyfileId == 0) {
-      psError(PS_ERR_UNEXPECTED_NULL, false, "Unable to find identifier for image %d", i);
-      return false;
-    }
-
-    if (!psFitsMoveExtName(fits, "SkyChip.psf")) {
-      psError(PS_ERR_IO, false, "Unable to move to HDU with detections");
-      return false;
-    }
-    int skyChipPsfVersion = ppMopsGetSkyChipPsfVersion(fits);
-    if (args->version == 0) {
-      psTrace("ppMops.read", 1, "Changing args->version to %d\n", skyChipPsfVersion);
-      args->version = (unsigned short) skyChipPsfVersion;
-    }
-    if (skyChipPsfVersion == 0) {
-      // Try to read with the user specified version?
-      skyChipPsfVersion = args->version;
-    }
-    /* Display a warning message if there are version
-       inconsistencies between the file and the flag (note that
-       those inconsistencies might be wanted) */
-    if (skyChipPsfVersion != args->version) {
-      if (skyChipPsfVersion > args->version) {
-        psWarning("The FITS data will be downgraded from PS1_DV%d to PS1_DV%d\n",
-                  skyChipPsfVersion, args->version);
-      } else { // Necessarily: skyChipPsfVersion > args->version
-        psWarning("The FITS data will be upgraded from PS1_DV%d to PS1_DV%d (new values set to default 0, NaN...)\n",
-                  skyChipPsfVersion, args->version);
-      }
-    }
-
-    long size = psFitsTableSize(fits); // Size of table
-    if (size <= 0) {
-      psErrorStackPrint(stderr, "Unable to determine size of table %d", i);
-      psErrorClear();
-      psWarning("Ignoring input %d", i);
-      psFree(header);
-      psFitsClose(fits);
-      continue;
-    }
-    ppMopsDetections *det = ppMopsDetectionsAlloc(size);
-    detections->data[i] = det;
-    det->component = psStringNCopy(name, strrchr(name, '.') - name); // Strip off extension
-    det->num = size;
-    det->diffSkyfileId = diffSkyfileId;
-
-    psTrace("ppMops.read", 3, "Reading %ld rows from %s\n", size, (const char*)inNames->data[i]);
-
-    det->raBoresight = psMemIncrRefCounter(psMetadataLookupStr(NULL, header, "FPA.RA"));
-    det->decBoresight = psMemIncrRefCounter(psMetadataLookupStr(NULL, header, "FPA.DEC"));
-    det->filter = psMemIncrRefCounter(psMetadataLookupStr(NULL, header, "FPA.FILTER"));
-    det->airmass = psMetadataLookupF32(NULL, header, "AIRMASS");
-    det->exptime = psMetadataLookupF32(NULL, header, "EXPTIME");
-    det->posangle = psMetadataLookupF64(NULL, header, "FPA.POSANGLE");
-    det->alt = psMetadataLookupF64(NULL, header, "FPA.ALT");
-    det->az = psMetadataLookupF64(NULL, header, "FPA.AZ");
-    det->mjd = psMetadataLookupF64(NULL, header, "MJD-OBS") + det->exptime / 2.0 / 3600 / 24;
-
-    det->seeing = (float) 0.5 * (psMetadataLookupF32(NULL, header, "FWHM_MAJ") +
-                                 psMetadataLookupF32(NULL, header, "FWHM_MIN"));
-
-    det->naxis1 = psMetadataLookupS32(NULL, header, "IMNAXIS1"); // Number of columns
-    det->naxis2 = psMetadataLookupS32(NULL, header, "IMNAXIS2"); // Number of rows
-
-    psFree(header);
-
-    psMetadata *table = psFitsReadTableAllColumns(fits); // Table of interest
-    if (!table) {
-      psError(PS_ERR_IO, false, "Unable to read table %d", i);
-      return NULL;
-    }
-    //Sort table by idet if necessary
-    psMetadata* sortedTable = sort_by_increasing_idet(table, det->num, det->num);
-    if (sortedTable != table) {
-      psFree(table);
-      table = sortedTable;
-    } // else det->table is sorted
-    det->table = table;
-    psTrace("ppMops.read", 9, "PSF table:\n%s\n", psMetadataConfigFormat(det->table));
-    psTrace("ppMops.read", 9, "End of PSF table\n");
-
-    //fittedTrail extension (Supported for releases 2 and above)
-    if (skyChipPsfVersion >= 2) {
-      if (!psFitsMoveExtName(fits, "SkyChip.xfit")) {
-        psTrace("ppMops.read", 3, "No fitted trails extension");
-      } else {
-        psTrace("ppMops.read", 3, "Fitted trails extension found\n");
-        det->fittedTrailsSize = psFitsTableSize(fits);
-        if (det->fittedTrailsSize <= 0) {
-          psErrorStackPrint(stderr, "Unable to determine size of fitted trails extension table %d", i);
-          psTrace("ppMops.read", 3, "No entry in fitted trails extension!!!!\n");
-          psErrorClear();
-          det->fittedTrails = NULL;
-        } else {
-          det->fittedTrails = psFitsReadTableAllColumns(fits); // Table of interest
-          if (!det->fittedTrails) {
-            psError(PS_ERR_IO, false, "Unable to read fittedTrails table %d", i);
-            return NULL;
-          }
-          psMetadata* sortedFittedTrails = sort_by_increasing_idet(det->fittedTrails, det->fittedTrailsSize, det->num);
-          if (sortedFittedTrails != det->fittedTrails) {
-            psFree(det->fittedTrails);
-            det->fittedTrails = sortedFittedTrails;
-          } // else det->fittedTrails is sorted
-          psTrace("ppMops.read", 9, "Fitted trail:\n%s\n", psMetadataConfigFormat(det->fittedTrails));
-          psTrace("ppMops.read", 9, "End of Fitted trail\n");
-        }
-      }
-    }
-
-    psFitsClose(fits);
-
-    if (args->version == 0) {
-      if (skyChipPsfVersion < 2) {
-        // XXX: TODO: Do we need to add dummy vectors for the missing columns?
-      }
-    }
-
-    psVector *ra = psMetadataLookupVector(NULL, table, "RA_PSF");
-    psVector *dec = psMetadataLookupVector(NULL, table, "DEC_PSF");
-
-    det->raErr = psVectorAlloc(size, PS_TYPE_F64);
-    det->decErr = psVectorAlloc(size, PS_TYPE_F64);
-    det->mask = psVectorAlloc(size, PS_TYPE_U8);
-
-    // convert ra and dec to radians for use in the purge duplicates function
-    det->ra = (psVector*)psBinaryOp(NULL, ra, "*", psScalarAlloc(DEG_TO_RAD(1.0), PS_TYPE_F64));
-    det->dec = (psVector*)psBinaryOp(NULL, dec, "*", psScalarAlloc(DEG_TO_RAD(1.0), PS_TYPE_F64));
-    det->x = psMemIncrRefCounter(psMetadataLookupVector(NULL, table, "X_PSF"));
-    det->y = psMemIncrRefCounter(psMetadataLookupVector(NULL, table, "Y_PSF"));
-    if (!det->ra || !det->dec || !det->x || !det->y) {
-      psError(PS_ERR_UNEXPECTED_NULL, true, "Unable to find all of RA, Dec, X and Y columns");
-      return NULL;
-    }
-
-    // Add our new vectors to the table so that duplicates and masked items may be purged
-    psMetadataAddVector(table, PS_LIST_HEAD, "DEC_ERR", 0, NULL, det->decErr);
-    psMetadataAddVector(table, PS_LIST_HEAD, "RA_ERR", 0, NULL, det->raErr);
-
-    psTrace("ppMops.read", 2, "Read %ld rows from %s\n", det->num, det->component);
-
-    psVector *mag    = psMetadataLookupVector(NULL, table, "PSF_INST_MAG");
-    psVector *magErr = psMetadataLookupVector(NULL, table, "PSF_INST_MAG_SIG");
-    psVector *xErrV  = psMetadataLookupVector(NULL, table, "X_PSF_SIG");
-    psVector *yErrV  = psMetadataLookupVector(NULL, table, "Y_PSF_SIG");
-    psVector *scaleV = psMetadataLookupVector(NULL, table, "PLTSCALE");
-    psVector *angleV = psMetadataLookupVector(NULL, table, "POSANGLE");
-    psVector *flagsV = psMetadataLookupVector(NULL, table, "FLAGS");
-
-    double plateScale = 0.0;        // Plate scale
-    long numGood = 0;               // Number of good rows
-    for (long row = 0; row < size; row++) {
-
-      psU32 flags = flagsV->data.U32[row]; // psFitsTableGetU32(NULL, table, row, "FLAGS");
-      if (flags & SOURCE_MASK) {
-        psTrace("ppMops.read", 10, "Discarding row %ld from input %d because of flags: %ud", row, i, flags);
-        det->mask->data.U8[row] = 0xFF;
-        continue;
-      }
-
-      // Calculate error in RA, Dec
-
-      double xErr = xErrV->data.F32[row];
-      double yErr = yErrV->data.F32[row];
-      double scale = scaleV->data.F32[row];
-      double angle = angleV->data.F32[row];
-
-      if (!isfinite(det->x->data.F32[row]) || !isfinite(det->y->data.F32[row]) ||
-          !isfinite(det->ra->data.F64[row]) || !isfinite(det->dec->data.F64[row]) ||
-          !isfinite(mag->data.F32[row]) || !isfinite(magErr->data.F32[row]) ||
-          !isfinite(xErr) || !isfinite(yErr) || !isfinite(scale) || !isfinite(angle)) {
-        psTrace("ppMops.read", 10,
-                "Discarding row %ld from input %d because of non-finite values: "
-                "%f %f %lf %lf %f %f %f %f %f %f",
-                row, i,
-                det->x->data.F32[row], det->y->data.F32[row],
-                det->ra->data.F64[row], det->dec->data.F64[row],
-                mag->data.F32[row], magErr->data.F32[row],
-                xErr, yErr, scale, angle);
-        det->mask->data.U8[row] = 0xFF;
-        continue;
-      }
-
-      // XXX Not at all sure I've got the angles around the right way here...
-      double cosAngle = cos(angle), sinAngle = sin(angle);
-      double cosAngle2 = PS_SQR(cosAngle), sinAngle2 = PS_SQR(sinAngle);
-      double xErr2 = PS_SQR(xErr), yErr2 = PS_SQR(yErr);
-      double errScale = scale / 3600.0;
-      det->raErr->data.F64[row] = errScale * sqrt(cosAngle2 * xErr2 + sinAngle2 * yErr2);
-      det->decErr->data.F64[row] = errScale * sqrt(sinAngle2 * xErr2 + cosAngle2 * yErr2);
-      det->mask->data.U8[row] = 0;
-      plateScale += scale;
-      numGood++;
-    }
-    det->seeing *= ((float) plateScale) / ((float) numGood);
-
-    // Are we using numGood for anything outside of this function?
-    det->numGood = numGood;
-
-    if (isfinite(args->zp) && numGood > 0) {
-      psBinaryOp(mag, mag, "+", psScalarAlloc(args->zp, PS_TYPE_F32));
-    }
-
-    psTrace("ppMops.read", 2, "Read %ld good rows from %s\n", numGood, (const char*)name);
-  }
-
-  psTrace("ppMops.read", 1, "Done reading input detections\n");
-
-  return detections;
+    psTrace("ppMops.read", 1, "Done reading input detections\n");
+
+    return detections;
 }
-
-static long* sorted_idets_indexes(psVector *idets);
-static int compare (const void *a, const void *b);
-static psVector* fill_F32_vector(psVector* in, int size, long* indexes, int indexes_length);
-
-// We want the entries in table/fittedTrails to be sorted by increasing idet
-static psMetadata* sort_by_increasing_idet(psMetadata *aTable, long length, long fillLength) {
-  if (aTable == NULL) {
-    return aTable;
-  }
-  psVector *idets = psMetadataLookupVector(NULL, aTable, "IPP_IDET");
-  if (!idets) {
-    psError(PS_ERR_PROGRAMMING, true, "Failed to find IPP_IDET column");
-    return NULL;
-  }
-  bool is_sorted = true;
-  int count = 0;
-  long previous_value = -1;
-  while (is_sorted && (count<length) ){
-    is_sorted = (previous_value < idets->data.S64[count]);
-    previous_value = idets->data.S64[count];
-    count += 1;
-  }
-  if (is_sorted && (length==fillLength)) {
-    psTrace("ppMops.read", 5, "Values are already sorted by idet and the table is complete\n");
-    return aTable;
-  }
-  //otherwise sort (even if it should not happen) and fill
-  long* sortedIndexVector = sorted_idets_indexes(idets);
-  psMetadata* sortedTable = psMetadataAlloc();
-  char* column_name = "X_EXT";
-  psMetadataAddVector( sortedTable, PS_LIST_TAIL, column_name, 0, NULL, 
-                       fill_F32_vector(psMetadataLookupVector(NULL, aTable, column_name), fillLength, sortedIndexVector, idets->n) );
-  column_name = "Y_EXT";
-  psMetadataAddVector( sortedTable, PS_LIST_TAIL, column_name, 0, NULL, 
-                       fill_F32_vector(psMetadataLookupVector(NULL, aTable, column_name), fillLength, sortedIndexVector, idets->n) );
-  column_name = "X_EXT_SIG";
-  psMetadataAddVector( sortedTable, PS_LIST_TAIL, column_name, 0, NULL, 
-                       fill_F32_vector(psMetadataLookupVector(NULL, aTable, column_name), fillLength, sortedIndexVector, idets->n) );
-  column_name = "Y_EXT_SIG";
-  psMetadataAddVector( sortedTable, PS_LIST_TAIL, column_name, 0, NULL, 
-                       fill_F32_vector(psMetadataLookupVector(NULL, aTable, column_name), fillLength, sortedIndexVector, idets->n) );
-  column_name = "EXT_INST_MAG";
-  psMetadataAddVector( sortedTable, PS_LIST_TAIL, column_name, 0, NULL, 
-                       fill_F32_vector(psMetadataLookupVector(NULL, aTable, column_name), fillLength, sortedIndexVector, idets->n) );
-  column_name = "EXT_INST_MAG_SIG";
-  psMetadataAddVector( sortedTable, PS_LIST_TAIL, column_name, 0, NULL, 
-                       fill_F32_vector(psMetadataLookupVector(NULL, aTable, column_name), fillLength, sortedIndexVector, idets->n) );
-  column_name = "NPARAMS";
-  psMetadataAddVector( sortedTable, PS_LIST_TAIL, column_name, 0, NULL, 
-                       fill_F32_vector(psMetadataLookupVector(NULL, aTable, column_name), fillLength, sortedIndexVector, idets->n) );
-  column_name = "EXT_WIDTH_MAJ";
-  psMetadataAddVector( sortedTable, PS_LIST_TAIL, column_name, 0, NULL, 
-                       fill_F32_vector(psMetadataLookupVector(NULL, aTable, column_name), fillLength, sortedIndexVector, idets->n) );
-  column_name = "EXT_WIDTH_MIN";
-  psMetadataAddVector( sortedTable, PS_LIST_TAIL, column_name, 0, NULL, 
-                       fill_F32_vector(psMetadataLookupVector(NULL, aTable, column_name), fillLength, sortedIndexVector, idets->n) );
-  column_name = "EXT_THETA";
-  psMetadataAddVector( sortedTable, PS_LIST_TAIL, column_name, 0, NULL, 
-                       fill_F32_vector(psMetadataLookupVector(NULL, aTable, column_name), fillLength, sortedIndexVector, idets->n) );
-  column_name = "EXT_WIDTH_MAJ_ERR";
-  psMetadataAddVector( sortedTable, PS_LIST_TAIL, column_name, 0, NULL, 
-                       fill_F32_vector(psMetadataLookupVector(NULL, aTable, column_name), fillLength, sortedIndexVector, idets->n) );
-  column_name = "EXT_WIDTH_MIN_ERR";
-  psMetadataAddVector( sortedTable, PS_LIST_TAIL, column_name, 0, NULL, 
-                       fill_F32_vector(psMetadataLookupVector(NULL, aTable, column_name), fillLength, sortedIndexVector, idets->n) );
-  column_name = "EXT_THETA_ERR";
-  psMetadataAddVector( sortedTable, PS_LIST_TAIL, column_name, 0, NULL, 
-                       fill_F32_vector(psMetadataLookupVector(NULL, aTable, column_name), fillLength, sortedIndexVector, idets->n) );
-  return sortedTable;
-}
-
-static psVector* fill_F32_vector(psVector* in, int length, long* indexes, int indexes_length) {
-  if (in == NULL) {
-    return NULL;
-  }
-  psVector* out = psVectorAlloc(length, in->type.type);
-  int current_index = 0;
-  for (int i = 0; i < indexes_length; i++) {
-    while ( (current_index < indexes[i]) && (current_index < length) ) {
-      out->data.F32[current_index] = NAN;
-      current_index += 1;
-    }
-    out->data.F32[current_index] = in->data.F32[i];
-    current_index += 1;
-  }
-  while ( current_index < length ) {
-    out->data.F32[current_index] = NAN;
-    current_index += 1;
-  }
-  return out;
-}
-
-static int compare (const void *a, const void *b) {
-  return ( *((long*) a)-*((long*) b));
-}
-
-static long* sorted_idets_indexes(psVector *idets) {
-  long* sortedIndexes = (long*) malloc(idets->n * sizeof(long));
-  for (int i=0;i<idets->n;i++) {
-    sortedIndexes[i] = idets->data.S64[i];
-  }
-  qsort(sortedIndexes, idets->n, sizeof(long), (void *) compare);
-  return sortedIndexes;
-}

tags/ipp-20120802/ppTranslate/src/ppMopsWrite.c

@@ -9 +9 @@
 #include "ppTranslateVersion.h"
 
-static bool addOutputColumn(psMetadata *table, const psArray *detections, long total, int extension, char *outColName, char *inColName, bool convertTo32);
+static bool addOutputColumn(psMetadata *table, const psArray *detections, long total, char *outColName, char *inColName, bool convertTo32);
 static bool addSkyfileIDColumn(psMetadata *table, const psArray *detections, long total, char *colName);
 
 bool ppMopsWrite(const psArray *detections, const ppMopsArguments *args)
 {
-  psFits *fits = psFitsOpen(args->output, "w"); // FITS file
-  if (!fits) {
-    psError(PS_ERR_IO, false, "Unable to open output file.");
-    return false;
-  }
-
-  psMetadata *header = psMetadataAlloc(); // Header to write
-  psString source = ppTranslateSource(), version = ppTranslateVersion();
-  psMetadataAddStr(header, PS_LIST_TAIL, "SWSOURCE", 0, "Software source", source);
-  psMetadataAddStr(header, PS_LIST_TAIL, "SWVERSN", 0, "Software version", version);
-  ppTranslateVersionHeader(header);
-  psFree(source);
-  psFree(version);
-
-  psMetadataAddStr(header, PS_LIST_TAIL, "EXP_NAME", 0, "Exposure name", args->exp_name);
-  psMetadataAddS64(header, PS_LIST_TAIL, "EXP_ID", 0, "Exposure identifier", args->exp_id);
-  psMetadataAddS64(header, PS_LIST_TAIL, "CHIP_ID", 0, "Chip stage identifier", args->chip_id);
-  psMetadataAddS64(header, PS_LIST_TAIL, "CAM_ID", 0, "Cam stage identifier", args->cam_id);
-  psMetadataAddS64(header, PS_LIST_TAIL, "FAKE_ID", 0, "Fake stage identifier", args->fake_id);
-  psMetadataAddS64(header, PS_LIST_TAIL, "WARP_ID", 0, "Warp stage identifier", args->warp_id);
-  psMetadataAddS64(header, PS_LIST_TAIL, "DIFF_ID", 0, "Diff stage identifier", args->diff_id);
-  psMetadataAddBool(header, PS_LIST_TAIL, "DIFF_POS", 0, "Positive subtraction?", args->positive);
-
-  // Get these header words from the first input non null input
-  ppMopsDetections *det = NULL;
-  for (int d = 0; d < detections->n && det == NULL; d++) {
-    det = detections->data[d];
-  }
-  if (det != NULL) {
-    psMetadataAddF64(header, PS_LIST_TAIL, "MJD-OBS", 0, "MJD of exposure midpoint", det->mjd);
-    psMetadataAddStr(header, PS_LIST_TAIL, "RA", 0, "Right Ascension of boresight", det->raBoresight);
-    psMetadataAddStr(header, PS_LIST_TAIL, "DEC", 0, "Declination of boresight", det->decBoresight);
-    psMetadataAddF64(header, PS_LIST_TAIL, "TEL_ALT", 0, "Telescope altitude", det->alt);
-    psMetadataAddF64(header, PS_LIST_TAIL, "TEL_AZ", 0, "Telescope azimuth", det->az);
-    psMetadataAddF64(header, PS_LIST_TAIL, "EXPTIME", 0, "Exposure time (sec)", det->exptime);
-    psMetadataAddF64(header, PS_LIST_TAIL, "ROTANGLE", 0, "Rotator position angle", det->posangle);
-    psMetadataAddStr(header, PS_LIST_TAIL, "FILTER", 0, "Filter name", det->filter);
-    psMetadataAddF32(header, PS_LIST_TAIL, "AIRMASS", 0, "Airmass of exposure", det->airmass);
-    psMetadataAddF32(header, PS_LIST_TAIL, "SEEING", 0, "Mean seeing", det->seeing);
-  } else {
-    psWarning("no inputs with surviving detections. output header will be incomplete");
-  }
-  psMetadataAddStr(header, PS_LIST_TAIL, "OBSCODE", 0, "IAU Observatory code", OBSERVATORY_CODE);
-  psMetadataAddF32(header, PS_LIST_TAIL, "MAGZP", 0, "Magnitude zero point", args->zp);
-  psMetadataAddF32(header, PS_LIST_TAIL, "MAGZPERR", 0, "Error in magnitude zero point", args->zpErr);
-  psMetadataAddF32(header, PS_LIST_TAIL, "ASTRORMS", 0, "RMS of astrometric fit", args->rmsAstrom);
-
-  //field in header that tells about the CMF version
-  char cmfVersion[8];
-  sprintf(cmfVersion, PS1_DV_FORMAT, args->version);
-  psMetadataAddStr(header, PS_LIST_TAIL, "CMFVERSION", 0, "CMF version", cmfVersion);
-
-  // Find the total number of detections
-
-  long total = 0;
-  for (long i=0; i<detections->n; i++) {
-    ppMopsDetections *det = detections->data[i];
-    if (!det) {
-      continue;
-    }
-    total += det->num;
-  }
-
-  psTrace("ppMops.write", 1, "Writing %ld rows to %s", total, args->output);
-
-  if (total == 0) {
-    // Write dummy table
-    psMetadata *row = psMetadataAlloc(); // Output row
-    psMetadataAddF64(row, PS_LIST_TAIL, "RA", 0, "Right ascension (degrees)", NAN);
-    psMetadataAddF64(row, PS_LIST_TAIL, "RA_ERR", 0, "Right ascension error (degrees)", NAN);
-    psMetadataAddF64(row, PS_LIST_TAIL, "DEC", 0, "Declination (degrees)", NAN);
-    psMetadataAddF64(row, PS_LIST_TAIL, "DEC_ERR", 0, "Declination error (degrees)", NAN);
-    psMetadataAddF32(row, PS_LIST_TAIL, "MAG", 0, "Magnitude", NAN);
-    psMetadataAddF32(row, PS_LIST_TAIL, "MAG_ERR", 0, "Magnitude error", NAN);
-    psMetadataAddF32(row, PS_LIST_TAIL, "PSF_CHI2", 0, "chi^2 of PSF fit", NAN);
-    psMetadataAddS32(row, PS_LIST_TAIL, "PSF_DOF", 0, "Degrees of freedom of PSF fit", 0);
-    psMetadataAddF32(row, PS_LIST_TAIL, "CR_SIGNIFICANCE", 0, "Significance of CR", NAN);
-    psMetadataAddF32(row, PS_LIST_TAIL, "EXT_SIGNIFICANCE", 0, "Significance of extendedness", NAN);
-    psMetadataAddF32(row, PS_LIST_TAIL, "PSF_MAJOR", 0, "PSF major axis (pixels)", NAN);
-    psMetadataAddF32(row, PS_LIST_TAIL, "PSF_MINOR", 0, "PSF minor axis (pixels)", NAN);
-    psMetadataAddF32(row, PS_LIST_TAIL, "PSF_THETA", 0, "PSF position angle (deg on chip)", NAN);
-    psMetadataAddF32(row, PS_LIST_TAIL, "PSF_QUALITY", 0, "PSF quality factor", NAN);
-    psMetadataAddS32(row, PS_LIST_TAIL, "PSF_NPIX", 0, "Number of pixels in PSF", 0);
-    psMetadataAddF32(row, PS_LIST_TAIL, "MOMENTS_XX", 0, "xx moment", NAN);
-    psMetadataAddF32(row, PS_LIST_TAIL, "MOMENTS_XY", 0, "xy moment", NAN);
-    psMetadataAddF32(row, PS_LIST_TAIL, "MOMENTS_YY", 0, "yy moment", NAN);
-    psMetadataAddU32(row, PS_LIST_TAIL, "FLAGS", 0, "Detection bit flags", 0);
-    psMetadataAddS64(row, PS_LIST_TAIL, "DIFF_SKYFILE_ID", 0, "Identifier for diff skyfile", 0);
-
-    psMetadataAddS32(row, PS_LIST_TAIL, "N_POS", 0, "Number of positive pixels", 0);
-    psMetadataAddF32(row, PS_LIST_TAIL, "F_POS", 0, "Fraction of positive pixels", NAN);
-    psMetadataAddF32(row, PS_LIST_TAIL, "RATIO_BAD", 0, "Ratio of positive pixels to negative", NAN);
-    psMetadataAddF32(row, PS_LIST_TAIL, "RATIO_MASK", 0, "Ratio of positive pixels to masked", NAN);
-    psMetadataAddF32(row, PS_LIST_TAIL, "RATIO_ALL", 0, "Ratio of positive pixels to all", NAN);
-
-    if (args->version == 2) {
-      // Write data of version 2 (see ICD)
-      psMetadataAdd (row, PS_LIST_TAIL, "IPP_IDET",         PS_DATA_U32, "IPP detection identifier index",
-                     NAN);
-      psMetadataAdd (row, PS_LIST_TAIL, "PSF_INST_FLUX",    PS_DATA_F32, "PSF fit instrumental magnitude",
-                     NAN);
-      psMetadataAdd (row, PS_LIST_TAIL, "PSF_INST_FLUX_SIG",PS_DATA_F32, "Sigma of PSF instrumental magnitude",
-                     NAN);
-      psMetadataAdd (row, PS_LIST_TAIL, "AP_MAG",           PS_DATA_F32, "magnitude in standard aperture",
-                     NAN);
-      psMetadataAdd (row, PS_LIST_TAIL, "AP_MAG_RAW",       PS_DATA_F32, "magnitude in real aperture",
-                     NAN);
-      psMetadataAdd (row, PS_LIST_TAIL, "AP_MAG_RADIUS",    PS_DATA_F32, "radius used for aperture mags",
-                     NAN);
-      psMetadataAdd (row, PS_LIST_TAIL, "AP_FLUX",          PS_DATA_F32, "instrumental flux in standard aperture",
-                     NAN);
-      psMetadataAdd (row, PS_LIST_TAIL, "AP_FLUX_SIG",      PS_DATA_F32, "aperture flux error",
-                     NAN);
-      psMetadataAdd (row, PS_LIST_TAIL, "PEAK_FLUX_AS_MAG", PS_DATA_F32, "Peak flux expressed as magnitude",
-                     NAN);
-      psMetadataAdd (row, PS_LIST_TAIL, "CAL_PSF_MAG",      PS_DATA_F32, "PSF Magnitude using supplied calibration",
-                     NAN);
-      psMetadataAdd (row, PS_LIST_TAIL, "CAL_PSF_MAG_SIG",  PS_DATA_F32, "measured scatter of zero point calibration",
-                     NAN);
-      psMetadataAdd (row, PS_LIST_TAIL, "SKY",              PS_DATA_F32, "Sky level",
-                     NAN);
-      psMetadataAdd (row, PS_LIST_TAIL, "SKY_SIGMA",        PS_DATA_F32, "Sigma of sky level",
-                     NAN);
-      psMetadataAdd (row, PS_LIST_TAIL, "PSF_QF_PERFECT",   PS_DATA_F32, "PSF coverage/quality factor (poor)",
-                     NAN);
-      psMetadataAdd (row, PS_LIST_TAIL, "MOMENTS_R1",       PS_DATA_F32, "first radial moment",
-                     NAN);
-      psMetadataAdd (row, PS_LIST_TAIL, "MOMENTS_RH",       PS_DATA_F32, "half radial moment",
-                     NAN);
-      psMetadataAdd (row, PS_LIST_TAIL, "KRON_FLUX",        PS_DATA_F32, "Kron Flux (in 2.5 R1)",
-                     NAN);
-      psMetadataAdd (row, PS_LIST_TAIL, "KRON_FLUX_ERR",    PS_DATA_F32, "Kron Flux Error",
-                     NAN);
-      psMetadataAdd (row, PS_LIST_TAIL, "KRON_FLUX_INNER",  PS_DATA_F32, "Kron Flux (in 1.0 R1)", 
-                     NAN);
-      psMetadataAdd (row, PS_LIST_TAIL, "KRON_FLUX_OUTER",  PS_DATA_F32, "Kron Flux (in 4.0 R1)",
-                     NAN);
-      psMetadataAdd (row, PS_LIST_TAIL, "DIFF_R_P",         PS_DATA_F32, "distance to positive match source",
-                     NAN);
-      psMetadataAdd (row, PS_LIST_TAIL, "DIFF_SN_P",        PS_DATA_F32, "signal-to-noise of pos match src",
-                     NAN);
-      psMetadataAdd (row, PS_LIST_TAIL, "DIFF_R_M",         PS_DATA_F32, "distance to negative match source",
-                     NAN);
-      psMetadataAdd (row, PS_LIST_TAIL, "DIFF_SN_M",        PS_DATA_F32, "signal-to-noise of neg match src",
-                     NAN);
-      psMetadataAdd (row, PS_LIST_TAIL, "FLAGS2",           PS_DATA_U32, "psphot analysis flags (group 2)",
-                     0);
-      psMetadataAdd (row, PS_LIST_TAIL, "N_FRAMES",         PS_DATA_U16, "Number of frames overlapping source center", 
-                     0);
-      psMetadataAdd (row, PS_LIST_TAIL, "PADDING",          PS_DATA_S16, "padding", 
-                     0);
-      if (args->version == 3) {
-        // TODO: Write data of version 3 (see ICD)
-      }
-    }
-
-    if (!psFitsWriteTableEmpty(fits, header, row, OUT_EXTNAME)) {
-      psErrorStackPrint(stderr, "Unable to write empty table.");
-      psFree(header);
-      psFree(row);
-      return false;
-    }
-    psFree(row);
-  } else {
-
-#define addColumn(_outName, _inName, _convertTo32)                      \
-    if (!addOutputColumn(table, detections, total, 0, _outName, _inName, _convertTo32)) { \
-      psError(PS_ERR_UNKNOWN, false, "Failed to add column %s", _outName); \
-      return false;                                                     \
-    }
-
-    // Allocate the output table
-    psMetadata *table = psMetadataAlloc();
-    addColumn("RA", "RA_PSF", 0);
-    addColumn("RA_ERR", NULL, 0);      // calculated from various parameters including X_PSF_SIG Y_PSF_SIG and POSANG
-    addColumn("DEC", "DEC_PSF", 0);
-    addColumn("DEC_ERR", NULL, 0);     // calculated from various parameters including X_PSF_SIG Y_PSF_SIG and POSANG
-    addColumn("MAG", "PSF_INST_MAG", 0);
-    addColumn("MAG_ERR", "PSF_INST_MAG_SIG", 0);
-    addColumn("PSF_CHI2", "PSF_CHISQ", 0);
-    addColumn("PSF_DOF", "PSF_NDOF", 1);
-    addColumn("CR_SIGNIFICANCE", "CR_NSIGMA", 0);
-    addColumn("EXT_SIGNIFICANCE", "EXT_NSIGMA", 0);
-    addColumn("PSF_MAJOR", NULL, 0);
-    addColumn("PSF_MINOR", NULL, 0);     
-    addColumn("PSF_THETA", NULL, 0);    
-    addColumn("PSF_QUALITY", "PSF_QF", 0); 
-    addColumn("PSF_NPIX", NULL, 1);        
-    addColumn("MOMENTS_XX", NULL, 0);
-    addColumn("MOMENTS_XY", NULL, 0);
-    addColumn("MOMENTS_YY", NULL, 0);
-    addColumn("N_POS", "DIFF_NPOS", 1);
-    addColumn("F_POS", "DIFF_FRATIO", 0);
-    addColumn("RATIO_BAD", "DIFF_NRATIO_BAD", 0);
-    addColumn("RATIO_MASK", "DIFF_NRATIO_MASK", 0);
-    addColumn("RATIO_ALL", "DIFF_NRATIO_ALL", 0);
-    addColumn("FLAGS", "FLAGS", 1);
-    addSkyfileIDColumn(table, detections, total, "DIFF_SKYFILE_ID");
-    if (args->version >= 2) {
-      addColumn("IPP_IDET", NULL, 1);
-      addColumn("PSF_INST_FLUX", NULL, 0);
-      addColumn("PSF_INST_FLUX_SIG", NULL, 0);
-      addColumn("AP_MAG", NULL, 0);
-      addColumn("AP_MAG_RAW", NULL, 0);
-      addColumn("AP_MAG_RADIUS", NULL, 0);
-      addColumn("AP_FLUX", NULL, 0);
-      addColumn("AP_FLUX_SIG", NULL, 0);
-      addColumn("PEAK_FLUX_AS_MAG", NULL, 0);
-      addColumn("CAL_PSF_MAG", NULL, 0);
-      addColumn("CAL_PSF_MAG_SIG", NULL, 0);
-      addColumn("SKY", NULL, 0);
-      addColumn("SKY_SIGMA", NULL, 0);
-      addColumn("PSF_QF_PERFECT", NULL, 0);
-      addColumn("MOMENTS_R1", NULL, 0);
-      addColumn("MOMENTS_RH", NULL, 0);
-      addColumn("KRON_FLUX", NULL, 0);
-      addColumn("KRON_FLUX_ERR", NULL, 0);
-      addColumn("KRON_FLUX_INNER", NULL, 0);
-      addColumn("KRON_FLUX_OUTER", NULL, 0);
-      addColumn("DIFF_R_P", NULL, 0);
-      addColumn("DIFF_SN_P", NULL, 0);
-      addColumn("DIFF_R_M", NULL, 0);
-      addColumn("DIFF_SN_M", NULL, 0);
-      addColumn("FLAGS2", NULL, 1);
-      addColumn("IPP_IDET", NULL, 0);
-      addColumn("N_FRAMES", NULL, 0);
-      addColumn("PADDING", NULL, 0);
-      if (det->fittedTrails != NULL) {
-#define addFittedTrailColumn(_outName, _inName, _convertTo32)                   \
-        if (!addOutputColumn(table, detections, total, 1, _outName, _inName, _convertTo32)) { \
-          psTrace("ppMops.write", 1, "Failed to add column %s -> Replaced with NAN", _outName); \
+
+    psFits *fits = psFitsOpen(args->output, "w"); // FITS file
+    if (!fits) {
+        psError(PS_ERR_IO, false, "Unable to open output file.");
+        return false;
+    }
+
+    psMetadata *header = psMetadataAlloc(); // Header to write
+    psString source = ppTranslateSource(), version = ppTranslateVersion();
+    psMetadataAddStr(header, PS_LIST_TAIL, "SWSOURCE", 0, "Software source", source);
+    psMetadataAddStr(header, PS_LIST_TAIL, "SWVERSN", 0, "Software version", version);
+    ppTranslateVersionHeader(header);
+    psFree(source);
+    psFree(version);
+
+    psMetadataAddStr(header, PS_LIST_TAIL, "EXP_NAME", 0, "Exposure name", args->exp_name);
+    psMetadataAddS64(header, PS_LIST_TAIL, "EXP_ID", 0, "Exposure identifier", args->exp_id);
+    psMetadataAddS64(header, PS_LIST_TAIL, "CHIP_ID", 0, "Chip stage identifier", args->chip_id);
+    psMetadataAddS64(header, PS_LIST_TAIL, "CAM_ID", 0, "Cam stage identifier", args->cam_id);
+    psMetadataAddS64(header, PS_LIST_TAIL, "FAKE_ID", 0, "Fake stage identifier", args->fake_id);
+    psMetadataAddS64(header, PS_LIST_TAIL, "WARP_ID", 0, "Warp stage identifier", args->warp_id);
+    psMetadataAddS64(header, PS_LIST_TAIL, "DIFF_ID", 0, "Diff stage identifier", args->diff_id);
+    psMetadataAddBool(header, PS_LIST_TAIL, "DIFF_POS", 0, "Positive subtraction?", args->positive);
+
+    // Get these header words from the first input non null input
+    ppMopsDetections *det = NULL;
+    for (int d = 0; d < detections->n && det == NULL; d++) {
+        det = detections->data[d];
+    }
+    if (det != NULL) {
+        psMetadataAddF64(header, PS_LIST_TAIL, "MJD-OBS", 0, "MJD of exposure midpoint", det->mjd);
+        psMetadataAddStr(header, PS_LIST_TAIL, "RA", 0, "Right Ascension of boresight", det->raBoresight);
+        psMetadataAddStr(header, PS_LIST_TAIL, "DEC", 0, "Declination of boresight", det->decBoresight);
+        psMetadataAddF64(header, PS_LIST_TAIL, "TEL_ALT", 0, "Telescope altitude", det->alt);
+        psMetadataAddF64(header, PS_LIST_TAIL, "TEL_AZ", 0, "Telescope azimuth", det->az);
+        psMetadataAddF64(header, PS_LIST_TAIL, "EXPTIME", 0, "Exposure time (sec)", det->exptime);
+        psMetadataAddF64(header, PS_LIST_TAIL, "ROTANGLE", 0, "Rotator position angle", det->posangle);
+        psMetadataAddStr(header, PS_LIST_TAIL, "FILTER", 0, "Filter name", det->filter);
+        psMetadataAddF32(header, PS_LIST_TAIL, "AIRMASS", 0, "Airmass of exposure", det->airmass);
+        psMetadataAddF32(header, PS_LIST_TAIL, "SEEING", 0, "Mean seeing", det->seeing);
+    } else {
+        psWarning("no inputs with surviving detections. output header will be incomplete");
+    }
+    psMetadataAddStr(header, PS_LIST_TAIL, "OBSCODE", 0, "IAU Observatory code", OBSERVATORY_CODE);
+    psMetadataAddF32(header, PS_LIST_TAIL, "MAGZP", 0, "Magnitude zero point", args->zp);
+    psMetadataAddF32(header, PS_LIST_TAIL, "MAGZPERR", 0, "Error in magnitude zero point", args->zpErr);
+    psMetadataAddF32(header, PS_LIST_TAIL, "ASTRORMS", 0, "RMS of astrometric fit", args->rmsAstrom);
+
+    //field in header that tells about the CMF version
+    char cmfVersion[8];
+    sprintf(cmfVersion, PS1_DV_FORMAT, args->version);
+    psMetadataAddStr(header, PS_LIST_TAIL, "CMFVERSION", 0, "CMF version", cmfVersion);
+
+    // Find the total number of detections
+
+    long total = 0;
+    for (long i=0; i<detections->n; i++) {
+        ppMopsDetections *det = detections->data[i];
+        if (!det) {
+            continue;
         }
-        addFittedTrailColumn("X_EXT", NULL, 0);
-        addFittedTrailColumn("Y_EXT", NULL, 0);
-        addFittedTrailColumn("X_EXT_SIG", NULL, 0);
-        addFittedTrailColumn("Y_EXT_SIG", NULL, 0);
-        addFittedTrailColumn("EXT_INST_MAG", NULL, 0);
-        addFittedTrailColumn("EXT_INST_MAG_SIG", NULL, 0);
-        addFittedTrailColumn("NPARAMS", NULL, 0);
-        addFittedTrailColumn("EXT_WIDTH_MAJ", NULL, 0);
-        addFittedTrailColumn("EXT_WIDTH_MIN", NULL, 0);
-        addFittedTrailColumn("EXT_THETA", NULL, 0);
-        addFittedTrailColumn("EXT_WIDTH_MAJ_ERR", NULL, 0);
-        addFittedTrailColumn("EXT_WIDTH_MIN_ERR", NULL, 0);
-        addFittedTrailColumn("EXT_THETA_ERR", NULL, 0);
-      }
-    }
-    if (!psFitsWriteTableAllColumns(fits, header, table, OUT_EXTNAME)) {
-      psError(psErrorCodeLast(), false, "Unable to write table");
-      return false;
-    }
-    psFree(table);
-  }
-
-  psFree(header);
-  psFitsClose(fits);
-
-  psTrace("ppMops.write", 1, "Done writing %ld rows to %s", det->num, args->output);
-
-  return true;
+        total += det->num;
+    }
+
+    psTrace("ppMops.write", 1, "Writing %ld rows to %s", total, args->output);
+
+    if (total == 0) {
+        // Write dummy table
+        psMetadata *row = psMetadataAlloc(); // Output row
+        psMetadataAddF64(row, PS_LIST_TAIL, "RA", 0, "Right ascension (degrees)", NAN);
+        psMetadataAddF64(row, PS_LIST_TAIL, "RA_ERR", 0, "Right ascension error (degrees)", NAN);
+        psMetadataAddF64(row, PS_LIST_TAIL, "DEC", 0, "Declination (degrees)", NAN);
+        psMetadataAddF64(row, PS_LIST_TAIL, "DEC_ERR", 0, "Declination error (degrees)", NAN);
+        psMetadataAddF32(row, PS_LIST_TAIL, "MAG", 0, "Magnitude", NAN);
+        psMetadataAddF32(row, PS_LIST_TAIL, "MAG_ERR", 0, "Magnitude error", NAN);
+        psMetadataAddF32(row, PS_LIST_TAIL, "PSF_CHI2", 0, "chi^2 of PSF fit", NAN);
+        psMetadataAddS32(row, PS_LIST_TAIL, "PSF_DOF", 0, "Degrees of freedom of PSF fit", 0);
+        psMetadataAddF32(row, PS_LIST_TAIL, "CR_SIGNIFICANCE", 0, "Significance of CR", NAN);
+        psMetadataAddF32(row, PS_LIST_TAIL, "EXT_SIGNIFICANCE", 0, "Significance of extendedness", NAN);
+        psMetadataAddF32(row, PS_LIST_TAIL, "PSF_MAJOR", 0, "PSF major axis (pixels)", NAN);
+        psMetadataAddF32(row, PS_LIST_TAIL, "PSF_MINOR", 0, "PSF minor axis (pixels)", NAN);
+        psMetadataAddF32(row, PS_LIST_TAIL, "PSF_THETA", 0, "PSF position angle (deg on chip)", NAN);
+        psMetadataAddF32(row, PS_LIST_TAIL, "PSF_QUALITY", 0, "PSF quality factor", NAN);
+        psMetadataAddS32(row, PS_LIST_TAIL, "PSF_NPIX", 0, "Number of pixels in PSF", 0);
+        psMetadataAddF32(row, PS_LIST_TAIL, "MOMENTS_XX", 0, "xx moment", NAN);
+        psMetadataAddF32(row, PS_LIST_TAIL, "MOMENTS_XY", 0, "xy moment", NAN);
+        psMetadataAddF32(row, PS_LIST_TAIL, "MOMENTS_YY", 0, "yy moment", NAN);
+        psMetadataAddU32(row, PS_LIST_TAIL, "FLAGS", 0, "Detection bit flags", 0);
+        psMetadataAddS64(row, PS_LIST_TAIL, "DIFF_SKYFILE_ID", 0, "Identifier for diff skyfile", 0);
+
+        psMetadataAddS32(row, PS_LIST_TAIL, "N_POS", 0, "Number of positive pixels", 0);
+        psMetadataAddF32(row, PS_LIST_TAIL, "F_POS", 0, "Fraction of positive pixels", NAN);
+        psMetadataAddF32(row, PS_LIST_TAIL, "RATIO_BAD", 0, "Ratio of positive pixels to negative", NAN);
+        psMetadataAddF32(row, PS_LIST_TAIL, "RATIO_MASK", 0, "Ratio of positive pixels to masked", NAN);
+        psMetadataAddF32(row, PS_LIST_TAIL, "RATIO_ALL", 0, "Ratio of positive pixels to all", NAN);
+
+        if (args->version == 2) {
+              // Write data of version 2 (see ICD)
+              psMetadataAdd (row, PS_LIST_TAIL, "IPP_IDET",         PS_DATA_U32, "IPP detection identifier index",
+                             NAN);
+              psMetadataAdd (row, PS_LIST_TAIL, "PSF_INST_FLUX",    PS_DATA_F32, "PSF fit instrumental magnitude",
+                             NAN);
+              psMetadataAdd (row, PS_LIST_TAIL, "PSF_INST_FLUX_SIG",PS_DATA_F32, "Sigma of PSF instrumental magnitude",
+                             NAN);
+              psMetadataAdd (row, PS_LIST_TAIL, "AP_MAG",           PS_DATA_F32, "magnitude in standard aperture",
+                             NAN);
+              psMetadataAdd (row, PS_LIST_TAIL, "AP_MAG_RAW",       PS_DATA_F32, "magnitude in real aperture",
+                             NAN);
+              psMetadataAdd (row, PS_LIST_TAIL, "AP_MAG_RADIUS",    PS_DATA_F32, "radius used for aperture mags",
+                             NAN);
+              psMetadataAdd (row, PS_LIST_TAIL, "AP_FLUX",          PS_DATA_F32, "instrumental flux in standard aperture",
+                             NAN);
+              psMetadataAdd (row, PS_LIST_TAIL, "AP_FLUX_SIG",      PS_DATA_F32, "aperture flux error",
+                             NAN);
+              psMetadataAdd (row, PS_LIST_TAIL, "PEAK_FLUX_AS_MAG", PS_DATA_F32, "Peak flux expressed as magnitude",
+                             NAN);
+              psMetadataAdd (row, PS_LIST_TAIL, "CAL_PSF_MAG",      PS_DATA_F32, "PSF Magnitude using supplied calibration",
+                             NAN);
+              psMetadataAdd (row, PS_LIST_TAIL, "CAL_PSF_MAG_SIG",  PS_DATA_F32, "measured scatter of zero point calibration",
+                             NAN);
+              psMetadataAdd (row, PS_LIST_TAIL, "SKY",              PS_DATA_F32, "Sky level",
+                             NAN);
+              psMetadataAdd (row, PS_LIST_TAIL, "SKY_SIGMA",        PS_DATA_F32, "Sigma of sky level",
+                             NAN);
+              psMetadataAdd (row, PS_LIST_TAIL, "PSF_QF_PERFECT",   PS_DATA_F32, "PSF coverage/quality factor (poor)",
+                             NAN);
+              psMetadataAdd (row, PS_LIST_TAIL, "MOMENTS_R1",       PS_DATA_F32, "first radial moment",
+                             NAN);
+              psMetadataAdd (row, PS_LIST_TAIL, "MOMENTS_RH",       PS_DATA_F32, "half radial moment",
+                             NAN);
+              psMetadataAdd (row, PS_LIST_TAIL, "KRON_FLUX",        PS_DATA_F32, "Kron Flux (in 2.5 R1)",
+                             NAN);
+              psMetadataAdd (row, PS_LIST_TAIL, "KRON_FLUX_ERR",    PS_DATA_F32, "Kron Flux Error",
+                             NAN);
+              psMetadataAdd (row, PS_LIST_TAIL, "KRON_FLUX_INNER",  PS_DATA_F32, "Kron Flux (in 1.0 R1)", 
+                             NAN);
+              psMetadataAdd (row, PS_LIST_TAIL, "KRON_FLUX_OUTER",  PS_DATA_F32, "Kron Flux (in 4.0 R1)",
+                             NAN);
+              psMetadataAdd (row, PS_LIST_TAIL, "DIFF_R_P",         PS_DATA_F32, "distance to positive match source",
+                             NAN);
+              psMetadataAdd (row, PS_LIST_TAIL, "DIFF_SN_P",        PS_DATA_F32, "signal-to-noise of pos match src",
+                             NAN);
+              psMetadataAdd (row, PS_LIST_TAIL, "DIFF_R_M",         PS_DATA_F32, "distance to negative match source",
+                             NAN);
+              psMetadataAdd (row, PS_LIST_TAIL, "DIFF_SN_M",        PS_DATA_F32, "signal-to-noise of neg match src",
+                             NAN);
+              psMetadataAdd (row, PS_LIST_TAIL, "FLAGS2",           PS_DATA_U32, "psphot analysis flags (group 2)",
+                             0);
+              psMetadataAdd (row, PS_LIST_TAIL, "N_FRAMES",         PS_DATA_U16, "Number of frames overlapping source center", 
+                             0);
+              psMetadataAdd (row, PS_LIST_TAIL, "PADDING",          PS_DATA_S16, "padding", 
+                             0);
+        }
+
+        if (!psFitsWriteTableEmpty(fits, header, row, OUT_EXTNAME)) {
+            psErrorStackPrint(stderr, "Unable to write empty table.");
+            psFree(header);
+            psFree(row);
+            return false;
+        }
+        psFree(row);
+    } else {
+
+#define addColumn(_outName, _inName, _convertTo32) \
+        if (!addOutputColumn(table, detections, total, _outName, _inName, _convertTo32)) { \
+            psError(PS_ERR_UNKNOWN, false, "Failed to add column %s", _outName); \
+            return false; \
+        }
+
+        // Allocate the output table
+        psMetadata *table = psMetadataAlloc();
+        addColumn("RA", "RA_PSF", 0);
+        addColumn("RA_ERR", NULL, 0);      // calculated from various parameters including X_PSF_SIG Y_PSF_SIG and POSANG
+        addColumn("DEC", "DEC_PSF", 0);
+        addColumn("DEC_ERR", NULL, 0);     // calculated from various parameters including X_PSF_SIG Y_PSF_SIG and POSANG
+        addColumn("MAG", "PSF_INST_MAG", 0);
+        addColumn("MAG_ERR", "PSF_INST_MAG_SIG", 0);
+        addColumn("PSF_CHI2", "PSF_CHISQ", 0);
+        addColumn("PSF_DOF", "PSF_NDOF", 1);
+        addColumn("CR_SIGNIFICANCE", "CR_NSIGMA", 0);
+        addColumn("EXT_SIGNIFICANCE", "EXT_NSIGMA", 0);
+        addColumn("PSF_MAJOR", NULL, 0);
+        addColumn("PSF_MINOR", NULL, 0);     
+        addColumn("PSF_THETA", NULL, 0);    
+        addColumn("PSF_QUALITY", "PSF_QF", 0); 
+        addColumn("PSF_NPIX", NULL, 1);        
+        addColumn("MOMENTS_XX", NULL, 0);
+        addColumn("MOMENTS_XY", NULL, 0);
+        addColumn("MOMENTS_YY", NULL, 0);
+        addColumn("N_POS", "DIFF_NPOS", 1);
+        addColumn("F_POS", "DIFF_FRATIO", 0);
+        addColumn("RATIO_BAD", "DIFF_NRATIO_BAD", 0);
+        addColumn("RATIO_MASK", "DIFF_NRATIO_MASK", 0);
+        addColumn("RATIO_ALL", "DIFF_NRATIO_ALL", 0);
+        addColumn("FLAGS", "FLAGS", 1);
+        addSkyfileIDColumn(table, detections, total, "DIFF_SKYFILE_ID");
+        if (args->version == 2) {
+            addColumn("IPP_IDET", NULL, 1);
+            addColumn("PSF_INST_FLUX", NULL, 0);
+            addColumn("PSF_INST_FLUX_SIG", NULL, 0);
+            addColumn("AP_MAG", NULL, 0);
+            addColumn("AP_MAG_RAW", NULL, 0);
+            addColumn("AP_MAG_RADIUS", NULL, 0);
+            addColumn("AP_FLUX", NULL, 0);
+            addColumn("AP_FLUX_SIG", NULL, 0);
+            addColumn("PEAK_FLUX_AS_MAG", NULL, 0);
+            addColumn("CAL_PSF_MAG", NULL, 0);
+            addColumn("CAL_PSF_MAG_SIG", NULL, 0);
+            addColumn("SKY", NULL, 0);
+            addColumn("SKY_SIGMA", NULL, 0);
+            addColumn("PSF_QF_PERFECT", NULL, 0);
+            addColumn("MOMENTS_R1", NULL, 0);
+            addColumn("MOMENTS_RH", NULL, 0);
+            addColumn("KRON_FLUX", NULL, 0);
+            addColumn("KRON_FLUX_ERR", NULL, 0);
+            addColumn("KRON_FLUX_INNER", NULL, 0);
+            addColumn("KRON_FLUX_OUTER", NULL, 0);
+            addColumn("DIFF_R_P", NULL, 0);
+            addColumn("DIFF_SN_P", NULL, 0);
+            addColumn("DIFF_R_M", NULL, 0);
+            addColumn("DIFF_SN_M", NULL, 0);
+            addColumn("FLAGS2", NULL, 1);
+            addColumn("IPP_IDET", NULL, 0);
+            addColumn("N_FRAMES", NULL, 0);
+            addColumn("PADDING", NULL, 0);
+        }
+        if (!psFitsWriteTableAllColumns(fits, header, table, OUT_EXTNAME)) {
+            psError(psErrorCodeLast(), false, "Unable to write table");
+            return false;
+        }
+        psFree(table);
+    }
+
+    psFree(header);
+    psFitsClose(fits);
+
+    psTrace("ppMops.write", 1, "Done writing %ld rows to %s", det->num, args->output);
+
+    return true;
 }
 
-// extension parameter values:
-// 0: SkyChip.psf
-// 1: SkyChip.xfit
-// Any other value is ignored
-static bool addOutputColumn(psMetadata *table, const psArray *detections, long outputSize, int extension, char *outColumnName, char *inColumnName, bool convertTo32)
+static bool addOutputColumn(psMetadata *table, const psArray *detections, long outputSize, char *outColumnName, char *inColumnName, bool convertTo32)
 {
-  if (inColumnName == NULL) {
-    inColumnName = outColumnName;
-  }
-  psVector *out = NULL;
-  if (convertTo32) {
-    // psFitsReadTableAllColumns reads columns of cfitsio type LONG and ULONG into a 64 bit integers
-    // We want to write 32 bits to the output.
-    int next = 0;
-    for (long i=0; i<detections->n; i++) {
-      ppMopsDetections *det = detections->data[i];
-      if (!det || det->num == 0) {
-        // no detections survived for this input
-        continue;
-      }
-      psVector *in = NULL;
-      if (extension == 0) {
-        in = psMetadataLookupVector(NULL, det->table, inColumnName);
-      } else if (extension == 1) {
-        in = psMetadataLookupVector(NULL, det->fittedTrails, inColumnName);
-      } else {
-        //Error?
-      }
-      if (!in) {
-        psError(PS_ERR_PROGRAMMING, true, "failed to find input column: %s (convertTo32 is true)", inColumnName);
-        return false;
-      }
-      if (in->type.type != PS_TYPE_S64 && in->type.type != PS_TYPE_U64) {
-        psError(PS_ERR_PROGRAMMING, true, "input column to convert is not S64 or U64: %s %d",
-                inColumnName, in->type.type);
-        return false;
-      }
-      if (out == NULL) {
-        // First time through set up the output vector and the copy parameters
-        if (in->type.type == PS_TYPE_S64) {
-          out = psVectorAlloc(outputSize, PS_TYPE_S32);
-        } else {
-          out = psVectorAlloc(outputSize, PS_TYPE_U32);
-        }
-      }
-      for (long d=0; d < det->num; d++) {
-        if (in->type.type == PS_TYPE_S64) {
-          out->data.S32[next++] = in->data.S64[d];
-        } else {
-          out->data.U32[next++] = in->data.U64[d];
-        }
-      }
-    }
-  } else {
-    void *next = NULL;
-    int elementSize = 0;    // size of elements in vector... We are making assumptions here about the organization of primitives in memory so we can use memcopy
-    for (long i=0; i<detections->n; i++) {
-      ppMopsDetections *det = detections->data[i];
-      if (!det || det->num == 0) {
-        // no detections survived for this input
-        continue;
-      }
-      psVector *in = NULL;
-      if (extension == 0) {
-        psTrace("ppMops.write", 1, "Using det->table for %s", inColumnName);
-        in = psMetadataLookupVector(NULL, det->table, inColumnName);
-      } else if (extension == 1) {
-        psTrace("ppMops.write", 1, "Using det->fittedTrails for %s", inColumnName);
-        in = psMetadataLookupVector(NULL, det->fittedTrails, inColumnName);
-      } else {
-        //Error?
-      }
-      if (!in) {
-        psError(PS_ERR_PROGRAMMING, true, "failed to find input column: %s (convertTo32 is false)", inColumnName);
-        out = psVectorAlloc(outputSize, PS_TYPE_F32);
-        psVectorInit(out, NAN);
-        psMetadataAddVector(table, PS_LIST_TAIL, outColumnName, 0, NULL, out);
-        return false;
-      }
-      if (out == NULL) {
-        // First time through set up the output vector and the copy parameters
-        out = psVectorAlloc(outputSize, in->type.type);
-        next = (void *) out->data.U8;
-        switch (in->type.type) {
-        case PS_TYPE_S8:
-          elementSize = sizeof(psS8);
-          break;
-        case PS_TYPE_U8:
-          elementSize = sizeof(psU8);
-          break;
-        case PS_TYPE_S16:
-          elementSize = sizeof(psS16);
-          break;
-        case PS_TYPE_U16:
-          elementSize = sizeof(psU16);
-          break;
-        case PS_TYPE_S32:
-          elementSize = sizeof(psS32);
-          break;
-        case PS_TYPE_U32:
-          elementSize = sizeof(psU32);
-          break;
-        case PS_TYPE_S64:
-          elementSize = sizeof(psS64);
-          break;
-        case PS_TYPE_U64:
-          elementSize = sizeof(psU64);
-          break;
-        case PS_TYPE_F32:
-          elementSize = sizeof(psF32);
-          break;
-        case PS_TYPE_F64:
-          elementSize = sizeof(psF64);
-          break;
-        default:
-          psError(PS_ERR_BAD_PARAMETER_TYPE, true, "Unknown vector type %d", in->type.type);
-          return false;
-        }
-      }
-      // We are doing nasty things here so we can use memcpy. 
-      // It would be safer to do a proper loop over the elements.
-      long toCopy = det->num * elementSize;
-      memcpy(next, in->data.U8, toCopy);
-      next += toCopy;
-    }
-  }
-  // Finally add the new column to the output table
-  psMetadataAddVector(table, PS_LIST_TAIL, outColumnName, 0, NULL, out);
-  psFree(out);    // drop reference
-  return true;
+    if (inColumnName == NULL) {
+        inColumnName = outColumnName;
+    }
+
+    psVector *out = NULL;
+    if (convertTo32) {
+        // psFitsReadTableAllColumns reads columns of cfitsio type LONG and ULONG into a 64 bit integers
+        // We want to write 32 bits to the output.
+        int next = 0;
+        for (long i=0; i<detections->n; i++) {
+            ppMopsDetections *det = detections->data[i];
+            if (!det || det->num == 0) {
+                // no detections survived for this input
+                continue;
+            }
+            psVector *in = psMetadataLookupVector(NULL, det->table, inColumnName);
+            if (!in) {
+                psError(PS_ERR_PROGRAMMING, true, "failed to find input column: %s", inColumnName);
+                return false;
+            }
+            if (in->type.type != PS_TYPE_S64 && in->type.type != PS_TYPE_U64) {
+                psError(PS_ERR_PROGRAMMING, true, "input column to convert is not S64 or U64: %s %d",
+                    inColumnName, in->type.type);
+                return false;
+            }
+            if (out == NULL) {
+                // First time through set up the output vector and the copy parameters
+                if (in->type.type == PS_TYPE_S64) {
+                    out = psVectorAlloc(outputSize, PS_TYPE_S32);
+                } else {
+                    out = psVectorAlloc(outputSize, PS_TYPE_U32);
+                }
+            }
+            for (long d=0; d < det->num; d++) {
+                if (in->type.type == PS_TYPE_S64) {
+                    out->data.S32[next++] = in->data.S64[d];
+                } else {
+                    out->data.U32[next++] = in->data.U64[d];
+                }
+            }
+        }
+    } else {
+        void *next = NULL;
+        int elementSize = 0;    // size of elements in vector... We are making assumptions here about the organization of primitives in memory so we can use memcopy
+        for (long i=0; i<detections->n; i++) {
+            ppMopsDetections *det = detections->data[i];
+            if (!det || det->num == 0) {
+                // no detections survived for this input
+                continue;
+            }
+            psVector *in = psMetadataLookupVector(NULL, det->table, inColumnName);
+            if (!in) {
+                psError(PS_ERR_PROGRAMMING, true, "failed to find input column: %s", inColumnName);
+                return false;
+            }
+            if (out == NULL) {
+                // First time through set up the output vector and the copy parameters
+                out = psVectorAlloc(outputSize, in->type.type);
+                next = (void *) out->data.U8;
+                switch (in->type.type) {
+                    case PS_TYPE_S8:
+                    elementSize = sizeof(psS8);
+                    break;
+                case PS_TYPE_U8:
+                    elementSize = sizeof(psU8);
+                    break;
+                case PS_TYPE_S16:
+                    elementSize = sizeof(psS16);
+                    break;
+                case PS_TYPE_U16:
+                    elementSize = sizeof(psU16);
+                    break;
+                case PS_TYPE_S32:
+                    elementSize = sizeof(psS32);
+                    break;
+                case PS_TYPE_U32:
+                    elementSize = sizeof(psU32);
+                    break;
+                case PS_TYPE_S64:
+                    elementSize = sizeof(psS64);
+                    break;
+                case PS_TYPE_U64:
+                    elementSize = sizeof(psU64);
+                    break;
+                case PS_TYPE_F32:
+                    elementSize = sizeof(psF32);
+                    break;
+                case PS_TYPE_F64:
+                    elementSize = sizeof(psF64);
+                    break;
+                default:
+                    psError(PS_ERR_BAD_PARAMETER_TYPE, true, "Unknown vector type %d", in->type.type);
+                    return false;
+
+                }
+            }
+            // We are doing nasty things here so we can use memcpy. 
+            // It would be safer to do a proper loop over the elements.
+            long toCopy = det->num * elementSize;
+            memcpy(next, in->data.U8, toCopy);
+            next += toCopy;
+        }
+    }
+
+    // Finally add the new column to the output table
+    psMetadataAddVector(table, PS_LIST_TAIL, outColumnName, 0, NULL, out);
+    psFree(out);    // drop reference
+
+    return true;
 }
-
 static bool addSkyfileIDColumn(psMetadata *table, const psArray *detections, long total, char *colName)
 {
-  psVector *out = psVectorAlloc(total, PS_TYPE_S64);
-  long next = 0;
-  for (long i = 0; i<detections->n; i++) {
-    ppMopsDetections *det = detections->data[i];
-    if (!det) {
-      continue;
-    }
-    psS64 diffSkyfileId = det->diffSkyfileId;
-    for (long j = 0; j < det->num; j++) {
-      out->data.S64[next++] = diffSkyfileId;
-    }
-  }
-  psMetadataAddVector(table, PS_LIST_TAIL, colName, 0, NULL, out);
-  psFree(out);
-  return true;
+    psVector *out = psVectorAlloc(total, PS_TYPE_S64);
+    long next = 0;
+    for (long i = 0; i<detections->n; i++) {
+        ppMopsDetections *det = detections->data[i];
+        if (!det) {
+            continue;
+        }
+        psS64 diffSkyfileId = det->diffSkyfileId;
+        for (long j = 0; j < det->num; j++) {
+            out->data.S64[next++] = diffSkyfileId;
+        }
+    }
+    psMetadataAddVector(table, PS_LIST_TAIL, colName, 0, NULL, out);
+    psFree(out);
+    return true;
 }