/[MITgcm]/MITgcm_contrib/verification_other/offline_cheapaml/input/gendata.m
ViewVC logotype

Annotation of /MITgcm_contrib/verification_other/offline_cheapaml/input/gendata.m

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


Revision 1.2 - (hide annotations) (download)
Tue Jun 11 01:23:52 2013 UTC (12 years, 1 month ago) by jmc
Branch: MAIN
CVS Tags: checkpoint64k, checkpoint64x, checkpoint64z, checkpoint64o, checkpoint64p, checkpoint64r, checkpoint64w, checkpoint64v, checkpoint66g, checkpoint66f, checkpoint66e, checkpoint66d, checkpoint66c, checkpoint66b, checkpoint66a, checkpoint66o, checkpoint66n, checkpoint66m, checkpoint66l, checkpoint66k, checkpoint66j, checkpoint66i, checkpoint66h, checkpoint64m, checkpoint65z, checkpoint65x, checkpoint65y, checkpoint65r, checkpoint65s, checkpoint65p, checkpoint65q, checkpoint65v, checkpoint65w, checkpoint65u, checkpoint65j, checkpoint67a, checkpoint67b, checkpoint65i, checkpoint67d, checkpoint65m, checkpoint65c, checkpoint65a, checkpoint65f, checkpoint65e, checkpoint64i, checkpoint65h, checkpoint65, checkpoint64j, checkpoint65n, HEAD
Changes since 1.1: +14 -4 lines
- new wind V-component wind (to produce small uniform convergence)
- rename windx.bin --> windx_10ms.bin ; add more diagnostics.

