Index: trunk/psLib/src/imageops/psImageCovariance.c
===================================================================
--- trunk/psLib/src/imageops/psImageCovariance.c	(revision 24820)
+++ trunk/psLib/src/imageops/psImageCovariance.c	(revision 24832)
@@ -282,4 +282,55 @@
 }
 
+psKernel *psImageCovarianceAverageWeighted(const psArray *array, const psVector *weights)
+{
+    PS_ASSERT_ARRAY_NON_NULL(array, NULL);
+    PS_ASSERT_ARRAY_NON_EMPTY(array, NULL);
+    if (!weights) {
+        return psImageCovarianceAverage(array);
+    }
+    PS_ASSERT_VECTOR_TYPE(weights, PS_TYPE_F32, NULL);
+
+    int xMin = INT_MAX, xMax = INT_MIN, yMin = INT_MAX, yMax = INT_MIN; // Range for covariance
+    double sumWeights = 0.0;            // Sum of weights
+    for (int i = 0; i < array->n; i++) {
+        psKernel *covar = array->data[i]; // Covariance matrix
+        if (!covar) {
+            continue;
+        }
+        xMin = PS_MIN(xMin, covar->xMin);
+        xMax = PS_MAX(xMax, covar->xMax);
+        yMin = PS_MIN(yMin, covar->yMin);
+        yMax = PS_MAX(yMax, covar->yMax);
+        sumWeights += weights->data.F32[i];
+    }
+    if (sumWeights == 0) {
+        psError(PS_ERR_BAD_PARAMETER_SIZE, true, "No covariance matrices supplied for summation");
+        return NULL;
+    }
+
+    psKernel *sum = psKernelAlloc(xMin, xMax, yMin, yMax); // Summed covariance
+    for (int i = 0; i < array->n; i++) {
+        psKernel *covar = array->data[i]; // Covariance matrix
+        if (!covar) {
+            continue;
+        }
+        for (int y = covar->yMin; y <= covar->yMax; y++) {
+            for (int x = covar->xMin; x <= covar->xMax; x++) {
+                if (!isfinite(covar->kernel[y][x])) {
+                    psError(PS_ERR_BAD_PARAMETER_VALUE, true,
+                            "Non-finite covariance matrix element at %d,%d for input %d",
+                            x, y, i);
+                    psFree(sum);
+                    return NULL;
+                }
+                sum->kernel[y][x] += weights->data.F32[i] * covar->kernel[y][x];
+            }
+        }
+    }
+    psBinaryOp(sum->image, sum->image, "/", psScalarAlloc((float)sumWeights, PS_TYPE_F32));
+
+    return sum;
+}
+
 
 psKernel *psImageCovarianceTruncate(const psKernel *covar, float frac)
Index: trunk/psLib/src/imageops/psImageCovariance.h
===================================================================
--- trunk/psLib/src/imageops/psImageCovariance.h	(revision 24820)
+++ trunk/psLib/src/imageops/psImageCovariance.h	(revision 24832)
@@ -57,4 +57,10 @@
     );
 
+/// Weighted average of multiple covariance pseudo-matrices
+psKernel *psImageCovarianceAverageWeighted(
+    const psArray *array,               ///< Array of covariance pseudo-matrices
+    const psVector *weights             ///< Weights for each (F32)
+    );
+
 /// Truncate covariance pseudo-matrix
 ///