七政四余三王的算法

目录

1，儒略数d的计算 2，行星轨道要素。

3， 根据儒略历计算太阳的平黄经。 4， 根据儒略历年计算月球的平黄经。 5， 天球坐标转换。 6, 恒星时和地方时间 7一个精确计算的例子

正文

1，儒略数d的计算

YYYY年MM月DD日TT时

d = 367×Y - (7×(Y + ((M+9)/12)))/4 + (275×M)/9 + D+TT/24 - 730530

2，行星轨道要素：

2－1）参数介绍 a （半长轴） e （偏心率）

I（轨道面和黄道交角），0N （升交点黄经）（大写Ω） 从春分点到升交点的角度。

w （升交点到近日点的角度--近日点黄经）ω) 升交点到近日点角度。 M（平近点角）

L（平均黄经度） = M +w+ N E（偏近点角） M = E - e × sin(E) V（真近点角）

r（向径）

2－2）太阳系主要天体的轨道要素

太阳轨道要素：

w = 282.9404_deg + 4.70935E-5_deg * d

a = 1.000000 e = 0.016709 - 1.151E-9 * d M = 356.0470_deg + 0.9856002585_deg * d

太阴轨道要素：

N = 125.1228_deg - 0.0529538083_deg * d i = 5.14_deg w = 318.0634_deg + 0.13573223_deg * d ( a = 60.2666

e = 0.0900

水星轨道要素：

N = 48.3313_deg + 3.24587E-5_deg * d i = 7.0047_deg + 5.00E-8_deg * d w = 29.1241_deg + 1.01444E-5_deg * d a = 0.387098 e = 0.205635 + 5.59E-10 * d M = 168.6562_deg + 4.0923344368_deg * d

金星轨道要素：

N = 76.6799_deg + 2.46590E-5_deg * d i = 3.3946_deg + 2.75E-8_deg * d w = .10_deg + 1.38374E-5_deg * d a = 0.723330

e = 0.006773 - 1.302E-9 * d M = 48.0052_deg + 1.6021302244_deg * d

火星轨道要素：

N = 49.5574_deg + 2.11081E-5_deg * d i = 1.8497_deg - 1.78E-8_deg * d w = 286.5016_deg + 2.92961E-5_deg * d a = 1.523688

e = 0.093405 + 2.516E-9 * d M = 18.6021_deg + 0.5240207766_deg * d

木星轨道要素：

N = 100.42_deg + 2.768E-5_deg * d i = 1.3030_deg - 1.557E-7_deg * d w = 273.8777_deg + 1.505E-5_deg * d a = 5.20256

e = 0.048498 + 4.469E-9 * d M = 19.50_deg + 0.0830853001_deg * d

土星轨道要素：

N = 113.6634_deg + 2.380E-5_deg * d i = 2.4886_deg - 1.081E-7_deg * d w = 339.3939_deg + 2.97661E-5_deg * d a = 9.575

e = 0.0556 - 9.499E-9 * d M = 316.9670_deg + 0.0334442282_deg * d

3，根据公元纪年计算太阳的平黄经及黄赤交角

即黄经的变化量与儒略数d的关系 PI = 3.141592653579323846

l w = 282.9404_deg （度）+ 4.70935E-5_deg × d a = 1.000000 e = 0.016709 - 1.151E-9 × d l M = 356.0470_deg + 0.9856002585_deg × d

l 黄赤交角 oblecl = 23.4393_deg - 3.563E-7_deg × d l L = M +w+ N

N （升交点黄经）（大写Ω） 从春分点到升交点的角度.因为对太阳来说，升交点就是春分点，故N=0 L=M+w

2005年6月17日23：00UT（UT：格林威治时间）

d=1995.35, 取整， d=1995

w = 282.9404 + 4.70935E-5 × d=283.03437 M = 356.0470 + 0.9856002585 × d =2322.6686

即:M=162.6686

太阳的平黄经L=M+ w=445.7029度 L=85.7029

岁差：3.82394E-5 × ( 365.2422 × ( 公元纪年 - 2000.0 ) - d )=0.2 L=85.7029+0.2=85.9

黄赤交角oblecl = 23.4393 - 3.563E-7 × d=23.4393度

4,根据公元纪年计算月球的平黄经及轨道角度。

月球的计算较复杂，因其轨道与黄道有一个交角(5.1415度),精确计算还要考虑太阳木星金星的摄动，其黄经黄纬的校正项多达100多项.限于篇幅此处就不介绍了，只讲求近似值。 轨道要素:

N = 125.1228 - 0.0529538083 × d (升交点经度)

i = 5.14度 (轨道倾角)

w = 318.0634 + 0.13573223 × d (Arg. of perigee(近地点角) a = 60.2666 (Mean distance) e = 0.0900 (Eccentricity) M = 115.36 + 13.09929509 × d (平近点角) l 试计算1990年4月19日0:00的月球黄经

4－1）,求得儒略数为: d=-33,代入轨道要素得:

N = 312.7381 i = 5.14

w = -2.26 修正为95.74 a = 60.2666 e = 0.0900

M = -46173.9046,用360×129加上修正为266.09

4－2）,根据开普勒第二(面积)定理推得

E=M+e×sin(E) 用迭代法求解此超越方程 先求E0, PI=3.1415926

E0≈M + (180/pi) × e × sin(M) × (1 + e × cos(M)) 度

E1 = E0 - (E0 - (180/pi) × e × sin(E0) - M) / (1 - e × cos(E0))

En+1=En-(En-(180/pi) × e × sin( En)-M)/(1- e × cos (En)) 直到En+1和En的差非常小为止 代入数值,得到:E0≈262.96度 迭代后解得E=262.9735

4－3），求轨道面上的直角坐标

x = r × cos(v) = a × (cos(E) - e)

y = r × sin(v) = a × sqrt(1 - e×e) × sin(E) 得到 x = -10.68095 y = -59.72377

4－4），求真近点角v及向径r

r = sqrt( x×x + y×y ) = 60.67134 地球半径

V== 180/3.1415926*atan(-59.72377/-10.68095)-180=-100.1395 加上360度的倍数,使之大于0,得到 V==259.8605_deg

4－5）,化为黄道直角坐标:

公式:

xeclip = r * ( cos(N) * cos(v+w) - sin(N) * sin(v+w) * cos(i) ) yeclip = r * ( sin(N) * cos(v+w) + cos(N) * sin(v+w) * cos(i) ) zeclip = r * sin(v+w) * sin(i) 代入数值,注意度要化弧度, N= 312.7381

r=60.67134 ,v=259.8605*PI/180 w=95.74*PI/180 I=5.14*PI/180

得到日心三维直角坐标：

xeclip = +37.65311 yeclip = -47.57180 zeclip = -0.41687

4－6）,化为黄道坐标系

6.1公式:

向径:r =sqrt((x*x+y*y+z*z) 黄经L: RA=atan(y/x)

黄纬Lat : Decl=atan(z/sqrt(x*x+y*y))

代入数值,注意求反正切时要用第4中判断角的象限:方法

得到:

ra = 308.3616_deg lat = -0.3937_deg r = 60.6713

7，恒星时与时角、地方时

SIDTIME = GMST0 + UT + LON/15

地方恒星时=格林尼治00：00恒星时+世界时(格林尼治时间)+经度（东经为正，西经为负）/15

GMST0（格林尼治0点恒星时）= ( L + 180_deg ) / 15 = L/15 + 12h

L：太阳平黄经， L=26.8388 : 经度为东经15度的地方恒星时为:

GMST0 = 26.8388_deg/15 + 12h = 13.725 hours

SIDTIME = GMST0 + UT + LON/15 = 13.725h + 0 + 15_deg/15 = 14.725 hours

SIDTIME = 14h 47m 21.3s

计算地平坐标需要知道时角HA:

HA = SIDTIME - RA

HA = 14.725h - 1.77720h = 13.01205h = 195.1808_deg

坐标转换，求出七政黄道后可以做出七政静盘，弧角星盘的算法，要用到坐标变换，公式如下，有心人可以自己研究推算，另有黄道十二宫的详细算法，因为笔者在其他网站已经发表过，此处就不写了

地平坐标系: 设天体P点的高度为a，是从天文地平沿通过P点和天文天顶的大圆向天文天顶(或天底)计量的角距离，向天顶为正(向天底为负)， P点和天文天顶(天底)重合，则a＝十90。(一90。)．高度的余数称为天项距，北天极的高度是观测者的天文纬度φ． P点的方位角A，是从南点向西沿着天文地平到通过P 点和z点的大圆与天文地平的交点所量度的角距离．

δ赤纬，h时角，设α为赤经，h=τ-α, τ恒星时，φ天文纬度，A是地平方位角 β黄纬，λ黄经，ε黄赤交角。a是地平高度 7．1，从地平坐标到赤道坐标

cosδ*cosh=cosφ*sina+sinφ*cosa*cosA cosδ*sinh= cosa*sinA

sinδ=sinφ*sina-cosφcosa*cosA 7．2从时角和赤纬坐标系到地平坐标系 cosa*cosA= sinφ*cosδ*cos h-cosφsinδ cosa*sinA= cosδ*sinh

sina=sinφ*sinδ+cosφcosδ*cosh 7．3 从黄道坐标系到赤道坐标系 cosδcosα=cosβcosλ

cosδsinα= cosεcosβsinλ-sinεsinβ sinδ=cosεsinβ+sinεcosβsinλ 7. 4 从赤道坐标系到黄道坐标系 cosβcosλ=cosδcosα

cosβsinλ=sinεsinδ+cosεcosδsinα sinβ=cosεsinδ-sinεcosδsinα

举一个精确计算的例子：

昆明市七政四余三王真躔

经度 纬度 102.7 25 年 月 日 时 分

北京时间 2007 4 28 12 0 2007 4 28 10 50.8

高 度 黄 经 赤 经 赤 纬 到 方

日 124.966 72.23571 37.8379 2.365037 14.12268044 月 261.4755 -.7374 241.0451 15.85572 -25.0010048 金 76.45095 37.21378 78.75123 5.171086 25.1310498 木 259.7284 -35.9382 259.0905 17.21336 -22.26662 水 140.4888 72.52965 31.81744 2.00087 11.14019674 火 227.2013 46.76688 346.78 23.2265 -6.53865875 土 59.57946 -20.1115 138.2662 9.408092 16.59498244 计都 224.7706 53.2332 353.2759 23.58842 -2.66942139 罗睺 44.77057 -53.2332 173.2759 11.58842 2.669421395 紫炁 218.2315 59.58156 0.575021 0.035172 0.228718286 月孛 38.23147 -59.5816 180.575 12.03517 -0.22871829 天王 227.6742 47.69747 347.42 23.25243 -5.63336956 海王 239.5127 23.85027 321.8535 21.6204 -14.4575283 冥王 256.2059 -29.5288 266.1162 17.71795 -23.3481804

2007年4月28日12：00，昆明的十二宫黄道如下 宫数 黄道

1/7宫 117.461917 2/8宫 142.7933 3/9宫 169.520136 4/10宫 200.996901 5/11宫 236.405226 6/12宫 269.524383 7/1宫 297.461917 8/2宫 322.7933 9/3宫 349.520136 10/4宫 20.996901 11/5宫 56.405226 12/6宫 .524383