IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Ignore:
Timestamp:
Oct 13, 2005, 2:07:37 PM (21 years ago)
Author:
drobbin
Message:

Fixed/Updated SphereRot functions and tests.

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/psLib/test/astro/tst_psSphereOps.c

    r5306 r5319  
    1 /** @file  tst_psCoord.c
     1/** @file  tst_psSphereOps.c
    22*
    3 *  @brief The code will ...
     3*  @brief The code will ..... Work ....
    44*
    55*
    6 *  @author GLG, MHPCC
     6*  @author d-Rob, MHPCC
    77*
    8 *  @version $Revision: 1.4 $ $Name: not supported by cvs2svn $
    9 *  @date $Date: 2005-10-13 20:23:57 $
     8*  @version $Revision: 1.5 $ $Name: not supported by cvs2svn $
     9*  @date $Date: 2005-10-14 00:07:37 $
    1010*
    1111*  Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii
     
    1515#include "pslib_strict.h"
    1616static psS32 testSphereRotAlloc(void);
     17static psS32 testSphereRotQuat(void);
    1718static psS32 testSphereRotApply1(void);
    18 static psS32 testSphereRotApply2(void);
    1919static psS32 testSphereRotApplyCelestial(void);
    2020static psS32 testSphereRotPrecess(void);
     
    2222testDescription tests[] = {
    2323                              {testSphereRotAlloc, 819, "psSphereRotAlloc()", 0, false},
    24                               {testSphereRotApply1, 820, "psSphereRotApply()", 0, false},
    25                               {testSphereRotApply2, 820, "psSphereRotApply()", 0, false},
    26                               {testSphereRotApplyCelestial, 820, "psSphereRotApply()", 0, false},
    27                               //                              {testSphereRotApplyCelestial, 821, "psSphereRotICRSToEcliptic()", true},
    28                               //                              {testSphereRotApplyCelestial, 822, "psSphereRotEclipticToICRS()", 0, true},
    29                               //                              {testSphereRotApplyCelestial, 824, "psSphereRotICRSToGalactic()", 0, true},
    30                               //                              {testSphereRotApplyCelestial, 823, "psSphereRotGalacticToICRS()", 0, true},
    31                               {testSphereRotPrecess, 825, "psSphereRotPrecess()", 0, false},
     24                              {testSphereRotQuat, 820, "psSphereRotQuat()", 0, false},
     25                              {testSphereRotApply1, 821, "psSphereRotApply()", 0, false},
     26                              {testSphereRotApplyCelestial, 822, "psSphereRotApplyCel()", 0, false},
     27                              {testSphereRotPrecess, 823, "psSphereRotPrecess()", 0, false},
     28                              //                           {testSphereRotApplyCelestial, 821, "psSphereRotICRSToEcliptic()", true},
     29                              //                         {testSphereRotApplyCelestial, 822, "psSphereRotEclipticToICRS()", 0, true},
     30                              //                         {testSphereRotApplyCelestial, 824, "psSphereRotICRSToGalactic()", 0, true},
     31                              //                         {testSphereRotApplyCelestial, 823, "psSphereRotGalacticToICRS()", 0, true},
    3232                              {NULL}
    3333                          };
     
    5050#define ALPHA_P M_PI/6
    5151#define DELTA_P M_PI/4
    52 #define PHI_P M_PI/2
     52#define PHI_P M_PI/3
    5353
    5454psS32 testSphereRotAlloc( void )
     
    7373    double q2 = sin(a2)*cos(DELTA_P/2);
    7474
    75     printf("\n q0=%lf, q1=%lf, q2=%lf, q3=%lf,\n", q0, q1, q2, q3);
    76     printf("myq0=%lf, myq1=%lf, myq2=%lf, myq3=%lf\n", myST->q0, myST->q1, myST->q2, myST->q3);
     75    //    printf("\n q0=%lf, q1=%lf, q2=%lf, q3=%lf,\n", q0, q1, q2, q3);
     76    //    printf("myq0=%lf, myq1=%lf, myq2=%lf, myq3=%lf\n", myST->q0, myST->q1, myST->q2, myST->q3);
    7777    if (FLT_EPSILON < fabs(q0 - myST->q0)) {
    7878        psError(PS_ERR_UNKNOWN,true,"myST->q0 is %lf, should be %lf\n", myST->q0, q0);
     
    9898}
    9999
     100psS32 testSphereRotQuat(void)
     101{
     102    double a0 = (ALPHA_P - PHI_P)/2.0;
     103    double a1 = (ALPHA_P - PHI_P)/2.0;
     104    double a2 = (ALPHA_P + PHI_P)/2.0;
     105    double a3 = (ALPHA_P + PHI_P)/2.0;
     106
     107    double q3 = cos(a3)*cos(DELTA_P/2);
     108    double q0 = sin(a0)*sin(DELTA_P/2);
     109    double q1 = cos(a1)*sin(DELTA_P/2);
     110    double q2 = sin(a2)*cos(DELTA_P/2);
     111
     112    psSphereRot *myST = psSphereRotQuat(q0*2.0, q1*2.0, q2*2.0, q3*2.0);
     113    // Verify null not returned
     114    if(myST == NULL) {
     115        psError(PS_ERR_UNKNOWN,true,"Returned NULL with valid parameters");
     116        return 1;
     117    }
     118
     119    //    printf("\n q0=%lf, q1=%lf, q2=%lf, q3=%lf,\n", q0, q1, q2, q3);
     120    //    printf("myq0=%lf, myq1=%lf, myq2=%lf, myq3=%lf\n", myST->q0, myST->q1, myST->q2, myST->q3);
     121    if (FLT_EPSILON < fabs(q0 - myST->q0)) {
     122        psError(PS_ERR_UNKNOWN,true,"myST->q0 is %lf, should be %lf\n", myST->q0, q0);
     123        return 2;
     124    }
     125    if (FLT_EPSILON < fabs(q1 - myST->q1)) {
     126        psError(PS_ERR_UNKNOWN,true,"myST->q1 is %f, should be %f\n", myST->q1, q1);
     127        return 3;
     128    }
     129    if (FLT_EPSILON < fabs(q2 - myST->q2)) {
     130        psError(PS_ERR_UNKNOWN,true,"myST->q2 is %f, should be %f\n", myST->q2, q2);
     131        return 4;
     132    }
     133    if (FLT_EPSILON < fabs(q3 - myST->q3)) {
     134        psError(PS_ERR_UNKNOWN,true,"myST->q0 is %f, should be %f\n", myST->q3, q3);
     135        return 5;
     136    }
     137
     138    // Free data structure
     139    psFree(myST);
     140
     141    return 0;
     142}
     143
     144
    100145// We do a simple identity transformation on a few RA, DEC pairs.
    101146psS32 testSphereRotApply1( void )
    102147{
    103     /*
    104         psSphere *in = psSphereAlloc();
    105         psSphere *out = psSphereAlloc();
    106         psSphere *temp = NULL;
    107         psSphere *rc = NULL;
    108         psSphereRot *myST = psSphereRotAlloc(0.0, 0.0, 0.0);
    109 
    110         for (float r=0.0;r<180.0;r+=DEG_INC) {
    111             for (float d=0.0;d<90.0;d+=DEG_INC) {
    112                 in->r = DEG_TO_RAD(r);
    113                 in->d = DEG_TO_RAD(d);
    114                 in->rErr = 0.0;
    115                 in->dErr = 0.0;
    116 
    117                 if(psSphereRotApply(out, myST, in) != out) {
    118                     psError(PS_ERR_UNKNOWN,true,"Did not return output pointer.");
    119                     return 1;
    120                 }
    121                 if (ERROR_TOL < fabs(out->r - in->r)) {
    122                     psError(PS_ERR_UNKNOWN,true,"out->r is %f, should be %f\n", out->r, in->r);
    123                     return 2;
    124                 }
    125                 if (ERROR_TOL < fabs(out->d - in->d)) {
    126                     psError(PS_ERR_UNKNOWN,true,"out->d is %f, should be %f\n", out->d, in->d);
    127                     return 3;
    128                 }
    129             }
    130         }
    131 
    132         // Verify new sphere object is created if out parameter NULL
    133         temp = psSphereRotApply(NULL, myST, in);
    134         if ( temp == NULL) {
    135             psError(PS_ERR_UNKNOWN,true,"Returned NULL when out parameter was null");
    136             return 4;
    137         }
    138         psFree(temp);
    139 
    140         // Verify NULL returned if transform structure null
    141         psLogMsg(__func__,PS_LOG_INFO,"Following should generate an error");
    142         rc = psSphereRotApply(NULL, NULL, in);
    143         if (rc != NULL) {
    144             psError(PS_ERR_UNKNOWN,true,"psSphereRotApply() did not return NULL.");
    145             return 5;
    146         }
    147 
    148         // Verify NULL returned when input sphere is NULL
    149         psLogMsg(__func__,PS_LOG_INFO,"Following should generate an error");
    150         rc = psSphereRotApply(NULL, myST, NULL);
    151         if (rc != NULL) {
    152             psError(PS_ERR_UNKNOWN,true,"psSphereRotApply() did not return NULL");
    153             return 6;
    154         }
    155 
    156         psFree(myST);
    157         psFree(out);
    158         psFree(in);
    159     */
     148
     149    psSphere *in = psSphereAlloc();
     150    psSphere *out = psSphereAlloc();
     151    psSphere *temp = NULL;
     152    psSphere *rc = NULL;
     153    //        psSphereRot *myST = psSphereRotAlloc(0.0, 0.0, 0.0);
     154    psSphereRot *myST = psSphereRotAlloc(ALPHA_P, DELTA_P, PHI_P);
     155
     156    for (float r=0.0;r<180.0;r+=DEG_INC) {
     157        for (float d=0.0;d<90.0;d+=DEG_INC) {
     158            in->r = DEG_TO_RAD(r);
     159            in->d = DEG_TO_RAD(d);
     160            in->rErr = 0.0;
     161            in->dErr = 0.0;
     162
     163            if(psSphereRotApply(out, myST, in) != out) {
     164                psError(PS_ERR_UNKNOWN,true,"Did not return output pointer.");
     165                return 1;
     166            }
     167            psSphereRotInvert(myST);
     168            psSphereRotApply(out, myST, out);
     169
     170            if (ERROR_TOL < fabs(out->r - in->r)) {
     171                psError(PS_ERR_UNKNOWN,true,"out->r is %f, should be %f\n", out->r, in->r);
     172                return 2;
     173            }
     174            if (ERROR_TOL < fabs(out->d - in->d)) {
     175                psError(PS_ERR_UNKNOWN,true,"out->d is %f, should be %f\n", out->d, in->d);
     176                return 3;
     177            }
     178        }
     179    }
     180    // Verify new sphere object is created if out parameter NULL
     181    temp = psSphereRotApply(NULL, myST, in);
     182    if ( temp == NULL) {
     183        psError(PS_ERR_UNKNOWN,true,"Returned NULL when out parameter was null");
     184        return 4;
     185    }
     186    psFree(temp);
     187
     188    // Verify NULL returned if transform structure null
     189    psLogMsg(__func__,PS_LOG_INFO,"Following should generate an error");
     190    rc = psSphereRotApply(NULL, NULL, in);
     191    if (rc != NULL) {
     192        psError(PS_ERR_UNKNOWN,true,"psSphereRotApply() did not return NULL.");
     193        return 5;
     194    }
     195
     196    // Verify NULL returned when input sphere is NULL
     197    psLogMsg(__func__,PS_LOG_INFO,"Following should generate an error");
     198    rc = psSphereRotApply(NULL, myST, NULL);
     199    if (rc != NULL) {
     200        psError(PS_ERR_UNKNOWN,true,"psSphereRotApply() did not return NULL");
     201        return 6;
     202    }
     203
     204    psFree(myST);
     205    psFree(out);
     206    psFree(in);
     207
    160208    return 0;
    161209}
    162210
    163 /******************************************************************************
    164 testSphereRotApply2(): This test verifies that psSphereRotApply()
    165 works properly.  We create two psSphereRots: a forward transform and a
    166 reverse transform (which is the mathematical inverse of the forward transform).
    167 We apply both transforms to several spherical coordinates and ensure that the
    168 original input coordinate is obtained after applying both transforms.
    169  
    170 XXX: We currently test the alpha and delta offsets independently.  Attempts to
    171 test them both concurrently failed.  Determine why this is.  Are the following
    172 spherical transforms not mathematical inverses?
    173     psSphereRotAlloc(X, Y, 0.0)
    174     psSphereRotAlloc(-X, -Y, 0.0)
    175  *****************************************************************************/
    176211#define ERROR_PERCENT 0.01
    177 psS32 testSphereRotApply2( void )
    178 {
    179     /*
    180         psS32 testStatus = 0;
    181         psSphere in;
    182         psSphere out;
    183         psSphere out2;
    184         psSphereRot *mySphereRotForward = NULL;
    185         psSphereRot *mySphereRotReverse = NULL;
    186 
    187 
    188         mySphereRotForward = psSphereRotAlloc(DEG_TO_RAD(22.0),
    189                                               0.0,
    190                                               0.0);
    191         mySphereRotReverse = psSphereRotAlloc(DEG_TO_RAD(-22.0),
    192                                               0.0,
    193                                               0.0);
    194 
    195         for (float r=0.1;r<180.0;r+=(DEG_INC/5.0)) {
    196             for (float d=0.1;d<90.0;d+=(DEG_INC/5.0)) {
    197                 in.r = DEG_TO_RAD(r);
    198                 in.d = DEG_TO_RAD(d);
    199                 in.rErr = 0.0;
    200                 in.dErr = 0.0;
    201 
    202                 psSphereRotApply(&out, mySphereRotForward, &in);
    203                 psSphereRotApply(&out2, mySphereRotReverse, &out);
    204 
    205                 if ((fabs((in.r - out2.r) / in.r) > ERROR_PERCENT) ||
    206                         (fabs((in.d - out2.d) / in.d) > ERROR_PERCENT)) {
    207                     printf("ERROR: \n");
    208                     printf("Input  coords (R, D) are (%f, %f)\n", in.r, in.d);
    209                     printf("Output coords (R, D) are (%f, %f)\n", out2.r, out2.d);
    210                     testStatus = 4;
    211                 }
    212             }
    213         }
    214         psFree(mySphereRotForward);
    215         psFree(mySphereRotReverse);
    216 
    217         mySphereRotForward = psSphereRotAlloc(0.0,
    218                                               DEG_TO_RAD(33.0),
    219                                               0.0);
    220         mySphereRotReverse = psSphereRotAlloc(0.0,
    221                                               DEG_TO_RAD(-33.0),
    222                                               0.0);
    223         for (float r=0.1;r<180.0;r+=(DEG_INC/5.0)) {
    224             for (float d=0.1;d<90.0;d+=(DEG_INC/5.0)) {
    225                 in.r = DEG_TO_RAD(r);
    226                 in.d = DEG_TO_RAD(d);
    227                 in.rErr = 0.0;
    228                 in.dErr = 0.0;
    229 
    230                 psSphereRotApply(&out, mySphereRotForward, &in);
    231                 psSphereRotApply(&out2, mySphereRotReverse, &out);
    232 
    233                 if ((fabs((in.r - out2.r) / in.r) > ERROR_PERCENT) ||
    234                         (fabs((in.d - out2.d) / in.d) > ERROR_PERCENT)) {
    235                     printf("ERROR: \n");
    236                     printf("Input  coords (R, D) are (%f, %f)\n", in.r, in.d);
    237                     printf("Output coords (R, D) are (%f, %f)\n", out2.r, out2.d);
    238                     testStatus = 4;
    239                 }
    240             }
    241         }
    242         psFree(mySphereRotForward);
    243         psFree(mySphereRotReverse);
    244 
    245         return(testStatus);
    246         */
    247     return 0;
    248 }
    249212
    250213psS32 testSphereRotApplyCelestial( void)
Note: See TracChangeset for help on using the changeset viewer.