- Thank you received: 0
Requiem for Relativity
15 years 9 months ago #23662
by Maurol
Replied by Maurol on topic Reply from Mauro Lacy
<blockquote id="quote"><font size="2" face="Verdana, Arial, Helvetica" id="quote">quote:<hr height="1" noshade id="quote"><i>Originally posted by nemesis</i>
Is there any particular reason the equinoxes are preferred to the solstices? It's like the lunar cycle - a lunar month is traditionally reckoned from new moon to new moon but could just as well go from full moon to full moon. That makes more sense to me because it's easier to see when the moon is full.
<hr height="1" noshade id="quote"></blockquote id="quote"></font id="quote">
en.wikipedia.org/wiki/Gregorian_calendar
Read it. All of it.
Just joking. The answer to your question is under "Gregorian reform" in that article. Besides, the complete article is well written and worth reading.
Is there any particular reason the equinoxes are preferred to the solstices? It's like the lunar cycle - a lunar month is traditionally reckoned from new moon to new moon but could just as well go from full moon to full moon. That makes more sense to me because it's easier to see when the moon is full.
<hr height="1" noshade id="quote"></blockquote id="quote"></font id="quote">
en.wikipedia.org/wiki/Gregorian_calendar
Read it. All of it.
Just joking. The answer to your question is under "Gregorian reform" in that article. Besides, the complete article is well written and worth reading.
Please Log in or Create an account to join the conversation.
15 years 9 months ago #15779
by Stoat
Replied by Stoat on topic Reply from Robert Turner
Hi Nemesis, hunter gatherers were hugely interested in the moon, then with the birth of a more sedentary agricultural society, the sun become more important. This was a major cultural crisis.
A bronze age wizard's hat was dug up in Germany, like all wizard's hats it had the moon and stars on it. The astonishing thing was that it showed the nineteen year cycle of the moon. That knowledge had to have taken a long time to establish. The new thinking about henges, is that they are constructed to show people that the setting sun and the rising moon align on a certain day. This restores harmony to the world.
Over the year, the sun rises north or south of east. So mark those two furthest points. In Oedipus Rex the shepherd explains when he found the dumped baby, by talking about the rising of Canupus. When that particular star rose just before the sun, he knew that it was time to bring the sheep down from the mountain pastures. Likewise the star Sirius was important in Egypt, as it brought with it the flooding of the Nile. A social contract with the gods.
People celebrated the equinox and the solstice and they knew the stars in a manner that modern man has totally lost. The newer religions had to fit their festivals into this world view. Problems arose with the Roman Empire. It was big and lasting! The calendar needed to be adjusted, by adding days to the year and setting up the leap year. Only the priests misinterpreted the new set up and had a leap year every three years. That was stopped and the julian calendar was adopted. However, July and August had to have the same number of days. Because Julius Caesar and Augustus Caesar had to have equal status. It worked though. Zoom on to the late middle ages, Easter is falling later and later. The church wants it close to the equinox. so the Gregorian calendar is brought in. Add ten days to the year and a new clause for working out leap years (a year divisible by four, except those that are divisible by 100 but not by 400; 1900 was not a leap year)
You've guessed it, it worked but again was rather messy. Britain and the colonies didn't join in until a hundred plus years later, and of course the general population thought it was a huge scam to rob them of their wages. The Russian October revolution wasn't in October, so rich Times readers had apoplexy over their brandy.
Anyway back to this twenty five thousand year cycle. There's something called the "argument of pericentre," if it is zero degrees, then the planet is closest to the sun in its orbit, at the ascending node of the orbit (for earth the spring equinox) if it's 180 degrees, then the planet is furthest away in its orbit at the ascending node. The argument for the earth is 95 degrees.
This doesn't change but if it did wander slightly over time, then at mid winter we could be a little further out in our slightly eccentric orbit. Let's say that the wander takes about 500 million years. Then we could think that the sun's output varies. Or we could think that the earth's orbit alters periodically in eccentricity. Or, that the tilt of the earth alters periodically. This does have relevance to the ice ages, as ice advanced in both hemispheres at the same time.
I once sent a letter to the New Musical Express, suggesting that we had a sliding rate of time. People's watches would be set by how much money they had. The richer they were the faster their watch would run. Pity it never caught on, as everyone would have learned calculus and no one would have a clue what the real time was.[][8D]
A bronze age wizard's hat was dug up in Germany, like all wizard's hats it had the moon and stars on it. The astonishing thing was that it showed the nineteen year cycle of the moon. That knowledge had to have taken a long time to establish. The new thinking about henges, is that they are constructed to show people that the setting sun and the rising moon align on a certain day. This restores harmony to the world.
Over the year, the sun rises north or south of east. So mark those two furthest points. In Oedipus Rex the shepherd explains when he found the dumped baby, by talking about the rising of Canupus. When that particular star rose just before the sun, he knew that it was time to bring the sheep down from the mountain pastures. Likewise the star Sirius was important in Egypt, as it brought with it the flooding of the Nile. A social contract with the gods.
People celebrated the equinox and the solstice and they knew the stars in a manner that modern man has totally lost. The newer religions had to fit their festivals into this world view. Problems arose with the Roman Empire. It was big and lasting! The calendar needed to be adjusted, by adding days to the year and setting up the leap year. Only the priests misinterpreted the new set up and had a leap year every three years. That was stopped and the julian calendar was adopted. However, July and August had to have the same number of days. Because Julius Caesar and Augustus Caesar had to have equal status. It worked though. Zoom on to the late middle ages, Easter is falling later and later. The church wants it close to the equinox. so the Gregorian calendar is brought in. Add ten days to the year and a new clause for working out leap years (a year divisible by four, except those that are divisible by 100 but not by 400; 1900 was not a leap year)
You've guessed it, it worked but again was rather messy. Britain and the colonies didn't join in until a hundred plus years later, and of course the general population thought it was a huge scam to rob them of their wages. The Russian October revolution wasn't in October, so rich Times readers had apoplexy over their brandy.
Anyway back to this twenty five thousand year cycle. There's something called the "argument of pericentre," if it is zero degrees, then the planet is closest to the sun in its orbit, at the ascending node of the orbit (for earth the spring equinox) if it's 180 degrees, then the planet is furthest away in its orbit at the ascending node. The argument for the earth is 95 degrees.
This doesn't change but if it did wander slightly over time, then at mid winter we could be a little further out in our slightly eccentric orbit. Let's say that the wander takes about 500 million years. Then we could think that the sun's output varies. Or we could think that the earth's orbit alters periodically in eccentricity. Or, that the tilt of the earth alters periodically. This does have relevance to the ice ages, as ice advanced in both hemispheres at the same time.
I once sent a letter to the New Musical Express, suggesting that we had a sliding rate of time. People's watches would be set by how much money they had. The richer they were the faster their watch would run. Pity it never caught on, as everyone would have learned calculus and no one would have a clue what the real time was.[][8D]
Please Log in or Create an account to join the conversation.
15 years 9 months ago #15780
by nemesis
Replied by nemesis on topic Reply from
For what it's worth, the Muslims use a pure lunar calendar and have no problem with their months rotating through the seasons. A lunisolar calendar had been used in pre-Islamic Arabia, but evidently Muhammad personally decreed an uncorrected lunar calendar.
Please Log in or Create an account to join the conversation.
15 years 9 months ago #15783
by Maurol
Replied by Maurol on topic Reply from Mauro Lacy
Hi,
I in no way want to turn this into a religious issue, but as Easter is the most important religious feast for the Christians, and as it is a moveable feast, that depends on the full moon and in the time of the Vernal equinox, for the Catholics and Christians it was important to fix the Vernal equinox.
They did't want the timing of a feast, which depends on astronomical events, to always advance with the passing of the years.
I in no way want to turn this into a religious issue, but as Easter is the most important religious feast for the Christians, and as it is a moveable feast, that depends on the full moon and in the time of the Vernal equinox, for the Catholics and Christians it was important to fix the Vernal equinox.
They did't want the timing of a feast, which depends on astronomical events, to always advance with the passing of the years.
Please Log in or Create an account to join the conversation.
- Joe Keller
- Offline
- Platinum Member
Less
More
- Thank you received: 0
15 years 9 months ago #15784
by Joe Keller
Replied by Joe Keller on topic Reply from
REM program to fit Barbarossa/Frey orbits to 3 data points, predict others
REM and check binary orbit for Kepler's 2nd law
REM initialize constants
PRINT : PRINT
REM GOTO 9000
pi# = 4 * ATN(1): pi180# = pi# / 180: pi2# = 2 * pi#: crit = 0
kk# = 0 - (pi2# / 365.25#) ^ 2 * 332447 / 332001: corr# = 1 / 1.01#
REM kk# is solar accel. const. in AU/day^2, including known planets
REM dimension variables; double precision throughout
REM positive integration weights for 1954, 1986, 1987:
DIM w(3) AS DOUBLE
w(1) = 1 / 4: w(2) = 3 / 8: w(3) = 3 / 8
DIM t(7) AS DOUBLE: DIM ra(7) AS DOUBLE: DIM decl(7) AS DOUBLE
DIM x(7) AS DOUBLE: DIM y(7) AS DOUBLE: DIM z(7) AS DOUBLE
DIM r(7) AS DOUBLE: DIM rpl(7) AS DOUBLE: DIM rmi(7) AS DOUBLE
DIM drdth(7) AS DOUBLE
DIM xe(7) AS DOUBLE: DIM ye(7) AS DOUBLE: DIM ze(7) AS DOUBLE
DIM re(7) AS DOUBLE
DIM rab(7) AS DOUBLE: DIM declb(7) AS DOUBLE: DIM raf(7) AS DOUBLE
DIM declf(7) AS DOUBLE
DIM nax(7) AS DOUBLE: DIM nay(7) AS DOUBLE: DIM naz(7) AS DOUBLE
DIM x0(7) AS DOUBLE: DIM y0(7) AS DOUBLE: DIM z0(7) AS DOUBLE
DIM theta(7, 7) AS DOUBLE: DIM the(7) AS DOUBLE
DIM rapred(7) AS DOUBLE: DIM declpred(7) AS DOUBLE
DIM cc(5, 6) AS DOUBLE: DIM co(5) AS DOUBLE: DIM area(7) AS DOUBLE
REM get times
GOSUB 1000
dt1# = 1 / (t(2) - t(1)): dt2# = 1 / (t(3) - t(2)): dt3# = 2 / (t(3) - t(1))
REM get coordinates
GOSUB 1100
REM get barycentric Earth positions
GOSUB 1200
REM PRINT "Done initializing; now searching..."
REM best mass ratio, radius, dr/dtheta
REM using work and general orbit are q0#=.8898#; r00#=205.2#
REM rp0# = -1.0# rpp0#=-23.3#
nn = 1
FOR i = 1 TO nn
q# = .8898# + 0 * .0001# * i: p# = 1 - q#: rel0# = 10 ^ 9
FOR r0# = 205.2# TO 205.2# STEP .1#
REM dr/dtheta
FOR rp# = -1 TO -1 STEP .1#
REM d2r/dtheta2
FOR rpp# = -23.3# TO -23.3# STEP .1#
REM get dr/dt
rt# = rp# * pi2# / 3080 / 365.25#
REM get d2r/dt2
rtt# = rpp# * (pi2# / 3080 / 365.25#) ^ 2
REM getting c.o.m. coords. by weighted ave. of Barbarossa & Frey coords.,
REM causes < 0.5" error
FOR j = 1 TO 6
ra(j) = rab(j) * q# + raf(j) * p#
decl(j) = declb(j) * q# + declf(j) * p#
zz# = SIN(decl(j)): cs# = COS(decl(j))
xx# = cs# * COS(ra(j)): yy# = cs# * SIN(ra(j))
csthe# = (xx# * xe(j) + yy# * ye(j) + zz# * ze(j)) / re(j)
r1# = r0# + rt# * (t(j) - t(1)) + .5# * rtt# * (t(j) - t(1)) ^ 2
REM cosine rule: rr# ^ 2 + 2 * rr# * re(j) * csthe# + re(j) ^ 2 - r1#^2 = 0
REM use quadratic equation for rr#
b# = 2 * re(j) * csthe#: c# = re(j) ^ 2 - r1# ^ 2
rr# = (0 - b# + SQR(b# ^ 2 - 4 * c#)) * .5#
x(j) = xe(j) + rr# * xx#
y(j) = ye(j) + rr# * yy#
z(j) = ze(j) + rr# * zz#: r(j) = r1#
NEXT j
REM one-time light time correction
IF crit = 0 THEN GOSUB 1050
crit = 1
REM Minimize difference between Newtonian accel.
REM and accel. of trial curve.
REM accel. of trial curve
GOSUB 2000
REM Newtonian accel.
GOSUB 2100
REM difference^2
GOSUB 2200
IF rel# < rel0# THEN GOSUB 3000
NEXT rpp#: NEXT rp#: NEXT r0#
REM IF i = 1 THEN PRINT i; "/"; nn; "done...";
REM IF i > 1 THEN PRINT i; "/"; nn; "...";
PRINT : PRINT "unexplained accel., AU/day^2:"; SQR(rel0#)
PRINT "best mass ratio, initial radius, dr/dtheta & d2r/dtheta2:"
PRINT q0#, r00#, rp0#, rpp0#
PRINT "t(3) radius, midrange of t(1) & t(3) radii";
rmid# = .5# * (r0(1) + r0(3))
PRINT r0(3), rmid#
PRINT "circular period at this distance"
PRINT rmid# ^ 1.5# * SQR(332001 / 332447)
FOR i = 2 TO 3
FOR j = 1 TO i - 1
cs# = (x0(i) * x0(j) + y0(i) * y0(j) + z0(i) * z0(j)) / r0(i) / r0(j)
theta(i, j) = ATN(SQR(1 - cs# ^ 2) / cs#)
PRINT "sector "; j; i; "in radians: ";
PRINT theta(i, j)
PRINT "radians/day"; " ";
om# = theta(i, j) / (t(i) - t(j))
PRINT om#
IF i = 3 AND j = 1 THEN GOSUB 3100
NEXT j: NEXT i
PRINT "excess travel"
PRINT theta(2, 1) + theta(3, 2) - theta(3, 1)
NEXT i
PRINT "first order elliptical extrapolation"
xx# = x0(3) / r0(3): yy# = y0(3) / r0(3): zz# = z0(3) / r0(3)
x0# = x0(1) / r0(1): y0# = y0(1) / r0(1): z0# = z0(1) / r0(1)
cs# = xx# * x0# + yy# * y0# + zz# * z0#
xx# = xx# - cs# * x0#: yy# = yy# - cs# * y0#: zz# = zz# - cs# * z0#
rr# = SQR(xx# ^ 2 + yy# ^ 2 + zz# ^ 2)
xx# = xx# / rr#: yy# = yy# / rr#: zz# = zz# / rr#
th0# = ATN(SQR(1 - cs# ^ 2) / cs#)
REM To get all of them, let i=4 to 7.
FOR i = 6 TO 6
th# = th0# * (t(i) - t(1)) / (t(3) - t(1))
ri# = r00# + .5# * rt0# * (t(i) - t(1)) + .1667# * rtt0# * (t(i) - t(1)) ^ 2
r3# = r00# + .5# * rt0# * (t(3) - t(1)) + .1667# * rtt0# * (t(3) - t(1)) ^ 2
th# = th# * (r3# / ri#) ^ 2
cs# = COS(th#): sn# = SIN(th#)
r# = r00# + rt0# * (t(i) - t(1)) + .5# * rtt0# * (t(i) - t(1)) ^ 2
xpred# = r# * (cs# * x0# + sn# * xx#)
ypred# = r# * (cs# * y0# + sn# * yy#)
zpred# = r# * (cs# * z0# + sn# * zz#)
REM PRINT "barycentric c.o.m. X, Y, Z for i-th time": PRINT
REM PRINT xpred#, ypred#, zpred#
REM PRINT "geocentric X, Y, Z"
xpred# = xpred# - xe(i): ypred# = ypred# - ye(i): zpred# = zpred# - ze(i)
REM PRINT xpred#, ypred#, zpred#
REM PRINT "J2000 celestial coords."
decl# = ATN(zpred# / SQR(xpred# ^ 2 + ypred# ^ 2)) / pi180#
ra# = ATN(ypred# / xpred#) / pi180# + 180
rapred(i) = ra# * pi180#: declpred(i) = decl# * pi180#
PRINT i; "-th point:"; ra#, decl#
PRINT "RA"; INT(ra# / 15); "h"; 4 * (ra# - INT(ra# / 15) * 15); "m"
PRINT "Decl"; -INT(-decl#); "deg"; 60 * (-decl# - INT(-decl#)); "arcmin"
PRINT
NEXT i
600 PRINT "PART 2: binary orbit - A, B, C; and Dec. 2008 U. of Iowa photo"
FOR i = 1 TO 3
y(i) = declf(i) - declb(i): th# = declf(i) * .1102# + declb(i) * .8898#
x(i) = (raf(i) - rab(i)) * COS(th#)
y(i) = y(i) * r0(i): x(i) = x(i) * r0(i)
NEXT i
y(4) = (declf(6) - declpred(6)) / .8898#
x(4) = (raf(6) - rapred(6)) / .8898# * COS(declpred(6))
r# = r00# + rt0# * (t(6) - t(1)) + .5# * rtt0# * (t(6) - t(1)) ^ 2
x(4) = x(4) * r#: y(4) = y(4) * r#
t(4) = t(6)
REM try to correct for precession
REM slope of 1986-1987 line is 54deg15'+180
dth# = .159#
REM precession correction estimated graphically is .157#
REM this is the aspect change expected from orbital motion
cs# = .81157#: sn# = .58425#
x(1) = x(1) + cs# * dth#
x(4) = x(4) - cs# * dth#
y(1) = y(1) + sn# * dth#
y(4) = y(4) - sn# * dth#
FOR i = 1 TO 4
PRINT x(i), y(i)
NEXT i
PRINT
PRINT (y(3) - y(2)) / (x(3) - x(2)), (y(4) - y(3)) / (x(4) - x(3))
PRINT
REM vary 5th point & fit ellipse
y(5) = 0: rel0# = 10 ^ 9
FOR x0# = .355# TO .355# STEP .001#
x(5) = -x0#
FOR i = 1 TO 5
cc(i, 1) = x(i) ^ 2: cc(i, 2) = x(i) * y(i): cc(i, 3) = y(i) ^ 2
cc(i, 4) = x(i): cc(i, 5) = y(i): cc(i, 6) = 1
NEXT i
REM solve for ellipse coeffs.
FOR i = 1 TO 5
FOR j = 1 TO 5
IF j = i GOTO 650
f# = cc(j, i) / cc(i, i)
FOR k = 1 TO 6
cc(j, k) = cc(j, k) - cc(i, k) * f#
NEXT k
650 NEXT j
NEXT i
FOR i = 1 TO 5
co(i) = cc(i, 6) / cc(i, i)
NEXT i
discrim# = co(2) ^ 2 - 4 * co(1) * co(3)
IF discrim# > 0 GOTO 690
REM PRINT x0#, discrim#
REM find relative standard deviation of areal speed
GOSUB 4000
690 NEXT x0#
PRINT "5th x-intercept, best relative std. dev. of areal speed"
PRINT -x00#, rel0#
990 END
REM Get Julian date - 2400000.
1000 t(1) = 34798.5# + 8 / 24 + 5 / 1440
t(2) = 46504.5# + 14 / 24 + 18 / 1440
t(3) = t(2) + 365 - 28 - 31 + 16 + 3 / 24 - 16 / 1440
t(4) = t(3) + 10 * 365 + 3 + 28 + 3 - 4 / 24 + 57 / 1440
tfrey3# = 54184.5# + 42 / 1440
tbarbarossa3# = 54191.5# + 7 / 24 + 39 / 1440
REM Interpolate t(5) using est. of mass ratio.
t(5) = .8743# * tbarbarossa3# + .1257# * tfrey3#
REM Guess that 12/22/08 epoch, t(6), is approx. time on meridian;
REM Barbarossa was near stationarity anyway.
t(6) = 54822.5# + 12 / 24 + 47 / 1440
t(7) = 54851.5# + 9 / 24 + 48.5# / 1440
RETURN
REM light time correction
1050 FOR j = 1 TO 7
t(j) = t(j) - rr# * 1.496 * 10 ^ 13 / (2.99793 * 10 ^ 10)
NEXT j
dt1# = 1 / (t(2) - t(1)): dt2# = 1 / (t(3) - t(2)): dt3# = 2 / (t(3) - t(1))
RETURN
REM get J2000 geocentric Barbarossa & Frey coordinates
REM "A" plate (1954); objects A2 & A
1100 raf(1) = pi180# * (11 * 15 + 2 / 4 + 25.16# / 240)
declf(1) = 0 - pi180# * (5 + 56 / 60 + 11.3# / 3600)
rab(1) = pi180# * (11 * 15 + 3 / 4 + 12.4# / 240)
declb(1) = 0 - pi180# * (5 + 58 / 60 + 9 / 3600)
REM "B" plate (1986); objects B3 & B
raf(2) = pi180# * (11 * 15 + 16 / 4 + 51.67# / 240)
declf(2) = 0 - pi180# * (7 + 49 / 60 + 40.4# / 3600)
rab(2) = pi180# * (11 * 15 + 16 / 4 + 55.78# / 240)
declb(2) = 0 - pi180# * (7 + 55 / 60 + 14# / 3600)
REM "C" plate (1987); objects C & C11
raf(3) = pi180# * (11 * 15 + 18 / 4 + 3.18# / 240)
declf(3) = 0 - pi180# * (7 + 58 / 60 + 46.1 / 3600)
rab(3) = pi180# * (11 * 15 + 18 / 4 + .41 / 240)
declb(3) = 0 - pi180# * (8 + 1 / 60 + 57.7 / 3600)
REM "D" plate (1997); objects uncertain
REM Genebriera March 25, 2007
raf(5) = pi180# * (11 * 15 + 26 / 4 + 22.2# / 240)
declf(5) = 0 - pi180# * (9 + 4 / 60 + 59 / 3600)
REM Riley April 1, 2007
rab(5) = pi180# * (11 * 15 + 26 / 4 + 25 / 240)
declb(5) = 0 - pi180# * (8 + 57 / 60 + 26 / 3600)
REM alternate coords. for Riley's "Barbarossa"
REM (see messageboard May 7, 2007; called "Frey" there)
REM rab(5) = pi180# * (11 * 15 + 26 / 4 + 24.6 / 240)
REM declb(5) = 0 - pi180# * (8 + 57 / 60 + 48.5 / 3600)
REM Genebriera's "Barbarossa"
REM (see messageboard May 15, 2007; called "Frey" there)
REM rab(5) = pi180# * (11 * 15 + 26 / 4 + 31.8 / 240)
REM declb(5) = 0 - pi180# * (9 + 0 / 60 + 11 / 3600)
REM Dec. 22, 2008 U. of Iowa photo
raf(6) = pi180# * (11 * 15 + 27 / 4 + 30.17# / 240)
declf(6) = 0 - pi180# * (9 + 21 / 60 + 48.6# / 3600)
REM (Barbarossa would have been off west edge)
REM objects identified for Dec. 22 in earlier theory
REM rab(6) = pi180# * (11 * 15 + 28 / 4 + 22.08# / 240)
REM declb(6) = 0 - pi180# * (9 + 16 / 60 + 6.4# / 3600)
REM raf(6) = pi180# * (11 * 15 + 29 / 4 + 4.66# / 240)
REM declf(6) = 0 - pi180# * (9 + 7 / 60 + 2.3# / 3600)
RETURN
REM Get J2000 barycentric Earth XYZ position (excludes Barbarossa).
REM In this program, "barycentric" = sun + traditional planets.
REM all but 1954, from Astronomical Almanac
1200 a# = 14 / 24 + 18 / 1440: i = 2
GOSUB 1220
a# = 17 / 24 + 2 / 1440: i = 3
GOSUB 1220
a# = 13 / 24 + 59 / 1440: i = 4
GOSUB 1220
REM as in subroutine 1000, weighted time for 2007
REM warning: need to consult Almanac to accurize
a# = 42 / 1440 * .1102# + (7 + 7 / 24 + 39 / 1440) * .8898#: i = 5
GOSUB 1220
REM as in subroutine 1000, guess that 12/2008 epoch was meridian
a# = 12 / 24 + 47 / 1440: i = 6
GOSUB 1220
a# = 9 / 24 + 48.5# / 1440: i = 7
GOSUB 1220
GOSUB 1250
RETURN
1220 b# = 1 - a#
READ xe1#, xe2#, ye1#, ye2#, ze1#, ze2#
xe(i) = b# * xe1# + a# * xe2#
ye(i) = b# * ye1# + a# * ye2#
ze(i) = b# * ze1# + a# * ze2#
re(i) = SQR(xe(i) ^ 2 + ye(i) ^ 2 + ze(i) ^ 2)
RETURN
REM find 1954 barycentric Earth position
1250 msun# = 332000: mjup# = 318.1#
msat# = 95.2#: mur# = 14.6#: mnep# = 17.2#
me# = 1: mtot# = msun# + mjup# + msat# + mur# + mnep# + me#
sx# = 0: sy# = 0: sz# = 0
REM interpolate sun coordinates
a# = 8 / 24 + 5 / 1440: b# = 1 - a#
raa# = 22 * 15 + 30 / 4 + 33.54# / 240
rab# = 22 * 15 + 34 / 4 + 20.69# / 240
ra# = (raa# * b# + rab# * a#) * pi180#
decla# = 9 + 22 / 60 + 8.300000000000001# / 3600
declb# = 8 + 59 / 60 + 52 / 3600
decl# = 0 - (decla# * b# + declb# * a#) * pi180#
ra# = .9899423000000001#: rb# = .9901818999999999#
r# = ra# * b# + rb# * a#: w# = msun#
GOSUB 1270
REM planets interpolated to nearest RA second, or Decl arcminute
REM Jupiter
ra# = pi180# * (5 * 15 + 2 / 4 + 41 / 240)
decl# = pi180# * (22 + 29 / 60 + 51 / 3600)
r# = 4.846#: w# = mjup#
GOSUB 1270
REM Saturn
ra# = pi180# * (14 * 15 + 30 / 4 + 57 / 240)
decl = 0 - pi180# * (12 + 9 / 60 + 22 / 3600)
r# = 9.308#: w# = msat#
GOSUB 1270
REM Uranus
ra# = pi180# * (7 * 15 + 24 / 4 + 19 / 240)
decl# = pi180# * (22 + 31 / 60 + 46 / 3600)
r# = 18.001#: w# = mur#
GOSUB 1270
REM Neptune
ra# = pi180# * (13 * 15 + 38 / 4 + 22 / 240)
decl# = 0 - pi180# * (8 + 22 / 60 + 9 / 3600)
r# = 29.662#: w# = mnep#
GOSUB 1270
REM convert to 1950 ecliptic coords.
theta# = 23.45# * pi180#
GOSUB 1290
REM convert to 2000 ecliptic coords.
theta# = pi180# * 50.2619# / 3600 * 50
cs# = COS(theta#): sn# = SIN(theta#): sx0# = sx#
sx# = cs# * sx# - sn# * sy#
sy# = cs# * sy# + sn# * sx0#
REM convert to 2000 celestial coords.
theta# = 0 - 23.45# * pi180#
GOSUB 1290
xe(1) = 0 - sx# / mtot#
ye(1) = 0 - sy# / mtot#
ze(1) = 0 - sz# / mtot#
re(1) = SQR(xe(1) ^ 2 + ye(1) ^ 2 + ze(1) ^ 2)
RETURN
1270 z# = r# * SIN(decl#): cs# = COS(decl#)
x# = r# * cs# * COS(ra#): y# = r# * cs# * SIN(ra#)
sx# = sx# + w# * x#
sy# = sy# + w# * y#
sz# = sz# + w# * z#
RETURN
1290 cs# = COS(theta#): sn# = SIN(theta#): sz0# = sz#
sz# = cs# * sz# - sn# * sy#
sy# = cs# * sy# + sn# * sz0#
RETURN
REM find the accel. of the parabolas fitting 1st thru 3rd points
2000 vx1# = (x(2) - x(1)) * dt1#
vy1# = (y(2) - y(1)) * dt1#
vz1# = (z(2) - z(1)) * dt1#
vx2# = (x(3) - x(2)) * dt2#
vy2# = (y(3) - y(2)) * dt2#
vz2# = (z(3) - z(2)) * dt2#
REM accel., corrected to fixed frame
ax# = (vx2# - vx1#) * dt3# * corr#
ay# = (vy2# - vy1#) * dt3# * corr#
az# = (vz2# - vz1#) * dt3# * corr#
RETURN
REM find Newtonian accel. at 1st thru 3rd points
2100 FOR j = 1 TO 3
dr# = kk# / r(j) ^ 3
nax(j) = x(j) * dr#
nay(j) = y(j) * dr#
naz(j) = z(j) * dr#
NEXT j
RETURN
REM weighted ave. Newtonian accel., and work done
2200 avenax# = 0: avenay# = 0: avenaz# = 0
FOR j = 1 TO 3
avenax# = avenax# + w(j) * nax(j)
avenay# = avenay# + w(j) * nay(j)
avenaz# = avenaz# + w(j) * naz(j)
NEXT j
rel# = (avenax# - ax#) ^ 2 + (avenay# - ay#) ^ 2 + (avenaz# - az#) ^ 2
RETURN
REM update best
3000 rel0# = rel#: i0 = i: q0# = q#
r00# = r0#: rp0# = rp#: rt0# = rt#: rpp0# = rpp#: rtt0# = rtt#
FOR m = 1 TO 6
x0(m) = x(m): y0(m) = y(m): z0(m) = z(m): r0(m) = r(m)
NEXT m
RETURN
3100 PRINT "apparent circular period"
PRINT pi2# / om# / 365.25#
RETURN
REM find areas of the four sectors
REM PRINT "radii"
4000 FOR i = 1 TO 4
the(i) = ATN(-y(i) / x(i))
IF x(i) > 0 THEN GOSUB 4100
IF the(i) < 0 THEN GOSUB 4120
r(i) = SQR(x(i) ^ 2 + y(i) ^ 2)
REM PRINT r(i),
NEXT i
REM PRINT
r(5) = r(1): the(5) = the(1) + pi2#
REM PRINT "drdth"
FOR i = 1 TO 4
u# = 2 * co(1) * x(i) + co(2) * y(i) + co(4)
w# = co(2) * x(i) + 2 * co(3) * y(i) + co(5)
sn# = 0 - (u# * y(i) - w# * x(i)) / SQR(u# ^ 2 + w# ^ 2)
drdth(i) = sn# / SQR(1 - sn# ^ 2)
REM PRINT drdth(i),
NEXT i
REM PRINT
drdth(5) = drdth(1)
PRINT "interpolated radii"
REM interpolate r quadratically
FOR i = 1 TO 5
IF i < 5 THEN dthpl# = the(i + 1) - the(i)
IF i < 5 THEN drpl# = r(i + 1) - r(i)
IF i > 1 THEN dthmi# = the(i) - the(i - 1)
IF i > 1 THEN drmi# = r(i) - r(i - 1)
IF i < 5 THEN rpl(i) = r(i) + 2 / 9 * drdth(i) * dthpl# + 1 / 9 * drpl#
IF i > 1 THEN rmi(i) = r(i) - 2 / 9 * drdth(i) * dthmi# - 1 / 9 * drmi#
PRINT rpl(i), rmi(i)
NEXT i
REM Cotes' 3/8 rule
REM PRINT "areas"
FOR i = 1 TO 4
r2m# = r(i) ^ 2 + r(i + 1) ^ 2 + 3 * (rpl(i) ^ 2 + rmi(i + 1) ^ 2)
area(i) = r2m# / 8 * (the(i + 1) - the(i)) / 2
REM PRINT area(i),
NEXT i
REM PRINT
are# = area(1) + area(2) + area(3) + area(4)
s# = 0: ss# = 0
PRINT "rates"
FOR i = 1 TO 3
IF i = 1 THEN lap# = 1
IF i = 2 THEN lap# = 0
IF i = 3 THEN lap# = 1
rate# = (area(i) + lap# * are#) / (t(i + 1) - t(i))
PRINT rate#; " ";
s# = s# + rate#: ss# = ss# + rate# ^ 2
NEXT i
rel# = SQR((ss# - s# ^ 2 / 3) / 2) / (s# / 3)
IF rel# < rel0# THEN GOSUB 4200
RETURN
4100 the(i) = the(i) + pi#
RETURN
4120 the(i) = the(i) + pi2#
RETURN
4200 rel0# = rel#: x00# = x0#
RETURN
REM Astronomical Almanac barycentric Earth XYZ '86
DATA -.992473#,-.99431#,.097062#,.081274#,.042041#,.035195#
REM '87
DATA -.647098#,-.660316#,.689355#,.678874#,.298895#,.29435#
REM '97
DATA -.951971#,-.957296#,.278791#,.263632#,.121036#,.114464#
REM '07
DATA -.992917#,-.991887#,-.057633#,-.073429#,-.025095#,-.031943#
REM '08
DATA -.008637#,-.026129#,.906672#,.906369#,.39304#,.392909#
REM '09
DATA -.493191#,-.508316#,.786488#,.778431#,.340941#,.337449#
REM (times of plates, Sun positions,
REM coords. of identified disappearing dots)
REM plate A (Red)
REM DATA 1954.154,.9900229, 22,34, 27.6,-8,-59,-7
REM object A2
REM DATA 11,2,25.47,-5,-56,-12
REM object A
REM DATA 11,3,12.4,-5,-58,-9
REM plate B (Red)
REM DATA 1986.201,.9945905,23,40,30.9,-2,-10,-49
REM object B3
REM DATA 11,16,51.67,-7,-49,-40.4
REM object B
REM DATA 11,16,55.78,-7,-55,-14
REM plate C (Red)
REM DATA 1987.08215,.9852006,20,53,6.9,-17,31,13
REM object C11 (named for easy recall, like "CII")
REM DATA 11,18,.41,-8,-1,-57.7
REM object C
REM DATA 11,18,3.18,-7,-58,-46.1
REM plate D (optical infrared)
REM DATA 1997.16711,.9914845,22,57,21.1,-6,-39,-54
REM DATA 11,22,32.9,-8,-26,-56
REM DATA 11,22,16.77,-8,-29,-32.6
REM DATA 11,22,1.33,-8,-34,-36.9
REM and check binary orbit for Kepler's 2nd law
REM initialize constants
PRINT : PRINT
REM GOTO 9000
pi# = 4 * ATN(1): pi180# = pi# / 180: pi2# = 2 * pi#: crit = 0
kk# = 0 - (pi2# / 365.25#) ^ 2 * 332447 / 332001: corr# = 1 / 1.01#
REM kk# is solar accel. const. in AU/day^2, including known planets
REM dimension variables; double precision throughout
REM positive integration weights for 1954, 1986, 1987:
DIM w(3) AS DOUBLE
w(1) = 1 / 4: w(2) = 3 / 8: w(3) = 3 / 8
DIM t(7) AS DOUBLE: DIM ra(7) AS DOUBLE: DIM decl(7) AS DOUBLE
DIM x(7) AS DOUBLE: DIM y(7) AS DOUBLE: DIM z(7) AS DOUBLE
DIM r(7) AS DOUBLE: DIM rpl(7) AS DOUBLE: DIM rmi(7) AS DOUBLE
DIM drdth(7) AS DOUBLE
DIM xe(7) AS DOUBLE: DIM ye(7) AS DOUBLE: DIM ze(7) AS DOUBLE
DIM re(7) AS DOUBLE
DIM rab(7) AS DOUBLE: DIM declb(7) AS DOUBLE: DIM raf(7) AS DOUBLE
DIM declf(7) AS DOUBLE
DIM nax(7) AS DOUBLE: DIM nay(7) AS DOUBLE: DIM naz(7) AS DOUBLE
DIM x0(7) AS DOUBLE: DIM y0(7) AS DOUBLE: DIM z0(7) AS DOUBLE
DIM theta(7, 7) AS DOUBLE: DIM the(7) AS DOUBLE
DIM rapred(7) AS DOUBLE: DIM declpred(7) AS DOUBLE
DIM cc(5, 6) AS DOUBLE: DIM co(5) AS DOUBLE: DIM area(7) AS DOUBLE
REM get times
GOSUB 1000
dt1# = 1 / (t(2) - t(1)): dt2# = 1 / (t(3) - t(2)): dt3# = 2 / (t(3) - t(1))
REM get coordinates
GOSUB 1100
REM get barycentric Earth positions
GOSUB 1200
REM PRINT "Done initializing; now searching..."
REM best mass ratio, radius, dr/dtheta
REM using work and general orbit are q0#=.8898#; r00#=205.2#
REM rp0# = -1.0# rpp0#=-23.3#
nn = 1
FOR i = 1 TO nn
q# = .8898# + 0 * .0001# * i: p# = 1 - q#: rel0# = 10 ^ 9
FOR r0# = 205.2# TO 205.2# STEP .1#
REM dr/dtheta
FOR rp# = -1 TO -1 STEP .1#
REM d2r/dtheta2
FOR rpp# = -23.3# TO -23.3# STEP .1#
REM get dr/dt
rt# = rp# * pi2# / 3080 / 365.25#
REM get d2r/dt2
rtt# = rpp# * (pi2# / 3080 / 365.25#) ^ 2
REM getting c.o.m. coords. by weighted ave. of Barbarossa & Frey coords.,
REM causes < 0.5" error
FOR j = 1 TO 6
ra(j) = rab(j) * q# + raf(j) * p#
decl(j) = declb(j) * q# + declf(j) * p#
zz# = SIN(decl(j)): cs# = COS(decl(j))
xx# = cs# * COS(ra(j)): yy# = cs# * SIN(ra(j))
csthe# = (xx# * xe(j) + yy# * ye(j) + zz# * ze(j)) / re(j)
r1# = r0# + rt# * (t(j) - t(1)) + .5# * rtt# * (t(j) - t(1)) ^ 2
REM cosine rule: rr# ^ 2 + 2 * rr# * re(j) * csthe# + re(j) ^ 2 - r1#^2 = 0
REM use quadratic equation for rr#
b# = 2 * re(j) * csthe#: c# = re(j) ^ 2 - r1# ^ 2
rr# = (0 - b# + SQR(b# ^ 2 - 4 * c#)) * .5#
x(j) = xe(j) + rr# * xx#
y(j) = ye(j) + rr# * yy#
z(j) = ze(j) + rr# * zz#: r(j) = r1#
NEXT j
REM one-time light time correction
IF crit = 0 THEN GOSUB 1050
crit = 1
REM Minimize difference between Newtonian accel.
REM and accel. of trial curve.
REM accel. of trial curve
GOSUB 2000
REM Newtonian accel.
GOSUB 2100
REM difference^2
GOSUB 2200
IF rel# < rel0# THEN GOSUB 3000
NEXT rpp#: NEXT rp#: NEXT r0#
REM IF i = 1 THEN PRINT i; "/"; nn; "done...";
REM IF i > 1 THEN PRINT i; "/"; nn; "...";
PRINT : PRINT "unexplained accel., AU/day^2:"; SQR(rel0#)
PRINT "best mass ratio, initial radius, dr/dtheta & d2r/dtheta2:"
PRINT q0#, r00#, rp0#, rpp0#
PRINT "t(3) radius, midrange of t(1) & t(3) radii";
rmid# = .5# * (r0(1) + r0(3))
PRINT r0(3), rmid#
PRINT "circular period at this distance"
PRINT rmid# ^ 1.5# * SQR(332001 / 332447)
FOR i = 2 TO 3
FOR j = 1 TO i - 1
cs# = (x0(i) * x0(j) + y0(i) * y0(j) + z0(i) * z0(j)) / r0(i) / r0(j)
theta(i, j) = ATN(SQR(1 - cs# ^ 2) / cs#)
PRINT "sector "; j; i; "in radians: ";
PRINT theta(i, j)
PRINT "radians/day"; " ";
om# = theta(i, j) / (t(i) - t(j))
PRINT om#
IF i = 3 AND j = 1 THEN GOSUB 3100
NEXT j: NEXT i
PRINT "excess travel"
PRINT theta(2, 1) + theta(3, 2) - theta(3, 1)
NEXT i
PRINT "first order elliptical extrapolation"
xx# = x0(3) / r0(3): yy# = y0(3) / r0(3): zz# = z0(3) / r0(3)
x0# = x0(1) / r0(1): y0# = y0(1) / r0(1): z0# = z0(1) / r0(1)
cs# = xx# * x0# + yy# * y0# + zz# * z0#
xx# = xx# - cs# * x0#: yy# = yy# - cs# * y0#: zz# = zz# - cs# * z0#
rr# = SQR(xx# ^ 2 + yy# ^ 2 + zz# ^ 2)
xx# = xx# / rr#: yy# = yy# / rr#: zz# = zz# / rr#
th0# = ATN(SQR(1 - cs# ^ 2) / cs#)
REM To get all of them, let i=4 to 7.
FOR i = 6 TO 6
th# = th0# * (t(i) - t(1)) / (t(3) - t(1))
ri# = r00# + .5# * rt0# * (t(i) - t(1)) + .1667# * rtt0# * (t(i) - t(1)) ^ 2
r3# = r00# + .5# * rt0# * (t(3) - t(1)) + .1667# * rtt0# * (t(3) - t(1)) ^ 2
th# = th# * (r3# / ri#) ^ 2
cs# = COS(th#): sn# = SIN(th#)
r# = r00# + rt0# * (t(i) - t(1)) + .5# * rtt0# * (t(i) - t(1)) ^ 2
xpred# = r# * (cs# * x0# + sn# * xx#)
ypred# = r# * (cs# * y0# + sn# * yy#)
zpred# = r# * (cs# * z0# + sn# * zz#)
REM PRINT "barycentric c.o.m. X, Y, Z for i-th time": PRINT
REM PRINT xpred#, ypred#, zpred#
REM PRINT "geocentric X, Y, Z"
xpred# = xpred# - xe(i): ypred# = ypred# - ye(i): zpred# = zpred# - ze(i)
REM PRINT xpred#, ypred#, zpred#
REM PRINT "J2000 celestial coords."
decl# = ATN(zpred# / SQR(xpred# ^ 2 + ypred# ^ 2)) / pi180#
ra# = ATN(ypred# / xpred#) / pi180# + 180
rapred(i) = ra# * pi180#: declpred(i) = decl# * pi180#
PRINT i; "-th point:"; ra#, decl#
PRINT "RA"; INT(ra# / 15); "h"; 4 * (ra# - INT(ra# / 15) * 15); "m"
PRINT "Decl"; -INT(-decl#); "deg"; 60 * (-decl# - INT(-decl#)); "arcmin"
NEXT i
600 PRINT "PART 2: binary orbit - A, B, C; and Dec. 2008 U. of Iowa photo"
FOR i = 1 TO 3
y(i) = declf(i) - declb(i): th# = declf(i) * .1102# + declb(i) * .8898#
x(i) = (raf(i) - rab(i)) * COS(th#)
y(i) = y(i) * r0(i): x(i) = x(i) * r0(i)
NEXT i
y(4) = (declf(6) - declpred(6)) / .8898#
x(4) = (raf(6) - rapred(6)) / .8898# * COS(declpred(6))
r# = r00# + rt0# * (t(6) - t(1)) + .5# * rtt0# * (t(6) - t(1)) ^ 2
x(4) = x(4) * r#: y(4) = y(4) * r#
t(4) = t(6)
REM try to correct for precession
REM slope of 1986-1987 line is 54deg15'+180
dth# = .159#
REM precession correction estimated graphically is .157#
REM this is the aspect change expected from orbital motion
cs# = .81157#: sn# = .58425#
x(1) = x(1) + cs# * dth#
x(4) = x(4) - cs# * dth#
y(1) = y(1) + sn# * dth#
y(4) = y(4) - sn# * dth#
FOR i = 1 TO 4
PRINT x(i), y(i)
NEXT i
PRINT (y(3) - y(2)) / (x(3) - x(2)), (y(4) - y(3)) / (x(4) - x(3))
REM vary 5th point & fit ellipse
y(5) = 0: rel0# = 10 ^ 9
FOR x0# = .355# TO .355# STEP .001#
x(5) = -x0#
FOR i = 1 TO 5
cc(i, 1) = x(i) ^ 2: cc(i, 2) = x(i) * y(i): cc(i, 3) = y(i) ^ 2
cc(i, 4) = x(i): cc(i, 5) = y(i): cc(i, 6) = 1
NEXT i
REM solve for ellipse coeffs.
FOR i = 1 TO 5
FOR j = 1 TO 5
IF j = i GOTO 650
f# = cc(j, i) / cc(i, i)
FOR k = 1 TO 6
cc(j, k) = cc(j, k) - cc(i, k) * f#
NEXT k
650 NEXT j
NEXT i
FOR i = 1 TO 5
co(i) = cc(i, 6) / cc(i, i)
NEXT i
discrim# = co(2) ^ 2 - 4 * co(1) * co(3)
IF discrim# > 0 GOTO 690
REM PRINT x0#, discrim#
REM find relative standard deviation of areal speed
GOSUB 4000
690 NEXT x0#
PRINT "5th x-intercept, best relative std. dev. of areal speed"
PRINT -x00#, rel0#
990 END
REM Get Julian date - 2400000.
1000 t(1) = 34798.5# + 8 / 24 + 5 / 1440
t(2) = 46504.5# + 14 / 24 + 18 / 1440
t(3) = t(2) + 365 - 28 - 31 + 16 + 3 / 24 - 16 / 1440
t(4) = t(3) + 10 * 365 + 3 + 28 + 3 - 4 / 24 + 57 / 1440
tfrey3# = 54184.5# + 42 / 1440
tbarbarossa3# = 54191.5# + 7 / 24 + 39 / 1440
REM Interpolate t(5) using est. of mass ratio.
t(5) = .8743# * tbarbarossa3# + .1257# * tfrey3#
REM Guess that 12/22/08 epoch, t(6), is approx. time on meridian;
REM Barbarossa was near stationarity anyway.
t(6) = 54822.5# + 12 / 24 + 47 / 1440
t(7) = 54851.5# + 9 / 24 + 48.5# / 1440
RETURN
REM light time correction
1050 FOR j = 1 TO 7
t(j) = t(j) - rr# * 1.496 * 10 ^ 13 / (2.99793 * 10 ^ 10)
NEXT j
dt1# = 1 / (t(2) - t(1)): dt2# = 1 / (t(3) - t(2)): dt3# = 2 / (t(3) - t(1))
RETURN
REM get J2000 geocentric Barbarossa & Frey coordinates
REM "A" plate (1954); objects A2 & A
1100 raf(1) = pi180# * (11 * 15 + 2 / 4 + 25.16# / 240)
declf(1) = 0 - pi180# * (5 + 56 / 60 + 11.3# / 3600)
rab(1) = pi180# * (11 * 15 + 3 / 4 + 12.4# / 240)
declb(1) = 0 - pi180# * (5 + 58 / 60 + 9 / 3600)
REM "B" plate (1986); objects B3 & B
raf(2) = pi180# * (11 * 15 + 16 / 4 + 51.67# / 240)
declf(2) = 0 - pi180# * (7 + 49 / 60 + 40.4# / 3600)
rab(2) = pi180# * (11 * 15 + 16 / 4 + 55.78# / 240)
declb(2) = 0 - pi180# * (7 + 55 / 60 + 14# / 3600)
REM "C" plate (1987); objects C & C11
raf(3) = pi180# * (11 * 15 + 18 / 4 + 3.18# / 240)
declf(3) = 0 - pi180# * (7 + 58 / 60 + 46.1 / 3600)
rab(3) = pi180# * (11 * 15 + 18 / 4 + .41 / 240)
declb(3) = 0 - pi180# * (8 + 1 / 60 + 57.7 / 3600)
REM "D" plate (1997); objects uncertain
REM Genebriera March 25, 2007
raf(5) = pi180# * (11 * 15 + 26 / 4 + 22.2# / 240)
declf(5) = 0 - pi180# * (9 + 4 / 60 + 59 / 3600)
REM Riley April 1, 2007
rab(5) = pi180# * (11 * 15 + 26 / 4 + 25 / 240)
declb(5) = 0 - pi180# * (8 + 57 / 60 + 26 / 3600)
REM alternate coords. for Riley's "Barbarossa"
REM (see messageboard May 7, 2007; called "Frey" there)
REM rab(5) = pi180# * (11 * 15 + 26 / 4 + 24.6 / 240)
REM declb(5) = 0 - pi180# * (8 + 57 / 60 + 48.5 / 3600)
REM Genebriera's "Barbarossa"
REM (see messageboard May 15, 2007; called "Frey" there)
REM rab(5) = pi180# * (11 * 15 + 26 / 4 + 31.8 / 240)
REM declb(5) = 0 - pi180# * (9 + 0 / 60 + 11 / 3600)
REM Dec. 22, 2008 U. of Iowa photo
raf(6) = pi180# * (11 * 15 + 27 / 4 + 30.17# / 240)
declf(6) = 0 - pi180# * (9 + 21 / 60 + 48.6# / 3600)
REM (Barbarossa would have been off west edge)
REM objects identified for Dec. 22 in earlier theory
REM rab(6) = pi180# * (11 * 15 + 28 / 4 + 22.08# / 240)
REM declb(6) = 0 - pi180# * (9 + 16 / 60 + 6.4# / 3600)
REM raf(6) = pi180# * (11 * 15 + 29 / 4 + 4.66# / 240)
REM declf(6) = 0 - pi180# * (9 + 7 / 60 + 2.3# / 3600)
RETURN
REM Get J2000 barycentric Earth XYZ position (excludes Barbarossa).
REM In this program, "barycentric" = sun + traditional planets.
REM all but 1954, from Astronomical Almanac
1200 a# = 14 / 24 + 18 / 1440: i = 2
GOSUB 1220
a# = 17 / 24 + 2 / 1440: i = 3
GOSUB 1220
a# = 13 / 24 + 59 / 1440: i = 4
GOSUB 1220
REM as in subroutine 1000, weighted time for 2007
REM warning: need to consult Almanac to accurize
a# = 42 / 1440 * .1102# + (7 + 7 / 24 + 39 / 1440) * .8898#: i = 5
GOSUB 1220
REM as in subroutine 1000, guess that 12/2008 epoch was meridian
a# = 12 / 24 + 47 / 1440: i = 6
GOSUB 1220
a# = 9 / 24 + 48.5# / 1440: i = 7
GOSUB 1220
GOSUB 1250
RETURN
1220 b# = 1 - a#
READ xe1#, xe2#, ye1#, ye2#, ze1#, ze2#
xe(i) = b# * xe1# + a# * xe2#
ye(i) = b# * ye1# + a# * ye2#
ze(i) = b# * ze1# + a# * ze2#
re(i) = SQR(xe(i) ^ 2 + ye(i) ^ 2 + ze(i) ^ 2)
RETURN
REM find 1954 barycentric Earth position
1250 msun# = 332000: mjup# = 318.1#
msat# = 95.2#: mur# = 14.6#: mnep# = 17.2#
me# = 1: mtot# = msun# + mjup# + msat# + mur# + mnep# + me#
sx# = 0: sy# = 0: sz# = 0
REM interpolate sun coordinates
a# = 8 / 24 + 5 / 1440: b# = 1 - a#
raa# = 22 * 15 + 30 / 4 + 33.54# / 240
rab# = 22 * 15 + 34 / 4 + 20.69# / 240
ra# = (raa# * b# + rab# * a#) * pi180#
decla# = 9 + 22 / 60 + 8.300000000000001# / 3600
declb# = 8 + 59 / 60 + 52 / 3600
decl# = 0 - (decla# * b# + declb# * a#) * pi180#
ra# = .9899423000000001#: rb# = .9901818999999999#
r# = ra# * b# + rb# * a#: w# = msun#
GOSUB 1270
REM planets interpolated to nearest RA second, or Decl arcminute
REM Jupiter
ra# = pi180# * (5 * 15 + 2 / 4 + 41 / 240)
decl# = pi180# * (22 + 29 / 60 + 51 / 3600)
r# = 4.846#: w# = mjup#
GOSUB 1270
REM Saturn
ra# = pi180# * (14 * 15 + 30 / 4 + 57 / 240)
decl = 0 - pi180# * (12 + 9 / 60 + 22 / 3600)
r# = 9.308#: w# = msat#
GOSUB 1270
REM Uranus
ra# = pi180# * (7 * 15 + 24 / 4 + 19 / 240)
decl# = pi180# * (22 + 31 / 60 + 46 / 3600)
r# = 18.001#: w# = mur#
GOSUB 1270
REM Neptune
ra# = pi180# * (13 * 15 + 38 / 4 + 22 / 240)
decl# = 0 - pi180# * (8 + 22 / 60 + 9 / 3600)
r# = 29.662#: w# = mnep#
GOSUB 1270
REM convert to 1950 ecliptic coords.
theta# = 23.45# * pi180#
GOSUB 1290
REM convert to 2000 ecliptic coords.
theta# = pi180# * 50.2619# / 3600 * 50
cs# = COS(theta#): sn# = SIN(theta#): sx0# = sx#
sx# = cs# * sx# - sn# * sy#
sy# = cs# * sy# + sn# * sx0#
REM convert to 2000 celestial coords.
theta# = 0 - 23.45# * pi180#
GOSUB 1290
xe(1) = 0 - sx# / mtot#
ye(1) = 0 - sy# / mtot#
ze(1) = 0 - sz# / mtot#
re(1) = SQR(xe(1) ^ 2 + ye(1) ^ 2 + ze(1) ^ 2)
RETURN
1270 z# = r# * SIN(decl#): cs# = COS(decl#)
x# = r# * cs# * COS(ra#): y# = r# * cs# * SIN(ra#)
sx# = sx# + w# * x#
sy# = sy# + w# * y#
sz# = sz# + w# * z#
RETURN
1290 cs# = COS(theta#): sn# = SIN(theta#): sz0# = sz#
sz# = cs# * sz# - sn# * sy#
sy# = cs# * sy# + sn# * sz0#
RETURN
REM find the accel. of the parabolas fitting 1st thru 3rd points
2000 vx1# = (x(2) - x(1)) * dt1#
vy1# = (y(2) - y(1)) * dt1#
vz1# = (z(2) - z(1)) * dt1#
vx2# = (x(3) - x(2)) * dt2#
vy2# = (y(3) - y(2)) * dt2#
vz2# = (z(3) - z(2)) * dt2#
REM accel., corrected to fixed frame
ax# = (vx2# - vx1#) * dt3# * corr#
ay# = (vy2# - vy1#) * dt3# * corr#
az# = (vz2# - vz1#) * dt3# * corr#
RETURN
REM find Newtonian accel. at 1st thru 3rd points
2100 FOR j = 1 TO 3
dr# = kk# / r(j) ^ 3
nax(j) = x(j) * dr#
nay(j) = y(j) * dr#
naz(j) = z(j) * dr#
NEXT j
RETURN
REM weighted ave. Newtonian accel., and work done
2200 avenax# = 0: avenay# = 0: avenaz# = 0
FOR j = 1 TO 3
avenax# = avenax# + w(j) * nax(j)
avenay# = avenay# + w(j) * nay(j)
avenaz# = avenaz# + w(j) * naz(j)
NEXT j
rel# = (avenax# - ax#) ^ 2 + (avenay# - ay#) ^ 2 + (avenaz# - az#) ^ 2
RETURN
REM update best
3000 rel0# = rel#: i0 = i: q0# = q#
r00# = r0#: rp0# = rp#: rt0# = rt#: rpp0# = rpp#: rtt0# = rtt#
FOR m = 1 TO 6
x0(m) = x(m): y0(m) = y(m): z0(m) = z(m): r0(m) = r(m)
NEXT m
RETURN
3100 PRINT "apparent circular period"
PRINT pi2# / om# / 365.25#
RETURN
REM find areas of the four sectors
REM PRINT "radii"
4000 FOR i = 1 TO 4
the(i) = ATN(-y(i) / x(i))
IF x(i) > 0 THEN GOSUB 4100
IF the(i) < 0 THEN GOSUB 4120
r(i) = SQR(x(i) ^ 2 + y(i) ^ 2)
REM PRINT r(i),
NEXT i
REM PRINT
r(5) = r(1): the(5) = the(1) + pi2#
REM PRINT "drdth"
FOR i = 1 TO 4
u# = 2 * co(1) * x(i) + co(2) * y(i) + co(4)
w# = co(2) * x(i) + 2 * co(3) * y(i) + co(5)
sn# = 0 - (u# * y(i) - w# * x(i)) / SQR(u# ^ 2 + w# ^ 2)
drdth(i) = sn# / SQR(1 - sn# ^ 2)
REM PRINT drdth(i),
NEXT i
REM PRINT
drdth(5) = drdth(1)
PRINT "interpolated radii"
REM interpolate r quadratically
FOR i = 1 TO 5
IF i < 5 THEN dthpl# = the(i + 1) - the(i)
IF i < 5 THEN drpl# = r(i + 1) - r(i)
IF i > 1 THEN dthmi# = the(i) - the(i - 1)
IF i > 1 THEN drmi# = r(i) - r(i - 1)
IF i < 5 THEN rpl(i) = r(i) + 2 / 9 * drdth(i) * dthpl# + 1 / 9 * drpl#
IF i > 1 THEN rmi(i) = r(i) - 2 / 9 * drdth(i) * dthmi# - 1 / 9 * drmi#
PRINT rpl(i), rmi(i)
NEXT i
REM Cotes' 3/8 rule
REM PRINT "areas"
FOR i = 1 TO 4
r2m# = r(i) ^ 2 + r(i + 1) ^ 2 + 3 * (rpl(i) ^ 2 + rmi(i + 1) ^ 2)
area(i) = r2m# / 8 * (the(i + 1) - the(i)) / 2
REM PRINT area(i),
NEXT i
REM PRINT
are# = area(1) + area(2) + area(3) + area(4)
s# = 0: ss# = 0
PRINT "rates"
FOR i = 1 TO 3
IF i = 1 THEN lap# = 1
IF i = 2 THEN lap# = 0
IF i = 3 THEN lap# = 1
rate# = (area(i) + lap# * are#) / (t(i + 1) - t(i))
PRINT rate#; " ";
s# = s# + rate#: ss# = ss# + rate# ^ 2
NEXT i
rel# = SQR((ss# - s# ^ 2 / 3) / 2) / (s# / 3)
IF rel# < rel0# THEN GOSUB 4200
RETURN
4100 the(i) = the(i) + pi#
RETURN
4120 the(i) = the(i) + pi2#
RETURN
4200 rel0# = rel#: x00# = x0#
RETURN
REM Astronomical Almanac barycentric Earth XYZ '86
DATA -.992473#,-.99431#,.097062#,.081274#,.042041#,.035195#
REM '87
DATA -.647098#,-.660316#,.689355#,.678874#,.298895#,.29435#
REM '97
DATA -.951971#,-.957296#,.278791#,.263632#,.121036#,.114464#
REM '07
DATA -.992917#,-.991887#,-.057633#,-.073429#,-.025095#,-.031943#
REM '08
DATA -.008637#,-.026129#,.906672#,.906369#,.39304#,.392909#
REM '09
DATA -.493191#,-.508316#,.786488#,.778431#,.340941#,.337449#
REM (times of plates, Sun positions,
REM coords. of identified disappearing dots)
REM plate A (Red)
REM DATA 1954.154,.9900229, 22,34, 27.6,-8,-59,-7
REM object A2
REM DATA 11,2,25.47,-5,-56,-12
REM object A
REM DATA 11,3,12.4,-5,-58,-9
REM plate B (Red)
REM DATA 1986.201,.9945905,23,40,30.9,-2,-10,-49
REM object B3
REM DATA 11,16,51.67,-7,-49,-40.4
REM object B
REM DATA 11,16,55.78,-7,-55,-14
REM plate C (Red)
REM DATA 1987.08215,.9852006,20,53,6.9,-17,31,13
REM object C11 (named for easy recall, like "CII")
REM DATA 11,18,.41,-8,-1,-57.7
REM object C
REM DATA 11,18,3.18,-7,-58,-46.1
REM plate D (optical infrared)
REM DATA 1997.16711,.9914845,22,57,21.1,-6,-39,-54
REM DATA 11,22,32.9,-8,-26,-56
REM DATA 11,22,16.77,-8,-29,-32.6
REM DATA 11,22,1.33,-8,-34,-36.9
Please Log in or Create an account to join the conversation.
- Joe Keller
- Offline
- Platinum Member
Less
More
- Thank you received: 0
15 years 9 months ago #15785
by Joe Keller
Replied by Joe Keller on topic Reply from
Hi *********, and Prof. *********!
Here's a guide star (too dim to be a real guide star)(at present, Jan. 30, Frey is about 6 arcminutes west of it, mainly due to change in Earth parallax since Dec. 22, and Barbarossa about 3.5 arcmin west of Frey; so, it would be best to ask those astronomers in Boston to center the photo (2*6+3.5)/2 = 7.75 arcminutes west of this star; i.e., 31 "seconds of Right Ascension" west of it) and 5.5*0.4 arcminutes north of it:
USNO-B 0806-0230017
USNO-B Red1 mag +19.61, USNO-B Red2 mag +19.21
Star's coordinates:
RA 11:27:32.00
Decl -9:21:48.3
Center of desired photo:
RA 11:27:01
decl -9:20
(Subtract about one second of RA, and add about 6" of Decl., that is, aim farther west and north, for each day past Jan. 30.)
I happen to know about this star, because its pixel appearance is very similar to a nearby "disappearing dot" (not in 1987 Red sky survey) not too near the west edge of the Dec. 22, 2008 U. of Iowa photo, at
RA 11:27:30.17
Decl -9:21:48.6
This object not only is very similar in appearance to the above Red mag +19.4 star, it also accords perfectly with my revised, corrected and accurized theory of the solar and binary Barbarossa/Frey orbits. This theory is completely contained within the BASIC program I posted to Dr. Van Flandern's messageboard about an hour ago.
There is now so little doubt in my mind of the reality of these objects, that I'm hardly in a rush to get photos. "They'll keep."
Your *******,
Joe
Here's a guide star (too dim to be a real guide star)(at present, Jan. 30, Frey is about 6 arcminutes west of it, mainly due to change in Earth parallax since Dec. 22, and Barbarossa about 3.5 arcmin west of Frey; so, it would be best to ask those astronomers in Boston to center the photo (2*6+3.5)/2 = 7.75 arcminutes west of this star; i.e., 31 "seconds of Right Ascension" west of it) and 5.5*0.4 arcminutes north of it:
USNO-B 0806-0230017
USNO-B Red1 mag +19.61, USNO-B Red2 mag +19.21
Star's coordinates:
RA 11:27:32.00
Decl -9:21:48.3
Center of desired photo:
RA 11:27:01
decl -9:20
(Subtract about one second of RA, and add about 6" of Decl., that is, aim farther west and north, for each day past Jan. 30.)
I happen to know about this star, because its pixel appearance is very similar to a nearby "disappearing dot" (not in 1987 Red sky survey) not too near the west edge of the Dec. 22, 2008 U. of Iowa photo, at
RA 11:27:30.17
Decl -9:21:48.6
This object not only is very similar in appearance to the above Red mag +19.4 star, it also accords perfectly with my revised, corrected and accurized theory of the solar and binary Barbarossa/Frey orbits. This theory is completely contained within the BASIC program I posted to Dr. Van Flandern's messageboard about an hour ago.
There is now so little doubt in my mind of the reality of these objects, that I'm hardly in a rush to get photos. "They'll keep."
Your *******,
Joe
Please Log in or Create an account to join the conversation.
Time to create page: 0.376 seconds