NavList:
A Community Devoted to the Preservation and Practice of Celestial Navigation and Other Methods of Traditional Wayfinding
Re: Angular Distance Between Stars By Camera and Sextant
From: Marcel Tschudin
Date: 2012 Sep 18, 02:45 +0300
From: Marcel Tschudin
Date: 2012 Sep 18, 02:45 +0300
Thanks, Paul, for your code. May I ask you where you obtained for the stars the high precission epoch 2000 coordinates from ?
Marcel
Marcel
On Tue, Sep 18, 2012 at 2:14 AM, Paul Hirose <cfuhb-acdgw@earthlink.net> wrote:
Greg Rudzinski wrote:My computations say:
The Vega Deneb distance was made at UT 2:25:00 9/16/2012 from Lat. 34* 10.4' N Lon. 119* 13.8' W
The Alioth Alkaid distance was made at UT 3:05:00 Temp. 72* F
Pressure 29.99 in.
48°03.90' +82°47.50' Vega unrefracted az el
48°03.90' +82°47.62' Vega refracted
57°21.56' +59°03.90' Deneb unrefracted az el
57°21.56' +59°04.46' Deneb refracted
23°50.82' Vega to Deneb (unrefracted)
23°50.37' Vega to Deneb (refracted)
310°14.96' +34°06.86' Alkaid unrefracted az el
310°14.96' +34°08.24' Alkaid refracted
320°26.60' +28°18.99' Alioth unrefracted az el
320°26.60' +28°20.73' Alioth refracted
10°27.66' Alkaid to Alioth (unrefracted)
10°27.34' Alkaid to Alioth (refracted)
C# source code to compute those values is included (one file, 20 k). For simplicity, the program does not interact with the user, except to display the output. All input data are embedded in the source code. The desired data are selected with preprocessor directives.
The calculation is precise, with corrections for proper motion,
parallax, radial velocity, aberration (including the diurnal component),
light deflection due to the Sun's gravity, and refraction (including
nonstandard pressure and temperature). In the source code are provisions
to input polar motion and deflection of the vertical, but I've set those
parameters to zero.
Microsoft Visual Studio will compile and run the program. It should also
run under the "Express" (free) version of Visual Studio, but I haven't
tested that.
http://www.microsoft.com/visualstudio/eng/downloads#d-2012-express
The .NET Framework must be installed on your computer. Recent Windows versions already have it. Otherwise, you can download it for free. Finally, the .NET Framework version of my free SofaJpl astronomy DLL must be installed.
http://home.earthlink.net/~s543t-24dst/sofajplNet/index.html
--
// SepAngle.cs
// Windows .NET Framework console app to calculate separation angle
// between stars.
// 2012 Sept 14 by Paul S. Hirose.
// This source code is in the public domain.
// To compile the program, the SofaJpl fundamental astronomy DLL must be
// installed on your system, and a reference to SofaJpl added to the Visual
// Studio project. To run the program, the SofaJpl leap second file is needed.
// It may be stored any convenient place. The DLL and leap second file are
// free at http://home.earthlink.net/~s543t-24dst/sofajplNet/index.html.
using System;
// The SofaJpl namespace is abbreviated to "SJ" for convenience.
using SJ = HirosePS.SofaJpl;
namespace SepAngle
{
class SepAngle
{
static void Main(string[] args)
{
// Load the leap second table.
string leapSecFile = @"G:\Documents and Settings\Hirose limited\" +
@"My Documents\astro\iers\leapSecTable.txt";
SJ.Utc.LoadTableFromFile(leapSecFile);
// date and time (UTC)
#if true
SJ.Utc utc = new SJ.Utc(2012, 9, 16, 2, 25, 0.0);
#else
SJ.Utc utc = new SJ.Utc(2012, 9, 16, 3, 5, 0.0);
#endif
// Delta T, computed from TT-TAI (always 32.184 s), TAI-UTC (35 s,
// until next leap second), and UT1-UTC.
// All quantities from IERS Bulletin A.
SJ.Duration deltaT = new SJ.Duration((32.184 + 35.0 - 0.38988)
/ 86400.0);
// TT and UT1
SJ.JulianDate tt = utc.TerrestrialTime;
SJ.JulianDate ut1 = tt - deltaT;
// polar motion (radians)
double xPole = SJ.Angle.DmsToRad(0, 0, 0);
double yPole = SJ.Angle.DmsToRad(0, 0, 0);
// Observer lat/lon (radians) and height (meters)
double lat = SJ.Angle.DmsToRad(34, 10.4, 0);
double lon = SJ.Angle.DmsToRad(-119, 13.8, 0);
double height = 0.0;
SJ.GeographicPosition observer =
new SJ.GeographicPosition(lat, lon, height);
// deflection of the vertical (radians)
double xi = SJ.Angle.DmsToRad(0, 0, 0);
double eta = SJ.Angle.DmsToRad(0, 0, 0);
// altimeter setting (millibars) and temperature (C)
double mb = 29.99 * SJ.Atmosphere.MillibarsPerInchHg;
double degC = (72.0 - 32.0) / 1.8;
// desired angle accuracy
double accuracy = SJ.Angle.DmsToRad(0, .01, 0);
// Sexagesimal display resolution. 1 = 1 degree, 60 = 1 minute,
// 600 = .1 minute, etc.
double resolution = 6000.0;
// constants
const double secsPerDay = 86400.0;
const double daysPerJY = 365.25;
const double kmPerAU = 1.495978707e8;
// Epoch of the star catalog.
SJ.JulianDate catalogEpoch = SJ.JulianDate.J2000Base;
// Catalog data, first star.
#if true
string star1Name = "Vega";
double star1Ra = SJ.Angle.HmsToRad(18, 36, 56.33635);
double star1Dec = SJ.Angle.DmsToRad(+38, 47, 01.2802);
double star1Parallax = SJ.Angle.DmsToRad(0, 0, 130.23 / 1000.0);
double star1RaPM = SJ.Angle.DmsToRad(0, 0, 200.94 / 1000.0)
/ daysPerJY;
double star1DecPM = SJ.Angle.DmsToRad(0, 0, 286.23 / 1000.0)
/ daysPerJY;
double star1RadVel = -13.9 * secsPerDay / kmPerAU;
#endif
#if false
string star1Name = "Alkaid";
double star1Ra = SJ.Angle.HmsToRad(13, 47, 32.43776);
double star1Dec = SJ.Angle.DmsToRad(+49, 18, 47.7602);
double star1Parallax = SJ.Angle.DmsToRad(0, 0, 31.38 / 1000.0);
double star1RaPM = SJ.Angle.DmsToRad(0, 0, -121.17 / 1000.0)
/ daysPerJY;
double star1DecPM = SJ.Angle.DmsToRad(0, 0, -14.91 / 1000.0)
/ daysPerJY;
double star1RadVel = -10.9 * secsPerDay / kmPerAU;
#endif
// Second star.
#if true
string star2Name = "Deneb";
double star2Ra = SJ.Angle.HmsToRad(20, 41, 25.91514);
double star2Dec = SJ.Angle.DmsToRad(+45, 16, 49.2197);
double star2Parallax = SJ.Angle.DmsToRad(0, 0, 2.31
/ 1000.0);
double star2RaPM = SJ.Angle.DmsToRad(0, 0, 2.01 / 1000.0)
/ daysPerJY;
double star2DecPM = SJ.Angle.DmsToRad(0, 0, 1.85 / 1000.0)
/ daysPerJY;
double star2RadVel = -4.5 * secsPerDay / kmPerAU;
#endif
#if false
string star2Name = "Alioth";
double star2Ra = SJ.Angle.HmsToRad(12, 54, 01.74959);
double star2Dec = SJ.Angle.DmsToRad(+55, 57, 35.3627);
double star2Parallax = SJ.Angle.DmsToRad(0, 0, 39.51 / 1000.0);
double star2RaPM = SJ.Angle.DmsToRad(0, 0, 111.91 / 1000.0)
/ daysPerJY;
double star2DecPM = SJ.Angle.DmsToRad(0, 0, -8.24 / 1000.0)
/ daysPerJY;
double star2RadVel = -9.3 * secsPerDay / kmPerAU;
#endif
// Sources of Sun and Earth position and velocity to correct for
// relativistic light deflection, parallax, aberration. The
// accuracy of the JPL ephemerides is not needed, so use the SOFA
// low accuracy routines.
SJ.SofaEphBody sun = new SJ.SofaEphBody(
SJ.Sofa.Body.Sun, accuracy);
SJ.SofaEphBody earth = new SJ.SofaEphBody(
SJ.Sofa.Body.Earth, accuracy);
// Create two Star objects.
SJ.Star star1 = new SJ.Star(catalogEpoch, star1Ra, star1Dec,
star1Parallax, star1RaPM, star1DecPM, star1RadVel, sun, earth);
SJ.Star star2 = new SJ.Star(catalogEpoch, star2Ra, star2Dec,
star2Parallax, star2RaPM, star2DecPM, star2RadVel, sun, earth);
// Rotation matrix, GCRS to the ITRS.
SJ.RMatrix crsToTrs = SJ.RMatrix.CrsToTrs06a(tt, ut1, xPole, yPole);
// GCRS position and velocity of observer.
SJ.PVVector obsPV = observer.ToPVVector(crsToTrs);
// GCRS vectors to the topocentric apparent places of the stars.
SJ.Vector body1Vec = star1.TopocentricApparent(tt, obsPV);
SJ.Vector body2Vec = star2.TopocentricApparent(tt, obsPV);
// Rotate the coordinate system from the GCRS to the ITRS...
body1Vec = crsToTrs * body1Vec;
body2Vec = crsToTrs * body2Vec;
// ... then from the ITRS to the local horizontal.
SJ.RMatrix trsToHor = SJ.RMatrix.TrsToHor(lat, lon, xi, eta);
body1Vec = trsToHor * body1Vec;
body2Vec = trsToHor * body2Vec;
// Convert coords from vector to spherical.
SJ.Spherical sph1Unref = new SJ.Spherical(body1Vec);
SJ.Spherical sph2Unref = new SJ.Spherical(body2Vec);
// Create an Atmosphere object, use it to refract the altitudes to
// the specified accuracy.
SJ.Atmosphere atmos = new SJ.Atmosphere(height, mb, degC);
double alt1Ref = atmos.Refract(sph1Unref.Alt, accuracy);
double alt2Ref = atmos.Refract(sph2Unref.Alt, accuracy);
// Re-assemble the refracted coordinates into new Spherical objects.
// The "Theta" property must be used because azimuth doesn't obey
// the usual conventions for spherical coordinates.
SJ.Spherical sph1Ref = new SJ.Spherical(sph1Unref.Theta, alt1Ref);
SJ.Spherical sph2Ref = new SJ.Spherical(sph2Unref.Theta, alt2Ref);
// Display unrefracted and refracted az/alt.
Console.WriteLine("{0} {1} unrefracted az alt",
azAltString(sph1Unref, resolution), star1Name);
Console.WriteLine("{0} {1} refracted az alt",
azAltString(sph1Ref, resolution), star1Name);
Console.WriteLine("\n{0} {1} unrefracted az alt",
azAltString(sph2Unref, resolution), star2Name);
Console.WriteLine("{0} {1} refracted az alt",
azAltString(sph2Ref, resolution), star2Name);
// Display unrefracted and refracted separation angles.
Console.WriteLine("\n{0} {1} to {2} (unrefracted)",
sepAngleString(sph1Unref, sph2Unref, resolution),
star1Name, star2Name);
Console.WriteLine("{0} {1} to {2} (refracted)\n",
sepAngleString(sph1Ref, sph2Ref, resolution),
star1Name, star2Name);
}
// Convert a Spherical object to sexagesimal az and alt in a
// string with specified resolution.
static string azAltString(SJ.Spherical sph, double resolution)
{
double az = SJ.Angle.RadiansToDegrees(sph.Az);
double alt = SJ.Angle.RadiansToDegrees(sph.El);
SJ.Sexagesimal azSex = new SJ.Sexagesimal(az, resolution);
SJ.Sexagesimal altSex = new SJ.Sexagesimal(alt, resolution);
return String.Format("{0:3b°′″} {1:+2b°′″}", azSex, altSex);
}
// Compute separation angle between two Spherical objects, convert to
// string in sexagesimal format with specified resolution.
static string sepAngleString(SJ.Spherical sph1, SJ.Spherical sph2,
double resolution)
{
// separation angle
double separation = sph1.SeparationAngle(sph2);
separation = SJ.Angle.RadiansToDegrees(separation);
// Convert degrees to sexagesimal and display.
SJ.Sexagesimal separationSex =
new SJ.Sexagesimal(separation, resolution);
return string.Format("{0:b°′″}", separationSex);
}
}
}