中华农历论坛知识讨论区历法知识 → 天文算法讨论


  共有209505人关注过本帖树形打印

主题:天文算法讨论

帅哥哟,离线,有人找我吗?
xjw01
  11楼 个性首页 | 信息 | 搜索 | 邮箱 | 主页 | UC


加好友 发短信
等级:蝙蝠侠 帖子:866 积分:3804 威望:3 精华:20 注册:2008/3/20 22:14:00
  发帖心情 Post By:2008/5/9 22:09:00


图片点击可在新窗口打开查看此主题相关图片如下:
图片点击可在新窗口打开查看

第(五)帖
  历元、摄动、自球自转轴、天北极、黄极、黄轴、黄经、黄纬、赤经、赤纬,这几个概念不再讲述,请在网络上找一些资料看看吧!
  春风点是黄道与赤道的交点!
  春风点是天球中的重要的参考点,遗憾的是,春风点会移动!
  在图中,有两条黄道、两道赤道。
  一条的历元2000.0时刻黄道,一条是当日(Date)黄道。
  一条的历元2000.0时刻赤道,一条是当日(Date)赤道。
  随时间推移,黄道为什么会发生变化?因为地球受至大行星的摄动。
  随时间推移,赤道也会发生移动,因为受月球影响,地球自转轴的方向发生移动,造成天北极围绕黄轴缓慢移动。找一个直铁棒,铁棒放置的方向与黄轴相同,用电焊机把铁棒焊在赤道面的中心。始终保持铁棒与黄轴平行,并转动铁棒,这时你将看到,北极围绕铁棒做圆周运动,同时春风点也同步移动,当北极运行一周(约26000年),春风点在黄道上也移动一周。
  春风点移动的后果:黄经、赤经的起算点是春风点,造成所有星体的黄经、赤经发生变化。假设,黄道不动(没有发生变化),赤道发生变化,引起春风点移动,如果春风点在黄道上向西移动了1度(在上图中,往左下方移动),那么,所有星体的黄经增加了1度!
  春风点的移动量称为岁差。岁差有很多种表达方法,最著名、最传统的岁差是指黄道上的岁差。
-----------------------
  上图中,有三个重要的交点:γo,γ1,γm,我们统称它为分点。γo称为“历元2000.0标准分点”,其实就是历元2000时刻的春风点,γm称为“Date分点”或“当日分点”或“当日春分点”。
  图中,两个黄道的交点是N,定义N到γo与N到γo'的距离相等,那么以γo'起算的黄经则与γo起算的黄经是相等的(除非天体过份靠近黄极),对于行星,它在黄道附近,对于大部分恒星,也是远离黄极的,所以γo'对于行星、月球、大部分恒星来说,是很重要的一个参考点。比如太阳,当它从γo出发,转动数周后达到了γo',所经历的时间是整倍数的太阳绕地球运转的周期。这种运行周期具有比较严格的力学意义。恒星在天空中的位置几乎是不变的,以它作为参考系,具有“惯性”系的特征,γo则是利用恒星确定的,也具有惯性特点。任意时刻,我们利用某一颗或几颗恒星的位置,轻松的标定出γo的位置(当然,首次建立FK5系统时,需要大量恒星标定出γo,以提高精度,所以工程巨大),当太阳运行到γo',可以认为太阳相对于惯性参考点γo移动了一周或数个1周,即相对于惯性参考点移动了数个周期,这就是力学意义上的周期。γo与γo'相差很小,所以太阳到了γo与太阳到了γo',太阳在星空的中相对位置几乎相同,因此,太阳经过这两点所经历的时间是整位数所恒星年(以恒星为参考,太阳在天球中运转一周的时间),由于这个原因,恒星年一般认为就是太阳(或地球)在惯性系中的运动周期,“在惯性系中的运动周期”常称为真周期。

   p,γo'到γm的角距离就是Date黄道上的岁差
  ψ,γo到γ1的角距离是J2000黄道上的岁差
  χ,γ1到γm是Date赤道上的岁差
  εo,是历元2000.0的黄赤交角
  ε, 是当日的黄赤交角
  ω, Date平赤道与J2000黄道的夹角
  π,是两个黄道之前的夹角
  П,是N到γo的角距离,它等于N到γo'的角距离
