Skip to content
Snippets Groups Projects
predictFromGps.m 3.5 KiB
Newer Older
function [gpsRafosDates, predictedToa, launchDateToa] = predictFromGps(pRfb, pSoundSource, pSoundspeed, pLeapsecondsMatrix)
%PREDICTFROMGPS Predicts TOA based on the given input parameter.
%   This is also known as ETA.

%% Collect Float data
gpsData = pRfb.SAT_DATA;
floatSchedule = pRfb.FLOAT.schedule;
floatPhasereftime = pRfb.FLOAT.phasereftime(1) + pRfb.FLOAT.phasereftime(2)/60;
floatWindowStart = pRfb.FLOAT.windowstart;
floatWindowLength = pRfb.FLOAT.windowlength;
floatLaunchTime = artoa.convert.dmy2rd( ...
    pRfb.FLOAT.launchtime(3), ...
    pRfb.FLOAT.launchtime(2), ...
    pRfb.FLOAT.launchtime(1) ...
) + pRfb.FLOAT.launchtime(4)/24 + pRfb.FLOAT.launchtime(5)/24/60;

%% Collect sound source data
soundSourceTime = pSoundSource.reftime(1) + pSoundSource.reftime(2)/60;
soundSourceEmissionBegin = artoa.convert.dmy2rd( ...
    [pSoundSource.begemis(3) pSoundSource.begemis(2) pSoundSource.begemis(1)] ...
) + pSoundSource.begemis(4)/24 + pSoundSource.begemis(5)/24/60;
soundSourcePosition = pSoundSource.position;

%% Collect generic variables
leapSeconds = artoa.data.calculateLeapseconds( ...
    floatLaunchTime, ...
    soundSourceEmissionBegin, ...
    pLeapsecondsMatrix ...
);

%% Calculations
loopLength = floor(24/floatSchedule);
floatTime = NaN(loopLength, 1);
for o = 1:loopLength
    floatTime(o) = mod(floatPhasereftime + (o-1) * floatSchedule, 24);
end
floatTime = sort(floatTime);

travelCorrection = soundSourceTime - floatTime;
maxTime = (floatWindowStart + floatWindowLength) / 60;
minTime = floatWindowStart / 60;
buffer = (maxTime - minTime) / 3;
indicesInRange = (travelCorrection <= maxTime + buffer) & (travelCorrection >= minTime - buffer);
windowStartTime = min(travelCorrection(indicesInRange));
clear maxTime minTime indicesInRange;

%% Predict toa for every gps position

gpsMeasurementCount = size(gpsData, 1);
predictedToa = NaN(gpsMeasurementCount, 1);
gpsRafosDates = NaN(gpsMeasurementCount, 1);
if (isempty(windowStartTime))
        return;
end

launchDateToa(1) = artoa.convert.dmy2rd( ...
    pRfb.FLOAT.launchtime(3), ...
    pRfb.FLOAT.launchtime(2), ...
    pRfb.FLOAT.launchtime(1) ...
    );
launchDateToa(2) = predictToa(pRfb.FLOAT.launchpos);
for o = 1:gpsMeasurementCount
    currentGpsPosition = [gpsData(o, pRfb.SAT_FORMAT.lat_sat), gpsData(o, pRfb.SAT_FORMAT.lon_sat)];
    % check if sat date is available
    if ~isnan(gpsData(o, pRfb.SAT_FORMAT.day_sat)) ...
            & ~isnan(gpsData(o, pRfb.SAT_FORMAT.month_sat)) ...
            & ~isnan(gpsData(o, pRfb.SAT_FORMAT.day_sat))
        gpsRafosDates(o) = artoa.convert.dmy2rd( ...
            gpsData(o, pRfb.SAT_FORMAT.day_sat), ...
            gpsData(o, pRfb.SAT_FORMAT.month_sat), ...
            gpsData(o, pRfb.SAT_FORMAT.year_sat) ...
        );
    else % get rtc date
        gpsRafosDates(o) = artoa.convert.dmy2rd( ...
            gpsData(o, pRfb.SAT_FORMAT.day_rtc), ...
            gpsData(o, pRfb.SAT_FORMAT.month_rtc), ...
            gpsData(o, pRfb.SAT_FORMAT.year_rtc) ...
        );
    end
    % calculate distance between sound source and gps
    function [predictedToa] = predictToa(pPosition)
        % calculate distance between sound source and gps
        distanceSoundSourceToFloat = artoa.data.calculateGeodist(pPosition, soundSourcePosition);

        predictedToa = windowStartTime * 3600 - leapSeconds ...
            + distanceSoundSourceToFloat / (pSoundspeed / 1000) - pRfb.FLOAT.windowstart * 60;
    end