AstroLib is a package of small generic routines useful above all in astronomical and astrophysical context, written in Julia.
Included are also translations of some IDL Astronomy User’s Library procedures, which are released under terms of BSD-2-Clause License. AstroLib’s functions are not drop-in replacement of those procedures, Julia standard data types are often used (e.g., DateTime type instead of generic string for dates) and the syntax may slightly differ.
An extensive error testing suite ensures old fixed bugs will not be brought back by future changes.
AstroLib is available for Julia 1.0 and later versions, and can be installed with Julia's built-in package manager. In a Julia session run the commands
AstroLib is a package of small generic routines useful above all in astronomical and astrophysical context, written in Julia.
Included are also translations of some IDL Astronomy User’s Library procedures, which are released under terms of BSD-2-Clause License. AstroLib’s functions are not drop-in replacement of those procedures, Julia standard data types are often used (e.g., DateTime type instead of generic string for dates) and the syntax may slightly differ.
An extensive error testing suite ensures old fixed bugs will not be brought back by future changes.
AstroLib is available for Julia 1.0 and later versions, and can be installed with Julia's built-in package manager. In a Julia session run the commands
Older versions are also available for Julia 0.4-0.6.
Note that, in order to work, a few functions require external files, which are automatically downloaded when building the package. Should these files be missing for some reason, you will be able to load the package but some functions may not work properly. You can manually build the package with
After installing the package, you can start using AstroLib with
using AstroLib
Many functions in AstroLib.jl are compatible with Measurements.jl package, which allows you to define quantities with uncertainty and propagate the error when performing calculations according to propagation of uncertainty rules. For example:
AstroLib.jl is developed on GitHub. You can contribute to the project in a number of ways: by translating more routines from IDL Astronomy User’s Library, or providing brand-new functions, or even improving existing ones (make them faster and more precise). Also bug reports are encouraged.
This project is a work-in-progress, only few procedures have been translated so far. In addition, function syntax may change from time to time. Check TODO.md out to see how you can help. Volunteers are welcome!
This is not the only effort to bundle astronomical functions written in Julia language. Other packages useful for more specific purposes are available at JuliaAstro.
Because of this, some of IDL AstroLib’s utilities are not provided in AstroLib.jl as they are already present in other Julia packages. Here is a list of such utilities:
readcol, use readdlm, part of Julia Base.DataFmt module. This is not a complete replacement for readcol but most of the time it does-the-right-thing even without using any option (it automatically identifies string and numerical columns) and you do not need to manually specify a variable for each column
AstroLib.jl is developed on GitHub. You can contribute to the project in a number of ways: by translating more routines from IDL Astronomy User’s Library, or providing brand-new functions, or even improving existing ones (make them faster and more precise). Also bug reports are encouraged.
This project is a work-in-progress, only few procedures have been translated so far. In addition, function syntax may change from time to time. Check TODO.md out to see how you can help. Volunteers are welcome!
This is not the only effort to bundle astronomical functions written in Julia language. Other packages useful for more specific purposes are available at JuliaAstro.
Because of this, some of IDL AstroLib’s utilities are not provided in AstroLib.jl as they are already present in other Julia packages. Here is a list of such utilities:
readcol, use readdlm, part of Julia Base.DataFmt module. This is not a complete replacement for readcol but most of the time it does-the-right-thing even without using any option (it automatically identifies string and numerical columns) and you do not need to manually specify a variable for each column
AstroLib.jl defines a new Observatory type. This can be used to define a new object holding information about an observing site. It is a composite type whose fields are
name (String type): the name of the site
latitude (Float64 type): North-ward latitude of the site in degrees
longitude (Float64 type): East-ward longitude of the site in degrees
altitude (Float64 type): altitude of the site in meters
tz (Float64 type): the number of hours of offset from UTC
The type constructor Observatory can be used to create a new Observatory object. Its syntax is
Observatory(name, lat, long, alt, tz)
name should be a string; lat, long, and tz should be anything that can be converted to a floating number with ten function; alt should be a real number.
A predefined list of some observing sites is provided with AstroLib.observatories constant. It is a dictionary whose keys are the abbreviated names of the observatories. For example, you can access information of the European Southern Observatory with
AstroLib.jl defines a new Observatory type. This can be used to define a new object holding information about an observing site. It is a composite type whose fields are
name (String type): the name of the site
latitude (Float64 type): North-ward latitude of the site in degrees
longitude (Float64 type): East-ward longitude of the site in degrees
altitude (Float64 type): altitude of the site in meters
tz (Float64 type): the number of hours of offset from UTC
The type constructor Observatory can be used to create a new Observatory object. Its syntax is
Observatory(name, lat, long, alt, tz)
name should be a string; lat, long, and tz should be anything that can be converted to a floating number with ten function; alt should be a real number.
A predefined list of some observing sites is provided with AstroLib.observatories constant. It is a dictionary whose keys are the abbreviated names of the observatories. For example, you can access information of the European Southern Observatory with
List of planets of the Solar System, from Mercury to Pluto. The elements of the list have Planet type.
Reference for most quantities is the Planetary Fact Sheet: http://nssdc.gsfc.nasa.gov/planetary/factsheet/index.html and the Keplerian Elements for Approximate Positions of the Major Planets: https://ssd.jpl.nasa.gov/txt/pelemt1.txt
List of planets of the Solar System, from Mercury to Pluto. The elements of the list have Planet type.
Reference for most quantities is the Planetary Fact Sheet: http://nssdc.gsfc.nasa.gov/planetary/factsheet/index.html and the Keplerian Elements for Approximate Positions of the Major Planets: https://ssd.jpl.nasa.gov/txt/pelemt1.txt
Returns right ascension and declination as string(s) in sexagesimal format.
Explanation
Takes right ascension and declination expressed in decimal format, converts them to sexagesimal and return a formatted string. The precision of right ascension and declination can be specified.
Arguments
Arguments of this function are:
ra: right ascension in decimal degrees. It is converted to hours before printing.
dec: declination in decimal degrees.
The function can be called in different ways:
Two numeric arguments: first is ra, the second is dec.
An iterable (array, tuple) of two elements: (ra, dec).
One numeric argument: it is assumed only dec is provided.
Optional keywords affecting the output format are always available:
precision (optional integer keyword): specifies the number of digits of declination seconds. The number of digits for right ascension seconds is always assumed to be one more precision. If the function is called with only dec as input, precision default to 1, in any other case defaults to 0.
truncate (optional boolean keyword): if true, then the last displayed digit in the output is truncated in precision rather than rounded. This option is useful if adstring is used to form an official IAU name (see http://vizier.u-strasbg.fr/Dic/iau-spec.htx) with coordinate specification.
Output
The function returns one string. The format of strings can be specified with precision and truncate keywords, see above.
Wavelengths are corrected for the index of refraction of air under standard conditions. Wavelength values below $2000 Å$ will not be altered, take care within $[1 Å, 2000 Å]$. Uses relation of Ciddor (1996).
Arguments
wave_air: the wavelength in air.
Output
Vacuum wavelength in angstroms.
Method
Uses relation of Ciddor (1996), Applied Optics 62, 958.
Example
If the air wavelength is w = 6056.125 (a Krypton line), then airtovac(w) yields a vacuum wavelength of 6057.8019.
Wavelengths are corrected for the index of refraction of air under standard conditions. Wavelength values below $2000 Å$ will not be altered, take care within $[1 Å, 2000 Å]$. Uses relation of Ciddor (1996).
Arguments
wave_air: the wavelength in air.
Output
Vacuum wavelength in angstroms.
Method
Uses relation of Ciddor (1996), Applied Optics 62, 958.
Example
If the air wavelength is w = 6056.125 (a Krypton line), then airtovac(w) yields a vacuum wavelength of 6057.8019.
julia> using AstroLib
julia> airtovac(6056.125)
-6057.801930991426
Notes
vactoair converts vacuum wavelengths to air wavelengths.
Code of this function is based on IDL Astronomy User's Library.
Convert longitude l and latitude b to (x, y) using an Aitoff projection.
Explanation
This function can be used to create an all-sky map in Galactic coordinates with an equal-area Aitoff projection. Output map coordinates are zero longitude centered.
Arguments
l: longitude, scalar or vector, in degrees.
b: latitude, number of elements as l, in degrees.
Coordinates can be given also as a 2-tuple (l, b).
Output
2-tuple (x, y).
x: x coordinate, same number of elements as l. x is normalized to be in $[-180, 180]$.
y: y coordinate, same number of elements as l. y is normalized to be in $[-90, 90]$.
Example
Get $(x ,y)$ Aitoff coordinates of Sirius, whose Galactic coordinates are $(227.23, -8.890)$.
julia> using AstroLib
+6057.801930991426
Notes
vactoair converts vacuum wavelengths to air wavelengths.
Code of this function is based on IDL Astronomy User's Library.
Convert longitude l and latitude b to (x, y) using an Aitoff projection.
Explanation
This function can be used to create an all-sky map in Galactic coordinates with an equal-area Aitoff projection. Output map coordinates are zero longitude centered.
Arguments
l: longitude, scalar or vector, in degrees.
b: latitude, number of elements as l, in degrees.
Coordinates can be given also as a 2-tuple (l, b).
Output
2-tuple (x, y).
x: x coordinate, same number of elements as l. x is normalized to be in $[-180, 180]$.
y: y coordinate, same number of elements as l. y is normalized to be in $[-90, 90]$.
Example
Get $(x ,y)$ Aitoff coordinates of Sirius, whose Galactic coordinates are $(227.23, -8.890)$.
julia> using AstroLib
julia> x, y = aitoff(227.23, -8.890)
-(-137.92196683723276, -11.772527357473054)
Notes
See AIPS memo No. 46, page 4, for details of the algorithm. This version of aitoff assumes the projection is centered at b=0 degrees.
Code of this function is based on IDL Astronomy User's Library.
Convert Horizon (Alt-Az) coordinates to Hour Angle and Declination.
Explanation
Can deal with the NCP singularity. Intended mainly to be used by program hor2eq.
Arguments
Input coordinates may be either a scalar or an array, of the same dimension.
alt: local apparent altitude, in degrees, scalar or array.
az: the local apparent azimuth, in degrees, scalar or vector, measured east of north!!! If you have measured azimuth west-of-south (like the book Meeus does), convert it to east of north via: az = (az + 180) % 360.
lat: the local geodetic latitude, in degrees, scalar or array.
alt and az can be given as a 2-tuple (alt, az).
Output
2-tuple (ha, dec)
ha: the local apparent hour angle, in degrees. The hour angle is the time that right ascension of 0 hours crosses the local meridian. It is unambiguously defined.
dec: the local apparent declination, in degrees.
The output coordinates are always floating points and have the same type (scalar or array) as the input coordinates.
Example
Arcturus is observed at an apparent altitude of 59d,05m,10s and an azimuth (measured east of north) of 133d,18m,29s while at the latitude of +43.07833 degrees. What are the local hour angle and declination of this object?
julia> using AstroLib
+(-137.92196683723276, -11.772527357473054)
Notes
See AIPS memo No. 46, page 4, for details of the algorithm. This version of aitoff assumes the projection is centered at b=0 degrees.
Code of this function is based on IDL Astronomy User's Library.
Convert Horizon (Alt-Az) coordinates to Hour Angle and Declination.
Explanation
Can deal with the NCP singularity. Intended mainly to be used by program hor2eq.
Arguments
Input coordinates may be either a scalar or an array, of the same dimension.
alt: local apparent altitude, in degrees, scalar or array.
az: the local apparent azimuth, in degrees, scalar or vector, measured east of north!!! If you have measured azimuth west-of-south (like the book Meeus does), convert it to east of north via: az = (az + 180) % 360.
lat: the local geodetic latitude, in degrees, scalar or array.
alt and az can be given as a 2-tuple (alt, az).
Output
2-tuple (ha, dec)
ha: the local apparent hour angle, in degrees. The hour angle is the time that right ascension of 0 hours crosses the local meridian. It is unambiguously defined.
dec: the local apparent declination, in degrees.
The output coordinates are always floating points and have the same type (scalar or array) as the input coordinates.
Example
Arcturus is observed at an apparent altitude of 59d,05m,10s and an azimuth (measured east of north) of 133d,18m,29s while at the latitude of +43.07833 degrees. What are the local hour angle and declination of this object?
julia> using AstroLib
julia> ha, dec = altaz2hadec(ten(59,05,10), ten(133,18,29), 43.07833)
(336.6828582472844, 19.182450965120402)
The widely available XEPHEM code gets:
Hour Angle = 336.683
-Declination = 19.1824
Notes
hadec2altaz converts Hour Angle and Declination to Horizon (Alt-Az) coordinates.
Code of this function is based on IDL Astronomy User's Library.
The 3-vectors outputs dvelh and dvelb are given in a right-handed coordinate system with the +X axis toward the Vernal Equinox, and +Z axis toward the celestial pole.
Code of this function is based on IDL Astronomy User's Library.
The 3-vectors outputs dvelh and dvelb are given in a right-handed coordinate system with the +X axis toward the Vernal Equinox, and +Z axis toward the celestial pole.
Code of this function is based on IDL Astronomy User's Library.
Precess positions from J2000.0 (FK5) to B1950.0 (FK4).
Explanation
Calculates the mean place of a star at B1950.0 on the FK4 system from the mean place at J2000.0 on the FK5 system.
bprecess function has two methods, one for each of the following cases:
the proper motion is known and non-zero
the proper motion is unknown or known to be exactly zero (i.e. extragalactic radio sources). Better precision can be achieved in this case by inputting the epoch of the original observations.
Arguments
The function has 2 methods. The common mandatory arguments are:
ra: input J2000 right ascension, in degrees.
dec: input J2000 declination, in degrees.
The two methods have a different third argument (see "Explanation" section for more details). It can be one of the following:
muradec: 2-element vector containing the proper motion in seconds of arc per tropical century in right ascension and declination.
epoch: scalar giving epoch of original observations.
If none of these two arguments is provided (so bprecess is fed only with right ascension and declination), it is assumed that proper motion is exactly zero and epoch = 2000.
If it is used the method involving muradec argument, the following keywords are available:
parallax (optional numerical keyword): stellar parallax, in seconds of arc.
radvel (optional numerical keyword): radial velocity in km/s.
Right ascension and declination can be passed as the 2-tuple (ra, dec). You can also pass ra, dec, parallax, and radvel as arrays, all of the same length N. In that case, muradec should be a matrix 2×N.
Output
The 2-tuple of right ascension and declination in 1950, in degrees, of input coordinates is returned. If ra and dec (and other possible optional arguments) are arrays, the 2-tuple of arrays (ra1950, dec1950) of the same length as the input coordinates is returned.
Method
The algorithm is taken from the Explanatory Supplement to the Astronomical Almanac 1992, page 186. See also Aoki et al (1983), A&A, 128, 263. URL: http://adsabs.harvard.edu/abs/1983A%26A...128..263A.
Example
The SAO2000 catalogue gives the J2000 position and proper motion for the star HD 119288. Find the B1950 position.
"When transferring individual observations, as opposed to catalog mean place, the safest method is to transform the observations back to the epoch of the observation, on the FK4 system (or in the system that was used to to produce the observed mean place), convert to the FK5 system, and transform to the the epoch and equinox of J2000.0" – from the Explanatory Supplement (1992), p. 180
jprecess performs the precession to J2000 coordinates.
Code of this function is based on IDL Astronomy User's Library.
Deredden a galaxy spectrum using the Calzetti et al. (2000) recipe.
Explanation
Calzetti et al. (2000, ApJ 533, 682; http://adsabs.harvard.edu/abs/2000ApJ...533..682C) developed a recipe for dereddening the spectra of galaxies where massive stars dominate the radiation output, valid between $0.12$ to $2.2$ microns. (calz_unred extrapolates between $0.12$ and $0.0912$ microns.)
Arguments
wave: wavelength (Angstroms)
flux: calibrated flux.
ebv: color excess E(B-V). If a negative ebv is supplied, then fluxes will be reddened rather than deredenned. Note that the supplied color excess should be that derived for the stellar continuum, EBV(stars), which is related to the reddening derived from the gas, EBV(gas), via the Balmer decrement by EBV(stars) = 0.44*EBV(gas).
r_v (optional): ratio of total to selective extinction, default is 4.05. Calzetti et al. (2000) estimate $r_v = 4.05 ± 0.80$ from optical-IR observations of 4 starbursts.
Output
Unreddened flux, same units as flux. Flux values will be left unchanged outside valid domain ($0.0912$ - $2.2$ microns).
Example
Estimate how a flat galaxy spectrum (in wavelength) between $1200 Å$ and $3200 Å$ is altered by a reddening of E(B-V) = 0.1.
"When transferring individual observations, as opposed to catalog mean place, the safest method is to transform the observations back to the epoch of the observation, on the FK4 system (or in the system that was used to to produce the observed mean place), convert to the FK5 system, and transform to the the epoch and equinox of J2000.0" – from the Explanatory Supplement (1992), p. 180
jprecess performs the precession to J2000 coordinates.
Code of this function is based on IDL Astronomy User's Library.
Deredden a galaxy spectrum using the Calzetti et al. (2000) recipe.
Explanation
Calzetti et al. (2000, ApJ 533, 682; http://adsabs.harvard.edu/abs/2000ApJ...533..682C) developed a recipe for dereddening the spectra of galaxies where massive stars dominate the radiation output, valid between $0.12$ to $2.2$ microns. (calz_unred extrapolates between $0.12$ and $0.0912$ microns.)
Arguments
wave: wavelength (Angstroms)
flux: calibrated flux.
ebv: color excess E(B-V). If a negative ebv is supplied, then fluxes will be reddened rather than deredenned. Note that the supplied color excess should be that derived for the stellar continuum, EBV(stars), which is related to the reddening derived from the gas, EBV(gas), via the Balmer decrement by EBV(stars) = 0.44*EBV(gas).
r_v (optional): ratio of total to selective extinction, default is 4.05. Calzetti et al. (2000) estimate $r_v = 4.05 ± 0.80$ from optical-IR observations of 4 starbursts.
Output
Unreddened flux, same units as flux. Flux values will be left unchanged outside valid domain ($0.0912$ - $2.2$ microns).
Example
Estimate how a flat galaxy spectrum (in wavelength) between $1200 Å$ and $3200 Å$ is altered by a reddening of E(B-V) = 0.1.
Code of this function is based on IDL Astronomy User's Library.
The output dra is not multiplied by cos(dec), so that apparentra = ra + d_ra/3600.
These formula are from Meeus, Chapters 23. Accuracy is much better than 1 arcsecond. The maximum deviation due to annual aberration is 20.49'' and occurs when the Earth's velocity is perpendicular to the direction of the star.
Calculate changes in RA and Dec due to nutation of the Earth's rotation
Explanation
Calculates necessary changes to ra and dec due to the nutation of the Earth's rotation axis, as described in Meeus, Chap 23. Uses formulae from Astronomical Almanac, 1984, and does the calculations in equatorial rectangular coordinates to avoid singularities at the celestial poles.
Arguments
jd: julian date, scalar or vector
ra: right ascension in degrees, scalar or vector
dec: declination in degrees, scalar or vector
Output
The 5-tuple (d_ra, d_dec, eps, d_psi, d_eps):
d_ra: correction to right ascension due to nutation, in degrees
d_dec: correction to declination due to nutation, in degrees
eps: the true obliquity of the ecliptic
d_psi: nutation in the longitude of the ecliptic
d_eps: nutation in the obliquity of the ecliptic
Example
Example 23a in Meeus: On 2028 Nov 13.19 TD the mean position of Theta Persei is 2h 46m 11.331s 49d 20' 54.54''. Determine the shift in position due to the Earth's nutation.
julia> using AstroLib
+(30.04404628365077, 6.699400463119431)
d_ra = 30.04404628365103'' (≈ 2.003s)
d_dec = 6.699400463118504''
Notes
Code of this function is based on IDL Astronomy User's Library.
The output dra is not multiplied by cos(dec), so that apparentra = ra + d_ra/3600.
These formula are from Meeus, Chapters 23. Accuracy is much better than 1 arcsecond. The maximum deviation due to annual aberration is 20.49'' and occurs when the Earth's velocity is perpendicular to the direction of the star.
Calculate changes in RA and Dec due to nutation of the Earth's rotation
Explanation
Calculates necessary changes to ra and dec due to the nutation of the Earth's rotation axis, as described in Meeus, Chap 23. Uses formulae from Astronomical Almanac, 1984, and does the calculations in equatorial rectangular coordinates to avoid singularities at the celestial poles.
Arguments
jd: julian date, scalar or vector
ra: right ascension in degrees, scalar or vector
dec: declination in degrees, scalar or vector
Output
The 5-tuple (d_ra, d_dec, eps, d_psi, d_eps):
d_ra: correction to right ascension due to nutation, in degrees
d_dec: correction to declination due to nutation, in degrees
eps: the true obliquity of the ecliptic
d_psi: nutation in the longitude of the ecliptic
d_eps: nutation in the obliquity of the ecliptic
Example
Example 23a in Meeus: On 2028 Nov 13.19 TD the mean position of Theta Persei is 2h 46m 11.331s 49d 20' 54.54''. Determine the shift in position due to the Earth's nutation.
Calculate correction to altitude due to atmospheric refraction.
Explanation
Because the index of refraction of air is not precisely 1.0, the atmosphere bends all incoming light, making a star or other celestial object appear at a slightly different altitude (or elevation) than it really is. It is important to understand the following definitions:
Observed Altitude: The altitude that a star is seen to be, with a telescope. This is where it appears in the sky. This is should be always greater than the apparent altitude.
Apparent Altitude: The altitude that a star would be at, if ~there were no atmosphere~ (sometimes called the "true" altitude). This is usually calculated from an object's celestial coordinates. Apparent altitude should always be smaller than the observed altitude.
Thus, for example, the Sun's apparent altitude when you see it right on the horizon is actually -34 arcminutes.
This program uses a couple of simple formulae to estimate the effect for most optical and radio wavelengths. Typically, you know your observed altitude (from an observation), and want the apparent altitude. To go the other way, this program uses an iterative approach.
Arguments
old_alt: observed altitude in degrees. If to_observe is set to true, this should be apparent altitude
altitude (optional): the height of the observing location, in meters. This is only used to determine an approximate temperature and pressure, if these are not specified separately. Default is 0 i.e. sea level
pressure (optional): the pressure at the observing location, in millibars. Default is NaN
temperature (optional): the temperature at the observing location, in Kelvins. Default is NaN
epsilon (optional): the accuracy to obtain, in arcseconds. If to_observe is true, then it will be calculated. Default is 0.25 arcseconds
to_observe (optional boolean keyword): if set to true, it is assumed that old_alt has apparent altitude as its input and the observed altitude will be found
Output
aout: apparent altitude, in degrees. Observed altitude is returned if to_observe is set to true
Example
The lower limb of the Sun is observed to have altitude of 0d 30'. Calculate the the true (i.e. apparent) altitude of the Sun's lower limb using mean conditions of air pressure and temperature.
julia> using AstroLib
julia> co_refract(0.5)
-0.02584736873098442
Notes
If altitude is set but the temperature or pressure is not, the program will make an intelligent guess for the temperature and pressure.
Wavelength Dependence
This correction is 0 at zenith, about 1 arcminute at 45 degrees, and 34 arcminutes at the horizon for optical wavelengths. The correction is non-negligible at all wavelengths, but is not very easily calculable. These formulae assume a wavelength of 550 nm, and will be accurate to about 4 arcseconds for all visible wavelengths, for elevations of 10 degrees and higher. Amazingly, they are also accurate for radio frequencies less than ~ 100 GHz.
References
Meeus, Astronomical Algorithms, Chapter 15.
Explanatory Supplement to the Astronomical Almanac, 1992.
Methods of Experimental Physics, Vol 12 Part B, Astrophysics, Radio Telescopes, Chapter 2.5, "Refraction Effects in the Neutral Atmosphere", by R.K. Crane.
Code of this function is based on IDL Astronomy User's Library.
If altitude is set but the temperature or pressure is not, the program will make an intelligent guess for the temperature and pressure.
Wavelength Dependence
This correction is 0 at zenith, about 1 arcminute at 45 degrees, and 34 arcminutes at the horizon for optical wavelengths. The correction is non-negligible at all wavelengths, but is not very easily calculable. These formulae assume a wavelength of 550 nm, and will be accurate to about 4 arcseconds for all visible wavelengths, for elevations of 10 degrees and higher. Amazingly, they are also accurate for radio frequencies less than ~ 100 GHz.
References
Meeus, Astronomical Algorithms, Chapter 15.
Explanatory Supplement to the Astronomical Almanac, 1992.
Methods of Experimental Physics, Vol 12 Part B, Astrophysics, Radio Telescopes, Chapter 2.5, "Refraction Effects in the Neutral Atmosphere", by R.K. Crane.
Code of this function is based on IDL Astronomy User's Library.
Convert from Local Civil Time to Local Mean Sidereal Time.
Arguments
The function can be called in two different ways. The only argument common to both methods is longitude:
longitude: the longitude in degrees (east of Greenwich) of the place for which the local sidereal time is desired. The Greenwich mean sidereal time (GMST) can be found by setting longitude = 0.
The civil date to be converted to mean sidereal time can be specified either by providing the Julian days:
jd: this is number of Julian days for the date to be converted.
or the time zone and the date:
tz: the time zone of the site in hours, positive East of the Greenwich meridian (ahead of GMT). Use this parameter to easily account for Daylight Savings time (e.g. -4=EDT, -5 = EST/CDT).
date: this is the local civil time with type DateTime.
Output
The local sidereal time for the date/time specified in hours.
Method
The Julian days of the day and time is question is used to determine the number of days to have passed since 2000-01-01. This is used in conjunction with the GST of that date to extrapolate to the current GST; this is then used to get the LST. See Astronomical Algorithms by Jean Meeus, p. 84 (Eq. 11-4) for the constants used.
Example
Find the Greenwich mean sidereal time (GMST) on 2008-07-30 at 15:53 in Baltimore, Maryland (longitude=-76.72 degrees). The timezone is EDT or tz=-4
julia> using AstroLib, Dates
julia> lst = ct2lst(-76.72, -4, DateTime(2008, 7, 30, 15, 53))
@@ -107,18 +107,18 @@
3-element StaticArraysCore.SVector{3, Float64} with indices SOneTo(3):
17.0
8.0
- 26.466615619137883
Notes
Code of this function is based on IDL Astronomy User's Library.
eci2geo(x, y, z, jd) -> latitude, longitude, altitude
Purpose
Convert Earth-centered inertial coordinates to geographic spherical coordinates.
Explanation
Converts from ECI (Earth-Centered Inertial) (x, y, z) rectangular coordinates to geographic spherical coordinates (latitude, longitude, altitude). Julian day is also needed as input.
ECI coordinates are in km from Earth center at the supplied time (True of Date). Geographic coordinates assume the Earth is a perfect sphere, with radius equal to its equatorial radius.
Arguments
x: ECI x coordinate at jd, in kilometers.
y: ECI y coordinate at jd, in kilometers.
z: ECI z coordinate at jd, in kilometers.
jd: Julian days.
The three coordinates can be passed as a 3-tuple (x, y, z). In addition, x, y, z, and jd can be given as arrays of the same length.
Output
The 3-tuple of geographical coordinate (latitude, longitude, altitude).
latitude: latitude, in degrees.
longitude: longitude, in degrees.
altitude: altitude, in kilometers.
If ECI coordinates are given as arrays, a 3-tuple of arrays of the same length is returned.
Example
Obtain the geographic direction of the vernal point on 2015-06-30T14:03:12.857, in geographic coordinates, at altitude 600 km. Note: equatorial radii of Solar System planets in meters are stored into AstroLib.planets dictionary.
julia> using AstroLib
+(-0.3, 1.165, 0.905, -0.665)
Notes
Code of this function is based on IDL Astronomy User's Library.
eci2geo(x, y, z, jd) -> latitude, longitude, altitude
Purpose
Convert Earth-centered inertial coordinates to geographic spherical coordinates.
Explanation
Converts from ECI (Earth-Centered Inertial) (x, y, z) rectangular coordinates to geographic spherical coordinates (latitude, longitude, altitude). Julian day is also needed as input.
ECI coordinates are in km from Earth center at the supplied time (True of Date). Geographic coordinates assume the Earth is a perfect sphere, with radius equal to its equatorial radius.
Arguments
x: ECI x coordinate at jd, in kilometers.
y: ECI y coordinate at jd, in kilometers.
z: ECI z coordinate at jd, in kilometers.
jd: Julian days.
The three coordinates can be passed as a 3-tuple (x, y, z). In addition, x, y, z, and jd can be given as arrays of the same length.
Output
The 3-tuple of geographical coordinate (latitude, longitude, altitude).
latitude: latitude, in degrees.
longitude: longitude, in degrees.
altitude: altitude, in kilometers.
If ECI coordinates are given as arrays, a 3-tuple of arrays of the same length is returned.
Example
Obtain the geographic direction of the vernal point on 2015-06-30T14:03:12.857, in geographic coordinates, at altitude 600 km. Note: equatorial radii of Solar System planets in meters are stored into AstroLib.planets dictionary.
julia> using AstroLib
julia> x = AstroLib.planets["earth"].eqradius*1e-3 + 600;
julia> lat, long, alt = eci2geo(x, 0, 0, jdcnv("2015-06-30T14:03:12.857"))
-(0.0, 230.87301833205856, 600.0)
These coordinates can be further transformed into geodetic coordinates using geo2geodetic or into geomagnetic coordinates using geo2mag.
Notes
geo2eci converts geographic spherical coordinates to Earth-centered inertial coordinates.
Code of this function is based on IDL Astronomy User's Library.
Convert right ascension $l$ and declination $b$ to coordinate $(x, y)$ using an equal-area polar projection.
Explanation
The output $x$ and $y$ coordinates are scaled to be in the range $[-90, 90]$ and to go from equator to pole to equator. Output map points can be centered on the north pole or south pole.
Arguments
l: longitude, scalar or vector, in degrees
b: latitude, same number of elements as right ascension, in degrees
southpole (optional boolean keyword): keyword to indicate that the plot is to be centered on the south pole instead of the north pole. Default is false.
Output
The 2-tuple $(x, y)$:
$x$ coordinate, same number of elements as right ascension, normalized to be in the range $[-90, 90]$.
$y$ coordinate, same number of elements as declination, normalized to be in the range $[-90, 90]$.
Example
julia> using AstroLib
+" 17 42 25.6 +16 25 26"
Notes
Code of this function is based on IDL Astronomy User's Library.
Convert right ascension $l$ and declination $b$ to coordinate $(x, y)$ using an equal-area polar projection.
Explanation
The output $x$ and $y$ coordinates are scaled to be in the range $[-90, 90]$ and to go from equator to pole to equator. Output map points can be centered on the north pole or south pole.
Arguments
l: longitude, scalar or vector, in degrees
b: latitude, same number of elements as right ascension, in degrees
southpole (optional boolean keyword): keyword to indicate that the plot is to be centered on the south pole instead of the north pole. Default is false.
Output
The 2-tuple $(x, y)$:
$x$ coordinate, same number of elements as right ascension, normalized to be in the range $[-90, 90]$.
$y$ coordinate, same number of elements as declination, normalized to be in the range $[-90, 90]$.
Transform between Galactic, celestial, and ecliptic coordinates.
Explanation
The function is used by the astro procedure.
Arguments
ai: input longitude, scalar or vector.
bi: input latitude, scalar or vector.
select : integer input specifying type of coordinate transformation. SELECT From To | SELECT From To 1 RA-Dec (2000) Galactic | 4 Ecliptic RA-Dec 2 Galactic RA-DEC | 5 Ecliptic Galactic 3 RA-Dec Ecliptic | 6 Galactic Ecliptic
FK4 (optional boolean keyword) : if this keyword is set to true, then input and output celestial and ecliptic coordinates should be given in equinox B1950. When false, by default, they should be given in equinox J2000.
radians (optional boolean keyword) : if this keyword is set to true, all input and output angles are in radians rather than degrees.
Output
a 2-tuple (ao, bo):
ao: output longitude in degrees.
bo: output latitude in degrees.
Example
Find the Galactic coordinates of Cyg X-1 (ra=299.590315, dec=35.201604)
julia> using AstroLib
+(72.78853915267848, 12.83458333897169)
Notes
Code of this function is based on IDL Astronomy User's Library.
Transform between Galactic, celestial, and ecliptic coordinates.
Explanation
The function is used by the astro procedure.
Arguments
ai: input longitude, scalar or vector.
bi: input latitude, scalar or vector.
select : integer input specifying type of coordinate transformation. SELECT From To | SELECT From To 1 RA-Dec (2000) Galactic | 4 Ecliptic RA-Dec 2 Galactic RA-DEC | 5 Ecliptic Galactic 3 RA-Dec Ecliptic | 6 Galactic Ecliptic
FK4 (optional boolean keyword) : if this keyword is set to true, then input and output celestial and ecliptic coordinates should be given in equinox B1950. When false, by default, they should be given in equinox J2000.
radians (optional boolean keyword) : if this keyword is set to true, all input and output angles are in radians rather than degrees.
Output
a 2-tuple (ao, bo):
ao: output longitude in degrees.
bo: output latitude in degrees.
Example
Find the Galactic coordinates of Cyg X-1 (ra=299.590315, dec=35.201604)
julia> using AstroLib
julia> euler(299.590315, 35.201604, 1)
-(71.33498957116959, 3.0668335310640984)
Notes
Code of this function is based on IDL Astronomy User's Library.
gal_uvw(ra, dec, pmra, pmdec, vrad, plx[, lsr=true]) -> u, v, w
Purpose
Calculate the Galactic space velocity $(u, v, w)$ of a star.
Explanation
Calculates the Galactic space velocity $(u, v, w)$ of a star given its (1) coordinates, (2) proper motion, (3) parallax, and (4) radial velocity.
Arguments
User must supply a position, proper motion, radial velocity and parallax. Either scalars or arrays all of the same length can be supplied.
Position:
ra: right ascension, in degrees
dec: declination, in degrees
Proper Motion
pmra: proper motion in right ascension in arc units (typically milli-arcseconds/yr). If given $μ_α$ – proper motion in seconds of time/year – then this is equal to $15 μ_α \cos(\text{dec})$.
pmdec: proper motion in declination (typically mas/yr).
Radial Velocity
vrad: velocity in km/s
Parallax
plx: parallax with same distance units as proper motion measurements typically milliarcseconds (mas)
If you know the distance in parsecs, then set plx to $1000/\text{distance}$, if proper motion measurements are given in milli-arcseconds/yr.
There is an additional optional keyword:
lsr (optional boolean keyword): if this keyword is set to true, then the output velocities will be corrected for the solar motion $(u, v, w)_⊙ = (-8.5, 13.38, 6.49)$ (Coşkunoǧlu et al. 2011 MNRAS, 412, 1237; DOI:10.1111/j.1365-2966.2010.17983.x) to the local standard of rest (LSR). Note that the value of the solar motion through the LSR remains poorly determined.
Output
The 3-tuple $(u, v, w)$
$u$: velocity (km/s) positive toward the Galactic anticenter
$v$: velocity (km/s) positive in the direction of Galactic rotation
$w$: velocity (km/s) positive toward the North Galactic Pole
Method
Follows the general outline of Johnson & Soderblom (1987, AJ, 93, 864; DOI:10.1086/114370) except that $u$ is positive outward toward the Galactic anticenter, and the J2000 transformation matrix to Galactic coordinates is taken from the introduction to the Hipparcos catalog.
Example
Compute the U,V,W coordinates for the halo star HD 6755. Use values from Hipparcos catalog, and correct to the LSR.
julia> using AstroLib
+27.423535345634598
Notes
Code of this function is based on IDL Astronomy User's Library.
gal_uvw(ra, dec, pmra, pmdec, vrad, plx[, lsr=true]) -> u, v, w
Purpose
Calculate the Galactic space velocity $(u, v, w)$ of a star.
Explanation
Calculates the Galactic space velocity $(u, v, w)$ of a star given its (1) coordinates, (2) proper motion, (3) parallax, and (4) radial velocity.
Arguments
User must supply a position, proper motion, radial velocity and parallax. Either scalars or arrays all of the same length can be supplied.
Position:
ra: right ascension, in degrees
dec: declination, in degrees
Proper Motion
pmra: proper motion in right ascension in arc units (typically milli-arcseconds/yr). If given $μ_α$ – proper motion in seconds of time/year – then this is equal to $15 μ_α \cos(\text{dec})$.
pmdec: proper motion in declination (typically mas/yr).
Radial Velocity
vrad: velocity in km/s
Parallax
plx: parallax with same distance units as proper motion measurements typically milliarcseconds (mas)
If you know the distance in parsecs, then set plx to $1000/\text{distance}$, if proper motion measurements are given in milli-arcseconds/yr.
There is an additional optional keyword:
lsr (optional boolean keyword): if this keyword is set to true, then the output velocities will be corrected for the solar motion $(u, v, w)_⊙ = (-8.5, 13.38, 6.49)$ (Coşkunoǧlu et al. 2011 MNRAS, 412, 1237; DOI:10.1111/j.1365-2966.2010.17983.x) to the local standard of rest (LSR). Note that the value of the solar motion through the LSR remains poorly determined.
Output
The 3-tuple $(u, v, w)$
$u$: velocity (km/s) positive toward the Galactic anticenter
$v$: velocity (km/s) positive in the direction of Galactic rotation
$w$: velocity (km/s) positive toward the North Galactic Pole
Method
Follows the general outline of Johnson & Soderblom (1987, AJ, 93, 864; DOI:10.1086/114370) except that $u$ is positive outward toward the Galactic anticenter, and the J2000 transformation matrix to Galactic coordinates is taken from the introduction to the Hipparcos catalog.
Example
Compute the U,V,W coordinates for the halo star HD 6755. Use values from Hipparcos catalog, and correct to the LSR.
julia> using AstroLib
julia> ra = ten(1,9,42.3)*15.; dec = ten(61,32,49.5);
@@ -157,13 +157,13 @@
julia> vrad = -321.4; dis = 129; # distance in parsecs
julia> u, v, w = gal_uvw(ra, dec, pmra, pmdec, vrad, 1e3/dis, lsr=true)
-(118.2110474553902, -466.4828898385057, 88.16573278565097)
Notes
This function does not take distance as input. See "Arguments" section above for how to provide it using parallax argument.
Code of this function is based on IDL Astronomy User's Library.
geo2eci(latitude, longitude, altitude, jd) -> x, y, z
Purpose
Convert geographic spherical coordinates to Earth-centered inertial coordinates.
Explanation
Converts from geographic spherical coordinates (latitude, longitude, altitude) to ECI (Earth-Centered Inertial) (x, y, z) rectangular coordinates. Julian days is also needed.
Geographic coordinates assume the Earth is a perfect sphere, with radius equal to its equatorial radius. ECI coordinates are in km from Earth center at epoch TOD (True of Date).
Arguments
latitude: geographic latitude, in degrees.
longitude: geographic longitude, in degrees.
altitude: geographic altitude, in kilometers.
jd: Julian days.
The three coordinates can be passed as a 3-tuple (latitude, longitude, altitude). In addition, latitude, longitude, altitude, and jd can be given as arrays of the same length.
Output
The 3-tuple of ECI (x, y, z) coordinates, in kilometers. The TOD epoch is the supplied jd time.
If geographical coordinates are given as arrays, a 3-tuple of arrays of the same length is returned.
Example
Obtain the ECI coordinates of the intersection of the equator and Greenwich's meridian on 2015-06-30T14:03:12.857
julia> using AstroLib
+1.590442261600714
Notes
The function sphdist provides an alternate method of computing a spherical distance.
The Haversine formula can give rounding errors for antipodal points.
Code of this function is based on IDL Astronomy User's Library.
geo2eci(latitude, longitude, altitude, jd) -> x, y, z
Purpose
Convert geographic spherical coordinates to Earth-centered inertial coordinates.
Explanation
Converts from geographic spherical coordinates (latitude, longitude, altitude) to ECI (Earth-Centered Inertial) (x, y, z) rectangular coordinates. Julian days is also needed.
Geographic coordinates assume the Earth is a perfect sphere, with radius equal to its equatorial radius. ECI coordinates are in km from Earth center at epoch TOD (True of Date).
Arguments
latitude: geographic latitude, in degrees.
longitude: geographic longitude, in degrees.
altitude: geographic altitude, in kilometers.
jd: Julian days.
The three coordinates can be passed as a 3-tuple (latitude, longitude, altitude). In addition, latitude, longitude, altitude, and jd can be given as arrays of the same length.
Output
The 3-tuple of ECI (x, y, z) coordinates, in kilometers. The TOD epoch is the supplied jd time.
If geographical coordinates are given as arrays, a 3-tuple of arrays of the same length is returned.
Example
Obtain the ECI coordinates of the intersection of the equator and Greenwich's meridian on 2015-06-30T14:03:12.857
Convert from geographic (or planetographic) to geodetic coordinates.
Explanation
Converts from geographic (latitude, longitude, altitude) to geodetic (latitude, longitude, altitude). In geographic coordinates, the Earth is assumed a perfect sphere with a radius equal to its equatorial radius. The geodetic (or ellipsoidal) coordinate system takes into account the Earth's oblateness.
Geographic and geodetic longitudes are identical. Geodetic latitude is the angle between local zenith and the equatorial plane. Geographic and geodetic altitudes are both the closest distance between the satellite and the ground.
Arguments
The function has two base methods. The arguments common to all methods and always mandatory are latitude, longitude, and altitude:
latitude: geographic latitude, in degrees.
longitude: geographic longitude, in degrees.
altitude: geographic altitude, in kilometers.
In order to convert to geodetic coordinates, you can either provide custom equatorial and polar radii of the planet or use the values of one of the planets of Solar System (Pluto included).
If you want to use the method with explicit equatorial and polar radii the additional mandatory arguments are:
equatorial_radius: value of the equatorial radius of the body, in kilometers.
polar_radius: value of the polar radius of the body, in kilometers.
Instead, if you want to use the method with the selection of a planet, the only additional argument is the planet name:
planet (optional string argument): string with the name of the Solar System planet, from "Mercury" to "Pluto". If omitted (so, when only latitude, longitude, and altitude are provided), the default is "Earth".
In all cases, the three coordinates can be passed as a 3-tuple (latitude, longitude, altitude). In addition, geographical latitude, longitude, and altitude can be given as arrays of the same length.
Output
The 3-tuple (latitude, longitude, altitude) in geodetic coordinates, for the body with specified equatorial and polar radii (Earth by default).
If geographical coordinates are given as arrays, a 3-tuple of arrays of the same length is returned.
Method
Stephen P. Keeler and Yves Nievergelt, "Computing geodetic coordinates", SIAM Rev. Vol. 40, No. 2, pp. 300-309, June 1998 (DOI:10.1137/S0036144597323921).
Planetary constants are from Planetary Fact Sheet (http://nssdc.gsfc.nasa.gov/planetary/factsheet/index.html).
Example
Locate the Earth geographic North pole (latitude: 90°, longitude: 0°, altitude 0 km), in geodetic coordinates:
Convert from geodetic (or planetodetic) to geographic coordinates.
Explanation
Converts from geodetic (latitude, longitude, altitude) to geographic (latitude, longitude, altitude). In geographic coordinates, the Earth is assumed a perfect sphere with a radius equal to its equatorial radius. The geodetic (or ellipsoidal) coordinate system takes into account the Earth's oblateness.
Geographic and geodetic longitudes are identical. Geodetic latitude is the angle between local zenith and the equatorial plane. Geographic and geodetic altitudes are both the closest distance between the satellite and the ground.
Arguments
The function has two base methods. The arguments common to all methods and always mandatory are latitude, longitude, and altitude:
latitude: geodetic latitude, in degrees.
longitude: geodetic longitude, in degrees.
altitude: geodetic altitude, in kilometers.
In order to convert to geographic coordinates, you can either provide custom equatorial and polar radii of the planet or use the values of one of the planets of Solar System (Pluto included).
If you want to use the method with explicit equatorial and polar radii the additional mandatory arguments are:
equatorial_radius: value of the equatorial radius of the body, in kilometers.
polar_radius: value of the polar radius of the body, in kilometers.
Instead, if you want to use the method with the selection of a planet, the only additional argument is the planet name:
planet (optional string argument): string with the name of the Solar System planet, from "Mercury" to "Pluto". If omitted (so, when only latitude, longitude, and altitude are provided), the default is "Earth".
In all cases, the three coordinates can be passed as a 3-tuple (latitude, longitude, altitude). In addition, geodetic latitude, longitude, and altitude can be given as arrays of the same length.
Output
The 3-tuple (latitude, longitude, altitude) in geographic coordinates, for the body with specified equatorial and polar radii (Earth by default).
If geodetic coordinates are given as arrays, a 3-tuple of arrays of the same length is returned.
Method
Stephen P. Keeler and Yves Nievergelt, "Computing geodetic coordinates", SIAM Rev. Vol. 40, No. 2, pp. 300-309, June 1998 (DOI:10.1137/S0036144597323921).
Planetary constants from "Allen's Astrophysical Quantities", Fourth Ed., (2000).
Example
Find geographic coordinates of geodetic North pole (latitude: 90°, longitude: 0°, altitude 0 km) of the Earth:
julia> using AstroLib
@@ -194,14 +194,14 @@
(90.0, 0.0, -4638.0)
Find geographic coordinates for point of geodetic coordinates (latitude, longitude, altitude) = (43.16°, -24.32°, 3.87 km) on a planet with equatorial radius 8724.32 km and polar radius 8619.19 km:
Returns the UTC date in "CCYY-MM-DD" format for FITS headers.
Explanation
This is the format required by the DATE and DATE-OBS keywords in a FITS header.
Argument
date (optional): the date in UTC standard. If omitted, defaults to the current UTC time. Each element can be either a DateTime type or anything that can be converted to that type. In the case of vectorial input, each element is considered as a date, so you cannot provide a date by parts.
old (optional boolean keyword): see below.
timetag (optional boolean keyword): see below.
Output
A string with the date formatted according to the given optional keywords.
When no optional keywords (timetag and old) are supplied, the format of the output string is "CCYY-MM-DD" (year-month-day part of the date), where CCYY represents a 4-digit calendar year, MM the 2-digit ordinal number of a calendar month within the calendar year, and DD the 2-digit ordinal number of a day within the calendar month.
If the boolean keyword old is true (default: false), the year-month-day part of date has "DD/MM/YY" format. This is the formerly (pre-1997) recommended for FITS. Note that this format is now deprecated because it uses only a 2-digit representation of the year.
If the boolean keyword timetag is true (default: false), "Thh:mm:ss" is appended to the year-month-day part of the date, where <hh> represents the hour in the day, <mm> the minutes, <ss> the seconds, and the literal 'T' the ISO 8601 time designator.
Note that old and timetag keywords can be used together, so that the output string will have "DD/MM/YYThh:mm:ss" format.
Example
julia> using AstroLib, Dates
+(42.46772711708433, -24.32, -44.52902080669082)
Notes
geo2geodetic converts from geographic (or planetographic) to geodetic coordinates.
Code of this function is based on IDL Astronomy User's Library.
Returns the UTC date in "CCYY-MM-DD" format for FITS headers.
Explanation
This is the format required by the DATE and DATE-OBS keywords in a FITS header.
Argument
date (optional): the date in UTC standard. If omitted, defaults to the current UTC time. Each element can be either a DateTime type or anything that can be converted to that type. In the case of vectorial input, each element is considered as a date, so you cannot provide a date by parts.
old (optional boolean keyword): see below.
timetag (optional boolean keyword): see below.
Output
A string with the date formatted according to the given optional keywords.
When no optional keywords (timetag and old) are supplied, the format of the output string is "CCYY-MM-DD" (year-month-day part of the date), where CCYY represents a 4-digit calendar year, MM the 2-digit ordinal number of a calendar month within the calendar year, and DD the 2-digit ordinal number of a day within the calendar month.
If the boolean keyword old is true (default: false), the year-month-day part of date has "DD/MM/YY" format. This is the formerly (pre-1997) recommended for FITS. Note that this format is now deprecated because it uses only a 2-digit representation of the year.
If the boolean keyword timetag is true (default: false), "Thh:mm:ss" is appended to the year-month-day part of the date, where <hh> represents the hour in the day, <mm> the minutes, <ss> the seconds, and the literal 'T' the ISO 8601 time designator.
Note that old and timetag keywords can be used together, so that the output string will have "DD/MM/YYThh:mm:ss" format.
A discussion of the DATExxx syntax in FITS headers can be found in http://www.cv.nrao.edu/fits/documents/standards/year2000.txt
Those who wish to use need further flexibility in their date formats (e.g. to use TAI time) should look at Bill Thompson's time routines in http://sohowww.nascom.nasa.gov/solarsoft/gen/idl/time
Convert Hour Angle and Declination to Horizon (Alt-Az) coordinates.
Explanation
Can deal with the NCP singularity. Intended mainly to be used by program eq2hor.
Arguments
Input coordinates may be either a scalar or an array, of the same dimension.
ha: the local apparent hour angle, in degrees. The hour angle is the time that right ascension of 0 hours crosses the local meridian. It is unambiguously defined.
dec: the local apparent declination, in degrees.
lat: the local geodetic latitude, in degrees, scalar or array.
ws (optional boolean keyword): if true, the output azimuth is measured West from South. The default is to measure azimuth East from North.
ha and dec can be given as a 2-tuple (ha, dec).
Output
2-tuple (alt, az)
alt: local apparent altitude, in degrees.
az: the local apparent azimuth, in degrees.
The output coordinates are always floating points and have the same type (scalar or array) as the input coordinates.
Example
Arcturus is observed at an apparent hour angle of 336.6829 and a declination of 19.1825 while at the latitude of +43° 4' 42''. What are the local altitude and azimuth of this object?
julia> using AstroLib
+"21937-05-30T09:59:00"
Notes
A discussion of the DATExxx syntax in FITS headers can be found in http://www.cv.nrao.edu/fits/documents/standards/year2000.txt
Those who wish to use need further flexibility in their date formats (e.g. to use TAI time) should look at Bill Thompson's time routines in http://sohowww.nascom.nasa.gov/solarsoft/gen/idl/time
Convert Hour Angle and Declination to Horizon (Alt-Az) coordinates.
Explanation
Can deal with the NCP singularity. Intended mainly to be used by program eq2hor.
Arguments
Input coordinates may be either a scalar or an array, of the same dimension.
ha: the local apparent hour angle, in degrees. The hour angle is the time that right ascension of 0 hours crosses the local meridian. It is unambiguously defined.
dec: the local apparent declination, in degrees.
lat: the local geodetic latitude, in degrees, scalar or array.
ws (optional boolean keyword): if true, the output azimuth is measured West from South. The default is to measure azimuth East from North.
ha and dec can be given as a 2-tuple (ha, dec).
Output
2-tuple (alt, az)
alt: local apparent altitude, in degrees.
az: the local apparent azimuth, in degrees.
The output coordinates are always floating points and have the same type (scalar or array) as the input coordinates.
Example
Arcturus is observed at an apparent hour angle of 336.6829 and a declination of 19.1825 while at the latitude of +43° 4' 42''. What are the local altitude and azimuth of this object?
julia> using AstroLib
julia> alt, az = hadec2altaz(336.6829, 19.1825, ten(43, 4, 42))
-(59.08617155005685, 133.3080693440254)
Notes
altaz2hadec converts Horizon (Alt-Az) coordinates to Hour Angle and Declination.
Code of this function is based on IDL Astronomy User's Library.
The mean orbital elements for epoch J2000 are used. These are derived from a 250 yr least squares fit of the DE 200 planetary ephemeris to a Keplerian orbit where each element is allowed to vary linearly with time. Useful mainly for dates between 1800 and 2050, this solution fits the terrestrial planet orbits to ~25'' or better, but achieves only ~600'' for Saturn.
The mean orbital elements for epoch J2000 are used. These are derived from a 250 yr least squares fit of the DE 200 planetary ephemeris to a Keplerian orbit where each element is allowed to vary linearly with time. Useful mainly for dates between 1800 and 2050, this solution fits the terrestrial planet orbits to ~25'' or better, but achieves only ~600'' for Saturn.
Convert geocentric (reduced) Julian date to heliocentric Julian date.
Explanation
This procedure corrects for the extra light travel time between the Earth and the Sun.
An online calculator for this quantity is available at http://www.physics.sfasu.edu/astro/javascript/hjd.html
Users requiring more precise calculations and documentation should look at the IDL code available at http://astroutils.astronomy.ohio-state.edu/time/
Arguments
date: reduced Julian date (= JD - 2400000). You can use juldate() to calculate the reduced Julian date.
ra and dec: right ascension and declination in degrees. Default equinox is J2000.
B1950 (optional boolean keyword): if set to true, then input coordinates are assumed to be in equinox B1950 coordinates. Default is false.
diff (optional boolean keyword): if set to true, the function returns the time difference (heliocentric JD - geocentric JD) in seconds. Default is false.
Output
The return value depends on the value of diff optional keywords:
if diff is false (default), then the heliocentric reduced Julian date is returned.
if diff is true, then the time difference in seconds between the geocentric and heliocentric Julian date is returned.
Example
What is the heliocentric Julian date of an observation of V402 Cygni (J2000: RA = 20 9 7.8, Dec = 37 09 07) taken on June 15, 2016 at 11:40 UT?
The user should ensure consistency with all time systems being used (i.e. jd and t should be in the same units and time system). Generally, users should reduce large time values by subtracting a large constant offset, which may improve numerical accuracy.
If using the the function juldate and helio_jd, the reduced HJD time system must be used throughtout.
Code of this function is based on IDL Astronomy User's Library.
The user should ensure consistency with all time systems being used (i.e. jd and t should be in the same units and time system). Generally, users should reduce large time values by subtracting a large constant offset, which may improve numerical accuracy.
If using the the function juldate and helio_jd, the reduced HJD time system must be used throughtout.
Code of this function is based on IDL Astronomy User's Library.
Compute an N-component power-law logarithmic initial mass function (IMF).
Explanation
The function is normalized so that the total mass distribution equals one solar mass.
Arguments
mass: mass in units of solar mass, vector.
expon: power law exponent, vector. The number of values in expon equals the number of different power-law components in the IMF.
mass_range: vector containing the mass upper and lower limits of the IMF and masses where the IMF exponent changes. The number of values in massrange should be one more than in expon. The values in massrange should be monotonically increasing and positive.
Output
psi: mass function, number of stars per unit logarithmic mass interval evaluated for supplied masses.
Example
Show the number of stars per unit mass interval at 3 Msun for a Salpeter (expon = -1.35) IMF, with a mass range from 0.1 MSun to 110 Msun.
julia> using AstroLib
+" 00 13 17.4 +15 11 26"
Notes
Code of this function is based on IDL Astronomy User's Library.
Compute an N-component power-law logarithmic initial mass function (IMF).
Explanation
The function is normalized so that the total mass distribution equals one solar mass.
Arguments
mass: mass in units of solar mass, vector.
expon: power law exponent, vector. The number of values in expon equals the number of different power-law components in the IMF.
mass_range: vector containing the mass upper and lower limits of the IMF and masses where the IMF exponent changes. The number of values in massrange should be one more than in expon. The values in massrange should be monotonically increasing and positive.
Output
psi: mass function, number of stars per unit logarithmic mass interval evaluated for supplied masses.
Example
Show the number of stars per unit mass interval at 3 Msun for a Salpeter (expon = -1.35) IMF, with a mass range from 0.1 MSun to 110 Msun.
ismeuv(wave, hcol[, he1col=hcol*0.1, he2col=0, fano=false]) -> tau
Purpose
Compute the continuum interstellar EUV optical depth
Explanation
The EUV optical depth is computed from the photoionization of hydrogen and helium.
Arguments
wave: wavelength value (in Angstroms). Useful range is 40 - 912 A; at shorter wavelength metal opacity should be considered, at longer wavelengths there is no photoionization.
hcol: interstellar hydrogen column density in cm-2.
he1col (optional): neutral helium column density in cm-2. Default is 0.1*hcol (10% of hydrogen column)
he2col (optional): ionized helium column density in cm-2 Default is 0.
fano (optional boolean keyword): If this keyword is true, then the 4 strongest auto-ionizing resonances of He I are included. The shape of these resonances is given by a Fano profile - see Rumph, Bowyer, & Vennes 1994, AJ, 107, 2108. If these resonances are included then the input wavelength vector should have a fine (>~0.01 A) grid between 190 A and 210 A, since the resonances are very narrow.
One has a model EUV spectrum with wavelength, w (in Angstroms). Find the EUV optical depth by 1e18 cm-2 of HI, with N(HeI)/N(HI) = N(HeII)/N(HI) = 0.05.
julia> using AstroLib
+ 0.01294143518151214
Notes
Code of this function is based on IDL Astronomy User's Library.
ismeuv(wave, hcol[, he1col=hcol*0.1, he2col=0, fano=false]) -> tau
Purpose
Compute the continuum interstellar EUV optical depth
Explanation
The EUV optical depth is computed from the photoionization of hydrogen and helium.
Arguments
wave: wavelength value (in Angstroms). Useful range is 40 - 912 A; at shorter wavelength metal opacity should be considered, at longer wavelengths there is no photoionization.
hcol: interstellar hydrogen column density in cm-2.
he1col (optional): neutral helium column density in cm-2. Default is 0.1*hcol (10% of hydrogen column)
he2col (optional): ionized helium column density in cm-2 Default is 0.
fano (optional boolean keyword): If this keyword is true, then the 4 strongest auto-ionizing resonances of He I are included. The shape of these resonances is given by a Fano profile - see Rumph, Bowyer, & Vennes 1994, AJ, 107, 2108. If these resonances are included then the input wavelength vector should have a fine (>~0.01 A) grid between 190 A and 210 A, since the resonances are very narrow.
One has a model EUV spectrum with wavelength, w (in Angstroms). Find the EUV optical depth by 1e18 cm-2 of HI, with N(HeI)/N(HI) = N(HeII)/N(HI) = 0.05.
Precess positions from B1950.0 (FK4) to J2000.0 (FK5).
Explanation
Calculate the mean place of a star at J2000.0 on the FK5 system from the mean place at B1950.0 on the FK4 system.
jprecess function has two methods, one for each of the following cases:
the proper motion is known and non-zero
the proper motion is unknown or known to be exactly zero (i.e. extragalactic radio sources). Better precision can be achieved in this case by inputting the epoch of the original observations.
Arguments
The function has 2 methods. The common mandatory arguments are:
ra: input B1950 right ascension, in degrees.
dec: input B1950 declination, in degrees.
The two methods have a different third argument (see "Explanation" section for more details). It can be one of the following:
muradec: 2-element vector containing the proper motion in seconds of arc per tropical century in right ascension and declination.
epoch: scalar giving epoch of original observations.
If none of these two arguments is provided (so jprecess is fed only with right ascension and declination), it is assumed that proper motion is exactly zero and epoch = 1950.
If it is used the method involving muradec argument, the following keywords are available:
parallax (optional numerical keyword): stellar parallax, in seconds of arc.
radvel (optional numerical keyword): radial velocity in km/s.
Right ascension and declination can be passed as the 2-tuple (ra, dec). You can also pass ra, dec, parallax, and radvel as arrays, all of the same length N. In that case, muradec should be a matrix 2×N.
Output
The 2-tuple of right ascension and declination in 2000, in degrees, of input coordinates is returned. If ra and dec (and other possible optional arguments) are arrays, the 2-tuple of arrays (ra2000, dec2000) of the same length as the input coordinates is returned.
Method
The algorithm is taken from the Explanatory Supplement to the Astronomical Almanac 1992, page 184. See also Aoki et al (1983), A&A, 128, 263. URL: http://adsabs.harvard.edu/abs/1983A%26A...128..263A.
Example
The SAO catalogue gives the B1950 position and proper motion for the star HD 119288. Find the J2000 position.
"When transferring individual observations, as opposed to catalog mean place, the safest method is to tranform the observations back to the epoch of the observation, on the FK4 system (or in the system that was used to to produce the observed mean place), convert to the FK5 system, and transform to the the epoch and equinox of J2000.0" – from the Explanatory Supplement (1992), p. 180
bprecess performs the precession to B1950 coordinates.
Code of this function is based on IDL Astronomy User's Library.
Julian Day Number is a count of days elapsed since Greenwich mean noon on 1 January 4713 B.C. Julian Days are the number of Julian days followed by the fraction of the day elapsed since the preceding noon.
This function takes the given date and returns the number of Julian calendar days since epoch 1858-11-16T12:00:00 (Reduced Julian Days = Julian Days - 2400000).
Argument
date: date in Julian Calendar, UTC standard. Each element can be given in DateTime type or anything that can be converted to that type.
Output
The number of Reduced Julian Days is returned.
Example
Get number of Reduced Julian Days at 2016-03-20T15:24:00.
"When transferring individual observations, as opposed to catalog mean place, the safest method is to tranform the observations back to the epoch of the observation, on the FK4 system (or in the system that was used to to produce the observed mean place), convert to the FK5 system, and transform to the the epoch and equinox of J2000.0" – from the Explanatory Supplement (1992), p. 180
bprecess performs the precession to B1950 coordinates.
Code of this function is based on IDL Astronomy User's Library.
Julian Day Number is a count of days elapsed since Greenwich mean noon on 1 January 4713 B.C. Julian Days are the number of Julian days followed by the fraction of the day elapsed since the preceding noon.
This function takes the given date and returns the number of Julian calendar days since epoch 1858-11-16T12:00:00 (Reduced Julian Days = Julian Days - 2400000).
Argument
date: date in Julian Calendar, UTC standard. Each element can be given in DateTime type or anything that can be converted to that type.
Output
The number of Reduced Julian Days is returned.
Example
Get number of Reduced Julian Days at 2016-03-20T15:24:00.
Julian Calendar is assumed, thus before 1582-10-15T00:00:00 this function is not the inverse of daycnv. For the conversion proleptic Gregorian date to number of Julian days, use jdcnv, which is the inverse of daycnv.
Code of this function is based on IDL Astronomy User's Library.
Solve Kepler's equation in the elliptic motion regime ($0 ≤ e ≤ 1$) and return eccentric anomaly $E$.
Explanation
In order to find the position of a body in elliptic motion (e.g., in the two-body problem) at a given time $t$, one has to solve the Kepler's equation
\[M(t) = E(t) - e \sin E(t)\]
where $M(t) = (t - t_0)/P$ is the mean anomaly, $E(t)$ the eccentric anomaly, $e$ the eccentricity of the orbit, $t_0$ is the time of periapsis passage, and $P$ is the period of the orbit. Usually the eccentricity is given and one wants to find the eccentric anomaly $E(t)$ at a specific time $t$, so that also the mean anomaly $M(t)$ is known.
Arguments
M: mean anomaly.
e: eccentricity, in the elliptic motion regime ($0 ≤ e ≤ 1$)
Output
The eccentric anomaly $E$, restricted to the range $[-π, π]$.
Method
Many different numerical methods exist to solve Kepler's equation. This function implements the algorithm proposed in Markley (1995) Celestial Mechanics and Dynamical Astronomy, 63, 101 (DOI:10.1007/BF00691917). This method is not iterative, requires only four transcendental function evaluations, and has been proved to be fast and efficient over the entire range of elliptic motion $0 ≤ e ≤ 1$.
Example
Find the eccentric anomaly for an orbit with eccentricity $e = 0.7$ and for $M(t) = 8π/3$.
julia> using AstroLib
+57468.14166666667
Notes
Julian Calendar is assumed, thus before 1582-10-15T00:00:00 this function is not the inverse of daycnv. For the conversion proleptic Gregorian date to number of Julian days, use jdcnv, which is the inverse of daycnv.
Code of this function is based on IDL Astronomy User's Library.
Solve Kepler's equation in the elliptic motion regime ($0 ≤ e ≤ 1$) and return eccentric anomaly $E$.
Explanation
In order to find the position of a body in elliptic motion (e.g., in the two-body problem) at a given time $t$, one has to solve the Kepler's equation
\[M(t) = E(t) - e \sin E(t)\]
where $M(t) = (t - t_0)/P$ is the mean anomaly, $E(t)$ the eccentric anomaly, $e$ the eccentricity of the orbit, $t_0$ is the time of periapsis passage, and $P$ is the period of the orbit. Usually the eccentricity is given and one wants to find the eccentric anomaly $E(t)$ at a specific time $t$, so that also the mean anomaly $M(t)$ is known.
Arguments
M: mean anomaly.
e: eccentricity, in the elliptic motion regime ($0 ≤ e ≤ 1$)
Output
The eccentric anomaly $E$, restricted to the range $[-π, π]$.
Method
Many different numerical methods exist to solve Kepler's equation. This function implements the algorithm proposed in Markley (1995) Celestial Mechanics and Dynamical Astronomy, 63, 101 (DOI:10.1007/BF00691917). This method is not iterative, requires only four transcendental function evaluations, and has been proved to be fast and efficient over the entire range of elliptic motion $0 ≤ e ≤ 1$.
Example
Find the eccentric anomaly for an orbit with eccentricity $e = 0.7$ and for $M(t) = 8π/3$.
julia> using AstroLib
julia> ecc = 0.7;
@@ -303,8 +303,8 @@
M = range(0, stop=2pi, length=1001)[1:end-1];
for ecc in (0, 0.5, 0.9)
plot(M, mod2pi.(kepler_solver.(M, ecc)))
-end
Notes
The true anomaly can be calculated with trueanom function.
Create a 1-d convolution kernel to broaden a spectrum from a rotating star.
Explanation
Can be used to derive the broadening effect (LSF, line spread function) due to rotation on a synthetic stellar spectrum. Assumes constant limb darkening across the disk.
Arguments
delta_v: the step increment (in km/s) in the output rotation kernel
v_sin_i: the rotational velocity projected along the line of sight (km/s)
epsilon (optional numeric argument): the limb-darkening coefficient, default = 0.6 which is typical for photospheric lines. The specific intensity $I$ at any angle $θ$ from the specific intensity $I_{\text{cen}}$ at the center of the disk is given by:
\[I = I_{\text{cen}} ⋅ (1 - ε ⋅ (1 - \cos(θ)))\]
Output
The 2-tuple (velocity_grid, lsf):
velocity_grid: vector of velocity grid with the same number of elements as lsf (see below)
lsf: the convolution kernel vector for the specified rotational velocity. The number of points in lsf will be always be odd (the kernel is symmetric) and equal to either ceil(2*v_sin_i/delta_v) or ceil(2*v_sin_i/delta_v) + 1, whichever number is odd. Elements of lsf will always be of type AbstractFloat. To actually compute the broadening, the spectrum should be convolved with the rotational lsf
Example
Plot the line spread function for a star rotating at 90 km/s in velocity space every 3 km/s. Use Plots.jl for plotting.
using Plots
-plot(lsf_rotate(3, 90)...)
Notes
Code of this function is based on IDL Astronomy User's Library.
Convert from magnitudes to flux expressed in erg/(s cm² Å).
Explanation
This is the reverse of flux2mag.
Arguments
mag: the magnitude to be converted in flux.
zero_point: the zero point level of the magnitude. If not supplied then defaults to 21.1 (Code et al 1976). Ignored if the ABwave keyword is supplied
ABwave (optional numeric keyword): wavelength, in Angstroms. If supplied, then the input mag is assumed to contain Oke AB magnitudes (Oke & Gunn 1983, ApJ, 266, 713; http://adsabs.harvard.edu/abs/1983ApJ...266..713O).
Output
The flux.
If the ABwave keyword is set, then the flux is given by the expression
Create a 1-d convolution kernel to broaden a spectrum from a rotating star.
Explanation
Can be used to derive the broadening effect (LSF, line spread function) due to rotation on a synthetic stellar spectrum. Assumes constant limb darkening across the disk.
Arguments
delta_v: the step increment (in km/s) in the output rotation kernel
v_sin_i: the rotational velocity projected along the line of sight (km/s)
epsilon (optional numeric argument): the limb-darkening coefficient, default = 0.6 which is typical for photospheric lines. The specific intensity $I$ at any angle $θ$ from the specific intensity $I_{\text{cen}}$ at the center of the disk is given by:
\[I = I_{\text{cen}} ⋅ (1 - ε ⋅ (1 - \cos(θ)))\]
Output
The 2-tuple (velocity_grid, lsf):
velocity_grid: vector of velocity grid with the same number of elements as lsf (see below)
lsf: the convolution kernel vector for the specified rotational velocity. The number of points in lsf will be always be odd (the kernel is symmetric) and equal to either ceil(2*v_sin_i/delta_v) or ceil(2*v_sin_i/delta_v) + 1, whichever number is odd. Elements of lsf will always be of type AbstractFloat. To actually compute the broadening, the spectrum should be convolved with the rotational lsf
Example
Plot the line spread function for a star rotating at 90 km/s in velocity space every 3 km/s. Use Plots.jl for plotting.
using Plots
+plot(lsf_rotate(3, 90)...)
Notes
Code of this function is based on IDL Astronomy User's Library.
Convert from magnitudes to flux expressed in erg/(s cm² Å).
Explanation
This is the reverse of flux2mag.
Arguments
mag: the magnitude to be converted in flux.
zero_point: the zero point level of the magnitude. If not supplied then defaults to 21.1 (Code et al 1976). Ignored if the ABwave keyword is supplied
ABwave (optional numeric keyword): wavelength, in Angstroms. If supplied, then the input mag is assumed to contain Oke AB magnitudes (Oke & Gunn 1983, ApJ, 266, 713; http://adsabs.harvard.edu/abs/1983ApJ...266..713O).
Output
The flux.
If the ABwave keyword is set, then the flux is given by the expression
julia> using AstroLib
julia> mean_obliquity(jdcnv(1978,01,7,11, 01))
-0.4091425159336512
Notes
The algorithm used to find the mean obliquity(eps0) is mentioned in USNO Circular 179, but the canonical reference for the IAU adoption is apparently Hilton et al., 2006, Celest.Mech.Dyn.Astron. 94, 351. 2000
The algorithm used to find the mean obliquity(eps0) is mentioned in USNO Circular 179, but the canonical reference for the IAU adoption is apparently Hilton et al., 2006, Celest.Mech.Dyn.Astron. 94, 351. 2000
Compute the right ascension and declination of the Moon at specified Julian date.
Arguments
jd: the Julian ephemeris date. It can be either a scalar or an array
radians (optional boolean keyword): if set to true, then all output angular quantities are given in radians rather than degrees. The default is false
Output
The 5-tuple (ra, dec, dis, geolong, geolat):
ra: apparent right ascension of the Moon in degrees, referred to the true equator of the specified date(s)
dec: the declination of the Moon in degrees
dis: the distance between the centre of the Earth and the centre of the Moon in kilometers
geolong: apparent longitude of the moon in degrees, referred to the ecliptic of the specified date(s)
geolat: apparent longitude of the moon in degrees, referred to the ecliptic of the specified date(s)
If jd is an array, then all output quantities are arrays of the same length as jd.
Method
Derived from the Chapront ELP2000/82 Lunar Theory (Chapront-Touze' and Chapront, 1983, 124, 50), as described by Jean Meeus in Chapter 47 of ``Astronomical Algorithms'' (Willmann-Bell, Richmond), 2nd edition, 1998. Meeus quotes an approximate accuracy of 10" in longitude and 4" in latitude, but he does not give the time range for this accuracy.
Comparison of the IDL procedure with the example in ``Astronomical Algorithms'' reveals a very small discrepancy (~1 km) in the distance computation, but no difference in the position calculation.
Compute the right ascension and declination of the Moon at specified Julian date.
Arguments
jd: the Julian ephemeris date. It can be either a scalar or an array
radians (optional boolean keyword): if set to true, then all output angular quantities are given in radians rather than degrees. The default is false
Output
The 5-tuple (ra, dec, dis, geolong, geolat):
ra: apparent right ascension of the Moon in degrees, referred to the true equator of the specified date(s)
dec: the declination of the Moon in degrees
dis: the distance between the centre of the Earth and the centre of the Moon in kilometers
geolong: apparent longitude of the moon in degrees, referred to the ecliptic of the specified date(s)
geolat: apparent longitude of the moon in degrees, referred to the ecliptic of the specified date(s)
If jd is an array, then all output quantities are arrays of the same length as jd.
Method
Derived from the Chapront ELP2000/82 Lunar Theory (Chapront-Touze' and Chapront, 1983, 124, 50), as described by Jean Meeus in Chapter 47 of ``Astronomical Algorithms'' (Willmann-Bell, Richmond), 2nd edition, 1998. Meeus quotes an approximate accuracy of 10" in longitude and 4" in latitude, but he does not give the time range for this accuracy.
Comparison of the IDL procedure with the example in ``Astronomical Algorithms'' reveals a very small discrepancy (~1 km) in the distance computation, but no difference in the position calculation.
Return the illuminated fraction of the Moon at given Julian date(s).
Arguments
jd: the Julian ephemeris date.
Output
The illuminated fraction $k$ of Moon's disk, with $0 ≤ k ≤ 1$. $k = 0$ indicates a new moon, while $k = 1$ stands for a full moon.
Method
Algorithm from Chapter 46 of "Astronomical Algorithms" by Jean Meeus (Willmann-Bell, Richmond) 1991. sunpos and moonpos are used to get positions of the Sun and the Moon, and the Moon distance. The selenocentric elongation of the Earth from the Sun (phase angle) is then computed, and used to determine the illuminated fraction.
Example
Plot the illuminated fraction of the Moon for every day in January 2018 with a hourly sampling. Use Plots.jl for plotting
using Dates
+plot(points, moonpos(jdcnv.(points))[3])
Notes
Code of this function is based on IDL Astronomy User's Library.
Return the illuminated fraction of the Moon at given Julian date(s).
Arguments
jd: the Julian ephemeris date.
Output
The illuminated fraction $k$ of Moon's disk, with $0 ≤ k ≤ 1$. $k = 0$ indicates a new moon, while $k = 1$ stands for a full moon.
Method
Algorithm from Chapter 46 of "Astronomical Algorithms" by Jean Meeus (Willmann-Bell, Richmond) 1991. sunpos and moonpos are used to get positions of the Sun and the Moon, and the Moon distance. The selenocentric elongation of the Earth from the Sun (phase angle) is then computed, and used to determine the illuminated fraction.
Example
Plot the illuminated fraction of the Moon for every day in January 2018 with a hourly sampling. Use Plots.jl for plotting
using Dates
using Plots
points = DateTime(2018,01,01):Dates.Hour(1):DateTime(2018,01,31,23,59,59);
-plot(points, mphase.(jdcnv.(points)))
Note that in this calendar month there are two full moons, this event is called blue moon.
Notes
Code of this function is based on IDL Astronomy User's Library.
Return the nutation in longitude and obliquity for a given Julian date.
Arguments
jd: Julian ephemeris date, it can be either a scalar or a vector
Output
The 2-tuple (long, obliq), where
long: the nutation in longitude
obl: the nutation in latitude
If jd is an array, long and obl are arrays of the same length.
Method
Uses the formula in Chapter 22 of "Astronomical Algorithms" by Jean Meeus (1998, 2nd ed.) which is based on the 1980 IAU Theory of Nutation and includes all terms larger than 0.0003".
Example
Find the nutation in longitude and obliquity 1987 on Apr 10 at 0h. This is example 22.a from Meeus
julia> using AstroLib
+plot(points, mphase.(jdcnv.(points)))
Note that in this calendar month there are two full moons, this event is called blue moon.
Notes
Code of this function is based on IDL Astronomy User's Library.
Return the nutation in longitude and obliquity for a given Julian date.
Arguments
jd: Julian ephemeris date, it can be either a scalar or a vector
Output
The 2-tuple (long, obliq), where
long: the nutation in longitude
obl: the nutation in latitude
If jd is an array, long and obl are arrays of the same length.
Method
Uses the formula in Chapter 22 of "Astronomical Algorithms" by Jean Meeus (1998, 2nd ed.) which is based on the 1980 IAU Theory of Nutation and includes all terms larger than 0.0003".
Example
Find the nutation in longitude and obliquity 1987 on Apr 10 at 0h. This is example 22.a from Meeus
julia> using AstroLib
julia> jd = jdcnv(1987, 4, 10);
@@ -352,7 +352,7 @@
using Plots
years = DateTime(2000):DateTime(2100);
long, obl = nutate(jdcnv.(years));
-plot(years, long); plot(years, obl)
You can see both the dominant large scale period of nutation, of 18.6 years, and smaller oscillations with shorter periods.
Notes
Code of this function is based on IDL Astronomy User's Library.
Find right ascension and declination for the planets when provided a date as input.
Explanation
This function uses the helio to get the heliocentric ecliptic coordinates of the planets at the given date which it then converts these to geocentric ecliptic coordinates ala "Astronomical Algorithms" by Jean Meeus (1991, p 209). These are then converted to right ascension and declination using euler.
The accuracy between the years 1800 and 2050 is better than 1 arcminute for the terrestial planets, but reaches 10 arcminutes for Saturn. Before 1850 or after 2050 the accuracy can get much worse.
Arguments
date: Can be either a single date or an array of dates. Each element can be either a DateTime type or Julian Date. It can be a scalar or vector.
num: integer denoting planet number, scalar or vector 1 = Mercury, 2 = Venus, ... 9 = Pluto. If not in that change, then the program will throw an error.
Output
ra: right ascension of planet(J2000), in degrees
dec: declination of the planet(J2000), in degrees
Example
Find the RA, Dec of Venus on 1992 Dec 20
julia> using AstroLib, Dates
+plot(wavelength, flux)
Notes
planck_freq calculates the flux of a black body per unit frequency.
Code of this function is based on IDL Astronomy User's Library.
Find right ascension and declination for the planets when provided a date as input.
Explanation
This function uses the helio to get the heliocentric ecliptic coordinates of the planets at the given date which it then converts these to geocentric ecliptic coordinates ala "Astronomical Algorithms" by Jean Meeus (1991, p 209). These are then converted to right ascension and declination using euler.
The accuracy between the years 1800 and 2050 is better than 1 arcminute for the terrestial planets, but reaches 10 arcminutes for Saturn. Before 1850 or after 2050 the accuracy can get much worse.
Arguments
date: Can be either a single date or an array of dates. Each element can be either a DateTime type or Julian Date. It can be a scalar or vector.
num: integer denoting planet number, scalar or vector 1 = Mercury, 2 = Venus, ... 9 = Pluto. If not in that change, then the program will throw an error.
Convert 2D polar coordinates to rectangular coordinates.
Explanation
This is the partial inverse function of recpol.
Arguments
radius: radial coordinate of the point. It may be a scalar or an array.
angle: the angular coordinate of the point. It may be a scalar or an array of the same lenth as radius.
degrees (optional boolean keyword): if true, the angle is assumed to be in degrees, otherwise in radians. It defaults to false.
Mandatory arguments can also be passed as the 2-tuple (radius, angle), so that it is possible to execute recpol(polrec(radius, angle)).
Output
A 2-tuple (x, y) with the rectangular coordinate of the input. If radius and angle are arrays, x and y are arrays of the same length as radius and angle.
Example
Get rectangular coordinates $(x, y)$ of the point with polar coordinates $(r, φ) = (1.7, 227)$, with angle $φ$ expressed in degrees.
julia> using AstroLib
+" 21 05 02.8 -18 51 41"
Notes
Code of this function is based on IDL Astronomy User's Library.
Convert 2D polar coordinates to rectangular coordinates.
Explanation
This is the partial inverse function of recpol.
Arguments
radius: radial coordinate of the point. It may be a scalar or an array.
angle: the angular coordinate of the point. It may be a scalar or an array of the same lenth as radius.
degrees (optional boolean keyword): if true, the angle is assumed to be in degrees, otherwise in radians. It defaults to false.
Mandatory arguments can also be passed as the 2-tuple (radius, angle), so that it is possible to execute recpol(polrec(radius, angle)).
Output
A 2-tuple (x, y) with the rectangular coordinate of the input. If radius and angle are arrays, x and y are arrays of the same length as radius and angle.
Example
Get rectangular coordinates $(x, y)$ of the point with polar coordinates $(r, φ) = (1.7, 227)$, with angle $φ$ expressed in degrees.
julia> using AstroLib
julia> x, y = polrec(1.7, 227, degrees=true)
-(-1.1593972121062475, -1.2433012927525897)
Compute rigorous position angle of point 2 relative to point 1.
Explanation
Computes the rigorous position angle of point 2 (with given right ascension and declination) using point 1 (with given right ascension and declination) as the center.
Arguments
units: integer, can be either 0, or 1, or 2. Describes units of inputs and output:
0: everything (input right ascensions and declinations, and output distance) is radians
1: right ascensions are in decimal hours, declinations in decimal degrees, output distance in degrees
2: right ascensions and declinations are in degrees, output distance in degrees
ra1: right ascension or longitude of point 1
dec1: declination or latitude of point 1
ra2: right ascension or longitude of point 2
dec2: declination or latitude of point 2
Both ra1 and dec1, and ra2 and dec2 can be given as 2-tuples (ra1, dec1) and (ra2, dec2).
Output
Angle of the great circle containing [ra2, dec2] from the meridian containing [ra1, dec1], in the sense north through east rotating about [ra1, dec1]. See units argument above for units.
Method
The "four-parts formula" from spherical trigonometry (p. 12 of Smart's Spherical Astronomy or p. 12 of Green' Spherical Astronomy).
Example
Mizar has coordinates (ra, dec) = (13h 23m 55.5s, +54° 55' 31''). Its companion, Alcor, has coordinates (ra, dec) = (13h 25m 13.5s, +54° 59' 17''). Find the position angle of Alcor with respect to Mizar.
julia> using AstroLib
+(-1.1593972121062475, -1.2433012927525897)
Compute rigorous position angle of point 2 relative to point 1.
Explanation
Computes the rigorous position angle of point 2 (with given right ascension and declination) using point 1 (with given right ascension and declination) as the center.
Arguments
units: integer, can be either 0, or 1, or 2. Describes units of inputs and output:
0: everything (input right ascensions and declinations, and output distance) is radians
1: right ascensions are in decimal hours, declinations in decimal degrees, output distance in degrees
2: right ascensions and declinations are in degrees, output distance in degrees
ra1: right ascension or longitude of point 1
dec1: declination or latitude of point 1
ra2: right ascension or longitude of point 2
dec2: declination or latitude of point 2
Both ra1 and dec1, and ra2 and dec2 can be given as 2-tuples (ra1, dec1) and (ra2, dec2).
Output
Angle of the great circle containing [ra2, dec2] from the meridian containing [ra1, dec1], in the sense north through east rotating about [ra1, dec1]. See units argument above for units.
Method
The "four-parts formula" from spherical trigonometry (p. 12 of Smart's Spherical Astronomy or p. 12 of Green' Spherical Astronomy).
Example
Mizar has coordinates (ra, dec) = (13h 23m 55.5s, +54° 55' 31''). Its companion, Alcor, has coordinates (ra, dec) = (13h 25m 13.5s, +54° 59' 17''). Find the position angle of Alcor with respect to Mizar.
The default (ra, dec) system is FK5 based on epoch J2000.0 but FK4 based on B1950.0 is available via the FK4 boolean keyword.
Arguments
ra: input right ascension, scalar or vector, in degrees, unless the radians keyword is set to true
dec: input declination, scalar or vector, in degrees, unless the radians keyword is set to true
equinox1: original equinox of coordinates, numeric scalar.
equinox2: equinox of precessed coordinates.
FK4 (optional boolean keyword): if this keyword is set to true, the FK4 (B1950.0) system precession angles are used to compute the precession matrix. When it is false, the default, use FK5 (J2000.0) precession angles.
radians (optional boolean keyword): if this keyword is set to true, then the input and output right ascension and declination vectors are in radians rather than degrees.
Output
The 2-tuple (ra, dec) of coordinates modified by precession.
Example
The Pole Star has J2000.0 coordinates (2h, 31m, 46.3s, 89d 15' 50.6"); compute its coordinates at J1985.0
julia> using AstroLib
+-108.46011246802047
Notes
The function sphdist provides an alternate method of computing a spherical
distance.
Note that posang is not commutative: the position angle between A and B is $θ$, then the position angle between B and A is $180 + θ$.
Code of this function is based on IDL Astronomy User's Library.
The default (ra, dec) system is FK5 based on epoch J2000.0 but FK4 based on B1950.0 is available via the FK4 boolean keyword.
Arguments
ra: input right ascension, scalar or vector, in degrees, unless the radians keyword is set to true
dec: input declination, scalar or vector, in degrees, unless the radians keyword is set to true
equinox1: original equinox of coordinates, numeric scalar.
equinox2: equinox of precessed coordinates.
FK4 (optional boolean keyword): if this keyword is set to true, the FK4 (B1950.0) system precession angles are used to compute the precession matrix. When it is false, the default, use FK5 (J2000.0) precession angles.
radians (optional boolean keyword): if this keyword is set to true, then the input and output right ascension and declination vectors are in radians rather than degrees.
Output
The 2-tuple (ra, dec) of coordinates modified by precession.
Example
The Pole Star has J2000.0 coordinates (2h, 31m, 46.3s, 89d 15' 50.6"); compute its coordinates at J1985.0
Algorithm from "Computational Spherical Astronomy" by Taff (1983), p. 24. (FK4). FK5 constants from "Explanatory Supplement To The Astronomical Almanac" 1992, page 104 Table 3.211.1 (https://archive.org/details/131123ExplanatorySupplementAstronomicalAlmanac).
Notes
Accuracy of precession decreases for declination values near 90 degrees. precess should not be used more than 2.5 centuries from 2000 on the FK5 system (1950.0 on the FK4 system). If you need better accuracy, use bprecess or jprecess as needed.
Code of this function is based on IDL Astronomy User's Library.
precess_cd(cd, epoch1, epoch2, crval_old, crval_new[, FK4=true]) -> cd
Purpose
Precess the coordinate description matrix.
Explanation
The coordinate matrix is precessed from epoch1 to epoch2.
Arguments
cd: 2×2 coordinate description matrix in degrees
epoch1: original equinox of coordinates, scalar
epoch2: equinox of precessed coordinates, scalar
crval_old: 2 element vector containing right ascension and declination in degrees of the reference pixel in the original equinox
crval_new: 2 element vector giving crval in the new equinox
FK4 (optional boolean keyword): if this keyword is set to true, then the precession constants are taken in the FK4 reference frame. When it is false, the default is the FK5 frame
Algorithm from "Computational Spherical Astronomy" by Taff (1983), p. 24. (FK4). FK5 constants from "Explanatory Supplement To The Astronomical Almanac" 1992, page 104 Table 3.211.1 (https://archive.org/details/131123ExplanatorySupplementAstronomicalAlmanac).
Notes
Accuracy of precession decreases for declination values near 90 degrees. precess should not be used more than 2.5 centuries from 2000 on the FK5 system (1950.0 on the FK4 system). If you need better accuracy, use bprecess or jprecess as needed.
Code of this function is based on IDL Astronomy User's Library.
precess_cd(cd, epoch1, epoch2, crval_old, crval_new[, FK4=true]) -> cd
Purpose
Precess the coordinate description matrix.
Explanation
The coordinate matrix is precessed from epoch1 to epoch2.
Arguments
cd: 2×2 coordinate description matrix in degrees
epoch1: original equinox of coordinates, scalar
epoch2: equinox of precessed coordinates, scalar
crval_old: 2 element vector containing right ascension and declination in degrees of the reference pixel in the original equinox
crval_new: 2 element vector giving crval in the new equinox
FK4 (optional boolean keyword): if this keyword is set to true, then the precession constants are taken in the FK4 reference frame. When it is false, the default is the FK5 frame
Code of this function is based on IDL Astronomy User's Library. This function should not be used for values more than 2.5 centuries from the year 1900. This function calls sec2rad, precess and bprecess.
x, y, z: scalars or vectors giving heliocentric rectangular coordinates.
equinox1: original equinox of coordinates, numeric scalar.
equinox2: equinox of precessed coordinates, numeric scalar.
Input coordinates can be given also a 3-tuple (x, y, z).
Output
The 3-tuple (x, y, z) of coordinates modified by precession.
Example
Precess 2000 equinox coordinates (1, 1, 1) to 2050.
julia> using AstroLib
+ 110.188 110.365
Notes
Code of this function is based on IDL Astronomy User's Library. This function should not be used for values more than 2.5 centuries from the year 1900. This function calls sec2rad, precess and bprecess.
The equatorial geocentric rectangular coordinates are converted to right ascension and declination, precessed in the normal way, then changed back to x, y and z using unit vectors.
Notes
Code of this function is based on IDL Astronomy User's Library.
Return the precession matrix needed to go from equinox1 to equinox2.
Explanation
This matrix is used by precess and baryvel to precess astronomical coordinates.
Arguments
equinox1: original equinox of coordinates.
equinox2: equinox of precessed coordinates.
FK4 (optional boolean keyword): if this keyword is set to true, the FK4 (B1950.0) system precession angles are used to compute the precession matrix. When it is false, the default, use FK5 (J2000.0) precession angles.
Output
A 3×3 matrix, used to precess equatorial rectangular coordinates.
Example
Return the precession matrix from 1950.0 to 1975.0 in the FK4 system
julia> using AstroLib
+(0.9838854500981734, 1.0110925876508692, 1.0048189888146941)
Method
The equatorial geocentric rectangular coordinates are converted to right ascension and declination, precessed in the normal way, then changed back to x, y and z using unit vectors.
Notes
Code of this function is based on IDL Astronomy User's Library.
Return the precession matrix needed to go from equinox1 to equinox2.
Explanation
This matrix is used by precess and baryvel to precess astronomical coordinates.
Arguments
equinox1: original equinox of coordinates.
equinox2: equinox of precessed coordinates.
FK4 (optional boolean keyword): if this keyword is set to true, the FK4 (B1950.0) system precession angles are used to compute the precession matrix. When it is false, the default, use FK5 (J2000.0) precession angles.
Output
A 3×3 matrix, used to precess equatorial rectangular coordinates.
Example
Return the precession matrix from 1950.0 to 1975.0 in the FK4 system
julia> using AstroLib
julia> premat(1950, 1975, FK4=true)
3×3 StaticArraysCore.SMatrix{3, 3, Float64, 9} with indices SOneTo(3)×SOneTo(3):
0.999981 -0.00558775 -0.00242909
0.00558775 0.999984 -6.78691e-6
- 0.00242909 -6.78633e-6 0.999997
Method
FK4 constants from "Computational Spherical Astronomy" by Taff (1983), p. 24. (FK4). FK5 constants from "Explanatory Supplement To The Astronomical Almanac" 1992, page 104 Table 3.211.1 (https://archive.org/details/131123ExplanatorySupplementAstronomicalAlmanac).
Notes
Code of this function is based on IDL Astronomy User's Library.
julia> using AstroLib
+ 0.00242909 -6.78633e-6 0.999997
Method
FK4 constants from "Computational Spherical Astronomy" by Taff (1983), p. 24. (FK4). FK5 constants from "Explanatory Supplement To The Astronomical Almanac" 1992, page 104 Table 3.211.1 (https://archive.org/details/131123ExplanatorySupplementAstronomicalAlmanac).
Notes
Code of this function is based on IDL Astronomy User's Library.
Convert 2D rectangular coordinates to polar coordinates.
Explanation
This is the partial inverse function of polrec.
Arguments
x: the abscissa coordinate of the point. It may be a scalar or an array.
y: the ordinate coordinate of the point. It may be a scalar or an array of the same lenth as x.
degrees (optional boolean keyword): if true, the output angle is given in degrees, otherwise in radians. It defaults to false.
Mandatory arguments may also be passed as the 2-tuple (x, y), so that it is possible to execute polrec(recpol(x, y)).
Output
A 2-tuple (radius, angle) with the polar coordinates of the input. The coordinate angle coordinate lies in the range $[-π, π]$ if degrees=false, or $[-180, 180]$ when degrees=true.
If x and y are arrays, radius and angle are arrays of the same length as radius and angle.
Example
Calculate polar coordinates $(r, φ)$ of point with rectangular coordinates $(x, y) = (2.24, -1.87)$.
julia> using AstroLib
+(6.0, 45.0, 9.0, -16.0, 42.0, 57.9600000000064)
Convert 2D rectangular coordinates to polar coordinates.
Explanation
This is the partial inverse function of polrec.
Arguments
x: the abscissa coordinate of the point. It may be a scalar or an array.
y: the ordinate coordinate of the point. It may be a scalar or an array of the same lenth as x.
degrees (optional boolean keyword): if true, the output angle is given in degrees, otherwise in radians. It defaults to false.
Mandatory arguments may also be passed as the 2-tuple (x, y), so that it is possible to execute polrec(recpol(x, y)).
Output
A 2-tuple (radius, angle) with the polar coordinates of the input. The coordinate angle coordinate lies in the range $[-π, π]$ if degrees=false, or $[-180, 180]$ when degrees=true.
If x and y are arrays, radius and angle are arrays of the same length as radius and angle.
Example
Calculate polar coordinates $(r, φ)$ of point with rectangular coordinates $(x, y) = (2.24, -1.87)$.
Calculate the separation and position angle of a binary star.
Explanation
This function will return the separation $ρ$ and position angle $θ$ of a visual binary star derived from its orbital elements. The algorithms described in the following book will be used: Meeus J., 1992, Astronomische Algorithmen, Barth. Compared to the examples given at page 400 and no discrepancy found.
Arguments
period: period [year]
periastro: time of periastron passage [year]
eccentricity: eccentricity of the orbit
semimajor_axis: semi-major axis [arc second]
inclination: inclination angle [degree]
omega: node [degree]
omega2: longitude of periastron [degree]
epoch: epoch of observation [year]
All input parameters have to be scalars.
Output
The 2-tuple $(ρ, θ)$, where
$ρ$ is separation [arc second], and
$θ$ is position angle (degree).
Example
Find the position of Eta Coronae Borealis at the epoch 2016
julia> using AstroLib
+(2.9179616172938263, -0.6956158538564537)
Calculate the separation and position angle of a binary star.
Explanation
This function will return the separation $ρ$ and position angle $θ$ of a visual binary star derived from its orbital elements. The algorithms described in the following book will be used: Meeus J., 1992, Astronomische Algorithmen, Barth. Compared to the examples given at page 400 and no discrepancy found.
Arguments
period: period [year]
periastro: time of periastron passage [year]
eccentricity: eccentricity of the orbit
semimajor_axis: semi-major axis [arc second]
inclination: inclination angle [degree]
omega: node [degree]
omega2: longitude of periastron [degree]
epoch: epoch of observation [year]
All input parameters have to be scalars.
Output
The 2-tuple $(ρ, θ)$, where
$ρ$ is separation [arc second], and
$θ$ is position angle (degree).
Example
Find the position of Eta Coronae Borealis at the epoch 2016
degrees (optional boolean keyword): if true, all angles, including the output distance, are assumed to be in degrees, otherwise they are all in radians. It defaults to false.
Output
Angular distance on a sphere between points 1 and 2, as an AbstractFloat. It is expressed in radians unless degrees keyword is set to true.
Example
julia> using AstroLib
+ 54.0
Notes
Code of this function is based on IDL Astronomy User's Library.
degrees (optional boolean keyword): if true, all angles, including the output distance, are assumed to be in degrees, otherwise they are all in radians. It defaults to false.
Output
Angular distance on a sphere between points 1 and 2, as an AbstractFloat. It is expressed in radians unless degrees keyword is set to true.
Example
julia> using AstroLib
julia> sphdist(120, -43, 175, +22)
-1.5904422616007134
Notes
gcirc function is similar to sphdist, but may be more suitable for astronomical applications.
Code of this function is based on IDL Astronomy User's Library.
Compute the right ascension and declination of the Sun at a given date.
Arguments
jd: the Julian date of when you want to calculate Sun position. It can be either a scalar or a vector. Use jdcnv to get the Julian date for a given date and time.
radians (optional boolean keyword): if set to true, all output quantities are given in radians. The default is false, so all quantities are given in degrees.
Output
The 4-tuple (ra, dec, elong, obliquity):
ra: the right ascension of the Sun at that date
dec: the declination of the Sun at that date
elong: ecliptic longitude of the Sun at that date
obliquity: the obliquity of the ecliptic
All quantities are given in degrees, unless radians keyword is set to true (see "Arguments" section). If jd is an array, arrays of the same given as jd are returned.
Method
Uses a truncated version of Newcomb's Sun. Adapted from the IDL routine SUN_POS by CD Pike, which was adapted from a FORTRAN routine by B. Emerson (RGO).
Example
Find the apparent right ascension and declination of the Sun on May 1, 1982
julia> using AstroLib
+1.5904422616007134
Notes
gcirc function is similar to sphdist, but may be more suitable for astronomical applications.
Code of this function is based on IDL Astronomy User's Library.
Compute the right ascension and declination of the Sun at a given date.
Arguments
jd: the Julian date of when you want to calculate Sun position. It can be either a scalar or a vector. Use jdcnv to get the Julian date for a given date and time.
radians (optional boolean keyword): if set to true, all output quantities are given in radians. The default is false, so all quantities are given in degrees.
Output
The 4-tuple (ra, dec, elong, obliquity):
ra: the right ascension of the Sun at that date
dec: the declination of the Sun at that date
elong: ecliptic longitude of the Sun at that date
obliquity: the obliquity of the ecliptic
All quantities are given in degrees, unless radians keyword is set to true (see "Arguments" section). If jd is an array, arrays of the same given as jd are returned.
Method
Uses a truncated version of Newcomb's Sun. Adapted from the IDL routine SUN_POS by CD Pike, which was adapted from a FORTRAN routine by B. Emerson (RGO).
Example
Find the apparent right ascension and declination of the Sun on May 1, 1982
The Astronomical Almanac gives 02 31 32.58 +14 54 34.9 so the error for this case is < 0.5".
Plot the apparent right ascension, in hours, and declination of the Sun, in degrees, for every day in 2016. Use Plots.jl for plotting.
using Plots
@@ -445,7 +445,7 @@
days = DateTime(2016):Day(1):DateTime(2016, 12, 31);
ra, declin = sunpos(jdcnv.(days));
-plot(days, ra/15); plot(days, declin)
Notes
Patrick Wallace (Rutherford Appleton Laboratory, UK) has tested the accuracy of a C adaptation of the present algorithm and found the following results. From 1900-2100 sunpos gave 7.3 arcsec maximum error, 2.6 arcsec RMS. Over the shorter interval 1950-2050 the figures were 6.4 arcsec max, 2.2 arcsec RMS.
The returned ra and dec are in the given date's equinox.
Code of this function is based on IDL Astronomy User's Library.
Patrick Wallace (Rutherford Appleton Laboratory, UK) has tested the accuracy of a C adaptation of the present algorithm and found the following results. From 1900-2100 sunpos gave 7.3 arcsec maximum error, 2.6 arcsec RMS. Over the shorter interval 1950-2050 the figures were 6.4 arcsec max, 2.2 arcsec RMS.
The returned ra and dec are in the given date's equinox.
Code of this function is based on IDL Astronomy User's Library.
Converts a sexagesimal number or string to decimal.
Explanation
ten is the inverse of the sixty function.
Arguments
ten takes as argument either three scalars (deg, min, sec) or a string. The string should have the form "deg:min:sec" or "deg min sec". Also any iterable like (deg, min, sec) or [deg, min, sec] is accepted as argument.
If minutes and seconds are not specified they default to zero.
Output
The decimal conversion of the sexagesimal numbers provided is returned.
These functions cannot deal with -0 (negative integer zero) in numeric input. If it is important to give sense to negative zero, you can either make sure to pass a floating point negative zero -0.0 (this is the best option), or use negative minutes and seconds, or non-integer negative degrees and minutes.
Determine the position of the first tic mark for astronomical images.
Explanation
For use in labelling images with right ascension and declination axes. This routine determines the position in pixels of the first tic.
Arguments
zmin: astronomical coordinate value at axis zero point (degrees or hours).
pixx: distance in pixels between tic marks (usually obtained from tics).
incr - increment in minutes for labels (usually an even number obtained from the procedure tics).
ra (optional boolean keyword): if true, incremental value being entered is in minutes of time, else it is assumed that value is in else it's in minutes of arc. Default is false.
Output
The 2 tuple (min2, tic1):
min2: astronomical coordinate value at first tic mark
tic1: position in pixels of first tic mark
Example
Suppose a declination axis has a value of 30.2345 degrees at its zero point. A tic mark is desired every 10 arc minutes, which corresponds to 12.74 pixels, with increment for labels being 10 minutes. Then
julia> using AstroLib
+-10.433333333333334
Notes
These functions cannot deal with -0 (negative integer zero) in numeric input. If it is important to give sense to negative zero, you can either make sure to pass a floating point negative zero -0.0 (this is the best option), or use negative minutes and seconds, or non-integer negative degrees and minutes.
Determine the position of the first tic mark for astronomical images.
Explanation
For use in labelling images with right ascension and declination axes. This routine determines the position in pixels of the first tic.
Arguments
zmin: astronomical coordinate value at axis zero point (degrees or hours).
pixx: distance in pixels between tic marks (usually obtained from tics).
incr - increment in minutes for labels (usually an even number obtained from the procedure tics).
ra (optional boolean keyword): if true, incremental value being entered is in minutes of time, else it is assumed that value is in else it's in minutes of arc. Default is false.
Output
The 2 tuple (min2, tic1):
min2: astronomical coordinate value at first tic mark
tic1: position in pixels of first tic mark
Example
Suppose a declination axis has a value of 30.2345 degrees at its zero point. A tic mark is desired every 10 arc minutes, which corresponds to 12.74 pixels, with increment for labels being 10 minutes. Then
julia> using AstroLib
julia> tic_one(30.2345, 12.74, 10)
-(30.333333333333332, 7.554820000000081)
yields values of min2 ≈ 30.333 and tic1 ≈ 7.55482, i.e. the first tic mark should be labeled 30 deg 20 minutes and be placed at pixel value 7.55482.
Notes
Code of this function is based on IDL Astronomy User's Library.
ticpos(deglen, pixlen, ticsize) -> ticsize, incr, units
Purpose
Specify distance between tic marks for astronomical coordinate overlays.
Explanation
User inputs number an approximate distance between tic marks, and the axis length in degrees. ticpos will return a distance between tic marks such that the separation is a round multiple in arc seconds, arc minutes, or degrees.
Arguments
deglen: length of axis in degrees, positive scalar
pixlen: length of axis in plotting units (pixels), postive scalar
ticsize: distance between tic marks (pixels). This value will be adjusted by ticpos such that the distance corresponds to a round multiple in the astronomical coordinate.
Output
The 3-tuple (ticsize, incr, units):
ticsize: distance between tic marks (pixels), positive scalar
incr: incremental value for tic marks in round units given by the units parameter
units: string giving units of ticsize, either 'Arc Seconds', 'Arc Minutes', or 'Degrees'
Example
Suppose a 512 × 512 image array corresponds to 0.2 × 0.2 degrees on the sky. A tic mark is desired in round angular units, approximately every 75 pixels. Then
julia> using AstroLib
+(30.333333333333332, 7.554820000000081)
yields values of min2 ≈ 30.333 and tic1 ≈ 7.55482, i.e. the first tic mark should be labeled 30 deg 20 minutes and be placed at pixel value 7.55482.
Notes
Code of this function is based on IDL Astronomy User's Library.
ticpos(deglen, pixlen, ticsize) -> ticsize, incr, units
Purpose
Specify distance between tic marks for astronomical coordinate overlays.
Explanation
User inputs number an approximate distance between tic marks, and the axis length in degrees. ticpos will return a distance between tic marks such that the separation is a round multiple in arc seconds, arc minutes, or degrees.
Arguments
deglen: length of axis in degrees, positive scalar
pixlen: length of axis in plotting units (pixels), postive scalar
ticsize: distance between tic marks (pixels). This value will be adjusted by ticpos such that the distance corresponds to a round multiple in the astronomical coordinate.
Output
The 3-tuple (ticsize, incr, units):
ticsize: distance between tic marks (pixels), positive scalar
incr: incremental value for tic marks in round units given by the units parameter
units: string giving units of ticsize, either 'Arc Seconds', 'Arc Minutes', or 'Degrees'
Example
Suppose a 512 × 512 image array corresponds to 0.2 × 0.2 degrees on the sky. A tic mark is desired in round angular units, approximately every 75 pixels. Then
Compute a nice increment between tic marks for astronomical images.
Explanation
For use in labelling a displayed image with right ascension or declination axes. An approximate distance between tic marks is input, and a new value is computed such that the distance between tic marks is in simple increments of the tic label values.
Arguements
radec_min : minimum axis value (degrees).
radec_min : maximum axis value (degrees).
numx : number of pixels in x direction.
ticsize : distance between tic marks (pixels).
ra (optional boolean keyword): if true, incremental value would be in minutes of time. Default is false.
Output
A 2-tuple (ticsize, incr):
ticsize : distance between tic marks (pixels).
incr : incremental value for tic labels. The format is dependent on the optional keyword. If true (i.e for right ascension), it's in minutes of time, else it's in minutes of arc (for declination).
Example
julia> using AstroLib
+(85.33333333333333, 2.0, "Arc Minutes")
i.e. a good tic mark spacing is every 2 arc minutes, corresponding to 85.333 pixels.
Notes
All the arguments taken as input are assumed to be positive in nature.
Code of this function is based on IDL Astronomy User's Library.
Compute a nice increment between tic marks for astronomical images.
Explanation
For use in labelling a displayed image with right ascension or declination axes. An approximate distance between tic marks is input, and a new value is computed such that the distance between tic marks is in simple increments of the tic label values.
Arguements
radec_min : minimum axis value (degrees).
radec_min : maximum axis value (degrees).
numx : number of pixels in x direction.
ticsize : distance between tic marks (pixels).
ra (optional boolean keyword): if true, incremental value would be in minutes of time. Default is false.
Output
A 2-tuple (ticsize, incr):
ticsize : distance between tic marks (pixels).
incr : incremental value for tic labels. The format is dependent on the optional keyword. If true (i.e for right ascension), it's in minutes of time, else it's in minutes of arc (for declination).
Calculate true anomaly for a particle in elliptic orbit with eccentric anomaly $E$ and eccentricity $e$.
Explanation
In the two-body problem, once that the Kepler's equation is solved and $E(t)$ is determined, the polar coordinates $(r(t), θ(t))$ of the body at time $t$ in the elliptic orbit are given by
Calculate true anomaly for a particle in elliptic orbit with eccentric anomaly $E$ and eccentricity $e$.
Explanation
In the two-body problem, once that the Kepler's equation is solved and $E(t)$ is determined, the polar coordinates $(r(t), θ(t))$ of the body at time $t$ in the elliptic orbit are given by
in which $a$ is the semi-major axis of the orbit, and $θ_0$ the value of angular coordinate at time $t = t_0$.
Arguments
E: eccentric anomaly.
e: eccentricity, in the elliptic motion regime ($0 ≤ e ≤ 1$)
Output
The true anomaly.
Example
Plot the true anomaly as a function of mean anomaly for eccentricity $e = 0, 0.5, 0.9$. Use Plots.jl for plotting.
using Plots
M = range(0, stop=2pi, length=1001)[1:end-1];
p = plot()
for ecc in (0, 0.5, 0.9)
plot!(p, M, mod2pi.(trueanom.(kepler_solver.(M, ecc), ecc)))
end
-p
Notes
The eccentric anomaly can be calculated with kepler_solver function.
hbeta (optional): H-beta line strength index. If it is not supplied, then by default its value will be NaN and the code will estimate a value based on by, m1,and c1. It is not used for stars in group 8.
eby_in (optional): specifies the E(b-y) color to use. If not supplied, then by default its value will be NaN and E(b-y) will be estimated from the Stromgren colors.
Output
te: approximate effective temperature
mv: absolute visible magnitude
eby: color excess E(b-y)
delm0: metallicity index, delta m0, may not be calculable for early B stars and so returns NaN.
radius: stellar radius (R/R(solar))
Example
Determine the stellar parameters of 5 stars given their Stromgren parameters
julia> using AstroLib
+p
Notes
The eccentric anomaly can be calculated with kepler_solver function.
hbeta (optional): H-beta line strength index. If it is not supplied, then by default its value will be NaN and the code will estimate a value based on by, m1,and c1. It is not used for stars in group 8.
eby_in (optional): specifies the E(b-y) color to use. If not supplied, then by default its value will be NaN and E(b-y) will be estimated from the Stromgren colors.
Output
te: approximate effective temperature
mv: absolute visible magnitude
eby: color excess E(b-y)
delm0: metallicity index, delta m0, may not be calculable for early B stars and so returns NaN.
radius: stellar radius (R/R(solar))
Example
Determine the stellar parameters of 5 stars given their Stromgren parameters
Corrects for the index of refraction of air under standard conditions. Wavelength values below $2000 Å$ will not be altered. Uses relation of Ciddor (1996).
Arguments
wave_vacuum: vacuum wavelength in angstroms. Wavelengths are corrected for the index of refraction of air under standard conditions. Wavelength values below $2000 Å$ will not be altered, take care within $[1 Å, 2000 Å]$.
Corrects for the index of refraction of air under standard conditions. Wavelength values below $2000 Å$ will not be altered. Uses relation of Ciddor (1996).
Arguments
wave_vacuum: vacuum wavelength in angstroms. Wavelengths are corrected for the index of refraction of air under standard conditions. Wavelength values below $2000 Å$ will not be altered, take care within $[1 Å, 2000 Å]$.
Calculate geocentric $x$, $y$, and $z$ and velocity coordinates of the Sun.
Explanation
Calculates geocentric $x$, $y$, and $z$ vectors and velocity coordinates ($dx$, $dy$ and $dz$) of the Sun. (The positive $x$ axis is directed towards the equinox, the $y$-axis, towards the point on the equator at right ascension 6h, and the $z$ axis toward the north pole of the equator). Typical position accuracy is $<10^{-4}$ AU (15000 km).
Arguments
jd: number of Reduced Julian Days for the wanted date. It can be either a scalar or a vector.
equinox (optional numeric argument): equinox of output. Default is 1950.
You can use juldate to get the number of Reduced Julian Days for the selected dates.
Output
The 6-tuple $(x, y, z, v_x, v_y, v_z)$, where
$x, y, z$: scalars or vectors giving heliocentric rectangular coordinates (in AU) for each date supplied. Note that $\sqrt{x^2 + y^2 + z^2}$ gives the Earth-Sun distance for the given date.
$v_x, v_y, v_z$: velocity vectors corresponding to $x, y$, and $z$.
Example
What were the rectangular coordinates and velocities of the Sun on 1999-01-22T00:00:00 (= JD 2451200.5) in J2000 coords? Note: Astronomical Almanac (AA) is in TDT, so add 64 seconds to UT to convert.
julia> using AstroLib, Dates
+1999.3526230448367
Notes
airtovac converts air wavelengths to vacuum wavelengths.
Code of this function is based on IDL Astronomy User's Library.
Calculate geocentric $x$, $y$, and $z$ and velocity coordinates of the Sun.
Explanation
Calculates geocentric $x$, $y$, and $z$ vectors and velocity coordinates ($dx$, $dy$ and $dz$) of the Sun. (The positive $x$ axis is directed towards the equinox, the $y$-axis, towards the point on the equator at right ascension 6h, and the $z$ axis toward the north pole of the equator). Typical position accuracy is $<10^{-4}$ AU (15000 km).
Arguments
jd: number of Reduced Julian Days for the wanted date. It can be either a scalar or a vector.
equinox (optional numeric argument): equinox of output. Default is 1950.
You can use juldate to get the number of Reduced Julian Days for the selected dates.
Output
The 6-tuple $(x, y, z, v_x, v_y, v_z)$, where
$x, y, z$: scalars or vectors giving heliocentric rectangular coordinates (in AU) for each date supplied. Note that $\sqrt{x^2 + y^2 + z^2}$ gives the Earth-Sun distance for the given date.
$v_x, v_y, v_z$: velocity vectors corresponding to $x, y$, and $z$.
Example
What were the rectangular coordinates and velocities of the Sun on 1999-01-22T00:00:00 (= JD 2451200.5) in J2000 coords? Note: Astronomical Almanac (AA) is in TDT, so add 64 seconds to UT to convert.
Return the zenith right ascension and declination in radians for a given Julian date or a local civil time and timezone.
Explanation
The local sidereal time is computed with the help of ct2lst, which is the right ascension of the zenith. This and the observatories latitude (corresponding to the declination) are converted to radians and returned as the zenith direction.
Arguments
The function can be called in two different ways. The arguments common to both methods are latitude and longitude:
latitude : latitude of the desired location.
longitude : longitude of the desired location.
The zenith direction can be computed either by providing the Julian date:
jd : the Julian date of the date and time for which the zenith position is desired.
or the time zone and the date:
tz: the time zone (in hours) of the desired location (e.g. 4 = EDT, 5 = EST)
date: the local civil time with type DateTime.
Output
A 2-tuple (ra, dec):
ra : the right ascension (in radians) of the zenith.