另一组常用的岁差计算参数
  90°-ζ,γo到Q的角距离,是一个常用的岁差计算参数
  90°+z, γm到Q的角距离,是一个常用的岁差计算参数
  θ是两个赤道面的夹角

  也许你会问,怎么会有这么多岁差参数,有必要吗?回答,为了坐标变换方便而已。通常,只需只道3个岁数参数,就可利用坐标变换方法及多项式插值法导出其它的岁差参数。
--------------------------
一组常用的岁差表达,取自IAU2000
  P:(     0.0,        41.99604,  19.39715, -0.22350, -0.01035,  0.00019,  0.0,     0.0)
  Q:(     0.0,      -468.09550,   5.10421,  0.52228, -0.00569, -0.00014,  0.00001, 0.0)
 π:(     0.0,       469.97560,  -3.35050, -0.12370,  0.00030,  0.0,      0.0,     0.0)
 II:(629543.988,   -8679.218,    15.342,    0.005,   -0.037,   -0.001,    0.0,     0.0)
  p:(     0.0,     50287.92262, 111.24406,  0.07699, -0.23479, -0.00178,  0.00018, 0.00001)
 θ:(     0.0,     20041.90936, -42.66980,-41.82364, -0.07291, -0.01127,  0.00036, 0.00009)
 ζ:(     2.72767, 23060.80472,  30.23262, 18.01752, -0.05708, -0.03040, -0.00013, 0.0)
 z :(    -2.72767, 23060.76070, 109.56768, 18.26676, -0.28276, -0.02486, -0.00005, 0.0)
 ε:( 84381.40880,  -468.36051,  -0.01667,  1.99911, -0.00523, -0.00248, -0.00003, 0.0)
 ω:( 84381.40880,    -0.26501,   5.12769, -7.72723, -0.00492,  0.03329, -0.00031,-0.00006)
 ψ:(     0.0,     50384.78750,-107.19530, -1.14366,  1.32832, -0.00940, -0.00350, 0.00017)
 χ:(     0.0,       105.57686,-238.13769, -1.21258,  1.70238, -0.00770, -0.00399, 0.00016)
  上表是关于t的多项式系数,t等于J2000起算的儒略千年数,儒略千年=364250天
  第一列对应t的0次方,第二列对应t的1次方,第三列对应t的2次方,其余类推。
  如π=0.0*t^0 + 469.97560*t^1 - 3.35050*t^2 + ...
  结果的单位是角秒。

-------------------------
  以上岁差表达是非常精确的,误差小于1毫角秒。
  不过t不能太大,t应小于1。
  如果t=10,那么,会有什么后果?比如计算θ,t的7次方项为0.00001*10000000=100角秒!实际上,表达式系数的舍入误本身就是0.00001,所以误差高达100角秒。
  如果由于某种需要,你非要计算t=10的情况,最好把t的5、6、7次方项忽略。这样,如果t=5,误差也可控制在1角秒之内。


