Определение ПСТ по TLE

[Sleepwalker]

Вот...
А у меня 3.4759871545186476 3.4863432329915702

Тут и нужно ковырять.

theta ( theta = atan2(Z,X) ) = -0.0099568520437508332
theta ( theta = atan2(Y,X) ) = 0.00039922642917146830

Среднее звездное время (Результат функции ThetaG(time)) = 2.7972413006171877

UPD: поправил оси

[art193]

 theta=-1,5705733533
 ThetaG(time)=0,000686873947513983

Странно... У вас используется atan2(Z,X), у меня вместо нее используется ArcTan, по описанию они должны были возвращать одинаковый результат. Функции atan2 в Делфи нет...

[Sleepwalker]

Так, временно отойдем от арктангенсов, т.к. звездное время уже не верное.
Функция среднего звездного времени:

double ThetaG(double jd) { //jd = 2454805.7385309837

double       temp, tut1;
tut1 = (jd - 2451545.0) / 36525.0;
temp = -6.2e-6* tut1 * tut1 * tut1 + 0.093104 * tut1 * tut1 + (876600.0*3600 + 8640184.812866) * tut1 + 67310.54841;  // sec
temp = fmod(temp*degToRad/240.0, Pi2); //360/86400 = 1/240, to deg, to rad
if (temp < 0.0)temp += Pi2;
}

[art193]

   Ого, по моему эта функция несколько отличается от той, которую вы давали ранее...

double CSGP4_SDP4::ThetaG(double jd)
{
// Reference:  The 1992 Astronomical Almanac, page B6.
   double   UT,TU,GMST, dummy;
  UT   = modf(jd + 0.5, &dummy);
  jd   -= UT;
  TU   = (jd - 2451545.0)/36525.0;
  GMST   = 24110.54841 + TU * (8640184.812866 + TU * (0.093104 - TU * 6.2E-6));
  GMST   = Modulus(GMST + secday*omega_E*UT,secday);
  return (2.0*PI*GMST/secday);
}

  У меня по ней считает...

[Sleepwalker]

А, может быть, секундочку...

[Sleepwalker]

Я по привел пример более "легкой" и менее точной реализации.
Проверил по старой, у меня сошлись значения.

Ок, тогда чему у нас равны:
omega_E
secday

[art193]

 Они константы, берутся из соответствующего модуля и равны
  secday   = 86400.0;
  omega_E  = 1.00273790934;

[Sleepwalker]

Хм.. где ж тогда расхождение
UT = 0.23853098368272185
TU = 0.089267624914442165
GMST = 38464.828992863419
результат 2.7972413005918377

[art193]

Вот...
UT = 0.23853098368272185
TU = 0.089267624914442165
GMST = 9,44519478

[art193]

 GMST   = Modulus(GMST + secday*omega_E*UT,secday); эквивалентно

 GMST   = (GMST + secday*omega_E*UT)/(secday)

[Sleepwalker]

Ну собственно все. где-то тут в этих двух строчках косяк

Код: [Выделить]

GMST   = 24110.54841 + TU * (8640184.812866 + TU * (0.093104 - TU * 6.2E-6));
GMST   = Modulus(GMST + secday*omega_E*UT,secday);

промежуточное значение GMST в первой строке : 795399.32621829608

[Sleepwalker]

GMST   = Modulus(GMST + secday*omega_E*UT,secday); эквивалентно

 GMST   = (GMST + secday*omega_E*UT)/(secday)

Неа
вот код

Код: [Выделить]

static double Modulus(double arg1, double arg2)
{
double modu;
modu = arg1 - (long)(arg1/arg2) * arg2;
if (modu >= 0.0) return modu;
else return modu + arg2;
}


[art193]

   Да, не эквивалентны, уже сам понял... все ThetaG(time) = 2,79724130059183 - не совсем, но думаю, что отличия в 10 знаке после запятой не сильно скажутся...
   А что делать с theta = atan2(Y,X) ?

[Ivan Mhitarov]

Функции atan2 в Делфи нет...

Гм. Никогда не говори никогда.

Calculates the arctangent angle and quadrant of a given number.

Unit

Math

Category

Trigonometry routines

Delphi syntax:

function ArcTan2(const Y, X: Extended): Extended;

Description

ArcTan2 calculates ArcTan(Y/X), and returns an angle in the correct quadrant.  The values of X and Y must be between –2^64 and 2^64.  In addition, the value of X can’t be 0.  The return value will fall in the range from -Pi to Pi radians.

[art193]

  Все верно, ArcTan2(Y,X) существуе и вычисляет значение ArcTan(Y/X). Я использую ArcTan(Y/X)... (тоже самое, только аргументы делятся, а не перечисляются), замена на ArcTan2 результат не меняет.

[Sleepwalker]

Не, нефига. Не просто делятся, но и выбирается квадрант.
А то что результат не меняется... это мы сейчас разберемся.

[art193]

 Получается:
ArcTan2(3,06777284853136,-0,000684032094947294)=1,5710193003061

[Sleepwalker]

Хм... а вот это что за числа 3,06777284853136,-0,000684032094947294 ?

Такое впечатление, что это скорость КА

[art193]

Результаты расчета Convert_Sat_State

Пример TLE:
1 26056U 00002A   08340.23853098  .00000088  00000-0  10000-3 0  4865
2 26056 000.5746 083.1707 0002318 142.0726 134.7604 00.99643563 32538

SDP4

X     42346,1467857657
Y     16,9057018686055
Z     -421,648252175922

после Convert_Sat_State

XDOT  -0,000684032094947279
YDOT  3,06777284853136
ZDOT  0,00366539951216375

У нас вроде эти данные совпадали... (Ответ 11)

[Sleepwalker]

Все правильно.

Только на вход ПСТ идет вектор положения, а не скоростей.
Соответственно арктан2 будет иметь свои аргументы Y,X : 16,9057018686055,  42346,1467857657

а чего метод Convert_Sat_State делает ?

UPD: Все, нашел.  Convert_Sat_State преобразует вектора положения и скоростей из относительных в км и км/c и вычисляет модуль.

XDOT,YDOT, ZDOT это элементы вектора скорости, а не результат преобразования координат XYZ