Index: trunk/psModules/src/objects/pmSourceMatch.c
===================================================================
--- trunk/psModules/src/objects/pmSourceMatch.c	(revision 23229)
+++ trunk/psModules/src/objects/pmSourceMatch.c	(revision 23233)
@@ -152,5 +152,5 @@
 
 
-psArray *pmSourceMatchSources(const psArray *sourceArrays, float radius)
+psArray *pmSourceMatchSources(const psArray *sourceArrays, float radius, bool cullSingles)
 {
     PS_ASSERT_ARRAY_NON_NULL(sourceArrays, NULL);
@@ -220,38 +220,38 @@
 
         // Match with the master list
-
-        psTree *tree = psTreePlant(2, SOURCES_MAX_LEAF, xMaster, yMaster); // kd Tree with sources
-        long numMatch = 0;              // Number of matches
-
-        long size = numMaster + numSources; // New size
-        xMaster = psVectorRealloc(xMaster, size);
-        yMaster = psVectorRealloc(yMaster, size);
-        matches = psArrayRealloc(matches, size);
-
-        psVector *coords = psVectorAlloc(2, PS_TYPE_F32); // Coordinates for tree lookup
-        for (int j = 0; j < numSources; j++) {
-            coords->data.F32[0] = xImage->data.F32[j];
-            coords->data.F32[1] = yImage->data.F32[j];
-            long index = psTreeNearestWithin(tree, coords, radius); // Match index
-            if (index >= 0) {
-                // Record the match
-                pmSourceMatch *match = matches->data[index]; // Match data
-                pmSourceMatchAdd(match, magImage->data.F32[j], magErrImage->data.F32[j], i, j);
-                numMatch++;
-            } else {
-                // Add to the master list
-                pmSourceMatch *match = pmSourceMatchAlloc(numImages); // Match data
-                pmSourceMatchAdd(match, magImage->data.F32[j], magErrImage->data.F32[j], i, j);
-                xMaster->data.F32[numMaster] = xImage->data.F32[j];
-                yMaster->data.F32[numMaster] = yImage->data.F32[j];
-                matches->data[numMaster] = match;
-                numMaster++;
-                xMaster->n = yMaster->n = matches->n = numMaster;
-            }
-
-        }
-        psFree(coords);
-        psFree(tree);
-
+        {
+            psTree *tree = psTreePlant(2, SOURCES_MAX_LEAF, xMaster, yMaster); // kd Tree with sources
+            long numMatch = 0;              // Number of matches
+
+            long size = numMaster + numSources; // New size
+            xMaster = psVectorRealloc(xMaster, size);
+            yMaster = psVectorRealloc(yMaster, size);
+            matches = psArrayRealloc(matches, size);
+
+            psVector *coords = psVectorAlloc(2, PS_TYPE_F32); // Coordinates for tree lookup
+            for (int j = 0; j < numSources; j++) {
+                coords->data.F32[0] = xImage->data.F32[j];
+                coords->data.F32[1] = yImage->data.F32[j];
+                long index = psTreeNearestWithin(tree, coords, radius); // Match index
+                if (index >= 0) {
+                    // Record the match
+                    pmSourceMatch *match = matches->data[index]; // Match data
+                    pmSourceMatchAdd(match, magImage->data.F32[j], magErrImage->data.F32[j], i, j);
+                    numMatch++;
+                } else {
+                    // Add to the master list
+                    pmSourceMatch *match = pmSourceMatchAlloc(numImages); // Match data
+                    pmSourceMatchAdd(match, magImage->data.F32[j], magErrImage->data.F32[j], i, j);
+                    xMaster->data.F32[numMaster] = xImage->data.F32[j];
+                    yMaster->data.F32[numMaster] = yImage->data.F32[j];
+                    matches->data[numMaster] = match;
+                    numMaster++;
+                    xMaster->n = yMaster->n = matches->n = numMaster;
+                }
+
+            }
+            psFree(coords);
+            psFree(tree);
+        }
 
     match_image_done:
@@ -263,26 +263,71 @@
     }
 
