Положение Солнца, уравнение времени: НЕПРАВИЛЬНО К 27 ДНЯМ

Question

Я пытаюсь улучшить алгоритм Военно-морской обсерватории США для расчета положения солнца , написав его более простым способом, где все используемые числа были правильно идентифицированы - чтобы нормальные люди могли понятьэто также.

Но так или иначе, уравнение времени отключено на 27 дней.Может быть, кто-то здесь может определить, что не так?

Тестовый прогон:

1) Failure:
test_equation_of_time(TestSolarCalculations):
Expected 2011-10-01 10:23:00 UTC, not 2011-10-28 10:59:31 UTC.

2) Failure:
test_suns_declination(TestSolarCalculations):
Expected -3.18, not -3.2087920449753007.

Новый алгоритм:

# Constants for J2000.0 / 1 Jan 2000 12:00Z:
#
EPOCHS_JULIAN_DATE = 2451545
ANOMALISTIC_YEAR_IN_DAYS = 365.259636
TROPICAL_YEAR_IN_DAYS = 365.2421897
SUNS_MEAN_ANOMALY_AT_EPOCH = degrees_to_radians 357.5291
SUNS_MEAN_LONGITUDE_AT_EPOCH = degrees_to_radians 280.459
SUNS_GEODETIC_PRECESSION = degrees_to_radians 1.915
EARTHS_ORBITAL_ECCENTRICITY = 0.020
EARTHS_ADJUSTED_AVERAGE_RADIUS_IN_AU = 1.00014
EARTHS_APPROXIMATE_ATMOSPHERIC_REFRACTION = degrees_to_radians 0.01671
EARTHS_ECLIPTIC_MEAN_OBLIQUITY = degrees_to_radians 23.439
EARTHS_ECLIPTIC_OBLIQUITY_CHANGE_RATE = degrees_to_radians 0.00000036

def days_from_epoch
  @todays_julian_date - EPOCHS_JULIAN_DATE
end

def suns_daily_mean_anomaly_change_rate
  degrees_to_radians(360 / ANOMALISTIC_YEAR_IN_DAYS)
end

def suns_mean_anomaly
  SUNS_MEAN_ANOMALY_AT_EPOCH + suns_daily_mean_anomaly_change_rate * days_from_epoch
end

def suns_daily_mean_longitude_change_rate
  degrees_to_radians(360 / TROPICAL_YEAR_IN_DAYS)
end

def suns_mean_longitude
  SUNS_MEAN_LONGITUDE_AT_EPOCH + suns_daily_mean_longitude_change_rate * days_from_epoch
end

def suns_apparent_ecliptic_longitude
  suns_mean_longitude + SUNS_GEODETIC_PRECESSION * sin(suns_mean_anomaly) + EARTHS_ORBITAL_ECCENTRICITY * sin(2 * suns_mean_anomaly)
end

def suns_distance_from_earth_in_au
  EARTHS_ADJUSTED_AVERAGE_RADIUS_IN_AU - EARTHS_APPROXIMATE_ATMOSPHERIC_REFRACTION * cos(suns_mean_anomaly) - (EARTHS_APPROXIMATE_ATMOSPHERIC_REFRACTION ^ 2 / 2) * cos(2 * suns_mean_anomaly)
end

def earths_ecliptic_mean_obliquity
  EARTHS_ECLIPTIC_MEAN_OBLIQUITY - EARTHS_ECLIPTIC_OBLIQUITY_CHANGE_RATE * days_from_epoch
end

def suns_right_ascension
  atan2(cos(suns_apparent_ecliptic_longitude) * sin(suns_apparent_ecliptic_longitude), cos(suns_apparent_ecliptic_longitude)) / 15
end

# Time.utc(Time.now.year) => 2011-01-01 00:00:00 UTC
#
def equation_of_time
  Time.utc(Time.now.year) + (radians_to_degrees(suns_mean_longitude) / 15 - suns_right_ascension) * 60 * 60 * 24
end

