Index: /trunk/ppSub/src/Makefile.am
===================================================================
--- /trunk/ppSub/src/Makefile.am	(revision 20610)
+++ /trunk/ppSub/src/Makefile.am	(revision 20611)
@@ -7,8 +7,9 @@
 	ppSub.c			\
 	ppSubArguments.c	\
+	ppSubBackground.c	\
 	ppSubCamera.c		\
 	ppSubLoop.c		\
 	ppSubReadout.c		\
-	ppSubVersion.c            
+	ppSubVersion.c
 
 ppSubKernel_CPPFLAGS = $(PSLIB_CFLAGS) $(PSMODULE_CFLAGS) $(PPSUB_CFLAGS)
Index: /trunk/ppSub/src/ppSub.h
===================================================================
--- /trunk/ppSub/src/ppSub.h	(revision 20610)
+++ /trunk/ppSub/src/ppSub.h	(revision 20611)
@@ -27,4 +27,10 @@
     );
 
+/// Higher-order background subtraction
+bool ppSubBackground(pmReadout *ro,     ///< Readout of interest
+                     pmConfig *config,  ///< Configuration
+                     const pmFPAview *view ///< View to readout
+    );
+
 /// Put the program version information into a metadata
 void ppSubVersionMetadata(psMetadata *metadata ///< Metadata to populate
Index: /trunk/ppSub/src/ppSubBackground.c
===================================================================
--- /trunk/ppSub/src/ppSubBackground.c	(revision 20611)
+++ /trunk/ppSub/src/ppSubBackground.c	(revision 20611)
@@ -0,0 +1,79 @@
+#ifdef HAVE_CONFIG_H
+#include <config.h>
+#endif
+
+#include <stdio.h>
+#include <pslib.h>
+#include <psmodules.h>
+#include <psphot.h>
+#include "ppSub.h"
+
+// Copied from ppImageSubtractBackground()
+bool ppSubBackground(pmReadout *ro, pmConfig *config, const pmFPAview *view)
+{
+    psAssert(config, "Need configuration");
+    psAssert(view, "Need view to chip");
+
+    // The background model file may be defined, though the model hasn't been built
+    // If so, we will build it below.
+    bool status;                        // Status of lookup
+    pmFPAfile *modelFile = psMetadataLookupPtr(&status, config->files, "PSPHOT.BACKMDL"); // Background model
+    if (!status) {
+        psError(PS_ERR_UNEXPECTED_NULL, true, "PSPHOT.BACKMDL file is not defined");
+        return false;
+    }
+
+    psMetadata *ppSubRecipe = psMetadataLookupPtr(NULL, config->recipes, PPSUB_RECIPE);
+    psAssert(ppSubRecipe, "Need PPSUB recipe");
+    psMetadata *psphotRecipe = psMetadataLookupPtr(NULL, config->recipes, PSPHOT_RECIPE);
+    psAssert(psphotRecipe, "Need PSPHOT recipe");
+
+    psString maskBadStr = psMetadataLookupStr(NULL, ppSubRecipe, "MASK.BAD"); // Name of bits to mask for bad
+    psMaskType maskBad = pmConfigMaskGet(maskBadStr, config); // Bits to mask for bad pixels
+
+    // user-defined masks to test for good/bad pixels (build from recipe list if not yet set)
+    psMetadataAddU8(psphotRecipe, PS_LIST_TAIL, "MASK.PSPHOT", PS_META_REPLACE, "user-defined mask", maskBad);
+
+    psImage *image = ro->image, *mask = ro->mask; // Image and mask of interest
+
+    pmReadout *modelRO = pmFPAviewThisReadout(view, modelFile->fpa); // Background model
+    if (!modelRO) {
+        // Create the background model
+        if (!psphotModelBackground(config, view, "PPIMAGE.CHIP")) {
+            psError(PS_ERR_UNKNOWN, false, "Unable to model background");
+            return false;
+        }
+        modelRO = (modelFile->mode == PM_FPA_MODE_INTERNAL) ? modelFile->readout :
+                   pmFPAviewThisReadout(view, modelFile->fpa);
+        if (!modelRO) {
+            psError(PS_ERR_UNEXPECTED_NULL, true, "Unable to find background model");
+            return false;
+        }
+    }
+    psImageBinning *binning = psMetadataLookupPtr(&status, psphotRecipe,
+                                                  "PSPHOT.BACKGROUND.BINNING"); // Binning for model
+    psImage *modelImage = modelRO->image; // Background model
+
+    // Do the background subtraction
+    int numCols = image->numCols, numRows = image->numRows; // Size of image
+    for (int y = 0; y < numRows; y++) {
+        for (int x = 0; x < numCols; x++) {
+            if (mask && mask->data.PS_TYPE_MASK_DATA[y][x] & maskBad) {
+                image->data.F32[y][x] = 0.0;
+            } else {
+                float value = psImageUnbinPixel_V2(x, y, modelImage, binning); // Background value
+                if (!isfinite(value)) {
+                    image->data.F32[y][x] = NAN;
+                    mask->data.PS_TYPE_MASK_DATA[y][x] |= maskBad;
+                } else {
+                    image->data.F32[y][x] -= value;
+                }
+            }
+        }
+    }
+
+    pmFPAfileDropInternal(config->files, "PSPHOT.BACKMDL");
+    pmFPAfileDropInternal(config->files, "PSPHOT.BACKMDL.STDEV");
+
+    return true;
+}
Index: /trunk/ppSub/src/ppSubReadout.c
===================================================================
--- /trunk/ppSub/src/ppSubReadout.c	(revision 20610)
+++ /trunk/ppSub/src/ppSubReadout.c	(revision 20611)
@@ -623,4 +623,11 @@
     }
 
+    // Higher order background subtraction using psphot
+    if (!ppSubBackground(outRO, config, view)) {
+        psError(PS_ERR_UNKNOWN, false, "Unable to subtract background.");
+        psFree(outRO);
+        return false;
+    }
+
     // Renormalising for pixels, because that's what magic desires
     if (psMetadataLookupBool(&mdok, recipe, "RENORM")) {