| 1 |
function []=example_MITprof(); |
| 2 |
% EXAMPLE_MITPROF illustrates the use of MITprof datasets by computing |
| 3 |
% Argo stats over a regional domain that includes the ACC core |
| 4 |
|
| 5 |
gcmfaces_global; |
| 6 |
|
| 7 |
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
| 8 |
if myenv.verbose>0; |
| 9 |
gcmfaces_msg('==============================================='); |
| 10 |
gcmfaces_msg(['*** entering example_MITprof: illustrate ' ... |
| 11 |
'the use of MITprof datasets by computing Argo ACC stats'],''); |
| 12 |
end; |
| 13 |
|
| 14 |
%%%%%%%%%%%%%%%%% |
| 15 |
%load grid and setup test case: |
| 16 |
%%%%%%%%%%%%%%%%% |
| 17 |
|
| 18 |
if isempty(mygrid); |
| 19 |
grid_load; |
| 20 |
end; |
| 21 |
|
| 22 |
myenv.nctilesdir=fullfile('release2_climatology',filesep,'nctiles_climatology',filesep); |
| 23 |
myenv.mitprofdir=fullfile('release2_climatology',filesep,'profiles',filesep); |
| 24 |
|
| 25 |
if ~isdir(myenv.nctilesdir); |
| 26 |
diags=[]; |
| 27 |
help gcmfaces_demo; |
| 28 |
warning(['skipping example_MITprof (missing ' myenv.nctilesdir ')']); |
| 29 |
return; |
| 30 |
end; |
| 31 |
|
| 32 |
if ~isdir(myenv.mitprofdir); |
| 33 |
diags=[]; |
| 34 |
help gcmfaces_demo; |
| 35 |
warning(['skipping example_MITprof (missing ' myenv.mitprofdir ')']); |
| 36 |
return; |
| 37 |
end; |
| 38 |
|
| 39 |
%%%%%%%%%%%%%%%%%%%%%%% |
| 40 |
%define regional mask |
| 41 |
|
| 42 |
if myenv.verbose>0; |
| 43 |
gcmfaces_msg('* call calc_barostream and define regional mask accordingly'); |
| 44 |
end; |
| 45 |
|
| 46 |
%compute mean streamfunction to define mask: |
| 47 |
fldU=mygrid.mskW.*mean(read_nctiles([myenv.nctilesdir 'UVELMASS/UVELMASS']),4); |
| 48 |
fldV=mygrid.mskS.*mean(read_nctiles([myenv.nctilesdir 'VVELMASS/VVELMASS']),4); |
| 49 |
[fldBAR]=calc_barostream(fldU,fldV); |
| 50 |
|
| 51 |
%define 2D mask: |
| 52 |
msk2D=fldBAR; |
| 53 |
msk2D(fldBAR>100)=0; msk2D(fldBAR<30)=0; msk2D(mygrid.YC>0)=0; |
| 54 |
msk2D(msk2D>0)=1; msk2D=msk2D.*mygrid.mskC(:,:,1); |
| 55 |
|
| 56 |
%max climatological MLD could aternatively be used to define mask: |
| 57 |
% maxMLD=mygrid.mskC(:,:,1).*max(read_nctiles('release2_climatology/nctiles_climatology/MXLDEPTH/MXLDEPTH'),[],3); |
| 58 |
|
| 59 |
%define 3D mask: |
| 60 |
msk3D=repmat(msk2D,[1 1 50]).*mygrid.mskC; |
| 61 |
msk3D(:,:,26:50)=0.*msk3D(:,:,26:50); |
| 62 |
|
| 63 |
|
| 64 |
%%%%%%%%%%%%%%%%%%%%%%% |
| 65 |
%read in MITprof data sets and matching climatology profiles |
| 66 |
|
| 67 |
if myenv.verbose>0; |
| 68 |
gcmfaces_msg('* read-in and select Argo profiles'); |
| 69 |
end; |
| 70 |
|
| 71 |
listData={'MITprof_mar2016_argo9506*','MITprof_mar2016_argo0708*','MITprof_mar2016_argo0910*',... |
| 72 |
'MITprof_mar2016_argo1112*','MITprof_mar2016_argo1314*','MITprof_mar2016_argo1515*'}; |
| 73 |
listVar={'prof_T'}; listV={'T'}; |
| 74 |
|
| 75 |
for vv=1:length(listVar); |
| 76 |
tmpprof=MITprof_stats_load(myenv.mitprofdir,listData,listV{vv},listVar{vv}); |
| 77 |
% |
| 78 |
loc_tile=gcmfaces_loc_tile(90,90,tmpprof.prof_lon,tmpprof.prof_lat); |
| 79 |
tmpprof.prof_msk=convert2array(msk2D); |
| 80 |
tmpprof.prof_msk=tmpprof.prof_msk(loc_tile.point); |
| 81 |
tmpprof=MITprof_subset(tmpprof,'msk',1); |
| 82 |
% |
| 83 |
if vv==1; myprof=tmpprof; end; |
| 84 |
myprof.(listVar{vv})=tmpprof.prof; |
| 85 |
% |
| 86 |
clear tmpprof; |
| 87 |
end; |
| 88 |
myprof=rmfield(myprof,{'prof','prof_msk'}); |
| 89 |
|
| 90 |
if myenv.verbose>0; |
| 91 |
gcmfaces_msg('* subsample monthly climatology to profile locations'); |
| 92 |
end; |
| 93 |
|
| 94 |
fldIn.fil=fullfile(myenv.nctilesdir,filesep,'THETA'); |
| 95 |
fldIn.tim='monclim'; fldIn.name='prof_Tclim'; |
| 96 |
myprof=MITprof_resample(myprof,fldIn); |
| 97 |
|
| 98 |
|
| 99 |
%%%%%%%%%%%%%%%%%%%%%%% |
| 100 |
%compute monthly averages |
| 101 |
|
| 102 |
if myenv.verbose>0; |
| 103 |
gcmfaces_msg('* call MITprof_wrapper to compute monthly climatology'); |
| 104 |
end; |
| 105 |
|
| 106 |
%specify operation to compute via MITprof_wrapper |
| 107 |
myop=[]; |
| 108 |
%myop.op_name='mean'; |
| 109 |
myop.op_name='cycle'; myop.op_tim=[0:7:365]; |
| 110 |
myop.op_vars={'prof_T','prof_Tclim'}; |
| 111 |
|
| 112 |
%compute cycle via MITprof_wrapper |
| 113 |
[cy_T,cy_Tclim]=MITprof_wrapper(myprof,myop); |
| 114 |
|
| 115 |
if myenv.verbose>0; |
| 116 |
gcmfaces_msg('* call MITprof_wrapper via bootstrp to compute std'); |
| 117 |
end; |
| 118 |
|
| 119 |
%store myprof and myop into global variable |
| 120 |
global myprofmyop; myprofmyop=myprof; |
| 121 |
for ii=fieldnames(myop)'; myprofmyop.(ii{1})=myop.(ii{1}); end; |
| 122 |
|
| 123 |
%store mask since it could be useful for plotting or other purposes |
| 124 |
myprofmyop.msk2D=msk2D; |
| 125 |
myprofmyop.msk3D=msk3D; |
| 126 |
|
| 127 |
%compute boostrap samples using MITprof_wrapper |
| 128 |
K=[1:myprofmyop.np]; N=100; |
| 129 |
myout = bootstrp(N, @MITprof_wrapper,K); |
| 130 |
|
| 131 |
%compute standard deviation |
| 132 |
bo_cy_T=reshape(cell2mat(myout(:,1)'),[size(cy_T) N]); |
| 133 |
bo_cy_Tclim=reshape(cell2mat(myout(:,2)'),[size(cy_T) N]); |
| 134 |
std_cy_T=std(bo_cy_T,[],3); |
| 135 |
std_cy_Tclim=std(bo_cy_Tclim,[],3); |
| 136 |
|
| 137 |
%%%%%%%%%%%%%%%%%%%%%%% |
| 138 |
%display results |
| 139 |
|
| 140 |
figure; |
| 141 |
subplot(2,2,1); imagescnan(cy_T'); colorbar; title('mean using Argo'); |
| 142 |
subplot(2,2,2); imagescnan(cy_Tclim'); colorbar; title('mean using climatology'); |
| 143 |
subplot(2,2,3); imagescnan(std_cy_T'); colorbar; title('std using Argo'); |
| 144 |
subplot(2,2,4); imagescnan(std_cy_Tclim'); colorbar; title('std using climatology'); |
| 145 |
|
| 146 |
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
| 147 |
if myenv.verbose>0; |
| 148 |
gcmfaces_msg('*** leaving example_MITprof'); |
| 149 |
gcmfaces_msg('==============================================='); |
| 150 |
end; |
| 151 |
|
Parent Directory
| 
Revision Log
| 
Revision Graph