支持(0中立(0反对(0单帖管理 | 引用 | 回复 回到顶部
帅哥哟,离线,有人找我吗?
xjw01
  12楼 个性首页 | 信息 | 搜索 | 邮箱 | 主页 | UC


加好友 发短信
等级:蝙蝠侠 帖子:866 积分:3804 威望:3 精华:20 注册:2008/3/20 22:14:00
  发帖心情 Post By:2008/5/9 22:45:00

第(六)帖
  顺便说一下,儒略千年数是一种时间计数方法。日常生活中,时间的计数法使用“年月日时分秒”,在天文计算过程,这种时间表达是很不方便的,原因是明显的,需要6个数才能表达一个时间。
  天文学中,采用2000年1月1日12时0分0秒时作为计时的起点,并连续计时,单位是天。1天日86400秒。1秒是多少?1秒时间长度由原子钟得到,现代天文学的时间基础已不再使用星体运动位置来确定,改用原子钟,原子钟的精度要好得多。
  为什么又称为“儒略”千年数。这还要谈到另一个问题,一年有多少天?你可能立刻回答到:365.2422日。然而,在现代天文学家眼里,一年的长度是变化的,并不是固定为365.2422日,何况地球运动时还受到各种干扰,365.2422只是一定时期内的平均值。因此,任意一年某时刻加上365.2422后,并不是下一年的同一时刻,这样,把365.2422看做1个单位意义不大,天文学家更原意把365.25看做1个单位,它与公历的置闰规则具有更好的“兼容”性,便于日历转换。365.25天是儒略年的长度。儒略1千年就是365250天。
  再次强调,1天=86400原子秒,与地球自转没有必然的关系。在历史上,曾使用某一年的地球公转与自转来定义秒长,之后,地球自转速度发生改变,1天不等于86400秒了。注意,此“天”非彼“天”,前面的天是时间单位,后面的天指平均自转一周所用的时间长度。然而,人类的日常生活与地球自转息息相关,日期与时间的计数必须尽可能与地球自转同步,这就造成,我们有两套基本的时间系统,一个是手表时(或称为协调世界时或称为UTC时),一个是力学时。力学时又时什么东西呢?力学时具有鲜明的动力学色彩,假如我们使用2000年地球公转一圈所用时间定义为1个单位,2万年后,公转周期可能发生了变化,如果不通过力学原理推算,我们将无法得到那时候一年的长度,但是,通过力学方法,我们可以得到2万年后公转的周期,比如是1.00001(乱写一个)个单位。简单的说,当我们定义了时间单位,我们就可以使用力学原理推算天体在某一位置所对应时间,这就够成时间基础(根据星体位置确定时间)。原子时就简单多了,它不是通过天体位置来确定时间的,而是通过原子钟内部的振荡次数来来确定时间。原子钟曾与力学时对准过,秒长是相同的,起点相差零点几秒钟。
  理论上,力学时、原子时是非常均匀的。手表时与地球自转(以恒星为参考点确定自转速度)同步,而地球自转时快时慢,所以不均匀。
  手表时(去除地区的时差后的手表时,即UT1)与原子时的差称为ΔT,当然手表时含有时区信息,如北京时间多了正好8个时。ΔT可由一组多项式表达式计算得到。
  手表时的秒长是原子秒,为什么手表时不是原子时呢?因为国际上有人动了手脚,在每年的最后一天做了“正或负跳秒(正负闰秒)”处理,实现手表时与地球自转同步。

//原子时与UT1的差
//部分参考了美国国家航空航天局网站里的算法,www.nasa.gov
//2050以后的外推参考了skymap10.5,在skymap的星历表中分段找几个数据,再用最小二乘法拟合
double deltatT(double y){ //输入年份(可带小数点),如2000.78
 static double DTxs[15][12]={//{范围下,范围上,偏移,时间单位,多项各系数}
   {-9999,-500,1820,100,-20,      0,        32},
   {-500, 500, 0,   100, 10583.6,-1014.41,  33.78311, -5.952053,  -0.1798452,   0.022174192, 0.0090316521},
   {500,  1600,1000,100, 1574.2, -556.01,   71.23472,  0.319781,  -0.8503463,  -0.005050998, 0.0083572073},
   {1600, 1700,1600,1,   120,    -0.9808,  -0.015320,  1.0/7129},
   {1700, 1800,1700,1,   8.83,    0.1603,  -0.0059285, 0.00013336,-1.0/1174000},
   {1800, 1860,1800,1,   13.72,  -0.332447, 0.0068612, 0.0041116, -0.00037436,  1.21272e-5, -1.699e-7, 8.75e-10},
   {1860, 1900,1860,1,   7.62,    0.5737,  -0.2517540, 0.01680668,-0.0004473624,1.0/233174},
   {1900, 1920,1900,1,  -2.79,    1.494119,-0.0598939, 0.0061966, -0.000197},
   {1920, 1941,1920,1,   21.20,   0.84493, -0.0761000, 0.0020936},
   {1941, 1961,1950,1,   29.07,   0.407,   -1.0/233,   1.0/2547},
   {1961, 1986,1975,1,   45.45,   1.067,   -1.0/260,  -1.0/718},
   {1986, 2005,2000,1,   63.86,   0.3345,  -0.060374,  0.0017275,  0.000651814, 0.00002373599},
   {2005, 2150,2000,100, 64.723,-11.050,    232.182,  -86.08},
   {2150, 9999,2000,100, 59.033, 103.4335,  28.980314, 0.000276}
 };
 //计算TD-UT,力学时与世界时之差
 double *p,dt=0,u;
 int i,j;
 for(i=0;i<15;i++){
   p=DTxs[i];
   if(y<p[0]||y>p[1]) continue;
   y=(y-p[2])/p[3];
   for(j=4,u=1;j<12;j++) dt+=u*p[j],u*=y;
   return dt;
 }
 return 0;
}


支持(0中立(0反对(0单帖管理 | 引用 | 回复 回到顶部
帅哥哟,离线,有人找我吗?
xjw01
  13楼 个性首页 | 信息 | 搜索 | 邮箱 | 主页 | UC


加好友 发短信
等级:蝙蝠侠 帖子:866 积分:3804 威望:3 精华:20 注册:2008/3/20 22:14:00
  发帖心情 Post By:2008/5/10 9:58:00

  第(七)帖 关于岁差与章动
  [岁差]:
  我们在第(五)贴已经谈到了赤道与黄道的变化将引起岁差。受到月球运动的影响,地球自转轴发生围绕黄轴的缓慢的长期的圆周运动,造成赤道发生变化,进而引起春风点变化。从上面的那个岁差示意图中不难发现,春风点的变化不单单是赤道变化引起的,黄道的变化也会引起春风点的变化。不过,图中的角π非常小,所以黄道变化引起和春风点移动是非常小的。
  从岁差示道意图中的几何关系:
    弧Nγ1 = 弧Nγm + χcos(ω),注意,这是近似计算,但精度仍非常高
  两边同时减去П得:
    Nγ1 - П = (Nγm - П)+ χcos(ω),即 ψ = p+χcos(ω)
  其实ψ就是:黄道不动,赤道变化引起的黄经岁差,我们称之是日月黄经岁差。
  其实χ就是:Date赤道上的行星岁差。
  如果把ψ = p+χcos(ω)改写为p = ψ-χcos(ω)
  那么 -χcos(ω) 就是行星黄经岁差。
  同样,从几何关系,还可轻松得到行星黄纬岁差:χsin(ω)
----------------
  χ的岁差速度是:106"/千年,加速度为-238"/平方千年
  χ岁差造成春风点东移:0".106*cos(23.5°)/年
  日月岁差(黄经)造成春风点西移:50.384.78750/年
  黄经总岁差速度约为二者之差:50".29/年
----------------
  [章动]:
  受到日月运动的影响,地球自转轴发生振动,这种振动称为章动。“章”字,从何讲起?地球自转的振周期主要与月球轨道升交点运动有关,升交点运动周期为18.6年,古人称之为一章。
  在第(五)贴的图中,已给出受到章动后的赤道位置。图中的Δψ就是黄经章动值。
  章动的后果:春风点在黄道上移动了Δψ(称为黄经章动),黄赤交角变化了Δε(称为交角章动)
  考虑章动以后的赤道称为真赤道,否则为平赤道。
  黄道的振动是非常微弱,所以一般不再区分“真”与“平”
  要得到真春点起算的黄经,必须加上Δψ,要继续旋转到真赤道坐标中,应使用真黄赤交角(即加上Δε)
-----------------
  Δψ与Δε的计算
  前面已发了一个关于IAU1980与IAU2000B章动序列比较的贴子,里面含有相关的所有计算,请阅读程序。
----------------
  精度控制
  完整的计算,计算量偏大,程序冗长。如果你只是为了计算农历,根本不需要那么高的精度。
  为什么呢?
  因为,现代的紫金山中国农历是基于北京手表时的,这种时间与地球自转有关,考虑到地球地转的不规律性,无法对几百年后的北京时间进行测算,误差0.5分钟是正常的,如果对5千年后的北京时间预测计算,误差可达几小时。0.5分钟的误差意味差什么呢?意味着太阳位置误差可达1".2,月亮位置误差可达15",因此没有必要把Δψ算得很准,0".1已经足够。
  以下javascript程序给出几千年内精度为0".05的黄经章动及交角章动。
<script language=javascript>
/*********************
 对IAU1980章动计算的简化处理
 许剑伟 2008年5月2日
*********************/
t=1; //J2000.0起算的儒略世纪数
var nuB=new Array(//章动表
  2.1824,   -33.75705, 36e-6,-171996,-1742,92025, 89,
  3.5069,  1256.66393, 11e-6, -13187,  -16, 5736,-31,
  1.3375, 16799.41822,-51e-6,  -2274,   -2,  977, -5,
  4.3649,   -67.51409, 72e-6,   2062,    2, -895,  5,
  0.0431,  -628.30196,  3e-6,  -1426,   34,   54, -1,
  2.3556,  8328.69143,155e-6,    712,    1,   -7,  0,
  3.4638,  1884.96589,  8e-6,   -517,   12,  224, -6,
  5.4383, 16833.17527,-87e-6,   -386,   -4,  200,  0,
  3.6931, 25128.10965,103e-6,   -301,    0,  129, -1,
  3.5501,   628.36198, 13e-6,    217,   -5,  -95,  3
);
var dL3=0,dE3=0;
for(i=0;i<nuB.length;i+=7){
  c=nuB[i] +nuB[i+1]*t +nuB[i+2]*t*t;
  dL3+=(nuB[i+3]+nuB[i+4]*t/10)*Math.sin(c); //黄经章动
  dE3+=(nuB[i+5]+nuB[i+6]*t/10)*Math.cos(c); //交角章动
}
dL3/=10000; //黄经章动
dE3/=10000; //交角章动
document.write("[精简的章动计算,精度0.05角秒,单位角秒]<br>");
document.write("黄经章动:" + dL3  + " 交角章动:" + dE3  + "<br>");

</script>

---------------------
  [问题]:程序中使用了“J2000.0起算的儒略世纪数”,你会计算它吗?
  [例题]:比如,2002年1月1日12:00:00 TD的儒略世纪数T是多少?TD表示力学时。
  [解]:
   J2000.0对应2000年1月1日12:00:00 TD
   由于2000年是闰年,含有366天,也就是说:
     2001年1月1日12:00:00加上366天后变为2001年1月1日12:00:00
   又因为2001年是平年,所以再加上365天就是2002年1月1日12:00:00
   所以,所求时刻相对于J2000.0相差366+365=731天
   那么,T = 731/36525 = 0.020013689

 


支持(0中立(0反对(0单帖管理 | 引用 | 回复 回到顶部
帅哥哟,离线,有人找我吗?
xjw01
  14楼 个性首页 | 信息 | 搜索 | 邮箱 | 主页 | UC


加好友 发短信
等级:蝙蝠侠 帖子:866 积分:3804 威望:3 精华:20 注册:2008/3/20 22:14:00
  发帖心情 Post By:2008/5/11 20:20:00

  第5帖中的“儒略千年=364250天”应改为儒略千年=365250天,“是整位数所恒星年”应改为“是整倍数的恒星年”

  使用VSOP行星理论及MPP02月球理论计算星历可以得到很高的精度,但是,多数情况下,我们并不需要高精度,比如业余天文观测、农历计算。
  注意,试图把未来几千年的农历精度提高到“1秒”,人类目前的知识是不可能的做到这一点的。
  为此,本文阐述VSOP87及MPP02的截断处理。

 下载信息  [文件大小:   下载次数: ]
图片点击可在新窗口打开查看点击浏览该文件:


 

[此贴子已经被作者于2008-5-12 7:23:01编辑过]

支持(0中立(0反对(0单帖管理 | 引用 | 回复 回到顶部
帅哥哟,离线,有人找我吗?
春光
  15楼 个性首页 | 信息 | 搜索 | 邮箱 | 主页 | UC


加好友 发短信
等级:版主 帖子:917 积分:7567 威望:10 精华:1 注册:2005/3/17 11:37:00
  发帖心情 Post By:2008/5/13 8:26:00

    本人也发一个关于任意经纬度的太阳,月亮升落的计算程序,摘自国外的网站,供大家研读。

程序作者说他的程序参考了 Montenbruck and Pfleger's “Astonomy on the Computer 2nd english ed” 的第3.8章。

可惜的是“Astonomy on the Computer 2nd english ed” 这本书我还没有找到,没有看到原文。

对于生活在地球上的人来说,一日有几个重要时刻:真夜半(太阳下中天)时刻,日出时刻,真正午(太阳上中天)时刻,日没时刻。

     还有民用晨光始时刻,航海晨光始时刻,天文晨光始时刻,民用昏影终时刻,航海昏影终时刻,天文昏影终时刻。

     月亮下中天时刻(和潮汐有关),月出时刻,月亮上中天时刻(和潮汐有关),月没时刻。

这个程序的精度由于天文资料较少,所以我没有测试,如果计算太阳和月亮上中天和下中天时刻还得需要增加一些程序段。

 

 下载信息  [文件大小:   下载次数: ]
图片点击可在新窗口打开查看点击浏览该文件:


支持(0中立(0反对(0单帖管理 | 引用 | 回复 回到顶部
帅哥哟,离线,有人找我吗?
xjw01
  16楼 个性首页 | 信息 | 搜索 | 邮箱 | 主页 | UC


加好友 发短信
等级:蝙蝠侠 帖子:866 积分:3804 威望:3 精华:20 注册:2008/3/20 22:14:00
  发帖心情 Post By:2008/5/13 9:09:00

好东西啊,有参考价值。

程序中真位置的计算精度较底,但这不影响“升、中天、落”的计算精度。

具体原因请参考《天文算法》中的“关于精度”一章的描述。


支持(0中立(0反对(0单帖管理 | 引用 | 回复 回到顶部
帅哥哟,离线,有人找我吗?
xjw01
  17楼 个性首页 | 信息 | 搜索 | 邮箱 | 主页 | UC


加好友 发短信
等级:蝙蝠侠 帖子:866 积分:3804 威望:3 精华:20 注册:2008/3/20 22:14:00
  发帖心情 Post By:2008/5/13 9:15:00

行星光行差

剑 2008-5-13日

行星光行差计算起来比较麻烦,也没有现成的简单公式可以应用,我们总是根据具体的行星问题展开计算。用下文所述的方法得到的光行差是非常准确的,最后得到的行星视位置与skymap的计算结果分毫不差。当然,日月计算也会涉及行星光行差问题,但这只以行星光行差中的最简单的一个特例。

 
图片点击可在新窗口打开查看此主题相关图片如下:
图片点击可在新窗口打开查看

t时刻,地球在A点,天体在B点。此时,我们在地球上看到的星光并不是并非t时刻天体发出来的光。由于光速是有限值,星光传到地球需要一定的时间τ,所以t时刻我们看到的星光是天体在t-τ时刻发出来的。在图中,O和S分别表示t-τ时刻地球和天体的直位置。我们称τ为光行时。

可以把τ看作一个时间单位,图中矢量的既可表示速度大小和方向,也可表示位移的大小和方向。

t时刻看到的光线(图中的SA矢量),可分角为SO矢和SC矢。显然SC矢与地球的速度(或位移)OA是相同的,这样,我们站在地心看光波,感觉不到SC分量,光的速度只剩下SO分量(这好比两个人并行走路,感觉对方是相对自已是静止的)。SO矢量与CA矢量是相同的,考虑到光线最终落到A点,所以星体的视觉方向为AC。因此,视方向与真方向的差为图中的角CAB。以上推理是基于伽俐略变换的,我们总结如下:

t时刻星体的视方向是由t-τ时刻星体发出的并传到观测者的光线决定的,在地心坐标中,视方向由该光线的SO分量决定。换句话说,在地心坐标中,t时刻星体的视方向是t-τ是刻星体的真位置方向。

我们定义:t时刻,行星光行差等于视方向与真方向的角度差。这个角度差就是图中的角CAB。通常,角度差(行星光行差)用经纬度修正量来表度。反过来说,真方向加上行星光行差就得到视位置。

我们已经知道,行星光行差与t-τ=时刻的真方向有关,由于光行时τ是未知的,这给计算带来一点小麻烦。不过,多数情况下,B点的地心坐标(经度U、纬度V、距离Δ)已经精确算出,那么我们可以用τ=Δ/c来估算光行时,式中c为光速,Δ的单位中千米,τ的单位是秒。

接下来我们有两种计算方法求解修正量:

1)使用低精度算法重新计算t-τ时刻与t时刻天体的地心坐标。然后算出这两个坐标的经纬度差值即可,这个差值就是行星光行差。应注意,用t-τ的坐标值减t时刻的坐标值。还应注意,应使用完全相同的算法算出这两个坐标,并且计算过程中请保留所有的数字(16位),不要做任何舍入,以保证能够分辨出差值。

刚才说到“低精度”算法,具体含义是什么?我们的目标是准确计算出位置的差值,实际上是差值的精度达到3至4位有效数字即可。因为光行星光行差通常为几十角秒,如果达到3位有效数字,误差已小于0.05角秒。

当使用VSOP87算法时,周期项的表达为P=A*(t^n)*COS(B+C*t),进行差值计算时,P对t的导数P'越大,对差值的影响就越大,每个项的影响量=P'*τ/(365250*86400)弧度。除了几个A的值比较大的中短周期项的P’比较大,其它的都非常小,实际应用中每个表取4项已足够,这相当于把行星运动看作椭圆运动并加上1两个主要摄动项。

如果你使用VSOP87方法进行星历计算,建议使用这种方法进行光行差修正。因为,你已经拥有了星体位置计算的函数,调用该函数计算两次位置即可轻松算出行星光行差。当然,要适当控制计算的项数,以免影响程序效率。

2)先计算出因地球运动引起的偏差,在天文学中,这个偏差称为“光行差”,对应图中的角ASO。地球运动造成光线传播期间地球发生位移(图中的OA),它对S点的张角(即角ASO)就是恒星周年光行差,除了使用上述的“速度或位移的分解法”,我们也可以用速度合成的方法去理解它,恒星周年光行差有现成的公式可能使用。同样,由于行星运动,也形成了一个张角SAB,恒星周年光行差减去这个张角就得到行星光行差。行星光行差加上真位置就得到视位置。

几个注意点:

1)虽然星体B的运动改变了光的方向(伽俐略原理),但对于观测者,只能接收到到来自SA路径的光线,所以观测者并没有觉察到运动方向的改变。

2)计算光行时τ的时候,我们使用τ=Δ/c计算,式中Δ=|AB|。严格的说,τ=|CA|/c,然而CA是未知的。不过,由于行星的运动速度远小于光速c,所以|CA|与|AB|相差很小。因此,用τ=Δ/c计算产生的误差最多为v/c,式中v是行星的地心速度,v/c只有万分之几,引起的误差仅是毫角秒数量级的。另外,如果使用伽利变换,在地心坐标中,行星发出的光线相对地球的速度最大可达|c±v|,而我们计算τ时使用c,这也将造成一定的误差。其实,如果考虑相对论效应,伽利略修正的误差与Δ、c取值不严格所造成的误差是同一数量级的,所以没有必要把τ计算得很严格。


支持(0中立(0反对(0单帖管理 | 引用 | 回复 回到顶部
帅哥哟,离线,有人找我吗?
春光
  18楼 个性首页 | 信息 | 搜索 | 邮箱 | 主页 | UC


加好友 发短信
等级:版主 帖子:917 积分:7567 威望:10 精华:1 注册:2005/3/17 11:37:00
  发帖心情 Post By:2008/5/23 20:53:00

感谢xjw01兄的奉献。

支持(0中立(0反对(0单帖管理 | 引用 | 回复 回到顶部
帅哥哟,离线,有人找我吗?
xjw01
  19楼 个性首页 | 信息 | 搜索 | 邮箱 | 主页 | UC


加好友 发短信
等级:蝙蝠侠 帖子:866 积分:3804 威望:3 精华:20 注册:2008/3/20 22:14:00
  发帖心情 Post By:2008/5/24 23:05:00

通过以下改正,可以大幅度提高VSOP87的精度
1、以下对VSOP87日心球面地球序列拟合到DE406:
   利用VSOP87计算出黄经L后,计算黄经的修正量dL如下:
      dL = 0.2 -0.22*t -0.054*t2 + 0.1*cos(t)*t - 0.11*cos(1.5*t)*t -0.13*cos(3*t)
   式中t是J2000.0起算的儒略千年数
   最后 L = L+dL
   修正前的最大误差:
     J2000+-3000年1.00"
          +-2000年0.67"
          +-1000年0.39"
          +- 500年0.28"
          +- 250年0.17"
   修正后的最大误差:
     J2000+-3000年0.6"
          +-2000年0.25"
          +-1000年0.12"
          +- 500年0.02"
          +- 250年0.012"
   距离与黄纬不再修正,其最大误差为:
     距离:+-250年3E-8AU, +-1000年1E-7AU, +-3000年1E-7AU
     黄纬:+-250年0.01" , +-1000年0.024", +-3000年0.08"
2、以下对VSOP87的date球面地球序列拟合到DE406:
   dL = 0.2 +2.95*t -0.15*t2 - 0.02*t3
        - 0.01*sin(5*t) -0.13*cos(3*t)- 0.1*cos(1*t)*t
   最后L = L + dL
   修正后的最大误差:
     J2000+-3000年0.8"
          +-2000年0.4"
          +-1000年0.1"
          +- 500年0.02"
          +- 250年0.014"
  黄经拟合时已考虑了最新的岁差速度订正: -2.9965"/儒略千年
  黄纬与距离不用修正。
3、使用J2000.0动力学分点,因此J2000黄道坐标中的黄经与DE405/406相差一个估定值0.04"。
   注意1,使用J2000.0分点(非常接近ICRF的赤经起标点)得到的星体坐标才能与最新的《中国天文年历》相同。
   注意2,公元3000以后的拟合了“瑞士星历表”
4、许剑伟,地震后第12天


[此贴子已经被作者于2008-5-25 8:09:22编辑过]


支持(0中立(0反对(0单帖管理 | 引用 | 回复 回到顶部
帅哥哟,离线,有人找我吗?
xjw01
  20楼 个性首页 | 信息 | 搜索 | 邮箱 | 主页 | UC


加好友 发短信
等级:蝙蝠侠 帖子:866 积分:3804 威望:3 精华:20 注册:2008/3/20 22:14:00
  发帖心情 Post By:2008/5/26 7:46:00

上文中的“ 手表时(去除地区的时差后的手表时,即UT1)与原子时的差称为ΔT”应改为

 手表时(去除地区的时差后的手表时,即通过跳秒与UT1同步的时间)与原子时的差称为ΔT

UT1同由自转与恒星(春风点)确定的,是天文意义的时间系统,严格说,UT1与时角有关,即与角度有关

原子时则与“振荡次数”有关

  手表时UTC的秒长与时角无关,它与原子时的秒长完全相同。不过,UTC通过年终的跳秒实现与UT1同步,

这样就造成UTC与原子时不同步。


支持(0中立(0反对(0单帖管理 | 引用 | 回复 回到顶部
总数 115 上一页 1 2 3 4 5 6 7 8 9 10 下一页 ..12

返回版面帖子列表

天文算法讨论








签名