Телескопы покупают здесь


A A A A Автор Тема: Расчёт даты григорианского календаря для большой отрицательной юлианской даты  (Прочитано 946 раз)

0 Пользователей и 1 Гость просматривают эту тему.

Оффлайн BerezkaАвтор темы

  • Новичок
  • *
  • Сообщений: 17
  • Благодарностей: 0
    • Сообщения от Berezka
Доброго времени суток. Мне понадобилось научиться определять дату григорианского календаря по очень большой отрицательной юлианской дате. Это, к примеру, JD = -3100015.5. Я работаю с эфемеридами NASA и в их 431-м выпуске это начальная дата. Согласно их данным, данное значение соответствует следующей дате: 15 августа -13200. Обычно для перевода я использовал алгоритм с википедии, однако он не подходит, так как не учитывает количество накопившихся суток для такого большого периода. И вот тут я попал в тупик. Пришлось перерыть кучу ресурсов с алгоритмами переводов, но ни один из них не подходил для такого большого значения.
  • Я смотрел алгоритм  в "David A. Vallado — Fundamentals of Astrodynamics and Applications", однако он там верен на очень маленьком периоде (от 1900-го до 2100-го года).
  • Также есть алгоритм в "P. Kenneth Seidelmann — Explanatory Supplement to the Astronomical Almanac, 1992", который использует IAU SOFA (http://www.iausofa.org/2018_0130_C/sofa/jd2cal.c). Однако он верен только для JD >= -68569.5).
  • Пытался найти алгоритм перевода у самих NASA (NAIF SPICE TOOLKIT), но там исходный код скомпилирован в библиотеки, так что ничего подсмотреть не удалось.
Но не всё так плохо, этот перевод возможен не только средствами NASA. Я нашёл онлайн-калькулятор для перевода юлианской даты в календарную на русскоязычном ресурсе: http://cosmos-online.myhomeinet.ru/jd.php#.XS81tegzaUk. И вот он определяет данную юлианскую дату верно, то есть, как её и представляют сами NASA, однако процесс вычисления там не описан.

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

Оффлайн xd

  • *****
  • Сообщений: 17 977
  • Благодарностей: 378
    • Skype - deimos.belastro.net
  • Награды Открытие комет, астероидов, сверхновых звезд, научно значимые исследования.
    • Сообщения от xd
    • Белорусская любительская астрономическая сеть
А какой смысл григорианского календаря в столь отдалённом прошлом?
У природы нет плохой погоды, у неё просто на нас аллергия.

Учение без размышления бесполезно, но и размышление без учения опасно /Конфуций/
Слово есть поступок. /Л. Толстой/

Оффлайн BerezkaАвтор темы

  • Новичок
  • *
  • Сообщений: 17
  • Благодарностей: 0
    • Сообщения от Berezka
А какой смысл григорианского календаря в столь отдалённом прошлом?

Для меня смысл только в том, как в него перевести Юлианскую Дату, то есть чисто практический.

Оффлайн xd

  • *****
  • Сообщений: 17 977
  • Благодарностей: 378
    • Skype - deimos.belastro.net
  • Награды Открытие комет, астероидов, сверхновых звезд, научно значимые исследования.
    • Сообщения от xd
    • Белорусская любительская астрономическая сеть
Может сразу работать в системе непрерывного счёта времени и не заморачиваться?
У природы нет плохой погоды, у неё просто на нас аллергия.

Учение без размышления бесполезно, но и размышление без учения опасно /Конфуций/
Слово есть поступок. /Л. Толстой/

Оффлайн BerezkaАвтор темы

  • Новичок
  • *
  • Сообщений: 17
  • Благодарностей: 0
    • Сообщения от Berezka
Может сразу работать в системе непрерывного счёта времени и не заморачиваться?
У меня просто сейчас стоит задача именно научиться переводить даты таких порядков. Не более, не менее.

Оффлайн Klapaucius

  • *****
  • Сообщений: 11 260
  • Благодарностей: 184
  • Илья
    • Сообщения от Klapaucius
В общем, надеюсь на вашу помощь в данном вопросе. Если у кого-то есть методика для вычисления юлианской даты на большом отрезке времени, то буду очень благодарен за предоставленную помощь.
Вопрос довольно простой: https//ru.wikipedia.org/wiki/Юлианская_дата , инфы более чем достаточно для вычисления. Всё остальное - тень на плетень.
Carthago restituenda est

Оффлайн gasha

  • ****
  • Сообщений: 350
  • Благодарностей: 17
  • 61:45:28,6 N 34:21:39.5 E (134 + 15 этаж) м
    • Сообщения от gasha
    • Siä Karjalassa
Из книги Монтенбрука и Пфлегера "Астрономия на персональном компьютере".

//------------------------------------------------------------------------------
//
// Ddd: Conversion of angular degrees, minutes and seconds of arc to decimal
//   representation of an angle
//
// Input:
//
//   D        Angular degrees
//   M        Minutes of arc
//   S        Seconds of arc
//
// <return>:  Angle in decimal representation
//
//------------------------------------------------------------------------------
double Ddd(int D, int M, double S)
{
  //
  // Variables
  //
  double sign;


  if ( (D<0) || (M<0) || (S<0) ) sign = -1.0; else sign = 1.0;
   
  return  sign * ( fabs((double)D)+fabs((double)M)/60.0+fabs(S)/3600.0 );
}

//------------------------------------------------------------------------------
//
// Mjd: Modified Julian Date from calendar date and time
//
// Input:
//
//   Year      Calendar date components
//   Month
//   Day
//   Hour      Time components (optional)
//   Min
//   Sec
//
// <return>:   Modified Julian Date
//
//------------------------------------------------------------------------------
double Mjd ( int Year, int Month, int Day,
             int Hour, int Min, double Sec )
{
  //
  // Variables
  //
  long    MjdMidnight;
  double  FracOfDay;
  int     b;


  if (Month<=2) { Month+=12; --Year;}
 
  if ( (10000L*Year+100L*Month+Day) <= 15821004L )
    b = -2 + ((Year+4716)/4) - 1179;     // Julian calendar
  else
    b = (Year/400)-(Year/100)+(Year/4);  // Gregorian calendar
   
  MjdMidnight = 365L*Year - 679004L + b + int(30.6001*(Month+1)) + Day;
  FracOfDay   = Ddd(Hour,Min,Sec) / 24.0;

  return MjdMidnight + FracOfDay;
}


// DE-0431LE-0431      -13200-AUG-15 00:00  JD  -3100015.50 to  17191-MAR-15 00:00  JD   8000016.50
cout << Mjd ( -13200, 8, 15, 0, 0, 0 ) + 2400000.5 << endl; // Ответ: -3100015.5
cout << Mjd (  17191, 3, 15, 0, 0, 0 ) + 2400000.5 << endl; // Ответ:   8000016.5
БПЦ 15х50, Nikon Aculon 7x50, Celestron Advanced VX 8" N, Sky-Watcher BK 909AZ3, ТАЛ-65, Таир-3ФС, Canon EOS 60D, Sony Alpha NEX-3.

Оффлайн gasha

  • ****
  • Сообщений: 350
  • Благодарностей: 17
  • 61:45:28,6 N 34:21:39.5 E (134 + 15 этаж) м
    • Сообщения от gasha
    • Siä Karjalassa
Вот старая статья Hatcher, D. A. Generalized Equations for Julian Day Numbers and Calendar Dates

Также PDF прикреплён к сообщению.

В ней приводится расширенный вариант расчёта Юлианской даты для разных календарей. О статье узнал из Трудов ИПА РАН, выпуск 10, 2004.

void Time::calcJD(bool cal) // https://web.archive.org/web/https://articles.adsabs.harvard.edu/pdf/1985QJRAS..26..151H
{
    // For Gregorian calendar
    int Ji = 1721426;
    int y = 4716;
    int j = 1401;
    int m = 3;
    int n = 12;
    int i = 1;
    int q = 0;
    int r = 4;
    int p = 1461;
    int v = 3;
    int u = 5;
    int s = 153;
    int t = 2;
    int w = 2;

    int Ys = year + y - (int)((n+m-i-month)/n);
    int g = (int)((int)((Ys+184)/100)*3/4)-38;
    int Ms = (month-m) % n;
    int Ds = day-i;
    jdn = (int)((p*Ys+q)/r) + (int)((s*Ms+t)/u) + Ds-j-g-0.5;
}

Из файла header.441
JPL Planetary Ephemeris DE441/LE441
Start Epoch: JED= -3100015.5-13200-AUG-15 00:00:00
Final Epoch: JED=  8000016.5 17191-MAR-15 00:00:00

На практике выходит
-13200-8-15 0:0:0.0, Gregorian JD: -3099915.5
 17191-3-15 0:0:0.0, Gregorian JD:  8000016.5

То есть в + работает нормально, а в - ошибка ровно в 100 дней. Хз)