def suns_declination
  radians_to_degrees(asin(sin(earths_ecliptic_mean_obliquity) * sin(suns_apparent_ecliptic_longitude)))
end

Спасибо!

Матс

Answer 1

Здесь:

def suns_apparent_ecliptic_longitude
  suns_mean_longitude 
  + SUNS_GEODETIC_PRECESSION * sin(suns_mean_anomaly) 
  + EARTHS_ORBITAL_ECCENTRICITY * sin(2 * suns_mean_anomaly)
end

Вы смешиваете единицы. suns_mean_longitude - это радианы, как и SUNS_GEODETIC_PRECESSION (и, следовательно, так же, как и второе слагаемое), но EARTHS_ORBITAL_ECCENTRICITY - это 0.020 градусов .

Дополнительно:

EARTHS_APPROXIMATE_ATMOSPHERIC_REFRACTION = degrees_to_radians 0.01671

должно быть просто

EARTHS_APPROXIMATE_ATMOSPHERIC_REFRACTION = 0.01671

Я думаю, поскольку это используется в формуле, которая определяет расстояние .

Дополнительно дополнительно

def suns_right_ascension
  atan2(cos(suns_apparent_ecliptic_longitude) * sin(suns_apparent_ecliptic_longitude), cos(suns_apparent_ecliptic_longitude)) / 15
end

Я не проверял, какие единицы возвращает эта atan2 функция - вы уверены, что это градусы, которые должны были бы быть для деления на 15, чтобы иметь смысл?

Могут быть и другие проблемы; это только те, которые я вижу.

Answer 2

Часть моего кода.

# From angles.rb:<br>
# eccentricity of elliptical Earth orbit around Sun
# Horner calculation method  
def eccentricity_Earth( ta = A2000 )
  ta = check_jct_zero( ta )      
  # 0.016708617 - ta[ 0 ] * ( 0.000042037 + ta[ 0 ] * 0.0000001235 )
  [-0.0000001235, -0.000042037, 0.016708617].inject(0.0) {|p, a| p * ta[0] + a}
end
alias_method :eccentricity_earth_orbit, :eccentricity_Earth

От: Уравнение времени Рубиновый камень

Где он используется?

# From angles.rb:<br>
# equation of centre
# added to mean anomaly to get true anomaly. 
def center( ta = A2000)
  ta = check_jct_zero( ta )    
  sine_1M = sin( 1.0 * deg_to_rad( @ma ) )
  sine_2M = sin( 2.0 * deg_to_rad( @ma ) )
  sine_3M = sin( 3.0 * deg_to_rad( @ma ) )
  sine_4M = sin( 4.0 * deg_to_rad( @ma ) )
  sine_5M = sin( 5.0 * deg_to_rad( @ma ) )
  e = eccentricity_Earth( ta )
  rad_to_deg( sine_1M * (     2.0  * e    - e**3/4.0 + 5/96.0 * e**5 ) +  
              sine_2M * (   5/4.0  * e**2 - 11/24.0 * e**4 )           + 
              sine_3M * ( 13/12.0  * e**3 - 43/64.0 * e**5 )           +
              sine_4M *  103/96.0  * e**4                              +
              sine_5M * 1097/960.0 * e**5                              )

  # sine_1M *( 1.914602 - ta[ 0 ] * ( 0.004817 + ta[ 0 ] * 0.000014 )) +
  # sine_2M *( 0.019993 - ta[ 0 ] * 0.000101 )                         +
  # sine_3M *  0.000289
end
alias_method :equation_of_center, :center

Answer 3

Может быть, это проблема точности чисел с плавающей запятой в Ruby, а не логика вашей программы? Я не использовал ruby, но, кажется, есть некоторые проблемы в операциях с плавающей запятой в ruby, согласно данным постам. Ошибка с точностью математических операций на Ruby

Согласно сообщению, вам может понадобиться использовать класс BigDecimal для операций с плавающей запятой. Я на самом деле не прошел через вашу программу, это всего лишь предположение.