MITgcm_contrib/mitgcm_tools/loadforcing.m

function [FORCING] = loadforcing(varargin)
%loadforcing()
%loadforcing(DIRECTORY)
%
%Reads MITgcm parameter input files ("data") file to create a FORCING structure
%If DIRECTORY is not specified the current working directory is used.
%
%e.g.
% >> FORCING=loadforcing
% FORCING = 
% >> FORCING2=loadforcing('/scratch/john/run2/');
%
%Written by adcroft@mit.edu, 2002
%$Header:

if nargin==0
 Dir='./';
elseif nargin==1
 Dir=[varargin{1} '/'];
else
 error('I don''t know what to do with the second argument');
end

% Extract names of forcing files from "data" file
datafile=[Dir 'data'];

fid=fopen(datafile,'r');
if fid==-1
 error(['Could not open file:' datafile ' for reading']);
end

tinitfile=grepparameter(datafile,'hydrogthetafile');
sinitfile=grepparameter(datafile,'hydrogsaltfile');
sstfile=grepparameter(datafile,'thetaclimfile');
sssfile=grepparameter(datafile,'saltclimfile');
tauxfile=grepparameter(datafile,'zonalwindfile');
tauyfile=grepparameter(datafile,'meridwindfile');
qfile=grepparameter(datafile,'surfqfile');
empmrfile=grepparameter(datafile,'empmrfile');

GRID=loadgrid(Dir);
nxy=size(GRID.rac);
nxyr=size(GRID.hfacc);
clear GRID

FORCING.fmt=grepparameter(datafile,'readbinaryprec');
if isempty(FORCING.fmt)
 FORCING.fmt=32;
end
if FORCING.fmt==64
 fmt='real*8';
else
 fmt='real*4';
end
FORCING.fmt=fmt;

FORCING.tinit=readforcing(Dir,tinitfile,fmt,nxy,nxyr);
FORCING.sinit=readforcing(Dir,sinitfile,fmt,nxy,nxyr);
FORCING.sst=readforcing(Dir,sstfile,fmt,nxy,nxyr);
FORCING.sss=readforcing(Dir,sssfile,fmt,nxy,nxyr);
FORCING.taux=readforcing(Dir,tauxfile,fmt,nxy,nxyr);
FORCING.tauy=readforcing(Dir,tauyfile,fmt,nxy,nxyr);
FORCING.q=readforcing(Dir,qfile,fmt,nxy,nxyr);
FORCING.empmr=readforcing(Dir,empmrfile,fmt,nxy,nxyr);

function [forcingdata] = readforcing(Dir,filename,fmt,nxy,nxyr)

forcingdata=[];
if isempty(filename)
 return;
end

fid=fopen([Dir filename],'r','b');
if fid==-1
  error(['Could not open file:' filename ' for reading']);
end
forcingdata=fread(fid,inf,fmt);
N=size(forcingdata,1);
if N>prod(nxyr)
 nshape=[nxyr size(forcingdata,1)/prod(nxyr)];
else
 nshape=[nxy N/prod(nxy)];
end
forcingdata=reshape(forcingdata,nshape);
1	adcroft	1.1	function [FORCING] = loadforcing(varargin)
2			%loadforcing()
3			%loadforcing(DIRECTORY)
4			%
5			%Reads MITgcm parameter input files ("data") file to create a FORCING structure
6			%If DIRECTORY is not specified the current working directory is used.
7			%
8			%e.g.
9			% >> FORCING=loadforcing
10			% FORCING =
11			% >> FORCING2=loadforcing('/scratch/john/run2/');
12			%
13			%Written by adcroft@mit.edu, 2002
14			%$Header:
15
16			if nargin==0
17			Dir='./';
18			elseif nargin==1
19			Dir=[varargin{1} '/'];
20			else
21			error('I don''t know what to do with the second argument');
22			end
23
24			% Extract names of forcing files from "data" file
25			datafile=[Dir 'data'];
26
27			fid=fopen(datafile,'r');
28			if fid==-1
29			error(['Could not open file:' datafile ' for reading']);
30			end
31
32			tinitfile=grepparameter(datafile,'hydrogthetafile');
33			sinitfile=grepparameter(datafile,'hydrogsaltfile');
34			sstfile=grepparameter(datafile,'thetaclimfile');
35			sssfile=grepparameter(datafile,'saltclimfile');
36			tauxfile=grepparameter(datafile,'zonalwindfile');
37			tauyfile=grepparameter(datafile,'meridwindfile');
38			qfile=grepparameter(datafile,'surfqfile');
39			empmrfile=grepparameter(datafile,'empmrfile');
40
41			GRID=loadgrid(Dir);
42			nxy=size(GRID.rac);
43			nxyr=size(GRID.hfacc);
44			clear GRID
45
46			FORCING.fmt=grepparameter(datafile,'readbinaryprec');
47			if isempty(FORCING.fmt)
48			FORCING.fmt=32;
49			end
50			if FORCING.fmt==64
51			fmt='real*8';
52			else
53			fmt='real*4';
54			end
55			FORCING.fmt=fmt;
56
57			FORCING.tinit=readforcing(Dir,tinitfile,fmt,nxy,nxyr);
58			FORCING.sinit=readforcing(Dir,sinitfile,fmt,nxy,nxyr);
59			FORCING.sst=readforcing(Dir,sstfile,fmt,nxy,nxyr);
60			FORCING.sss=readforcing(Dir,sssfile,fmt,nxy,nxyr);
61			FORCING.taux=readforcing(Dir,tauxfile,fmt,nxy,nxyr);
62			FORCING.tauy=readforcing(Dir,tauyfile,fmt,nxy,nxyr);
63			FORCING.q=readforcing(Dir,qfile,fmt,nxy,nxyr);
64			FORCING.empmr=readforcing(Dir,empmrfile,fmt,nxy,nxyr);
65
66			function [forcingdata] = readforcing(Dir,filename,fmt,nxy,nxyr)
67
68			forcingdata=[];
69			if isempty(filename)
70			return;
71			end
72
73			fid=fopen([Dir filename],'r','b');
74			if fid==-1
75			error(['Could not open file:' filename ' for reading']);
76			end
77			forcingdata=fread(fid,inf,fmt);
78			N=size(forcingdata,1);
79			if N>prod(nxyr)
80			nshape=[nxyr size(forcingdata,1)/prod(nxyr)];
81			else
82			nshape=[nxy N/prod(nxy)];
83			end
84			forcingdata=reshape(forcingdata,nshape);