/[MITgcm]/MITgcm_contrib/gael/matlab_class/ecco_v4/cost_altimeter.m
ViewVC logotype

Annotation of /MITgcm_contrib/gael/matlab_class/ecco_v4/cost_altimeter.m

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph


Revision 1.1 - (hide annotations) (download)
Mon Feb 7 01:07:34 2011 UTC (14 years, 5 months ago) by gforget
Branch: MAIN
- compute and displays sea level cost functions.

1 gforget 1.1
2     doComp=0;
3     doSave=0;
4     doPlot=1;
5     doPrint=1;
6    
7     %%%%%%%%%%%
8     %load grid:
9     %%%%%%%%%%%
10    
11     dirGrid='/net/weddell/raid3/gforget/gcmfaces/sample_input/GRIDv4/';
12     global mygrid;
13     if isempty(mygrid);
14     gcmfaces_path; grid_load(dirGrid,5);
15     LATS=[-89:89]'; LATS_MASKS=line_zonal_TUV_MASKS(LATS);
16     % SECTIONS_MASKS=line_greatC_TUV_MASKS_v4;
17     end;
18    
19     %%%%%%%%%%%%%%%
20     %define pathes:
21     %%%%%%%%%%%%%%%
22    
23     dirData='/net/altix3700/raid4/gforget/mysetups/ecco_v4/INPUTS_90x50/';
24     dirIn='./';
25     dirOut=[dirIn 'matlabDiagsCost/']; if isempty(dir(dirOut)); eval(['mkdir ' dirOut ';']); end;
26     runName=pwd; tmp1=strfind(runName,'/'); runName=runName(tmp1(end)+1:end);
27    
28     %%%%%%%%%%%%%%%%%
29     %do computations:
30     %%%%%%%%%%%%%%%%%
31    
32     if doComp==0;
33     eval(['load ' dirOut 'cost_altimeter.mat myflds;']);
34     else;
35    
36     tic;
37     %mdt cost function term (misfit plot)
38     dif_mdt=rdmds2gcmfaces('barfiles/mdtdiff_smooth');
39     sig_mdt=v4_read_bin([dirData 'sigma_MDT_glob_eccollc.bin'],1,0);
40     %store:
41     myflds.dif_mdt=dif_mdt;
42     myflds.sig_mdt=sig_mdt;
43    
44     toc; tic;
45     %lsc cost function term:
46     sladiff_smooth=rdmds2gcmfaces('barfiles/sladiff_smooth');
47     %skip blanks:
48     sladiff_smooth=sladiff_smooth(:,:,1:7:6147); sladiff_smooth(sladiff_smooth==0)=NaN;
49     %compute rms:
50     rms_sladiff_smooth=sqrt(nanmean(sladiff_smooth.^2,3));
51     %get weight:
52     sig_sladiff_smooth=v4_read_bin([dirData 'sigma_SLA_smooth_eccollc.bin'],1,0);
53     %store:
54     myflds.rms_sladiff_smooth=rms_sladiff_smooth;
55     myflds.sig_sladiff_smooth=sig_sladiff_smooth;
56    
57     toc; tic;
58     %pointwise/35days cost function term:
59     sladiff_35d=rdmds2gcmfaces('barfiles/sladiff_raw');
60     %skip blanks:
61     sladiff_35d=sladiff_35d(:,:,1:7:6147); sladiff_35d(sladiff_35d==0)=NaN;
62     %compute rms:
63     rms_sladiff_35d=sqrt(nanmean(sladiff_35d.^2,3));
64     %store:
65     myflds.rms_sladiff_35d=rms_sladiff_35d;
66    
67     toc; tic;
68     %pointwise/1day terms:
69     sum_all=0*mygrid.XC; msk_all=sum_all;
70     for ii=1:3;
71     if ii==1; myset='tp'; elseif ii==2; myset='gfo'; else; myset='ers'; end;
72     %topex pointwise misfits:
73     sladiff_point=rdmds2gcmfaces(['barfiles/sladiff_' myset '_raw']);
74     %compute rms:
75     msk_tmp=1*(sladiff_point~=0);
76     msk_tmp=sum(msk_tmp,3); sum_tmp=sum(sladiff_point.^2,3);
77     sum_all=sum_all+sum_tmp; msk_all=msk_all+msk_tmp;
78     eval(['rms_' myset '=sum_tmp./msk_tmp;']);
79     end;
80     %compute overall rms:
81     rms_sladiff_point=sqrt(sum_all./msk_all);
82     %fill blanks:
83     rms_sladiff_point=diffsmooth2D_extrap_inv(rms_sladiff_point,mygrid.mskC(:,:,1));
84     %get weight:
85     sig_sladiff_point=v4_read_bin([dirData 'sigma_SLA_PWS07r2_glob_eccollc.bin'],1,0);
86     %store:
87     myflds.rms_tp=rms_tp; myflds.rms_gfo=rms_gfo; myflds.rms_ers=rms_ers;
88     myflds.rms_sladiff_point=rms_sladiff_point;
89     myflds.sig_sladiff_point=sig_sladiff_point;
90    
91     %compute zonal mean/median:
92     for ii=1:3;
93     switch ii;
94     case 1; tmp1='mdt'; cost_fld=(mygrid.mskC(:,:,1).*myflds.dif_mdt./myflds.sig_mdt).^2;
95     case 2; tmp1='lsc'; cost_fld=(mygrid.mskC(:,:,1).*myflds.rms_sladiff_smooth./myflds.sig_sladiff_smooth).^2;
96     case 3; tmp1='point'; cost_fld=(mygrid.mskC(:,:,1).*myflds.rms_sladiff_point./myflds.sig_sladiff_point).^2;
97     end;
98     cost_zmean=calc_zonmean_T(cost_fld,LATS_MASKS); eval(['mycosts_mean.' tmp1 '=cost_zmean;']);
99     cost_zmedian=calc_zonmedian_T(cost_fld,LATS_MASKS); eval(['mycosts_median.' tmp1 '=cost_zmedian;']);
100     end;
101    
102     toc; %write to disk:
103     if doSave; eval(['save ' dirOut 'cost_altimeter.mat myflds mycosts_mean mycosts_median;']); end;
104    
105     end;%if doComp
106    
107     if doPlot;
108     figure; m_map_gcmfaces(myflds.dif_mdt,0,[-0.2:0.02:0.2]); drawnow;
109     figure; m_map_gcmfaces(myflds.rms_sladiff_smooth,0,[0:0.01:0.1]); drawnow;
110     figure; m_map_gcmfaces(myflds.rms_sladiff_35d,0,[0:0.01:0.1]); drawnow;
111     figure; m_map_gcmfaces(myflds.rms_sladiff_point,0,[0:0.02:0.2]); drawnow;
112     end;
113    
114     if doPlot&doPrint;
115     dirFig='../figs/altimeter/'; ff0=gcf-4;
116     for ff=1:4;
117     figure(ff+ff0); saveas(gcf,[dirFig runName '_fig' num2str(ff)],'fig');
118     eval(['print -depsc ' dirFig runName '_fig' num2str(ff) '.eps;']);
119     eval(['print -djpeg90 ' dirFig runName '_fig' num2str(ff) '.jpg;']);
120     end;
121     end;
122    

  ViewVC Help
Powered by ViewVC 1.1.22