Changeset 15756 for trunk/psModules/src/camera/pmReadoutFake.c

Timestamp:

Dec 6, 2007, 3:57:15 PM (18 years ago)

Author:

Paul Price

Message:

Implementing dual-convolution. This required reworking those
functions involved with calculating and solving the equations; moved
these into a separate file and made various other cleanups (e.g.,
assertions for pmSubtractionKernels). Removed the normalisation
(central pixel) term from all kernels, because this shouldn't be in
the second kernel (there's only one normalisation term between the two
kernels); the normalisation is treated explicitly in the equations,
along with the background (still only a constant background is
supported for now, but there's a lot of support for a polynomial
background for when I get around to putting it in the equations).
Tested the single convolution, and it's working; same results as
before, apparently. Haven't tested dual convolution yet.

File:

• trunk/psModules/src/camera/pmReadoutFake.c (modified) (4 diffs)

trunk/psModules/src/camera/pmReadoutFake.c

@@ -20 +20 @@
 #include "pmSource.h"
 #include "pmSourceUtils.h"
+#include "pmModelUtils.h"
 
 #define MODEL_TYPE "PS_MODEL_RGAUSS"     // Type of model to use
 
 
-pmReadout *pmReadoutFakeFromSources(int numCols, int numRows, const psArray *sources,
-                                    float fwhm, float minFlux)
+bool pmReadoutFakeFromSources(pmReadout *readout, int numCols, int numRows, const psArray *sources,
+                              const psVector *xOffset, const psVector *yOffset, pmPSF *psf, float minFlux,
+                              int radius)
 {
-    PS_ASSERT_INT_LARGER_THAN(numCols, 0, NULL);
-    PS_ASSERT_INT_LARGER_THAN(numRows, 0, NULL);
-    PS_ASSERT_ARRAY_NON_NULL(sources, NULL);
-    PS_ASSERT_FLOAT_LARGER_THAN(fwhm, 0.0, NULL);
-    PS_ASSERT_FLOAT_LARGER_THAN(minFlux, 0.0, NULL);
+    PS_ASSERT_PTR_NON_NULL(readout, false);
+    PS_ASSERT_INT_LARGER_THAN(numCols, 0, false);
+    PS_ASSERT_INT_LARGER_THAN(numRows, 0, false);
+    PS_ASSERT_ARRAY_NON_NULL(sources, false);
+    if (xOffset || yOffset) {
+        PS_ASSERT_VECTOR_NON_NULL(xOffset, false);
+        PS_ASSERT_VECTOR_NON_NULL(yOffset, false);
+        PS_ASSERT_VECTORS_SIZE_EQUAL(xOffset, yOffset, false);
+        PS_ASSERT_VECTOR_TYPE(xOffset, PS_TYPE_S32, false);
+        PS_ASSERT_VECTOR_TYPE_EQUAL(xOffset, yOffset, false);
+        if (xOffset->n != sources->n) {
+            psError(PS_ERR_BAD_PARAMETER_SIZE, true,
+                    "Number of offset vectors (%ld) and sources (%ld) doesn't match",
+                    xOffset->n, sources->n);
+            return false;
+        }
+    }
+    PS_ASSERT_PTR_NON_NULL(psf, false);
+    if (radius > 0 && isfinite(minFlux) && minFlux > 0.0) {
+        psError(PS_ERR_BAD_PARAMETER_VALUE, true, "Cannot define both minimum flux and fixed radius.");
+        return false;
+    }
 
-    pmReadout *readout = pmReadoutAlloc(NULL); // Output readout
-    readout->image = psImageAlloc(numCols, numRows, PS_TYPE_F32);
+
+    readout->image = psImageRecycle(readout->image, numCols, numRows, PS_TYPE_F32);
     psImageInit(readout->image, 0);
 
     int numSources = sources->n;          // Number of stars
 
+#if 0
     pmModelType modelType = pmModelClassGetType(MODEL_TYPE); // Type of PSF model
     assert(modelType >= 0);
@@ -65 +85 @@
         psAbort("Unsupported model type: %d", modelType);
     }
+#endif
 
+    pmModel *fakeModel = pmModelFromPSFforXY(psf, (float)numCols / 2.0, (float)numRows / 2.0,
+                                             1.0); // Fake model, with central intensity of 1.0
     float flux0 = fakeModel->modelFlux(fakeModel->params); // Flux for central intensity of 1.0
-
 
     for (int i = 0; i < numSources; i++) {
@@ -93 +115 @@
         pmSource *fakeSource = pmSourceAlloc(); // Fake source to generate
         fakeSource->peak = pmPeakAlloc(x, y, fakeModel->params->data.F32[PM_PAR_I0], PM_PEAK_LONE);
-        float radius = fakeModel->modelRadius(fakeModel->params, minFlux); // Radius of interest for source
+        float fakeRadius = radius > 0 ? radius : fakeModel->modelRadius(fakeModel->params, minFlux); // Radius
 
-        if (!pmSourceDefinePixels(fakeSource, readout, x, y, radius)) {
-            psError(PS_ERR_UNKNOWN, false, "Unable to define pixels for source.");
-            psFree(readout);
-            psFree(fakeModel);
-            return NULL;
-        }
-
-        if (!pmModelAdd(fakeSource->pixels, NULL, fakeModel, PM_MODEL_OP_FULL, 0)) {
-            psError(PS_ERR_UNKNOWN, false, "Unable to add model of source to image.");
-            psFree(readout);
-            psFree(fakeModel);
-            return NULL;
+        if (xOffset) {
+            if (!pmSourceDefinePixels(fakeSource, readout, x + xOffset->data.S32[i],
+                                      y + yOffset->data.S32[i], fakeRadius)) {
+                psError(PS_ERR_UNKNOWN, false, "Unable to define pixels for source.");
+                psFree(readout);
+                psFree(fakeModel);
+                return false;
+            }
+            if (!pmModelAddWithOffset(fakeSource->pixels, NULL, fakeModel, PM_MODEL_OP_FULL, 0,
+                                      - xOffset->data.S32[i], - yOffset->data.S32[i])) {
+                psError(PS_ERR_UNKNOWN, false, "Unable to add model of source to image.");
+                psFree(readout);
+                psFree(fakeModel);
+                return false;
+            }
+        } else {
+            if (!pmSourceDefinePixels(fakeSource, readout, x, y, fakeRadius)) {
+                psError(PS_ERR_UNKNOWN, false, "Unable to define pixels for source.");
+                psFree(readout);
+                psFree(fakeModel);
+                return false;
+            }
+            if (!pmModelAdd(fakeSource->pixels, NULL, fakeModel, PM_MODEL_OP_FULL, 0)) {
+                psError(PS_ERR_UNKNOWN, false, "Unable to add model of source to image.");
+                psFree(readout);
+                psFree(fakeModel);
+                return false;
+            }
         }
         psFree(fakeSource);
@@ -113 +151 @@
     psFree(fakeModel);
 
-    return readout;
+    return true;
 }