Загрузите дорожку судна

Загрузите astKonTikiData.mat файл. Это содержит траекторию поставки, скорость и курс для этого примера. Это хранилища файлов, широта и долгота для различной дорожки указывает в переменных "lat" и "долго", соответственно. Переменные содержат достаточно данных для одной точки дорожки в день от порта Кальяо к атоллу Рэроия. Кроме того, этот файл также хранит значения для скорости судна каждого дня в узлах в день "V" и курс в градусе "T".

load astKonTikiData

Создайте структуру наблюдения

Навигационный процесс сокращения является серией шагов, которые навигатор выполняет, чтобы определить широту и долготу его судна. Это основано на теории, описанной в американском Практическом Навигаторе [1], Навигационный Альманах [2], и Объяснительное Дополнение к Астрономическому Альманаху [3]. Процесс использует наблюдательные данные, полученные из секстанта, часов, компаса и навигационных графиков. Это возвращает прерывание (p) и истинный азимут (Z) для каждого из наблюдаемых объектов. Этот пример использует массив структур наблюдения, obs, чтобы содержать наблюдательные данные. Поля для массива структур:

Для простоты примите, что все измерения проведены в том же местоположении в лодке, с тем же секстантом, при той же температуре окружающей среды и давлении:

obs.h = 4;
obs.IC = 0;
obs.P = 982;
obs.T = 15;

Из экспедиции отбывают 28-го апреля 1947. В результате инициализируйте структуру для наблюдения для этой даты:

obs.year = 1947;
obs.month = 4;
obs.day = 28;

Инициализируйте процесс точного расчета для навигации

Чтобы запустить процесс точного расчета, задайте начальные условия для положения судна. Сохраните широту для фиксированного решения для широты и долготы в latFix и longFix переменных, соответственно. В этом примере, для первого местоположения фиксации, используют широту и долготу Кальяо, Перу:

longFix = zeros(size(long));
latFix = zeros(size(lat));
longFix(1) = long(1);
latFix(1) = lat(1);

Ежедневный точный расчет

В данном примере примите, что фиксация ежедневно получается с помощью наблюдательных данных. Поэтому этот пример использует "цикл for" в каждом наблюдении. Переменная "m" действует как счетчик, представляющий каждый прошедший день начиная с отклонения от порта:

for m = 1:size(lat,1)-1

Постепенно увеличьте день и внесите дневные корректировки в течение месяцев июня и апреля, оба из которых имеют только 30 дней:

    obs.day = obs.day + 1;
    [obs.month,obs.day] = astHelperDayCheck(obs.year,obs.month,obs.day);

Фактическая широта и долгота

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

    longActual = long(m+1);
    latActual = lat(m+1);

Выбор планеты

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

    if longActual>-90
        obs.object = {'Saturn';'Neptune'};
    elseif longActual<=-90 && longActual>-95
        obs.object = {'Saturn';'Neptune';'Jupiter'};
    elseif longActual<=-95
        obs.object = {'Neptune';'Jupiter'};
    end

Вычисление времени UTC

Настройте местное время к UTC в зависимости от принятой долготы. В данном примере примите, что все наблюдения берутся одновременно каждый день в 20:00 по местному времени.

    obs.hour = 20;

Для процесса точного расчета обновите структуру наблюдения с оценкой текущего местоположения. В этом случае местоположение оценивается с помощью предыдущей фиксации, скорость судна V и курс T.

    obs.longitude = longFix(m)+(1/60)*V(m)*sind(T(m))/cosd(latFix(m));
    obs.latitude = latFix(m)+(1/60)*V(m)*cosd(T(m));

Настройте местное время к UTC с помощью функции помощника astHelperLongitudeHour. Эта функция настраивает время наблюдения UTC в зависимости от предполагаемой долготы судна.

    obs.UTC = astHelperLongitudeHour(obs);

Высотный расчет секстанта

Для каждой из планет astHelperNauticalCalculation функция помощника вычисляет измерение секстанта, которое измерила бы команда в Коне-Тики. Это модели функции фактическое поведение планет, при компенсации локальных условий. Эта функция использует планетарную эфемериду и ECI к матрице преобразования ECEF. Анализ не включает планетарную аберрацию, гравитационное легкое отклонение и аберрацию легких явлений.

    obs.Hs = astHelperNauticalCalculation(obs,latActual,longActual);

Вычислите положение

