Legacy Documentation

Comets

http://www.wolfram.com/applications/astronomer/index.html

setup

<<Astronomer`HomeSite`

SetLocation[GeoLongitude145. Degree,
GeoLatitude-37.8 Degree,
TimeZone11];

Comets

This notebook deals with how to add comets and then analyze the circumstances of their apparition.

By far the most famous comet is Comet Halley. For the past few thousand years its every appearance has been recorded in history.

Every now and then, other great comets can appear. For example, the next great hope for a good naked-eye visible comet is comet Hale-Bopp, first discovered in mid-1995. Hopefully it should be bright when it comes closer in March, 1997.

Comet Halley: Orbital Elements and Orbit

Comet Halley orbits the Sun roughly every 76 years. Its closest approach to the Sun is about 0.587 AU, and because its orbital eccentricity is 0.967 the semi-major axis of its elliptical orbit is 0.587/(1-0.967) AU, or about 17.9 AU. Thus its farthest distance from the Sun is about twice that, at 35.3 AU, which is just past the orbit of Neptune.

Those numbers above are just some of the orbital elements of comet Halley. If you know the full orbital elements then you can add the object to the package.

Seven numbers are needed to define the Orbital Elements of an object.

For example the Orbital Elements of comet Halley are as follows. First it orbits the Sun in an elliptical orbit with an eccentricity of 0.967, a semi-major axis of 17.9 AU and a period of about 76 years as mentioned above. The orbit itself is inclined 162 degrees to the plane of the Earth orbit, and the orbit cuts that plane upwards at an ascending longitude of 58.8 degrees. The perigee (or point of closest approach to the Sun) occurs 111.8 degrees further along the orbit.

The final piece of information required to define the orbital elements is the position of the object at some date. In the case of comet Halley, it is known that a perigee occurs on 1986 February 9. At perigee the so-called mean anomaly of an object is 0 degrees. (In this context an anomaly is the angle of an object away from the perigee point. And the mean anomaly is the angle of a fictional object that moves uniformly around the ellipse.)

Putting all these pieces together you can set the orbital elements of Comet Halley as follows.

SetOrbitalElements[Halley,
ViewPoint -> Sun,
Date -> {1986, 2, 9.4589},
TimeZone -> 0,
OrbitalSemiMajorAxis ->
(0.587104*AU)/(1 - 0.967277),
OrbitalEccentricity -> 0.967277,
OrbitalInclination -> 162.2422*Degree,
OrbitalPeriod ->
(0.587104/(1 - 0.967277))^(3/2)*
365.25*Day,
MeanAnomaly -> 0*Degree,
PerigeeArgument -> 111.8657*Degree,
AscendingLongitude -> 58.8601*Degree];

You can now use Halley throughout the package in any context where you could use a planet, asteroid, moon or star etc.

In particular you could use it in the OrbitPlot function to show how its orbit compares with the outer planets.

OrbitPlot[{Jupiter,Saturn,Uranus,
Neptune,Pluto,Halley},
Distance48 AU];

The above shows that comet Halley is in a highly elliptic orbit that passes just outside the orbit of Neptune. The direction of the major axis, as measured counter-clockwise from the x-axis is about 135 degrees or equivalently 9 hours of right ascension which is in the constellation of Cancer.

One thing that is not clear from the above is that Halley orbits in a clockwise sense, unlike all the planets which orbit anticlockwise.

Hale-Bopp: Orbital Elements and Magnitude

