本来是球面天文学的作业，不过感觉挺有意思，于是打算发出来。

JPL，指的是 Jet Propulsion Laboratory。喷气推进实验室。DE，指的是 Development Ephemeris，星历。使用较为广泛的为DE405，儒略历范围为2305424.5至2525008.5，即为公历1599.12. 09至2201.2.20。包含了八大行星及太阳、月亮、冥王星、日地质心、地月质心、平动章动，可选择为质心或目标从而得出数据。

以Linux系统为例，首先通过open命令连接到ssd.jpl.nasa.gov，采用anonymous匿名登陆。

通过cd与ls命令移动至

/pub/eph/planets路径下

本次实践采用de405星表与Fortran语言，需要的文件路径为：

采用get或mget命令将其下载至本地。

以Linux系统为例，将下载的de405中的文件合并，每一个ascp文件代表20年的星表数据，例如ascp2000.405为2000年至2020年的数据。

首先尝试制作2000-2040年星表。

将header.405放在前，生成合并ascii文件。

需要对fortran文件夹下asc2eph.f进行编辑，确认NRECL值。用vi编辑器打开。

到文件74、75行，按 i 进入编辑模式，

将其中一项句首用于注释的c删除

然后按Esc退出编辑模式，

依次输入冒号，wq，enter保存退出。

编译asc2eph.f文件，可通过-o选项指定生成文件名。

将先前生成的生成的文件jpleph00-40用其打开：

若正常则会返回常数与历表年份等值。

即可发现文件夹下生成了历表文件JPLEPH。

若该过程发生了错误，可能是是NRECL参数设置错误，将其更改后再次重复上述步骤。

Windows系统下过程差别不大，生成制作部分为使用copy命令：

下载文件中/fortran/testeph.f中包含了许多实用的子程序，可将其单独分离出来保存。

首先对testeph.f中未输入的参数进行补全，分别为251、316、388、417行，

前三项为FSIZER1、FSIZER2、FSIZER3的NRECL参数，推荐填写与先前asc2eph.f一样的数值。

417行处填写与采用的星表有关，对于de405填写2036

此外，在939-941行，选择FSIZER三种其中一项，通常Linux采用2，其余采用3。将句首C删除即可。保存。

打开testeph.f文件，219行后为子程序即SUBROUTINE部分，将其单独保存，取名为jplsub.f

将前面部分保存命名为testeph.f

随后对其进行测试。

若能读出天文常数则运行成功。

若出现错误，尝试修改CALL FSIZER的选择，或者检查NAMFIL一项是否路径正确的指向了JPLEPH文件。

在接下来的操作中，只需自己编写fortran文件，CALL命令调用PLEPH，由userguide.txt文件可知，通过输入儒略历时间，参考系中心，目标天体，则可返回目标天体的三维位置与速度共六个值。

例如编写了venus.f文件后，

则可返回结果。

以地球为中心，得出同一时刻太阳与金星坐标。

凌日时应该满足条件为

为避免太阳遮挡金星，太阳距离为金星距离一倍以上

由于太阳视直径约为0.5度，两者间角距离小于0.25度

由于为弧度制，所以设定为0.0045

限定时间为2000.1.1至2013.1.1

代码如下

计算结果为：

可得发生金星凌日

2453164 即为2004.6.7

2456084 即为2012.6.5

制作2000-2200星表，并在代码中修改限定时间2000.1.1至2219.12.31

对代码进行了部分改进。

计算结果为：

可得下一次发生金星凌日日期：

2494622 即为2117.12.10

2497542 即为2125.12.8

取用第三次凌日结果，用stellarium进行验证：

起始时：

结束时：

符合精度，较为满意。

首先利用余弦公式计算出水星与太阳间的夹角。

达到水星西大距时需要同时满足条件

设定从2020年3月1日开始计算，即为ET=2458910，代码如下

运行结果为

可得第一个时间节点为

2458932.5830037016 即为2020,3,24,上午10点左右

设定结束日期为：2760073

可见每次水星西大距的间隔约为115天左右。

JPL星表采用的是ICRS，与经典赤道直角坐标系略有偏差，在此忽略不计。三维坐标xyz满足如下关系：

其中$\delta$为赤纬，$\alpha$为赤经。

计算过程中会产生截断误差，在此忽略不计。参考马文章老师《球面天文学》天体出没的时刻与方位角一章，得到赤经赤纬后可以计算太阳与水星升起时的时角与方位角。

需要注意的是，计算赤经赤纬应用反正弦函数，若用反余弦可能会产生正负号错误。

假设太阳刚刚升起时为水星最佳观测时间。

计算太阳升起时的时角

$t$会有正负两个值，对应升起与落下时角。

代入太阳赤经得到太阳升起时的恒星时。

再代入水星赤经得到此时水星的时角。

对于地平坐标与时角坐标转换有：

在此为避免代码堆在一起显示比较混乱，单独编写一个计算观测方位的代码：

可得结果

可得地平高度为$3^\circ17’$，方位角为$106^\circ27’$

用stellarium进行验证，对齐恒星时$-98.865303$，即为$17h24m32s$左右。

可以发现符合得很好！

金星凌日发生日期:

2453164 即为2004.6.7

2456084 即为2012.6.5

2494622 即为2117.12.10

2497542 即为2125.12.8

最近一次水星西大距：

对于判断思路可以更加提高，觉得队友“本来无一物”的方案就很好，而且试验过后效果显著，其判断思路为：

对于这一部分感觉思路还行，暂无改进思路

步骤的关键为确定恒星时作为标准，在这里采用了马文章老师确定太阳升时方法，此刻太阳尚未升起，如图。

但此刻天空即将亮起，所以此刻可为适宜观测时机。

对于稍有误差可能是由于星表采用的是ICRS，与经典赤道直角坐标系略有偏差，但总体符合得非常好。

在刚才提交的作业版本中，对该部分计算错误。

正确值应为$1.58534$重新代入计算。可得结果：

软件显示为：

可见这次计算才是正确的，太阳即将升起。但是水星却因为天亮了，看不清。

联系晨光始概念：

在此选用民用。

所以可以采用的其他解题思路为：

利用地平坐标与时角坐标变换，设太阳地平高度为$-6^\circ$，代入太阳赤纬与地理纬度，得到此时太阳时角，这里由于太阳还在地平面下，就不用考虑太阳圆面与大气作用了。

利用太阳时角与赤纬，得到此刻恒星时。

由此刻恒星时与水星赤纬，得到水星时角。

由水星时角与地平坐标时角坐标关系，得到水星此刻方位角与地平高度。

利用软件，

当太阳地平高度为$-6^\circ$时，水星地平高度$10^\circ50’$，方位角$107^\circ34’$

此外，还有改进的地方是需要在代码中添加计算由当地恒星时转换到当地平太阳时的部分。