1 jmc 1.1
2     kwr=0;
3     %kwr=-1;
4     nx=80; ny=42; nr=3; nt=1;
5    
6     xc=[1:nx]; xc=xc-mean(xc);
7 jmc 1.2 yc=[1:ny]-.5; ymid=mean(yc);
8     yv=yc-0.5;
9 jmc 1.1
10     %------------------------------------------------------
11    
12     windx=10.;
13     H0=-100.;
14    
15     namf='channel.bin';
16     depth=H0*ones(nx,ny); depth(:,1)=0.; depth(:,ny)=0.;
17     if kwr > 0,
18     fprintf('write to file: %s\n',namf);
19     fid=fopen(namf,'w','b'); fwrite(fid,depth,'real*8'); fclose(fid);
20     end
21    
22 jmc 1.2 namf=['windx_',int2str(windx),'ms.bin'];
23     uwind=windx*ones(nx,ny,nt);
24 jmc 1.1 if kwr > 0,
25     fprintf('write to file: %s\n',namf);
26 jmc 1.2 fid=fopen(namf,'w','b'); fwrite(fid,uwind,'real*8'); fclose(fid);
27     end
28    
29     namf='windy_conv.bin';
30     dvdy=-1.e-6*5.e+3; %- uniform convergence: wWind = 10^-6 m/s
31     yy=yv-ymid; vwind=dvdy*yy; vwind(1)=0;
32     fld=ones(nx,1)*vwind;
33     if kwr > 0,
34     fprintf('write to file: %s\n',namf);
35     fid=fopen(namf,'w','b'); fwrite(fid,fld,'real*8'); fclose(fid);
36 jmc 1.1 end
37    
38     %- file name convention: "const_{xx}.bin" <-> uniform value = xx (in percent)
39     namf='const_00.bin';
40     fld=0*ones(nx,ny,nt);
41     if kwr > 0,
42     fprintf('write to file: %s\n',namf);
43     fid=fopen(namf,'w','b'); fwrite(fid,fld,'real*8'); fclose(fid);
44     end
45    
46     namf='const100.bin'; w0=1.;
47     %var=w0*ones(nx,ny);
48     %if kwr > 0,
49     % fprintf('write to file: %s\n',namf);
50     % fid=fopen(namf,'w','b'); fwrite(fid,var,'real*8'); fclose(fid);
51     %end
52    
53     namf='const+20.bin'; w0=0.2;
54     var=w0*ones(nx,ny);
55     if kwr > 0,
56     fprintf('write to file: %s\n',namf);
57     fid=fopen(namf,'w','b'); fwrite(fid,var,'real*8'); fclose(fid);
58     end
59    
60     namf='ice0_area.bin'; iceC0=1.;
61     iceConc=zeros(nx,ny);
62     iceConc(21:60,:)=iceC0;
63     iceConc(:,1)=0; iceConc(:,ny)=0;
64     if kwr > 0,
65     fprintf('write to file: %s\n',namf);
66     fid=fopen(namf,'w','b'); fwrite(fid,iceConc,'real*8'); fclose(fid);
67     end
68    
69     namf='ice0_heff.bin'; iceH0=0.2;
70     %iceVol=iceConc*iceH0;
71     %if kwr > 0,
72     % fprintf('write to file: %s\n',namf);
73     % fid=fopen(namf,'w','b'); fwrite(fid,iceVol,'real*8'); fclose(fid);
74     %end
75    
76     %------------------------------------------------------
77    
78     dsw0=70; dswA=30;
79     ymid=(yc(2)+yc(ny))/2;
80     ysn=yc-ymid ; ysn=ysn/(yc(ny)-yc(2)); ysn=sin(pi*ysn);
81    
82     namf=['dsw_',int2str(dsw0),'y.bin'];
83     dsw_y=dsw0+dswA*ysn;
84     fld=ones(nx,1)*dsw_y;
85     if kwr > 0,
86     fprintf('write to file: %s\n',namf);
87     fid=fopen(namf,'w','b'); fwrite(fid,fld,'real*8'); fclose(fid);
88     end
89    
90     dlw0=270; dlwA=-20;
91     namf=['dlw_',int2str(dlw0),'y.bin'];
92     dlw_y=dlw0+dlwA*ysn;
93     fld=ones(nx,1)*dlw_y;
94     if kwr > 0,
95     fprintf('write to file: %s\n',namf);
96     fid=fopen(namf,'w','b'); fwrite(fid,fld,'real*8'); fclose(fid);
97     end
98    
99     cel2K=273.15; taC=-10;
100     %ta0=cel2K+taC;
101     %- cheapAML works in deg.C
102     ta0=taC;
103     ta=ta0*ones(nx,ny);
104     namf=['tair_',int2str(taC),'.bin'];
105     if kwr > 0,
106     fprintf('write to file: %s\n',namf);
107     fid=fopen(namf,'w','b'); fwrite(fid,ta,'real*8'); fclose(fid);
108     end;
109    
110     cvapor_fac = 640380.000 ;
111     cvapor_exp = 5107.400 ;
112     atmrho = 1.200 ;
113     rh=70; %- specific humid <--> 70.% relative humid
114     taK=cel2K+taC;
115     tmpbulk = cvapor_fac*exp(-cvapor_exp./taK);
116     qa0= (rh/100.)*tmpbulk/atmrho ;
117     qa=qa0*ones(nx,ny);
118     namf=['qa',int2str(rh),'_',int2str(taC),'.bin'];
119     if kwr > 0,
120     fprintf('write to file: %s\n',namf);
121     fid=fopen(namf,'w','b'); fwrite(fid,qa,'real*8'); fclose(fid);
122     end;
123    
124     %- salinity
125     sCst=30;
126     so=sCst*ones(nx,ny,nt);
127     namf='socn.bin';
128     %if kwr > 0,
129     % fprintf('write to file: %s\n',namf);
130     % fid=fopen(namf,'w','b'); fwrite(fid,so,'real*8'); fclose(fid);
131     %end;
132    
133     muTf = 5.4e-2;
134     tfreeze=-muTf*sCst;
135     fprintf('T-freeze = %10.6f\n',tfreeze);
136     dtx=1; %- dtx = amplitude of relaxing ocean temp variations in X-dir
137     %to_x=tfreeze + dtx*sin(pi*(1+2*xc'/nx)) + dtx;
138     to_x=tfreeze + dtx*sin(pi*(2*xc'/nx)) + dtx;
139     mnV=min(to_x); MxV=max(to_x); Avr=mean(to_x);
140     fprintf(' SST* av,mn,Mx: %9.6f , %9.6f , %9.6f , %9.6f\n',Avr,mnV,MxV,MxV-mnV);
141     to=repmat(to_x,[1 ny nt]);
142     namf=['tocn_',int2str(dtx),'x.bin'];
143     if kwr > 0,
144     fprintf('write to file: %s\n',namf);
145     fid=fopen(namf,'w','b'); fwrite(fid,to,'real*8'); fclose(fid);
146     end;
147    
148     %-- make some plots to check: ----------------
149    
150     hScal=[-1.1 0.1]*abs(H0);
151     figure(1); clf;
152     subplot(211);
153     var=depth;
154     imagesc(xc,yc,var'); set(gca,'YDir','normal');
155     %caxis(hScal);
156     %change_colmap(-1);
157     colorbar;
158     title('Depth [m]');
159    
160     subplot(413);
161     var=depth;
162     j1=2;
163     j2=ny/2;
164     j3=j2+1;
165     plot(xc,var(:,j1),'k-')
166     hold on; j=j+1;
167     plot(xc,var(:,j2),'ro-')
168     plot(xc,var(:,j3),'b-')
169     hold off;
170     axis([-nx/2 nx/2 hScal]);
171     grid
172     legend(int2str(j1),int2str(j2),int2str(j3));
173     title('Depth @ j= cst');
174    
175     subplot(414);
176     i=nx/2;
177     plot(yc,var(i,:),'k-')
178     axis([0 ny H0*1.1 -H0*.1]);
179     grid
180     title(['Depth @ i=',int2str(i)]);
181    
182     %--
183    
184     figure(2);clf;
185     subplot(211)
186     plot(xc,to_x,'r-'); hold on;
187     %plot(xc,dewPt-cel2K,'b-');
188     plot(xc,tfreeze*ones(nx,1),'k-');
189     hold off;
190     AA=axis; axis([-nx/2 nx/2 AA(3:4)]);
191     legend('To','Tfrz');
192     grid
193     title('Ocean Temp ^oC');
194     %title(['del-Temp-X= ',int2str(dtx),' ; RH= ',int2str(rh),' ; Air Temp (^oC)']);
195     subplot(212)
196     plot(yc,dsw_y,'r-'); hold on;
197     plot(yc,dlw_y-200,'b-');
198     hold off;
199     %AA=axis; axis([0 ny AA(3:4)]);
200     axis([0 ny 30 110]);
201     grid
202     legend('sw','lw','Location','South');
203     title('Downward SW and LW (-200) [W/m^2]');
204    
205     %--
206    
207     figure(3);clf;
208     subplot(311)
209     var=iceConc; ccB=[-1 12]/10;
210     imagesc(xc,yc,var'); set(gca,'YDir','normal');
211     caxis(ccB);
212     %change_colmap(-1);
213     colorbar;
214     title('Ice Concentration in Channel');
215    
216     iceVol=iceH0*iceConc;
217     subplot(312)
218     var=iceVol; ccB=[-1 12]/50;
219     imagesc(xc,yc,var'); set(gca,'YDir','normal');
220     caxis(ccB);
221     %change_colmap(-1);
222     colorbar;
223     title('Effective ice thickness in Channel');
224    
225     subplot(313)
226     var=iceConc(1,:);
227     %plot(yc,var,'b-x'); hold on;
228     semilogy(yc,var,'b-x'); hold on;
229     var=iceVol(1,:);
230     %plot(yc,var,'r-x'); hold off;
231     semilogy(yc,var,'r-x'); hold on;
232     AA=axis; axis([0 ny [0 2]*iceC0]);
233     grid
234     legend('iceC','hEff','Location','South');
235     title('Initial ice in Channel : y-section');
236     %--

  ViewVC Help
Powered by ViewVC 1.1.22