function [data, flags] = flagdata(data, lim, varargin)
%  1.  NAME
%
%  flagdata - Flag data
%
%  2.  SYNOPSIS
%
%  [data, flags] = flagdata(data, lim, options )
%
%  3.  DESCRIPTION
%
%  This function parses the input data and flags them as bad according to
%  specified limits. The synthax of this vector is defined in stat(5m).
%  This function returns a list of flags for each data being tested.
%  These have the following meaning:
%
%     -Inf
%        Corresponding data was less than minimal value,and limit mode
%        was inside.
%
%     +Inf
%        Corresponding data was more than maximal value,and limit mode
%        was inside.
%
%     1  Corresponding data passed the limit test.
%
%     0  Corresponding data was inside intervall,but limit mode was set
%        to outside.
%
%     NaN
%        Corresponding data was already set to NaN on input.
%
%  4.  PARAMETERS
%
%     data
%        Array (MxN) of M variables,with N different values that should
%        be flagged.
%
%     lim
%        Array (Mx3) specifying the limits for the M variables.
%
%     options keyword , value
%        The rest of the arguments is a list of options. See below for
%        implemented options.
%
%  5.  OPTIONS
%
%  Implemented options are:
%
%     -hdg vector
%        Non-NaN values in the specified vector correspond to column
%        numbers in data containing heading values. These are converted
%        within 0 and 360 before applying the limit tests.
%
%  6.  RETURNS
%
%     data
%        Array (MxN),with all data not passing the limit test set to NaN.
%
%     flags
%        Array (MxN) with the flags valeus resulting from the limit tests
%        " as described above.
%
%  7.  SEE ALSO
%
%  flagdata(3m)stat(5m)statcalc(3m)statinit(3m)statlim(3m)stat_str(3m)statupd(3m)
%
%  8.  AUTHORS
%
%  Pierre Jaccard, Geophysical Institute, University of Bergen, 1999
%

%--- Defaults:
nvars = size(data,1);
if size(lim,1) ~= nvars,
  error([ 'Limits and data do not have same number of variables' ]);
end
hdg = ones(1,nvars) * NaN;
d = data;

%--- Options:
i = 1;
while i < length(varargin),
  if isstr(varargin{i}),
	 switch(varargin{i})
		case '-hdg',
		  hdg = varargin{i+1};
		  i = i+2;
		otherwise,
		  error([ 'Option ' varargin{i} ' not implemented' ]);
	 end
  end
end

%--- Convert for heading:
ihdg = find(~isnan(hdg));
if ~isempty(ihdg),
  d(ihdg,:) = mod(d(ihdg,:), 360);
end

%--- Indexes to limit array:
min_ = 1; max_ = 2; mode_ = 3;

%--- Flag data:
ii = 1:nvars;
flags = ones(size(d));
for i=ii,
  %--- Get NaNs:
  bad = find(isnan(d(i,:)));
  if ~isempty(bad),
	 flags(i,bad) = NaN;
  end
  %--- Inside checks:
  if lim(i,mode_) == 1,
	 bad = find(d(i,:)>lim(i,max_));
	 if ~isempty(bad),
		flags(i,bad) = Inf;
	 end
	 bad = find(d(i,:)<lim(i,min_));
	 if ~isempty(bad),
		flags(i,bad) = -Inf;
	 end
  end
  %--- Outside checks:
  if lim(i,mode_) == -1,
	 bad = find(d(i,:)<lim(i,max_));
	 if ~isempty(bad),
		flags(i,bad) = 0;
	 end
	 bad = find(d(i,:)>lim(i,min_));
	 if ~isempty(bad),
		flags(i,bad) = 0;
	 end
  end
end

%--- Flag data:
for i=ii,
  bad = find(flags(i,:) ~= 1);
  if ~isempty(bad),
	 data(i,bad) = NaN;
  end
end

return;


  
