PyEphem：通过解析字符串表示来转换负角度的问题

Question

我一直在使用PyEphem来解决行星相对于地球某个位置的运动。然而，我有时会注意到，我从PyEphem得到的结果对于Declination而言似乎是不连续的，而Right Ascension看起来是连续的。每隔30分钟在这些图中拍摄RA和Dec坐标。

我希望恒星的身体能够以连续的方式运动，但是当赤纬是负的时候，它看起来是不连续的。

任何想法为什么会这样？我可以发布我的脚本，如果这也有帮助。

码：

from ephem import *
import datetime
import re

def cSTI(number):
    degrees = int(number[0])
    minutes = int(number[1])
    seconds = float(number[2])
    full = degrees + (minutes/60) + (seconds/60/60)

    return round(full,2)

planets = [Sun(),Mercury(),Venus(),Moon(),Mars(),Jupiter(),Saturn(),Uranus(),Neptune(),Pluto()]
masses =  [1.989*10**30,.330*10**24,4.87*10**24,0.073*10**24,0.642*10**24,1898*10**24,568*10**24,86.8*10**24,102*10**24,0.0146*10**24]
earthMass = 5.97*10**24
gravitational_constant = 6.67408 * 10**-11
theYear = '2018'
theMonth = '9'
theDay = '8'
revere = Observer()
revere.lon = '-42.4084'
revere.lat = '71.0120'
currentDate = datetime.datetime(2018,11,30,12,0,0)
lowerDateLimit = datetime.datetime(2018,9,1,12,0,0)
revere.date = currentDate
print('DATE SUN(RA) SUN(DEC)    SUN(AZI)    SUN(GFORCE) MERCURY(RA) MERCURY(DEC)    MERCURY(AZI)    MERCURY(GFORCE) VENUS(RA)   VENUS(DEC)  VENUS(AZI)  VENUS(GFORCE)   MOON(RA)    MOON(DEC)   MOON(AZI)   MOON(GFORCE)    MARS(RA)    MARS(DEC)   MARS(AZI)   MARS(GFORCE)    JUPITER(RA) JUPITER(DEC)    JUPITER(AZI)    JUPITER(GFORCE) SATURN(RA)  SATURN(DEC) SATURN(AZI) SATURN(GFORCE)  URANUS(RA)  URANUS(DEC) URANUS(AZI) URANUS(GFORCE)  NEPTUNE(RA) NEPTUNE(DEC)    NEPTUNE(AZI)    NEPTUNE(GFORCE) PLUTO(RA)   PLUTO(DEC)  PLUTO(AZI)  PLUTO(GFORCE)   ')

while (currentDate> lowerDateLimit):
    print('%s   ' % (revere.date),end = ' ')
    planetIndex = 0;
    for planet in planets:
        planet.compute(revere)

        rightascension = str(planet.ra)
        declination = str(planet.dec)
        azimuth = str(planet.az)

        rightascension = re.split('[:]',rightascension)
        declination = re.split('[:]',declination)
        azimuth = re.split('[:]',azimuth )

        rightascension = cSTI(rightascension);
        declination = cSTI(declination);
        azimuth = cSTI(azimuth);
        GFORCE = gravitational_constant * ((masses[planetIndex]*earthMass)/(planet.earth_distance**2))
        print('%s   %s  %s  %s  ' % (rightascension,declination,azimuth,GFORCE),end = ' ')
        planetIndex+=1
    print()
    currentDate += datetime.timedelta(minutes=-30)
    revere.date = currentDate

Answer 1

我相信问题在于你从ephem.Angle()（ra，azi等的对象）到float的手动转换。

---

EDIT

特别是，问题出现是因为在你的cSTI()函数中，当值为负时，你应该减去（而不是添加）不同的值。

更正后的实施方式如下：

import math

def cSTI(number):
    degrees = int(number[0])
    minutes = int(number[1])
    seconds = float(number[2])
    full = degrees + \
        math.copysign((minutes / 60), degrees) + \
        math.copysign((seconds / 60 / 60), degrees)
    return round(full, 2)

