Index: /trunk/psLib/src/astro/psSphereOps.c =================================================================== --- /trunk/psLib/src/astro/psSphereOps.c (revision 5436) +++ /trunk/psLib/src/astro/psSphereOps.c (revision 5437) @@ -8,6 +8,6 @@ * @author Dave Robbins, MHPCC * - * @version $Revision: 1.6 $ $Name: not supported by cvs2svn $ - * @date $Date: 2005-10-18 22:14:47 $ + * @version $Revision: 1.7 $ $Name: not supported by cvs2svn $ + * @date $Date: 2005-10-21 02:14:02 $ * * Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii @@ -62,6 +62,7 @@ // calculate t*s*r. - psSphereRot* result = psSphereRotCombine(NULL,&t,&s); - psSphereRotCombine(result,result,&r); + psSphereRot* temp = psSphereRotCombine(NULL,&t,&s); + psSphereRot* result = psSphereRotCombine(NULL, temp, &r); + psFree(temp); return result; @@ -89,4 +90,23 @@ return rot; +} + +psSphereRot* psSphereRotConjugate(psSphereRot *out, const psSphereRot *in) +{ + if (in == NULL) { + psError(PS_ERR_BAD_PARAMETER_NULL, true, + "psSphereRot input cannot be NULL.\n"); + return NULL; + } + if (out == NULL) { + out = (psSphereRot* ) psAlloc(sizeof(psSphereRot)); + psMemSetDeallocator(out, (psFreeFunc)sphereRotFree); + } + out->q0 = -in->q0; + out->q1 = -in->q1; + out->q2 = -in->q2; + out->q3 = in->q3; + + return out; } @@ -115,4 +135,26 @@ 0.0); + // Inserted by PAP + psSphereRot *conjugate = psSphereRotConjugate(NULL, transform); + psSphereRot *temp = psSphereRotCombine(NULL, transform, coordQuat); + psSphereRot *result = psSphereRotCombine(NULL, temp, conjugate); + out->r = atan2(result->q1, result->q0); + // out->r = atan2(result->q1, result->q0); + out->d = asin(result->q2); + out->rErr = 0.0; + out->dErr = 0.0; + + if (out->r < -0.0001) { + out->r += 2.0 * M_PI; + } + + psFree(conjugate); + psFree(temp); + psFree(result); + psFree(coordQuat); + return out; + // Pau. + + #if 0 // psSphereRot* coordQuatConjugate = psSphereRotQuat( // coordQuat->q0, coordQuat->q1, coordQuat->q2, coordQuat->q3); @@ -158,4 +200,5 @@ return out; + #endif } @@ -194,21 +237,43 @@ } -psSphereRot *psSphereRotInvert(psSphereRot *rot) -{ - if (rot == NULL) { - return NULL; - // XXX: Error? - } - - rot->q0 = -rot->q0; - rot->q1 = -rot->q1; - rot->q2 = -rot->q2; - - // rot->q3 = rot->q3; - return rot; -} - +psSphereRot *psSphereRotInvert(double alphaP, + double deltaP, + double phiP) +{ + return (psSphereRotAlloc(-phiP, -deltaP, -alphaP)); +} psSphereRot* psSphereRotEclipticToICRS(const psTime *time) +{ + psF64 T; + + // Check for null parameter + PS_ASSERT_PTR_NON_NULL(time, NULL); + + // Convert psTime to MJD + psF64 MJD = psTimeToMJD(time); + + // Check the specified MJD is greater than 1900 + if ( MJD < MJD_1900 ) { + psError(PS_ERR_BAD_PARAMETER_TYPE,true,PS_ERRORTEXT_psCoord_INVALID_MJD); + return NULL; + } + + // Calculate number of Julian centuries since 1900 + T = ( MJD - MJD_1900 ) / JULIAN_CENTURY; + + psF64 phiP = - DEG_TO_RAD(270.0); + psF64 deltaP = - (DEG_TO_RAD(23.0) + + MIN_TO_RAD(27.0) + + SEC_TO_RAD(8.26) - + (SEC_TO_RAD(46.845) * T) - + (SEC_TO_RAD(0.0059) * T * T) + + (SEC_TO_RAD(0.00181) * T * T * T)); + psF64 alphaP = - DEG_TO_RAD(90.0); + + return (psSphereRotAlloc(alphaP, deltaP, phiP)); +} + +psSphereRot* psSphereRotICRSToEcliptic(const psTime *time) { psF64 T; @@ -241,22 +306,24 @@ } -psSphereRot* psSphereRotICRSToEcliptic(const psTime *time) -{ - return psSphereRotInvert(psSphereRotEclipticToICRS(time)); -} - // XXX: This is bug 245: alphaP swaps with phiP from psSphereTransformGalacticToICRS() psSphereRot* psSphereRotGalacticToICRS(void) { + /* psF64 phiP = - DEG_TO_RAD(180.0-192.85948); + psF64 deltaP = - DEG_TO_RAD(90.0 - 27.12825); + psF64 alphaP = - DEG_TO_RAD(90.0+32.93192); + */ psF64 alphaP = DEG_TO_RAD(180.0-192.85948); - psF64 deltaP = DEG_TO_RAD(90.0-62.87175); + psF64 deltaP = DEG_TO_RAD(90.0 - 27.12825); psF64 phiP = DEG_TO_RAD(90.0+32.93192); + return (psSphereRotAlloc(-phiP,-deltaP,-alphaP)); +} + +psSphereRot* psSphereRotICRSToGalactic(void) +{ + psF64 alphaP = DEG_TO_RAD(180.0-192.85948); + psF64 deltaP = DEG_TO_RAD(90.0 - 27.12825); + psF64 phiP = DEG_TO_RAD(90.0+32.93192); return (psSphereRotAlloc(alphaP, deltaP, phiP)); -} - -psSphereRot* psSphereRotICRSToGalactic(void) -{ - return psSphereRotInvert(psSphereRotGalacticToICRS()); } @@ -479,5 +546,5 @@ // Calculate conversion constants // psF64 alphaP = DEG_TO_RAD(90.0) - ((DEG_TO_RAD(0.6406161) * T) + - psF64 alphaP = DEG_TO_RAD(180.0) - ((DEG_TO_RAD(0.6406161) * T) + + psF64 alphaP = DEG_TO_RAD(180.0) + ((DEG_TO_RAD(0.6406161) * T) + (DEG_TO_RAD(0.0000839) * T * T) + (DEG_TO_RAD(0.000005) * T * T * T)); @@ -497,4 +564,7 @@ // Apply transform to coordinates psSphere *out = psSphereRotApply(NULL, tmpST, coords); + if (out->r < -0.0001) { + out->r += 2.0 * M_PI; + } psFree(tmpST); Index: /trunk/psLib/src/astro/psSphereOps.h =================================================================== --- /trunk/psLib/src/astro/psSphereOps.h (revision 5436) +++ /trunk/psLib/src/astro/psSphereOps.h (revision 5437) @@ -7,6 +7,6 @@ * @author Robert DeSonia, MHPCC * - * @version $Revision: 1.5 $ $Name: not supported by cvs2svn $ - * @date $Date: 2005-10-14 00:07:37 $ + * @version $Revision: 1.6 $ $Name: not supported by cvs2svn $ + * @date $Date: 2005-10-21 02:14:02 $ * * Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii @@ -122,10 +122,27 @@ ); +/** Returns the conjugate of a specified rotation. + * + * Stores the conjugate of a specified rotation in an existing psSphereRot* or creates a + * new psSphereRot* if NULL. The conjugate of a rotation given in quaternions is -q0, + * -q1, -q2, q3. + * + * @return psSphereRot* the Conjugate of the specified psSphereRot, in. + */ +psSphereRot* psSphereRotConjugate( + psSphereRot *out, ///< a psSphereRot to recycle or NULL + const psSphereRot *in ///< the psSphereRot from which to obtain the conjugate +); + /** Inverts the rotation. * + * A given rotation is inverted by creating a psSphereRot using -phiP, -deltaP, -alphaP. + * * @return psSphereRot* Inverted input psSphereRot */ psSphereRot* psSphereRotInvert( - psSphereRot *rot + double alphaP, ///< north pole latitude + double deltaP, ///< north pole longitude + double phiP ///< defines the longitude in the input system of the equatorial intersection between the two systems (e.g, the first point of Ares). ); Index: /trunk/psLib/test/astro/tst_psSphereOps.c =================================================================== --- /trunk/psLib/test/astro/tst_psSphereOps.c (revision 5436) +++ /trunk/psLib/test/astro/tst_psSphereOps.c (revision 5437) @@ -6,6 +6,6 @@ * @author d-Rob, MHPCC * -* @version $Revision: 1.7 $ $Name: not supported by cvs2svn $ -* @date $Date: 2005-10-18 22:14:47 $ +* @version $Revision: 1.8 $ $Name: not supported by cvs2svn $ +* @date $Date: 2005-10-21 02:14:02 $ * * Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii @@ -19,4 +19,6 @@ static psS32 testSphereRotApplyCelestial(void); static psS32 testSphereRotPrecess(void); +static psS32 testSphereToFromCube(void); +static psS32 testSphereOffset(void); testDescription tests[] = { @@ -26,4 +28,6 @@ {testSphereRotApplyCelestial, 822, "psSphereRotApplyCel()", 0, false}, {testSphereRotPrecess, 823, "psSphereRotPrecess()", 0, false}, + {testSphereToFromCube, 824, "testSphereToFromCube()", 0, false}, + {testSphereOffset, 825, "testSphereOffset()", 0, false}, {NULL} }; @@ -147,6 +151,8 @@ psSphere *temp = NULL; psSphere *rc = NULL; + psSphere *temp2 = psSphereAlloc(); // psSphereRot *myST = psSphereRotAlloc(0.0, 0.0, 0.0); psSphereRot *myST = psSphereRotAlloc(ALPHA_P, DELTA_P, PHI_P); + psSphereRot *yourST = psSphereRotInvert(ALPHA_P, DELTA_P, PHI_P); for (float r=0.0;r<180.0;r+=DEG_INC) { @@ -157,10 +163,6 @@ in->dErr = 0.0; - if(psSphereRotApply(out, myST, in) != out) { - psError(PS_ERR_UNKNOWN,true,"Did not return output pointer."); - return 1; - } - psSphereRotInvert(myST); - psSphereRotApply(out, myST, out); + temp2 = psSphereRotApply(temp2, myST, in); + out = psSphereRotApply(out, yourST, temp2); if (ERROR_TOL < fabs(out->r - in->r)) { @@ -199,4 +201,6 @@ psFree(myST); + psFree(yourST); + psFree(temp2); psFree(out); psFree(in); @@ -209,5 +213,4 @@ psS32 testSphereRotApplyCelestial( void) { - int numTestPoints = 3; // ICRS coordinates @@ -220,10 +223,6 @@ double l[] = { 96.337272, 122.93192, 195.639488}; double b[] = {-60.188553, 27.12825, 78.353806}; - double t[] = { MJD_2000, MJD_2000, MJD_2100}; - double TOLERANCE = 0.001; - - for (int x = 0; x < numTestPoints; x++) { @@ -286,5 +285,4 @@ // return 4; } - psFree(galactic); psFree(icrsFromGalactic); @@ -294,7 +292,5 @@ psFree(toGalactic); psFree(fromGalactic); - - } - + } return 0; } @@ -432,2 +428,13 @@ return 0; } + +psS32 testSphereToFromCube(void) +{ + return 0; +} + +psS32 testSphereOffset(void) +{ + return 0; +} + Index: /trunk/psLib/test/astro/verified/tst_psSphereOps.stderr =================================================================== --- /trunk/psLib/test/astro/verified/tst_psSphereOps.stderr (revision 5436) +++ /trunk/psLib/test/astro/verified/tst_psSphereOps.stderr (revision 5437) @@ -40,4 +40,6 @@ \**********************************************************************************/ +