Расчет положения солнца и луны

[Sleepwalker]

Какой размер невязок?

[StormBS]

+- 0.5 град.

[Sleepwalker]

странно, у меня все красиво поулчилось

[StormBS]

Возможно орбитрон считает по другому  .

[Sleepwalker]

Напишите для теста входные параметры.

[StormBS]

Точка наблюдения 0,0; часовой пояс (zona) +4; 2010.10.17 14:32:06
_________________________________________________ __

      Dist   409427.10192282428   
      LatLune   -10.728224427141967   
           LonLune    139.36119864118135 
      H   -48.802665821805221   
      Az   106.21991237721139   

[Tau]

но к сожалению большая погрешность в расстоянии от точки наблюдения до луны (сравнивал с орбитроном), а также в расчете пр. восхождения и склонения

Какой размер невязок?

+- 0.5 град.

В сравнении с Орбитроном по прямому восхождению отклонение 7', по склонению 2'. Если вам нужно точнее, то используйте более точное разложение эклиптических параметров для Луны. Ну и заодно прецессию учтите, чтобы к видимой эпохе перейти.
Но дистанция как в Орбитроне все равно не получится, дистанция там удивительная.

[StormBS]

Более точно не нужно, просто интересно почему такая разница.

[Tau]

1) Полная аналитическая теория движения Луны содержит несколько тысяч членов. Есть более точные разложения, есть менее точные. Здесь использовалась самое короткое разложение, оно само по себе точность лучше 10' не дает.
2) Аргументом времени при вычислении элементов Луны является динамическое барицентрическое время в юлианских столетиях. Здесь использовалась юлианская дата в юлианских столетиях, патамушта некоторые дамы в принципе не любят вводить в свои игрушки изменяющиеся поправки, а разница между всемирным временем и динамическим - величина изменяющаяся. Сейчас разница равна 80 секунд, кажется. Но чтобы почувствовать влияние этих 80 секунд, нужно использовать более точное разложение для координат Луны.
3) Разложения для координат Луны обычно даются на какую-то эпоху, современные - на эпоху 2000 года. Из-за прецессии система координат смещается и нужно переводить координаты на видимую эпоху.

Если вы пишете какую-то программку с заявкой на серьезность, то нужно эти пункты учесть, убрать все упрощения.

[Sleepwalker]

      Dist   409427.10192282428   
      LatLune   -10.728224427141967   
           LonLune    139.36119864118135 
      H   -48.802665821805221   
      Az   106.21991237721139

Странно, если я подставляю входные данные в LunePos то даже таких точных значений не получается
А если заменяю блок получения геоц. координат Луны то выходит нечто более близкое к орбитрону:

Dist 409285.914
LonLune 139.273
LatLune -10.743
H -48.716
Az  106.215

Луну считал вот этим

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