-    // Now cull the matches that contain only a single star
-    int numGood = 0;                    // Number of good matches
+    if (cullSingles) {
+        // Now cull the matches that contain only a single star
+        int numGood = 0;                    // Number of good matches
+        for (int i = 0; i < matches->n; i++) {
+            pmSourceMatch *match = matches->data[i]; // Match of interest
+            if (match->num > 1) {
+                if (i != numGood) {
+                    psFree(matches->data[numGood]);
+                    matches->data[numGood] = psMemIncrRefCounter(match);
+                }
+                numGood++;
+            }
+        }
+        matches->n = numGood;
+        for (int i = numGood; i < numMaster; i++) {
+            psFree(matches->data[i]);
+            matches->data[i] = NULL;
+        }
+    }
+
+    return matches;
+}
+
+
+psArray *pmSourceMatchMerge(psArray *sourceArrays, float radius)
+{
+    PS_ASSERT_ARRAY_NON_NULL(sourceArrays, NULL);
+    PS_ASSERT_FLOAT_LARGER_THAN(radius, 0.0, NULL);
+
+    psArray *matches = pmSourceMatchSources(sourceArrays, radius, false); // Source matches
+    if (!matches) {
+        psError(PS_ERR_UNKNOWN, false, "Unable to match source lists.");
+        return NULL;
+    }
+
+    int index = 0;                      // Index to current position on list
     for (int i = 0; i < matches->n; i++) {
+        pmSource *source = NULL;        // Source to put in merged list
         pmSourceMatch *match = matches->data[i]; // Match of interest
-        if (match->num > 1) {
-            if (i != numGood) {
-                psFree(matches->data[numGood]);
-                matches->data[numGood] = psMemIncrRefCounter(match);
-                //                psFree(match);
-            }
-            numGood++;
-        }
-    }
-    matches->n = numGood;
-    for (int i = numGood; i < numMaster; i++) {
-        psFree(matches->data[i]);
-        matches->data[i] = NULL;
-    }
+        for (int j = 0; j < match->num && !source; j++) {
+            if (!isfinite(match->mag->data.F32[j]) || !isfinite(match->magErr->data.F32[j])) {
+                continue;
+            }
+            int image = match->image->data.S32[j]; // Index of image
+            int index = match->index->data.S32[j]; // Index of source for image
+            psArray *list = sourceArrays->data[image]; // List of interest
+            source = list->data[index];
+            break;
+        }
+
+        if (source) {
+            psFree(matches->data[index]);
+            matches->data[index] = psMemIncrRefCounter(source);
+            index++;
+        }
+    }
+
+    // Clear out the rest of the list
+    int num = index;                    // Number of good sources
+    for (; index < matches->n; index++) {
+        psFree(matches->data[index]);
+        matches->data[index] = NULL;
+    }
+    matches->n = num;
 
     return matches;
 }
-
 
 // Iterate on the star magnitudes and image transparencies
Index: trunk/psModules/src/objects/pmSourceMatch.h
===================================================================
--- trunk/psModules/src/objects/pmSourceMatch.h	(revision 23229)
+++ trunk/psModules/src/objects/pmSourceMatch.h	(revision 23233)
@@ -48,7 +48,15 @@
 ///
 /// Returns an array of psSourceMatch
-psArray *pmSourceMatchSources(const psArray *sourceArrays, // Array of arrays of sources on each image
-                              float radius // Correlation radius
+psArray *pmSourceMatchSources(const psArray *sourceArrays, ///< Array of arrays of sources on each image
+                              float radius, ///< Matching radius
+                              bool cullSingles ///< Cull "matches" with only a single source?
                               );
+
+/// Merge two source lists
+///
+/// Sources are pulled from the lists into a new list, with no effort made to adjust them.
+psArray *pmSourceMatchMerge(psArray *sourceArrays, ///< Array of arrays of sources on each image
+                            float radius ///< Matching radius
+    );
 
 /// Perform relative photometry to calibrate images