matlab_class/gcmfaces_IO/grid_load_native.m

function []=grid_load_native(dirGrid,nFaces,doWarn);
%object:    load NATIVE FORMAT grid information, convert it to gcmfaces format
%           and encapsulate it in the global mygrid structure.
%inputs: EITHER
%           dirGrid is the directory where the grid files (gcm output) can be found.
%           nFaces is the number of faces in this gcm set-up of current interest.
%           doWarn issue warning if grid is not complete (e.g. hFacC is missing)
%        OR
%           [] if mygrid has already read by grid_load.m from XC.data, YC.data, etc.
%           from which grid_load_native can get dirGrid and nFaces


global mygrid;

alsoDoTheOldWay=0;
if alsoDoTheOldWay;
    global MM Nfaces;
    global dyG dxG dxF dyF dxC dyC dyU dxV rA rAw rAs rAz xC yC xG yG xS yS xW yW;
end;

if isempty(whos('doWarn')); doWarn=1; end;

if nargin>0;
  mygrid=[];
  mygrid.dirGrid=dirGrid;
  mygrid.nFaces=nFaces;
  mygrid.fileFormat='native';
end;

%search for native grid files
files=dir([mygrid.dirGrid 'grid_cs32*bin']);
if isempty(files); files=dir([mygrid.dirGrid 'tile*.mitgrid']); end;
if isempty(files); files=dir([dirGrid '*bin']); end;
while length(files)~=mygrid.nFaces;
  fprintf('dirGrid does not contain the expected number of ''*.bin'' grid files\n');
  nm=input('please specify file name pattern for the native grid such as ''llc*.bin''\n');
  files=dir([dirGrid nm]);
end;

%discard intermediate step grid files:
tmp1=[];
for ii=1:length(files); 
    if isempty(strfind(files(ii).name,'FM')); tmp1=[tmp1;ii]; end; 
end;
files=files(tmp1);

%list fields of interest and their sizes
list_fields2={'XC','YC','DXF','DYF','RAC','XG','YG','DXV','DYU','RAZ',...
    'DXC','DYC','RAW','RAS','DXG','DYG'};
list_fields={'xC','yC','dxF','dyF','rA','xG','yG','dxV','dyU','rAz',...
    'dxC','dyC','rAw','rAs','dxG','dyG'};
list_x={'xC','xC','xC','xC','xC','xG','xG','xG','xG','xG',...
    'xW','xS','xW','xS','xS','xW'};
list_y={'yC','yC','yC','yC','yC','yG','yG','yG','yG','yG',...
    'yW','yS','yW','yS','yS','yW'};
list_ni={'ni','ni','ni','ni','ni','ni+1','ni+1','ni+1','ni+1','ni+1',...
    'ni+1','ni','ni+1','ni','ni','ni+1'};
list_nj={'nj','nj','nj','nj','nj','nj+1','nj+1','nj+1','nj+1','nj+1',...
    'nj','nj+1','nj','nj+1','nj+1','nj'};

if mygrid.nFaces~=length(files); 
    error('wrong specification of nFaces');
else;
    Nfaces=length(files);
end;

