1 |
dimitri |
1.1 |
krd=2; kgr=1; kwr=0; Tprt=0; |
2 |
|
|
|
3 |
|
|
% $Header: /u/gcmpack/MITgcm/utils/matlab/cs_grid/bk_line/use_psiLine.m,v 1.1 2005/09/15 16:46:28 jmc Exp $ |
4 |
|
|
% $Name: $ |
5 |
|
|
|
6 |
|
|
if krd > 0, |
7 |
|
|
Rac='grid_files/'; |
8 |
|
|
drF=squeeze(rdmds([Rac,'DRF']))'; nr=length(drF); |
9 |
|
|
|
10 |
|
|
%- to make a plot with grph_CS: |
11 |
|
|
fprintf(' load xc,yc,xg,yg file ...'); |
12 |
|
|
xcs=rdmds([Rac,'XC']); |
13 |
|
|
ycs=rdmds([Rac,'YC']); |
14 |
|
|
xcg=rdmds([Rac,'XG']); |
15 |
|
|
ycg=rdmds([Rac,'YG']); |
16 |
|
|
fprintf(' done\n'); |
17 |
|
|
nc=size(xcs,2); nPg=nc*nc*6; |
18 |
|
|
%- add the 2 missing corners: |
19 |
|
|
xcg=reshape(xcg,[nPg 1]); ycg=reshape(ycg,[nPg 1]); |
20 |
|
|
xcg(nPg+1)=xcg(1); ycg(nPg+1)=ycg(1+2*nc); |
21 |
|
|
xcg(nPg+2)=xcg(1+3*nc); ycg(nPg+2)=ycg(1); |
22 |
|
|
|
23 |
|
|
psiLineF=[Rac,'psiLine_N2S_cs',int2str(nc)]; |
24 |
|
|
%psiLineF=['psiLine_N2S_cs',int2str(nc)]; |
25 |
|
|
|
26 |
|
|
fprintf([' load psiLine from file=',psiLineF,' ...']); |
27 |
|
|
%- load : 'psi_N','psi_C','psi_P','psiUV','psi_T' : |
28 |
|
|
load(psiLineF); |
29 |
|
|
fprintf(' done\n'); |
30 |
|
|
%- load the grid dx,dy , convert to 10^6 m : |
31 |
|
|
%dxg=rdmds([Rac,'DXG']); dxg=dxg*1.e-6; |
32 |
|
|
%dyg=rdmds([Rac,'DYG']); dyg=dyg*1.e-6; |
33 |
|
|
|
34 |
|
|
end |
35 |
|
|
|
36 |
|
|
if krd > 1, |
37 |
|
|
%-- load Vertically Integrated transport : |
38 |
|
|
rac='output_files/'; |
39 |
|
|
%nit=72020; |
40 |
|
|
%uu=rdmds([rac,'hUtave'],nit); |
41 |
|
|
%vv=rdmds([rac,'hVtave'],nit); |
42 |
|
|
namf=[rac,'UintTrsYr.bin']; |
43 |
|
|
fprintf([' read file: ',namf,' ...']); |
44 |
|
|
ut=rdda(namf,[nc*6*nc 1],1,'real*8','b'); |
45 |
|
|
fprintf(' done\n'); |
46 |
|
|
namf=[rac,'VintTrsYr.bin']; |
47 |
|
|
fprintf([' read file: ',namf,' ...']); |
48 |
|
|
vt=rdda(namf,[nc*6*nc 1],1,'real*8','b'); |
49 |
|
|
fprintf(' done\n'); |
50 |
|
|
end |
51 |
|
|
|
52 |
|
|
%- compute detph integrated flow : |
53 |
|
|
%delR=[50 70 100 140 190 240 290 340 390 440 490 540 590 640 690]; % = delR in m |
54 |
|
|
%ut=reshape(uu,nPg,nr); vt=reshape(vv,nPg,nr); |
55 |
|
|
%dxg=reshape(dxg,nPg,1); dyg=reshape(dyg,nPg,1); |
56 |
|
|
%var=dyg*delR; ut=var.*ut; ut=sum(ut,2); |
57 |
|
|
%var=dxg*delR; vt=var.*vt; vt=sum(vt,2); |
58 |
|
|
|
59 |
|
|
%- compute Barotropic Stream Function : |
60 |
|
|
psiNx=size(psi_C,1);psiNy=size(psi_C,2); nPc2=nPg+2; |
61 |
|
|
psiBa=zeros(nPc2,1); |
62 |
|
|
ufac=rem(psi_T,2) ; vfac=fix(psi_T/2) ; |
63 |
|
|
fprintf(' start computing stream-function ...'); |
64 |
|
|
if Tprt, TimeT0=clock; end |
65 |
|
|
for j=2:psiNy, |
66 |
|
|
for i=1:psi_N(j), |
67 |
|
|
i1=psi_C(i,j); i0=psi_P(i,j); i2=psiUV(i,j); |
68 |
|
|
psiBa(i1)=psiBa(i0)+ufac(i,j)*ut(i2)+vfac(i,j)*vt(i2); |
69 |
|
|
end |
70 |
|
|
end |
71 |
|
|
if Tprt, TimeT1=clock; end |
72 |
|
|
fprintf(' done\n'); |
73 |
|
|
if Tprt, fprintf(' Total time: %9.6f \n',etime(TimeT1,TimeT0) ); end |
74 |
|
|
|
75 |
|
|
if kwr==1, |
76 |
|
|
namfil='Psi_Horiz.bin'; |
77 |
|
|
fprintf([' write to file: ',namfil,' ...']); |
78 |
|
|
fid=fopen(namfil,'w','b'); fwrite(fid,psiBa,'real*8'); fclose(fid); |
79 |
|
|
fprintf(' done\n'); |
80 |
|
|
end |
81 |
|
|
|
82 |
|
|
if kgr > 0, |
83 |
|
|
shift=0; ccB=[0 0]; cbV=3; AxBx=[-180 180 -90 90]; kEnv=1; |
84 |
|
|
ccB=[-240 240]; |
85 |
|
|
%AxBx=[-210 210 30 90]; |
86 |
|
|
grph_CSz(psiBa,xcs,ycs,xcg,ycg,ccB(1),ccB(2),shift,cbV,AxBx,kEnv) |
87 |
|
|
end |