static void SunMoon( double T, double *Sun, double *Moon )
{
static double AU = 1.49597870E+11;   
static double ABERTN = 99.3664E-06;  /* constant of aberation in rad */
static double T0 = 2415020.0;        /* epoch for this fourier series */
static double REarth = 6378140.0;

double DT;
double MEARTH, COSME, SINME, NURA;
double LONSUN, COSL, SINL, RS, E, OB, SOB, COB;
double LONPRT, NMOON, MR, MRi, LMMean, F, MM;
double cosmm, cosf, sinmm, sinf, sin2d;
double cos2d, sin2f, cos2f, sin2mm, cos2mm;
double SL, CL, S, cosd, sind, LATMN, LONGMN;
double MVENUS, MMARS, MJUP, DMNSUN;
double myTrue[3];

/* convert Julian time to ephemeris time wrt to epoch 1900 Jan 0.5 ET*/
DT = (T - T0 + 6.31E-04 + 2.85E-08 * (T-2446066.5))/36525.0;
/* compute earth mean anomaly, eqn 14 */
MEARTH = fmod(6.25658358 + 628.30194572 * DT, M2PI );
COSME = cos( MEARTH ); SINME = sin( MEARTH );
/* Mean anomalies for perturbing bodies */
MVENUS = fmod(3.70795199 + DT * 1021.32292286, M2PI);
MMARS  = fmod(5.57772322 + DT *  334.05561740, M2PI);
MJUP   = fmod(3.93187774 + DT *   52.96346477, M2PI);
/* mean elongation of the moon from the sun, eqn 15 */
DMNSUN = fmod(6.12152394 + DT * 7771.37719393, M2PI);
/* celestial longitude of the mean asc node of the moon's orbit, eqn 16 */
NMOON  = fmod(4.52360151 - DT *   33.757146246, M2PI);
NURA = -83.55E-06 * sin( NMOON );   /* nutation of ra ?  eqn 25 */
/* perturbation expressions are from Newcomb, Tables of the sun
A.P.A.E. vol. 6, pp.14-18 (1898)*/
LONPRT =  23.4553E-06 * cos( 5.220309 + MVENUS - MEARTH ) /*eqn 13 */
+ 26.7908E-06 * cos(2.588556 + 2.0 * (MVENUS - MEARTH))
+ 12.1058E-06 * cos(5.514251 + 2.0 * MVENUS - 3.0 * MEARTH)
+  9.9047E-06 * cos(6.001983 + 2.0 * (MEARTH - MMARS))
+ 34.9454E-06 * cos(3.133418 + MEARTH - MJUP)
+ 12.6051E-06 * cos(4.593997 - MJUP)
+ 13.2402E-06 * cos(1.520967 + 2.0 * (MEARTH - MJUP))
+ 31.0281E-06 * cos(4.035026 + DT * 0.352556)
+ 31.2898E-06 * sin( DMNSUN ) + NURA - ABERTN;
LONSUN = fmod(4.88162793 + DT * (628.33195099 + DT * 5.2844E-06), M2PI)
+ SINME * (33495.79E-06 - 83.58E-06 * DT
+ COSME * (701.41E-06 - 3.50E-06 * DT + COSME * 20.44E-06))
+ LONPRT;
COSL = cos( LONSUN );
SINL = sin( LONSUN );
RS = AU * (1.0002806 - COSME * (1675.104E-05 - 4.180E-05 * DT
+ COSME * 28.060E-05));
/* true obliquity of the ecliptic, eqn 1 & 26 */
OB = 0.40931970 - 2.27111E-04 * DT + 4.465E-05 * cos ( NMOON );
COB = cos( OB ); SOB = sin( OB );
E = NURA * COB;                   /* eqn of equinoxes, eqn 2 */
/* convert to cartesian coordinates, eqns. 6-10 scaled to meters */
Sun[0] = RS * (COSL + SINL * COB * E ); 
Sun[1] = RS * (SINL * COB - COSL * E );
Sun[2] = RS * SINL * SOB;

/* lunar ephemeris is based on Lrown's Lunar Theory Ammended.
The expressions are adapted from Woolard, theory of the rotation the
Earth, A.P.A.E Vol.15, Part 1 (1953), Tables 4, 6, and 7.

The moon vector is computed at the Earth's center in the geometric
direction of the moon.  Neglecting the light travel time of 1.3 sec
between Earth and moon gives a negligible error of 0.0002 deg for
the moons's direction.
*/

LMMean = fmod( 4.71996657 + DT * (8399.70914492 - DT * 1.9780E-05), M2PI);
F      = fmod( 0.19636505 + DT * (8433.46629117 - DT * 5.6044E-05), M2PI);
MM     = fmod( 5.16800034 + DT * (8328.69110367 + DT * 1.60425E-04),M2PI);
DMNSUN -= (2.5065E-05 * DT * DT);   /* eqn 33 */
NMOON  += (3.6264E-05 * DT * DT);   
/* fundamental frequency sines and cosines */
cosmm = cos(MM);     sinmm = sin(MM);
cosd  = cos(DMNSUN); sind  = sin(DMNSUN);
cosf  = cos(F);      sinf  = sin(F);
/* double frequency sines and cosines */
sin2d   = 2.0 * sind * cosd;    cos2d  =  cosd * cosd  -  sind * sind;
sin2f   = 2.0 * sinf * cosf;    cos2f  =  cosf * cosf  -  sinf * sinf;
sin2mm  = 2.0 * sinmm * cosmm;  cos2mm = cosmm * cosmm - sinmm * sinmm;
/* lunar longitude */
LONGMN = sind * (-607.0 + 84.0 * COSME - 131.0 * cosmm )
+ sin2d * (11490.0 + 134.0 * cos2d + 239.0 * cos2f
+ 683.0 * COSME
+  cosmm * (23098.0 + 392.0 * cos2d + 918.0 * COSME
- 80.0 * cos2f + 138.0 * cos2mm)
+  cos2mm * (1096.0 + 300.0 * cos2d))
+  sin2f * (-1996.0 - 295.0 * cos2d
+ cosmm * (-411.0 - 20.0 * cos2d))
+  SINME * (-3238.0 + 90.0 * cosd - 919.0 * cos2d - 72.0 * COSME
+  cosmm * (-1248.0 - 1222.0 * cos2d) - 84.0 * cos2mm)
+ sin2mm * ( 3876.0 + 10.0 * COSME
+ cos2d * ( -956.0 - 296.0 * cos2d ) )
+ (sin2mm * cosmm + cos2mm * sinmm ) * ( 175.0 -59.0 * cos2d )
+ 68.9 * sin( 6.042967 + 2.318888 * DT ) - 49.0 * sin( NMOON )
+ sinmm * ( sin2d * ( -1052.0 * SINME -72.0 * sin2f )
+ COSME * ( 184.0 - 804.0 * cos2d )
+ cos2d * ( -21305.0 - 352.0 * cos2d )
+ cos2f * (  -27.0 + 8.0 * cos2d )
+ 49.0 * cosd + 109936.0 );
LONGMN = LONGMN * 1.0E-06 + LMMean;
/* lunar latitude */
LATMN = cosf * ( -49.0E-06 * sind
+ sin2d * ( 3592.0E-06 +  119.0E-06 * COSME + 2008.0E-06 * cosmm
+   85.0E-06 * cos2mm )
+ SINME  * (  -55.0E-06 -   24.0E-06 * cos2d - 112.0E-06 * cosmm )
+ sinmm  * ( 9744.0E-06 - 1540.0E-06 * cos2d +  10.0E-06 * COSME )
+ sin2mm * (  455.0E-06 -   49.0E-06 * cos2d ) )
+ sinf * (-3.0E-06 * cosd - 2454.0E-06 * cos2d - 60.0E-06 * cos2f
+ 6.0E-06 * sinmm * SINME
+ sin2d * (-157.0E-06 * SINME+248.0E-06 * sinmm-79.0E-06 * sin2mm )
+ COSME * (  -7.0E-06 -   53.0E-06 * cos2d )
+ cosmm * (  50.0E-06 + 72.0E-06 * cos2d + 4.0E-06 * COSME )
+ cos2mm * (147.0E-06 - 87.0E-06 * cos2d )
+ 89473.0E-06 );
SL = sin( LATMN ); CL = cos( LATMN );
/* the following numerator is the inverse of the constant of sine
parallax for the moon, 3422.451 arc-sec, but iin radians. */
MRi = 9999237.0 - 2858.0 * cosd
+ sin2d * ( 6481.0 * SINME + sin2mm * (-1715.0 + 2144.0 * cos2d))
+ sinmm * (-2083.0 * sin2f +  SINME * ( 6142.0 - 4067.0 * cos2d )
+ sin2d * (91230.0 + 3260.0 * cos2d + 3027.0 * COSME))
+ cos2d * (  82488.0 + 1526.0 * cos2d)
+ COSME  * ( -1169.0 + 4727.0 * cos2d +   67.0 * cos2mm)
+ cosmm  * (541309.0 - 2089.0 * cos2f + 3634.0 * cos2mm
+ cos2d  * (109264.0         + 3764.0 * cos2d)
+ SINME  * (  5691.0 * sin2d + 1350.0 * sinmm)
+ COSME  * (   596.0         + 4091.0 * cos2d ))
+ cos2mm * ( 28598.0 + cos2d * ( -61.0 + 2204.0 * cos2d ));
MR = REarth * 60.26816E+07 / MRi;
myTrue[0] =  CL * cos( LONGMN );
S = CL * sin( LONGMN );
myTrue[1] =  ( S * COB - SL * SOB);
Moon[2] = MR *  ( S * SOB + SL * COB );
/* convert x and y components from true equinox to mean equinox of date
so components are based on true equator and mean equinox. */
Moon[0] = MR * ( myTrue[0] + myTrue[1] * NURA );
Moon[1] = MR * ( myTrue[1] - myTrue[0] * NURA );
return;
}


