Eigenized reference frames.
Chris Laurel
parent
f81650c69a
commit
5a241fe28c
|
|
@ -361,8 +361,7 @@ Vector3d
|
|||
VelocityVectorArrow::getDirection(double tdb) const
|
||||
{
|
||||
const TimelinePhase* phase = body.getTimeline()->findPhase(tdb);
|
||||
return toEigen(phase->orbit()->velocityAtTime(tdb) *
|
||||
phase->orbitFrame()->getOrientation(tdb).toMatrix3());
|
||||
return phase->orbitFrame()->getOrientation(tdb).conjugate() * toEigen(phase->orbit()->velocityAtTime(tdb));
|
||||
}
|
||||
|
||||
|
||||
|
|
@ -414,12 +413,7 @@ Vector3d
|
|||
SpinVectorArrow::getDirection(double tdb) const
|
||||
{
|
||||
const TimelinePhase* phase = body.getTimeline()->findPhase(tdb);
|
||||
return toEigen(phase->rotationModel()->angularVelocityAtTime(tdb) *
|
||||
phase->bodyFrame()->getOrientation(tdb).toMatrix3());
|
||||
#if 0
|
||||
return body.getRotationModel(tdb)->angularVelocityAtTime(tdb) *
|
||||
body.getEclipticToFrame(tdb).toMatrix3();
|
||||
#endif
|
||||
return phase->bodyFrame()->getOrientation(tdb).conjugate() * toEigen(phase->rotationModel()->angularVelocityAtTime(tdb));
|
||||
}
|
||||
|
||||
|
||||
|
|
|
|||
|
|
@ -476,12 +476,12 @@ UniversalCoord Body::getPosition(double tdb) const
|
|||
while (frame->getCenter().getType() == Selection::Type_Body)
|
||||
{
|
||||
phase = frame->getCenter().body()->timeline->findPhase(tdb);
|
||||
position += toEigen(frame->getOrientation(tdb)).conjugate() * p;
|
||||
position += frame->getOrientation(tdb).conjugate() * p;
|
||||
p = toEigen(phase->orbit()->positionAtTime(tdb));
|
||||
frame = phase->orbitFrame();
|
||||
}
|
||||
|
||||
position += toEigen(frame->getOrientation(tdb)).conjugate() * p;
|
||||
position += frame->getOrientation(tdb).conjugate() * p;
|
||||
|
||||
#if 0
|
||||
if (frame->getCenter().star())
|
||||
|
|
@ -504,8 +504,8 @@ Quaterniond Body::getOrientation(double tdb) const
|
|||
// different things. The other method should be renamed to
|
||||
// getModelOrientation().
|
||||
const TimelinePhase* phase = timeline->findPhase(tdb);
|
||||
return toEigen(phase->rotationModel()->orientationAtTime(tdb) *
|
||||
phase->bodyFrame()->getOrientation(tdb));
|
||||
return toEigen(phase->rotationModel()->orientationAtTime(tdb)) *
|
||||
phase->bodyFrame()->getOrientation(tdb);
|
||||
}
|
||||
|
||||
|
||||
|
|
@ -518,7 +518,7 @@ Vector3d Body::getVelocity(double tdb) const
|
|||
ReferenceFrame* orbitFrame = phase->orbitFrame();
|
||||
|
||||
Vector3d v = toEigen(phase->orbit()->velocityAtTime(tdb));
|
||||
v = toEigen(orbitFrame->getOrientation(tdb)).conjugate() * v + orbitFrame->getCenter().getVelocity(tdb);
|
||||
v = orbitFrame->getOrientation(tdb).conjugate() * v + orbitFrame->getCenter().getVelocity(tdb);
|
||||
|
||||
if (!orbitFrame->isInertial())
|
||||
{
|
||||
|
|
@ -528,7 +528,7 @@ Vector3d Body::getVelocity(double tdb) const
|
|||
v += cross(orbitFrame->getAngularVelocity(tdb), r);
|
||||
#endif
|
||||
Vector3d r = Selection(const_cast<Body*>(this)).getPosition(tdb).offsetFromKm(orbitFrame->getCenter().getPosition(tdb));
|
||||
v += toEigen(orbitFrame->getAngularVelocity(tdb)).cross(r);
|
||||
v += orbitFrame->getAngularVelocity(tdb).cross(r);
|
||||
}
|
||||
|
||||
return v;
|
||||
|
|
@ -544,10 +544,10 @@ Vector3d Body::getAngularVelocity(double tdb) const
|
|||
Vector3d v = toEigen(phase->rotationModel()->angularVelocityAtTime(tdb));
|
||||
|
||||
ReferenceFrame* bodyFrame = phase->bodyFrame();
|
||||
v = toEigen(bodyFrame->getOrientation(tdb)).conjugate() * v;
|
||||
v = bodyFrame->getOrientation(tdb).conjugate() * v;
|
||||
if (!bodyFrame->isInertial())
|
||||
{
|
||||
v += toEigen(bodyFrame->getAngularVelocity(tdb));
|
||||
v += bodyFrame->getAngularVelocity(tdb);
|
||||
}
|
||||
|
||||
return v;
|
||||
|
|
@ -563,7 +563,7 @@ Vector3d Body::getAngularVelocity(double tdb) const
|
|||
Matrix4d Body::getLocalToAstrocentric(double tdb) const
|
||||
{
|
||||
const TimelinePhase* phase = timeline->findPhase(tdb);
|
||||
Vector3d p = toEigen(phase->orbitFrame()->convertToAstrocentric(phase->orbit()->positionAtTime(tdb), tdb));
|
||||
Vector3d p = phase->orbitFrame()->convertToAstrocentric(toEigen(phase->orbit()->positionAtTime(tdb)), tdb);
|
||||
return Transform3d(Translation3d(p)).matrix();
|
||||
}
|
||||
|
||||
|
|
@ -574,7 +574,7 @@ Vector3d Body::getAstrocentricPosition(double tdb) const
|
|||
{
|
||||
// TODO: Switch the iterative method used in getPosition
|
||||
const TimelinePhase* phase = timeline->findPhase(tdb);
|
||||
return toEigen(phase->orbitFrame()->convertToAstrocentric(phase->orbit()->positionAtTime(tdb), tdb));
|
||||
return phase->orbitFrame()->convertToAstrocentric(toEigen(phase->orbit()->positionAtTime(tdb)), tdb);
|
||||
}
|
||||
|
||||
|
||||
|
|
@ -583,7 +583,7 @@ Vector3d Body::getAstrocentricPosition(double tdb) const
|
|||
Quaterniond Body::getEclipticToFrame(double tdb) const
|
||||
{
|
||||
const TimelinePhase* phase = timeline->findPhase(tdb);
|
||||
return toEigen(phase->bodyFrame()->getOrientation(tdb));
|
||||
return phase->bodyFrame()->getOrientation(tdb);
|
||||
}
|
||||
|
||||
|
||||
|
|
@ -594,7 +594,7 @@ Quaterniond Body::getEclipticToEquatorial(double tdb) const
|
|||
{
|
||||
const TimelinePhase* phase = timeline->findPhase(tdb);
|
||||
return toEigen(phase->rotationModel()->equatorOrientationAtTime(tdb)) *
|
||||
toEigen(phase->bodyFrame()->getOrientation(tdb));
|
||||
phase->bodyFrame()->getOrientation(tdb);
|
||||
}
|
||||
|
||||
|
||||
|
|
@ -605,7 +605,7 @@ Quaterniond Body::getEclipticToBodyFixed(double tdb) const
|
|||
{
|
||||
const TimelinePhase* phase = timeline->findPhase(tdb);
|
||||
return toEigen(phase->rotationModel()->orientationAtTime(tdb)) *
|
||||
toEigen(phase->bodyFrame()->getOrientation(tdb));
|
||||
phase->bodyFrame()->getOrientation(tdb);
|
||||
}
|
||||
|
||||
|
||||
|
|
|
|||
|
|
@ -53,14 +53,14 @@ ReferenceFrame::release() const
|
|||
|
||||
// High-precision rotation using 64.64 fixed point path. Rotate uc by
|
||||
// the rotation specified by unit quaternion q.
|
||||
static UniversalCoord rotate(const UniversalCoord& uc, const Quatd& q)
|
||||
static UniversalCoord rotate(const UniversalCoord& uc, const Quaterniond& q)
|
||||
{
|
||||
Mat3d r = q.toMatrix3();
|
||||
Matrix3d r = q.toRotationMatrix();
|
||||
UniversalCoord uc1;
|
||||
|
||||
uc1.x = uc.x * BigFix(r[0].x) + uc.y * BigFix(r[1].x) + uc.z * BigFix(r[2].x);
|
||||
uc1.y = uc.x * BigFix(r[0].y) + uc.y * BigFix(r[1].y) + uc.z * BigFix(r[2].y);
|
||||
uc1.z = uc.x * BigFix(r[0].z) + uc.y * BigFix(r[1].z) + uc.z * BigFix(r[2].z);
|
||||
uc1.x = uc.x * BigFix(r(0, 0)) + uc.y * BigFix(r(1, 0)) + uc.z * BigFix(r(2, 0));
|
||||
uc1.y = uc.x * BigFix(r(0, 1)) + uc.y * BigFix(r(1, 1)) + uc.z * BigFix(r(2, 1));
|
||||
uc1.z = uc.x * BigFix(r(0, 2)) + uc.y * BigFix(r(1, 2)) + uc.z * BigFix(r(2, 2));
|
||||
|
||||
return uc1;
|
||||
}
|
||||
|
|
@ -76,14 +76,14 @@ UniversalCoord
|
|||
ReferenceFrame::convertFromUniversal(const UniversalCoord& uc, double tjd) const
|
||||
{
|
||||
UniversalCoord uc1 = uc.difference(centerObject.getPosition(tjd));
|
||||
return rotate(uc1, conjugate(getOrientation(tjd)));
|
||||
return rotate(uc1, getOrientation(tjd).conjugate());
|
||||
}
|
||||
|
||||
|
||||
Quatd
|
||||
ReferenceFrame::convertFromUniversal(const Quatd& q, double tjd) const
|
||||
Quaterniond
|
||||
ReferenceFrame::convertFromUniversal(const Quaterniond& q, double tjd) const
|
||||
{
|
||||
return q * conjugate(getOrientation(tjd));
|
||||
return q * getOrientation(tjd).conjugate();
|
||||
}
|
||||
|
||||
|
||||
|
|
@ -105,53 +105,53 @@ ReferenceFrame::convertToUniversal(const UniversalCoord& uc, double tjd) const
|
|||
}
|
||||
|
||||
|
||||
Quatd
|
||||
ReferenceFrame::convertToUniversal(const Quatd& q, double tjd) const
|
||||
Quaterniond
|
||||
ReferenceFrame::convertToUniversal(const Quaterniond& q, double tjd) const
|
||||
{
|
||||
return q * getOrientation(tjd);
|
||||
}
|
||||
|
||||
|
||||
Point3d
|
||||
ReferenceFrame::convertFromAstrocentric(const Point3d& p, double tjd) const
|
||||
Vector3d
|
||||
ReferenceFrame::convertFromAstrocentric(const Vector3d& p, double tjd) const
|
||||
{
|
||||
if (centerObject.getType() == Selection::Type_Body)
|
||||
{
|
||||
Point3d center = ptFromEigen(centerObject.body()->getAstrocentricPosition(tjd));
|
||||
return Point3d(0.0, 0.0, 0.0) + (p - center) * conjugate(getOrientation(tjd)).toMatrix3();
|
||||
Vector3d center = centerObject.body()->getAstrocentricPosition(tjd);
|
||||
return getOrientation(tjd) * (p - center);
|
||||
}
|
||||
else if (centerObject.getType() == Selection::Type_Star)
|
||||
{
|
||||
return p * conjugate(getOrientation(tjd)).toMatrix3();
|
||||
return getOrientation(tjd) * p;
|
||||
}
|
||||
else
|
||||
{
|
||||
// TODO:
|
||||
// bad if the center object is a galaxy
|
||||
// what about locations?
|
||||
return Point3d(0.0, 0.0, 0.0);
|
||||
return Vector3d::Zero();
|
||||
}
|
||||
}
|
||||
|
||||
|
||||
Point3d
|
||||
ReferenceFrame::convertToAstrocentric(const Point3d& p, double tjd) const
|
||||
Vector3d
|
||||
ReferenceFrame::convertToAstrocentric(const Vector3d& p, double tjd) const
|
||||
{
|
||||
if (centerObject.getType() == Selection::Type_Body)
|
||||
{
|
||||
Point3d center = ptFromEigen(centerObject.body()->getAstrocentricPosition(tjd));
|
||||
return center + p * getOrientation(tjd).toMatrix3();
|
||||
Vector3d center = centerObject.body()->getAstrocentricPosition(tjd);
|
||||
return center + getOrientation(tjd).conjugate() * p;
|
||||
}
|
||||
else if (centerObject.getType() == Selection::Type_Star)
|
||||
{
|
||||
return p * getOrientation(tjd).toMatrix3();
|
||||
return getOrientation(tjd).conjugate() * p;
|
||||
}
|
||||
else
|
||||
{
|
||||
// TODO:
|
||||
// bad if the center object is a galaxy
|
||||
// what about locations?
|
||||
return Point3d(0.0, 0.0, 0.0);
|
||||
return Vector3d::Zero();
|
||||
}
|
||||
}
|
||||
|
||||
|
|
@ -165,19 +165,22 @@ ReferenceFrame::getCenter() const
|
|||
}
|
||||
|
||||
|
||||
Vec3d
|
||||
Vector3d
|
||||
ReferenceFrame::getAngularVelocity(double tjd) const
|
||||
{
|
||||
Quatd q0 = getOrientation(tjd);
|
||||
Quatd q1 = getOrientation(tjd + ANGULAR_VELOCITY_DIFF_DELTA);
|
||||
Quatd dq = ~q0 * q1;
|
||||
Quaterniond q0 = getOrientation(tjd);
|
||||
Quaterniond q1 = getOrientation(tjd + ANGULAR_VELOCITY_DIFF_DELTA);
|
||||
Quaterniond dq = q0.conjugate() * q1;
|
||||
|
||||
if (fabs(dq.w) > 0.99999999)
|
||||
return Vec3d(0.0, 0.0, 0.0);
|
||||
|
||||
Vec3d v(dq.x, dq.y, dq.z);
|
||||
if (std::abs(dq.w()) > 0.99999999)
|
||||
return Vector3d::Zero();
|
||||
else
|
||||
return dq.vec().normalized() * (2.0 * acos(dq.w()) / ANGULAR_VELOCITY_DIFF_DELTA);
|
||||
#if CELVEC
|
||||
Vector3d v(dq.x, dq.y, dq.z);
|
||||
v.normalize();
|
||||
return v * (2.0 * acos(dq.w) / ANGULAR_VELOCITY_DIFF_DELTA);
|
||||
#endif
|
||||
}
|
||||
|
||||
|
||||
|
|
@ -256,10 +259,10 @@ J2000EquatorFrame::J2000EquatorFrame(Selection center) :
|
|||
}
|
||||
|
||||
|
||||
Quatd
|
||||
Quaterniond
|
||||
J2000EquatorFrame::getOrientation(double /* tjd */) const
|
||||
{
|
||||
return Quatd::xrotation(astro::J2000Obliquity);
|
||||
return Quaterniond(AngleAxis<double>(astro::J2000Obliquity, Vector3d::UnitX()));
|
||||
}
|
||||
|
||||
|
||||
|
|
@ -288,7 +291,7 @@ BodyFixedFrame::BodyFixedFrame(Selection center, Selection obj) :
|
|||
}
|
||||
|
||||
|
||||
Quatd
|
||||
Quaterniond
|
||||
BodyFixedFrame::getOrientation(double tjd) const
|
||||
{
|
||||
// Rotation of 180 degrees about the y axis is required
|
||||
|
|
@ -298,26 +301,26 @@ BodyFixedFrame::getOrientation(double tjd) const
|
|||
switch (fixObject.getType())
|
||||
{
|
||||
case Selection::Type_Body:
|
||||
return fromEigen(yrot180 * fixObject.body()->getEclipticToBodyFixed(tjd));
|
||||
return yrot180 * fixObject.body()->getEclipticToBodyFixed(tjd);
|
||||
case Selection::Type_Star:
|
||||
return fromEigen(yrot180 * toEigen(fixObject.star()->getRotationModel()->orientationAtTime(tjd)));
|
||||
return yrot180 * toEigen(fixObject.star()->getRotationModel()->orientationAtTime(tjd));
|
||||
default:
|
||||
return fromEigen(yrot180);
|
||||
return yrot180;
|
||||
}
|
||||
}
|
||||
|
||||
|
||||
Vec3d
|
||||
Vector3d
|
||||
BodyFixedFrame::getAngularVelocity(double tjd) const
|
||||
{
|
||||
switch (fixObject.getType())
|
||||
{
|
||||
case Selection::Type_Body:
|
||||
return fromEigen(fixObject.body()->getAngularVelocity(tjd));
|
||||
return fixObject.body()->getAngularVelocity(tjd);
|
||||
case Selection::Type_Star:
|
||||
return fixObject.star()->getRotationModel()->angularVelocityAtTime(tjd);
|
||||
return toEigen(fixObject.star()->getRotationModel()->angularVelocityAtTime(tjd));
|
||||
default:
|
||||
return Vec3d(0.0, 0.0, 0.0);
|
||||
return Vector3d::Zero();
|
||||
}
|
||||
}
|
||||
|
||||
|
|
@ -370,7 +373,7 @@ BodyMeanEquatorFrame::BodyMeanEquatorFrame(Selection center,
|
|||
}
|
||||
|
||||
|
||||
Quatd
|
||||
Quaterniond
|
||||
BodyMeanEquatorFrame::getOrientation(double tjd) const
|
||||
{
|
||||
double t = isFrozen ? freezeEpoch : tjd;
|
||||
|
|
@ -378,21 +381,21 @@ BodyMeanEquatorFrame::getOrientation(double tjd) const
|
|||
switch (equatorObject.getType())
|
||||
{
|
||||
case Selection::Type_Body:
|
||||
return fromEigen(equatorObject.body()->getEclipticToEquatorial(t));
|
||||
return equatorObject.body()->getEclipticToEquatorial(t);
|
||||
case Selection::Type_Star:
|
||||
return equatorObject.star()->getRotationModel()->equatorOrientationAtTime(t);
|
||||
return toEigen(equatorObject.star()->getRotationModel()->equatorOrientationAtTime(t));
|
||||
default:
|
||||
return Quatd(1.0);
|
||||
return Quaterniond::Identity();
|
||||
}
|
||||
}
|
||||
|
||||
|
||||
Vec3d
|
||||
Vector3d
|
||||
BodyMeanEquatorFrame::getAngularVelocity(double tjd) const
|
||||
{
|
||||
if (isFrozen)
|
||||
{
|
||||
return Vec3d(0.0, 0.0, 0.0);
|
||||
return Vector3d::Zero();
|
||||
}
|
||||
else
|
||||
{
|
||||
|
|
@ -402,7 +405,7 @@ BodyMeanEquatorFrame::getAngularVelocity(double tjd) const
|
|||
}
|
||||
else
|
||||
{
|
||||
return Vec3d(0.0, 0.0, 0.0);
|
||||
return Vector3d::Zero();
|
||||
}
|
||||
}
|
||||
}
|
||||
|
|
@ -457,7 +460,7 @@ BodyMeanEquatorFrame::nestingDepth(unsigned int depth,
|
|||
CachingFrame::CachingFrame(Selection _center) :
|
||||
ReferenceFrame(_center),
|
||||
lastTime(-1.0e50),
|
||||
lastOrientation(1.0),
|
||||
lastOrientation(Quaterniond::Identity()),
|
||||
lastAngularVelocity(0.0, 0.0, 0.0),
|
||||
orientationCacheValid(false),
|
||||
angularVelocityCacheValid(false)
|
||||
|
|
@ -465,7 +468,7 @@ CachingFrame::CachingFrame(Selection _center) :
|
|||
}
|
||||
|
||||
|
||||
Quatd
|
||||
Quaterniond
|
||||
CachingFrame::getOrientation(double tjd) const
|
||||
{
|
||||
if (tjd != lastTime)
|
||||
|
|
@ -485,7 +488,7 @@ CachingFrame::getOrientation(double tjd) const
|
|||
}
|
||||
|
||||
|
||||
Vec3d CachingFrame::getAngularVelocity(double tjd) const
|
||||
Vector3d CachingFrame::getAngularVelocity(double tjd) const
|
||||
{
|
||||
if (tjd != lastTime)
|
||||
{
|
||||
|
|
@ -508,24 +511,26 @@ Vec3d CachingFrame::getAngularVelocity(double tjd) const
|
|||
* radians / Julian day.) The default implementation just
|
||||
* differentiates the orientation.
|
||||
*/
|
||||
Vec3d CachingFrame::computeAngularVelocity(double tjd) const
|
||||
Vector3d CachingFrame::computeAngularVelocity(double tjd) const
|
||||
{
|
||||
Quatd q0 = getOrientation(tjd);
|
||||
Quaterniond q0 = getOrientation(tjd);
|
||||
|
||||
// Call computeOrientation() instead of getOrientation() so that we
|
||||
// don't affect the cached value.
|
||||
// TODO: check the valid ranges of the frame to make sure that
|
||||
// jd+dt is still in range.
|
||||
Quatd q1 = computeOrientation(tjd + ANGULAR_VELOCITY_DIFF_DELTA);
|
||||
Quaterniond q1 = computeOrientation(tjd + ANGULAR_VELOCITY_DIFF_DELTA);
|
||||
|
||||
Quatd dq = ~q0 * q1;
|
||||
Quaterniond dq = q0.conjugate() * q1;
|
||||
|
||||
if (fabs(dq.w) > 0.99999999)
|
||||
return Vec3d(0.0, 0.0, 0.0);
|
||||
|
||||
Vec3d v(dq.x, dq.y, dq.z);
|
||||
v.normalize();
|
||||
return v * (2.0 * acos(dq.w) / ANGULAR_VELOCITY_DIFF_DELTA);
|
||||
if (std::abs(dq.w()) > 0.99999999)
|
||||
{
|
||||
return Vector3d::Zero();
|
||||
}
|
||||
else
|
||||
{
|
||||
return dq.vec().normalized() * (2.0 * acos(dq.w()) / ANGULAR_VELOCITY_DIFF_DELTA);
|
||||
}
|
||||
}
|
||||
|
||||
|
||||
|
|
@ -568,11 +573,11 @@ TwoVectorFrame::TwoVectorFrame(Selection center,
|
|||
}
|
||||
|
||||
|
||||
Quatd
|
||||
Quaterniond
|
||||
TwoVectorFrame::computeOrientation(double tjd) const
|
||||
{
|
||||
Vec3d v0 = primaryVector.direction(tjd);
|
||||
Vec3d v1 = secondaryVector.direction(tjd);
|
||||
Vector3d v0 = primaryVector.direction(tjd);
|
||||
Vector3d v1 = secondaryVector.direction(tjd);
|
||||
|
||||
// TODO: verify that v0 and v1 aren't zero length
|
||||
v0.normalize();
|
||||
|
|
@ -583,16 +588,16 @@ TwoVectorFrame::computeOrientation(double tjd) const
|
|||
if (secondaryAxis < 0)
|
||||
v1 = -v1;
|
||||
|
||||
Vec3d v2 = v0 ^ v1;
|
||||
Vector3d v2 = v0.cross(v1);
|
||||
|
||||
// Check for degenerate case when the primary and secondary vectors
|
||||
// are collinear. A well-chosen two vector frame should never have this
|
||||
// problem.
|
||||
double length = v2.length();
|
||||
double length = v2.norm();
|
||||
if (length < Tolerance)
|
||||
{
|
||||
// Just return identity . . .
|
||||
return Quatd(1.0);
|
||||
return Quaterniond::Identity();
|
||||
}
|
||||
else
|
||||
{
|
||||
|
|
@ -605,20 +610,21 @@ TwoVectorFrame::computeOrientation(double tjd) const
|
|||
rhAxis = 1;
|
||||
bool rhOrder = rhAxis == abs(secondaryAxis);
|
||||
|
||||
#if CELVEC
|
||||
// Set the rotation matrix axes
|
||||
Vec3d v[3];
|
||||
Vector3d v[3];
|
||||
v[abs(primaryAxis) - 1] = v0;
|
||||
|
||||
// Reverse the cross products if the axes are not in right
|
||||
// hand order.
|
||||
if (rhOrder)
|
||||
{
|
||||
v[abs(secondaryAxis) - 1] = v2 ^ v0;
|
||||
v[abs(secondaryAxis) - 1] = v2.cross(v0);
|
||||
v[abs(tertiaryAxis) - 1] = v2;
|
||||
}
|
||||
else
|
||||
{
|
||||
v[abs(secondaryAxis) - 1] = v0 ^ -v2;
|
||||
v[abs(secondaryAxis) - 1] = v0.cross(-v2);
|
||||
v[abs(tertiaryAxis) - 1] = -v2;
|
||||
}
|
||||
|
||||
|
|
@ -626,6 +632,28 @@ TwoVectorFrame::computeOrientation(double tjd) const
|
|||
// method must return the quaternion representation of the
|
||||
// orientation, so convert the rotation matrix to a quaternion now.
|
||||
Quatd q = Quatd::matrixToQuaternion(Mat3d(v[0], v[1], v[2]));
|
||||
#endif
|
||||
// Set the rotation matrix axes
|
||||
Matrix3d m;
|
||||
m.row(abs(primaryAxis) - 1) = v0;
|
||||
|
||||
// Reverse the cross products if the axes are not in right
|
||||
// hand order.
|
||||
if (rhOrder)
|
||||
{
|
||||
m.row(abs(secondaryAxis) - 1) = v2.cross(v0);
|
||||
m.row(abs(tertiaryAxis) - 1) = v2;
|
||||
}
|
||||
else
|
||||
{
|
||||
m.row(abs(secondaryAxis) - 1) = v0.cross(-v2);
|
||||
m.row(abs(tertiaryAxis) - 1) = -v2;
|
||||
}
|
||||
|
||||
// The axes are the rows of a rotation matrix. The getOrientation
|
||||
// method must return the quaternion representation of the
|
||||
// orientation, so convert the rotation matrix to a quaternion now.
|
||||
Quaterniond q(m);
|
||||
|
||||
// A rotation matrix will have a determinant of 1; if the matrix also
|
||||
// includes a reflection, the determinant will be -1, indicating that
|
||||
|
|
@ -745,7 +773,7 @@ FrameVector::createRelativeVelocityVector(const Selection& _observer,
|
|||
|
||||
|
||||
FrameVector
|
||||
FrameVector::createConstantVector(const Vec3d& _vec,
|
||||
FrameVector::createConstantVector(const Vector3d& _vec,
|
||||
const ReferenceFrame* _frame)
|
||||
{
|
||||
FrameVector fv(ConstantVector);
|
||||
|
|
@ -757,7 +785,7 @@ FrameVector::createConstantVector(const Vec3d& _vec,
|
|||
}
|
||||
|
||||
|
||||
Vec3d
|
||||
Vector3d
|
||||
FrameVector::direction(double tjd) const
|
||||
{
|
||||
Vector3d v;
|
||||
|
|
@ -778,9 +806,9 @@ FrameVector::direction(double tjd) const
|
|||
|
||||
case ConstantVector:
|
||||
if (frame == NULL)
|
||||
v = toEigen(vec);
|
||||
v = vec;
|
||||
else
|
||||
v = toEigen(vec * frame->getOrientation(tjd).toMatrix3());
|
||||
v = frame->getOrientation(tjd).conjugate() * vec;
|
||||
break;
|
||||
|
||||
default:
|
||||
|
|
@ -789,7 +817,7 @@ FrameVector::direction(double tjd) const
|
|||
break;
|
||||
}
|
||||
|
||||
return fromEigen(v);
|
||||
return v;
|
||||
}
|
||||
|
||||
|
||||
|
|
|
|||
|
|
@ -1,6 +1,7 @@
|
|||
// frame.h
|
||||
//
|
||||
// Copyright (C) 2003, Chris Laurel <claurel@shatters.net>
|
||||
// Copyright (C) 2003-2009, the Celestia Development Team
|
||||
// Original version by Chris Laurel <claurel@gmail.com>
|
||||
//
|
||||
// This program is free software; you can redistribute it and/or
|
||||
// modify it under the terms of the GNU General Public License
|
||||
|
|
@ -10,10 +11,10 @@
|
|||
#ifndef _CELENGINE_FRAME_H_
|
||||
#define _CELENGINE_FRAME_H_
|
||||
|
||||
#include <celmath/vecmath.h>
|
||||
#include <celmath/quaternion.h>
|
||||
#include <celengine/astro.h>
|
||||
#include <celengine/selection.h>
|
||||
#include <Eigen/Core>
|
||||
#include <Eigen/Geometry>
|
||||
|
||||
|
||||
/*! A ReferenceFrame object has a center and set of orthogonal axes.
|
||||
|
|
@ -33,16 +34,16 @@ class ReferenceFrame
|
|||
|
||||
UniversalCoord convertFromUniversal(const UniversalCoord& uc, double tjd) const;
|
||||
UniversalCoord convertToUniversal(const UniversalCoord& uc, double tjd) const;
|
||||
Quatd convertFromUniversal(const Quatd& q, double tjd) const;
|
||||
Quatd convertToUniversal(const Quatd& q, double tjd) const;
|
||||
Eigen::Quaterniond convertFromUniversal(const Eigen::Quaterniond& q, double tjd) const;
|
||||
Eigen::Quaterniond convertToUniversal(const Eigen::Quaterniond& q, double tjd) const;
|
||||
|
||||
Point3d convertFromAstrocentric(const Point3d& p, double tjd) const;
|
||||
Point3d convertToAstrocentric(const Point3d& p, double tjd) const;
|
||||
Eigen::Vector3d convertFromAstrocentric(const Eigen::Vector3d& p, double tjd) const;
|
||||
Eigen::Vector3d convertToAstrocentric(const Eigen::Vector3d& p, double tjd) const;
|
||||
|
||||
Selection getCenter() const;
|
||||
|
||||
virtual Quatd getOrientation(double tjd) const = 0;
|
||||
virtual Vec3d getAngularVelocity(double tdb) const;
|
||||
virtual Eigen::Quaterniond getOrientation(double tjd) const = 0;
|
||||
virtual Eigen::Vector3d getAngularVelocity(double tdb) const;
|
||||
|
||||
virtual bool isInertial() const = 0;
|
||||
|
||||
|
|
@ -73,15 +74,15 @@ class CachingFrame : public ReferenceFrame
|
|||
CachingFrame(Selection _center);
|
||||
virtual ~CachingFrame() {};
|
||||
|
||||
Quatd getOrientation(double tjd) const;
|
||||
Vec3d getAngularVelocity(double tjd) const;
|
||||
virtual Quatd computeOrientation(double tjd) const = 0;
|
||||
virtual Vec3d computeAngularVelocity(double tjd) const;
|
||||
Eigen::Quaterniond getOrientation(double tjd) const;
|
||||
Eigen::Vector3d getAngularVelocity(double tjd) const;
|
||||
virtual Eigen::Quaterniond computeOrientation(double tjd) const = 0;
|
||||
virtual Eigen::Vector3d computeAngularVelocity(double tjd) const;
|
||||
|
||||
private:
|
||||
mutable double lastTime;
|
||||
mutable Quatd lastOrientation;
|
||||
mutable Vec3d lastAngularVelocity;
|
||||
mutable Eigen::Quaterniond lastOrientation;
|
||||
mutable Eigen::Vector3d lastAngularVelocity;
|
||||
mutable bool orientationCacheValid;
|
||||
mutable bool angularVelocityCacheValid;
|
||||
};
|
||||
|
|
@ -94,9 +95,9 @@ class J2000EclipticFrame : public ReferenceFrame
|
|||
J2000EclipticFrame(Selection center);
|
||||
virtual ~J2000EclipticFrame() {};
|
||||
|
||||
Quatd getOrientation(double /* tjd */) const
|
||||
Eigen::Quaterniond getOrientation(double /* tjd */) const
|
||||
{
|
||||
return Quatd(1.0);
|
||||
return Eigen::Quaterniond::Identity();
|
||||
}
|
||||
|
||||
virtual bool isInertial() const;
|
||||
|
|
@ -113,7 +114,7 @@ class J2000EquatorFrame : public ReferenceFrame
|
|||
public:
|
||||
J2000EquatorFrame(Selection center);
|
||||
virtual ~J2000EquatorFrame() {};
|
||||
Quatd getOrientation(double tjd) const;
|
||||
Eigen::Quaterniond getOrientation(double tjd) const;
|
||||
virtual bool isInertial() const;
|
||||
virtual unsigned int nestingDepth(unsigned int depth,
|
||||
unsigned int maxDepth,
|
||||
|
|
@ -132,8 +133,8 @@ class BodyFixedFrame : public ReferenceFrame
|
|||
public:
|
||||
BodyFixedFrame(Selection center, Selection obj);
|
||||
virtual ~BodyFixedFrame() {};
|
||||
Quatd getOrientation(double tjd) const;
|
||||
virtual Vec3d getAngularVelocity(double tjd) const;
|
||||
Eigen::Quaterniond getOrientation(double tjd) const;
|
||||
virtual Eigen::Vector3d getAngularVelocity(double tjd) const;
|
||||
virtual bool isInertial() const;
|
||||
virtual unsigned int nestingDepth(unsigned int depth,
|
||||
unsigned int maxDepth,
|
||||
|
|
@ -150,8 +151,8 @@ class BodyMeanEquatorFrame : public ReferenceFrame
|
|||
BodyMeanEquatorFrame(Selection center, Selection obj, double freeze);
|
||||
BodyMeanEquatorFrame(Selection center, Selection obj);
|
||||
virtual ~BodyMeanEquatorFrame() {};
|
||||
Quatd getOrientation(double tjd) const;
|
||||
virtual Vec3d getAngularVelocity(double tjd) const;
|
||||
Eigen::Quaterniond getOrientation(double tjd) const;
|
||||
virtual Eigen::Vector3d getAngularVelocity(double tjd) const;
|
||||
virtual bool isInertial() const;
|
||||
virtual unsigned int nestingDepth(unsigned int depth,
|
||||
unsigned int maxDepth,
|
||||
|
|
@ -169,11 +170,13 @@ class BodyMeanEquatorFrame : public ReferenceFrame
|
|||
class FrameVector
|
||||
{
|
||||
public:
|
||||
EIGEN_MAKE_ALIGNED_OPERATOR_NEW
|
||||
|
||||
FrameVector(const FrameVector& fv);
|
||||
~FrameVector();
|
||||
FrameVector& operator=(const FrameVector&);
|
||||
|
||||
Vec3d direction(double tjd) const;
|
||||
Eigen::Vector3d direction(double tjd) const;
|
||||
|
||||
/*! Frames can be defined in reference to other frames; this method
|
||||
* counts the depth of such nesting, up to some specified maximum
|
||||
|
|
@ -193,7 +196,7 @@ class FrameVector
|
|||
const Selection& _target);
|
||||
static FrameVector createRelativeVelocityVector(const Selection& _observer,
|
||||
const Selection& _target);
|
||||
static FrameVector createConstantVector(const Vec3d& _vec,
|
||||
static FrameVector createConstantVector(const Eigen::Vector3d& _vec,
|
||||
const ReferenceFrame* _frame);
|
||||
|
||||
private:
|
||||
|
|
@ -205,7 +208,7 @@ class FrameVector
|
|||
FrameVectorType vecType;
|
||||
Selection observer;
|
||||
Selection target;
|
||||
Vec3d vec; // constant vector
|
||||
Eigen::Vector3d vec; // constant vector
|
||||
const ReferenceFrame* frame; // frame for constant vector
|
||||
};
|
||||
|
||||
|
|
@ -231,7 +234,7 @@ class TwoVectorFrame : public CachingFrame
|
|||
int secAxis);
|
||||
virtual ~TwoVectorFrame() {};
|
||||
|
||||
Quatd computeOrientation(double tjd) const;
|
||||
Eigen::Quaterniond computeOrientation(double tjd) const;
|
||||
virtual bool isInertial() const;
|
||||
virtual unsigned int nestingDepth(unsigned int depth,
|
||||
unsigned int maxDepth,
|
||||
|
|
@ -248,5 +251,4 @@ class TwoVectorFrame : public CachingFrame
|
|||
int tertiaryAxis;
|
||||
};
|
||||
|
||||
|
||||
#endif // _CELENGINE_FRAME_H_
|
||||
|
|
|
|||
|
|
@ -16,6 +16,7 @@
|
|||
static const double maximumSimTime = 730486721060.00073; // 2000000000 Jan 01 12:00:00 UTC
|
||||
static const double minimumSimTime = -730498278941.99951; // -2000000000 Jan 01 12:00:00 UTC
|
||||
|
||||
using namespace Eigen;
|
||||
using namespace std;
|
||||
|
||||
|
||||
|
|
@ -603,7 +604,7 @@ void Observer::computeGotoParameters(const Selection& destination,
|
|||
else
|
||||
{
|
||||
ObserverFrame offsetFrame(offsetCoordSys, destination);
|
||||
offset = offset * offsetFrame.getFrame()->getOrientation(getTime()).toMatrix3();
|
||||
offset = offset * fromEigen(offsetFrame.getFrame()->getOrientation(getTime())).toMatrix3();
|
||||
}
|
||||
jparams.to = targetPosition + offset;
|
||||
|
||||
|
|
@ -615,7 +616,7 @@ void Observer::computeGotoParameters(const Selection& destination,
|
|||
else
|
||||
{
|
||||
ObserverFrame upFrame(upCoordSys, destination);
|
||||
upd = upd * upFrame.getFrame()->getOrientation(getTime()).toMatrix3();
|
||||
upd = upd * fromEigen(upFrame.getFrame()->getOrientation(getTime())).toMatrix3();
|
||||
}
|
||||
|
||||
jparams.initialOrientation = getOrientation();
|
||||
|
|
@ -667,7 +668,7 @@ void Observer::computeGotoParametersGC(const Selection& destination,
|
|||
jparams.from = getPosition();
|
||||
|
||||
ObserverFrame offsetFrame(offsetCoordSys, destination);
|
||||
offset = offset * offsetFrame.getFrame()->getOrientation(getTime()).toMatrix3();
|
||||
offset = offset * fromEigen(offsetFrame.getFrame()->getOrientation(getTime())).toMatrix3();
|
||||
|
||||
jparams.to = targetPosition + offset;
|
||||
|
||||
|
|
@ -679,7 +680,7 @@ void Observer::computeGotoParametersGC(const Selection& destination,
|
|||
else
|
||||
{
|
||||
ObserverFrame upFrame(upCoordSys, destination);
|
||||
upd = upd * upFrame.getFrame()->getOrientation(getTime()).toMatrix3();
|
||||
upd = upd * fromEigen(upFrame.getFrame()->getOrientation(getTime())).toMatrix3();
|
||||
}
|
||||
|
||||
jparams.initialOrientation = getOrientation();
|
||||
|
|
@ -1557,14 +1558,14 @@ ObserverFrame::convertToUniversal(const UniversalCoord& uc, double tjd) const
|
|||
Quatd
|
||||
ObserverFrame::convertFromUniversal(const Quatd& q, double tjd) const
|
||||
{
|
||||
return frame->convertFromUniversal(q, tjd);
|
||||
return fromEigen(frame->convertFromUniversal(toEigen(q), tjd));
|
||||
}
|
||||
|
||||
|
||||
Quatd
|
||||
ObserverFrame::convertToUniversal(const Quatd& q, double tjd) const
|
||||
{
|
||||
return frame->convertToUniversal(q, tjd);
|
||||
return fromEigen(frame->convertToUniversal(toEigen(q), tjd));
|
||||
}
|
||||
|
||||
|
||||
|
|
@ -1630,7 +1631,7 @@ ObserverFrame::createFrame(CoordinateSystem _coordSys,
|
|||
|
||||
case PhaseLock_Old:
|
||||
{
|
||||
FrameVector rotAxis(FrameVector::createConstantVector(Vec3d(0, 1, 0),
|
||||
FrameVector rotAxis(FrameVector::createConstantVector(Vector3d::UnitY(),
|
||||
new BodyMeanEquatorFrame(_refObject, _refObject)));
|
||||
return new TwoVectorFrame(_refObject,
|
||||
FrameVector::createRelativePositionVector(_refObject, _targetObject), 3,
|
||||
|
|
@ -1639,7 +1640,7 @@ ObserverFrame::createFrame(CoordinateSystem _coordSys,
|
|||
|
||||
case Chase_Old:
|
||||
{
|
||||
FrameVector rotAxis(FrameVector::createConstantVector(Vec3d(0, 1, 0),
|
||||
FrameVector rotAxis(FrameVector::createConstantVector(Vector3d::UnitY(),
|
||||
new BodyMeanEquatorFrame(_refObject, _refObject)));
|
||||
|
||||
return new TwoVectorFrame(_refObject,
|
||||
|
|
|
|||
|
|
@ -1528,7 +1528,7 @@ CreateFrameVector(const Universe& universe,
|
|||
return NULL;
|
||||
}
|
||||
|
||||
return new FrameVector(FrameVector::createConstantVector(vec, f));
|
||||
return new FrameVector(FrameVector::createConstantVector(toEigen(vec), f));
|
||||
}
|
||||
else
|
||||
{
|
||||
|
|
@ -1652,7 +1652,7 @@ CreateTopocentricFrame(const Selection& center,
|
|||
const Selection& observer)
|
||||
{
|
||||
BodyMeanEquatorFrame* eqFrame = new BodyMeanEquatorFrame(target, target);
|
||||
FrameVector north = FrameVector::createConstantVector(Vec3d(0.0, 1.0, 0.0), eqFrame);
|
||||
FrameVector north = FrameVector::createConstantVector(Vector3d::UnitY(), eqFrame);
|
||||
FrameVector up = FrameVector::createRelativePositionVector(observer, target);
|
||||
|
||||
return new TwoVectorFrame(center, up, -2, north, -3);
|
||||
|
|
|
|||
|
|
@ -2245,7 +2245,7 @@ void Renderer::renderOrbit(const OrbitPathListEntry& orbitPath,
|
|||
Quatd orientation(1.0);
|
||||
if (body != NULL)
|
||||
{
|
||||
orientation = body->getOrbitFrame(t)->getOrientation(t);
|
||||
orientation = fromEigen(body->getOrbitFrame(t)->getOrientation(t));
|
||||
}
|
||||
|
||||
// Equivalent to:
|
||||
|
|
@ -8850,7 +8850,7 @@ void Renderer::buildRenderLists(const Point3d& astrocentricObserverPos,
|
|||
// Get the position of the body relative to the sun.
|
||||
Point3d p = phase->orbit()->positionAtTime(now);
|
||||
ReferenceFrame* frame = phase->orbitFrame();
|
||||
Vec3d pos_s = frameCenter + Vec3d(p.x, p.y, p.z) * frame->getOrientation(now).toMatrix3();
|
||||
Vec3d pos_s = frameCenter + Vec3d(p.x, p.y, p.z) * fromEigen(frame->getOrientation(now)).toMatrix3();
|
||||
|
||||
// We now have the positions of the observer and the planet relative
|
||||
// to the sun. From these, compute the position of the body
|
||||
|
|
@ -9249,7 +9249,7 @@ void Renderer::buildLabelLists(const Frustum& viewFrustum,
|
|||
if (primary != lastPrimary)
|
||||
{
|
||||
Point3d p = phase->orbit()->positionAtTime(now) *
|
||||
phase->orbitFrame()->getOrientation(now).toMatrix3();
|
||||
fromEigen(phase->orbitFrame()->getOrientation(now)).toMatrix3();
|
||||
Vec3d v(iter->position.x - p.x, iter->position.y - p.y, iter->position.z - p.z);
|
||||
|
||||
primarySphere = Sphered(Point3d(v.x, v.y, v.z),
|
||||
|
|
|
|||
Loading…
Reference in New Issue