Index: /trunk/psphot/src/psphotImageMedian.c
===================================================================
--- /trunk/psphot/src/psphotImageMedian.c	(revision 10171)
+++ /trunk/psphot/src/psphotImageMedian.c	(revision 10172)
@@ -10,4 +10,5 @@
     psRegion region;
     int MAX_SAMPLE_PIXELS;
+    static char *defaultStatsName = "FITTED_MEAN";
 
     pmReadout *model = NULL;
@@ -35,4 +36,14 @@
     if (!status) {
         SKY_BIAS = 0;
+    }
+
+    // supply the sky background statistic name:
+    char *statsName = psMetadataLookupStr (&status, recipe, "SKY_STAT");
+    if (statsName == NULL) {
+	statsName = defaultStatsName;
+    }
+    psStatsOptions statsOption = psStatsOptionFromString (statsName);
+    if (!(statsOption & (PS_STAT_SAMPLE_MEAN | PS_STAT_SAMPLE_MEDIAN | PS_STAT_ROBUST_MEDIAN | PS_STAT_ROBUST_QUARTILE | PS_STAT_CLIPPED_MEAN | PS_STAT_FITTED_MEAN))) {
+	statsOption = PS_STAT_FITTED_MEAN;
     }
 
@@ -105,6 +116,10 @@
 
             // XXX the value of the upper and lower cuts probably should be studied...
-            psStats *stats = psImageBackground(subset, submask, 0xff, 0.25, 0.75, MAX_SAMPLE_PIXELS, rng);
-            modelData[iy][ix] = stats->robustMedian + SKY_BIAS;
+            psStats *stats = psImageBackground(subset, submask, 0xff, 0.25, 0.75, MAX_SAMPLE_PIXELS, statsOption, rng);
+	    if (stats->options & PS_STAT_ROBUST_QUARTILE) {
+		modelData[iy][ix] = stats->robustMedian + SKY_BIAS;
+	    } else {	
+		modelData[iy][ix] = psStatsGetValue (stats, statsOption) + SKY_BIAS;
+	    }
 
             psFree (stats);