[judge]

В книгах Астрономические алгоритмы Jean Meeus есть разложения для долготы, широты, расстояния солнца, планет и тд.
Где их можно взять в электронном виде (киньте кто нибудь мне) или подскажете более новые разложения.

[Kirill]

В книгах Астрономические алгоритмы Jean Meeus есть разложения для долготы, широты, расстояния солнца, планет и тд.
Где их можно взять в электронном виде (киньте кто нибудь мне) или подскажете более новые разложения.

Возьмите DE405, получите координаты и по ним получите расстояния. Полный ряд для Луны это для гурманов. Не полный - думаю уже заметили, люди мучительно решают чей калькулятор правее и авторитетней

Имейте ввиду. Исходник с хорошей точностью сразу счастья не даст. Там весело будет уже начиная со времени. Для начала нужно вести базу с секундами (TAI - UTC), далее регулярно сливать файлы с EOP и по точкам интерполировать. А без хорошего времени не будет и полной точности.

[judge]

Возьмите DE405, получите координаты и по ним получите расстояния. Полный ряд для Луны это для гурманов. Не полный - думаю уже заметили, люди мучительно решают чей калькулятор правее и авторитетней

Имейте ввиду. Исходник с хорошей точностью сразу счастья не даст. Там весело будет уже начиная со времени. Для начала нужно вести базу с секундами (TAI - UTC), далее регулярно сливать файлы с EOP и по точкам интерполировать. А без хорошего времени не будет и полной точности.