P.S. Сейчас на некоторые зарубежные сайты вход закрыт. Один вариантов обхода это использование https://web.archive.org/web/. После web/ добавляем свою ссылку. Как раз верхняя ссылка так и сформирована. Да, ссылка будет устаревшей, но для доступа к архивам старых статей самое то.

P.P.S. Разница 100 дней, это разница между Юлианским и Григорианским календарями на тот момент. Если добавить условие,
if ((10000L*year+100L*month+day) > 15821004L) jdn -= g;  // Gregorian calendar
то всё встаёт на свои места.

Итоговый вариант:
void Time::calcJD(bool cal) // https://web.archive.org/web/https://articles.adsabs.harvard.edu/pdf/1985QJRAS..26..151H
{
    // For Gregorian calendar
    int Ji = 1721426;
    int y = 4716;
    int j = 1401;
    int m = 3;
    int n = 12;
    int i = 1;
    int q = 0;
    int r = 4;
    int p = 1461;
    int v = 3;
    int u = 5;
    int s = 153;
    int t = 2;
    int w = 2;

    int Ys = year + y - (int)((n+m-i-month)/n);
    int g = (int)((int)((Ys+184)/100)*3/4)-38;
    int Ms = (month-m) % n;
    int Ds = day-i;
    jdn = (int)((p*Ys+q)/r) + (int)((s*Ms+t)/u) + Ds-j-0.5;

    if ((10000L*year+100L*month+day) > 15821004L) jdn -= g;  // Gregorian calendar
}

