% make_ptflx
%
% This script interpolates the provided PE_Forcing flux file to the
% indicated point.

%-------------------------------------------------------------------------------
% Set run parameters.
%
%    FLXIN......Name of PE_forcing netCDF file.
%    GRDIN......Name of corresponding GRIDS netCDF file.
%    PTLON......Longitude of interpolation point.
%    PTLAT......Latitude of interpolation point.
%    WINDFILE...Name of output Wind Stress file.
%    HEATFILE...Name of output Net Heat Flux file.
%    EMPFILE....Name of output Water Flux file.
%-------------------------------------------------------------------------------

flxhome = '../PE_forcing/Plume/';
flxin   = [flxhome,'pe_disp_2003073100_2003091000plume.nc'];
clear flxhome

grdhome = '../Grids/Plume/';
grdin   = [grdhome,'grids_Plume_cnd05.nc'];
clear grdhome

ptlon = -122.0;
ptlat = 36.0;

windfile = '../windstress.dat';
heatfile = '../netheatflux.dat';
empfile  = '../waterflux.dat';

%-------------------------------------------------------------------------------
% Transform interpolation point to grid units.
%-------------------------------------------------------------------------------

%----------------------------------------------
%--- Read coordinate definition parameters. ---
%----------------------------------------------

[coord,nx,ny,lonc,latc,dx,dy,delx,dely,theta] = grdprm (grdin);

if isempty(theta)
   disp (' ');
   disp ('***Error:  MAKE_PTFLX - unable to read GRIDS netCDF file:');
   disp (['          ',setstr(34),grdin,setstr(34)]);
   disp (' ');
   return;
end;

clear grdin

%--------------------------------------
%--- Transform interpolation point. ---
%--------------------------------------

[ptxT,ptyT] = ll2xy (ptlon,ptlat,coord,nx,ny,dx,dy,lonc,latc,delx,dely,theta);

if ( (ptxT<1) | (ptxT>nx) | (ptyT<1) | (ptyT>ny) )
   disp (' ');
   disp ('***Error:  MAKE_PTFLX - interpolation point outside of tracer grid:');
   disp (['           ptxT = ',num2str(ptxT),'  valid = [1, ',num2str(nx),']']);
   disp (['           ptyT = ',num2str(ptyT),'  valid = [1, ',num2str(ny),']']);
   disp (' ');
   return;
end;

ptxV = ptxT - 0.5;
ptyV = ptyT - 0.5;

if ( (ptxV<1) | (ptxV>nx) | (ptyV<1) | (ptyV>ny) )
   disp (' ');
   disp (['***Error:  MAKE_PTFLX - interpolation point outside of ', ...
          'velocity grid:']);
   disp (['           ptxV = ',num2str(ptxV),'  valid = [1, ',num2str(nx),']']);
   disp (['           ptyV = ',num2str(ptyV),'  valid = [1, ',num2str(ny),']']);
   disp (' ');
   return;
end;

clear ptlon ptlat coord nx ny dx dy lonc latc delx dely theta

%-------------------------------------------------------------------------------
% Open PE_forcing netCDF file.
%-------------------------------------------------------------------------------

ncflx = netcdf (flxin);

if isempty(ncflx)
   disp (' ');
   disp ('***Error:  MAKE_PTFLX - unable to open PE_forcing netCDF file:');
   disp (['          ',setstr(34),flxin,setstr(34)]);
   disp (' ');
   return;
end;

%-------------------------------------------------------------------------------
% Interpolate wind stress.
%-------------------------------------------------------------------------------

disp ('   interpolating wind stress:');

vec_intp (ncflx,flxin,'smflux','smf_time',ptxV,ptyV,windfile);

clear ptxV ptyV windfile

%-------------------------------------------------------------------------------
% Interpolate net heat flux.
%-------------------------------------------------------------------------------

disp ('   interpolating net heat flux:');

scl_intp (ncflx,flxin,'shflux','shf_time',ptxT,ptyT,heatfile);

clear heatfile

%-------------------------------------------------------------------------------
% Interpolate water flux.
%-------------------------------------------------------------------------------

disp ('   interpolating water flux:');

scl_intp (ncflx,flxin,'swflux','swf_time',ptxT,ptyT,empfile);

clear ptxT ptyT empfile

%-------------------------------------------------------------------------------
% Close PE_forcing netCDF file.
%-------------------------------------------------------------------------------

status = close (ncflx);

if ~isempty(status)
   disp (' ');
   disp ('***Error:  MAKE_PTFLX - unable to close PE_forcing netCDF file:');
   disp (['          ',setstr(34),flxin,setstr(34)]);
   disp (' ');
   return;
end;

clear flxin ncflx status
