Орбитальные элементы Солнечной системы с использованием Pyephem

Я использую астрометрический модуль Python Pyephem и хочу получить орбитальные (кеплеровские) элементы для планет Солнечной системы.

Единственные значения, которые я нашел, это гелиоцентрическая широта, долгота и расстояние до солнца. Есть ли способ вычислить параметры орбиты на основе этих значений? Я пропустил функцию в Pyephem?

На данный момент у вас есть только вектор положения. Чтобы получить полный набор кеплеровских элементов, вам также нужен вектор скорости.
Хотя возможно использовать положение и скорость для расчета элементов орбиты, это работает в обратном направлении. pyephem использует элементы Кеплера для вычисления положения. Орбитальные элементы фиксированы, и их можно прочитать на сайте Horizons.
Вы можете попробовать посмотреть на небесное поле, но, похоже, даже оно не имеет непосредственно элементов орбиты. Я знаю, что в файлах эфемерид тоже нет элементов орбит, но вы можете попробовать взглянуть на VSOP 2013. Элементы орбит довольно постоянны, но изменяются медленно из-за возмущений.

Ответы (1)

Существует веб-сайт: http://orbitsimulator.com/formulas/OrbitalElements.html , на котором есть программа javascript для преобразования векторов состояния в элементы орбиты и обратно.

Источник этого веб-сайта можно преобразовать в python: вот Bodyкласс, в котором есть метод расчета элементов орбиты на основе этого веб-сайта. Метод принимает один аргумент, телом является Руководитель (т.е. Солнце)

G = 0.0002946        # in units of seconds, AU and solar masses.
class Body:
    """A body has attributes r and v, which are its position and 
velocity in cartesian coordinates and a mass. implied units are 
solar masses, AU and seconds.""" 
    def __init__(self,r,v,mass):
        self.r = np.array(r,dtype="float")
        self.v = np.array(v,dtype="float")
        self.mass = mass
        self.GM = self.mass*G

    def orbital_elements(self,principal):
        '''view-source:http://orbitsimulator.com/formulas/OrbitalElements.html'''
        mu = G*(principal.mass+self.mass)
        # calculate relative position,velocity
        r = self.r - principal.r
        v = self.v - principal.v
        try: #catch division by zero
            R = np.linalg.norm(r)
            V = np.linalg.norm(v)
            a = 1/(2/R - V**2/mu)  # semi major axis

            h = np.cross(r,v)
            H = np.linalg.norm(h)

            P = H**2/mu
            Q = np.dot(r,v)

            E= np.sqrt(1-P/a)  #eccentricity

            e = [1-R/a,Q/np.sqrt(a*mu)]
            i = np.arccos(h[2]/H)
            Omega = 0
            if i!=0: 
                Omega = np.arctan2(h[0],-h[1]) #Longitude of acending node

            ta = [H**2/(R*mu) -1,H*Q/(R*mu)]
            TA = np.arctan2(ta[1],ta[0])
            Cw = (r[0]*np.cos(Omega)+r[1]*np.sin(Omega))/R


            if i==0 or i==np.pi:
                Sw = (r[1]*np.cos(Omega) - r[0]*np.sin(Omega))/R
            else:
                Sw = r[2]/(R*np.sin(i))

            omega = np.arctan2(Sw,Cw) - TA  #argument of periapsis
            if omega<0: omega += 2*np.pi

            u = np.arctan2(e[1],e[0]) 
            M = u - E*np.sin(u) # mean anomaly

            return(a,E,omega,Omega,i,M)
        except ZeroDivisionError:
            #meaningless, but avoids crash
            return(0,0,0,0,0,0)
Мне понравился этот ответ, но меня смущает уменьшенная масса. В этой формулировке заданный входной вектор состояния ( р , в ) дает разные ответы в зависимости от обеих масс. Относительно какой системы координат предполагается измерять вектор состояния - системы главного или системы центра масс? К какой орбите относятся возвращенные элементы - к орбите Тела в системе координат главного, или к центру масс, или это однотельная орбита в центральном поле?