Ошибки в математике для вычисления поправок для параллакса небесного объекта в Java

Question

Выполняю расчеты книги Practical Astronomy with your Calculator or Spreadsheet. До сих пор мои расчеты дают тот же результат, что и примерные расчеты в книге.

Однако, достигнув §39 «Вычисление поправок для параллакса», я сталкиваюсь с различием, которое я не понимаю.

Задача под рукой описывается следующим образом:

В качестве примера вычислим кажущуюся прямое восхождение и склонение Луны 26 февраля 1979 года в 16ч 45м UT при наблюдении из места 60 метров над уровнем моря по долготе 100 ° з.д. и широте 50 ° Н. Геоцентрические координаты были α = 22h 35m 19s и δ = -7 ° 41 '13' ', а экваториальный горизонтальный параллакс Луны составлял 1 ° 01' 09 ' '.

В книге описана последовательность расчета следующим образом:

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

Я мог бы жить с этим, но там также является отличием на шаге 10. Мой результат -8,570634, результат книги -8.538165, довольно большая разница. Я снова и снова прочитал шаг 10, чтобы узнать, есть ли ошибка в моем коде, но я этого не вижу.

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

Моего Java-код для этой функции выглядит следующим образом:

static EquatorialCoordinate parallax(EquatorialCoordinate body, ObserverLocation observer, ZonedDateTime zdt, double P) { 
    double Hd = 15d * raha(body.α, zdt, observer.λ); 
    step("α", body.α); 
    step("δ", body.δ); 
    step("φ", observer.φ); 
    step("λ", observer.λ); 
    step("h", observer.h); 
    step("H", Hd); 

    double H = toRadians(Hd); 
    Parallax ρ = parallax(observer.φ, observer.h); 

    step("P", P); 
    P = toRadians(P); 

    double δ = toRadians(body.δ); 
    double r = 1d/sin(P); 
    step("r", r); 
    double ρsinφ = ρ.sin; 
    double ρcosφ = ρ.cos; 
    step("ρcosφ'",ρcosφ); 
    step("ρsinφ'",ρsinφ); 
    double Δ = atan((ρcosφ * sin(H))/((r * cos(δ)) - (ρcosφ * cos(H)))); 
    step("Δ", toDegrees(Δ)); 
    H += Δ; 
    step("H'", toDegrees(H)); 
    Δ = toDegrees(Δ); 
    double α$ = body.α - (Δ/15d); 
    step("α'", α$); 

    double divident = (r * sin(δ)) - ρsinφ; 
    double divisor = (r * cos(δ) * cos(H)) - ρcosφ; 
    double δ$ = atan(cos(H) * (divident/divisor)); 
    δ$ = toDegrees(δ$); 
    step("δ'", δ$); 
    return new EquatorialCoordinate(α$, δ$); 
} 

Функция «step» выполняет простой форматированный printf. Выход этой программы:

• α 22,588611
• δ -7,686944
• φ 50,000000
• λ -100,000000
• ч 60,000000
• Н -31,642500
• Р 1,019167
• r 56.221228
• ρcosφ '0.644060
• ρsinφ' 0,762422
• Δ -0,350915
• H '-31,993414
• α' 22,612005
• δ '-8,570634

Полученный δ' -8 ° 34' 14.28 "вместо -8 ° 32 '17"

Я заменил вычисленное значение H' значением книги, чтобы увидеть, содержит ли книга ошибку, но даже если я делаю это, это значение неверно.

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

(Edit :) Класс аннотируется с помощью strictfp, используя java.util.StrictMath.

Answer 1

Пишет

H += Δ; 

который изменяет значение H.

Тогда вы пишете

double δ$ = atan(cos(H) * (divident/divisor)); 

который использует новую версию H, когда она должна использовать старое значение.

Answer 2

Благодаря @svasa я обнаружил, что делитель шага 10 должен содержать H, а не H '. Правильный код:

static EquatorialCoordinate parallax(EquatorialCoordinate body, ObserverLocation observer, ZonedDateTime zdt, double P) { 
    double H = toRadians(15d * raha(body.α, zdt, observer.λ)); 
    P = toRadians(P); 
    Parallax ρ = parallax(observer.φ, observer.h); 
    double δ = toRadians(body.δ); 
    double r = 1d/sin(P); 
    double Δ = atan((ρ.cosφ * sin(H))/((r * cos(δ)) - (ρ.cosφ * cos(H)))); 
    double H$ = H + Δ; 
    double α$ = body.α - (toDegrees(Δ)/15d); 
    double δ$ = toDegrees(atan(cos(H$) * ((r * sin(δ) - ρ.sinφ)/(r * cos(δ) * cos(H) - ρ.cosφ)))); 

    return new EquatorialCoordinate(α$, δ$); 
}