计算任意时刻格林尼治视恒星时角

近期学习卫星轨道方面的一些知识，遇到计算任意时刻格林尼治视恒星时角的问题，在网上搜了好久也没有一个完整的解决方案，后来通过，网上的一些零碎的信息，终于完成了计算格林尼治视恒星时角的程序，先整理如下。

计算格林尼治视恒星时角，首先需要计算当前时间的儒略日，计算方法如下

设Y为给定的年份，M为给定的月份，D为给定的日期。运算符INT为取所给数据的整数部分，若M大于2则Y与M不变，否则M加上12，Y减去1，换句话说，如果M为1月或二月，则被看做是前一年的13月和14月。

对于格里高利历（1582年10月15日以后），有

A=INT（Y/100）

B=2-A+INT（A/4）。

对于儒略历（即1582年10月15日以前），B=0。

所求儒略日即为

JD=INT（365.25（Y+4716））+INT（30.6001*（M+1））+D+B-1524.5

实现代码如下：

            const long IGREG = (15 + 31 * (10 + 12 * 1582));int jul;int ja, jy = iyyy, jm;if (jy == 0){MessageBox.Show("ERROR in subroutine julday: there is no year zero!");return 0;}if (jy < 0)++jy;if (mm > 2){jm = mm + 1;}else{--jy;jm = mm + 13;}jul = (int)(Math.Floor(365.25 * (jy+4716)) + Math.Floor(30.6001 * jm) + id - 1524.5);if (id + 31 * (mm + 12 * iyyy) >= IGREG)//判断是否为格里高利历日IGREG{ja = (int)(0.01 * jy);jul += 2 - ja + (int)(0.25 * ja);//加百年闰}return jul;

然后计算任意时刻的格林尼治视恒星时角，需要首先计算格林尼治平恒星时Θm，计算公式如下

Θm=100.46061837 + 36000.770053608T + 0.000387933T*T - T*T*T / 38710000

其中T =(JD - 2451545)/36525 ，JD为相应的儒略日，实现代码如下

double calGMST(double jd，double st)
{double D = jd - 2451545.0;//Compute the number of days since Jully2000double T = D / 36525.0;double GMST = 6.697374558 + 2400.051336*T + 0.000025862*T*T+st / 3600 *1.00273790935GMST = GMST % 24 * 15;return GMST;}

然后计算格林尼治视恒星时角Θ，计算公式如下

Θ = Θm + ΔΨ cos(εm + Δε)

其中εm为黄道倾角，计算公式如下：

εm=23°26′21.″448 - 46.″8150T - 0.″00059T*T + 0.″001813T*T*T

倾角章动和经度设计一下三个三角参数

L = 280.4665 + 36000.7698T
L′ = 218.3165 + 481267.8813T
Ω = 125.04452 - 1934.136261T

计算倾角章动用如下公式

ΔΨ= -17.20 sinΩ - 1.32 sin2 L - 0.23 sin2 L′ + 0.21sin2 Ω

Δε= 9.20 cos Ω + 0.57cos2 L + 0.10 cos 2 L′ - 0.09 cos 2 Ω

这些参数均以弧秒为单位

实现代码如下

double calGAST(double day,double yr,double mo,double st)double fjd = fulljd(day, yr, mo);double TJD = fjd - 2451545.0;//2451545为2000 1月1.5的儒略日double T0 = TJD / 36525.0;double THETAm = calGMST(fjd,st);double EPSILONm = 23.4392911111111111 - 0.0130041666667* T0 - 1.638e-07 * T0 * T0 + 5.0361e-07 * T0 * T0 * T0;double L = 280.4665 + 36000.7698 * T0;double dL = 218.3165 + 481267.8813 * T0;double OMEGA = 125.04452 - 1934.136261 * T0;double dPSI = -17.20 * Math.Sin(OMEGA) - 1.32 * Math.Sin(2 * L) - 0.23 * Math.Sin(2 * dL) + 0.21 * Math.Sin(2 * OMEGA);double dEPSILON = 9.20 * Math.Cos(OMEGA) + 0.57 * Math.Cos(2 * L) + 0.10 * Math.Cos(2 * dL) - 0.09 * Math.Cos(2 * OMEGA);//Convert the units from arc-seconds to degreesdPSI /= 3600;dEPSILON /= 3600;double GAST = (THETAm + dPSI * Math.Cos(EPSILONm + dEPSILON)) % 360;return GAST;
}

计算任意时刻格林尼治视恒星时角相关推荐

计算任意时刻椭圆轨道上行星的位置的插值方法
目录需求搜集到的方法和问题直接用牛顿力学公式进行模拟的问题物理常数调整问题初始化参量不直观实时求解动力学微分方程插值方法问题简化求解 θ ( t ) \theta(t) θ(t) ( ...
T-SQL 2 Tips: 1.计算任意两日期之间的周一到周日分别各有几个! 2.根据出生日期计算精确年龄!...
这两个小技巧,不写不知道,一写吓一跳! 都是看似简单,实际做起来就懵,得仔细想一想,才能写对! 凡是有日期运算的程序都要细心哦! 先说第二个: 2.根据出生日期精确计算年龄! 所谓计算精确年龄就是 ...
[Android] 任意时刻从子线程切换到主线程的实现
[Android] 任意时刻从子线程切换到主线程的实现 - Android移动开发技术文章_手机开发 - 红黑联盟转载于:https://www.cnblogs.com/melons/p/57919 ...
Python：用生成器的方式计算任意起止范围内质数的和。
用生成器的方式计算任意起止范围内质数的和. (质数,又称素数,是大于 1 的自然数,除了 1 和它本身外,不能被其他自然数整除.) def judge(num):if num == 1:return ...
NSURLSessionDataTask与NSOperationQueue实现多文件断点下载（任意时刻终止进程，重启应用，自动重启下载）...
效果展示 gif有点大,直接连接:7qnbrb.com1.z0.glb.clouddn.com/download.gi- 知识要点 NSOperationQueue线程队列的管理 NSURLSessi ...
凸多边形面积_C++计算任意多边形的面积
任意多边形的面积计算_拾忆楓灵的博客-CSDN博客blog.csdn.net 计算任意多边形的面积 - tenos - 博客园www.cnblogs.com 完美解决计算3D空间任意多边形面积_S ...
FineReport:任意时刻只允许在一个客户端登陆账号的插件
在使用FineReport报表系统中,处于账户安全考虑,有些企业希望同一账号在任意时刻智能在统一客户端登录.那么当A用户在C1客户端登陆后,该账号又在另外一个C2客户端登陆,服务器如何取判断呢? 开发 ...
java计算任意2个日期内的工作日_计算任意2个日期内的工作日(摘抄)
思路不错,对时间的操作挺全面,参考一下.代码如下: package test; import java.text.SimpleDateFormat; import java.util.Calendar ...
C语言编程--根据麦克劳林公式计算任意角的正弦余弦
C语言编程–根据麦克劳林公式计算任意角的正弦余弦正弦余弦麦克劳林公式问题输入格式 x ϵ 注:x 为角(单位:弧度),ϵ 为计算精度. 输出格式 y 正弦输入样例注:y 为 x 的正弦(余弦 ...

计算任意时刻格林尼治视恒星时角

计算任意时刻格林尼治视恒星时角相关推荐

最新文章

热门文章