Do you have your own computer? And do you wish to use it for predicting astronomical phenomena such as the time of moonrise or the circumstances of the next eclipse? If so, this book may help you. It provides a kit of parts, a set of building blocks, to write complex programs in astronomy tailored to your own needs. Each building block is in the form of a subroutine, written in BASIC, which performs a well-defined function. There are 33 such routines in all which can link together in any combination. They cover most of the problems likely to be tackled by the amateur astronomer, ranging from very simple tasks like converting between decimal hours and hours, minutes, and seconds, to complex procedures for deducing the set of orbital elements which is consistent with the observed movement of a minor planet or comet through the heavens. It is possible that your particular requirement may be covered directly by a subroutine as it stands, for example, the calculation of the position of the Moon. You can then enter it from the keyboard or copy it from disk+ and run it with its supplied handling program. Or it may be that you have a much more demanding task, needing the results of several routines together to get the answer. Suppose, for instance, that you wished to predict the time of rising of a newly-discovered comet in Hawaii on the 27th August this year.

Routine RISET (3100) calculates the local sidereal times and azimuths of a celestial body given its right ascension and declination as input parameters. We may combine this routine with SUN (3400) and EQECL (2000) to calculate the right ascension and declination of the Sun on any given date, and hence the circumstances of sunrise and sunset. The procedure is complicated, however, by the fact that the Sun is not stationary in the sky with respect to the background of stars but moves by about 1 degree of ecliptic longitude per day. The calculated circumstances of rising and setting would be correct for an object at the Sun's position at the given instant, but would not correspond exactly with the observed circumstances of sunrise or sunset unless we had chosen the correct time already. We therefore have to adopt an iterative procedure, refining an initial crude estimate to achieve the required accuracy.

Routine SUNRS first calculates the Sun's right ascension and declination at local midday of the date in question, and hence the local sidereal times of rising and setting corresponding to the Sun's position at midday. The times are then corrected for the observer's longitude to find the Greenwich sidereal times, and converted to universal times using routine TIME. These times are already fair approximations to the actual times of rising and setting.

We saw in an earlier section how to calculate the position of any of the major planets using routine PLANS (4500) and its handling program. That routine itself called PELMENT (3800) to obtain the mean elliptical orbital elements for each of the planets, and then it applied lots of small corrections for the perturbing influences of other members of the Solar System. It was able to do this for each of the planets (except for Pluto) because dynamic theories of their motions have been worked out. It was therefore just a matter of plugging into the given formulae.

Most members of the Solar System, such as the periodic comets, minor planets, and Pluto, have not had dynamical theories worked out for their motions. Their orbital elements can be calculated on the basis of their observed orbits at any instant, but those elements change with time and will not be appropriate for a later instant. Elements which hold only over a short period near a given date are called osculating elements. If you use them to calculate positions within the period for which they are valid, you will get accurate results. The further you move outside that period, however, the more inaccurate the results become.

Osculating elements are given for the minor planets (and indeed for the major ones as well) in The Astronomical Almanac, correct for a specified instant. Having obtained them, you can use routine ELOSC and its handling program to calculate orbital positions. ELOSC is very similar to PLANS in many ways.

A point in the sky may most easily be fixed by an observer on the Earth with reference to his horizon. In the horizon coordinate system (see Figure 2), the position of the point is specified by its azimuth, the angle round from the northf point of the horizon (in the sense NESW) and by its altitude, the angle up from the horizon (positive if above the horizon, negative if below it). The positions of heavenly bodies, on the other hand, are very often described in the equatorial coordinate system (see Figure 3). Here the plane of the Earth's equator, extended to cut the celestial sphere, is used instead of the horizon, with the first point of Aires or vernal equinox taking the place of the north point of the horizon. A star's position is then given by the angle round from the vernal equinox along the equator (in the opposite sense to that in which it appears to move throughout the day) and the angle up from the equator (positive if to the north, negative if to the south). These coordinates are called the right ascension and declination respectively. Related to the right ascension is the hour angle which describes the angle along the equator from the observer's meridian. As the Earth rotates, so the star and the vernal equinox appear to move at the same rate, making the right ascension and the declination constants (but see routines PRCESS1 and PRCESS2). The hour angle, however, increases uniformly throughout the day from zero when the star crosses the meridian (moving westerly).

The microchip revolution has made mathematical positional astronomy easy. Computing power, such as would have been the envy of even the most prestigious institutions several decades ago, is now commonplace in the home. Amateur astronomers can perform complex calculations with comparative ease, making predictions of forthcoming celestial events or demonstrating what took place in the past, with a precision that is almost unlimited.

