IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Ignore:
Timestamp:
Oct 2, 2009, 12:11:34 PM (17 years ago)
Author:
eugene
Message:

merge changes from trunk

File:
1 edited

Legend:

Unmodified
Added
Removed

  • branches/eam_branches/20090715/psModules/src/objects/models/pmModel_PGAUSS.c

    r25683 r25745  
    1919 *****************************************************************************/
    2020
     21#include <stdio.h>
     22#include <pslib.h>
     23
     24#include "pmMoments.h"
     25#include "pmPeaks.h"
     26#include "pmSource.h"
     27#include "pmModel.h"
     28#include "pmModel_PGAUSS.h"
     29
    2130# define PM_MODEL_FUNC            pmModelFunc_PGAUSS
    2231# define PM_MODEL_FLUX            pmModelFlux_PGAUSS
    … …  
    2736# define PM_MODEL_PARAMS_FROM_PSF pmModelParamsFromPSF_PGAUSS
    2837# define PM_MODEL_FIT_STATUS      pmModelFitStatus_PGAUSS
     38# define PM_MODEL_SET_LIMITS      pmModelSetLimits_PGAUSS
     39
     40// Lax parameter limits
     41static float paramsMinLax[] = { -1.0e3, 1.0e-2, -100, -100, 0.5, 0.5, -1.0 };
     42static float paramsMaxLax[] = { 1.0e5, 1.0e8, 1.0e4, 1.0e4, 100, 100, 1.0 };
     43
     44// Strict parameter limits
     45static float *paramsMinStrict = paramsMinLax;
     46static float *paramsMaxStrict = paramsMaxLax;
     47
     48// Parameter limits to use
     49static float *paramsMinUse = paramsMinLax;
     50static float *paramsMaxUse = paramsMaxLax;
     51static float betaUse[] = { 1000, 3e6, 5, 5, 2.0, 2.0, 0.5 };
     52
     53static bool limitsApply = true;         // Apply limits?
    2954
    3055// the model is a function of the pixel coordinate (pixcoord[0,1] = x,y)
    … …  
    7196bool PM_MODEL_LIMITS (psMinConstraintMode mode, int nParam, float *params, float *beta)
    7297{
    73     float beta_lim = 0, params_min = 0, params_max = 0;
    74     float f1 = 0, f2 = 0, q1 = 0, q2 = 0;
     98    if (!limitsApply) {
     99        return true;
     100    }
     101    psAssert(nParam >= 0 && nParam <= PM_PAR_7, "Parameter index is out of bounds");
    75102
    76103    // we need to calculate the limits for SXY specially
     104    float q2 = NAN;
    77105    if (nParam == PM_PAR_SXY) {
    78         f1 = 1.0 / PS_SQR(params[PM_PAR_SYY]) + 1.0 / PS_SQR(params[PM_PAR_SXX]);
    79         f2 = 1.0 / PS_SQR(params[PM_PAR_SYY]) - 1.0 / PS_SQR(params[PM_PAR_SXX]);
    80         q1 = PS_SQR(f1)*AR_RATIO - PS_SQR(f2);
    81         q1 = PS_MAX (0.0, q1);
     106        float f1 = 1.0 / PS_SQR(params[PM_PAR_SYY]) + 1.0 / PS_SQR(params[PM_PAR_SXX]);
     107        float f2 = 1.0 / PS_SQR(params[PM_PAR_SYY]) - 1.0 / PS_SQR(params[PM_PAR_SXX]);
     108        float q1 = PS_SQR(f1)*AR_RATIO - PS_SQR(f2);
     109        q1 = (q1 < 0.0) ? 0.0 : q1;
    82110        // if q1 < 0.0, f2 ~ f1, we have a very large axis ratio near 45deg..  Saturate at that
    83111        // angle and let f2,f1 fight it out
    84         q2  = 0.5*sqrt (q1);
     112        q2 = 0.5*sqrtf(q1);
    85113    }
    86114
    87115    switch (mode) {
    88     case PS_MINIMIZE_BETA_LIMIT:
    89         switch (nParam) {
    90         case PM_PAR_SKY:
    91             beta_lim = 1000;
    92             break;
    93         case PM_PAR_I0:
    94             beta_lim = 3e6;
    95             break;
    96         case PM_PAR_XPOS:
    97             beta_lim = 5;
    98             break;
    99         case PM_PAR_YPOS:
    100             beta_lim = 5;
    101             break;
    102         case PM_PAR_SXX:
    103             beta_lim = 2.0;
    104             break;
    105         case PM_PAR_SYY:
    106             beta_lim = 2.0;
    107             break;
    108         case PM_PAR_SXY:
    109             beta_lim =  0.5*q2;
    110             break;
    111         default:
    112             psAbort("invalid parameter %d for beta test", nParam);
    113         }
    114         if (fabs(beta[nParam]) > fabs(beta_lim)) {
    115             beta[nParam] = (beta[nParam] > 0) ? fabs(beta_lim) : -fabs(beta_lim);
    116             psTrace ("psModules.objects", 5, "|beta[nParam==%d]| > |beta_lim|; %g v. %g",
    117                      nParam, beta[nParam], beta_lim);
    118             return false;
    119         }
    120         return true;
    121     case PS_MINIMIZE_PARAM_MIN:
    122         switch (nParam) {
    123         case PM_PAR_SKY:
    124             params_min = -1000;
    125             break;
    126         case PM_PAR_I0:
    127             params_min =  0.01;
    128             break;
    129         case PM_PAR_XPOS:
    130             params_min =  -100;
    131             break;
    132         case PM_PAR_YPOS:
    133             params_min =  -100;
    134             break;
    135         case PM_PAR_SXX:
    136             params_min =   0.5;
    137             break;
    138         case PM_PAR_SYY:
    139             params_min =   0.5;
    140             break;
    141         case PM_PAR_SXY:
    142             params_min =   -q2;
    143             break;
    144         default:
    145             psAbort("invalid parameter %d for param min test", nParam);
    146         }
    147         if (params[nParam] < params_min) {
    148             params[nParam] = params_min;
    149             psTrace ("psModules.objects", 5, "params[nParam==%d] < params_min; %g v. %g",
    150                      nParam, params[nParam], params_min);
    151             return false;
    152         }
    153         return true;
    154     case PS_MINIMIZE_PARAM_MAX:
    155         switch (nParam) {
    156         case PM_PAR_SKY:
    157             params_max =   1e5;
    158             break;
    159         case PM_PAR_I0:
    160             params_max =   1e8;
    161             break;
    162         case PM_PAR_XPOS:
    163             params_max =   1e4;
    164             break;
    165         case PM_PAR_YPOS:
    166             params_max =   1e4;
    167             break;
    168         case PM_PAR_SXX:
    169             params_max =   100;
    170             break;
    171         case PM_PAR_SYY:
    172             params_max =   100;
    173             break;
    174         case PM_PAR_SXY:
    175             params_max =   +q2;
    176             break;
    177         default:
    178             psAbort("invalid parameter %d for param max test", nParam);
    179         }
    180         if (params[nParam] > params_max) {
    181             params[nParam] = params_max;
    182             psTrace ("psModules.objects", 5, "params[nParam==%d] > params_max; %g v. %g",
    183                      nParam, params[nParam], params_max);
    184             return false;
    185         }
    186         return true;
    187     default:
     116      case PS_MINIMIZE_BETA_LIMIT: {
     117          psAssert(beta, "Require beta to limit beta");
     118          float limit = betaUse[nParam];
     119          if (nParam == PM_PAR_SXY) {
     120              limit *= q2;
     121          }
     122          if (fabs(beta[nParam]) > fabs(limit)) {
     123              beta[nParam] = (beta[nParam] > 0) ? fabs(limit) : -fabs(limit);
     124              psTrace("psModules.objects", 5, "|beta[nParam==%d]| > |beta_lim|; %g v. %g",
     125                      nParam, beta[nParam], limit);
     126              return false;
     127          }
     128          return true;
     129      }
     130      case PS_MINIMIZE_PARAM_MIN: {
     131          psAssert(params, "Require parameters to limit parameters");
     132          psAssert(paramsMinUse, "Require parameter limits to limit parameters");
     133          float limit = paramsMinUse[nParam];
     134          if (nParam == PM_PAR_SXY) {
     135              limit *= q2;
     136          }
     137          if (params[nParam] < limit) {
     138              params[nParam] = limit;
     139              psTrace("psModules.objects", 5, "params[nParam==%d] < params_min; %g v. %g",
     140                      nParam, params[nParam], limit);
     141              return false;
     142          }
     143          return true;
     144      }
     145      case PS_MINIMIZE_PARAM_MAX: {
     146          psAssert(params, "Require parameters to limit parameters");
     147          psAssert(paramsMaxUse, "Require parameter limits to limit parameters");
     148          float limit = paramsMaxUse[nParam];
     149          if (nParam == PM_PAR_SXY) {
     150              limit *= q2;
     151          }
     152          if (params[nParam] > limit) {
     153              params[nParam] = limit;
     154              psTrace("psModules.objects", 5, "params[nParam==%d] > params_max; %g v. %g",
     155                      nParam, params[nParam], limit);
     156              return false;
     157          }
     158          return true;
     159      }
     160      default:
    188161        psAbort("invalid choice for limits");
    189162    }
    … …  
    191164    return false;
    192165}
     166
    193167
    194168// make an initial guess for parameters
    … …  
    432406}
    433407
     408
     409void PM_MODEL_SET_LIMITS(pmModelLimitsType type)
     410{
     411    switch (type) {
     412      case PM_MODEL_LIMITS_NONE:
     413        paramsMinUse = NULL;
     414        paramsMaxUse = NULL;
     415        limitsApply = true;
     416        break;
     417      case PM_MODEL_LIMITS_IGNORE:
     418        paramsMinUse = NULL;
     419        paramsMaxUse = NULL;
     420        limitsApply = false;
     421        break;
     422      case PM_MODEL_LIMITS_LAX:
     423        paramsMinUse = paramsMinLax;
     424        paramsMaxUse = paramsMaxLax;
     425        limitsApply = true;
     426        break;
     427      case PM_MODEL_LIMITS_STRICT:
     428        paramsMinUse = paramsMinStrict;
     429        paramsMaxUse = paramsMaxStrict;
     430        limitsApply = true;
     431        break;
     432      default:
     433        psAbort("Unrecognised model limits type: %x", type);
     434    }
     435    return;
     436}
     437
    434438# undef PM_MODEL_FUNC
    435439# undef PM_MODEL_FLUX
    … …  
    440444# undef PM_MODEL_PARAMS_FROM_PSF
    441445# undef PM_MODEL_FIT_STATUS
     446# undef PM_MODEL_SET_LIMITS
Note: See TracChangeset for help on using the changeset viewer.