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	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);