IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 2650


Ignore:
Timestamp:
Dec 7, 2004, 12:28:29 PM (22 years ago)
Author:
gusciora
Message:

...

Location:
trunk/psLib/src
Files:
6 edited

Legend:

Unmodified
Added
Removed

  • trunk/psLib/src/astro/psCoord.c

    r2635 r2650  
    1010*  @author GLG, MHPCC
    1111*
    12 *  @version $Revision: 1.39 $ $Name: not supported by cvs2svn $
    13 *  @date $Date: 2004-12-06 19:49:36 $
     12*  @version $Revision: 1.40 $ $Name: not supported by cvs2svn $
     13*  @date $Date: 2004-12-07 22:28:28 $
    1414*
    1515*  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
    … …  
    189189
    190190/******************************************************************************
    191 XXX: It's not clear if this code correctly implements what the equations in
    192 the ADD imply.  Bug 243 will, hopefully, decide that.
    193  
    194 XXX: Once this code is verified, delete Gene's code.
    195191 *****************************************************************************/
    196192psSphere* psSphereTransformApply(psSphere* out,
    … …  
    224220
    225221    return(out);
    226 
    227     /*
    228     This algorithm comes from an email from Gene.  I assume (x,y) corresponds to
    229     (r,d) in the sphere coordinates.
    230      
    231     XXX: In Gene's email, there are different variables with similar names (y
    232     and Y) and in the code, there's cos_y, cos_Y, and cos(y).  Verify that there
    233     are no typo's.
    234      
    235         psF64 sinY = 0.0;
    236         psF64 cosY = 0.0;
    237         psF64 sinX = 0.0;
    238         psF64 cosX = 0.0;
    239         psF64 x = 0.0;
    240         psF64 y = 0.0;
    241         psF64 dx = 0.0;
    242      
    243         x = coord->r;
    244         y = coord->d;
    245         dx = x - transform->phiP;
    246         sinY = cos(y) * sin(dx) * transform->sinDeltaP +
    247                sin(y) * transform->cosDeltaP;
    248         cosY = sqrt(1.0 - sinY * sinY);
    249         sinX = (cos(y) * sin(dx) * transform->cosDeltaP -
    250                 sin(y) * transform->sinDeltaP) / cos(y);
    251         cosX = cos(y) * cos(dx) / cos(y);
    252      
    253         out->r = atan2(sinX, cosX) + transform->alphaP;
    254         out->d = atan2(sinY, cosY);
    255     */
    256 }
    257 
    258 // XXX: the ADD changed.  This code no longer conforms to the ADD.  In fact,
    259 // there is no algorithm for this in the ADD.
    260 // XXX: This is bug 244
    261 psSphereTransform* psSphereTransformICRSToEcliptic(psTime time)
    262 {
    263     PS_PTR_CHECK_NULL(&time, NULL);
    264 
    265     struct tm *tmTime = psTimeToTM(&time);
    266     psF64 year = (psF64)(1900 + tmTime->tm_year);
    267     psF64 T = year / 100.0;
    268     psF64 phi = DEG_TO_RAD(-23.452294) +
    269                 (DEG_TO_RAD(0.013013) * T) +
    270                 (DEG_TO_RAD(0.000001639) * T * T) -
    271                 (DEG_TO_RAD(0.000000503) * T * T * T);
    272     psF64 alphaP = 0.0;
    273     psF64 phiP = 0.0;
    274 
    275     psFree(tmTime);
    276     return (psSphereTransformAlloc(phi, alphaP, phiP));
    277 }
    278 
    279 // XXX: the ADD changed.  This code no longer conforms to the old ADD.
    280 // XXX: Waiting for IfA to sepcify if T should include days/months/etc.
    281 psSphereTransform* psSphereTransformEclipticToICRS(psTime time)
    282 {
    283     PS_PTR_CHECK_NULL(&time, NULL);
    284 
    285     struct tm *tmTime = psTimeToTM(&time);
    286     psF64 year = (psF64)(1900 + tmTime->tm_year);
     222}
     223
     224psSphereTransform* psSphereTransformICRSToEcliptic(psTime *time)
     225{
     226    PS_PTR_CHECK_NULL(time, NULL);
     227    psF64 MJD = psTimeToMJD(time);
     228    struct tm *tmTime = psTimeToTM(time);
     229    psF64 year = (psF64)(1900 - MJD + tmTime->tm_year);
    287230    psF64 T = year / 100.0;
    288231    psF64 alphaP = 0.0;
    … …  
    296239
    297240    psFree(tmTime);
     241
     242    // Don't neglect the minus sign on deltaP (bug 244):
     243    return (psSphereTransformAlloc(alphaP, -deltaP, phiP));
     244}
     245
     246
     247psSphereTransform* psSphereTransformEclipticToICRS(psTime *time)
     248{
     249    PS_PTR_CHECK_NULL(time, NULL);
     250    psF64 MJD = psTimeToMJD(time);
     251    struct tm *tmTime = psTimeToTM(time);
     252    psF64 year = (psF64)(1900 - MJD + tmTime->tm_year);
     253    psF64 T = year / 100.0;
     254    psF64 alphaP = 0.0;
     255    psF64 deltaP = DEG_TO_RAD(23.0) +
     256                   SEC_TO_RAD(27.0) +
     257                   MIN_TO_RAD(8.0) -
     258                   (SEC_TO_RAD(46.845) * T) -
     259                   (SEC_TO_RAD(0.0059) * T * T) +
     260                   (SEC_TO_RAD(0.00181) * T * T * T);
     261    psF64 phiP = 0.0;
     262
     263    psFree(tmTime);
    298264    return (psSphereTransformAlloc(alphaP, deltaP, phiP));
    299265}
    300266
    301 // XXX: Should these be in radians, not degrees?
    302 // XXX: the ADD changed.  This code no longer conforms to the ADD.  In fact,
    303 // there is no algorithm for this in the ADD.
    304 // XXX: This is bug 245
     267// XXX: This is bug 245: alphaP swaps with phiP from psSphereTransformGalaticToICRS()
    305268psSphereTransform* psSphereTransformICRSToGalatic(void)
    306269{
    307     psF64 alphaP = DEG_TO_RAD(-62.6);
    308     psF64 deltaP = DEG_TO_RAD(33.0);
    309     psF64 phiP = DEG_TO_RAD(282.25);
     270    psF64 alphaP = DEG_TO_RAD(32.93192);
     271    psF64 deltaP = DEG_TO_RAD(-62.87175);
     272    psF64 phiP = DEG_TO_RAD(282.85948);
    310273
    311274    return (psSphereTransformAlloc(alphaP, deltaP, phiP));
    312275}
    313276
    314 // XXX: Should these be in radians, not degrees?
    315 // XXX: the ADD changed.  This code no longer conforms to the old ADD.
    316277psSphereTransform* psSphereTransformGalaticToICRS(void)
    317278{
    … …  
    588549
    589550
     551
    590552/******************************************************************************
    591553psSpherePrecess(coords, fromTime, toTime):
    … …  
    599561                          const psTime *toTime)
    600562{
    601     psF64 fromMJD = fromTime->sec/86400.0 +
    602                     fromTime->usec/86400000000.0 +
    603                     40587.0;
    604     psF64 toMJD = toTime->sec/86400.0 +
    605                   toTime->usec/86400000000.0 +
    606                   40587.0;
     563    psF64 fromMJD = psTimeToMJD(fromTime);
     564    psF64 toMJD = psTimeToMJD(toTime);
    607565    psF64 T = (toMJD - fromMJD) / 36525.0;
    608566

  • trunk/psLib/src/astro/psCoord.h

    r2635 r2650  
    1 /** @file  psCoord.h
     1a/** @file  psCoord.h
    22*
    33*  @brief Contains basic coordinate transformation definitions and operations
    … …  
    1010*  @author GLG, MHPCC
    1111*
    12 *  @version $Revision: 1.21 $ $Name: not supported by cvs2svn $
    13 *  @date $Date: 2004-12-06 19:49:36 $
     12*  @version $Revision: 1.22 $ $Name: not supported by cvs2svn $
     13*  @date $Date: 2004-12-07 22:28:28 $
    1414*
    1515*  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
    … …  
    239239 */
    240240psSphereTransform* psSphereTransformICRSToEcliptic(
    241     psTime time                        ///< the time for which the resulting transform will be valid
     241    psTime *time                        ///< the time for which the resulting transform will be valid
    242242);
    243243
    … …  
    248248 */
    249249psSphereTransform* psSphereTransformEclipticToICRS(
    250     psTime time                        ///< the time for which the resulting transform will be valid
     250    psTime *time                        ///< the time for which the resulting transform will be valid
    251251);
    252252

  • trunk/psLib/src/astronomy/psCoord.c

    r2635 r2650  
    1010*  @author GLG, MHPCC
    1111*
    12 *  @version $Revision: 1.39 $ $Name: not supported by cvs2svn $
    13 *  @date $Date: 2004-12-06 19:49:36 $
     12*  @version $Revision: 1.40 $ $Name: not supported by cvs2svn $
     13*  @date $Date: 2004-12-07 22:28:28 $
    1414*
    1515*  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
    … …  
    189189
    190190/******************************************************************************
    191 XXX: It's not clear if this code correctly implements what the equations in
    192 the ADD imply.  Bug 243 will, hopefully, decide that.
    193  
    194 XXX: Once this code is verified, delete Gene's code.
    195191 *****************************************************************************/
    196192psSphere* psSphereTransformApply(psSphere* out,
    … …  
    224220
    225221    return(out);
    226 
    227     /*
    228     This algorithm comes from an email from Gene.  I assume (x,y) corresponds to
    229     (r,d) in the sphere coordinates.
    230      
    231     XXX: In Gene's email, there are different variables with similar names (y
    232     and Y) and in the code, there's cos_y, cos_Y, and cos(y).  Verify that there
    233     are no typo's.
    234      
    235         psF64 sinY = 0.0;
    236         psF64 cosY = 0.0;
    237         psF64 sinX = 0.0;
    238         psF64 cosX = 0.0;
    239         psF64 x = 0.0;
    240         psF64 y = 0.0;
    241         psF64 dx = 0.0;
    242      
    243         x = coord->r;
    244         y = coord->d;
    245         dx = x - transform->phiP;
    246         sinY = cos(y) * sin(dx) * transform->sinDeltaP +
    247                sin(y) * transform->cosDeltaP;
    248         cosY = sqrt(1.0 - sinY * sinY);
    249         sinX = (cos(y) * sin(dx) * transform->cosDeltaP -
    250                 sin(y) * transform->sinDeltaP) / cos(y);
    251         cosX = cos(y) * cos(dx) / cos(y);
    252      
    253         out->r = atan2(sinX, cosX) + transform->alphaP;
    254         out->d = atan2(sinY, cosY);
    255     */
    256 }
    257 
    258 // XXX: the ADD changed.  This code no longer conforms to the ADD.  In fact,
    259 // there is no algorithm for this in the ADD.
    260 // XXX: This is bug 244
    261 psSphereTransform* psSphereTransformICRSToEcliptic(psTime time)
    262 {
    263     PS_PTR_CHECK_NULL(&time, NULL);
    264 
    265     struct tm *tmTime = psTimeToTM(&time);
    266     psF64 year = (psF64)(1900 + tmTime->tm_year);
    267     psF64 T = year / 100.0;
    268     psF64 phi = DEG_TO_RAD(-23.452294) +
    269                 (DEG_TO_RAD(0.013013) * T) +
    270                 (DEG_TO_RAD(0.000001639) * T * T) -
    271                 (DEG_TO_RAD(0.000000503) * T * T * T);
    272     psF64 alphaP = 0.0;
    273     psF64 phiP = 0.0;
    274 
    275     psFree(tmTime);
    276     return (psSphereTransformAlloc(phi, alphaP, phiP));
    277 }
    278 
    279 // XXX: the ADD changed.  This code no longer conforms to the old ADD.
    280 // XXX: Waiting for IfA to sepcify if T should include days/months/etc.
    281 psSphereTransform* psSphereTransformEclipticToICRS(psTime time)
    282 {
    283     PS_PTR_CHECK_NULL(&time, NULL);
    284 
    285     struct tm *tmTime = psTimeToTM(&time);
    286     psF64 year = (psF64)(1900 + tmTime->tm_year);
     222}
     223
     224psSphereTransform* psSphereTransformICRSToEcliptic(psTime *time)
     225{
     226    PS_PTR_CHECK_NULL(time, NULL);
     227    psF64 MJD = psTimeToMJD(time);
     228    struct tm *tmTime = psTimeToTM(time);
     229    psF64 year = (psF64)(1900 - MJD + tmTime->tm_year);
    287230    psF64 T = year / 100.0;
    288231    psF64 alphaP = 0.0;
    … …  
    296239
    297240    psFree(tmTime);
     241
     242    // Don't neglect the minus sign on deltaP (bug 244):
     243    return (psSphereTransformAlloc(alphaP, -deltaP, phiP));
     244}
     245
     246
     247psSphereTransform* psSphereTransformEclipticToICRS(psTime *time)
     248{
     249    PS_PTR_CHECK_NULL(time, NULL);
     250    psF64 MJD = psTimeToMJD(time);
     251    struct tm *tmTime = psTimeToTM(time);
     252    psF64 year = (psF64)(1900 - MJD + tmTime->tm_year);
     253    psF64 T = year / 100.0;
     254    psF64 alphaP = 0.0;
     255    psF64 deltaP = DEG_TO_RAD(23.0) +
     256                   SEC_TO_RAD(27.0) +
     257                   MIN_TO_RAD(8.0) -
     258                   (SEC_TO_RAD(46.845) * T) -
     259                   (SEC_TO_RAD(0.0059) * T * T) +
     260                   (SEC_TO_RAD(0.00181) * T * T * T);
     261    psF64 phiP = 0.0;
     262
     263    psFree(tmTime);
    298264    return (psSphereTransformAlloc(alphaP, deltaP, phiP));
    299265}
    300266
    301 // XXX: Should these be in radians, not degrees?
    302 // XXX: the ADD changed.  This code no longer conforms to the ADD.  In fact,
    303 // there is no algorithm for this in the ADD.
    304 // XXX: This is bug 245
     267// XXX: This is bug 245: alphaP swaps with phiP from psSphereTransformGalaticToICRS()
    305268psSphereTransform* psSphereTransformICRSToGalatic(void)
    306269{
    307     psF64 alphaP = DEG_TO_RAD(-62.6);
    308     psF64 deltaP = DEG_TO_RAD(33.0);
    309     psF64 phiP = DEG_TO_RAD(282.25);
     270    psF64 alphaP = DEG_TO_RAD(32.93192);
     271    psF64 deltaP = DEG_TO_RAD(-62.87175);
     272    psF64 phiP = DEG_TO_RAD(282.85948);
    310273
    311274    return (psSphereTransformAlloc(alphaP, deltaP, phiP));
    312275}
    313276
    314 // XXX: Should these be in radians, not degrees?
    315 // XXX: the ADD changed.  This code no longer conforms to the old ADD.
    316277psSphereTransform* psSphereTransformGalaticToICRS(void)
    317278{
    … …  
    588549
    589550
     551
    590552/******************************************************************************
    591553psSpherePrecess(coords, fromTime, toTime):
    … …  
    599561                          const psTime *toTime)
    600562{
    601     psF64 fromMJD = fromTime->sec/86400.0 +
    602                     fromTime->usec/86400000000.0 +
    603                     40587.0;
    604     psF64 toMJD = toTime->sec/86400.0 +
    605                   toTime->usec/86400000000.0 +
    606                   40587.0;
     563    psF64 fromMJD = psTimeToMJD(fromTime);
     564    psF64 toMJD = psTimeToMJD(toTime);
    607565    psF64 T = (toMJD - fromMJD) / 36525.0;
    608566

  • trunk/psLib/src/astronomy/psCoord.h

    r2635 r2650  
    1 /** @file  psCoord.h
     1a/** @file  psCoord.h
    22*
    33*  @brief Contains basic coordinate transformation definitions and operations
    … …  
    1010*  @author GLG, MHPCC
    1111*
    12 *  @version $Revision: 1.21 $ $Name: not supported by cvs2svn $
    13 *  @date $Date: 2004-12-06 19:49:36 $
     12*  @version $Revision: 1.22 $ $Name: not supported by cvs2svn $
     13*  @date $Date: 2004-12-07 22:28:28 $
    1414*
    1515*  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
    … …  
    239239 */
    240240psSphereTransform* psSphereTransformICRSToEcliptic(
    241     psTime time                        ///< the time for which the resulting transform will be valid
     241    psTime *time                        ///< the time for which the resulting transform will be valid
    242242);
    243243
    … …  
    248248 */
    249249psSphereTransform* psSphereTransformEclipticToICRS(
    250     psTime time                        ///< the time for which the resulting transform will be valid
     250    psTime *time                        ///< the time for which the resulting transform will be valid
    251251);
    252252

  • trunk/psLib/src/dataManip/psMinimize.c

    r2583 r2650  
    99 *  @author GLG, MHPCC
    1010 *
    11  *  @version $Revision: 1.92 $ $Name: not supported by cvs2svn $
    12  *  @date $Date: 2004-12-01 19:56:06 $
     11 *  @version $Revision: 1.93 $ $Name: not supported by cvs2svn $
     12 *  @date $Date: 2004-12-07 22:28:29 $
    1313 *
    1414 *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
    … …  
    455455}
    456456
     457/*
     458XXX: from bug 230:
     459 
     460We first perform a rotation:
     461u = - (x-x0)*cos(theta) + (y-y0)*sin(theta)
     462v = (x-x0)*cos(theta) + (y-y0)*sin(theta)
     463 
     464Here u is the major axis, and v is the minor axis, x0,y0 is the centre, and
     465theta is the position angle.
     466 
     467Then the flux is
     468 
     469flux = norm * exp(-( u*u/2.0/sigmau/sigmau + v*v/2.0/sigmav/sigmav)
     470)/2.0/pi/sigmau/sigmav
     471 
     472Here sigmau and sigmav are the widths of the major and minor axes.
     473 
     474The "norm" parameter in the equation above corresponds to the normalisation.
     475 
     476Suggest order:
     477 
     478norm
     479x0
     480y0
     481sigma_u
     482sigma_v
     483theta
     484*/
    457485
    458486psVector *psMinimizeLMChi2Gauss2D(psImage *deriv,

  • trunk/psLib/src/math/psMinimize.c

    r2583 r2650  
    99 *  @author GLG, MHPCC
    1010 *
    11  *  @version $Revision: 1.92 $ $Name: not supported by cvs2svn $
    12  *  @date $Date: 2004-12-01 19:56:06 $
     11 *  @version $Revision: 1.93 $ $Name: not supported by cvs2svn $
     12 *  @date $Date: 2004-12-07 22:28:29 $
    1313 *
    1414 *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
    … …  
    455455}
    456456
     457/*
     458XXX: from bug 230:
     459 
     460We first perform a rotation:
     461u = - (x-x0)*cos(theta) + (y-y0)*sin(theta)
     462v = (x-x0)*cos(theta) + (y-y0)*sin(theta)
     463 
     464Here u is the major axis, and v is the minor axis, x0,y0 is the centre, and
     465theta is the position angle.
     466 
     467Then the flux is
     468 
     469flux = norm * exp(-( u*u/2.0/sigmau/sigmau + v*v/2.0/sigmav/sigmav)
     470)/2.0/pi/sigmau/sigmav
     471 
     472Here sigmau and sigmav are the widths of the major and minor axes.
     473 
     474The "norm" parameter in the equation above corresponds to the normalisation.
     475 
     476Suggest order:
     477 
     478norm
     479x0
     480y0
     481sigma_u
     482sigma_v
     483theta
     484*/
    457485
    458486psVector *psMinimizeLMChi2Gauss2D(psImage *deriv,
Note: See TracChangeset for help on using the changeset viewer.