Many books have already been published which list the steps (more or less clearly) to be followed to the desired answer.Some are for calculators, but others exploit the special powers of the microcomputer. Unfortunately, those in the former category do not provide programs which can be used immediately on the computer, while those in the latter category are often too specifically biased towards particular makes of machine and particular calculations. While it is always possible tomodify a program to work on another model or to make a different calculation, it is not always easy to do so, especially when detailed notes about its workings are not supplied.

This book is designed for the person who wishes to use his or her personal computer to make astronomical calculations with the minimum of fuss. It is not specific to any make of machine, neither are you confined to specific calculations. Rather, you are presented with a set of linking subroutines, each carrying out a specific task, which you can combine to form a complex program designed to your own requirements. You need only have a broad understanding of the problem; the subroutines themselves take care of the details.

The positions of members of our Solar System may best be specified in the ecliptic coordinate system (see Figure 4) since this uses the plane of the ecliptic as the fundamental plane and most members of the Solar System move in orbits close to it. The reference direction is that of the vernal equinox, as in the equatorial system. This is an obvious choice because it lies along the line of intersection of the planes of the Earth's equator and the ecliptic. The ecliptic longitude is the angle measured in the ecliptic round from the vernal equinox in the same sense as that in which the right ascension is measured, and the ecliptic latitude is the angle up from the ecliptic, positive if to the north, or ‘above’ the ecliptic, and negative if to the south, or ‘below’ the ecliptic.

The subroutine EQECL converts coordinates between the ecliptic and equatorial systems. The angles to be converted are input via the variables X (ecliptic longitude or right ascension) and Y (ecliptic latitude or declination), with the results output via P (right ascension or ecliptic longitude) and Q (declination or ecliptic latitude). The formulae are not quite symmetrical in the two pairs of coordinates so that the direction of conversion needs to be specified with a switch.

It is sometimes necessary to calculate the date and time of new or full moon corresponding to a particular calendar month. For example, since a solar eclipse can only occur at new moon, and a lunar eclipse at full moon, we need those times to predict the occurrence of an eclipse. This information could be found using the routine MOON (6000) and its handling program, but it would be a tedious exercise as the phase would have to be calculated several times in an iterative procedure of successive approximation. Routine MOONNF adopts the algorithm given by Meeus in his excellent book Astronomical Formulae for Calculators (second edition) to calculate the instants of new and full moon to within an accuracy of about 2 minutes.

The date is input as usual via the variables DY, MN, and YR. The instants of new and full moon nearest to the date are calculated in terms of the number of Julian days since 1900 January 0.5, being returned by the pairs of variables FI, FF for full moon and NI, NF for new moon. The first variable in each pair indicates the integer number of Julian days of the event since 1900 January 0.5, while the second variable indicates the fraction of the day. The Julian date of new moon is therefore

JD = 2415020 + NI + NF

since 1900 January 0.5 corresponds to Julian date 2415 020.0. Note that Julian days begin at noon, a fact which must be considered when converting to calendar date and time.

DEFAULT

The handling programs in this book, whose names all begin with the letter H (e.g. HEQECL, HALTAZ, etc.) require the user to supply answers to questions such as ‘Today's date?’ or ‘Geographical latitude?’. Often it is convenient to make use of default settings for the answers, especially if they are many characters in length. Routine DEFAULT allows you to do this very simply. To use it, you supply the question in the string variable Q$, the number of variables to be entered in N (= 1, 2, or 3) and the default settings for those variables in X, Y, and Z. The routine then displays the question on the screen together with the default settings in square brackets [], and waits for a response. You may enter new numbers for the variables as usual if you wish. Alternatively, you can enter a comma (,) in place of every variable, and the routine uses the default values instead. The new values, or default settings, are returned in the same variables X, Y, and Z.

Routine DEFAULT is used later in the book in the handling programs HELOSC, HPCOMET, HPFIT, and HEFIT. In HELOSC, for example, lines 160 and 165 get the date of the epoch (D,M,Y) by setting Q$=“Epoch (D,M,Y)”, X, Y, and Z to the default values of the day, month, and year respectively, and N=3 (since there are three variables to assign). Suppose that X, Y, and Z were given the default values 12, 3, 1986 respectively.

The result of an astronomical calculation is sometimes an instant of time, expressed as a number of days since a fundamental epoch. In our convention, that epoch is 1900 January 0.5, and the number of days elapsed since then is represented by the real variable DJ. For example, we may have calculated the time of the next lunar eclipse as DJ = 30 815.69, but unless we are familiar with the interpretation of such numbers, the result will not mean very much to us. We need a conversion routine to convert a Julian day number into a calendar date, the inverse function of JULDAY.

