Changeset 36228

Timestamp:

Oct 22, 2013, 6:27:11 AM (13 years ago)

Author:

eugene

Message:

renaming psImageSmooth_PreAlloc stuff to psImageSmoothCache; setting up for arbitrary shapes

Location:

branches/eam_branches/ipp-20130904

Files:

• psLib/src/imageops/psImageConvolve.c (modified) (4 diffs)
• psLib/src/imageops/psImageConvolve.h (modified) (3 diffs)
• psModules/src/objects/pmPCM_MinimizeChisq.c (modified) (1 diff)
• psModules/src/objects/pmPCMdata.c (modified) (4 diffs)

branches/eam_branches/ipp-20130904/psLib/src/imageops/psImageConvolve.c

@@ -780 +780 @@
 }
 
-void psImageSmooth_PreAlloc_DataFree (psImageSmooth_PreAlloc_Data *smdata) {
+void psImageSmoothCacheDataFree (psImageSmoothCacheData *smdata) {
     psFree (smdata->resultX);
     psFree (smdata->resultY);
@@ -786 +786 @@
 }
 
-psImageSmooth_PreAlloc_Data *psImageSmooth_PreAlloc_DataAlloc (psImage *image, double sigma, double Nsigma) {
+// allocate the psImageSmoothCache data structure, but do not define the kernel
+psImageSmoothCacheData *psImageSmoothCacheAlloc (psImage *image, double sigma, double Nsigma) {
 
     psImageSmooth_PreAlloc_Data *smdata = psAlloc(sizeof(psImageSmooth_PreAlloc_Data));
@@ -807 +808 @@
     smdata->Ny = image->numRows;            // Number of rows
 
-    IMAGE_SMOOTH_GAUSS(smdata->kernel, smdata->Nrange, sigma, F32);
+    // IMAGE_SMOOTH_GAUSS(smdata->kernel, smdata->Nrange, sigma, F32);
        
     // use a temp running buffer for X and Y directions.
@@ -819 +820 @@
 }
 
+// generate a Gaussian smoothing kernel for supplied sigma.  sigma here does not need to match
+// that used to allocate the structure, but it is recommended
+bool psImageSmoothCacheKernel_Gauss (psImageSmoothCacheData *smdata, float sigma) {
+    // check for NULL structure elements?
+    IMAGE_SMOOTH_GAUSS(smdata->kernel, smdata->Nrange, sigma, F32);
+    return true;
+}
+
 // we can use the same DATA structure on multiple images of the same size
-bool psImageSmooth_PreAlloc_F32(psImage *image, psImageSmooth_PreAlloc_Data *smdata)
+bool psImageSmoothCache_F32(psImage *image, psImageSmoothCacheData *smdata)
 {
     PS_ASSERT_IMAGE_NON_NULL(image, false);
+    PS_ASSERT_KERNEL_NON_NULL(smdata->kernel, false);
     // assert on data type
 

branches/eam_branches/ipp-20130904/psLib/src/imageops/psImageConvolve.h

@@ -26 +26 @@
 #define PS_TYPE_KERNEL_NAME "psF32"    ///< the data type for kernel as a string */
 
-/// a structure to contain data related to image smoothing with a 1D gauss kernel
+/// a structure to contain data related to image smoothing with a pre-cached 1D gauss kernel
 typedef struct {
     int Nx;
@@ -34 +34 @@
     psF32 *resultY;
     psVector *kernel;
-} psImageSmooth_PreAlloc_Data;
+} psImageSmoothCacheData;
 
 /// A convolution kernel
@@ -306 +306 @@
 );
 
-psImageSmooth_PreAlloc_Data *psImageSmooth_PreAlloc_DataAlloc (psImage *image, double sigma, double Nsigma);
-bool psImageSmooth_PreAlloc_F32(psImage *image, psImageSmooth_PreAlloc_Data *smdata);
+psImageSmoothCacheData *psImageSmoothCacheAlloc (psImage *image, double sigma, double Nsigma);
+bool psImageSmoothCache_F32(psImage *image, psImageSmoothCacheData *smdata);
 
 /// Control threading for image convolution functions

branches/eam_branches/ipp-20130904/psModules/src/objects/pmPCM_MinimizeChisq.c

@@ -581 +581 @@
 
             float ymodel  = pcm->modelConvFlux->data.F32[i][j];
+
+            // XXXX note this point here:::
             float yweight = 1.0 / source->variance->data.F32[i][j];
             float delta = ymodel - source->pixels->data.F32[i][j];

branches/eam_branches/ipp-20130904/psModules/src/objects/pmPCMdata.c

@@ -248 +248 @@
 }
 
