Changeset 36375 for trunk/psModules/src/objects/pmPCM_MinimizeChisq.c

Timestamp:

Dec 10, 2013, 2:55:11 PM (12 years ago)

Author:

eugene

Message:

merge changes from eam_branches/ipp-20130904

Location:

trunk

Files:

• . (modified) (1 prop)
• psModules (modified) (1 prop)
• psModules/src/objects/pmPCM_MinimizeChisq.c (modified) (8 diffs)

trunk/psModules/src/objects/pmPCM_MinimizeChisq.c

@@ -135 +135 @@
         }
 
+        if (min->isInteractive) {
+            fprintf (stderr, "%d : ", min->iter);
+            for (int ti = 0; ti < params->n; ti++) {
+                fprintf (stderr, "%f  ", params->data.F32[ti]);
+            }
+            fprintf (stderr, " : %f\n", min->value);
+        }
+
         char key[10]; // used for interactive responses
         bool testValue = false;
@@ -140 +148 @@
         // set a new guess for Alpha, Beta, Params
         if (!psMinLM_GuessABP(Alpha, Beta, Params, alpha, beta, params, paramMask, checkLimits, lambda, &dLinear)) {
-            if (min->isInteractive) {
+            if (false && min->isInteractive) {
                 fprintf (stdout, "guess failed (singular matrix or NaN values), continue? [Y,n] ");
                 if (!fgets(key, 8, stdin)) {
@@ -167 +175 @@
         }
 
-        if (min->isInteractive) {
+        if (false && min->isInteractive) {
             p_psVectorPrint(psTraceGetDestination(), Params, "current parameters: ");
             fprintf (stdout, "last chisq : %f\n", min->value);
@@ -473 +481 @@
 # else
     if (pcm->use1Dgauss) {
-        // do not use the threaded, mask-aware version of this code (psImageSmoothMaskPixelsThread):
-        // * the model flux is not masked
-        // * threading takes place above this level
-        pcm->modelConvFlux = psImageCopy (pcm->modelConvFlux, pcm->modelFlux, pcm->modelFlux->type.type);
-        psImageSmooth_PreAlloc_F32 (pcm->modelConvFlux, pcm->smdata);
-        // psImageSmooth (pcm->modelConvFlux, pcm->sigma, pcm->nsigma);
+
+        if (USE_1D_CACHE) {
+            // do not use the threaded, mask-aware version of this code (psImageSmoothMaskPixelsThread):
+            // * the model flux is not masked
+            // * threading takes place above this level
+            pcm->modelConvFlux = psImageCopy (pcm->modelConvFlux, pcm->modelFlux, pcm->modelFlux->type.type);
+            psImageSmoothCache_F32 (pcm->modelConvFlux, pcm->smdata);
+        } else {
+            pcm->modelConvFlux = psImageCopy (pcm->modelConvFlux, pcm->modelFlux, pcm->modelFlux->type.type);
+            psImageSmooth2dCache_F32 (pcm->modelConvFlux, pcm->smdata2d);
+        }
     } else {
         psImageConvolveKernel (pcm->modelConvFlux, pcm->modelFlux, NULL, 0, pcm->psfFFT);
@@ -493 +506 @@
 # else
         if (pcm->use1Dgauss) {
-            // do not use the threaded, mask-aware version of this code (psImageSmoothMaskPixelsThread):
-            // * the model flux is not masked
-            // * threading takes place above this level
-            dmodelConv = psImageCopy (dmodelConv, dmodel, dmodel->type.type);
-            psImageSmooth_PreAlloc_F32 (dmodelConv, pcm->smdata);
-            // psImageSmooth (dmodelConv, pcm->sigma, pcm->nsigma);
+            if (USE_1D_CACHE) {
+                // do not use the threaded, mask-aware version of this code (psImageSmoothMaskPixelsThread):
+                // * the model flux is not masked
+                // * threading takes place above this level
+                dmodelConv = psImageCopy (dmodelConv, dmodel, dmodel->type.type);
+                psImageSmoothCache_F32 (dmodelConv, pcm->smdata);
+            } else {
+                dmodelConv = psImageCopy (dmodelConv, dmodel, dmodel->type.type);
+                psImageSmooth2dCache_F32 (dmodelConv, pcm->smdata2d);
+            }
         } else {
             psImageConvolveKernel (dmodelConv, dmodel, NULL, 0, pcm->psfFFT);
@@ -514 +531 @@
 
         if (pcm->use1Dgauss) {
-            // do not use the threaded, mask-aware version of this code (psImageSmoothMaskPixelsThread):
-            // * the model flux is not masked
-            // * threading takes place above this level
-            dmodelConv = psImageCopy (dmodelConv, dmodel, dmodel->type.type);
-            psImageSmooth_PreAlloc_F32 (dmodelConv, pcm->smdata);
-            // psImageSmooth (dmodelConv, pcm->sigma, pcm->nsigma);
+            if (USE_1D_CACHE) {
+                // do not use the threaded, mask-aware version of this code (psImageSmoothMaskPixelsThread):
+                // * the model flux is not masked
+                // * threading takes place above this level
+                dmodelConv = psImageCopy (dmodelConv, dmodel, dmodel->type.type);
+                psImageSmoothCache_F32 (dmodelConv, pcm->smdata);
+            } else {
+                dmodelConv = psImageCopy (dmodelConv, dmodel, dmodel->type.type);
+                psImageSmooth2dCache_F32 (dmodelConv, pcm->smdata2d);
+            }
         } else {
             psImageConvolveFFT (dmodelConv, dmodel, NULL, 0, pcm->psf);
@@ -541 +562 @@
     static int Npass = 0;
     char name[128]; 
-    snprintf (name, 128, "psf.%03d.fits", Npass); psphotSaveImage (NULL, pcm->psf->image, name);
+    if (!pcm->use1Dgauss) {
+      snprintf (name, 128, "psf.%03d.fits", Npass); psphotSaveImage (NULL, pcm->psf->image, name);
+    }
     snprintf (name, 128, "mod.%03d.fits", Npass); psphotSaveImage (NULL, pcm->modelFlux, name);
     snprintf (name, 128, "cnv.%03d.fits", Npass); psphotSaveImage (NULL, pcm->modelConvFlux, name);
@@ -579 +602 @@
 
             float ymodel  = pcm->modelConvFlux->data.F32[i][j];
-            float yweight = 1.0 / source->variance->data.F32[i][j];
+
+            // XXXX note this point here:::
+            float yweight = pcm->poissonErrors ? 1.0 / source->variance->data.F32[i][j] : 1.0;
             float delta = ymodel - source->pixels->data.F32[i][j];