matlab_class/ecco_v4/cost_xx.m

function []=cost_xx(dirModel,dirMat,doComp,dirTex);
%object:        compute cost function term for atmospheric controls
%inputs:        dirModel is the model directory
%               dirMat is the directory where diagnozed .mat files will be saved
%                     -> set it to '' to use the default [dirModel 'mat/']
%               doComp is a switch (1->compute; 0->display)
%optional:      dirTex is the directory where tex and figures files are created
%                 (if not specified then display all results to screen instead)

if isempty(dirMat); dirMat=[dirModel 'mat/']; else; dirMat=[dirMat '/']; end;
if isempty(dir(dirMat)); eval(['mkdir ' dirMat ';']); end;

%determine if and where to create tex and figures files
dirMat=[dirMat '/'];
if isempty(who('dirTex'));
  addToTex=0;
else;
  if ~ischar(dirTex); error('mis-specified dirTex'); end;
  addToTex=1; fileTex=[dirTex 'myPlots.tex'];
end;

for ii=1:6;
switch ii;
case 1; xxName='atemp'; sigName='cap_sigma_tmp2m_degC_eccollc.bin'; cc=2; uni='K';
case 2; xxName='aqh'; sigName='cap_sigma_spfh2m_eccollc.bin'; cc=2; uni='g/kg';
case 3; xxName='tauu'; sigName='cap_sigma_ustr_eccollc.bin'; cc=0.04; uni='N/m2';
case 4; xxName='tauv'; sigName='cap_sigma_vstr_eccollc.bin'; cc=0.04; uni='N/m2';
case 5; xxName='lwdown'; sigName='cap_sigma_dlw_eccollc.bin'; cc=20; uni='W/m2';
case 6; xxName='swdown'; sigName='cap_sigma_dsw_eccollc.bin'; cc=40; uni='W/m2';
end;

if doComp;

%load grid
gcmfaces_global;
if ~isfield(mygrid,'XC'); grid_load('./GRID/',5,'compact'); end;
if ~isfield(mygrid,'LATS_MASKS'); gcmfaces_lines_zonal; end;

dirSig='/net/weddell/raid3/gforget/ecco_v4/input_files/';

%read model cost output
fld_xx=rdmds2gcmfaces([dirModel 'ADXXfiles/xx_' xxName '.effective.*']);
fld_sig=v4_read_bin([dirSig sigName],1,0);

if strcmp(xxName,'aqh'); fld_xx=fld_xx*1000; fld_sig=fld_sig*1000; end;

%compute xx stats
fld_rms=sqrt(mean(fld_xx.^2,3));
fld_mean=mean(fld_xx,3);
fld_std=std(fld_xx,[],3);

%mask
fld_rms=fld_rms.*mygrid.mskC(:,:,1); 
fld_sig=fld_sig.*mygrid.mskC(:,:,1);
fld_mean=fld_mean.*mygrid.mskC(:,:,1);
fld_std=fld_std.*mygrid.mskC(:,:,1);

clear fld_xx;

eval(['save ' dirMat '/cost_xx_' xxName '.mat fld_* cc uni;']);

else;%display previously computed results

global mygrid;

eval(['load ' dirMat '/cost_xx_' xxName '.mat;']);

figure; 
m_map_gcmfaces(fld_sig,0,{'myCaxis',[0:0.05:0.5 0.6:0.1:1 1.25]*cc});
myCaption={['prior uncertainty -- ' xxName ' (' uni ')']};
if addToTex; write2tex(fileTex,2,myCaption,gcf); end;

figure; 
m_map_gcmfaces(fld_rms,0,{'myCaxis',[0:0.05:0.5 0.6:0.1:1 1.25]*cc});
myCaption={['rms adjustment -- ' xxName ' (' uni ')']};
if addToTex; write2tex(fileTex,2,myCaption,gcf); end;

figure;
m_map_gcmfaces(fld_std,0,{'myCaxis',[0:0.05:0.5 0.6:0.1:1 1.25]*cc});
myCaption={['std adjustment -- ' xxName ' (' uni ')']};
if addToTex; write2tex(fileTex,2,myCaption,gcf); end;

figure;
m_map_gcmfaces(fld_mean,0,{'myCaxis',[-0.5:0.05:0.5]*cc});
myCaption={['mean adjustment -- ' xxName ' (' uni ')']};
if addToTex; write2tex(fileTex,2,myCaption,gcf); end;

end;

end;%for ii=1:6;