+// generate a Gaussian smoothing kernel for supplied sigma.  sigma here does not need to match
+// that used to allocate the structure, but it is recommended
+bool psImageSmoothCacheKernel_PS1_V1 (psImageSmoothCacheData *smdata, float sigma, float kappa) {
+    // check for NULL structure elements?
+
+    int size = smdata->Nrange + 1;
+
+    smdata->kernel = psVectorAlloc(2 * smdata->Nrange + 1, PS_TYPE_F32);
+
+    double sum = 0.0;                   // Sum of Gaussian, for normalization
+    double factor = 1.0 / (sigma * M_SQRT2);    // Multiplier for i -> z
+
+    // PS1_V1 is a power-law with fitted linear term:
+    // 1 / (1 + kappa z + z^1.666)  where z = (r/sigma)^2
+
+    // generate the kernel (not normalized)
+    for (int i = -size, j = 0; i <= size; i++, j++) {
+        float z = i / sigma / M_SQRT2;
+        sum += smdata->kernel->data.F32[j] = 1.0 / (1 + kappa * z + pow(z,1.666));
+    }
+
+    // renormalize kernel to integral of 1.0
+    for (int i = 0; i < 2 * size + 1; i++) {
+        smdata->kernel->data.F32[i] /= sum;
+    }
+
+    return true;
+}
+
 pmPCMdata *pmPCMinit(pmSource *source, pmSourceFitOptions *fitOptions, pmModel *model, psImageMaskType maskVal, float psfSize) {
 
@@ -312 +341 @@
     pcm->nsigma = 2.0;
 
-    pcm->smdata = psImageSmooth_PreAlloc_DataAlloc (source->pixels, pcm->sigma, pcm->nsigma);
+    // psImageSmoothCacheAlloc generates a structure but does not assign the smoothing vector
+    pcm->smdata = psImageSmoothCacheAlloc (source->pixels, pcm->sigma, pcm->nsigma);
+
+    float v1sigma = 0.5 * (axes.major + axes.minor);
+    psImageSmoothCacheKernel_Gauss (pcm->smdata, v1sigma, PAR[PM_PAR_7]);
+
 # else
     // make sure we save a cached copy of the psf flux
@@ -396 +430 @@
         psFree(pcm->smdata);
         pcm->smdata = psImageSmooth_PreAlloc_DataAlloc (source->pixels, pcm->sigma, pcm->nsigma);
+
+    // psImageSmoothCacheAlloc generates a structure but does not assign the smoothing vector
+    pcm->smdata = psImageSmoothCacheAlloc (source->pixels, pcm->sigma, pcm->nsigma);
+
+    float v1sigma = 0.5 * (axes.major + axes.minor);
+    psImageSmoothCacheKernel_Gauss (pcm->smdata, v1sigma, PAR[PM_PAR_7]);
     }
 
@@ -439 +479 @@
 
         psImageSmooth (source->modelFlux, sigma, nsigma);
+
+        // psImageSmoothCacheAlloc generates a structure but does not assign the smoothing vector
+        pcm->smdata = psImageSmoothCacheAlloc (source->pixels, pcm->sigma, pcm->nsigma);
+        
+        float v1sigma = 0.5 * (axes.major + axes.minor);
+        psImageSmoothCacheKernel_Gauss (pcm->smdata, v1sigma, PAR[PM_PAR_7]);
+
+        psImageSmoothCache_F32 (source->modelFlux, pcm->smdata);
     } else {
         // make sure we save a cached copy of the psf flux