Следующие вычисления заменяют использование навигационного альманаха. Они включают использование планетарного ephemerides и ECI к матрице преобразования ECEF.

Инициализируйте наклон, Гринвичский угол часа (GHA) и высоту для наблюдаемого объекта:

    obs.declination = zeros(size(obs.Hs));
    obs.GHA = zeros(size(obs.Hs));
    obs.altitude = zeros(size(obs.Hs));

Вычислите модифицированную дату Джулиана на время измерения:

    mjd = mjuliandate(obs.UTC);

Вычислите различие между UT1 и UTC:

    dUT1 = deltaUT1(mjd,'Action','None');

Вычислите ECI к матрице преобразования ECEF с помощью значений TAI-UTC (DAT) из американской Военно-морской Обсерватории:

    dAT = 1.4228180;
    TM = dcmeci2ecef('IAU-76/FK5',obs.UTC,dAT,dUT1);

Вычислите Юлианскую дату на Наземное Время, чтобы аппроксимировать Барицентрическое Динамическое Время:

    jdTT = juliandate(obs.UTC)+(dAT+32.184)/86400;

Вычислите наклон, Гринвичский Угол Часа и высоту для каждого из астрономических объектов:

    for k =1:length(obs.object)

Вычислите положение ECI для каждой планеты:

        posECI = planetEphemeris(jdTT,'Earth',obs.object{k},'405','km');

Вычислите положение ECEF каждой планеты:

        posECEF = TM*posECI';

Вычислите Гринвичский угол часа (GHA) и наклон с помощью положения ECEF:

        obs.GHA(k) = -atan2d(posECEF(2),posECEF(1));
        obs.declination(k) = atan2d(posECEF(3),sqrt(posECEF(1)^2 + ...
            posECEF(2)^2));

Вычислите расстояние от поверхности Земли к центру планеты с помощью ECEF для функции преобразования LLA:

        posLLA = ecef2lla(1000*posECEF');
        obs.altitude(k) = posLLA(3);

    end

Сокращение вида для планетарных объектов

Уменьшайте вид для каждой из планет, заданных в массиве структур наблюдения:

    [p,Z] = astHelperNauticalReduction(obs);

Вычислите шаг к широте и долготе от последней фиксации для текущей фиксации с помощью следующих уравнений. Эти уравнения основаны на Навигационном Альманахе.

    Ap = sum(cosd(Z).^2);
    Bp = sum(cosd(Z).*sind(Z));
    Cp = sum(sind(Z).^2);
    Dp = sum(p.*cosd(Z));
    Ep = sum(p.*sind(Z));
    Gp = Ap*Cp-Bp^2;

Вычислите шаг для широты и долготы согласно сокращению:

    deltaLongFix = (Ap*Ep-Bp*Dp)/(Gp*cosd(latFix(m)));
    deltaLatFix = (Cp*Dp-Bp*Ep)/Gp;

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

    longFix(m+1) = obs.longitude + deltaLongFix;
    latFix(m+1) = obs.latitude + deltaLatFix;

end

Визуализация решения для навигации

Следующий рисунок показывает фактическую дорожку и решения для сокращения вида:

astHelperVisualization(long,lat,longFix,latFix,'Plot')

Можно использовать Mapping Toolbox, чтобы получить более подробный график, изображающий решения для сокращения вида с американским континентом и Французской Полинезией.

astHelperVisualization(long,lat,longFix,latFix,'Map')

Относительная погрешность долготы и широты накапливается, когда судно приплывает от Кальяо до Raroia. Эта ошибка происходит из-за небольших ошибок в измерении высоты секстанта. На 9-е июня метод сокращения вычисляет истинный азимут (Z) для Нептуна и Юпитер. Истинные азимуты для Нептуна и Юпитер близко к на расстоянии в 180 градусов. Это различие вызывает небольшой пик в относительной погрешности. Эта ошибка, однако, все еще в ошибочных контурах метода сокращения.

Ссылки

[1] Bowditch, N. Американский практический навигатор. Национальная геопространственная спецслужба, 2012.

[2] Гидрографическая служба Соединенного Королевства. Навигационный альманах 2 012 коммерческих выпусков. Paradise Publications, Inc. 2011.

[3] Городской, Шон Э. и П. Кеннет Сейделманн. Объяснительное дополнение к астрономическому альманаху. 3-й Эд., университетские книги по науке, 2013.

[4] Военно-морская Обсерватория Соединенных Штатов. https://www.usno.navy.mil.