PyEphem：通过解析字符串表示转换负角的问题
PyEphem: Problem converting negative angles by parsing the string representation
我一直在使用 PyEphem 来研究行星相对于地球上某个位置的运动。但是，有时我注意到我从 PyEphem 获得的结果对于赤纬似乎是不连续的，而赤经似乎是连续的。 RA 和 Dec 坐标每 30 分钟在这些图中获取一次。
我希望星体以连续的方式移动，但是当赤纬为负时，它看起来是不连续的。
知道为什么会这样吗？如果这也有帮助，我可以 post 我的脚本。
代码：
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
我认为问题出在您从 ephem.Angle()（ra、azi 等的对象）到 float 的手动转换中。
编辑
特别是，问题的出现是因为在您的 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)])
此外，这里有一些绘图代码，表明 negative 值的问题已经解决：
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()
产生：
我一直在使用 PyEphem 来研究行星相对于地球上某个位置的运动。但是，有时我注意到我从 PyEphem 获得的结果对于赤纬似乎是不连续的，而赤经似乎是连续的。 RA 和 Dec 坐标每 30 分钟在这些图中获取一次。
我希望星体以连续的方式移动，但是当赤纬为负时，它看起来是不连续的。
知道为什么会这样吗？如果这也有帮助，我可以 post 我的脚本。
代码：
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
我认为问题出在您从 ephem.Angle()（ra、azi 等的对象）到 float 的手动转换中。
编辑
特别是，问题的出现是因为在您的 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)])
此外，这里有一些绘图代码，表明 negative 值的问题已经解决：
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()
产生：