Warning, /utils/matlab/cs_grid/use_psiLine.m is written in an unsupported language. File is not indexed.
view on githubraw file Latest commit 4ec37fd8 on 2023-10-03 22:00:32 UTC
4ec37fd829 Jean*0001 % input: krd, kfac, kgr
0002 % kfac = 0 : read-in 1 diagnostics output file containing both UVELMASS & VVELMASS
0003 % kfac > 0 : read-in 2 Time-ave output files ; kfac=1 with hFac in ; kfac=2 : add hFac here.
0004 krd=2; kfac=0; kgr=1;
0005
0006 if krd > 0,
0007 %-- load velocity :
0008 gDir='grid_files/'; rDir='output_files/'; nit=72000;
0009 %rDir='../res_s1y/'; nit=72360;
0010 %gDir=rDir;
0011 titexp=strrep(rDir(1:end-1),'_','\_'); tit0=['; it=',int2str(nit)];
0012 if kfac > 0,
0013 %-- time-ave output:
0014 uu=rdmds([rDir,'hUtave'],nit);
0015 vv=rdmds([rDir,'hVtave'],nit);
0016 if kfac==2,
0017 %- without included hFac:
0018 hFacW=rdmds([rDir,'hFacW']);
0019 hFacS=rdmds([rDir,'hFacS']);
0020 uu=hFacW.*uu; vv=hFacS.*vv;
0021 end
0022 else
0023 %-- diagnostics output:
0024 namFil='dynDiag';
0025 [v4d,dum,M]=rdmds([rDir,namFil],nit);
0026 eval(M);
0027 tmp=strcmp(fldList,'UVELMASS'); ju=find(tmp==1);
0028 tmp=strcmp(fldList,'VVELMASS'); jv=find(tmp==1);
0029 if isempty(ju) | isempty(jv),
0030 error(['At least 1 velocity Comp. missing from file: ',namFil]);
0031 end
0032 uu=squeeze(v4d(:,:,:,ju));
0033 vv=squeeze(v4d(:,:,:,jv));
0034 end
0035 n1h=size(uu,1); n2h=size(uu,2);
0036 if n1h == 6*n2h, nc=n2h;
0037 elseif n1h*6 == n2h, nc=n1h;
0038 else
0039 error([' var size: ',int2str(n1h),' x ',int2str(n2h),' does not fit regular cube !']);
0040 end
0041 nr=size(uu,3); nPg=nc*nc*6; nPp2=nPg+2;
0042 end
0043
0044 if krd > 1,
0045 lDir='grid_files/';
0046 psiLineF=[lDir,'psiLine_N2S_cs32'];
0047
0048 %- load : 'psi_N','psi_C','psi_P','psiUV','psi_T' :
0049 load(psiLineF);
0050 %- set ncdf=1 to load MNC (NetCDF) grid-files ;
0051 % or ncdf=0 to load MDS (binary) grid-files :
0052 ncdf=0;
0053 G=load_grid(gDir,ncdf);
0054 xcs=G.xC; ycs=G.yC; xcg=G.xG; ycg=G.yG;
0055 %- load the grid dx,dy , convert to 10^6 m :
0056 dxg=G.dXg; dxg=dxg*1.e-6;
0057 dyg=G.dYg; dyg=dyg*1.e-6;
0058
0059 fprintf([' load psiLine from file=',psiLineF, ...
0060 ' AND dxg,dyg files \n']);
0061 end
0062
0063 %- compute detph integrated flow :
0064 %delR=[50 70 100 140 190 240 290 340 390 440 490 540 590 640 690]; % = delR in m
0065 delR=G.dRf';
0066 ut=reshape(uu,nPg,nr); vt=reshape(vv,nPg,nr);
0067 dxg=reshape(dxg,nPg,1); dyg=reshape(dyg,nPg,1);
0068 var=dyg*delR; ut=var.*ut; ut=sum(ut,2);
0069 var=dxg*delR; vt=var.*vt; vt=sum(vt,2);
0070
0071 if n2h == nc,
0072 %- put long vectors ut,vt in "compact-compatible" format:
0073 ut=reshape(permute(reshape(ut,[nc 6 nc]),[1 3 2]),[nPg 1]);
0074 vt=reshape(permute(reshape(vt,[nc 6 nc]),[1 3 2]),[nPg 1]);
0075 end
0076
0077 %- compute Barotropic Stream Function :
0078 psiNx=size(psi_C,1);psiNy=size(psi_C,2);
0079 psiH=zeros(nPp2,1);
0080 ufac=rem(psi_T,2) ; vfac=fix(psi_T/2) ;
0081 for j=2:psiNy,
0082 for i=1:psi_N(j),
0083 i1=psi_C(i,j); i0=psi_P(i,j); i2=psiUV(i,j);
0084 psiH(i1)=psiH(i0)+ufac(i,j)*ut(i2)+vfac(i,j)*vt(i2);
0085 end
0086 end
0087
0088 if kgr > 0,
0089 %-
0090 shift=-1; ccB=[0 0]; cbV=1; AxBx=[-180 180 -90 90] ;
0091 ccB=[-140 140];
0092 %AxBx=[-210 210 30 90];
0093 grph_CSz(psiH,xcs,ycs,xcg,ycg,ccB(1),ccB(2),shift,cbV,AxBx)
0094 title(titexp);
0095 end
