Warning, /verification/tutorial_global_oce_latlon/diags_matlab/mit_barostream.m is written in an unsupported language. File is not indexed.
view on githubraw file Latest commit add29e06 on 2018-01-31 20:35:05 UTC
051ee7f715 Jean*0001 function [psi, psimask] = mit_barostream(u,varargin)
0002 %function [psi, psimask] = mit_barostream(u,umask,dy,dz);
0003 % or
0004 %function [psi, psimask] = mit_barostream(u,gridinformation);
0005 % global barotropic stream function for mitgcm model, time slabs are
0006 % handled as cell objects, integration from north to the south (by
0007 % convention).
0008
0009 if nargin == 2
0010 g = varargin{1};
0011 umask = g.umask;
0012 dy = g.dyg;
0013 dz = g.dz;
0014 elseif nargin == 4
0015 umask = varargin{1};
0016 dy = varargin{2};
0017 dz = varargin{3};
0018 else
0019 error('need 2 (one of which is the grid structure) or 4 arguments')
0020 end
0021
0022 [nx ny nz] = size(umask);
0023 for kz=1:nz
0024 dydzs(:,:,kz) = (umask(:,:,kz).*dy)*dz(kz);
0025 end
0026
0027 % mask for stream function
0028 pmask = change(squeeze(umask(:,:,1)),'==',NaN,0);
0029 % add psi-point to the north of all wet points
0030 pmask(1:nx,2:ny) = pmask(1:nx,2:ny)+pmask(1:nx,1:ny-1);
0031 pmask = change(pmask,'==',0,NaN);
0032 pmask = change(pmask,'~=',NaN,1);
0033 % integrate from the north to the south (by convention), change
0034 % integration direction by flipping the array ubar (transposed because
0035 % of MITgcm conventions)
0036 if iscell(u)
0037 nt = length(u);
0038 psi = cell(size(u));
0039 for k=1:nt
0040 udxdz = change(u{k}.*dydzs,'==',NaN,0);
0041 ubar = squeeze(sum(udxdz,3));
0042 psi{k} = fliplr(cumsum(fliplr(ubar),2)).*pmask;
0043 end
0044 else
0045 nt = size(u,4);
0046 psi = repmat(NaN,[nx ny nt]);
0047 udxdz = change(u.*repmat(dydzs,[1 1 1 nt]),'==',NaN,0);
0048 ubar = squeeze(sum(udxdz,3));
0049 for kt = 1:nt
0050 psi(:,:,kt) = fliplr(cumsum(fliplr(squeeze(ubar(:,:,kt))),2)).*pmask;
0051 end
0052 end
0053 if nargout == 2
0054 psimask = pmask;
0055 end
0056
0057 return
