Оцените аналемму Sun Используя планетарный Ephemerides и ECI к преобразованию AER

Открыть скрипт

В этом примере показано, как оценить аналемму Sun. Аналемма является кривой, которая представляет изменение углового смещения Sun от его средней позиции по астрономической сфере относительно определенной геолокации на Наземной поверхности. В этом примере аналемма оценивается относительно Королевской Обсерватории в Гринвиче, Соединенное Королевство. После оценки пример строит аналемму.

Этот пример использует данные, что можно загрузить использование aeroDataPackage команды.

Идентифицируйте Даты года, за Который можно Вычислить Аналемму Sun

Задайте даты, для которых можно вычислить аналемму. В этом примере эти даты располагаются в течение года 2014 с 1-го января до 31-го декабря в полдень UTC.

dv = datetime(2014,1,1:365,12,0,0);
dvUTC = [dv.Year' dv.Month' dv.Day' dv.Hour' dv.Minute' dv.Second'];

Вычислите позицию Sun

Используйте функцию planetEphemeris, чтобы вычислить позицию Sun. В этом примере:

Функция tdbjuliandate вычисляет дату Джулиана на динамическое барицентрическое время (TDB).
Функция tdbjuliandate требует наземного времени (TT).

Вычисление наземного времени в секундах от UTC требует различия во Всемирное координированное время (UTC) и Международное атомное время (TAI).

Для 2 014, этим различием (DAT) составляют 35 секунд.
Аппроксимированным наземным временем (secTT) является DAT + 32,184 секунды.
Наземное время в году, месяце, дне, часе, минутах и секундах содержится в dvTT массиве.

dAT = 35;
secTT = dAT + 32.184;
dvTT = dv + secTT/86400;

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

jdTDB = tdbjuliandate([dvTT.Year' dvTT.Month' dvTT.Day' dvTT.Hour' dvTT.Minute' dvTT.Second']);

Определите позицию Sun для этих дат:

posSun = planetEphemeris(jdTDB,'Earth','Sun')*1000;

Вычислите различие между UTC и основное всемирное время (UT1)

Вычислите различие между UTC и UT1, deltaUT1, с помощью модифицированных дат Джулиана UTC.

mjdUTC = mjuliandate(dvUTC);
dUT1 = deltaUT1(mjdUTC);

Вычислите полярное движение и смещение Астрономического промежуточного полюса (CIP)

Вычислите полярное движение и смещение CIP с помощью модифицированных дат Джулиана UTC.

PM = polarMotion(mjdUTC);
dCIP = deltaCIP(mjdUTC);

Задайте геопотенциальное расположение королевской обсерватории в Гринвиче, Соединенное Королевство

Задайте геопотенциальное местоположение для положения, по которому можно оценить аналемму. В этом примере это местоположение является широтой, долготой и высотой для Королевской Обсерватории в Гринвиче (51,48 градуса на север, 0,0015 градуса на запад, 0-метровая высота).

LLAGreenwich = [51.48,-0.0015,0];
aer = eci2aer(posSun,dvUTC,repmat(LLAGreenwich,length(jdTDB),1),...
              'deltaAT',dAT*ones(length(jdTDB),1),'deltaUT1',dUT1,...
              'PolarMotion',PM,'dCIP',dCIP);

Задайте дни в течение года после аналеммы, которую вы хотите построить

На аналемме можно построить дни интереса в течение года после аналеммы. Этот пример графики:

Первый день каждого месяца в 2 014.
Летние и зимние солнцестояния.
Равноденствия пружины и осени.

Получить первый день каждого месяца в 2 014:

aerFirstMonth = aer(dvUTC(:,3)==1,:);

Чтобы получить солнцестояния и равноденствия (для 2 014 3/20, 6/21, 9/22, 12/21):

solsticeEquinox = [ aer(and(dvUTC(:,2)==3,dvUTC(:,3)==20),1) aer(and(dvUTC(:,2)==3,dvUTC(:,3)==20),2); ...
                    aer(and(dvUTC(:,2)==6,dvUTC(:,3)==21),1) aer(and(dvUTC(:,2)==6,dvUTC(:,3)==21),2); ...
                    aer(and(dvUTC(:,2)==9,dvUTC(:,3)==22),1) aer(and(dvUTC(:,2)==9,dvUTC(:,3)==22),2); ...
                    aer(and(dvUTC(:,2)==12,dvUTC(:,3)==21),1) aer(and(dvUTC(:,2)==12,dvUTC(:,3)==21),2)];

Постройте результаты

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

for ii = 12:-1:1
	firstDays{ii} = [num2str(ii) '/' num2str(1)];
end

f = figure;
plot(aer(:,1),aer(:,2),'.',...
    solsticeEquinox(:,1),solsticeEquinox(:,2),'ks',...
    aerFirstMonth(:,1),aerFirstMonth(:,2),'ko',...
    'MarkerSize',8,'MarkerFaceColor','k');
title('Analemma observed at Greenwich Observatory');
xlabel('Azimuth [deg]');
ylabel('Elevation [deg]');
axis([176,185,10,70])
text(aerFirstMonth(:,1)+.1, aerFirstMonth(:,2)+1.2, firstDays, 'Color', 'k','HorizontalAlignment', 'Left');
text(solsticeEquinox(1,1)+.2, solsticeEquinox(1,2)-1.5, 'Spring Equinox', 'Color', 'k','HorizontalAlignment', 'Left');
text(solsticeEquinox(2,1), solsticeEquinox(2,2)+2.5, 'Summer Solstice', 'Color', 'k','HorizontalAlignment', 'Left');
text(solsticeEquinox(3,1)+.1, solsticeEquinox(3,2)+1.2, 'Fall Equinox', 'Color', 'k','HorizontalAlignment', 'Left');
text(solsticeEquinox(4,1)+.1, solsticeEquinox(4,2)-2.5, 'Winter Solstice', 'Color', 'k','HorizontalAlignment', 'Left');

Документация