matlab_class/gcmfaces_IO/process2interp.m

function [listInterp,listNot]=process2interp(listDo);
% [listInterp,listNot]=PROCESS2INTERP(listDo);
% 1) computes listInterp and listNot (if listDo is not provided)
% 2) interpolates and ouput fields in listDo (= precomputed listInterp)
%
% Usage: [listInterp,listNot]=process2interp;
%        process2interp({listInterp{1:3}});
%        process2interp({listInterp{4:25}});
%        process2interp({listInterp{26:47}});


dirIn=[pwd filesep 'nctiles/'];
dirTim=[pwd '/diags/STATE/']; filTim='state_2d_set1';
dirOut=[pwd filesep 'diags_interp_tmp/'];
filAvailDiag=[pwd filesep 'available_diagnostics.log'];
filReadme=[pwd filesep 'README'];

%% ======== PART 1 =======

if isempty(who('listDo'));
  listDirs=dir(dirIn);
  listInterp={};
  listNot={};
  for ii=1:length(listDirs);
    %get units and long name from available_diagnostics.log
    [avail_diag]=read_avail_diag(filAvailDiag,listDirs(ii).name);
    if ~isempty(avail_diag);
    if strcmp(avail_diag.loc_h,'C');
      ndiags=length(listInterp)+1;
      listInterp={listInterp{:},listDirs(ii).name};
    else;
      listNot={listNot{:},listDirs(ii).name};
    end;
    end;
  end;
  return;
end;

%% ======== PART 2 =======

lon=[-179.75:0.5:179.75]; lat=[-89.75:0.5:89.75];
[lat,lon] = meshgrid(lat,lon);
interp=gcmfaces_interp_coeffs(lon(:),lat(:));

if ~isdir(dirOut); mkdir(dirOut); end;

for ii=1:length(listDo);

tic;

nameDiag=listDo{ii};
myDiag=read_nctiles([dirIn nameDiag '/' nameDiag]);
[listTimes]=diags_list_times({dirTim},{filTim});

is3D=length(size(myDiag{1}))==4;

%loop over months and output result
for tt=1:240;
  filOut=sprintf('%s.%010i',nameDiag,listTimes(tt)); 
  if is3D; fldOut=myDiag(:,:,:,tt);
  else; fldOut=myDiag(:,:,tt);
  end;
  %interpolate one field
  tmp1=convert2vector(fldOut,'new');
  tmp0=1*~isnan(tmp1);
  tmp1(isnan(tmp1))=0;
  siz=[size(lon) size(tmp1,2)];
  tmp0=interp.SPM*tmp0;
  tmp1=interp.SPM*tmp1;
  fldOut=reshape(tmp1./tmp0,siz);
  sizOut=size(fldOut);
  %create subdirectory
  if ~isdir([dirOut nameDiag '/']); mkdir([dirOut nameDiag '/']); end;
  %write binary field (masked)
  write2file([dirOut nameDiag '/' filOut '.data'],fldOut,32,0);
  %create meta file
  write2meta([dirOut nameDiag '/' filOut '.data'],sizOut,32,{nameDiag});
end;

fprintf(['DONE: ' nameDiag ' (in ' num2str(toc) 's)\n']);
end;

%% ======== FUNCTIONS =======

function [avail_diag]=read_avail_diag(filAvailDiag,nameDiag);

gcmfaces_global;

avail_diag=[];

fid=fopen(filAvailDiag,'rt');
while ~feof(fid);
    tline = fgetl(fid);
    tmp1=8-length(nameDiag); tmp1=repmat(' ',[1 tmp1]);
    tname = ['|' sprintf('%s',nameDiag) tmp1 '|'];
    if ~isempty(strfind(tline,tname));
        %e.g. tline='   235 |SIatmQnt|  1 |       |SM      U1|W/m^2           |Net atmospheric heat flux, >0 decreases theta';
        %
        tmp1=strfind(tline,'|'); tmp1=tmp1(end-1:end);
        avail_diag.units=strtrim(tline(tmp1(1)+1:tmp1(2)-1));
        avail_diag.longNameDiag=tline(tmp1(2)+1:end);
        %
        tmp1=strfind(tline,'|'); tmp1=tmp1(4:5);
        pars=tline(tmp1(1)+1:tmp1(2)-1);
        %
        if strcmp(pars(2),'M'); avail_diag.loc_h='C';
        elseif strcmp(pars(2),'U'); avail_diag.loc_h='W';
        elseif strcmp(pars(2),'V'); avail_diag.loc_h='S';
        end;
        %
        avail_diag.loc_z=pars(9);
        %
        if strcmp(pars(10),'1'); avail_diag.nr=1;
        else; avail_diag.nr=length(mygrid.RC);
        end;
    end;
end;
fclose(fid);

