/[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.3 - (hide annotations) (download)
Fri Jun 24 17:02:03 2011 UTC (14 years, 1 month ago) by gforget
Branch: MAIN
Changes since 1.2: +3 -4 lines 
- add headers to gecmfaces_calc routines
- move line_greatC_TUV_mask.m, line_zonal_TUV_MASKS.m, line_zonal_TUV_mask.m to gcmfaces_legacy
- add replacements gcmfaces_calc/gcmfaces_lines_transp.m and gcmfaces_calc/gcmfaces_lines_zonal.m
- now the zonal and transport 'lines' masks are added to mygrid as LATS_MASKS and LINES_MASKS
- accordingly, the old LATS_MASKS is removed from arguments list in calc_overturn.m etc.
  (calc_transports and basic_diags_display_transport still have LINES_MASKS args. for now)
  and basic_diags_compute_v3_or_v4.m etc. now call gcmfaces_lines_zonal and gcmfaces_lines_transp.
- line_greatC_TUV_MASKS_v3.m and line_greatC_TUV_MASKS_v4.m now are in sample_analysis
  and they simply output the 'lines' definitions (rather than the masks off gcmfaces_lines_transp)
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 gforget 1.2 gcmfaces_path; grid_load(dirGrid,5,'compact');
15 gforget 1.3 gcmfaces_lines_zonal;
16 gforget 1.1 end;
17    
18     %%%%%%%%%%%%%%%
19     %define pathes:
20     %%%%%%%%%%%%%%%
21    
22     dirData='/net/altix3700/raid4/gforget/mysetups/ecco_v4/INPUTS_90x50/';
23     dirIn='./';
24     dirOut=[dirIn 'matlabDiagsCost/']; if isempty(dir(dirOut)); eval(['mkdir ' dirOut ';']); end;
25     runName=pwd; tmp1=strfind(runName,'/'); runName=runName(tmp1(end)+1:end);
26    
27     %%%%%%%%%%%%%%%%%
28     %do computations:
29     %%%%%%%%%%%%%%%%%
30    
31     if doComp==0;
32     eval(['load ' dirOut 'cost_altimeter.mat myflds;']);
33     else;
34    
35     tic;
36     %mdt cost function term (misfit plot)
37     dif_mdt=rdmds2gcmfaces('barfiles/mdtdiff_smooth');
38     sig_mdt=v4_read_bin([dirData 'sigma_MDT_glob_eccollc.bin'],1,0);
39     %store:
40     myflds.dif_mdt=dif_mdt;
41     myflds.sig_mdt=sig_mdt;
42    
43     toc; tic;
44     %lsc cost function term:
45     sladiff_smooth=rdmds2gcmfaces('barfiles/sladiff_smooth');
46     %skip blanks:
47     sladiff_smooth=sladiff_smooth(:,:,1:7:6147); sladiff_smooth(sladiff_smooth==0)=NaN;
48     %compute rms:
49     rms_sladiff_smooth=sqrt(nanmean(sladiff_smooth.^2,3));
50     %get weight:
51     sig_sladiff_smooth=v4_read_bin([dirData 'sigma_SLA_smooth_eccollc.bin'],1,0);
52     %store:
53     myflds.rms_sladiff_smooth=rms_sladiff_smooth;
54     myflds.sig_sladiff_smooth=sig_sladiff_smooth;
55    
56     toc; tic;
57     %pointwise/35days cost function term:
58     sladiff_35d=rdmds2gcmfaces('barfiles/sladiff_raw');
59     %skip blanks:
60     sladiff_35d=sladiff_35d(:,:,1:7:6147); sladiff_35d(sladiff_35d==0)=NaN;
61     %compute rms:
62     rms_sladiff_35d=sqrt(nanmean(sladiff_35d.^2,3));
63     %store:
64     myflds.rms_sladiff_35d=rms_sladiff_35d;
65    
66     toc; tic;
67     %pointwise/1day terms:
68     sum_all=0*mygrid.XC; msk_all=sum_all;
69     for ii=1:3;
70     if ii==1; myset='tp'; elseif ii==2; myset='gfo'; else; myset='ers'; end;
71     %topex pointwise misfits:
72     sladiff_point=rdmds2gcmfaces(['barfiles/sladiff_' myset '_raw']);
73     %compute rms:
74     msk_tmp=1*(sladiff_point~=0);
75     msk_tmp=sum(msk_tmp,3); sum_tmp=sum(sladiff_point.^2,3);
76     sum_all=sum_all+sum_tmp; msk_all=msk_all+msk_tmp;
77     eval(['rms_' myset '=sum_tmp./msk_tmp;']);
78     end;
79     %compute overall rms:
80     rms_sladiff_point=sqrt(sum_all./msk_all);
81     %fill blanks:
82     rms_sladiff_point=diffsmooth2D_extrap_inv(rms_sladiff_point,mygrid.mskC(:,:,1));
83     %get weight:
84     sig_sladiff_point=v4_read_bin([dirData 'sigma_SLA_PWS07r2_glob_eccollc.bin'],1,0);
85     %store:
86     myflds.rms_tp=rms_tp; myflds.rms_gfo=rms_gfo; myflds.rms_ers=rms_ers;
87     myflds.rms_sladiff_point=rms_sladiff_point;
88     myflds.sig_sladiff_point=sig_sladiff_point;
89    
90     %compute zonal mean/median:
91     for ii=1:3;
92     switch ii;
93     case 1; tmp1='mdt'; cost_fld=(mygrid.mskC(:,:,1).*myflds.dif_mdt./myflds.sig_mdt).^2;
94     case 2; tmp1='lsc'; cost_fld=(mygrid.mskC(:,:,1).*myflds.rms_sladiff_smooth./myflds.sig_sladiff_smooth).^2;
95     case 3; tmp1='point'; cost_fld=(mygrid.mskC(:,:,1).*myflds.rms_sladiff_point./myflds.sig_sladiff_point).^2;
96     end;
97 gforget 1.3 cost_zmean=calc_zonmean_T(cost_fld); eval(['mycosts_mean.' tmp1 '=cost_zmean;']);
98     cost_zmedian=calc_zonmedian_T(cost_fld); eval(['mycosts_median.' tmp1 '=cost_zmedian;']);
99 gforget 1.1 end;
100    
101     toc; %write to disk:
102     if doSave; eval(['save ' dirOut 'cost_altimeter.mat myflds mycosts_mean mycosts_median;']); end;
103    
104     end;%if doComp
105    
106     if doPlot;
107     figure; m_map_gcmfaces(myflds.dif_mdt,0,[-0.2:0.02:0.2]); drawnow;
108     figure; m_map_gcmfaces(myflds.rms_sladiff_smooth,0,[0:0.01:0.1]); drawnow;
109     figure; m_map_gcmfaces(myflds.rms_sladiff_35d,0,[0:0.01:0.1]); drawnow;
110     figure; m_map_gcmfaces(myflds.rms_sladiff_point,0,[0:0.02:0.2]); drawnow;
111     end;
112    
113     if doPlot&doPrint;
114     dirFig='../figs/altimeter/'; ff0=gcf-4;
115     for ff=1:4;
116     figure(ff+ff0); saveas(gcf,[dirFig runName '_fig' num2str(ff)],'fig');
117     eval(['print -depsc ' dirFig runName '_fig' num2str(ff) '.eps;']);
118     eval(['print -djpeg90 ' dirFig runName '_fig' num2str(ff) '.jpg;']);
119     end;
120     end;
121    
  ViewVC Help
Powered by ViewVC 1.1.22