Ive created a method to calculate the position of a celestial body in the sky given the observers location on a planet and the time.
Im not sure how accurate it is because I don't know what to trust more (celestia and my calcs - possibly with precision probs) or other sky map programs.
This would be useful to people trying to control telescopes with celestia.
Code: Select all
void Body::getDirectionOfBodyFromSurface(const Vec3f& lonLatAlt,Body * body,double when,
double& Dist,double& Azimuth,double&Altitude)
{
//get the bodies position in geographic coodinates
Point3d bodyPointd = body->getHeliocentricPosition(when)
*getGeographicToHeliocentric(when).inverse();
//get the position of the observer in geographic coords
Vec3f observerVector = planetocentricToCartesian(lonLatAlt);
//create a plane tangential to the earths surface where the
//observer is
Vec3f FloorN = observerVector;
FloorN.normalize();
Planef FloorPlane(FloorN,observerVector.length());
//find a point north of the observer by finding the floor planes
//intersection with the geographical y axis
Point3f northPoint = Planef::intersection(FloorPlane,
Planef(Vec3f(1,0,0),0),Planef(Vec3f(0,0,1),0));
//Do some precision conversions
Point3f bodyPoint(bodyPointd.x,bodyPointd.y,bodyPointd.z);
Vec3f bodyVector(bodyPoint.x,bodyPoint.y,bodyPoint.z);
Vec3f northVector(northPoint.x,northPoint.y,northPoint.z);
//calculate the vector from the observer to north
northVector -= observerVector;
//create a plane that splits the planet north south
//from the point of view of the observer
Vec3f NSPlaneN = northVector^observerVector;
NSPlaneN.normalize();
Planef NorthSouthPlane(NSPlaneN,0);
//Do the same for east west
Vec3f EWPlaneN = northVector;
EWPlaneN.normalize();
Planef EastWestPlane(EWPlaneN,0);
//Calculate the Azimuth
float bodyDistanceFromNSPlane = NorthSouthPlane.distanceTo(bodyPoint);
float bodyDistanceFromEWPlane = EastWestPlane.distanceTo(bodyPoint);
Azimuth = atan(bodyDistanceFromNSPlane/bodyDistanceFromEWPlane);
if (bodyDistanceFromNSPlane > 0) {
if (bodyDistanceFromEWPlane > 0) {
//do nothing
}else{
Azimuth += PI;
}
}else{
if (bodyDistanceFromEWPlane > 0) {
Azimuth = 2*PI + Azimuth;
}else{
Azimuth += PI;
}
}
//Calculate the Altitude
Altitude = asin(FloorPlane.distanceTo(bodyPoint)/Dist);
//work out the distance from the observer to the body
bodyVector -= observerVector;
Dist = bodyVector.length();
Dist -= body->getRadius();
}
Alex