1	gforget	1.1	function [listInterp,listNot]=process2interp(listDo);
2			% [listInterp,listNot]=PROCESS2INTERP(listDo);
3			% 1) computes listInterp and listNot (if listDo is not provided)
4			% 2) interpolates and ouput fields in listDo (= precomputed listInterp)
5			%
6			% Usage: [listInterp,listNot]=process2interp;
7			% process2interp({listInterp{1:3}});
8			% process2interp({listInterp{4:25}});
9			% process2interp({listInterp{26:47}});
10
11
12			dirIn=[pwd filesep 'nctiles/'];
13			dirTim=[pwd '/diags/STATE/']; filTim='state_2d_set1';
14			dirOut=[pwd filesep 'diags_interp_tmp/'];
15			filAvailDiag=[pwd filesep 'available_diagnostics.log'];
16			filReadme=[pwd filesep 'README'];
17
18			%% ======== PART 1 =======
19
20			if isempty(who('listDo'));
21			listDirs=dir(dirIn);
22			listInterp={};
23			listNot={};
24			for ii=1:length(listDirs);
25			%get units and long name from available_diagnostics.log
26			[avail_diag]=read_avail_diag(filAvailDiag,listDirs(ii).name);
27			if ~isempty(avail_diag);
28			if strcmp(avail_diag.loc_h,'C');
29			ndiags=length(listInterp)+1;
30			listInterp={listInterp{:},listDirs(ii).name};
31			else;
32			listNot={listNot{:},listDirs(ii).name};
33			end;
34			end;
35			end;
36			return;
37			end;
38
39			%% ======== PART 2 =======
40
41			lon=[-179.75:0.5:179.75]; lat=[-89.75:0.5:89.75];
42			[lat,lon] = meshgrid(lat,lon);
43			interp=gcmfaces_interp_coeffs(lon(:),lat(:));
44
45			if ~isdir(dirOut); mkdir(dirOut); end;
46
47			for ii=1:length(listDo);
48
49			tic;
50
51			nameDiag=listDo{ii};
52			myDiag=read_nctiles([dirIn nameDiag '/' nameDiag]);
53			[listTimes]=diags_list_times({dirTim},{filTim});
54
55			is3D=length(size(myDiag{1}))==4;
56
57			%loop over months and output result
58			for tt=1:240;
59			filOut=sprintf('%s.%010i',nameDiag,listTimes(tt));
60			if is3D; fldOut=myDiag(:,:,:,tt);
61			else; fldOut=myDiag(:,:,tt);
62			end;
63			%interpolate one field
64			tmp1=convert2vector(fldOut,'new');
65			tmp0=1*~isnan(tmp1);
66			tmp1(isnan(tmp1))=0;
67			siz=[size(lon) size(tmp1,2)];
68			tmp0=interp.SPM*tmp0;
69			tmp1=interp.SPM*tmp1;
70			fldOut=reshape(tmp1./tmp0,siz);
71			sizOut=size(fldOut);
72			%create subdirectory
73			if ~isdir([dirOut nameDiag '/']); mkdir([dirOut nameDiag '/']); end;
74			%write binary field (masked)
75			write2file([dirOut nameDiag '/' filOut '.data'],fldOut,32,0);
76			%create meta file
77			write2meta([dirOut nameDiag '/' filOut '.data'],sizOut,32,{nameDiag});
78			end;
79
80			fprintf(['DONE: ' nameDiag ' (in ' num2str(toc) 's)\n']);
81			end;
82
83			%% ======== FUNCTIONS =======
84
85			function [avail_diag]=read_avail_diag(filAvailDiag,nameDiag);
86
87			gcmfaces_global;
88
89			avail_diag=[];
90
91			fid=fopen(filAvailDiag,'rt');
92			while ~feof(fid);
93			tline = fgetl(fid);
94			tmp1=8-length(nameDiag); tmp1=repmat(' ',[1 tmp1]);
95			tname = ['\|' sprintf('%s',nameDiag) tmp1 '\|'];
96			if ~isempty(strfind(tline,tname));
97			%e.g. tline=' 235 \|SIatmQnt\| 1 \| \|SM U1\|W/m^2 \|Net atmospheric heat flux, >0 decreases theta';
98			%
99			tmp1=strfind(tline,'\|'); tmp1=tmp1(end-1:end);
100			avail_diag.units=strtrim(tline(tmp1(1)+1:tmp1(2)-1));
101			avail_diag.longNameDiag=tline(tmp1(2)+1:end);
102			%
103			tmp1=strfind(tline,'\|'); tmp1=tmp1(4:5);
104			pars=tline(tmp1(1)+1:tmp1(2)-1);
105			%
106			if strcmp(pars(2),'M'); avail_diag.loc_h='C';
107			elseif strcmp(pars(2),'U'); avail_diag.loc_h='W';
108			elseif strcmp(pars(2),'V'); avail_diag.loc_h='S';
109			end;
110			%
111			avail_diag.loc_z=pars(9);
112			%
113			if strcmp(pars(10),'1'); avail_diag.nr=1;
114			else; avail_diag.nr=length(mygrid.RC);
115			end;
116			end;
117			end;
118			fclose(fid);
119