Changeset 27838 for branches/tap_branches/psModules/src/objects/pmPSFtryModel.c

Timestamp:

May 3, 2010, 8:41:49 AM (16 years ago)

Author:

eugene

Message:

updates from trunk

Location:

branches/tap_branches

Files:

• . (modified) (1 prop)
• psModules (modified) (1 prop)
• psModules/src/objects/pmPSFtryModel.c (modified) (4 diffs)

branches/tap_branches/psModules/src/objects/pmPSFtryModel.c

@@ -72 +72 @@
     }
 
-    // XXX set the min number of needed source more carefully (depends on psfTrendMode?)
-    int orderMax = PS_MAX (options->psfTrendNx, options->psfTrendNy);
-    if ((sources->n < 15) && (orderMax >= 3)) orderMax = 2;
-    if ((sources->n < 11) && (orderMax >= 2)) orderMax = 1;
-    if ((sources->n <  8) && (orderMax >= 1)) orderMax = 0;
+    // hard limit on minimum number of stars
     if ((sources->n <  3)) {
         psError (PS_ERR_UNKNOWN, true, "failed to determine PSF parameters");
@@ -82 +78 @@
     }
 
-    int orderMin;
+    // this is a bit tricky, because we have two cases (MAP vs POLY), and they have a different
+    // definition for 'order' (order_MAP = order_POLY + 1).  in addition, we have a
+    // user-specified MAX order, which we should respect, regardless of the mode
+
+    // set the max order (0 = constant) which the number of psf stars can support:
+    // rule of thumb: require 3 stars per 'cell' (order+1)^2
+    int MaxOrderForStars = 0;
+    if (sources->n >= 12) MaxOrderForStars = 1; // 4 cells
+    if (sources->n >= 27) MaxOrderForStars = 2; // 9 cells
+    if (sources->n >= 48) MaxOrderForStars = 3; // 16 cells
+    if (sources->n >  75) MaxOrderForStars = 4; // 25 cells
+
+    int orderMax = PS_MAX (options->psfTrendNx, options->psfTrendNy);
+    int orderMin = 0;
     if (options->psfTrendMode == PM_TREND_MAP) {
-        orderMin = 1;
-        orderMax = PS_MAX(orderMax, 1);
-    } else {
-        orderMin = 0;
-    }
+        MaxOrderForStars ++;
+        orderMin ++;
+    }
+    orderMax = PS_MIN (orderMax, MaxOrderForStars);
 
     // save the raw source mask (generated by pmPSFtryFitEXT)
@@ -99 +107 @@
     // as we loop over orders, we need to refer to the initial selection, but we modify the
     // option values to match the current guess: save the max values here:
-    int Nx = options->psfTrendNx;
-    int Ny = options->psfTrendNy;
+    int Nx = (options->psfTrendMode == PM_TREND_MAP) ? options->psfTrendNx : options->psfTrendNx + 1;
+    int Ny = (options->psfTrendMode == PM_TREND_MAP) ? options->psfTrendNy : options->psfTrendNy + 1;
     for (int i = orderMin; i <= orderMax; i++) {
 
@@ -196 +204 @@
                   psStatsGetValue(stats, meanStat), psStatsGetValue(stats, stdevStat));
 
-        psFree(flux);
-        psFree(mask);
-        psFree(chisq);
-
         if (!result) {
             psError(PS_ERR_UNKNOWN, false, "Failed to fit psf->ChiTrend");
             psFree(psfTry);
-            return NULL;
-        }
-    }
+            psFree(flux);
+            psFree(mask);
+            psFree(chisq);
+            return NULL;
+        }
+    }
+
+    psFree(flux);
+    psFree(mask);
+    psFree(chisq);
 
     for (int i = 0; i < psfTry->psf->ChiTrend->nX + 1; i++) {