function []=cmpnav(raw, sm, plot_title, varargin)
%  1.  NAME
%
%  cmpnav - Compare ship velocities
%
%  2.  SYNOPSIS
%
%  cmpnav(raw, sm, plot_title, year_base, n )
%
%  3.  DESCRIPTION
%
%  This Matlab function generates plots of compared ship velocities. The
%  of input files must be the one used by putnav, i.e. ASCII data with
%  the following columns:
%
%  1. Time in decimal days
%
%  2. Not used
%
%  3. Not used
%
%  4. U ship in m/s
%
%  5. V ship in m/s
%
%  6. Latitude in degree
%
%  7. Longitude in degree
%
%  4.  PARAMETERS
%
%     raw
%        Name of file containing raw ship velocity data.
%
%     sm Name of file containing smoothed ship velocity data.
%
%     plot_title
%        Title for plot.
%
%     year_base optional
%        Base year to print together with the x-label.
%
%     n optional
%        Number of points per plots (default is 50).
%
%  5.  RETURNS
%
%  This function does not return anything.
%
%  6.  AUTHORS
%
%  Pierre Jaccard, Geophysical Institute, University of Bergen, 1999
%

% 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