Comet Hale-Bopp is currently in a 5,000 year orbit around the Sun. Its closest approach is 0.918 AU (just inside the Earth's orbit), but it travels out to about 300 AU or about 10 times further out than Neptune.

The great interest in comet Hale-Bopp is that it will probably be a bright naked eye visible object at closest approach at the end of March, 1997.

Current orbital elements for Hale-Bopp are:

SetOrbitalElements[HaleBopp,
ViewPoint -> Sun,
Date -> {1997, 3, 31.8154},
TimeZone -> 0,
OrbitalSemiMajorAxis ->
(0.918046*AU)/(1 - 0.997071),
OrbitalEccentricity -> 0.997071,
OrbitalInclination -> 88.9858*Degree,
OrbitalPeriod ->
(0.918046/(1 - 0.997071))^(3/2)*
365.25*Day,
MeanAnomaly -> 0*Degree,
PerigeeArgument -> 130.3508*Degree,
AscendingLongitude ->
282.4674*Degree];

Magnitude estimates

The estimated apparent magnitude for a comet like Hale-Bopp is given by the following formula.

0.0 + 5.0*Log[10,earthDistance] + 7.5*Log[10,solarDistance]

The leading 0.0 is related to the intrinsic magnitude of the comet, and the second term involving 5.0 is simply how the magnitude of any object will change as a function of its distance from the Earth viewing location. The final term involving 7.5 reflects how the comet might brighten as it gets closer to the Sun.

Using the formula above you can compute Hale-Bopp's apparent magnitude during the period say April, 1996 to April, 1998 as follows:

mag = Table[{m, 0. + 5.*
Log[10, Distance/AU /.
EquatorCoordinates[HaleBopp,
{1996, 0, 30.5*m}]] +
7.5*Log[10, Distance/AU /.
EquatorCoordinates[HaleBopp,
{1996, 0, 30.5*m},
ViewPoint -> Sun]]
}, {m, 4, 28, 0.5}];

And then plot it as follows:

ListPlot[mag,
PlotJoinedTrue,
AxesOrigin{13,0},
PlotRange{{4,27},{0,8}},
Ticks{{{ 4,"A"},{ 5,"M"},{ 6,"J"},
{ 7,"J"},{ 8,"A"},{ 9,"S"},
{10,"O"},{11,"N"},{12,"D"},
{13,"J"},{14,"F"},{15,"M"},
{16,"A"},{17,"M"},{18,"J"},
{19,"J"},{20,"A"},{21,"S"},
{22,"O"},{23,"N"},{24,"D"},
{25,"J"},{26,"F"},{27,"M"}},
Automatic}];1;

The above shows that Hale-Bopp should be brightest around the end of March, 1997, when it might reach magnitude 0.4. Between about January, 1997 and July, 1997 the comet should be brighter than magnitude 4.

Similar to the above, you can plot the expected magnitude for the critical period at the end of March, 1997.

mag1 = Table[{d, 0. + 5.*
Log[10, Distance/AU /.
EquatorCoordinates[HaleBopp,
{1997, 3, d}]] +
7.5*Log[10, Distance/AU /.
EquatorCoordinates[HaleBopp,
{1997, 3, d},
ViewPoint -> Sun]]},
{d, 16, 31, 1}];

ListPlot[mag1,PlotJoinedTrue];1;

In the above the brightnest magnitude occurs around 1997 March 28.

Hale-Bopp: Trajectory April 1996 to April 1998

Blue numbers, marked along the red trajectory of the comet, represent the number of months after the beginning of year 1996 in the following graphic. The corresponding position of the Sun is marked in orange along the ecliptic line.

StarChart[
MagnitudeScale0.7,
MagnitudeRange4.0,
Epilog{RGBColor[1,0,0],
OrbitTrack[HaleBopp,
{1996,4,1},{1998,4,1},
PlotPoints200],
RGBColor[0,0,1],
Thickness[0],
Table[OrbitMark[HaleBopp,
{1996,m,1},
PlotLabelm],
{m,6,26}],
RGBColor[1,0.5,0],
Thickness[0],
Table[OrbitMark[Sun,
{1996,m,1},
PlotLabelm],
{m,14,20}]}];1;

Because the comet is in a highly elliptic orbit with a long period, the start and stop points of the trajectory are roughly fixed points in the sky.

From about months 12 to 13 (December, 1996 to January, 1997) the comet will brighten rapidly, and will be visible predominantly to observers in the Northern hemisphere.

The comet moves fastest through the sky from about months 15 to 17 (from March to May, 1997), when it will be closest to the Earth. It will reach its brightest around 1997 March 26 through to 1997 April 12.

By the end of month 16 (late April, 1997) the comet will move closer in the sky to the Sun, and will also begin to fade.

By month 18 (June, 1997) it will have passed the Sun in the sky and will be visible only in the morning sky. Northern hemisphere observers will eventually lose the comet as it rapidly moves south and continues to fade.

After month 20 (October, 1997) it will be faint and only visible to Southern hemisphere observers with binoculars.

Hale-Bopp: Trajectory March 1997 to May 1997

Comet Hale-Bopp will be brightest in April, 1997. From about March to May it should be visible in the Northern hemisphere during evening hours just after dusk and roughly in the direction of the constellation of Andromeda.

The best dates should be about March 26 until April 12 when it will be visible in the north-west sky just after dusk and about 20 degrees above the horizon for most Northern hemisphere viewers.

setup annotate

RadialStarChart[Andromeda,
RadialAngle50 Degree,
MagnitudeScale0.4,
MagnitudeRange4.0,
ConstellationLabelsTrue,
Epilog{RGBColor[1,0,0],
OrbitTrack[HaleBopp,
{1997,3,1},{1997,5,1},
PlotPoints200],
annotate[]}];1;

The numbers along the trajectory are the day of the month in March and then April. The blue lines away from the red comet trajectory indicate the rough direction of the comet tail pointing away from the Sun.

The corresponding position of the Sun is shown as orange dots near the bottom of the graphic.

Hale-Bopp: 3D Stereographic Animation

To help visualise the motion of the comet you can use the Coordinates function to help plot the relative position of the Sun, Earth and comet in three-dimensional space. To make visualisation easier you can produce a pair of three-dimensional plots with slightly different viewpoints, so as to make a stereographic pair.

The following cell contains the setup for all that.

setup haleStereo

dateLabel[y_, m_] :=
Block[{day, mon, yea}, day = If[Abs[Round[m] - m] < 0.2, " 1 ", "16 "];
mon = {"Jan", "Feb", "Mar", "Apr", "May", "Jun", "Jul", "Aug", "Sep",
"Oct", "Nov", "Dec"}[[Mod[Floor[m] - 1, 12] + 1]];
yea = ToString[y + Floor[(m - 1)/12]]; StringJoin[day, mon, "\n ", yea]]

haleStereo[y_]:=Do[Block[{cs,ce,cc,d,p1,p2,sun,ear,com,text},d={y,1,m 30.5};cs={0,0,0};ce=Coordinates[Earth,d,ViewPointSun];cc=Coordinates[HaleBopp,d,ViewPointSun];sun={RGBColor[1,1,0],PointSize[0.05],Point[cs]};ear={RGBColor[0,0,1],PointSize[0.03],Point[ce]};com={RGBColor[1,0,0],PointSize[0.02],Point[cc]};text={Text[dateLabel[y,m],{-1.5,-3,0.0}]};p1=Graphics3D[{sun,ear,com,text},ViewPoint{1.3,-2.4,2} 2,SphericalRegionTrue,PlotRange{{-1,1},{-2,1},{-2,1}} 2];p2=Graphics3D[{sun,ear,com,text},ViewPoint{1.5,-2.4,2} 2,SphericalRegionTrue,PlotRange{{-1,1},{-2,1},{-2,1}} 2];Show[GraphicsArray[{p1,p2}],GraphicsSpacing0];],{m,4.0,27.99,0.5}]

If you can successfully cross your eyes, then you should be able to see the following animation in full 3D.

haleStereo[1996]

The fixed yellow dot represents the position of the Sun, and the blue dot rotating around it represents the Earth. The red dot that appears later is the comet Hale-Bopp.

If you can successfully cross your eyes to make the above two images merge into a full 3D image, then when you animate the cells you should see the comet sweep over the top of the Sun and then back and down behind it.

You can also get a feel for the relative positions of the Sun, Earth and comet during the interesting period of December, 1996 to April, 1997. Use the arrow keys on your front end to step forward and back through those frames.

Hale-Bopp: Rise/Set Charts

This section computes some rise/set charts for Hale-Bopp from the beginning of 1996 to the end of 1999.

Compute the charts

no = Table[RiseSetChart[HaleBopp, year,
GeoLatitude -> 37*Degree,
AspectRatio -> 1,
Mesh -> True,
Text -> False,
DisplayFunction -> Identity],
{year, 1996,1999}];

so = Table[RiseSetChart[HaleBopp, year,
GeoLatitude -> -37*Degree,
AspectRatio -> 1,
Mesh -> True,
Text -> False,
DisplayFunction -> Identity],
{year, 1996,1999}];

Following graphic shows RiseSetChart's for HaleBopp from 1996 through to 1998. On the left are the northern hemisphere charts for latitude 37 degrees. On the right are the southern hemisphere charts for latitude -37 degrees.

Remember the green areas are when you can see the object. The red line is when it rises. The blue line when it sets.

Note after about November 1997, Hale-Bopp can't be seen from the northern hemisphere - it is too far south. Southern hemisphere viewer can see it (with telescopes) whenever it is dark.

However, for the middle of 1998, Hale-Bopp dips below the horizon for southern hemisphere viewers. It would then be above the horizon in the northern hemisphere, but it happens in daytime (because it is close to the Sun).

Show[GraphicsArray[
Transpose[{Reverse[no],
Reverse[so]}],
GraphicsSpacing -> {0.1,0}],
     DisplayFunction ->
     $DisplayFunction];

Comet Hyakutake: Orbital Elements

When a new comet is discovered, and its orbital elements computed, the elements are often expressed in a special format. For instance for Comet Hyakutake, the elements were published as follows.

Comet 1996 B2 Hyakutake:

Epoch 1996 May 3.00000

e .999669 Eccentricity

q .230040 Perihelion passage distance (AU)

Node 188.0430 Longitude of ascending node (deg)

w 130.2092 Argument of perihelion (deg)

i 124.9101 Inclination (deg)

Tp 1996 May 1.39647 Perihelion passage time (TDB)

These elements can easily be mapped into the values required by SetOrbitalElements.

SetOrbitalElements[Hyakutake,
     ViewPointSun,
     Date{1996,5,3}, TimeZone0,
     OrbitalEccentricity0.999669,
     OrbitalPerigee0.230040 AU,
     AscendingLongitude188.0430 Degree,
     PerigeeArgument130.2092 Degree,
     OrbitalInclination124.9101 Degree,
     PerigeeDate{1996,5,1.39647}];

One point to remember is that the Date and PerigeeDate options used inside SetOrbitalElements are in local time and not UT. Thus a time zone value is added.

Other Regular Comets

cometdata={{Borrelly,3.6173472,0.6227398,30.27178,
353.36804,74.74414,
294.45626,{9.0,6.0}},
{Brooks2,3.6182664,0.4904944,5.54820,
197.99375,176.24794,
303.27207,{13.5,4.0}},
{DArrest,3.4833082,0.6145437,19.55004,
178.09251,138.37230,
250.02961,{16.5,2.0}},
{Encke,2.2093729,0.8502243,11.93531,
186.27588,334.02432,
302.23501,{14.5,2.0}},
{Faye,3.7987298,0.5777371,9.12066,
204.52051,198.78541,
82.98305,{8.0,6.0}},
{Finlay,3.6153369,0.7079640,3.68814,
322.72252,41.42866,
267.96095,{17.0,2.0}},
{GiacobiniZinner,
3.5228942,0.7061816,31.82992,
172.55570,194.67709,
70.78304,{15.5,2.0}},
{GriggSkjellerup,
2.9640167,0.6643389,21.10557,
359.29906,212.62543,
72.39923,{17.5,2.0}},
{Halley,17.8416126,0.9672619,162.25779,
111.34217,57.75577,
35.73643,{8.5,4.0}},
{Kohoutek,3.5435114,0.4964229,5.90660,
175.86579,268.93530,
310.81035,{14.5,2.0}},
{Kopff,3.4644633,0.5432679,4.72287,
162.75310,120.27729,
196.95778,{13.5,4.0}},
{PonsWinnecke,
3.4361969,0.6339418,22.29255,
172.26599,92.73658,
223.11103,{16.0,2.0}},
{SchwassmannWachmann1,
6.0250481,0.0428167,9.37093,
48.06751,312.10476,
93.30431,{7.5,4.0}},
{SchwassmannWachmann2,
3.4438341,0.3988348,3.75720,
358.13104,125.62449,
332.87576,{14.0,4.0}},
{SwiftTuttle,
26.2657927,0.9635175,113.43126,
152.99457,138.74259,
1.69636,{4.5,6.0}},
{Temple2,3.1077485,0.5225869,11.97741,
194.81958,117.58944,
319.02377,{15.0,2.0}},
{Wolf,4.0827747,0.4054527,27.48451,
162.25824,203.44110,
40.37950,{16.0,2.0}}};

comets=Apply[SetOrbitalElements[#1,
ViewPointSun,
Date{1993,8,1}, TimeZone0,
OrbitalSemiMajorAxis#2 AU,
OrbitalEccentricity#3,
OrbitalInclination#4 Degree,
PerigeeArgument#5 Degree,
AscendingLongitude#6 Degree,
MeanAnomaly#7 Degree]&,
cometdata,{1}]

OrbitPlot[comets];