Skip to content
Snippets Groups Projects
predictFromGps.m 3.98 KiB
Newer Older
function [gpsRafosDates, predictedToa] = predictFromGps(pRfb, pSoundSource, pSoundSpeedParameter, pLeapsecondsMatrix)
%ESTIMATETOAFROMGPS Summary of this function goes here
%   Detailed explanation goes here

%% 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;
% floatSurfaceEnd = artoa.convert.dmy2rd( ...
%     pRfb.FLOAT.cycle(13), pRfb.FLOAT.cycle(12), pRfb.FLOAT.cycle(11) ...
% ) + pRfb.FLOAT.cycle(14)/24 + pRfb.FLOAT.cycle(15)/24/60;

%% Collect sound source data
soundSourceTime = pSoundSource.reftime(1) + pSoundSource.reftime(2)/60;
%soundSourceSchedule = pSoundSource.schedule;
%soundSourceLaunchDate = [pSoundSource.launchtime(3) pSoundSource.launchtime(2) pSoundSource.launchtime(1)];
soundSourceEmissionBegin = artoa.convert.dmy2rd( ...
    [pSoundSource.begemis(3) pSoundSource.begemis(2) pSoundSource.begemis(1)] ...
) + pSoundSource.begemis(4)/24 + pSoundSource.begemis(5)/24/60;
%soundSourceOffset = 0;
% if (~isempty(pSoundSource.offset))
%     soundSourceOffset = pSoundSource.offset(4);
% end
soundSourceDrift = 0;
if ~isempty(pSoundSource.drift) && ~isnan(pSoundSource.drift)
    soundSourceDrift = pSoundSource.drift;
end
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;

% sound xmit time (seconds) = (hours*60 plus minutes)*60 -leapseconds - sosooffset
%soundXmitTime = windowStartTime * 3600 - leapSeconds + soundSourceOffset + pAdditionalOffset;

soundSourceDriftSinceLaunch = soundSourceDrift ...
    * (floatLaunchTime - soundSourceEmissionBegin);
%soundSourceDriftSinceSurfaceEnd = soundSourceDrift ...
%    * (floatSurfaceEnd - artoa.convert.dmy2rd(soundSourceEmissionBegin));

%predictedToaAtLaunch =  soundXmitTime + soundSourceDriftSinceLaunch;
%predictedToaAtSurface = soundXmitTime + soundSourceDriftSinceSurfaceEnd;


%% Predict toa for every gps position

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

for o = 1:gpsMeasurementCount
    currentGpsPosition = [gpsData(o, pRfb.SAT_FORMAT.lat_sat), gpsData(o, pRfb.SAT_FORMAT.lon_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) ...
    );
    currentSoundsourceDrift = soundSourceDrift * (gpsRafosDates(o) - soundSourceEmissionBegin);
    % calculate distance between sound source and gps
    distanceSoundSourceToFloat = artoa.data.calculateGeodist(currentGpsPosition, soundSourcePosition);
    
    predictedToa(o) = windowStartTime * 3600 - leapSeconds + currentSoundsourceDrift ...
        + distanceSoundSourceToFloat / (artoa.data.calculateSoundVelocity( ...
            pSoundSpeedParameter.temperature, pSoundSpeedParameter.pressure, ...
            pSoundSpeedParameter.method ...
        )/1000) - pRfb.FLOAT.windowstart * 60;