Index: /branches/eam_branches/ipp-20100823/psModules/test/objects/tap_pmSourceMoments.c
===================================================================
--- /branches/eam_branches/ipp-20100823/psModules/test/objects/tap_pmSourceMoments.c	(revision 29444)
+++ /branches/eam_branches/ipp-20100823/psModules/test/objects/tap_pmSourceMoments.c	(revision 29444)
@@ -0,0 +1,60 @@
+#include <stdio.h>
+#include <string.h>
+#include <pslib.h>
+#include <psmodules.h>
+#include "tap.h"
+#include "pstap.h"
+
+#define VERBOSE                 0
+#define ERR_TRACE_LEVEL         0
+
+float Gaussian(float Io, float sigma, float Xc, float Yc, int ix, int iy);
+
+int main(int argc, char* argv[])
+{
+    psLogSetFormat("HLNM");
+    psLogSetLevel(PS_LOG_INFO);
+    psTraceSetLevel("err", ERR_TRACE_LEVEL);
+    plan_tests(70);
+
+    // ----------------------------------------------------------------------
+    // pmSourceMomentsGetCentroid() tests
+    {
+        psMemId id = psMemGetId();
+
+	// generate a sample source (no noise)
+	pmSource *source = pmSourceAlloc();
+        ok(source, "source allocated");
+	
+	// need to have: peak, pixels, variance, 
+	int Nx = 100;
+	int Ny = 100;
+	float Xc = 52.0;
+	float Yc = 48.0;
+	float Io = 1000.0;
+
+	source->peak = pmPeakAlloc(Xc, Yc, Io, PS_PEAK_LONE);
+	source->moments = pmMomentsAlloc();
+	source->pixels = psImageAlloc(Nx, Ny, PS_TYPE_F32);
+	
+	// populate the pixels with an object (Gaussian with Io, sigma)
+	for (int iy = 0; iy < source->pixels->numRows; iy++) {
+	  for (int ix = 0; ix < source->pixels->numCols; ix++) {
+	    source->pixels->data.F32[iy][ix] = Gaussian(Io, sigma, Xc, Yc, ix, iy);
+	  }
+	}
+
+        psFree(growthCurve);
+        ok(!psMemCheckLeaks (id, NULL, NULL, false), "no memory leaks");
+    }
+
+}
+
+float Gaussian(float Io, float sigma, float Xc, float Yc, int ix, int iy) {
+  
+  float radius2 = PS_SQR(ix + 0.5 - Xc) + PS_SQR(iy + 0.5 - Yc);
+  float exponent = -0.5 * radius2 / PS_SQR(sigma);
+
+  float value = Io * exp(exponent);
+  return value;
+}