P.P.S. Условие для григорианского календаря я взял из Монтенбрука. Он ссылается на Sky & Telescope, April 1981. Там в статье сказано, что формулы взяты из Jean Meeus Astronomical formulae for calculators
« Последнее редактирование: 28 Дек 2022 [04:27:53] от gasha »
БПЦ 15х50, Nikon Aculon 7x50, Celestron Advanced VX 8" N, Sky-Watcher BK 909AZ3, ТАЛ-65, Таир-3ФС, Canon EOS 60D, Sony Alpha NEX-3.

Оффлайн gasha

  • ****
  • Сообщений: 350
  • Благодарностей: 17
  • 61:45:28,6 N 34:21:39.5 E (134 + 15 этаж) м
    • Сообщения от gasha
    • Siä Karjalassa
Вот старая статья Hatcher, D. A. Generalized Equations for Julian Day Numbers and Calendar Dates

Также PDF прикреплён к сообщению.

В ней приводится расширенный вариант расчёта Юлианской даты для разных календарей. О статье узнал из Трудов ИПА РАН, выпуск 10, 2004.

void Time::calcJD(bool cal) // https://web.archive.org/web/https://articles.adsabs.harvard.edu/pdf/1985QJRAS..26..151H
{
    // For Gregorian calendar
    int Ji = 1721426;
    int y = 4716;
    int j = 1401;
    int m = 3;
    int n = 12;
    int i = 1;
    int q = 0;
    int r = 4;
    int p = 1461;
    int v = 3;
    int u = 5;
    int s = 153;
    int t = 2;
    int w = 2;

    int Ys = year + y - (int)((n+m-i-month)/n);
    int g = (int)((int)((Ys+184)/100)*3/4)-38;
    int Ms = (month-m) % n;
    int Ds = day-i;
    jdn = (int)((p*Ys+q)/r) + (int)((s*Ms+t)/u) + Ds-j-g-0.5;
}

Из файла header.441
JPL Planetary Ephemeris DE441/LE441
Start Epoch: JED= -3100015.5-13200-AUG-15 00:00:00
Final Epoch: JED=  8000016.5 17191-MAR-15 00:00:00

На практике выходит
-13200-8-15 0:0:0.0, Gregorian JD: -3099915.5
 17191-3-15 0:0:0.0, Gregorian JD:  8000016.5

То есть в + работает нормально, а в - ошибка ровно в 100 дней. Хз)

P.S. Сейчас на некоторые зарубежные сайты вход закрыт. Один вариантов обхода это использование https://web.archive.org/web/. После web/ добавляем свою ссылку. Как раз верхняя ссылка так и сформирована. Да, ссылка будет устаревшей, но для доступа к архивам старых статей самое то.

P.P.S. Разница 100 дней, это разница между Юлианским и Григорианским календарями на тот момент. Если добавить условие,
if ((10000L*year+100L*month+day) > 15821004L) jdn -= g;  // Gregorian calendar
то всё встаёт на свои места.

