This page provides details of a method for precessing the
coordinates of objects from the B1950.0 epoch to the now standard
J2000.0 epoch, ignoring the proper motion of the stars. The algorithm
used can handle coordinates with declinations up to +90 or -90 degrees
(unlike many approximate methods) and the QBASIC program
provided will accept input and provide output in the notation used in
Burnham's Celestial Handbook, which uses B1950.0
coordinates. I have provided a DOS version of the program
compiled with the excellent FirstBasic compiler, for those who need to
convert a number of positions.
The celestial sphere can be mapped by using two main systems of coordinates, each referred to a reference plane. The celestial equator is the reference plane for equatorial coordinates (right ascension and declination), and the ecliptic plane is the basis of ecliptic coordinates (usually heliocentric longitude and latitude). The ecliptic and equator are represented as two great circles in the celestial sphere, and they meet at two points on the sphere. One of these points is used as the zero of right ascension. Declination has a natural zero at the celestial equator.
The Earth spins like a top on its own axis as it orbits the Sun. The Moon and the Sun exert tidal forces on the equatorial bulge of the Earth, and there are small effects from the rest of the planets in the solar system. This complex set of forces results in a movement in the celestial equator with time, which causes the zero point for RA to shift back along the ecliptic, and causes the tilt of the equator to change. This combined movement has short and long period components which are ususally broken down into 'corrections' as follows;
The two precession effects are lumped together into the general precession. The motion of the equator means that you must quote a date with the RA and DEC for an object. If you take account of nutation, then your positions are apparent positions. If you leave out nutation, then your positions are referred to the mean equator of the date you have chosen. All observations are referred to the equinox of the date, or epoch, chosen. Three main epochs are in common use;
Star charts, catalogues and observing lists use one of the epochs B1950.0 or J2000.0, and catalogues may list both. Star chart programs can produce coordinates refered to any epoch, but you should make sure you know which epoch the software is working in. The Yale Bright Star catalogue lists B1900.0 and J2000.0 positions. The 'B' and the 'J' refer to the use of the Besselian year and Julian year to count the epochs - the Besselian years are based on a solar year and are not now used. The Julian years are 365.25 days long, and can be related to the Julian day number easily.
The stars are not really fixed on the inside of a large sphere (even though it seems like that on a good dark night). They orbit the galactic centre and drift in their own local groups and clusters. On time scales which we can comprehend, this movement of the nearer stars can be modelled as the proper motion and is given in arcseconds per year (for both RA and DEC). The star with the largest proper motion is 61-cyg, at 4.1 and 3.2 arcseconds per year in RA and DEC. If you ignore the proper motion of this star when converting from B1950 to J2000, then you will have a J2000 RA which is about 0.2 minutes of time wrong, and a DEC which is about 3 arcminutes wrong. Most of the stars in the Bright Star Atlas have proper motions much less than this (about 17 more than 1 arcsec per year in RA out of 9000 odd stars). Deep space objects do not have proper motions that we can measure!
If you have an observing guide or comet ephemeris in B1950.0 coordinates, and you want to plot the positions of the object on a J2000.0 star atlas, you must correct the B1950 coordinates to J2000.0 by calculating the effects of precession. A rough correction, valid for shorter time periods, is often found in astronomical reference books - here is the recipe for B1950.0 to J2000.0;
1. change B1950.0 RA and DEC to decimal degrees
2. work out the following formulas
RA2000 = RA + 0.640265 + 0.278369 * sin(RA) * tan(DEC)
DEC2000 = DEC + 0.278369 * cos(RA)
3. convert RA back to hours, and format as Hrs, min, sec as
preferred.
The main problem with this recipe is that tan(DEC)
will grow very large for stars with a high declination. You can't
precess the coordinates of the Pole Star using the recipe above. For
objects which stay close to the ecliptic (like planets and the stars
in the 'occultation band') the recipe works well. A minor problem is
that stars with RA's close to 24h in 1950 coordinates may produce a
J2000 result higher than 24h, although that is easily taken care of.
A better method for use with programmable calculators or any form of computer program or spreadsheet is based on the idea of rotations around rectangular axes, so that no overflow can occur. This second method is developed in the next section.
As a practical application, I have written some format conversion
routines which can use the notation found in Burnham's Celestial
Handbook, a much loved observation guide for variable and
double stars and deep space objects. Burnham's uses
B1950.0 coordinates and does not include comprehensive star maps. If
you want to plot the objects mentioned on a J2000.0 star atlas, then
you may find the QBASIC program here useful.
The method here is found in Duffett-Smith, and is also explained in Meeus and (in much more detail) in the Explanatory Supplement.
As a concrete example, take the variable star R-lyr for which Burnham's gives the B1950 coordinates as 18538n4353. Translating from this compact notation, this is RA 18h 53.8 minutes and DEC 43 deg 53 minutes North. Wil Tirion's The Cambridge Star Atlas gives the J2000.0 coordinates of this star as 18 h 55.3 min and +43 57, so we can check the calculation to moderate accuracy. Other variable stars in the Burnham's list for Lyra are not shown on Tirion's map, so we might need to plot them from our converted coordinates. The basic steps in the algorithm are;
ra = (hrs + min/60) * 15
dec = degrees + minutes/60,
dec is negative if 's' is in the position as given in
Burnham's and positive if 'n'.
For our example;
Burnham's gives 18538n4353 as the 1950 coordinates, so
ra = (18 + 53.8/60) * 15 = 18.896667 *15 = 283.45 degrees
dec = 43 + 53/60 = 43.883333 degrees
The system of rectangular coordinates used has;
These are the formulas for converting to rectangular coordinates;
X = cos(ra) * cos(dec)
Y = sin(ra) * cos(dec)
Z = sin(dec)
You can derive these formulas using nothing more than simple trigonometry applied to two right angled triangles. An explanation is available.
Applying these formulas to our example we have;
X = cos(283.45) * cos(43.883333) = 0.1676447
Y = sin(283.45) * cos(43.883333) = -0.7009849
Z = sin(43.883333) = 0.6931922
To correct for the effects of precession, we need to rotate the coordinate axes. The star has remained in the same place, except for any proper motion. The rotations are applied as follows;
Rotations can be modelled using matrix algebra, each rotation is represented by a 3 by 3 matrix, and the three matrices are multiplied to calculate a final 'precession matrix'. Duffett-Smith has a good explanation of this method. Meeus approaches the whole task using formulas involving the coordinates. The point about precession is that the 'precession matrix' depends only on the difference in time between the two epochs.
Below are the entries for the precession matrix for conversion from B1950.0 to J2000.0;
NOTE: 1999 July 15. Eric Flesch has pointed out that the
middle column bottom row value should be -0.000027 not
0.000027. This will change subsequent results. Correction
needs to be made to the QBASIC code as well.
-- --
| 0.999926 -0.011179 -0.004859 |
| 0.011179 0.999938 -0.000027 |
| 0.004859 0.000027 0.999988 |
-- --
You can multiply the column vector representing the B1950 components by this precession matrix to generate the column vector of J2000 components. As neither QBASIC or my TI-80 calculator have matrix operations, below are the formulas for each component written out longhand;
x' = 0.999926*x -0.011179*y -0.004859*z y' = 0.011179*x +0.999938*y -0.000027*z z' = 0.004859*x +0.000027*y +0.999988*z
Applying this formula to our example;
The B1950 coordinates are; x = 0.1676447 y = -0.7009849 z = 0.6931922 Remember to do the multiplications before adding the terms, and that positive and negative numbers multiply to give negative numbers! x' = (0.999926 * 0.1676447) - (0.011179 * -0.7009849) - (0.004859 * 0.6931922) y' = (0.011179 * 0.1676447) + (0.999938 * -0.7009849) - (0.000027 * 0.6931922) z' = (0.004859 * 0.1676447) + (0.000027 * -0.7009849) + (0.999988 * 0.6931922) The result of each product is shown below x' = 0.1676323 + 0.0078363 - 0.0033682 = 0.1721004 y' = 0.0018741 - 0.7009414 - 0.0000187 = -0.6990860 z' = 0.0008146 - 0.0000189 + 0.6931839 = 0.6939796 and we have the J2000 components; x' = 0.1721004 y' = -0.6990860 z' = 0.6939796
The new values are close to the old values, precession is a small effect, so we can be reasonably confident that the calculations are correct. An error with a sign usually makes a large impact in the result!
You convert the rectangular components back to angular measure using the following formulas;
r = arctan(y/x)
if x < 0 then ra2000 = r + 180
if y < 0 and x > 0 then ra2000 = r + 360
dec2000 = arcsin(z)
Applying these formulas to our example data;
x' = 0.1721004
y' = -0.6990860
z' = 0.6939796
r = atan(-0.6990860/0.1721004)
= atan(-4.0620824)
= -76.169982 degrees
y' < 0 and x > 0, so add 360;
ra2000 = 283.83002 degs
= 18.922001 hrs
dec2000 = arcsin(0.6939796)
= 43.94596 degs
We can write this in various forms, including the notation used in
Burnham's and in degrees, minutes and seconds form.
Notations are important, and software is notorious for giving the wrong notation when you need results quickly! To convert the J2000.0 (or any) RA and DEC to degrees minutes and seconds form, use the following recipe;
RA = 18.922001 hrs
subtract the whole number of hours;
18.922001 - 18 = 0.922001 hrs
multiply the fractional part by 60;
0.922001 * 60 = 55.32006
the number of minutes is 55. Subtract this from what remains;
55.32006 - 55 = 0.32006
multiply by 60 again;
0.32006 * 60 = 19.20 seconds
The RA is thus, 18;55;19.2
The same method on the degrees gives;
43.94596 degs = 43;56;45
RA is meaured in hours of time, each of which is equal to 15 degrees. For this reason, RA is usually given to slightly more accuracy than DEC.
Many star atlases and observing guides list the positions to 0.1 minute of time for RA and 1 arcminute for DEC. Burnham's uses a very compact notation for this level of accuracy as explained earlier. To convert our example results into Burnham's notation, use the following recipe;
RA = 18.922001 hrs
DEC = 43.94596 degs
subtract the hours from the RA, and multiply remainder by 600
to get the tenths of minutes;
18.922001 - 18 = 0.922001
0.922001 * 600 = 553.2006
round number to the nearest whole unit, and join the hrs
figure to the minutes figure;
18553
Now subtract the whole degrees from the DEC figure and multiply
the remainder by 60, and round the number to the nearest whole figure;
43.94596 - 43 = 0.94596
0.94596 * 60 = 56.7576 = 57
join the degrees and arcminute figures together;
4357
Now join both figures together, but with the letter 'n' in between
for positive DEC, and 's' for negative dec;
18553n4357
It is now easy to see that the J2000.0 position for R-lyr coincides with the position given in the The Cambridge Star Atlas.
At the end of this section, you will find a full QBASIC
program which will prompt the user for a position in B1950.0
coordinates using the notation in Burnham's Celestial
Handbook. The actual calculation is done using the lines below;
' Corrected precession coefficients. ' form direction cosines (i.e. rectangular coordinates assuming ' radius vector of 1 unit) x = COS(ra) * COS(dec) y = SIN(ra) * COS(dec) z = SIN(dec) ' apply the precession matrix x2 = .999925708# * x - .0111789372# * y - .0048590035# * z y2 = .0111789372# * x + .9999375134# * y - .0000271626# * z z2 = .0048590036# * x - .0000271579# * y + .9999881946# * z ' convert the new direction cosines to RA, DEC ra2 = FNatn2(y2, x2) dec2 = FNasin(z2)
This should look similar to the formulas presented in the last section. Some programmable calculators (and spreadsheets such as Microsoft Excel and Lotus 1-2-3) have built-in matrix multiplication routines, so a more compact code may be used.
The rest of the code is concerned with defining some useful
functions, and processing input and output strings. I have used the
multiline version of the DEF FN statement to define the
functions. By doing this the code will compile using the FirstBasic
shareware compiler available from
PowerBasic.
The easiest way to transfer the listing into QBASIC would be to;
*******
All files (*.*), and
using the file extension .BAS
QBASIC or the FirstBas basic
compiler.
'*********************************************************************
CLS
' Make all variables double precision (about 16 significant figs)
DEFDBL A-Z
' Define some useful constants in radians
pi = 4 * ATN(1)
tpi = 2 * pi
degs = 180 / pi
rads = pi / 180
' Returns arctan in correct quadrant
DEF FNatn2 (y, x)
a = ATN(y / x)
IF x < 0 THEN a = a + pi
IF (y < 0) AND (x > 0) THEN a = a + tpi
FNatn2 = a
END DEF
' Returns arcsine, which QBASIC does not have
DEF FNasin (x)
c = SQR(1 - x * x)
FNasin = ATN(x / c)
END DEF
' returns a string which represents the RA as
' hhmmm, where the last figure is 1/10 minute of time
DEF FNra$ (r)
a = INT(r)
b = 600 * (r - a)
r$ = STR$(INT(a * 1000 + b + .5))
r$ = RIGHT$(r$, LEN(r$) - 1)
WHILE LEN(r$) < 5
r$ = "0" + r$
WEND
FNra$ = r$
END DEF
' returns a string which represents the DEC as ddmm
DEF FNdec$ (d)
IF d < 0 THEN
ds$ = "s"
ELSE
ds$ = "n"
END IF
c = ABS(d)
a = INT(c)
b = 60 * (c - a)
d$ = STR$(INT(a * 100 + b + .5))
d$ = RIGHT$(d$, LEN(d$) - 1)
WHILE LEN(d$) < 4
d$ = "0" + d$
WEND
FNdec$ = ds$ + d$
END DEF
' Start of main program
PRINT
PRINT "Precess B1950 to J2000"
PRINT "----------------------"
PRINT
PRINT "type 'q' to quit"
PRINT
10 PRINT 'start of loop
' Get input string from user
INPUT "B1950 : ", inpos$
' some simple error checking
IF inpos$ = "q" THEN 20
IF LEN(inpos$) < 10 THEN
GOTO 10
END IF
' get the numerical values of RA and DEC
ra = VAL(LEFT$(inpos$, 5)) / 1000
dec = VAL(RIGHT$(inpos$, 4)) / 100
' get sign of the DEC
SELECT CASE MID$(inpos$, 6, 1)
CASE "s"
decsgn = -1
CASE "n"
decsgn = 1
CASE ELSE 'detects wrong input format
GOTO 10
END SELECT
' extract decimal hrs/degs values from input format
rahrs = INT(ra)
ramin = (ra - rahrs) * 1000
decdeg = INT(dec)
decmin = (dec - decdeg) * 100
ra = rahrs + ramin / 600
dec = decsgn * (decdeg + decmin / 60)
' convert to radians
ra = ra * rads * 15
dec = dec * rads
' form direction cosines (i.e. rectangular coordinates assuming
' radius vector of 1 unit)
x = COS(ra) * COS(dec)
y = SIN(ra) * COS(dec)
z = SIN(dec)
' apply the precession matrix
x2 = .999925708# * x - .0111789372# * y - .0048590035# * z
y2 = .0111789372# * x + .9999375134# * y - .0000271626# * z
z2 = .0048590036# * x - .0000271579# * y + .9999881946# * z
' convert the new direction cosines to RA, DEC
ra2 = FNatn2(y2, x2)
dec2 = FNasin(z2)
' call functions to generate the Burnham notation and print
PRINT "J2000 : "; FNra$(ra2 * degs / 15) + FNdec$(dec2 * degs)
' loop round for next input
GOTO 10
PRINT
20 PRINT : PRINT "bye"
END
'********************************************************************
The output from a typical session is shown below, I added the comments
later! The position 01487n8902 means 1h 48.7 min RA and
+89 degs 2 min DEC, and is the 1950 position of the Pole Star!
B1950 : 01487n8902 'pole star J2000 : 02318n8916 B1950 : 02470n5541 'eta perseus J2000 : 02507n5553 B1950 : 04330n1625 'alpha taurus J2000 : 04359n1631 B1950 : 10057n1212 'alpha leo J2000 : 10084n1157 B1950 : 15573s2229 'delta scorpius J2000 : 16003s2238 B1950 : 01359s5730 'alpha eridanus J2000 : 01377s5715
The 1950 positions are apparent positions, including nutation, taken from Chris Mariott's Skymap v2.1 star chart program. I set the program to the date 1949 December 31st 2209 UT, which corresponds to the instant of B1950. The positions produced by this program coincide with the J2000.0 positions of the stars very closely. There is the odd difference here and there, almost certainly due to the nutation correction.
Duffett-Smith, Peter
Practical Astronomy with your calculator
Cambridge University Press
3rd edition 1988
ISBN 0-521-35699-7
Straight forward methods, easy to understand, helped me through my University astronomy module! Cheap and recommended if you are working with a programmable calculator or fancy writing some simple Basic programs.
Meeus, Jean
Astronomical Algorithms
Willmann-Bell
1st English edition, 1991
ISBN 0-943396-35-2
Comprehensive cookbook of astronomical calculation. Presents formulas and algorithms - you have to make your own program code, although early chapters do explain some common pitfalls and techniques. Not much in the way of derivation or linking of methods to dynamics. Each copy has a serial number, and you can buy a disc supplement with program code in Pascal, QuickBasic or C.
Seidelmann, P. Kenneth (ed)
Explanatory Supplement to the Astronomical Almanac
University Science Books
1992, completely revised
ISBN 0-935702-68-7
Definitive reference on all aspects of the ephemeris and
associated calculations. No hostages taken, no example calculations
and modern vectoral notation used throughout. Very thorough treatment
of high precision precession calculations in section 3.21.
Burnham, Robert Jr
The original hitchiker's guide to the Galaxy: detailed descriptions
of objects observable in amateur telescopes, and references to
mythical stories and historical observations.
Hoffleit E.D., Warren Jr. W.H.
Machine readable version is available for download from the NASA
Astronomical Data Center. Occupies about 2 Mb of disc space
in text form (800K download), and provides basic data on the 'naked
eye stars'.
Tirion, Wil
Nice clear maps based on the Bright Star Catalog. Covers the sky
in 20 maps, each taking one complete page. Opposite page lists details
of notable objects including their J2000 coordinates. Cheaper than
Norton's Star Atas in the UK, and the maps don't get
spread over the page fold.
[Root]
Last Modified 15th July 1999
Burnham's Celestial Handbook: An observer's guide
to the universe beyond the solar system
Dover Publications Inc
1977, revised and expanded edition, three volumes
ISBN 0-486-23567-X (Vol 1, paper)
ISBN 0-486-23568-8 (Vol 2, paper)
ISBN 0-486-23673-0 (Vol 3, paper)
The Bright Star Catalog
Yale University, National Space Science Data Center
1991, 5th Revised Ed.
The Cambridge Star Atlas
Cambridge University Press
2nd edition 1996
ISBN 0-521-56098-5
Keith Burnett