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

2024-10-03 11:23:47 发布

您现在位置:Python中文网/ 问答频道 /正文

我一直在用PyEphem来研究行星相对于地球上某个位置的运动。然而,我有时注意到,我从PyEphem得到的结果对于赤纬是不连续的,而赤经则是连续的。RA和Dec坐标每30分钟采集一次。在

我认为恒星体会以连续的方式运动,然而当赤纬为负时,它看起来是不连续的。在

你知道为什么会这样吗?我可以张贴我的脚本,如果这也有帮助。在

DecMa

代码:

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

Tags: renumberdatetimedecraprintplanetazimuth
1条回答
网友
1楼 · 发布于 2024-10-03 11:23:47

我相信问题出在您从ephem.Angle()(它是raazi等的对象)到{}的手动转换。在


编辑

特别是,这个问题是因为在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()对象开始,例如:

^{pr2}$

这可以直接调用到代码中的planet.ra或{}。在

但我不会用手工的。见下文。在


但好消息是,您不需要自己手动计算,您可以将其转换为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(使用radecazgforce的浮点值),可以这样做:

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()

产生:

Plots

相关问题 更多 >