Итоговый вариант:
void Time::calcJD(bool cal) // https://web.archive.org/web/https://articles.adsabs.harvard.edu/pdf/1985QJRAS..26..151H
{
    // For Gregorian calendar
    int Ji = 1721426;
    int y = 4716;
    int j = 1401;
    int m = 3;
    int n = 12;
    int i = 1;
    int q = 0;
    int r = 4;
    int p = 1461;
    int v = 3;
    int u = 5;
    int s = 153;
    int t = 2;
    int w = 2;

    int Ys = year + y - (int)((n+m-i-month)/n);
    int g = (int)((int)((Ys+184)/100)*3/4)-38;
    int Ms = (month-m) % n;
    int Ds = day-i;
    jdn = (int)((p*Ys+q)/r) + (int)((s*Ms+t)/u) + Ds-j-0.5;

    if ((10000L*year+100L*month+day) > 15821004L) jdn -= g;  // Gregorian calendar
}

P.P.S. Условие для григорианского календаря я взял из Монтенбрука. Он ссылается на Sky & Telescope, April 1981. Там в статье сказано, что формулы взяты из Jean Meeus Astronomical formulae for calculators

P.P.P.S. Оба метода дают, наверное, верный расчёт юлинаской даты в отрицательную сторону. Но оба они обратно конвертировать в календарную (при заявленных условиях) не могут.
БПЦ 15х50, Nikon Aculon 7x50, Celestron Advanced VX 8" N, Sky-Watcher BK 909AZ3, ТАЛ-65, Таир-3ФС, Canon EOS 60D, Sony Alpha NEX-3.

Оффлайн Александр Вольф

  • *****
  • Сообщений: 3 325
  • Благодарностей: 103
  • Stellarium Developer
    • Skype - alex.v.wolf
    • Jabber - alex.wolf@jabber.ru
    • DeepSkyHosting: alexwolf
    • Сообщения от Александр Вольф
    • 47 Tucanae
Ну так не удивительно, для таких больших значений чисел типа int при расчетах уже недостаточно. Если используется long для получения JD из календарной, то он же должен использоваться и для обратного преобразования.
С уважением, Александр
Астротоп | Stellarium: donate | KStars
SW ED80/SW AllView GOTO | Celestron 15x70 | Celestron 25-125x80 | Veber 25x100

Оффлайн gasha

  • ****
  • Сообщений: 350
  • Благодарностей: 17
  • 61:45:28,6 N 34:21:39.5 E (134 + 15 этаж) м
    • Сообщения от gasha
    • Siä Karjalassa
Итак. Найдена золотая середина. Алгоритм, работающий для любых дат и календарей (Юлианский обычный и пролептический, Григорианский обычный и пролептический). Туда и обратно.
Источник алгоритма.
В документе на странице 199 находим алгоритм 199 (ALGORITHM 199 CONVERSIONS BETWEEN CALENDAR DATE AND JULIAN DAY NUMBER ROBERT G. TANTZEN).
Он, правда, написан на Алголе, но не беда. В целом всё понятно.

Также в сети найдена его реализация на Python. Про лицензию ничего не сказано, будем считать, что Public domain.
Я преобразовал его в C++.

Вот пример расчёта:

/Astronomy/Ephemeris/GashaSoft/libMeneldil/JD/Tantzen$ ./a.out
Julian->Gregorian: 2022-12-25 -> 2023-1-7 Diff = 13
Gregorian->Julian: 2023-1-7 -> 2022-12-25 Diff = -13

Julian->JD -13200-8-15 -> -3100015.500000
JD->Julian -3100015.500000 -> -13200-8-15

Gregorian->JD 17191-3-15 -> 8000016.500000
JD->Gregorian 8000016.500000 -> 17191-3-15

Результаты можно сравнивать с JPL Калькулятором, но он автоматически переключается между календарями в 1582 году. Пределы -9999 - +9999
Или же с формой автора кода на Python. Ссылка на калькулятор. Ограничений нет. Наверное  :)

Во вложении исходный код функций и простенький пример. Файл julian_date.cpp.txt

P.S. В Stellarium переделанный код из книги Numerical recipes. Который на основе алгоритмов Hatcher & Richardson (кто из них первый не вникал). И код из Stellarium также работает без видимых ограничений.

БПЦ 15х50, Nikon Aculon 7x50, Celestron Advanced VX 8" N, Sky-Watcher BK 909AZ3, ТАЛ-65, Таир-3ФС, Canon EOS 60D, Sony Alpha NEX-3.

Оффлайн gasha

  • ****
  • Сообщений: 350
  • Благодарностей: 17
  • 61:45:28,6 N 34:21:39.5 E (134 + 15 этаж) м
    • Сообщения от gasha
    • Siä Karjalassa
P.S. Данный пост отменяет предыдущий пост.

Зацепила тема)

Источник данных.

