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

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

File:

• branches/eam_branches/ipp-20130904/psModules/src/objects/pmPCMdata.c (modified) (4 diffs)

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