for iFile=1:Nfaces;
    tmp1=files(iFile).name;
    %get face dimensions
    if ~isempty(strfind(mygrid.dirGrid,'cs32_tutorial_held_suarez_cs'))|...
            ~isempty(strfind(mygrid.dirGrid,'GRIDcube'));%special case of cs32
        ni=32; nj=32;
    elseif strcmp(tmp1(end-2:end),'bin');%get face dimension from file name
        tmp2=strfind(tmp1,'_');
        ni=str2num(tmp1(tmp2(2)+1:tmp2(3)-1));
        nj=str2num(tmp1(tmp2(3)+1:end-4));
    elseif isfield(mygrid,'XC');%get face dimention from previous grid load
        [ni,nj]=size(mygrid.XC{iFile});
    else;
        error('could not determine face size');
    end;
    if iFile==1&alsoDoTheOldWay; MM=ni; end;
    fid=fopen([mygrid.dirGrid files(iFile).name],'r','b');
    for iFld=1:length(list_fields);
        eval(['nni=' list_ni{iFld} ';']);
        eval(['nnj=' list_nj{iFld} ';']);
        tmp1=fread(fid,[ni+1 nj+1],'float64');
        if alsoDoTheOldWay;
           eval([list_fields{iFld} '{' num2str(iFile) '}.vals=tmp1(1:nni,1:nnj);']);
           eval([list_fields{iFld} '{' num2str(iFile) '}.x=''' list_x{iFld} ''';']);
        end;
        if ~isfield(mygrid,list_fields2{iFld}); eval(['mygrid.' list_fields2{iFld} '=gcmfaces;']); end;
        eval(['mygrid.' list_fields2{iFld} '{iFile}=tmp1(1:ni,1:nj);']);
    end;
    fclose(fid);
    if alsoDoTheOldWay;
       xS{iFile}.vals=(xG{iFile}.vals(2:end,:)+xG{iFile}.vals(1:end-1,:))/2;
       yS{iFile}.vals=(yG{iFile}.vals(2:end,:)+yG{iFile}.vals(1:end-1,:))/2;
       xW{iFile}.vals=(xG{iFile}.vals(:,2:end)+xG{iFile}.vals(:,1:end-1))/2;
       yW{iFile}.vals=(yG{iFile}.vals(:,2:end)+yG{iFile}.vals(:,1:end-1))/2;
   end;
end;

if doWarn;
  list0={'hFacC','hFacS','hFacW','Depth','AngleCS','AngleSN'};
  nWarn=0;
  for ff=1:length(list0); 
    if ~isfield(mygrid,list0{ff}); 
      warning(['The following variable is missing in mygrid: ' list0{ff}]);
      nWarn=nWarn+1;
    end;
  end;
  if nWarn>0; 
    warning('Missing these variables may restrict diagnostic possibilities.');
    nWarn
  end;
end;

1	gforget	1.3	function []=grid_load_native(dirGrid,nFaces,doWarn);
2	gforget	1.2	%object: load NATIVE FORMAT grid information, convert it to gcmfaces format
3			% and encapsulate it in the global mygrid structure.
4	gforget	1.4	%inputs: EITHER
5			% dirGrid is the directory where the grid files (gcm output) can be found.
6	gforget	1.2	% nFaces is the number of faces in this gcm set-up of current interest.
7	gforget	1.4	% doWarn issue warning if grid is not complete (e.g. hFacC is missing)
8			% OR
9			% [] if mygrid has already read by grid_load.m from XC.data, YC.data, etc.
10			% from which grid_load_native can get dirGrid and nFaces
11	gforget	1.1
12	gforget	1.4
13			global mygrid;
14
15			alsoDoTheOldWay=0;
16			if alsoDoTheOldWay;
17	gforget	1.2	global MM Nfaces;
18			global dyG dxG dxF dyF dxC dyC dyU dxV rA rAw rAs rAz xC yC xG yG xS yS xW yW;
19			end;
20	gforget	1.1
21	gforget	1.3	if isempty(whos('doWarn')); doWarn=1; end;
22
23	gforget	1.4	if nargin>0;
24			mygrid=[];
25			mygrid.dirGrid=dirGrid;
26			mygrid.nFaces=nFaces;
27			mygrid.fileFormat='native';
28			end;
29	gforget	1.1
30	gforget	1.3	%search for native grid files
31	gforget	1.4	files=dir([mygrid.dirGrid 'grid_cs32*bin']);
32			if isempty(files); files=dir([mygrid.dirGrid 'tile*.mitgrid']); end;
33			if isempty(files); files=dir([dirGrid '*bin']); end;
34			while length(files)~=mygrid.nFaces;
35	gforget	1.3	fprintf('dirGrid does not contain the expected number of ''*.bin'' grid files\n');
36			nm=input('please specify file name pattern for the native grid such as ''llc*.bin''\n');
37			files=dir([dirGrid nm]);
38			end;
39	gforget	1.2
40			%discard intermediate step grid files:
41			tmp1=[];
42	gforget	1.1	for ii=1:length(files);
43			if isempty(strfind(files(ii).name,'FM')); tmp1=[tmp1;ii]; end;
44			end;
45			files=files(tmp1);
46
47	gforget	1.2	%list fields of interest and their sizes
48	gforget	1.1	list_fields2={'XC','YC','DXF','DYF','RAC','XG','YG','DXV','DYU','RAZ',...
49			'DXC','DYC','RAW','RAS','DXG','DYG'};
50			list_fields={'xC','yC','dxF','dyF','rA','xG','yG','dxV','dyU','rAz',...
51			'dxC','dyC','rAw','rAs','dxG','dyG'};
52			list_x={'xC','xC','xC','xC','xC','xG','xG','xG','xG','xG',...
53			'xW','xS','xW','xS','xS','xW'};
54			list_y={'yC','yC','yC','yC','yC','yG','yG','yG','yG','yG',...
55			'yW','yS','yW','yS','yS','yW'};
56			list_ni={'ni','ni','ni','ni','ni','ni+1','ni+1','ni+1','ni+1','ni+1',...
57			'ni+1','ni','ni+1','ni','ni','ni+1'};
58			list_nj={'nj','nj','nj','nj','nj','nj+1','nj+1','nj+1','nj+1','nj+1',...
59			'nj','nj+1','nj','nj+1','nj+1','nj'};
60
61	gforget	1.4	if mygrid.nFaces~=length(files);
62	gforget	1.2	error('wrong specification of nFaces');
63			else;
64			Nfaces=length(files);
65			end;
66	gforget	1.1
67			for iFile=1:Nfaces;
68			tmp1=files(iFile).name;
69	gforget	1.4	%get face dimensions
70			if ~isempty(strfind(mygrid.dirGrid,'cs32_tutorial_held_suarez_cs'))\|...
71			~isempty(strfind(mygrid.dirGrid,'GRIDcube'));%special case of cs32
72	gforget	1.1	ni=32; nj=32;
73	gforget	1.4	elseif strcmp(tmp1(end-2:end),'bin');%get face dimension from file name
74	gforget	1.1	tmp2=strfind(tmp1,'_');
75			ni=str2num(tmp1(tmp2(2)+1:tmp2(3)-1));
76			nj=str2num(tmp1(tmp2(3)+1:end-4));
77	gforget	1.4	elseif isfield(mygrid,'XC');%get face dimention from previous grid load
78			[ni,nj]=size(mygrid.XC{iFile});
79			else;
80			error('could not determine face size');
81	gforget	1.1	end;
82	gforget	1.4	if iFile==1&alsoDoTheOldWay; MM=ni; end;
83			fid=fopen([mygrid.dirGrid files(iFile).name],'r','b');
84	gforget	1.1	for iFld=1:length(list_fields);
85			eval(['nni=' list_ni{iFld} ';']);
86			eval(['nnj=' list_nj{iFld} ';']);
87			tmp1=fread(fid,[ni+1 nj+1],'float64');
88	gforget	1.4	if alsoDoTheOldWay;
89	gforget	1.1	eval([list_fields{iFld} '{' num2str(iFile) '}.vals=tmp1(1:nni,1:nnj);']);
90			eval([list_fields{iFld} '{' num2str(iFile) '}.x=''' list_x{iFld} ''';']);
91			end;
92	gforget	1.4	if ~isfield(mygrid,list_fields2{iFld}); eval(['mygrid.' list_fields2{iFld} '=gcmfaces;']); end;
93	gforget	1.2	eval(['mygrid.' list_fields2{iFld} '{iFile}=tmp1(1:ni,1:nj);']);
94	gforget	1.1	end;
95			fclose(fid);
96	gforget	1.4	if alsoDoTheOldWay;
97	gforget	1.1	xS{iFile}.vals=(xG{iFile}.vals(2:end,:)+xG{iFile}.vals(1:end-1,:))/2;
98			yS{iFile}.vals=(yG{iFile}.vals(2:end,:)+yG{iFile}.vals(1:end-1,:))/2;
99			xW{iFile}.vals=(xG{iFile}.vals(:,2:end)+xG{iFile}.vals(:,1:end-1))/2;
100			yW{iFile}.vals=(yG{iFile}.vals(:,2:end)+yG{iFile}.vals(:,1:end-1))/2;
101			end;
102			end;
103	gforget	1.2
104	gforget	1.3	if doWarn;
105			list0={'hFacC','hFacS','hFacW','Depth','AngleCS','AngleSN'};
106	gforget	1.4	nWarn=0;
107			for ff=1:length(list0);
108			if ~isfield(mygrid,list0{ff});
109			warning(['The following variable is missing in mygrid: ' list0{ff}]);
110			nWarn=nWarn+1;
111			end;
112			end;
113			if nWarn>0;
114			warning('Missing these variables may restrict diagnostic possibilities.');
115			nWarn
116			end;
117	gforget	1.3	end;
118	gforget	1.2