Subroutine CALDAY performs this function. It converts DJ, the number of Julian days (and fraction of a day) elapsed since 1900 January 0.5, into the calendar date represented by DY, the day of the month including the fraction of the day, MN, the month number, and YR, the year number. The fraction of the day is also returned by FD. Thus DJ = 30815.69 gives DY = 15.19, MN = 5, YR = 1984 and FD = 0.19, from which we conclude that the lunar eclipse was in progress on 1984 15th May at 4 : 33 UT (Universal Time). The integer part of the day number is returned in the variable ID, and a string representation of the entire date in the string variable DT$.

The Sun's celestial coordinates are determined by its position in its apparent orbit of the Earth. Although the Earth actually orbits the Sun, it seems to us that it is the Sun which is in orbit about the Earth, and for the purposes of calculating the Sun's position it is convenient to regard this as the case. We can then compute the ecliptic longitude of the Sun for a given instant and, since the Sun-Earth orbit defines the plane of the ecliptic, we can take the ecliptic latitude to be 0. The right ascension and declination can be found by using routine EQECL (2000) followed by EQHOR (1500) if we wish to know the Sun's azimuth and altitude as well.

Subroutine SUN calculates the true geocentric ecliptic longitude of the Sun, SR, and its distance from the Earth, RR, at the instant given by UT (Universal Time, hours), DY (integer days), MN (months), and YR (years). The time may also be expressed as a fraction of a day and added to the value of DY, in which case UT must be set to 0. The true ecliptic latitude may be assumed to be 0 (always correct to within 1.2 arcseconds). SR is the longitude referred to the mean equinox of date and is the quantity required in the calculation of geocentric planetary positions.

As the Earth moves around the Sun, it sweeps out the plane of the ecliptic,making one complete circuit in one year. The Earth also rotates about its ownNorth-South polar axis, making one revolution in 24 hours. Since the Earth acts like a huge gyroscope, the direction of its polar axis remains relativelyfixed although precession (PRCESS1, 2500; PRCESS2, 2600) and nutation (NUTAT, 1800) conspire to cause a wobbling motion. This direction is inclined at about 23.5 degrees to the perpendicular to the plane of the ecliptic, and hence the Earth's equator is tilted at the same angle to the plane of the ecliptic. The angle is called the obliquity of the ecliptic and it is needed whenever we convert between ecliptic and equatorial coordinates (see EQECL, 2000). The angle changes slowly with time so it must be determined separately for calculations made at different epochs.

The subroutine OBLIQ returns the value of the obliquity of the ecliptic via the variable OB (degrees). The dateis input in the usual way via DY (days), MN (months), and YR (years). Execution is controlled by the flag FL(5), designed to save time if the calculation has already been made for the current date. When calling the routine for the first time, or for a new date, set FL(5) = 0. It is set to 1 by the routine (line 1705). Subsequent calls to OBLIQ with FL(5) =1 result in no calculation being performed, but control is returned immediately to the calling program.

Observations from the Earth's surface have to be made through the atmosphere. Although the air is tenuous, it has a refractive index which is significantly greater than 1, the value for a vacuum. When we observe a celestial object at an angle from the vertical the rays are bent slightly away from the path they would have followed if the atmosphere had not been present. The atmosphere clings to the curve of the Earth and it acts rather like a weak lens, making stars appear to be nearer to the vertical than they really are. The largest effect is seen at rising and setting where the horizontal refraction is about 34 minutes of arc. Since the effect is always to increase the Sun's apparent altitude, the length of the day is increased because of it.

Subroutine REFRACT assumes a model atmosphere whose refractive properties depend only on pressure, PR</b.>, and temperature, TR, at the place of observation. The amount of refraction actually observed depends on the detailed disposition of the atmosphere along the line of sight to the celestial object in question, and significant deviations from the values calculated here must be expected in some instances, especially near the horizon. Nevertheless, REFRACT will provide a very good approximation which will be more than adequate in most circumstances.

The altitude of the object is input via the parameter Y, corrected by the routine, and returned via the parameter Q. The direction of correction is indicated by the switch SW(6).

During the course of one day, the stars and other celestial objects appear to move in circles about a single point in the sky. This point is the north celestial pole for observers in the northern hemisphere or the south celestial pole in the southern hemisphere, and it marks the intersection of the Earth's spin axis with the celestial sphere. Stars which are close enough to the pole remain above the horizon all day and are called circumpolar stars. Stars which are far enough away from the pole spend part of the day below the horizon. As they cross the horizon on the way down they set, and as they cross on the way up they rise. Stars further still from the pole spend all day below the horizon and are never seen in that location.

Routine RISET finds the circumstances of rising and setting for a celestial object whose right ascension and declination are known. The local sidereal times are calculated for a given geographical latitude, GP, and are returned via the variables LU (rising) and LD (setting). The azimuths are also given by AU and AD. the local sidereal times may be converted to Greenwich sidereal times, Universal Times, and local civil times using routine TIME (1300).

The coordinates input to the routine by X and Y are usually the true geocentric coordinates.