Index: trunk/psphot/src/psphotOutput.c
===================================================================
--- trunk/psphot/src/psphotOutput.c	(revision 5837)
+++ trunk/psphot/src/psphotOutput.c	(revision 5980)
@@ -65,5 +65,5 @@
     }
 
-    float RADIUS = psMetadataLookupF32 (&status, config, "PSF_FIT_RADIUS");
+    float RADIUS = psMetadataLookupF32 (&status, config, "AP_RADIUS");
 
     // write sources with models 
@@ -116,5 +116,5 @@
     }
 
-    float RADIUS = psMetadataLookupF32 (&status, config, "PSF_FIT_RADIUS");
+    float RADIUS = psMetadataLookupF32 (&status, config, "AP_RADIUS");
 
     // write sources with models 
@@ -162,5 +162,5 @@
 
     // find config information for output header
-    float RADIUS = psMetadataLookupF32 (&status, config, "PSF_FIT_RADIUS");
+    float RADIUS = psMetadataLookupF32 (&status, config, "AP_RADIUS");
     float ZERO_POINT = psMetadataLookupF32 (&status, config, "ZERO_POINT");
     if (!status) ZERO_POINT = 25.0;
@@ -171,6 +171,9 @@
     psMetadataAdd (imdata->header, PS_LIST_TAIL, "ZERO_PT",  PS_DATA_F32 | PS_META_REPLACE, "zero point",          ZERO_POINT);
     psMetadataAdd (imdata->header, PS_LIST_TAIL, "APMIFIT",  PS_DATA_F32 | PS_META_REPLACE, "aperture residual",   psf->ApResid);
-    psMetadataAdd (imdata->header, PS_LIST_TAIL, "dAPMIFIT", PS_DATA_F32 | PS_META_REPLACE, "ap residual scatter", psf->dApResid);
+    psMetadataAdd (imdata->header, PS_LIST_TAIL, "DAPMIFIT", PS_DATA_F32 | PS_META_REPLACE, "ap residual scatter", psf->dApResid);
+    psMetadataAdd (imdata->header, PS_LIST_TAIL, "NAPMIFIT", PS_DATA_S32 | PS_META_REPLACE, "ap residual scatter", psf->nApResid);
+    psMetadataAdd (imdata->header, PS_LIST_TAIL, "NPSFSTAR", PS_DATA_S32 | PS_META_REPLACE, "ap residual scatter", psf->nPSFstars);
     psMetadataAdd (imdata->header, PS_LIST_TAIL, "SKYBIAS",  PS_DATA_F32 | PS_META_REPLACE, "aperture sky bias",   psf->skyBias);
+    psMetadataAdd (imdata->header, PS_LIST_TAIL, "SKYSAT",   PS_DATA_F32 | PS_META_REPLACE, "aperture sky bias",   psf->skySat);
     psMetadataAdd (imdata->header, PS_LIST_TAIL, "PHOTCODE", PS_DATA_STRING | PS_META_REPLACE, "photometry code",     PHOTCODE);
     psMetadataAdd (imdata->header, PS_LIST_TAIL, "FWHM_X",   PS_DATA_F32 | PS_META_REPLACE, "PSF FWHM X",          0.0);
@@ -250,5 +253,5 @@
 
     // find config information for output header
-    float RADIUS = psMetadataLookupF32 (&status, config, "PSF_FIT_RADIUS");
+    float RADIUS = psMetadataLookupF32 (&status, config, "AP_RADIUS");
     float ZERO_POINT = psMetadataLookupF32 (&status, config, "ZERO_POINT");
     char *PHOTCODE = psMetadataLookupPtr (&status, config, "PHOTCODE");
@@ -367,5 +370,5 @@
     bool status;
 
-    float RADIUS = psMetadataLookupF32 (&status, config, "PSF_FIT_RADIUS");
+    float RADIUS = psMetadataLookupF32 (&status, config, "AP_RADIUS");
 
     f = fopen (filename, "w");
@@ -400,6 +403,6 @@
 	    fprintf (f, "%9.6f ", dPAR[j]);
 	}