Если вас не затруднит скиньте мне DE405 пожалуйста, а немного запутался в поисках на файловых серверах...
На счет Луны, это отдельная тема и сейчас она не обязательна, но потом поспорю где-же все таки правее...

Суперская точность  конечно не нужна, но лучше чем преобразованиями Кеплера да...

[Kirill]

Если вас не затруднит скиньте мне DE405 пожалуйста, а немного запутался в поисках на файловых серверах...

ftp://ssd.jpl.nasa.gov/pub/eph/planets/

В папке ascii лежат данные по годам. Там же найдете два исходника на Яве и Си для доступа к этим данным.

Могу отсканировать страницы из Альманаха 2006 года (английский естественно). Там описывается последовательность расчетов с примером на дату. Правда нет полной серии ни для прецессии ни для нутации. В примере расчета используется интерполяция табличных значений и отдельно есть приблизительные формулы.

Все известные ошибки я аккуратно распечатал и заклеил поверх страницы поправками. Помоему уж лучше проверять свою реализацию по их примеру, чем гадать где я нахожусь по отношению к пяти выбранным планетариям.

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

Для Луны есть программа Виртуальный Лунный Атлас как то так называется. Исходники на Паскале, найдете. Написано хреново, но там можно выдрать две версии расчетов: с эфемеридами и без них.

Суперская точность  конечно не нужна, но лучше чем преобразованиями Кеплера да...

Дело не в том где лучше/хуже. Под конкретную задачу выбирается. Это же просто интересно знать как мы летим на шарике

Язык реализации алгоритмов вам какой роднее?

[judge]

Kirill
что касается языка то паскаль  , но язык для меня не главное, а главное само разложения с коэффициентами  и алгоритмом расчета.
Если вас не сильно затруднит , от сканируйте  страницы из Альманаха...

Дело не в том где лучше/хуже. Под конкретную задачу выбирается. Это же просто интересно знать как мы летим на шарике
Язык реализации алгоритмов вам какой роднее?

Согласен полностью с вами) Интересно

[Kirill]

Kirill
что касается языка то паскаль  , но язык для меня не главное, а главное само разложения с коэффициентами  и алгоритмом расчета.
Если вас не сильно затруднит , от сканируйте  страницы из Альманаха...

Сделаю, как только днем буду рядом со сканером.

[Kirill]

что касается языка то паскаль 

Лучшее что вспомнил на Паскале это Astronomy on the Personal Computer второе издание, там исходники прилагаются.

Отправил в приват.