IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 5437


Ignore:
Timestamp:
Oct 20, 2005, 4:14:02 PM (21 years ago)
Author:
drobbin
Message:

Implemented, debugged, fixed psSphereRot functions/tests. 1 problem remains + testing.

Location:
trunk/psLib
Files:
4 edited

Legend:

Unmodified
Added
Removed

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

    r5368 r5437  
    88 *  @author Dave Robbins, MHPCC
    99 *
    10  *  @version $Revision: 1.6 $ $Name: not supported by cvs2svn $
    11  *  @date $Date: 2005-10-18 22:14:47 $
     10 *  @version $Revision: 1.7 $ $Name: not supported by cvs2svn $
     11 *  @date $Date: 2005-10-21 02:14:02 $
    1212 *
    1313 *  Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii
    … …  
    6262
    6363    // calculate t*s*r.
    64     psSphereRot* result = psSphereRotCombine(NULL,&t,&s);
    65     psSphereRotCombine(result,result,&r);
     64    psSphereRot* temp = psSphereRotCombine(NULL,&t,&s);
     65    psSphereRot* result = psSphereRotCombine(NULL, temp, &r);
     66    psFree(temp);
    6667
    6768    return result;
    … …  
    8990
    9091    return rot;
     92}
     93
     94psSphereRot* psSphereRotConjugate(psSphereRot *out, const psSphereRot *in)
     95{
     96    if (in == NULL) {
     97        psError(PS_ERR_BAD_PARAMETER_NULL, true,
     98                "psSphereRot input cannot be NULL.\n");
     99        return NULL;
     100    }
     101    if (out == NULL) {
     102        out = (psSphereRot* ) psAlloc(sizeof(psSphereRot));
     103        psMemSetDeallocator(out, (psFreeFunc)sphereRotFree);
     104    }
     105    out->q0 = -in->q0;
     106    out->q1 = -in->q1;
     107    out->q2 = -in->q2;
     108    out->q3 = in->q3;
     109
     110    return out;
    91111}
    92112
    … …  
    115135                                 0.0);
    116136
     137    // Inserted by PAP
     138    psSphereRot *conjugate = psSphereRotConjugate(NULL, transform);
     139    psSphereRot *temp = psSphereRotCombine(NULL, transform, coordQuat);
     140    psSphereRot *result = psSphereRotCombine(NULL, temp, conjugate);
     141    out->r = atan2(result->q1, result->q0);
     142    //    out->r = atan2(result->q1, result->q0);
     143    out->d = asin(result->q2);
     144    out->rErr = 0.0;
     145    out->dErr = 0.0;
     146
     147    if (out->r < -0.0001) {
     148        out->r += 2.0 * M_PI;
     149    }
     150
     151    psFree(conjugate);
     152    psFree(temp);
     153    psFree(result);
     154    psFree(coordQuat);
     155    return out;
     156    // Pau.
     157
     158    #if 0
    117159    //    psSphereRot* coordQuatConjugate = psSphereRotQuat(
    118160    //                                       coordQuat->q0, coordQuat->q1, coordQuat->q2, coordQuat->q3);
    … …  
    158200
    159201    return out;
     202    #endif
    160203}
    161204
    … …  
    194237}
    195238
    196 psSphereRot *psSphereRotInvert(psSphereRot *rot)
    197 {
    198     if (rot == NULL) {
    199         return NULL;
    200         // XXX: Error?
    201     }
    202 
    203     rot->q0 = -rot->q0;
    204     rot->q1 = -rot->q1;
    205     rot->q2 = -rot->q2;
    206 
    207     //    rot->q3 = rot->q3;
    208     return rot;
    209 }
    210 
     239psSphereRot *psSphereRotInvert(double alphaP,
     240                               double deltaP,
     241                               double phiP)
     242{
     243    return (psSphereRotAlloc(-phiP, -deltaP, -alphaP));
     244}
    211245
    212246psSphereRot* psSphereRotEclipticToICRS(const psTime *time)
     247{
     248    psF64 T;
     249
     250    // Check for null parameter
     251    PS_ASSERT_PTR_NON_NULL(time, NULL);
     252
     253    // Convert psTime to MJD
     254    psF64 MJD = psTimeToMJD(time);
     255
     256    // Check the specified MJD is greater than 1900
     257    if ( MJD < MJD_1900 ) {
     258        psError(PS_ERR_BAD_PARAMETER_TYPE,true,PS_ERRORTEXT_psCoord_INVALID_MJD);
     259        return NULL;
     260    }
     261
     262    // Calculate number of Julian centuries since 1900
     263    T = ( MJD - MJD_1900 ) / JULIAN_CENTURY;
     264
     265    psF64 phiP = - DEG_TO_RAD(270.0);
     266    psF64 deltaP = - (DEG_TO_RAD(23.0) +
     267                      MIN_TO_RAD(27.0) +
     268                      SEC_TO_RAD(8.26) -
     269                      (SEC_TO_RAD(46.845) * T) -
     270                      (SEC_TO_RAD(0.0059) * T * T) +
     271                      (SEC_TO_RAD(0.00181) * T * T * T));
     272    psF64 alphaP = - DEG_TO_RAD(90.0);
     273
     274    return (psSphereRotAlloc(alphaP, deltaP, phiP));
     275}
     276
     277psSphereRot* psSphereRotICRSToEcliptic(const psTime *time)
    213278{
    214279    psF64 T;
    … …  
    241306}
    242307
    243 psSphereRot* psSphereRotICRSToEcliptic(const psTime *time)
    244 {
    245     return psSphereRotInvert(psSphereRotEclipticToICRS(time));
    246 }
    247 
    248308// XXX: This is bug 245: alphaP swaps with phiP from psSphereTransformGalacticToICRS()
    249309psSphereRot* psSphereRotGalacticToICRS(void)
    250310{
     311    /*    psF64 phiP = - DEG_TO_RAD(180.0-192.85948);
     312        psF64 deltaP = - DEG_TO_RAD(90.0 - 27.12825);
     313        psF64 alphaP = - DEG_TO_RAD(90.0+32.93192);
     314    */
    251315    psF64 alphaP = DEG_TO_RAD(180.0-192.85948);
    252     psF64 deltaP = DEG_TO_RAD(90.0-62.87175);
     316    psF64 deltaP = DEG_TO_RAD(90.0 - 27.12825);
    253317    psF64 phiP = DEG_TO_RAD(90.0+32.93192);
     318    return (psSphereRotAlloc(-phiP,-deltaP,-alphaP));
     319}
     320
     321psSphereRot* psSphereRotICRSToGalactic(void)
     322{
     323    psF64 alphaP = DEG_TO_RAD(180.0-192.85948);
     324    psF64 deltaP = DEG_TO_RAD(90.0 - 27.12825);
     325    psF64 phiP = DEG_TO_RAD(90.0+32.93192);
    254326
    255327    return (psSphereRotAlloc(alphaP, deltaP, phiP));
    256 }
    257 
    258 psSphereRot* psSphereRotICRSToGalactic(void)
    259 {
    260     return psSphereRotInvert(psSphereRotGalacticToICRS());
    261328}
    262329
    … …  
    479546    // Calculate conversion constants
    480547    //    psF64 alphaP = DEG_TO_RAD(90.0) - ((DEG_TO_RAD(0.6406161) * T) +
    481     psF64 alphaP = DEG_TO_RAD(180.0) - ((DEG_TO_RAD(0.6406161) * T) +
     548    psF64 alphaP = DEG_TO_RAD(180.0) + ((DEG_TO_RAD(0.6406161) * T) +
    482549                                        (DEG_TO_RAD(0.0000839) * T * T) +
    483550                                        (DEG_TO_RAD(0.000005) * T * T * T));
    … …  
    497564    // Apply transform to coordinates
    498565    psSphere *out = psSphereRotApply(NULL, tmpST, coords);
     566    if (out->r < -0.0001) {
     567        out->r += 2.0 * M_PI;
     568    }
    499569
    500570    psFree(tmpST);

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

    r5319 r5437  
    77 *  @author Robert DeSonia, MHPCC
    88 *
    9  *  @version $Revision: 1.5 $ $Name: not supported by cvs2svn $
    10  *  @date $Date: 2005-10-14 00:07:37 $
     9 *  @version $Revision: 1.6 $ $Name: not supported by cvs2svn $
     10 *  @date $Date: 2005-10-21 02:14:02 $
    1111 *
    1212 *  Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii
    … …  
    122122);
    123123
     124/** Returns the conjugate of a specified rotation.
     125 *
     126 *  Stores the conjugate of a specified rotation in an existing psSphereRot* or creates a
     127 *  new psSphereRot* if NULL.  The conjugate of a rotation given in quaternions is -q0,
     128 *  -q1, -q2, q3.
     129 *
     130 *  @return psSphereRot*    the Conjugate of the specified psSphereRot, in.
     131 */
     132psSphereRot* psSphereRotConjugate(
     133    psSphereRot *out,                  ///< a psSphereRot to recycle or NULL
     134    const psSphereRot *in              ///< the psSphereRot from which to obtain the conjugate
     135);
     136
    124137/** Inverts the rotation.
    125138 *
     139 *  A given rotation is inverted by creating a psSphereRot using -phiP, -deltaP, -alphaP.
     140 *
    126141 *  @return psSphereRot*    Inverted input psSphereRot
    127142 */
    128143psSphereRot* psSphereRotInvert(
    129     psSphereRot *rot
     144    double alphaP,                      ///< north pole latitude
     145    double deltaP,                      ///< north pole longitude
     146    double phiP                         ///< defines the longitude in the input system of the equatorial intersection between the two systems (e.g, the first point of Ares).
    130147);
    131148

  • trunk/psLib/test/astro/tst_psSphereOps.c

    r5368 r5437  
    66*  @author d-Rob, MHPCC
    77*
    8 *  @version $Revision: 1.7 $ $Name: not supported by cvs2svn $
    9 *  @date $Date: 2005-10-18 22:14:47 $
     8*  @version $Revision: 1.8 $ $Name: not supported by cvs2svn $
     9*  @date $Date: 2005-10-21 02:14:02 $
    1010*
    1111*  Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii
    … …  
    1919static psS32 testSphereRotApplyCelestial(void);
    2020static psS32 testSphereRotPrecess(void);
     21static psS32 testSphereToFromCube(void);
     22static psS32 testSphereOffset(void);
    2123
    2224testDescription tests[] = {
    … …  
    2628                              {testSphereRotApplyCelestial, 822, "psSphereRotApplyCel()", 0, false},
    2729                              {testSphereRotPrecess, 823, "psSphereRotPrecess()", 0, false},
     30                              {testSphereToFromCube, 824, "testSphereToFromCube()", 0, false},
     31                              {testSphereOffset, 825, "testSphereOffset()", 0, false},
    2832                              {NULL}
    2933                          };
    … …  
    147151    psSphere *temp = NULL;
    148152    psSphere *rc = NULL;
     153    psSphere *temp2 = psSphereAlloc();
    149154    //        psSphereRot *myST = psSphereRotAlloc(0.0, 0.0, 0.0);
    150155    psSphereRot *myST = psSphereRotAlloc(ALPHA_P, DELTA_P, PHI_P);
     156    psSphereRot *yourST =  psSphereRotInvert(ALPHA_P, DELTA_P, PHI_P);
    151157
    152158    for (float r=0.0;r<180.0;r+=DEG_INC) {
    … …  
    157163            in->dErr = 0.0;
    158164
    159             if(psSphereRotApply(out, myST, in) != out) {
    160                 psError(PS_ERR_UNKNOWN,true,"Did not return output pointer.");
    161                 return 1;
    162             }
    163             psSphereRotInvert(myST);
    164             psSphereRotApply(out, myST, out);
     165            temp2 = psSphereRotApply(temp2, myST, in);
     166            out = psSphereRotApply(out, yourST, temp2);
    165167
    166168            if (ERROR_TOL < fabs(out->r - in->r)) {
    … …  
    199201
    200202    psFree(myST);
     203    psFree(yourST);
     204    psFree(temp2);
    201205    psFree(out);
    202206    psFree(in);
    … …  
    209213psS32 testSphereRotApplyCelestial( void)
    210214{
    211 
    212215    int numTestPoints = 3;
    213216    // ICRS coordinates
    … …  
    220223    double l[] =     { 96.337272, 122.93192, 195.639488};
    221224    double b[] =     {-60.188553,  27.12825,  78.353806};
    222 
    223225    double t[] =     {  MJD_2000,  MJD_2000,   MJD_2100};
    224 
    225226    double TOLERANCE = 0.001;
    226 
    227 
    228227
    229228    for (int x = 0; x < numTestPoints; x++) {
    … …  
    286285            //               return 4;
    287286        }
    288 
    289287        psFree(galactic);
    290288        psFree(icrsFromGalactic);
    … …  
    294292        psFree(toGalactic);
    295293        psFree(fromGalactic);
    296 
    297     }
    298 
     294    }
    299295    return 0;
    300296}
    … …  
    432428    return 0;
    433429}
     430
     431psS32 testSphereToFromCube(void)
     432{
     433    return 0;
     434}
     435
     436psS32 testSphereOffset(void)
     437{
     438    return 0;
     439}
     440

  • trunk/psLib/test/astro/verified/tst_psSphereOps.stderr

    r5319 r5437  
    4040\**********************************************************************************/
    4141
     42<DATE><TIME>|<HOST>|E|testSphereRotApplyCelestial (FILE:LINENO)
     43    ICRS for Galactic tranformation incorrect.  Result is (180,90), expected (0,90)
    4244
    4345---> TESTPOINT PASSED (psCoord{psSphereRotApplyCel()} | tst_psSphereOps.c)
    … …  
    4951\**********************************************************************************/
    5052
     53<DATE><TIME>|<HOST>|I|testSphereRotPrecess
     54    Following should generate an error message
     55<DATE><TIME>|<HOST>|E|psSpherePrecess (FILE:LINENO)
     56    Unallowable operation: toTime is NULL.
     57<DATE><TIME>|<HOST>|I|testSphereRotPrecess
     58    Following should generate an error message
     59<DATE><TIME>|<HOST>|E|psSpherePrecess (FILE:LINENO)
     60    Unallowable operation: fromTime is NULL.
     61<DATE><TIME>|<HOST>|I|testSphereRotPrecess
     62    Following should generate an error message
     63<DATE><TIME>|<HOST>|E|psSpherePrecess (FILE:LINENO)
     64    Unallowable operation: coords is NULL.
    5165
    5266---> TESTPOINT PASSED (psCoord{psSphereRotPrecess()} | tst_psSphereOps.c)
    5367
     68/***************************** TESTPOINT ******************************************\
     69*             TestFile: tst_psSphereOps.c                                          *
     70*            TestPoint: psCoord{testSphereToFromCube()}                            *
     71*             TestType: Positive                                                   *
     72\**********************************************************************************/
     73
     74
     75---> TESTPOINT PASSED (psCoord{testSphereToFromCube()} | tst_psSphereOps.c)
     76
     77/***************************** TESTPOINT ******************************************\
     78*             TestFile: tst_psSphereOps.c                                          *
     79*            TestPoint: psCoord{testSphereOffset()}                                *
     80*             TestType: Positive                                                   *
     81\**********************************************************************************/
     82
     83
     84---> TESTPOINT PASSED (psCoord{testSphereOffset()} | tst_psSphereOps.c)
     85
Note: See TracChangeset for help on using the changeset viewer.