IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 3514


Ignore:
Timestamp:
Mar 25, 2005, 5:44:31 PM (21 years ago)
Author:
gusciora
Message:

...

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/psModules/src/pmObjects.c

    r3512 r3514  
    55 *  @author GLG, MHPCC
    66 *
    7  *  @version $Revision: 1.9 $ $Name: not supported by cvs2svn $
    8  *  @date $Date: 2005-03-25 23:23:29 $
     7 *  @version $Revision: 1.10 $ $Name: not supported by cvs2svn $
     8 *  @date $Date: 2005-03-26 03:44:31 $
    99 *
    1010 *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
     
    10821082XXX: This should probably be a public function.
    10831083 
    1084 XXX: We don't care about the model derivatives here.  So, for the individual
    1085 model evaluation functions, we should probably separate the derivative
    1086 calculation.
    1087  
    1088 XXX: Get rid of deriv.
    1089  
    10901084XXX: Use static vectors for x.
    10911085 
     
    10961090                psU32 col)
    10971091{
    1098     // XXX: You won't need this once you separate the function evaluation.
    1099     // from the derivative calculation.
    1100     psVector *deriv = psVectorAlloc(src->models->Nparams, PS_TYPE_F32);
    11011092    // XXX: The following step will not be necessary if the models->params
    11021093    // member is a psVector.
     
    11131104    switch (src->models->type) {
    11141105    case PS_MODEL_GAUSS:
    1115         tmpF = pmMinLM_Gauss2D(deriv, params, x);
     1106        tmpF = pmMinLM_Gauss2D(NULL, params, x);
    11161107        break;
    11171108    case PS_MODEL_PGAUSS:
    1118         tmpF = pmMinLM_PsuedoGauss2D(deriv, params, x);
     1109        tmpF = pmMinLM_PsuedoGauss2D(NULL, params, x);
    11191110        break;
    11201111    case PS_MODEL_TWIST_GAUSS:
    1121         tmpF = pmMinLM_TwistGauss2D(deriv, params, x);
     1112        tmpF = pmMinLM_TwistGauss2D(NULL, params, x);
    11221113        break;
    11231114    case PS_MODEL_WAUSS:
    1124         tmpF = pmMinLM_Wauss2D(deriv, params, x);
     1115        tmpF = pmMinLM_Wauss2D(NULL, params, x);
    11251116        break;
    11261117    case PS_MODEL_SERSIC:
    1127         tmpF = pmMinLM_Sersic(deriv, params, x);
     1118        tmpF = pmMinLM_Sersic(NULL, params, x);
    11281119        break;
    11291120    case PS_MODEL_SERSIC_CORE:
    1130         tmpF = pmMinLM_SersicCore(deriv, params, x);
     1121        tmpF = pmMinLM_SersicCore(NULL, params, x);
    11311122        break;
    11321123    default:
     
    11351126    }
    11361127
    1137     psFree(deriv);
    11381128    psFree(params);
    11391129    psFree(x);
     
    13881378    #define MAX_PARAMS 100
    13891379
    1390     psVector *deriv = psVectorAlloc(MAX_PARAMS, PS_TYPE_F32);
    13911380    psVector *params = psVectorAlloc(MAX_PARAMS, PS_TYPE_F32);
    13921381    psVector *x = psVectorAlloc(2, PS_TYPE_F32);
     
    14101399            switch (src->models->type) {
    14111400            case PS_MODEL_GAUSS:
    1412                 pixelValue = pmMinLM_Gauss2D(deriv, params, x);
     1401                pixelValue = pmMinLM_Gauss2D(NULL, params, x);
    14131402                break;
    14141403            case PS_MODEL_PGAUSS:
    1415                 pixelValue = pmMinLM_PsuedoGauss2D(deriv, params, x);
     1404                pixelValue = pmMinLM_PsuedoGauss2D(NULL, params, x);
    14161405                break;
    14171406            case PS_MODEL_TWIST_GAUSS:
    1418                 pixelValue = pmMinLM_TwistGauss2D(deriv, params, x);
     1407                pixelValue = pmMinLM_TwistGauss2D(NULL, params, x);
    14191408                break;
    14201409            case PS_MODEL_WAUSS:
    1421                 pixelValue = pmMinLM_Wauss2D(deriv, params, x);
     1410                pixelValue = pmMinLM_Wauss2D(NULL, params, x);
    14221411                break;
    14231412            case PS_MODEL_SERSIC:
    1424                 pixelValue = pmMinLM_Sersic(deriv, params, x);
     1413                pixelValue = pmMinLM_Sersic(NULL, params, x);
    14251414                break;
    14261415            case PS_MODEL_SERSIC_CORE:
    1427                 pixelValue = pmMinLM_SersicCore(deriv, params, x);
     1416                pixelValue = pmMinLM_SersicCore(NULL, params, x);
    14281417                break;
    14291418            default:
    14301419                psError(PS_ERR_UNKNOWN, true, "Undefined psModelType");
    14311420                psFree(x);
    1432                 psFree(deriv);
    14331421                psFree(params);
    14341422                return(false);
     
    14471435    }
    14481436    psFree(x);
    1449     psFree(deriv);
    14501437    psFree(params);
    14511438    return(true);
     
    14991486XXX: Consider getting rid of the parameter checks since this might consume
    15001487a significant fraction of this function CPU time.
     1488 
     1489XXX: I added the following.  Must conform with IfA.  If deriv==NULL, then
     1490we simply don't perform the derivative calculations.
    15011491 *****************************************************************************/
    15021492psF32 pmMinLM_Gauss2D(psVector *deriv,
     
    15041494                      psVector *x)
    15051495{
    1506     PS_VECTOR_CHECK_NULL(deriv, NAN);
    1507     PS_VECTOR_CHECK_TYPE(deriv, PS_TYPE_F32, NAN);
    1508     PS_VECTOR_CHECK_SIZE(deriv, 7, NAN);
    15091496    PS_VECTOR_CHECK_NULL(params, NAN);
    15101497    PS_VECTOR_CHECK_TYPE(params, PS_TYPE_F32, NAN);
     
    15221509    psF32 f = params->data.F32[1]*r + params->data.F32[0];
    15231510
    1524     psF32 q = params->data.F32[1]*r;
    1525     deriv->data.F32[0] = +1.0;
    1526     deriv->data.F32[1] = +r;
    1527     deriv->data.F32[2] = q*(2*px*params->data.F32[4] + params->data.F32[6]*Y);
    1528     deriv->data.F32[3] = q*(2*py*params->data.F32[5] + params->data.F32[6]*X);
    1529     deriv->data.F32[4] = -2.0*q*px*X;
    1530     deriv->data.F32[5] = -2.0*q*py*Y;
    1531     deriv->data.F32[6] = -q*X*Y;
     1511    if (deriv != NULL) {
     1512        PS_VECTOR_CHECK_TYPE(deriv, PS_TYPE_F32, NAN);
     1513        PS_VECTOR_CHECK_SIZE(deriv, 7, NAN);
     1514
     1515        psF32 q = params->data.F32[1]*r;
     1516        deriv->data.F32[0] = +1.0;
     1517        deriv->data.F32[1] = +r;
     1518        deriv->data.F32[2] = q*(2*px*params->data.F32[4] + params->data.F32[6]*Y);
     1519        deriv->data.F32[3] = q*(2*py*params->data.F32[5] + params->data.F32[6]*X);
     1520        deriv->data.F32[4] = -2.0*q*px*X;
     1521        deriv->data.F32[5] = -2.0*q*py*Y;
     1522        deriv->data.F32[6] = -q*X*Y;
     1523    }
    15321524
    15331525    return(f);
     
    15881580                            psVector *x)
    15891581{
    1590     PS_VECTOR_CHECK_NULL(deriv, NAN);
    1591     PS_VECTOR_CHECK_TYPE(deriv, PS_TYPE_F32, NAN);
    1592     PS_VECTOR_CHECK_SIZE(deriv, 7, NAN);
    15931582    PS_VECTOR_CHECK_NULL(params, NAN);
    15941583    PS_VECTOR_CHECK_TYPE(params, PS_TYPE_F32, NAN);
     
    16071596    psF32 f = params->data.F32[1]*r + params->data.F32[0];
    16081597
    1609     //
    1610     // note difference from a pure gaussian: q = params->data.F32[1]*r
    1611     //
    1612 
    1613     psF32 q = params->data.F32[1]*r*r*t;
    1614     deriv->data.F32[0] = +1.0;
    1615     deriv->data.F32[1] = +r;
    1616     deriv->data.F32[2] = q*(2.0*px*params->data.F32[4] + params->data.F32[6]*Y);
    1617     deriv->data.F32[3] = q *
    1618                          (2.0*py*params->data.F32[5] + params->data.F32[6]*X);
    1619     deriv->data.F32[4] = -2.0*q*px*X;
    1620     deriv->data.F32[5] = -2.0*q*py*Y;
    1621     deriv->data.F32[6] = -q*X*Y;
     1598    if (deriv != NULL) {
     1599        PS_VECTOR_CHECK_TYPE(deriv, PS_TYPE_F32, NAN);
     1600        PS_VECTOR_CHECK_SIZE(deriv, 7, NAN);
     1601
     1602        //
     1603        // note difference from a pure gaussian: q = params->data.F32[1]*r
     1604        //
     1605        psF32 q = params->data.F32[1]*r*r*t;
     1606        deriv->data.F32[0] = +1.0;
     1607        deriv->data.F32[1] = +r;
     1608        deriv->data.F32[2] = q*(2.0*px*params->data.F32[4] + params->data.F32[6]*Y);
     1609        deriv->data.F32[3] = q *
     1610                             (2.0*py*params->data.F32[5] + params->data.F32[6]*X);
     1611        deriv->data.F32[4] = -2.0*q*px*X;
     1612        deriv->data.F32[5] = -2.0*q*py*Y;
     1613        deriv->data.F32[6] = -q*X*Y;
     1614    }
    16221615
    16231616    return(f);
     
    16821675                      psVector *x)
    16831676{
    1684     PS_VECTOR_CHECK_NULL(deriv, NAN);
    1685     PS_VECTOR_CHECK_TYPE(deriv, PS_TYPE_F32, NAN);
    1686     PS_VECTOR_CHECK_SIZE(deriv, 9, NAN);
    16871677    PS_VECTOR_CHECK_NULL(params, NAN);
    16881678    PS_VECTOR_CHECK_TYPE(params, PS_TYPE_F32, NAN);
     
    17011691    psF32 f = params->data.F32[1]*r + params->data.F32[0];
    17021692
    1703     //
    1704     // note difference from gaussian: q = params->data.F32[1]*r
    1705     //
    1706     psF32 q = params->data.F32[1]*r*r*(1.0 + params->data.F32[7]*z*(1.0 + params->data.F32[8]*z/2.0));
    1707     deriv->data.F32[0] = +1.0;
    1708     deriv->data.F32[1] = +r;
    1709     deriv->data.F32[2] = q*(2.0*px*params->data.F32[4] + params->data.F32[6]*Y);
    1710     deriv->data.F32[3] = q*(2.0*py*params->data.F32[5] + params->data.F32[6]*X);
    1711     deriv->data.F32[4] = -2.0*q*px*X;
    1712     deriv->data.F32[5] = -2.0*q*py*Y;
    1713     deriv->data.F32[6] = -q*X*Y;
    1714     deriv->data.F32[7] = -100.0*params->data.F32[1]*r*r*t;
    1715     deriv->data.F32[8] = -100.0*params->data.F32[1]*r*r*params->data.F32[7]*(z*z*z)/6.0;
    1716     //
    1717     // The values of 100 dampen the swing of params->data.F32[7,8] */
    1718     //
     1693    if (deriv != NULL) {
     1694        PS_VECTOR_CHECK_TYPE(deriv, PS_TYPE_F32, NAN);
     1695        PS_VECTOR_CHECK_SIZE(deriv, 9, NAN);
     1696        //
     1697        // note difference from gaussian: q = params->data.F32[1]*r
     1698        //
     1699
     1700        psF32 q = params->data.F32[1]*r*r*(1.0 + params->data.F32[7]*z*(1.0 + params->data.F32[8]*z/2.0));
     1701        deriv->data.F32[0] = +1.0;
     1702        deriv->data.F32[1] = +r;
     1703        deriv->data.F32[2] = q*(2.0*px*params->data.F32[4] + params->data.F32[6]*Y);
     1704        deriv->data.F32[3] = q*(2.0*py*params->data.F32[5] + params->data.F32[6]*X);
     1705        deriv->data.F32[4] = -2.0*q*px*X;
     1706        deriv->data.F32[5] = -2.0*q*py*Y;
     1707        deriv->data.F32[6] = -q*X*Y;
     1708        deriv->data.F32[7] = -100.0*params->data.F32[1]*r*r*t;
     1709        deriv->data.F32[8] = -100.0*params->data.F32[1]*r*r*params->data.F32[7]*(z*z*z)/6.0;
     1710        //
     1711        // The values of 100 dampen the swing of params->data.F32[7,8] */
     1712        //
     1713    }
    17191714
    17201715    return(f);
     
    17861781                           psVector *x)
    17871782{
    1788     PS_VECTOR_CHECK_NULL(deriv, NAN);
    1789     PS_VECTOR_CHECK_TYPE(deriv, PS_TYPE_F32, NAN);
    1790     PS_VECTOR_CHECK_SIZE(deriv, 11, NAN);
    17911783    PS_VECTOR_CHECK_NULL(params, NAN);
    17921784    PS_VECTOR_CHECK_TYPE(params, PS_TYPE_F32, NAN);
     
    18081800
    18091801    psF32 f = params->data.F32[5]*r + params->data.F32[6];
    1810     psF32 q1 = params->data.F32[5]*PS_SQR(r);
    1811     psF32 q2 = params->data.F32[5]*PS_SQR(r)*params->data.F32[10]*pow(z2,(params->data.F32[10]-1.0));
    1812     deriv->data.F32[0] = +1.0;
    1813     deriv->data.F32[1] = +r;
    1814     deriv->data.F32[2] = q1*(2.0*px1*params->data.F32[4] + params->data.F32[6]*Y) + q2*(2*px2*params->data.F32[7] + params->data.F32[9]*Y);
    1815     deriv->data.F32[3] = q1*(2.0*py1*params->data.F32[5] + params->data.F32[6]*X) + q2*(2*py2*params->data.F32[8] + params->data.F32[9]*X);
    1816 
    1817     //
    1818     // These fudge factors impede the growth of params->data.F32[4] beyond
    1819     // params->data.F32[7].
    1820     //
    1821     psF32 f1 = fabs(params->data.F32[7]) / fabs(params->data.F32[4]);
    1822     psF32 f2 = (f1 < FSCALE) ? 1.0 : FFACTOR*(f1 - FSCALE) + 1.0;
    1823     deriv->data.F32[4] = -2.0*q1*px1*X*f2;
    1824 
    1825     //
    1826     // These fudge factors impede the growth of params->data.F32[5] beyond
    1827     // params->data.F32[8].
    1828     //
    1829     f1 = fabs(params->data.F32[8]) / fabs(params->data.F32[5]);
    1830     f2 = (f1 < FSCALE) ? 1.0 : FFACTOR*(f1 - FSCALE) + 1.0;
    1831     deriv->data.F32[5] = -2.0*q1*py1*Y*f2;
    1832     deriv->data.F32[6] = -q1*X*Y;
    1833     deriv->data.F32[7] = -2.0*q2*px2*X;
    1834     deriv->data.F32[8] = -2.0*q2*py2*Y;
    1835     deriv->data.F32[9] = -q2*X*Y;
    1836     deriv->data.F32[10] = -q1*log(z2);
     1802
     1803    if (deriv != NULL) {
     1804        PS_VECTOR_CHECK_TYPE(deriv, PS_TYPE_F32, NAN);
     1805        PS_VECTOR_CHECK_SIZE(deriv, 11, NAN);
     1806
     1807        psF32 q1 = params->data.F32[5]*PS_SQR(r);
     1808        psF32 q2 = params->data.F32[5]*PS_SQR(r)*params->data.F32[10]*pow(z2,(params->data.F32[10]-1.0));
     1809        deriv->data.F32[0] = +1.0;
     1810        deriv->data.F32[1] = +r;
     1811        deriv->data.F32[2] = q1*(2.0*px1*params->data.F32[4] + params->data.F32[6]*Y) + q2*(2*px2*params->data.F32[7] + params->data.F32[9]*Y);
     1812        deriv->data.F32[3] = q1*(2.0*py1*params->data.F32[5] + params->data.F32[6]*X) + q2*(2*py2*params->data.F32[8] + params->data.F32[9]*X);
     1813
     1814        //
     1815        // These fudge factors impede the growth of params->data.F32[4] beyond
     1816        // params->data.F32[7].
     1817        //
     1818        psF32 f1 = fabs(params->data.F32[7]) / fabs(params->data.F32[4]);
     1819        psF32 f2 = (f1 < FSCALE) ? 1.0 : FFACTOR*(f1 - FSCALE) + 1.0;
     1820        deriv->data.F32[4] = -2.0*q1*px1*X*f2;
     1821
     1822        //
     1823        // These fudge factors impede the growth of params->data.F32[5] beyond
     1824        // params->data.F32[8].
     1825        //
     1826        f1 = fabs(params->data.F32[8]) / fabs(params->data.F32[5]);
     1827        f2 = (f1 < FSCALE) ? 1.0 : FFACTOR*(f1 - FSCALE) + 1.0;
     1828        deriv->data.F32[5] = -2.0*q1*py1*Y*f2;
     1829        deriv->data.F32[6] = -q1*X*Y;
     1830        deriv->data.F32[7] = -2.0*q2*px2*X;
     1831        deriv->data.F32[8] = -2.0*q2*py2*Y;
     1832        deriv->data.F32[9] = -q2*X*Y;
     1833        deriv->data.F32[10] = -q1*log(z2);
     1834    }
    18371835
    18381836    return(f);
     
    18961894                     psVector *x)
    18971895{
    1898     PS_VECTOR_CHECK_NULL(deriv, NAN);
    1899     PS_VECTOR_CHECK_TYPE(deriv, PS_TYPE_F32, NAN);
    1900     PS_VECTOR_CHECK_SIZE(deriv, 8, NAN);
    19011896    PS_VECTOR_CHECK_NULL(params, NAN);
    19021897    PS_VECTOR_CHECK_TYPE(params, PS_TYPE_F32, NAN);
     
    19051900    PS_VECTOR_CHECK_TYPE(x, PS_TYPE_F32, NAN);
    19061901    PS_VECTOR_CHECK_SIZE(x, 2, NAN);
     1902
     1903    if (deriv != NULL) {
     1904        PS_VECTOR_CHECK_TYPE(deriv, PS_TYPE_F32, NAN);
     1905        PS_VECTOR_CHECK_SIZE(deriv, 8, NAN);
     1906    }
    19071907
    19081908    psError(PS_ERR_UNKNOWN, true, "This function is not implemented yet.");
     
    19681968                         psVector *x)
    19691969{
    1970     PS_VECTOR_CHECK_NULL(deriv, NAN);
    1971     PS_VECTOR_CHECK_TYPE(deriv, PS_TYPE_F32, NAN);
    1972     PS_VECTOR_CHECK_SIZE(deriv, 12, NAN);
    19731970    PS_VECTOR_CHECK_NULL(params, NAN);
    19741971    PS_VECTOR_CHECK_TYPE(params, PS_TYPE_F32, NAN);
     
    19771974    PS_VECTOR_CHECK_TYPE(x, PS_TYPE_F32, NAN);
    19781975    PS_VECTOR_CHECK_SIZE(x, 2, NAN);
     1976
     1977    if (deriv != NULL) {
     1978        PS_VECTOR_CHECK_TYPE(deriv, PS_TYPE_F32, NAN);
     1979        PS_VECTOR_CHECK_SIZE(deriv, 12, NAN);
     1980    }
    19791981
    19801982    psError(PS_ERR_UNKNOWN, true, "This function is not implemented yet.");
Note: See TracChangeset for help on using the changeset viewer.