C: Added Astronomy_Rotation_EQD_ECL, ECL_EQD, HOR_ECL, ECL_HOR.

These are the final 4 rotation functions to complete every
possible coordinate transform.

Finished unit tests of verifying that all triangular
cycles of transitive rotation are consistent.
This commit is contained in:
Don Cross
2019-12-09 16:53:02 -05:00
parent 6876aee963
commit 55d4cc3dde
5 changed files with 354 additions and 17 deletions

View File

@@ -1775,7 +1775,7 @@ static int CheckUnitVector(int lnum, const char *name, astro_rotation_t r, int i
sum += x*x;
}
x = fabs(sum - 1.0);
if (x > 5e-16)
if (x > 1.0e-15)
{
fprintf(stderr, "ERROR(%s line %d): unit error = %lg for i0=%d, j0=%d, di=%d, dj=%d\n", name, lnum, x, i0, j0, di, dj);
return 1;
@@ -2349,6 +2349,24 @@ fail:
#define CHECK_INVERSE(a,b) CHECK(CheckInverse(#a, #b, a, b))
static int CheckCycle(
const char *cyclename,
astro_rotation_t arot, astro_rotation_t brot, astro_rotation_t crot)
{
int error;
astro_rotation_t xrot;
xrot = Astronomy_CombineRotation(arot, brot);
CHECK_ROTMAT(xrot);
xrot = Astronomy_InverseRotation(xrot);
CHECK(CompareMatrices(cyclename, crot, xrot, 2.0e-15));
error = 0;
fail:
return error;
}
#define CHECK_CYCLE(a,b,c) CHECK(CheckCycle(("CheckCycle: " #a ", " #b ", " #c), a, b, c))
static int Test_RotRoundTrip(void)
{
int error;
@@ -2358,18 +2376,15 @@ static int Test_RotRoundTrip(void)
astro_rotation_t eqj_hor, hor_eqj;
astro_rotation_t eqj_eqd, eqd_eqj;
astro_rotation_t hor_eqd, eqd_hor;
//astro_rotation_t eqd_ecl, ecl_eqd;
//astro_rotation_t hor_ecl, ecl_hor;
astro_rotation_t eqd_ecl, ecl_eqd;
astro_rotation_t hor_ecl, ecl_hor;
time = Astronomy_MakeTime(2067, 5, 30, 14, 45, 0.0);
observer = Astronomy_MakeObserver(+28.0, -82.0, 0.0);
/*
In each round trip, we calculate a forward rotation and a backward rotation.
In each round trip, calculate a forward rotation and a backward rotation.
Verify the two are inverse matrices.
Then verify that combining different sequences of rotations result
in the expected combination.
For example, (EQJ ==> HOR ==> ECL) must be the same matrix as (EQJ ==> ECL).
*/
/* Round trip #1: EQJ <==> EQD. */
@@ -2393,14 +2408,28 @@ static int Test_RotRoundTrip(void)
CHECK_INVERSE(eqd_hor, hor_eqd);
/* Round trip #5: EQD <==> ECL. */
//eqd_ecl = Astronomy_Rotation_EQD_ECL(time, observer);
//ecl_eqd = Astronomy_Rotation_ECL_EQD(time, observer);
//CHECK_INVERSE(eqd_ecl, ecl_eqd);
eqd_ecl = Astronomy_Rotation_EQD_ECL(time);
ecl_eqd = Astronomy_Rotation_ECL_EQD(time);
CHECK_INVERSE(eqd_ecl, ecl_eqd);
/* Round trip #6: HOR <==> ECL. */
//hor_ecl = Astronomy_Rotation_HOR_ECL(time, observer);
//ecl_hor = Astronomy_Rotation_ECL_HOR(time, observer);
//CHECK_INVERSE(hor_ecl, ecl_hor);
hor_ecl = Astronomy_Rotation_HOR_ECL(time, observer);
ecl_hor = Astronomy_Rotation_ECL_HOR(time, observer);
CHECK_INVERSE(hor_ecl, ecl_hor);
/*
Verify that combining different sequences of rotations result
in the expected combination.
For example, (EQJ ==> HOR ==> ECL) must be the same matrix as (EQJ ==> ECL).
Each of these is a "triangle" of relationships between 3 orientations.
There are 4 possible ways to pick 3 orientations from the 4 to form a triangle.
Because we have just proved that each transformation is reversible,
we only need to verify the triangle in one cyclic direction.
*/
CHECK_CYCLE(eqj_ecl, ecl_eqd, eqd_eqj); /* excluded corner = HOR */
CHECK_CYCLE(eqj_hor, hor_ecl, ecl_eqj); /* excluded corner = EQD */
CHECK_CYCLE(eqj_hor, hor_eqd, eqd_eqj); /* excluded corner = ECL */
CHECK_CYCLE(ecl_eqd, eqd_hor, hor_ecl); /* excluded corner = EQJ */
printf("Test_RotRoundTrip: PASS\n");
error = 0;

View File

@@ -4617,6 +4617,110 @@ astro_rotation_t Astronomy_Rotation_EQJ_HOR(astro_time_t time, astro_observer_t
}
/**
* @brief
* Calculates a rotation matrix from equatorial of-date (EQD) to ecliptic J2000 (ECL).
*
* This is one of the family of functions that returns a rotation matrix
* for converting from one orientation to another.
* Source: EQD = equatorial system, using equator of date.
* Target: ECL = ecliptic system, using equator at J2000 epoch.
*
* @param time
* The date and time of the source equator.
*
* @return
* A rotation matrix that converts EQD to ECL.
*/
astro_rotation_t Astronomy_Rotation_EQD_ECL(astro_time_t time)
{
astro_rotation_t eqd_eqj;
astro_rotation_t eqj_ecl;
eqd_eqj = Astronomy_Rotation_EQD_EQJ(time);
eqj_ecl = Astronomy_Rotation_EQJ_ECL();
return Astronomy_CombineRotation(eqd_eqj, eqj_ecl);
}
/**
* @brief
* Calculates a rotation matrix from ecliptic J2000 (ECL) to equatorial of-date (EQD).
*
* This is one of the family of functions that returns a rotation matrix
* for converting from one orientation to another.
* Source: ECL = ecliptic system, using equator at J2000 epoch.
* Target: EQD = equatorial system, using equator of date.
*
* @param time
* The date and time of the desired equator.
*
* @return
* A rotation matrix that converts ECL to EQD.
*/
astro_rotation_t Astronomy_Rotation_ECL_EQD(astro_time_t time)
{
astro_rotation_t rot = Astronomy_Rotation_EQD_ECL(time);
return Astronomy_InverseRotation(rot);
}
/**
* @brief
* Calculates a rotation matrix from ecliptic J2000 (ECL) to horizontal (HOR).
*
* This is one of the family of functions that returns a rotation matrix
* for converting from one orientation to another.
* Source: ECL = ecliptic system, using equator at J2000 epoch.
* Target: HOR = horizontal system.
*
* Use #Astronomy_HorizonFromVector to convert the return value
* to a traditional altitude/azimuth pair.
*
* @param time
* The date and time of the desired horizontal orientation.
*
* @param observer
* A location near the Earth's mean sea level that defines the observer's horizon.
*
* @return
* A rotation matrix that converts ECL to HOR at `time` and for `observer`.
* The components of the horizontal vector are:
* x = north, y = west, z = zenith (straight up from the observer).
* These components are chosen so that the "right-hand rule" works for the vector
* and so that north represents the direction where azimuth = 0.
*/
astro_rotation_t Astronomy_Rotation_ECL_HOR(astro_time_t time, astro_observer_t observer)
{
astro_rotation_t ecl_eqd = Astronomy_Rotation_ECL_EQD(time);
astro_rotation_t eqd_hor = Astronomy_Rotation_EQD_HOR(time, observer);
return Astronomy_CombineRotation(ecl_eqd, eqd_hor);
}
/**
* @brief
* Calculates a rotation matrix from horizontal (HOR) to ecliptic J2000 (ECL).
*
* This is one of the family of functions that returns a rotation matrix
* for converting from one orientation to another.
* Source: HOR = horizontal system.
* Target: ECL = ecliptic system, using equator at J2000 epoch.
*
* @param time
* The date and time of the horizontal observation.
*
* @param observer
* The location of the horizontal observer.
*
* @return
* A rotation matrix that converts HOR to ECL.
*/
astro_rotation_t Astronomy_Rotation_HOR_ECL(astro_time_t time, astro_observer_t observer)
{
astro_rotation_t rot = Astronomy_Rotation_ECL_HOR(time, observer);
return Astronomy_InverseRotation(rot);
}
#ifdef __cplusplus
}
#endif

View File

@@ -760,6 +760,30 @@ This function transforms a vector in one orientation to a vector in another orie
---
<a name="Astronomy_Rotation_ECL_EQD"></a>
### Astronomy_Rotation_ECL_EQD(time) &#8658; [`astro_rotation_t`](#astro_rotation_t)
**Calculates a rotation matrix from ecliptic J2000 (ECL) to equatorial of-date (EQD).**
This is one of the family of functions that returns a rotation matrix for converting from one orientation to another. Source: ECL = ecliptic system, using equator at J2000 epoch. Target: EQD = equatorial system, using equator of date.
**Returns:** A rotation matrix that converts ECL to EQD.
| Type | Parameter | Description |
| --- | --- | --- |
| [`astro_time_t`](#astro_time_t) | `time` | The date and time of the desired equator. |
---
<a name="Astronomy_Rotation_ECL_EQJ"></a>
@@ -777,6 +801,57 @@ This is one of the family of functions that returns a rotation matrix for conver
---
<a name="Astronomy_Rotation_ECL_HOR"></a>
### Astronomy_Rotation_ECL_HOR(time, observer) &#8658; [`astro_rotation_t`](#astro_rotation_t)
**Calculates a rotation matrix from ecliptic J2000 (ECL) to horizontal (HOR).**
This is one of the family of functions that returns a rotation matrix for converting from one orientation to another. Source: ECL = ecliptic system, using equator at J2000 epoch. Target: HOR = horizontal system.
Use [`Astronomy_HorizonFromVector`](#Astronomy_HorizonFromVector) to convert the return value to a traditional altitude/azimuth pair.
**Returns:** A rotation matrix that converts ECL to HOR at `time` and for `observer`. The components of the horizontal vector are: x = north, y = west, z = zenith (straight up from the observer). These components are chosen so that the "right-hand rule" works for the vector and so that north represents the direction where azimuth = 0.
| Type | Parameter | Description |
| --- | --- | --- |
| [`astro_time_t`](#astro_time_t) | `time` | The date and time of the desired horizontal orientation. |
| [`astro_observer_t`](#astro_observer_t) | `observer` | A location near the Earth's mean sea level that defines the observer's horizon. |
---
<a name="Astronomy_Rotation_EQD_ECL"></a>
### Astronomy_Rotation_EQD_ECL(time) &#8658; [`astro_rotation_t`](#astro_rotation_t)
**Calculates a rotation matrix from equatorial of-date (EQD) to ecliptic J2000 (ECL).**
This is one of the family of functions that returns a rotation matrix for converting from one orientation to another. Source: EQD = equatorial system, using equator of date. Target: ECL = ecliptic system, using equator at J2000 epoch.
**Returns:** A rotation matrix that converts EQD to ECL.
| Type | Parameter | Description |
| --- | --- | --- |
| [`astro_time_t`](#astro_time_t) | `time` | The date and time of the source equator. |
---
<a name="Astronomy_Rotation_EQD_EQJ"></a>
@@ -896,6 +971,31 @@ Use [`Astronomy_HorizonFromVector`](#Astronomy_HorizonFromVector) to convert the
---
<a name="Astronomy_Rotation_HOR_ECL"></a>
### Astronomy_Rotation_HOR_ECL(time, observer) &#8658; [`astro_rotation_t`](#astro_rotation_t)
**Calculates a rotation matrix from horizontal (HOR) to ecliptic J2000 (ECL).**
This is one of the family of functions that returns a rotation matrix for converting from one orientation to another. Source: HOR = horizontal system. Target: ECL = ecliptic system, using equator at J2000 epoch.
**Returns:** A rotation matrix that converts HOR to ECL.
| Type | Parameter | Description |
| --- | --- | --- |
| [`astro_time_t`](#astro_time_t) | `time` | The date and time of the horizontal observation. |
| [`astro_observer_t`](#astro_observer_t) | `observer` | The location of the horizontal observer. |
---
<a name="Astronomy_Rotation_HOR_EQD"></a>

View File

@@ -5856,6 +5856,110 @@ astro_rotation_t Astronomy_Rotation_EQJ_HOR(astro_time_t time, astro_observer_t
}
/**
* @brief
* Calculates a rotation matrix from equatorial of-date (EQD) to ecliptic J2000 (ECL).
*
* This is one of the family of functions that returns a rotation matrix
* for converting from one orientation to another.
* Source: EQD = equatorial system, using equator of date.
* Target: ECL = ecliptic system, using equator at J2000 epoch.
*
* @param time
* The date and time of the source equator.
*
* @return
* A rotation matrix that converts EQD to ECL.
*/
astro_rotation_t Astronomy_Rotation_EQD_ECL(astro_time_t time)
{
astro_rotation_t eqd_eqj;
astro_rotation_t eqj_ecl;
eqd_eqj = Astronomy_Rotation_EQD_EQJ(time);
eqj_ecl = Astronomy_Rotation_EQJ_ECL();
return Astronomy_CombineRotation(eqd_eqj, eqj_ecl);
}
/**
* @brief
* Calculates a rotation matrix from ecliptic J2000 (ECL) to equatorial of-date (EQD).
*
* This is one of the family of functions that returns a rotation matrix
* for converting from one orientation to another.
* Source: ECL = ecliptic system, using equator at J2000 epoch.
* Target: EQD = equatorial system, using equator of date.
*
* @param time
* The date and time of the desired equator.
*
* @return
* A rotation matrix that converts ECL to EQD.
*/
astro_rotation_t Astronomy_Rotation_ECL_EQD(astro_time_t time)
{
astro_rotation_t rot = Astronomy_Rotation_EQD_ECL(time);
return Astronomy_InverseRotation(rot);
}
/**
* @brief
* Calculates a rotation matrix from ecliptic J2000 (ECL) to horizontal (HOR).
*
* This is one of the family of functions that returns a rotation matrix
* for converting from one orientation to another.
* Source: ECL = ecliptic system, using equator at J2000 epoch.
* Target: HOR = horizontal system.
*
* Use #Astronomy_HorizonFromVector to convert the return value
* to a traditional altitude/azimuth pair.
*
* @param time
* The date and time of the desired horizontal orientation.
*
* @param observer
* A location near the Earth's mean sea level that defines the observer's horizon.
*
* @return
* A rotation matrix that converts ECL to HOR at `time` and for `observer`.
* The components of the horizontal vector are:
* x = north, y = west, z = zenith (straight up from the observer).
* These components are chosen so that the "right-hand rule" works for the vector
* and so that north represents the direction where azimuth = 0.
*/
astro_rotation_t Astronomy_Rotation_ECL_HOR(astro_time_t time, astro_observer_t observer)
{
astro_rotation_t ecl_eqd = Astronomy_Rotation_ECL_EQD(time);
astro_rotation_t eqd_hor = Astronomy_Rotation_EQD_HOR(time, observer);
return Astronomy_CombineRotation(ecl_eqd, eqd_hor);
}
/**
* @brief
* Calculates a rotation matrix from horizontal (HOR) to ecliptic J2000 (ECL).
*
* This is one of the family of functions that returns a rotation matrix
* for converting from one orientation to another.
* Source: HOR = horizontal system.
* Target: ECL = ecliptic system, using equator at J2000 epoch.
*
* @param time
* The date and time of the horizontal observation.
*
* @param observer
* The location of the horizontal observer.
*
* @return
* A rotation matrix that converts HOR to ECL.
*/
astro_rotation_t Astronomy_Rotation_HOR_ECL(astro_time_t time, astro_observer_t observer)
{
astro_rotation_t rot = Astronomy_Rotation_ECL_HOR(time, observer);
return Astronomy_InverseRotation(rot);
}
#ifdef __cplusplus
}
#endif

View File

@@ -620,17 +620,17 @@ astro_spherical_t Astronomy_HorizonFromVector(astro_vector_t vector, astro_refra
astro_vector_t Astronomy_RotateVector(astro_rotation_t rotation, astro_vector_t vector);
astro_rotation_t Astronomy_Rotation_EQD_EQJ(astro_time_t time);
//astro_rotation_t Astronomy_Rotation_EQD_ECL(astro_time_t time);
astro_rotation_t Astronomy_Rotation_EQD_ECL(astro_time_t time);
astro_rotation_t Astronomy_Rotation_EQD_HOR(astro_time_t time, astro_observer_t observer);
astro_rotation_t Astronomy_Rotation_EQJ_EQD(astro_time_t time);
astro_rotation_t Astronomy_Rotation_EQJ_ECL(void);
astro_rotation_t Astronomy_Rotation_EQJ_HOR(astro_time_t time, astro_observer_t observer);
//astro_rotation_t Astronomy_Rotation_ECL_EQD(astro_time_t time);
astro_rotation_t Astronomy_Rotation_ECL_EQD(astro_time_t time);
astro_rotation_t Astronomy_Rotation_ECL_EQJ(void);
//astro_rotation_t Astronomy_Rotation_ECL_HOR(astro_time_t time, astro_observer_t observer);
astro_rotation_t Astronomy_Rotation_ECL_HOR(astro_time_t time, astro_observer_t observer);
astro_rotation_t Astronomy_Rotation_HOR_EQD(astro_time_t time, astro_observer_t observer);
astro_rotation_t Astronomy_Rotation_HOR_EQJ(astro_time_t time, astro_observer_t observer);
//astro_rotation_t Astronomy_Rotation_HOR_ECL(astro_time_t time, astro_observer_t observer);
astro_rotation_t Astronomy_Rotation_HOR_ECL(astro_time_t time, astro_observer_t observer);
#ifdef __cplusplus
}