1	gforget	1.3	function []=cost_xx(dirModel,dirMat,doComp,dirTex);
2	gforget	1.1	%object: compute cost function term for atmospheric controls
3			%inputs: dirModel is the model directory
4	gforget	1.3	% dirMat is the directory where diagnozed .mat files will be saved
5			% -> set it to '' to use the default [dirModel 'mat/']
6	gforget	1.1	% doComp is a switch (1->compute; 0->display)
7	gforget	1.3	%optional: dirTex is the directory where tex and figures files are created
8			% (if not specified then display all results to screen instead)
9	gforget	1.1
10	gforget	1.3	if isempty(dirMat); dirMat=[dirModel 'mat/']; else; dirMat=[dirMat '/']; end;
11			if isempty(dir(dirMat)); eval(['mkdir ' dirMat ';']); end;
12
13			%determine if and where to create tex and figures files
14			dirMat=[dirMat '/'];
15			if isempty(who('dirTex'));
16			addToTex=0;
17			else;
18			if ~ischar(dirTex); error('mis-specified dirTex'); end;
19			addToTex=1; fileTex=[dirTex 'myPlots.tex'];
20			end;
21	gforget	1.1
22			for ii=1:6;
23			switch ii;
24			case 1; xxName='atemp'; sigName='cap_sigma_tmp2m_degC_eccollc.bin'; cc=2; uni='K';
25			case 2; xxName='aqh'; sigName='cap_sigma_spfh2m_eccollc.bin'; cc=2; uni='g/kg';
26			case 3; xxName='tauu'; sigName='cap_sigma_ustr_eccollc.bin'; cc=0.04; uni='N/m2';
27			case 4; xxName='tauv'; sigName='cap_sigma_vstr_eccollc.bin'; cc=0.04; uni='N/m2';
28			case 5; xxName='lwdown'; sigName='cap_sigma_dlw_eccollc.bin'; cc=20; uni='W/m2';
29			case 6; xxName='swdown'; sigName='cap_sigma_dsw_eccollc.bin'; cc=40; uni='W/m2';
30			end;
31
32			if doComp;
33
34			%load grid
35	gforget	1.2	gcmfaces_global;
36			if ~isfield(mygrid,'XC'); grid_load('./GRID/',5,'compact'); end;
37			if ~isfield(mygrid,'LATS_MASKS'); gcmfaces_lines_zonal; end;
38
39	gforget	1.3	dirSig='/net/weddell/raid3/gforget/ecco_v4/input_files/';
40	gforget	1.1
41			%read model cost output
42			fld_xx=rdmds2gcmfaces([dirModel 'ADXXfiles/xx_' xxName '.effective.*']);
43			fld_sig=v4_read_bin([dirSig sigName],1,0);
44
45			if strcmp(xxName,'aqh'); fld_xx=fld_xx1000; fld_sig=fld_sig1000; end;
46
47			%compute xx stats
48			fld_rms=sqrt(mean(fld_xx.^2,3));
49			fld_mean=mean(fld_xx,3);
50			fld_std=std(fld_xx,[],3);
51
52			%mask
53			fld_rms=fld_rms.*mygrid.mskC(:,:,1);
54			fld_sig=fld_sig.*mygrid.mskC(:,:,1);
55			fld_mean=fld_mean.*mygrid.mskC(:,:,1);
56			fld_std=fld_std.*mygrid.mskC(:,:,1);
57
58			clear fld_xx;
59
60	gforget	1.3	eval(['save ' dirMat '/cost_xx_' xxName '.mat fld_* cc uni;']);
61	gforget	1.1
62			else;%display previously computed results
63
64			global mygrid;
65
66	gforget	1.3	eval(['load ' dirMat '/cost_xx_' xxName '.mat;']);
67	gforget	1.1
68			figure;
69			m_map_gcmfaces(fld_sig,0,{'myCaxis',[0:0.05:0.5 0.6:0.1:1 1.25]*cc});
70			myCaption={['prior uncertainty -- ' xxName ' (' uni ')']};
71	gforget	1.3	if addToTex; write2tex(fileTex,2,myCaption,gcf); end;
72	gforget	1.1
73			figure;
74			m_map_gcmfaces(fld_rms,0,{'myCaxis',[0:0.05:0.5 0.6:0.1:1 1.25]*cc});
75			myCaption={['rms adjustment -- ' xxName ' (' uni ')']};
76	gforget	1.3	if addToTex; write2tex(fileTex,2,myCaption,gcf); end;
77	gforget	1.1
78			figure;
79			m_map_gcmfaces(fld_std,0,{'myCaxis',[0:0.05:0.5 0.6:0.1:1 1.25]*cc});
80			myCaption={['std adjustment -- ' xxName ' (' uni ')']};
81	gforget	1.3	if addToTex; write2tex(fileTex,2,myCaption,gcf); end;
82	gforget	1.1
83			figure;
84			m_map_gcmfaces(fld_mean,0,{'myCaxis',[-0.5:0.05:0.5]*cc});
85			myCaption={['mean adjustment -- ' xxName ' (' uni ')']};
86	gforget	1.3	if addToTex; write2tex(fileTex,2,myCaption,gcf); end;
87	gforget	1.1
88			end;
89
90			end;%for ii=1:6;
91