请注意，这是对代码的一些最小修改，以使其工作。我建议你从PEP8开始写一些关于如何编写更好代码的文档。此外，你应该避免在你的代码中使用像野生60这样的神奇数字。

如果我手动执行此操作，我也会避免使用不必要的正则表达式，并且实际上将从ephem.Angle()对象开始，如：

import math

MIN_IN_DEG = SEC_IN_MIN = 60

def ephem_angle_to_float(angle):
    degrees, minutes, seconds = [float(s) for s in str(angle).split(':')]
    value = abs(degrees) + \
        minutes / MIN_IN_DEG + \
        seconds / SEC_IN_MIN / MIN_IN_DEG
    return math.copysign(value, degrees)

这可以在你的代码中直接调用planet.ra或planet.dec。

但我不会手动这样做。见下文。

---

但好消息是没有必要手动计算你自己，你可以把它转换为float的弧度值，如official documentation所示。

以下代码是您的精美版本，它将表格数据生成为列表列表（可以使用Python标准库中的csv轻松转换为CSV）。已经进行了一些修改：

• 天体和群众现在更明确
• ephem对象是动态创建的
• G常数现在取自scipy.constants（从最新的CODATA获取）
• 标签和数据是动态创建的
• 尽可能避免显式索引
• 开始和结束日期有些颠倒，因为它与问题无关

请注意，此处不执行转换，因为这会使我们失去信息。

这是代码：

import numpy as np
import scipy as sp

import ephem
import scipy.constants
import datetime

celestial_masses = {
    'sun': 1.989e30,
    'mercury': 0.330e24,
    'venus': 4.870e24,
    'earth': 5.970e24,
    'moon': 0.073e24,
    'mars': 0.642e24,
    'jupiter': 1898e24,
    'saturn': 568e24,
    'uranus': 86.8e24,
    'neptune': 102e24,
    'pluto': 0.0146e24, }
celestials = {
    name: (eval('ephem.{}()'.format(name.title())), mass)
    for name, mass in celestial_masses.items() if name != 'earth'}
gg = sp.constants.gravitational_constant

observer = ephem.Observer()
# Revere, Massachusetts, USA
observer.lon = '-42.4084'
observer.lat = '71.0120'

start_datetime = datetime.datetime(2018, 9, 1, 12, 0, 0)
end_datetime = datetime.datetime(2018, 11, 30, 12, 0, 0)

values = 'ra', 'dec', 'azi', 'gforce'
labels = ('date',) + tuple(
    '{}({})'.format(name, value)
    for name in celestials.keys()
    for value in values)
data = []
observer.date = start_datetime
delta_datetime = datetime.timedelta(minutes=30)
while (observer.date.datetime() < end_datetime):
    row = [observer.date]
    for name, (body, mass) in celestials.items():
        body.compute(observer)
        row.extend([
            body.ra, body.dec, body.az,
            gg * ((mass * celestial_masses['earth']) / (body.earth_distance ** 2))])
    data.append(row)
    observer.date = observer.date.datetime() + delta_datetime

要转换为CSV（ra，dec，az和gforce的浮点值），可以这样做：

import csv

filepath = 'celestial_bodies_revere_MA.csv'
with open(filepath, 'w') as file_obj:
    csv_obj = csv.writer(file_obj)
    csv_obj.writerow(labels)
    for row in data:
        # : use default string conversion
        # csv_obj.writerow(row)

        # : a possible conversion to float for all but the `date`
        csv_obj.writerow([
            float(x) if i != labels.index('date') else x
            for i, x in enumerate(row)])

此外，这里有一些图表的代码，表明负值的问题已经消失：

import matplotlib as mpl
import matplotlib.pyplot as plt

x = [row[labels.index('date')].datetime() for row in data]
fig, axs = plt.subplots(4, 1, figsize=(16, 26))
for i, value in enumerate(values):
    pos = i
    for name in celestials.keys():
        y = [row[labels.index('{}({})'.format(name, value))] for row in data]
        axs[pos].plot(x, y, label=name)
        axs[pos].set_title(value)
        axs[pos].legend(loc='center left', bbox_to_anchor=(1, 0.5))
fig.tight_layout()

产生：