Далее привожу алгоритмы, которые, вроде как, валидны в интервале от минус 1400000 (1 миллион четыреста тысяч) до плюс 1400000 года.
Чёткая конвертация туда и обратно, включая на миллионы лет. Что более, чем достаточно.

Код на C++. В стандартной библиотеке нет подходящей функции для определения неполного частного и остатка.
Поэтому, я определил её так:

void divmod(int a, int b, int& q, int& r)
{
    /*
     * [url]https://skysmart.ru/articles/mathematic/delenie-chisel-s-ostatkom[/url]
     * Теорема
     * a = b · q + r, где
     * a — делимое,
     * b — делитель,
     * q — неполное частное,
     * r — остаток. 0 ⩽ r < |b|
     */

    if (a%b==0) // Если делится без остатка
    {
        q = a/b; // Частное // quotient
        r = 0; // Остаток // remainder
    }
    else
    {
        if (a>0 && b>0)
        {
            q = static_cast<int>(a/b);
            r = a-b*q;
            if (r<0 && r>=abs(b)) std::cerr << " Error with rule r>=0 && r<abs(b)\n";
            if (a!=(b*q+r)) std::cerr << " Error with a>0 && b>0\n";
        }
        else if (a<0 && b>0)
        {
            q = -1*static_cast<int>(abs(a)/abs(b))-1;
            r = a-b*q;
            if (r<0 && r>=abs(b)) std::cerr << " Error with rule r>=0 && r<abs(b)\n";
            if (a!=(b*q+r)) std::cerr << " Error with a<0 && b>0";
        }
        else if (a>0 && b<0)
        {
            q = -1*static_cast<int>(abs(a)/abs(b));
            r = a-b*q;
            if (r<0 && r>=abs(b)) std::cerr << " Error with rule r>=0 && r<abs(b)\n";
            if (a!=(b*q+r)) std::cerr << " Error with a>0 && b<0";
        }
        else if (a<0 && b<0)
        {
            q = -1*static_cast<int>(abs(a)/abs(b))+1;
            r = a-b*q;
            if (r<0 && r>=abs(b)) std::cerr << " Error with rule r>=0 && r<abs(b)\n";
            if (a!=(b*q+r)) std::cerr << " Error with a<0 && b<0";
        }
    }
}

Далее непосредственно расчёт.

// From Gregorian Date to CJDN
    int a,m,a1,m1,c1,a2,d,J;
    a = -1000000;
    m = 1;
    d = 1;
    divmod((m-3), 12, a1, m1);
    divmod((a+a1), 100, c1, a2);
    J = floor(146097*c1/4.0);
    J += floor(36525*a2/100);
    J += floor((153*m1+2)/5.0);
    J += d+1721119;
    std::cout << std::fixed << " J= " << J << "\n";

    // From CJDN to Gregorian Date
    int e1, e2, e3, m2;
    divmod((4*J-6884477), 146097, c1, e1);
    divmod((100*floor(e1/4)+99), 36525, a1, e2);
    divmod((5*floor(e2/100)+2), 153, m1, e3);
    divmod((m1+2), 12, a2, m2);
    a = 100*c1+a1+a2;
    m=m2+1;
    d = floor(e3/5.0)+1;
    std::cout << a << "-" << m << "-" << d << "\n";

    // From Julian Date to CJDN
    divmod((m-3), 12, a1, m1);
    J = floor(1461*(a+a1)/4.0);
    J += floor((153*m1+2)/5.0);
    J += d + 1721117;
    std::cout << std::fixed << " J= " << J << "\n";

    // From CJDN to Julian Date
    int m0, alpha1;
    divmod((4*J-6884469), 1461, a1, e1);
    divmod((5*floor(e1/4.0)+2), 153, m1, e2);
    divmod((m1+2), 12, alpha1, m0);
    a = a1+alpha1;
    m = m0+1;
    d = floor(e2/5.0)+1;
    std::cout << a << "-" << m << "-" << d << "\n";

Далее видим, что минус миллионный год пролептического григорианского календаря конвертируется в юлианскую дату и обратно. Аналогично для пролептического юлианского календаря.
Проверка сайта [url]https://www.aa.quae.nl/en/reken/juliaansedag.html[/url]
 J= -363521440
-1000000-1-1
 J= -363528942
-1000000-1-1

Вроде бы, косяков нет. Надеюсь, что вопрос закрыт).
БПЦ 15х50, Nikon Aculon 7x50, Celestron Advanced VX 8" N, Sky-Watcher BK 909AZ3, ТАЛ-65, Таир-3ФС, Canon EOS 60D, Sony Alpha NEX-3.