function [average]=error_avg(OBSname,AVGname);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Copyright (c) 2004 Rutgers University.                                    %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%                                                                           %
% function [average]=error_avg(OBSname,AVGname)                             %
%                                                                           %
% This function over-writes observation error with global inverse time      %
% average squared.                                                          %
%                                                                           %
% On Input:                                                                 %
%                                                                           %
%    OBSname     Observation NetCDF file name.                              %
%    AVGname     Time-average file name.                                    %    
%                                                                           %
% On Output:                                                                %
%                                                                           %
%    average     global squared average.                                    %
%                                                                           %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%  Set state variable type indices.

isFsur=1;  got.Fsur=0;
isUbar=2;  got.Ubar=0;
isVbar=3;  got.Vbar=0;
isUvel=4;  got.Uvel=0;
isVvel=5;  got.Vvel=0;
isTemp=6;  got.Temp=0;
isSalt=7;  got.Salt=0;

%  Inquire NetCDF varaible names.

[vname,nvars]=nc_vname(AVGname);
for n=1:nvars,
  name=deblank(vname(n,:));
  switch (name),
    case 'zeta',
      got.Fsur=1;
    case 'ubar',
      got.Ubar=1;
    case 'vbar',
      got.Vbar=1;
    case 'u',
      got.Uvel=1;
    case 'v',
      got.Vvel=1;
    case 'temp',
      got.Temp=1;
    case 'salt',
      got.Salt=1;
  end,
end,

%  Get observation type and error.

type=nc_read(OBSname,'obs_type');
error=nc_read(OBSname,'obs_error');
variance=nc_read(OBSname,'obs_variance');

%  Compute global (time and space) squared average for each state variable.

if (got.Fsur),
  A2=avg2_obs(AVGname,'zeta');
  A2=A2(:);
  average(isFsur)=mean(A2);
  ind=find(type == isFsur);
  if (~isempty(ind)),
    error(ind)=1/average(isFsur);
  end,
end,

if (got.Ubar),
  A2=avg2_obs(AVGname,'ubar');
  A2=A2(:);
  average(isUbar)=mean(A2);
  ind=find(type == isUbar);
  if (~isempty(ind)),
    error(ind)=1/average(isUbar);
  end,
end,

if (got.Vbar),
  A2=avg2_obs(AVGname,'vbar');
  A2=A2(:);
  average(isVbar)=mean(A2);
  ind=find(type == isVbar);
  if (~isempty(ind)),
    error(ind)=1/average(isVbar);
  end,
end,

if (got.Uvel),
  A2=avg2_obs(AVGname,'u');
  A2=A2(:);
  average(isUvel)=mean(A2);
  ind=find(type == isUvel);
  if (~isempty(ind)),
    error(ind)=1/average(isUvel);
  end,
end,

if (got.Vvel),
  A2=avg2_obs(AVGname,'v');
  A2=A2(:);
  average(isVvel)=mean(A2);
  ind=find(type == isVvel);
  if (~isempty(ind)),
    error(ind)=1/average(isVvel);
  end,
end,

if (got.Temp),
  A2=avg2_obs(AVGname,'temp');
  A2=A2(:);
  average(isTemp)=mean(A2);
  ind=find(type == isTemp);
  if (~isempty(ind)),
    error(ind)=1/average(isTemp);
  end,
end,

if (got.Salt),
  A2=avg2_obs(AVGname,'salt');
  A2=A2(:);
  average(isSalt)=mean(A2);
  ind=find(type == isSalt);
  if (~isempty(ind)),
    error(ind)=1/average(isSalt);
  end,
end,

%  Write variance and new observation error.

status=nc_write(OBSname,'obs_variance',average)

status=nc_write(OBSname,'obs_error',error)

return
