function []=cmpnav(raw, sm, plot_title, varargin)
% function []=cmpnav(raw, sm, plot_title [, year_base, n ])
%
% Make plots of compared ship velocities. The file format must be the one used
% by 'putnav', i.e. ascii data with the following columns: 
%       
%   Column  
%      1     time in decimal days
%      4     u_ship in m/s
%      5     v_ship in m/s
%      6     latitude in degree
%      7     longitude in degree
%
% INPUT : raw       : name of file containing raw ship velocity data
%         sm        : name of file containing smoothed ship velocity data
%         plot_title: string with title for plot
%         
%         OPTIONS:
%   
%         year_base: year_base to display on the x-axis label (no default)
%         n        : number of points per plot (default is 50)
%
% OUTPUT: file 'cmpnav.ps' file with all plots.
%
%Time:<Last updated on 98/04/02 at 12:05:26 by jaccard@pc34.gfi.uib.no>
%File:</home1/jaccard/data1/codas3/gfi/matlab/cmpnav.m>

% Load raw data:
R = read_asc(raw);
R = tonan(R);

% Load smooth data:
S = read_asc(sm);
S = tonan(S);

% Definitions:
t = 1;
u = 4;
v = 5;
x = 6;
y = 7;

% Default optionnal arguments:
yb_str = '';
nx_max = 50;

% Process additionnal input arguments:
if nargin > 3,
  yb_str = sprintf('(%d)', varargin{1});
end;
if nargin > 4,
  nx_max = varargin{1};
end

% Divide the plot x-axis:
xmin = max( [ min(R(:,t)), min(S(:,t)) ] );
xmax = min( [ max(R(:,t)), max(S(:,t)) ] );
nx = min( [ size(R,1), size(S,1) ] );

% Determine number of plots:
np = round(nx/nx_max);
dn = round(nx/np);


% Create plots:
n0 = 1;
for i=1:np,

  figure;
  
  % Define plot indice:
  n1 = n0 + dn;
  if n1 > nx,
	 n1 = nx;
  end
  
  xmin = max( [ min(R(n0:n1,t)), min(S(n0:n1,t)) ] );
  xmax = min( [ max(R(n0:n1,t)), max(S(n0:n1,t)) ] );
  d = (xmax -xmin);
  xmin = xmin - d/50;
  xmax = xmax + d/50;

  % Plot u:
  subplot(411);
  hold on;
  plot(R(n0:n1,t), R(n0:n1,u), '+', 'Markersize', 4);
  plot(S(n0:n1,t), S(n0:n1,u), '-o', 'Markersize', 2);
  a = axis;
  axis([ xmin, xmax, a(3), a(4)]);
  ylabel('U (m/s)');
  
  % Plot title:
  title(sprintf('%s: Ship Speed [+:raw, -:smoothed]', plot_title));

  % Plot u:
  subplot(412);
  hold on;
  plot(R(n0:n1,t), R(n0:n1,v), '+', 'Markersize', 4);
  plot(S(n0:n1,t), S(n0:n1,v), '-o', 'Markersize', 2);
  a = axis;
  axis([ xmin, xmax, a(3), a(4)]);
  ylabel('V (m/s)');
  
  % Plot u:
  subplot(413);
  hold on;
  plot(R(n0:n1,t), R(n0:n1,x), '+', 'Markersize', 4);
  plot(S(n0:n1,t), S(n0:n1,x), '-');
  a = axis;
  axis([ xmin, xmax, a(3), a(4)]);
  ylabel('lon (deg)');
  
  % Plot u:
  subplot(414);
  hold on;
  plot(R(n0:n1,t), R(n0:n1,y), '+', 'Markersize', 4);
  plot(S(n0:n1,t), S(n0:n1,y), '-');
  a = axis;
  axis([ xmin, xmax, a(3), a(4)]);
  ylabel('lat (deg)');
  
  % x-axis label:
  xlabel(sprintf('decimal days %s', yb_str));
  
  % print the plot:
  if i==1,
	 print -dps cmpnav.ps
  else
	 print -dps -append cmpnav.ps
  end
  
  
  % Update plot indice:
  n0 = n1;
end


