1 |
gforget |
1.2 |
function []=cost_seaicearea(dirModel,dirMat,doComp,dirTex,nameTex); |
2 |
heimbach |
1.1 |
%object: compute cost function term for sea ice data |
3 |
|
|
%inputs: dimodel is the model directory |
4 |
|
|
% dirMat is the directory where diagnozed .mat files will be saved |
5 |
|
|
% -> set it to '' to use the default [dirModel 'mat/'] |
6 |
|
|
% doComp is a switch (1->compute; 0->display) |
7 |
|
|
%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.2 |
% nameTex is the tex file name (default : 'myPlots') |
10 |
heimbach |
1.1 |
|
11 |
|
|
if isempty(dirMat); dirMat=[dirModel 'mat/']; else; dirMat=[dirMat '/']; end; |
12 |
|
|
if isempty(dir(dirMat)); eval(['!mkdir ' dirMat ';']); end; |
13 |
|
|
|
14 |
|
|
%determine if and where to create tex and figures files |
15 |
|
|
dirMat=[dirMat '/']; |
16 |
|
|
if isempty(who('dirTex')); |
17 |
|
|
addToTex=0; |
18 |
|
|
else; |
19 |
|
|
if ~ischar(dirTex); error('mis-specified dirTex'); end; |
20 |
gforget |
1.2 |
addToTex=1; |
21 |
|
|
if isempty(who('nameTex')); nameTex='myPlots'; end; |
22 |
|
|
fileTex=[dirTex nameTex '.tex']; |
23 |
heimbach |
1.1 |
end; |
24 |
|
|
|
25 |
|
|
if doComp; |
26 |
|
|
|
27 |
gforget |
1.5 |
%grid, params and inputs |
28 |
|
|
|
29 |
|
|
gcmfaces_global; global myparms; |
30 |
heimbach |
1.1 |
if ~isfield(mygrid,'XC'); grid_load('./GRID/',5,'compact'); end; |
31 |
|
|
if ~isfield(mygrid,'LATS_MASKS'); gcmfaces_lines_zonal; end; |
32 |
gforget |
1.5 |
if isfield(myparms,'yearsta'); yearsta=myparms.yearsta; yearend=myparms.yearend; |
33 |
|
|
else; yearsta=1992; yearend=2011; |
34 |
|
|
end; |
35 |
heimbach |
1.1 |
|
36 |
|
|
fld_err=ones(90,1170); |
37 |
|
|
fld_err=convert2gcmfaces(fld_err); |
38 |
|
|
fld_w=fld_err.^-2; |
39 |
|
|
|
40 |
|
|
%%%dirData='/net/weddell/raid3/gforget/ecco_v4/input_files/'; |
41 |
heimbach |
1.3 |
dirData='/net/nares/raid11/ecco-shared/ecco-version-4/input/'; |
42 |
heimbach |
1.1 |
|
43 |
|
|
fileModel=dir([dirModel 'barfiles/gbar_area*data']); |
44 |
|
|
fileModel=['barfiles/' fileModel.name]; |
45 |
|
|
|
46 |
heimbach |
1.3 |
nyears=yearend-yearsta+1; |
47 |
|
|
nmonths=12*nyears; |
48 |
|
|
|
49 |
|
|
fld_dif=convert2gcmfaces(NaN*ones(90,90*13,nmonths)); |
50 |
|
|
|
51 |
gforget |
1.4 |
IceAreaNorthMod=NaN*zeros(1,nmonths); |
52 |
|
|
IceAreaNorthObs=NaN*zeros(1,nmonths); |
53 |
|
|
|
54 |
|
|
IceAreaSouthMod=NaN*zeros(6,nmonths); |
55 |
|
|
IceAreaSouthObs=NaN*zeros(6,nmonths); |
56 |
heimbach |
1.3 |
|
57 |
gforget |
1.5 |
[lonPairs,latPairs,names] = line_greatC_TUV_MASKS_core2_antarctic; |
58 |
|
|
lonLims=[lonPairs(1:5,1);lonPairs(1,1)]; |
59 |
heimbach |
1.1 |
|
60 |
heimbach |
1.3 |
for ycur=yearsta:yearend; |
61 |
heimbach |
1.1 |
tic; |
62 |
|
|
for mcur=1:12; |
63 |
heimbach |
1.3 |
mm=(ycur-yearsta)*12+mcur; |
64 |
|
|
|
65 |
heimbach |
1.1 |
%%% fld_dat=v4_read_bin([dirData 'gael_quick_interp/g_SST_monthly_r2_' num2str(ycur)],mcur,0); |
66 |
|
|
%%% fld_dat=v4_read_bin([dirData 'input_nsidc_postproc/NSIDC0079_v4_gael_' num2str(ycur)],mcur,0); |
67 |
|
|
fld_dat=v4_read_bin([dirData 'input_nsidc_all/nsidc79_monthly_' num2str(ycur)],mcur,0); |
68 |
|
|
|
69 |
gforget |
1.6 |
fld_dat=fld_dat.*mygrid.mskC(:,:,1);%land mask |
70 |
|
|
fld_dat(find(fld_dat<-99))=NaN;%missing data |
71 |
|
|
msk=1+0*fld_dat;%combined mask |
72 |
heimbach |
1.1 |
|
73 |
gforget |
1.6 |
fld_mod=v4_read_bin([dirModel fileModel],mm,0); |
74 |
|
|
fld_mod=fld_mod.*msk;%mask consistent with fld_dat |
75 |
heimbach |
1.1 |
|
76 |
|
|
fld_dif(:,:,mm)=fld_mod-fld_dat; |
77 |
|
|
|
78 |
|
|
fld=fld_mod.*mygrid.RAC.*(mygrid.YC>0); IceAreaNorthMod(mm)=nansum(fld); |
79 |
|
|
fld=fld_dat.*mygrid.RAC.*(mygrid.YC>0); IceAreaNorthObs(mm)=nansum(fld); |
80 |
gforget |
1.5 |
|
81 |
|
|
fld=fld_mod.*mygrid.RAC.*(mygrid.YC<0); IceAreaSouthMod(1,mm)=nansum(fld); |
82 |
|
|
fld=fld_dat.*mygrid.RAC.*(mygrid.YC<0); IceAreaSouthObs(1,mm)=nansum(fld); |
83 |
heimbach |
1.3 |
|
84 |
|
|
for kk=1:5; |
85 |
|
|
tmpmsk=0.*mygrid.XC; |
86 |
gforget |
1.5 |
if lonLims(kk+1) > lonLims(kk) |
87 |
|
|
tmpmsk(find(mygrid.XC >= lonLims(kk) & mygrid.XC < lonLims(kk+1)))=1.; |
88 |
heimbach |
1.3 |
else |
89 |
gforget |
1.5 |
tmpmsk(find(mygrid.XC >= lonLims(kk) & mygrid.XC <= 180.))=1.; |
90 |
|
|
tmpmsk(find(mygrid.XC >= -180. & mygrid.XC < lonLims(kk+1)))=1.; |
91 |
heimbach |
1.3 |
end |
92 |
|
|
tmpmsk=tmpmsk.*(mygrid.YC<0); |
93 |
|
|
% |
94 |
|
|
fld=fld_mod.*mygrid.RAC.*tmpmsk; |
95 |
gforget |
1.5 |
IceAreaSouthMod(kk+1,mm)=nansum(fld); |
96 |
heimbach |
1.3 |
% |
97 |
|
|
fld=fld_dat.*mygrid.RAC.*tmpmsk; |
98 |
gforget |
1.5 |
IceAreaSouthObs(kk+1,mm)=nansum(fld); |
99 |
heimbach |
1.3 |
end |
100 |
heimbach |
1.1 |
|
101 |
|
|
end; |
102 |
|
|
toc; |
103 |
|
|
end; |
104 |
|
|
|
105 |
|
|
fld_rms=sqrt(nanmean(fld_dif.^2,3)); |
106 |
|
|
fld_cost=fld_rms.^2.*fld_w; |
107 |
|
|
|
108 |
gforget |
1.5 |
eval(['save ' dirMat '/cost_seaicearea.mat fld_err fld_rms fld_cost IceArea*;']); |
109 |
heimbach |
1.1 |
|
110 |
|
|
else;%display previously computed results |
111 |
|
|
|
112 |
|
|
global mygrid; |
113 |
|
|
|
114 |
|
|
eval(['load ' dirMat '/cost_seaicearea.mat;']); |
115 |
|
|
|
116 |
|
|
%figure; m_map_gcmfaces(fld_cost,0,{'myCaxis',[0:0.2:1.2 1.5:0.5:3 4:1:6 8 10]}); |
117 |
|
|
%figure; m_map_gcmfaces(fld_err,0,{'myCaxis',[0:0.2:1.2 1.5:0.5:3 4:1:6 8 10]/2}); |
118 |
|
|
|
119 |
|
|
figure; m_map_gcmfaces(fld_rms,0,{'myCaxis',[0:0.1:1.]}); |
120 |
|
|
myCaption={'modeled-observed rms -- sea ice area (K)'}; |
121 |
|
|
|
122 |
|
|
if addToTex; write2tex(fileTex,2,myCaption,gcf); end; |
123 |
|
|
|
124 |
|
|
end; |
125 |
|
|
|
126 |
|
|
|