// SofaJplDemo.cs // Demonstration program for SofaJpl DLL v. 2.1.2. // Paul S. Hirose, 2018 August 26 // DOCUMENTATION // The Visual Studio Object Browser displays the main SofaJpl documentation. After adding a // reference to SofaJpl to your project, open Object Browser and expand the HirosePS. SofaJpl // namespace. // revision history #if false 2018-09-07 Added a compatibility mode to duplicate the equatorial, ecliptical, and az/el of JPL Horizons. The first two agree to 1 mas in my tests. Az/el agrees .0001 deg, which is the Horizons output precision. #endif // Select the input time scale. Exactly one must be #defined, and the others commented out. #define UTC // #define UT1 // #define TT // If #define is un-commented, convert a JPL ASCII ephemeris to binary. // After the binary is created, this #define should be commented out. // #define MAKE_EPHEMERIS // If #define is un-commented, use the delta T you provide. Otherwise use the SofaJpl delta T model. // Has no effect if time scale is UTC. // #define MANUAL_DELTA_T // If un-commented, use a precession nutation model compatible with JPL Horizons: IAU 1976 / 80 // precession / nutation, no frame bias, apply IERS pole offsets. // NOTE: To obtain coordinates in the horizontal system or the ITRS, compatible with Horizons, it is // also necessary to activate the MANUAL_DELTA_T option and supply the correct delta T of date. // In addition, for compatible coordinates in the horizontal system set the deflection of the // vertical parameters (xi and eta) to zero. // #define HORIZONS using System; // Everything in SofaJpl is in the HirosePS.SofaJpl namespace. using SofaJpl = HirosePS.SofaJpl; static class SofaJplDemo { // The leap second table. One is included in the SofaJpl distribution .zip file. static string _leapsecTable = @"C:\Users\Stan\Documents\astro\SofaJpl_2_1_2\leapSecTable.txt"; #if MAKE_EPHEMERIS // JPL ASCII ephemeris and header. Must be downloaded from JPL. Not needed after they are // converted to a binary ephemeris. static string _asciiEphemeris = @"C:\Users\Stan\Documents\astro\jpl\ascp2000.422"; static string _ephemerisHeader = @"C:\Users\Stan\Documents\astro\jpl\header.422"; #endif // The binary JPL ephemeris to create and use. static string _binaryEphemeris = @"binp2000_2099.422"; // mathematical constants const double _degPerHour = 15.0; // an enumeration of the three possible angle output formats enum angleFormat { D, DM, DMS }; // output format desired by the user static angleFormat _format; // These control the output resolution of angle and time. The values are computed automatically // from the angle accuracy set by the user. static double _angleResolution; // resolution units per degree static double _timeResolution; // resolution units per hour // Number of decimal places in the floating point 'f' format to achieve unit vector rectangular // coordinate precision comparable to the angles. static int _vectorResolution; static void Main(string[] args) { // Set the display format and angle accuracy. _format = angleFormat.DMS; double angleAccuracy = SofaJpl.Angle.DmsToRad(0, 0, 0.1); // Set the epoch. int year = 2018; int month = 6; int day = 6; int hour = 4; int minute = 0; double second = 0.0; bool julian = false; // true if date is in Julian calendar SofaJpl.Duration deltaT; #if UTC SofaJpl.Duration ut1MinusUtc = SofaJpl.Duration.FromSeconds(0.07522); #else #if MANUAL_DELTA_T deltaT = SofaJpl.Duration.FromSeconds(69.109); #else // Get delta T from the SofaJpl model. deltaT = SofaJpl.DeltaT.GetDeltaT( new SofaJpl.JulianDate(year, month, day, hour, minute, second, julian)); #endif #endif // Observer position. string obsName = "Kitt Peak Observatory"; // name of topocenter (optional) double lon = SofaJpl.Angle.DmsToRad(-111, 36, 00.0); // east longitude double lat = SofaJpl.Angle.DmsToRad(31, 57, 48.0); // north latitude float height = 2100f; // meters above ellipsoid // This instance of the Topocenter class will supply the observer position and velocity // with respect to the ICRS. SofaJpl.Topocenter obs = new SofaJpl.Topocenter(lat, lon, height); // Construct an Atmosphere object to apply refraction. The parameters required by the // constructor are type "float", thus the f suffix on the numbers. float altimiterSetting = 29.9f * SofaJpl.Atmosphere.MillibarsPerInchHg; // millibars float degC = 12f; // Celsius float dewPointC = 0f; // Celsius. SofaJpl.Atmosphere atm = new SofaJpl.Atmosphere(height, degC, altimiterSetting, dewPointC, false); // polar motion (radians) double poleX = SofaJpl.Angle.DmsToRad(0, 0, 0.0); double poleY = SofaJpl.Angle.DmsToRad(0, 0, 0.0); // deflection of the vertical at the observer (radians) double xi = SofaJpl.Angle.DmsToRad(0, 0, 0); double eta = SofaJpl.Angle.DmsToRad(0, 0, 0); // Create a binary ephemeris from JPL ASCII files. #if MAKE_EPHEMERIS Console.WriteLine("\nCREATING BINARY JPL EPHEMERIS.\n"); SofaJpl.JplEphemeris.AsciiToBinary(_ephemerisHeader, _asciiEphemeris, _binaryEphemeris); #endif // Open the binary JPL ephemeris. This loads the entire ephemeris into RAM. Even if the // target body is a star, a solar system ephemeris is required for parallax, aberration, // and light deflection due to solar gravitation. SofaJpl.JplEphemeris eph = new SofaJpl.JplEphemeris(_binaryEphemeris); // There are several ways to specify the target body. // Use the SofaJpl star catalog, a subset of the Hipparcos catalog complete to mag. 3. #if false // Create a HipparcosCatalog object. SofaJpl.HipparcosCatalog catalog = new SofaJpl.HipparcosCatalog(); // Look up the star name in the catalog we just constructed, and construct a Star object. // The current SofaJpl implementation requires a full match (not substring) to a designation // in the name dictionary. It may be necessary to examine the dictionary (star_names.vot) // with a text editor to find the correct name of a star. SofaJpl.Body body = catalog.GetStar("name Vega", eph); #endif // Star, from data you supply. The first parameter is the epoch of the data, J2000.0 in // this case. #if false SofaJpl.Body body = new SofaJpl.Star(SofaJpl.JulianDate.J2000Base, 293.0899579 / 15.0, +69.6611767, 173.77, 597.482, -1738.313, 26.78, eph, "sig Dra"); #endif // Solar system body in the JPL ephemeris. The correction for light time iterates until the // solution converges to the specified angle accuracy. #if true SofaJpl.Body body = new SofaJpl.JplBody(SofaJpl.JplEphemeris.Body.Venus, eph, angleAccuracy); #endif // Finished setting parameters. Calculate and display results. // Compute the angle, time, and rectangular coordinate resolution consistent with the angle // accuracy (radians) specified by user. These variables control output formatting. _angleResolution = 1.0 / SofaJpl.Angle.RadiansToDegrees(angleAccuracy); _timeResolution = _angleResolution * _degPerHour; _vectorResolution = (int)(-Math.Log10(angleAccuracy) + 1.0); if (_vectorResolution < 1) _vectorResolution = 1; // Always display at least 1 decimal place. // Load the leap second table. SofaJpl.Utc.LoadTableFromFile(_leapsecTable); // Calculate UT1 and TT. UTC is not calculated (yet) unless it's the input time scale. SofaJpl.JulianDate ut1, tt; SofaJpl.Utc utc; #if UTC utc = new SofaJpl.Utc(year, month, day, hour, minute, second, julian); tt = utc.TerrestrialTime; ut1 = new SofaJpl.JulianDate(year, month, day, hour, minute, second, julian) + ut1MinusUtc; deltaT = tt - ut1; #endif #if UT1 ut1 = new SofaJpl.JulianDate(year, month, day, hour, minute, second, julian); tt = ut1 + deltaT; #endif #if TT tt = new SofaJpl.JulianDate(year, month, day, hour, minute, second, julian); ut1 = tt - deltaT; #endif // Display date and time in UTC. // The TimeFields class breaks down a JulianDate into year, month, etc. SofaJpl.TimeFields tf1; #if UTC tf1 = utc.ToTimeFields(_timeResolution, julian); Console.WriteLine("{0} UTC", tf1); #else // The time input from the user was UT1 or TT, so it's not safe to assume conversion to UTC // is possible. It must be compared to the UTC table boundaries, which are in terms of TAI. SofaJpl.JulianDate tai = tt - SofaJpl.Duration.TTMinusTai; if (tai >= SofaJpl.Utc.DefaultTable.FirstTai && tai <= SofaJpl.Utc.DefaultTable.LastTai) { utc = new SofaJpl.Utc(tt); tf1 = utc.ToTimeFields(_timeResolution, julian); Console.WriteLine("{0} UTC", tf1); } #endif // Display UT1. tf1 = ut1.ToTimeFields(_timeResolution, julian); Console.WriteLine("{0} UT1", tf1); // Display TT as date and time, and also as 2-part Julian date. tf1 = tt.ToTimeFields(_timeResolution, julian); Console.WriteLine("{0} TT", tf1); Console.WriteLine("Dates are {0} calendar.", julian ? "Julian" : "Gregorian"); Console.WriteLine("JD {0} TT", tt); // Display delta T. SofaJpl.Sexagesimal sex1 = new SofaJpl.Sexagesimal(deltaT.ToHours(), _timeResolution); Console.WriteLine("{0:+fhms} delta T", sex1); // Display polar motion. Console.WriteLine("\npolar motion"); displayPolarMotion(poleX, poleY, "x, y"); // Display topocenter geodetic and rectangular coordinates. Console.WriteLine(); Console.WriteLine(obsName); // name of observatory displayGeodetic(lon, lat, "E lon, N lat"); Console.WriteLine("{0:f1} meters above ellipsoid", height); SofaJpl.Vector obsVec = obs.ToVector(); displayXyz(obsVec, "ITRS unit vector", "modulus (km)"); // Display deflection of the vertical. The same format as polar motion is appropriate. Console.WriteLine("\ndeflection of the vertical"); displayPolarMotion(xi, eta, "xi, eta"); // Display atmosphere conditions. Console.WriteLine("\n{0:f0} C ({1:f0} F) at observer", atm.StationTemperature, atm.StationTemperature * 1.8 + 32.0); Console.WriteLine("{0,6:f1} mb ({1,5:f2}\" Hg) altimeter setting", atm.AltimeterSetting, atm.AltimeterSetting / SofaJpl.Atmosphere.MillibarsPerInchHg); Console.WriteLine("{0,6:f1} mb ({1,5:f2}\" Hg) station pressure", atm.StationPressure, atm.StationPressure / SofaJpl.Atmosphere.MillibarsPerInchHg); if (atm.HumidityIsRelative) Console.WriteLine("{0:f0}% relative humidity", atm.Humidity); else Console.WriteLine("{0:f0} C ({1:f0} F) dew point", atm.Humidity, atm.Humidity * 1.8 + 32.0); // Create rotation matrices for the coordinate transformations. // ITRS to horizontal (east, north, zenith) system, including deflection of the vertical SofaJpl.RMatrix itrstoHor = SofaJpl.RMatrix.ItrsToHor(lat, lon, xi, eta); // polar motion matrix (terrestrial intermediate system to ITRS) SofaJpl.RMatrix tirsToItrs = SofaJpl.RMatrix.TirsToItrs(tt, poleX, poleY); // GCRS to true equator & equinox (IAU 2006 precession and 2000B nutation). SofaJpl.RMatrix gcrsToMean06 = SofaJpl.RMatrix.Precess06(tt); double eps06 = SofaJpl.RMatrix.MeanObliq06(tt); // mean obliquity double dPsi00, dEps00; // nutation in longitude and obliquity SofaJpl.RMatrix.NutationAngles00b(tt, out dPsi00, out dEps00); SofaJpl.RMatrix gcrsToTrue06 = SofaJpl.RMatrix.Nutate(eps06, dPsi00, dEps00) * gcrsToMean06; #if HORIZONS // Generate coordinates compatible with JPL Horizons (1976/80 precession/nutation). // In SofaJpl, 1976 precession includes frame bias, which must be removed for // Horizons compatibility. SofaJpl.RMatrix gcrsToMean = SofaJpl.RMatrix.Precess76(tt) * SofaJpl.RMatrix.IcrsToJ2000.Transpose(); // Calculate the GCRS to ecliptic (mean equinox) rotation matrix. double eps = SofaJpl.RMatrix.MeanObliq80(tt); // mean obl. SofaJpl.RMatrix gcrsToEclipticMean = SofaJpl.RMatrix.GcrsToEclip(gcrsToMean, eps); // GCRS vector to the CIP, based on the IAU 2006 precession & 2000 nutation models. SofaJpl.Vector cipGcrs06 = gcrsToTrue06.Row(3); // Transform it to spherical coords in the 1976/80 ecliptic and mean equinox system. SofaJpl.Vector cipEclip06 = gcrsToEclipticMean * cipGcrs06; SofaJpl.Spherical cipSph = new SofaJpl.Spherical(cipEclip06); // Derive and apply the nutation angles to obtain the GCRS to true equator/equinox matrix. double dPsi, dEps; dPsi = SofaJpl.Angle.HalfPi - cipSph.LonEast; dEps = cipSph.NPD - eps; SofaJpl.RMatrix gcrsToTrue = SofaJpl.RMatrix.Nutate(eps, dPsi, dEps) * gcrsToMean; // Compute the GCRS to terrestrial intermediate matrix double gast = SofaJpl.Angle.Gast94(ut1); SofaJpl.RMatrix gcrsToTirs = SofaJpl.RMatrix.GcrsToTirs(gcrsToTrue, gast); #else // Generate coordinates compatible with IAU 2006/00 precession/nutation. // Mean obliquity and nutation in obliquity have already been computed. double eps = eps06; double dEps = dEps00; // GCRS to true equator/equinox matrix has already been computed. SofaJpl.RMatrix gcrsToTrue = gcrsToTrue06; // Get X and Y of the celestial intermediate pole. double cipX, cipY; gcrsToTrue.CipXY(out cipX, out cipY); // Form the GCRS to celestial intermediate matrix. SofaJpl.RMatrix gcrsToCirs = SofaJpl.RMatrix.GcrsToCirs(cipX, cipY, SofaJpl.RMatrix.S06(tt, cipX, cipY)); // Greenwich apparent sidereal time. double gast = SofaJpl.Angle.Gast06b(ut1, tt); // Earth rotation angle double era = SofaJpl.Angle.Era00(ut1); // Compute the GCRS to terrestrial intermediate matrix SofaJpl.RMatrix gcrsToTirs = SofaJpl.RMatrix.GcrsToTirs(gcrsToCirs, era); #endif // End code compatible with IAU 2006/00 precession/nutation. The remaining rotation // matrix computations are common to IAU 2006/00 and JPL Horizons modes. // GCRS to ecliptic and true equinox SofaJpl.RMatrix gcrsToEcliptic = SofaJpl.RMatrix.GcrsToEclip(gcrsToTrue, eps + dEps); // GCRS to ITRS, including polar motion SofaJpl.RMatrix gcrsToItrs = tirsToItrs * gcrsToTirs; // GCRS to horizontal, including deflection of the vertical. SofaJpl.RMatrix gcrsToHor = itrstoHor * gcrsToItrs; // Display barycentric coordinates of the body. Console.WriteLine("\n{0} barycentric position & velocity", body.Name); SofaJpl.PVVector pv1 = body.Barycentric(tt); displayRaDec(pv1.Position, "RA, dec (ICRS)"); displayXyz(pv1.Position, "unit vector (ICRS)", "distance (km)"); displayXyz(pv1.Velocity, "velocity unit vector", "km/day"); // heliocentric coordinates Console.WriteLine("\n{0} heliocentric position & velocity", body.Name); pv1 = body.Heliocentric(tt); displayRaDec(pv1.Position, "RA, dec (ICRS)"); displayXyz(pv1.Position, "unit vector (ICRS)", "distance (km)"); displayXyz(pv1.Velocity, "velocity unit vector", "km/day"); // geocentric coordinates Console.WriteLine("\n{0} geocentric geometric position & velocity", body.Name); pv1 = body.GeocentricGeometric(tt); displayRaDec(pv1.Position, "RA, dec (ICRS)"); displayXyz(pv1.Position, "unit vector (ICRS)", "geometric distance (km)"); displayXyz(pv1.Velocity, "velocity unit vector", "km/day"); Console.WriteLine("\n{0} geocentric astrometric place", body.Name); SofaJpl.Vector vec1 = body.GeocentricAstrometric(tt); displayRaDec(vec1, "RA, dec (ICRS)"); displayModulus(vec1, "astrometric distance (km)"); Console.WriteLine("\n{0} geocentric apparent place", body.Name); vec1 = body.GeocentricApparent(tt); displayRaDec(vec1, "RA, dec (ICRS)"); displayRaDec(gcrsToTrue * vec1, "equinox RA, dec"); #if ! HORIZONS displayRaDec(gcrsToCirs * vec1, "intermediate RA, dec"); #endif displayEcliptical(gcrsToEcliptic * vec1, "ecliptic true lon, lat"); // Geographic position of the body in the ITRS. If the user provided the parameters, this // includes polar motion. SofaJpl.Spherical sph1 = new SofaJpl.Spherical(gcrsToItrs * vec1); displayGeodetic(sph1.LonEast, sph1.Lat, "E lon, N lat (ITRS)"); displayGeodetic(sph1.LonWest, sph1.Lat, "W lon, N lat (ITRS)"); displayHourAngle(SofaJpl.Angle.NormPlus(sph1.LonWest + lon), "LHA"); // Semidiameter. All bodies have a property to give radius in km. Stars default to zero. // Solar system bodies are automatically initialized to their adopted IAU radii. This // property can be modified by the user. double sd = SofaJpl.Angle.Semidiameter(vec1.Modulus(), body.Radius); displaySemidiameter(sd, "geocentric semidiameter"); // Greenwich apparent sidereal time Console.WriteLine(); displaySiderealTime(gast, "Greenwich apparent sidereal time"); #if ! HORIZONS // Earth rotation angle. displayHourAngle(era, "Earth rotation angle"); #endif // topocentric coordinates // Get the GCRS position and velocity of the topocenter in a PVVector (position and velocity // vector). The parameter passed to ToGcrsPV() must be the ITRS to GCRS rotation matrix, but // what we calculated earlier does the opposite transformation. Thus it's transposed to // reverse its sense. SofaJpl.PVVector obsPV = obs.ToGcrsPV(gcrsToItrs.Transpose()); Console.WriteLine("\n{0} topocentric geometric place", body.Name); vec1 = body.TopocentricGeometric(tt, obsPV).Position; displayRaDec(vec1, "RA, dec (ICRS)"); displayModulus(vec1, "km"); Console.WriteLine("\n{0} topocentric astrometric place", body.Name); vec1 = body.TopocentricAstrometric(tt, obsPV.Position); displayRaDec(vec1, "RA, dec (ICRS)"); displayModulus(vec1, "km"); Console.WriteLine("\n{0} topocentric apparent place", body.Name); vec1 = body.TopocentricApparent(tt, obsPV); displayRaDec(vec1, "RA, dec (ICRS)"); // Transform apparent place from the GCRS to the horizontal system. If the user supplied // the parameters, this includes polar motion and deflection of the vertical. vec1 = gcrsToHor * vec1; // Convert to spherical coordinates. SofaJpl.Spherical sphUnref = new SofaJpl.Spherical(vec1); // Use the Atmosphere object created earlier to generate refracted spherical coordinates. // The application of refraction is iterative, so the desired accuracy must be passed // to Refract(). SofaJpl.Spherical sphRefr = atm.Refract(sphUnref, angleAccuracy); Console.WriteLine("\n{0} azimuth, zenith distance, elevation", body.Name); displayAzZd(sphUnref, "az, unrefracted ZD"); displayAzEl(sphUnref, "az, unrefracted el"); displayAzEl(sphRefr, "az, refracted el"); // Semidiameter. All bodies have a property to give radius in km. Stars default to zero, but // solar system bodies are automatically initialized to their adopted IAU radii. sd = SofaJpl.Angle.Semidiameter(vec1.Modulus(), body.Radius); displaySemidiameter(sd, "topocentric semidiameter"); // If applicable, display phase angle: the separation angle, at the body, between // vectors directed to the Sun and the observer. These are the negatives of vectors to the // body's heliocentric and topocentric positions in the ICRS. If we omit the negations, // both vectors are off by 180 degrees and thus the angle between them is still correct. if (body.IsSolarSystemBody && body.IsSun == false) { SofaJpl.JplBody earth = new SofaJpl.JplBody( SofaJpl.JplEphemeris.Body.Earth, eph, angleAccuracy); SofaJpl.JplBody sun = new SofaJpl.JplBody( SofaJpl.JplEphemeris.Body.Sun, eph, angleAccuracy); SofaJpl.Vector earthVec = body.TopocentricGeometric(tt, obsPV).Position; SofaJpl.Vector sunVec = body.Heliocentric(tt).Position; double pa = earthVec.SeparationAngle(sunVec); Console.WriteLine("{0:f0}° phase angle (0 = full, 180 = new)", SofaJpl.Angle.RadiansToDegrees(pa)); } } // end Main() // helper methods to display data with the format and precision selected by user /// <summary> /// Display the modulus of a vector. /// </summary> /// <param name="vec">the vector</param> /// <param name="label">string to display after the modulus</param> static void displayModulus(SofaJpl.Vector vec, string label) { System.Text.StringBuilder sb1 = new System.Text.StringBuilder("{0:e", 8); // The _vectorResolution value is appropriate for the 'f' format. But the 'e' format gives // one more significant digit, so subtract 1. sb1.Append(_vectorResolution - 1); sb1.Append("} "); Console.WriteLine(sb1.ToString() + label, vec.Modulus()); } /// <summary> /// Display hour angle as a sexagesimal. /// </summary> /// <param name="angle">hour angle (radians) </param> /// <param name="str">string to display after the angle</param> static void displayHourAngle(double angle, string str) { string formatString; switch (_format) { case angleFormat.D: formatString = "{0:3c°'\"} {1}"; break; case angleFormat.DM: formatString = "{0:3b°'\"} {1}"; break; default: formatString = "{0:3a°'\"} {1}"; break; } SofaJpl.Sexagesimal sex = new SofaJpl.Sexagesimal( SofaJpl.Angle.RadiansToDegrees(angle), _angleResolution, 360); Console.WriteLine(formatString, sex, str); } /// <summary> /// Display sidereal time as a sexagesimal. /// </summary> /// <param name="gast">sidereal time (radians)</param> /// <param name="label">string to display after the angle</param> /// <remarks>Display degrees (not hours) unless format is DMS.</remarks> static void displaySiderealTime(double gast, string label) { string formatString; switch (_format) { case angleFormat.D: formatString = "{0:3c°'\"} {1}"; break; case angleFormat.DM: formatString = "{0:3b°'\"} {1}"; break; default: formatString = "{0:2ahms} {1}"; break; } SofaJpl.Sexagesimal sex; if (_format == angleFormat.DMS) sex = new SofaJpl.Sexagesimal(SofaJpl.Angle.RadiansToHours(gast), _timeResolution, 24); else sex = new SofaJpl.Sexagesimal(SofaJpl.Angle.RadiansToDegrees(gast), _angleResolution, 360); Console.WriteLine(formatString, sex, label); } /// <summary> /// Display semidiameter as a sexagesimal. /// </summary> /// <param name="sd">semidiameter (radians)</param> /// <param name="label">string to display after the angle</param> static void displaySemidiameter(double sd, string label) { string formatString; switch (_format) { case angleFormat.D: formatString = "{0:c°'\"} {1}"; break; case angleFormat.DM: formatString = "{0:e°'\"} {1}"; break; default: formatString = "{0:f°'\"} {1}"; break; } sd = SofaJpl.Angle.RadiansToDegrees(sd); SofaJpl.Sexagesimal sex = new SofaJpl.Sexagesimal(sd, _angleResolution); Console.WriteLine(formatString, sex, label); } /// <summary> /// Display polar motion angles as sexagesimals. /// </summary> /// <param name="poleX">pole X (radians)</param> /// <param name="poleY">pole Y (radians)</param> /// <param name="label">string to display after the angles</param> static void displayPolarMotion(double poleX, double poleY, string label) { string formatString; switch (_format) { case angleFormat.D: formatString = "{0:+c°'\"} {1:+c°'\"} {2}"; break; case angleFormat.DM: formatString = "{0:+e°'\"} {1:+e°'\"} {2}"; break; default: formatString = "{0:+f°'\"} {1:+f°'\"} {2}"; break; } double xDeg = SofaJpl.Angle.RadiansToDegrees(poleX); double yDeg = SofaJpl.Angle.RadiansToDegrees(poleY); SofaJpl.Sexagesimal sexX = new SofaJpl.Sexagesimal(xDeg, _angleResolution); SofaJpl.Sexagesimal sexY = new SofaJpl.Sexagesimal(yDeg, _angleResolution); Console.WriteLine(formatString, sexX, sexY, label); } /// <summary> /// Display a vector as right ascension and declination. /// </summary> /// <param name="vec">vector to the body with respect to the equatorial system</param> /// <param name="label">string to display after the angles</param> /// <remarks>Unless the format is DMS, display RA in degrees not hours.</remarks> static void displayRaDec(SofaJpl.Vector vec, string label) { string formatString; switch (_format) { case angleFormat.D: formatString = "{0:3c°'\"} {1:+2c°'\"} {2}"; break; case angleFormat.DM: formatString = "{0:3b°'\"} {1:+2b°'\"} {2}"; break; default: formatString = "{0:2ahms} {1:+2a°'\"} {2}"; break; } SofaJpl.Spherical sph = new SofaJpl.Spherical(vec); SofaJpl.Sexagesimal raSex; if (_format == angleFormat.DMS) raSex = new SofaJpl.Sexagesimal(SofaJpl.Angle.RadiansToHours(sph.RA), _timeResolution, 24); else raSex = new SofaJpl.Sexagesimal(SofaJpl.Angle.RadiansToDegrees(sph.RA), _angleResolution, 360); SofaJpl.Sexagesimal decSex = new SofaJpl.Sexagesimal( SofaJpl.Angle.RadiansToDegrees(sph.Dec), _angleResolution); Console.WriteLine(formatString, raSex, decSex, label); } /// <summary> /// Display a vector as azimuth and elevation. /// </summary> /// <param name="vec">vector to the body with respect to the horizontal system</param> /// <param name="label">string to display after the angles</param> static void displayAzEl(SofaJpl.Vector vec, string label) { string formatString; switch (_format) { case angleFormat.D: formatString = "{0:3c°'\"} {1:+2c°'\"} {2}"; break; case angleFormat.DM: formatString = "{0:3b°'\"} {1:+2b°'\"} {2}"; break; default: formatString = "{0:3a°'\"} {1:+2a°'\"} {2}"; break; } SofaJpl.Spherical sph = new SofaJpl.Spherical(vec); SofaJpl.Sexagesimal azSex = new SofaJpl.Sexagesimal( SofaJpl.Angle.RadiansToDegrees(sph.Az), _angleResolution, 360); SofaJpl.Sexagesimal elSex = new SofaJpl.Sexagesimal( SofaJpl.Angle.RadiansToDegrees(sph.El), _angleResolution); Console.WriteLine(formatString, azSex, elSex, label); } /// <summary> /// Display a vector as ecliptic longitude and latitude. /// </summary> /// <param name="vec">vector to the body with respect to the ecliptic system</param> /// <param name="label">string to display after the angles</param> static void displayEcliptical(SofaJpl.Vector vec, string label) { string formatString; switch (_format) { case angleFormat.D: formatString = "{0:3c°'\"} {1:+2c°'\"} {2}"; break; case angleFormat.DM: formatString = "{0:3b°'\"} {1:+2b°'\"} {2}"; break; default: formatString = "{0:3a°'\"} {1:+2a°'\"} {2}"; break; } SofaJpl.Spherical sph = new SofaJpl.Spherical(vec); SofaJpl.Sexagesimal lonSex = new SofaJpl.Sexagesimal( SofaJpl.Angle.RadiansToDegrees(sph.LonEast), _angleResolution, 360); SofaJpl.Sexagesimal latSex = new SofaJpl.Sexagesimal( SofaJpl.Angle.RadiansToDegrees(sph.Lat), _angleResolution); Console.WriteLine(formatString, lonSex, latSex, label); } /// <summary> /// Display geodetic longitude and latitude. /// </summary> /// <param name="lon">east longitude (radians)</param> /// <param name="lat">north latitude (radians)</param> /// <param name="label">string to display after the angles</param> static void displayGeodetic(double lon, double lat, string label) { string formatString; switch (_format) { case angleFormat.D: formatString = "{0:3c°'\"} {1:+2c°'\"} {2}"; break; case angleFormat.DM: formatString = "{0:3b°'\"} {1:+2b°'\"} {2}"; break; default: formatString = "{0:3a°'\"} {1:+2a°'\"} {2}"; break; } SofaJpl.Sexagesimal lonSex = new SofaJpl.Sexagesimal( SofaJpl.Angle.RadiansToDegrees(lon), _angleResolution, 360); SofaJpl.Sexagesimal latSex = new SofaJpl.Sexagesimal( SofaJpl.Angle.RadiansToDegrees(lat), _angleResolution); Console.WriteLine(formatString, lonSex, latSex, label); } /// <summary> /// Display a vector as xyz components of a unit vector and modulus. /// </summary> /// <param name="v">vector</param> /// <param name="xyzLabel">label for the unit vector</param> /// <param name="modulusLabel">label for the modulus</param> static void displayXyz(SofaJpl.Vector v, string xyzLabel, string modulusLabel) { System.Text.StringBuilder sb1 = new System.Text.StringBuilder("{0:f", 23); sb1.Append(_vectorResolution); sb1.Append("} {1:f"); sb1.Append(_vectorResolution); sb1.Append("} {2:f"); sb1.Append(_vectorResolution); sb1.Append("} "); SofaJpl.Vector uv = v.Unit(); Console.WriteLine(sb1.ToString() + xyzLabel, uv.X, uv.Y, uv.Z); displayModulus(v, modulusLabel); } /// <summary> /// Given a Spherical object, display azimuth and elevation as sexagesimals. /// </summary> /// <param name="sph">spherical coordinates in the horizontal system</param> /// <param name="label">string to display after the angles</param> static void displayAzEl(SofaJpl.Spherical sph, string label) { string formatString; switch (_format) { case angleFormat.D: formatString = "{0:3c°'\"} {1: 2c°'\"} {2}"; break; case angleFormat.DM: formatString = "{0:3b°'\"} {1: 2b°'\"} {2}"; break; default: formatString = "{0:3a°'\"} {1: 2a°'\"} {2}"; break; } double az = SofaJpl.Angle.RadiansToDegrees(sph.Az); double el = SofaJpl.Angle.RadiansToDegrees(sph.El); SofaJpl.Sexagesimal azSex = new SofaJpl.Sexagesimal(az, _angleResolution, 360); SofaJpl.Sexagesimal elSex = new SofaJpl.Sexagesimal(el, _angleResolution); Console.WriteLine(formatString, azSex, elSex, label); } /// <summary> /// Given a Spherical object, display azimuth and zenith distance as sexagesimals. /// </summary> /// <param name="sph">spherical coordinates</param> /// <param name="label">string to display after the angles</param> static void displayAzZd(SofaJpl.Spherical sph, string label) { string formatString; switch (_format) { case angleFormat.D: // Zenith distance needs a different format than elevation. There's never a negative // sign, and there can be up to 3 digits before the decimal point. formatString = "{0:3c°'\"} {1:3c°'\"} {2}"; break; case angleFormat.DM: formatString = "{0:3b°'\"} {1:3b°'\"} {2}"; break; default: formatString = "{0:3a°'\"} {1:3a°'\"} {2}"; break; } double az = SofaJpl.Angle.RadiansToDegrees(sph.Az); double zd = SofaJpl.Angle.RadiansToDegrees(sph.ZenithDistance); SofaJpl.Sexagesimal azSex = new SofaJpl.Sexagesimal(az, _angleResolution, 360); SofaJpl.Sexagesimal zdSex = new SofaJpl.Sexagesimal(zd, _angleResolution); Console.WriteLine(formatString, azSex, zdSex, label); } } // end class
[back]
(last modified 2018-09-09)