-	fprintf (f, ": %2d %#5x %7.3f %7.1f %7.2f %4d %2d\n", 
-		 source[0].type, source[0].mode, 
+	fprintf (f, ": %7.4f %2d %#5x %7.3f %7.1f %7.2f %4d %2d\n", 
+		 source[0].apMag, source[0].type, source[0].mode, 
 		 log10(model[0].chisq/model[0].nDOF), 
 		 source[0].moments->SN, 
@@ -453,5 +456,6 @@
 	    fprintf (f, "%9.6f ", dPAR[j]);
 	}
-	fprintf (f, ": %2d %#5x %7.3f %7.1f %7.2f %4d %2d\n", 
+	fprintf (f, ": %7.4f  %2d %#5x %7.3f %7.1f %7.2f %4d %2d\n", 
+		 source->apMag, 
 		 source[0].type, source[0].mode,
 		 log10(model[0].chisq/model[0].nDOF), 
@@ -694,4 +698,35 @@
 }
 
+psPolynomial4D *psPolynomial4DfromMD (psMetadata *folder) {
+
+    bool status;
+    char keyword[80];
+
+    // get polynomial orders
+    // XXX add status failures tests
+    int nXorder = psMetadataLookupS32 (&status, folder, "NORDER_X");
+    int nYorder = psMetadataLookupS32 (&status, folder, "NORDER_Y");
+    int nZorder = psMetadataLookupS32 (&status, folder, "NORDER_Z");
+    int nTorder = psMetadataLookupS32 (&status, folder, "NORDER_T");
+
+    psPolynomial4D *poly = psPolynomial4DAlloc (nXorder, nYorder, nZorder, nTorder, PS_POLYNOMIAL_ORD);
+
+    for (int nx = 0; nx < poly->nX + 1; nx++) {
+	for (int ny = 0; ny < poly->nY + 1; ny++) {
+	    for (int nz = 0; nz < poly->nZ + 1; nz++) {
+		for (int nt = 0; nt < poly->nT + 1; nt++) {
+		    sprintf (keyword, "VAL_X%02d_Y%02d_Z%02d_T%02d", nx, ny, nz, nt);
+		    poly->coeff[nx][ny][nz][nt] = psMetadataLookupF32 (&status, folder, keyword);
+		    if (!status) poly->mask[nx][ny][nz][nt] = 1;
+		    
+		    sprintf (keyword, "ERR_X%02d_Y%02d_Z%02d_T%02d", nx, ny, nz, nt);
+		    poly->coeffErr[nx][ny][nz][nt] = psMetadataLookupF32 (&status, folder, keyword);
+		}
+	    }
+	}
+    }
+    return (poly);
+}
+
 // XXX : these may need F64, or %g format for output
 bool psPolynomial2DtoMD (psMetadata *md, psPolynomial2D *poly, char *format, ...) {
@@ -774,4 +809,48 @@
 }
 
+bool psPolynomial4DtoMD (psMetadata *md, psPolynomial4D *poly, char *format, ...) {
+
+    int Nbyte;
+    char tmp;
+    char *root;
+    va_list argp;  
+
+    va_start (argp, format);
+    Nbyte = vsnprintf (&tmp, 0, format, argp);
+    va_end (argp);
+
+    if (!Nbyte) return false;
+
+    va_start (argp, format);
+    root = (char *) psAlloc (Nbyte + 1);
+    memset (root, 0, Nbyte + 1);
+    vsnprintf (root, Nbyte + 1, format, argp);
+    va_end (argp);
+
+    psMetadata *folder = psMetadataAlloc ();
+    psMetadataAdd (md, PS_LIST_TAIL, root, PS_DATA_METADATA, "folder for 4D polynomial", folder);
+    psFree (root);
+
+    // specify the polynomial orders
+    psMetadataAdd (folder, PS_LIST_TAIL, "NORDER_X", PS_DATA_S32, "number of x orders", poly->nX);
+    psMetadataAdd (folder, PS_LIST_TAIL, "NORDER_Y", PS_DATA_S32, "number of y orders", poly->nY);
+    psMetadataAdd (folder, PS_LIST_TAIL, "NORDER_Z", PS_DATA_S32, "number of z orders", poly->nZ);
+    psMetadataAdd (folder, PS_LIST_TAIL, "NORDER_T", PS_DATA_S32, "number of z orders", poly->nT);
+
+    // place polynomial entries on folder
+    for (int nx = 0; nx < poly->nX + 1; nx++) {
+	for (int ny = 0; ny < poly->nY + 1; ny++) {
+	    for (int nz = 0; nz < poly->nZ + 1; nz++) {
+		for (int nt = 0; nt < poly->nT + 1; nt++) {
+		    if (poly->mask[nx][ny][nz][nt]) continue;
+		    psMetadataAdd (folder, PS_LIST_TAIL, "VAL_X%02d_Y%02d_Z%02d_T%02d", PS_DATA_F32, "polynomial coefficient", poly->coeff[nx][ny][nz][nt], nx, ny, nz, nt);
+		    psMetadataAdd (folder, PS_LIST_TAIL, "ERR_X%02d_Y%02d_Z%02d_T%02d", PS_DATA_F32, "polynomial coeffficient error", poly->coeffErr[nx][ny][nz][nt], nx, ny, nz, nt);
+		}
+	    }
+	}
+    }
+    return true;
+}
+
 bool psphotWritePSF (pmPSF *psf, char *filename) {
 
@@ -788,5 +867,5 @@
 	psPolynomial2DtoMD (psfdata, poly, "PSF_PAR%02d", i);
     }
-    psPolynomial3DtoMD (psfdata, psf->ApTrend, "APTREND");
+    psPolynomial4DtoMD (psfdata, psf->ApTrend, "APTREND");
 
     psMetadataAdd (psfdata, PS_LIST_TAIL, "PSF_AP_RESID", PS_DATA_F32, "aperture residual", psf->ApResid);
@@ -828,5 +907,5 @@
     sprintf (keyword, "APTREND");
     psMetadata *folder = psMetadataLookupPtr (&status, psfdata, keyword);
-    psPolynomial3D *poly = psPolynomial3DfromMD (folder);
+    psPolynomial4D *poly = psPolynomial4DfromMD (folder);
     psFree (psf->ApTrend);
     psf->ApTrend = poly;
