/[MITgcm]/MITgcm_contrib/high_res_cube/matlab/use_psiLine.m
ViewVC logotype

Annotation of /MITgcm_contrib/high_res_cube/matlab/use_psiLine.m

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph


Revision 1.1 - (hide annotations) (download)
Wed Mar 28 18:07:27 2007 UTC (18 years, 3 months ago) by dimitri
Branch: MAIN
CVS Tags: HEAD
Files contributed by JM for computation of barotropic streamfunction on cs510.

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

  ViewVC Help
Powered by ViewVC 1.1.22