tr_map.m 4.76 KB
function tr_map
% TR_MAP - Display and save a map showing drifter trajectories specified
% in the input file.
%
% Syntax:  tr_map
% 
% Inputs:  The function reads an ascii input file called infile, which
%          must be written according to the following rules:
%          - The first line indicates the number of trajectories or GPS
%            devices to plot.
%          - The second line indicates the name of the device (e.g. spot15).
%          - The third line indicates the name of the mat file containing
%            the trajectory data. This file must be located in a directory
%            called mat/, which is located where infile is.
%          - Subsequent lines indicate, in the same order, the name of the
%            device and the name the mat file. The number of devices must 
%            equal the number specified at the first line.
%
%
% Example of infile:
%
%          3
%          spot15
%          spot15_20140910a
%          spot17
%          spot17_20140912a
%          spot19
%          spot19_20140910a
%
% 
% Other m-files required: Image Processing Toolbox
% Subfunctions: none
% MAT-files required: none
% 
% See also: TR_MERGE, GPX2MAT
% 
% Author: Dany Dumont
% UQAR/ISMER
% email: dany_dumont@uqar.ca
% Website: http://www.ismer.ca/dumont-dany
% September 2014
% ______________________________________________________________________


curdir = pwd;
matdir = [curdir,'/mat'];

if ~exist(matdir)
    disp(' ERROR : The /mat directory does not exist. Make')
    disp('         sure mat files exist and that they     ')
    disp('         are placed in a directory named mat.   ')
    return
end

if exist('infile')
    fid = fopen('infile');
else
    disp(' ERROR : You need an input file named infile    ')
    disp('         which must be written according to the ')
    disp('         following rules:                       ')
    disp('         3                % # of trajectories   ')
    disp('         spot15           % device 1            ')
    disp('         spot15_20140910a % mat file of device 1')
    disp('         spot17           % device 2            ')
    disp('         spot17_20140912a % mat file of device 2')
    disp('         spot19           % device 3            ')
    disp('         spot19_20140910a % mat file of device 3')
    return
end

ndstr = fgetl(fid);
nd    = str2double(ndstr);

lonmin =  360;
lonmax = -360;
latmin =  90;
latmax = -90;

for n = 1:nd
    A(n).dev  = fgetl(fid);
    A(n).file = fgetl(fid);
    load([matdir,'/',A(n).file]);
    
    [~,ds] = legs(data(:,3),data(:,2));         % distance (nm)
    dt     = data(2:end,1) - data(1:end-1,1);   % time difference (days)
    v      = 1000.*nm2km(ds)./(dt.*86400);      % average speed (m s-1)
    travel = nm2km(sum(ds));                    % km
    
    A(n).data   = data;
    A(n).ds     = ds;
    A(n).dt     = dt;
    A(n).v      = v;
    A(n).travel = travel;
    
    lonmin = min(lonmin,min(A(n).data(:,2)));
    lonmax = max(lonmax,max(A(n).data(:,2)));
    latmin = min(latmin,min(A(n).data(:,3)));
    latmax = max(latmax,max(A(n).data(:,3)));
end

%% Compute and display relevant quantities

disp(' -------------- Summary -------------- ')
disp(' ')
for n = 1:nd
    disp([' ',A(n).dev]);
    disp(' ')
    disp(['   Number of points : ',num2str(length(A(n).data(:,2)))]);
    disp(['   Total travel     : ',num2str(A(n).travel),' km']);
    disp(['   Initial position : ',num2str(A(n).data(1,2)),'E, ',num2str(A(n).data(1,3)),'N']);
    disp(['   Final position   : ',num2str(A(n).data(end,2)),'E, ',num2str(A(n).data(end,3)),'N']);
    disp(['   Deployment time  : ',datestr(A(n).data(1,1))]);
    disp(['   Beach time       : ',datestr(A(n).data(end,1))]);
    disp(['   Total time       : ',num2str(A(n).data(end,1) - A(n).data(1,1)),' days']);
    disp(['   Average speed    : ',num2str(1000.*A(n).travel./(86400.*(A(n).data(end,1) - A(n).data(1,1)))),' m/s']);
    disp(' ')
end
disp(' -------------------------------------- ')


%% Create the map
% Track and land colors
trk(1,:,1) = [0.8 0.2 0.2];
trk(1,:,2) = [0.2 0.8 0.2];
trk(1,:,3) = [0.2 0.2 0.8];
trk(1,:,4) = [0.4 0.2 0.8];
trk(1,:,5) = [0.2 0.4 0.8];

%land_color = [241 235 144]./255;
land_color = [0.7 0.7 0.7];

% Close-up
tol = (lonmax - lonmin)./10;
lonmin = lonmin - tol;
lonmax = lonmax + tol;
latmin = latmin - tol;
latmax = latmax + tol;

figure(1)
m_proj('mercator', ...
    'longitudes',[lonmin lonmax], ...
    'latitudes',[latmin latmax]);

m_gshhs_f('patch',land_color);
hold on

for n = 1:nd
    h(n) = m_plot(A(n).data(:,2),A(n).data(:,3),'-o','LineWidth',2, ...
        'Color',trk(1,:,n), ...
        'MarkerSize',1, ...
        'MarkerFaceColor',trk(1,:,n));
end

m_grid('box','fancy','linestyle','none','fontsize',13,'fontname','Verdana');
hold off

print(gcf,'map.png','-dpng','-painters');
print(gcf,'map.eps','-depsc2','-painters');