/[MITgcm]/MITgcm/verification/offline_exf_seaice/input/gendata.m
ViewVC logotype

Annotation of /MITgcm/verification/offline_exf_seaice/input/gendata.m

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


Revision 1.5 - (hide annotations) (download)
Fri Mar 7 22:31:46 2014 UTC (10 years, 3 months ago) by jmc
Branch: MAIN
CVS Tags: checkpoint64y, checkpoint64x, checkpoint64z, checkpoint64w, checkpoint64v, checkpoint66g, checkpoint66f, checkpoint66e, checkpoint66d, checkpoint66c, checkpoint66b, checkpoint66a, checkpoint66o, checkpoint66n, checkpoint66m, checkpoint66l, checkpoint66k, checkpoint66j, checkpoint66i, checkpoint66h, checkpoint65z, checkpoint65x, checkpoint65y, checkpoint65r, checkpoint65s, checkpoint65p, checkpoint65q, checkpoint65v, checkpoint65w, checkpoint65t, checkpoint65u, checkpoint65j, checkpoint65k, checkpoint65h, checkpoint65i, checkpoint65n, checkpoint65o, checkpoint65l, checkpoint65m, checkpoint65b, checkpoint65c, checkpoint65a, checkpoint65f, checkpoint65g, checkpoint65d, checkpoint65e, checkpoint65, HEAD
Changes since 1.4: +0 -1 lines
forgot a switch

1 jmc 1.1
2 jmc 1.4 kwr=1; kprt=0;
3 jmc 1.3 nx=80; ny=42; nr=3; nt=1;
4 jmc 1.1
5     xc=[1:nx]; xc=xc-mean(xc);
6     yc=[1:ny]-.5;
7    
8     %------------------------------------------------------
9    
10     windx=10.;
11     H0=-100.;
12    
13     namf='channel.bin';
14     depth=H0*ones(nx,ny); depth(:,1)=0.;
15 gforget 1.2 if kwr > 0,
16     fprintf('write to file: %s\n',namf);
17     fid=fopen(namf,'w','b'); fwrite(fid,depth,'real*8'); fclose(fid);
18     end
19 jmc 1.1
20     namf='bathy_3c.bin';
21     msk=abs(xc)'*ones(1,ny)+ones(nx,1)*yc;
22     depth=H0*ones(nx,ny); depth(:,1)=0.;
23     depth(find(msk < 24))=0.;
24     y2d=ones(nx,1)*yc;
25     depth(find(y2d > ny/2))=H0;
26     if kwr > 0,
27     fprintf('write to file: %s\n',namf);
28     fid=fopen(namf,'w','b'); fwrite(fid,depth,'real*8'); fclose(fid);
29     end
30    
31     namf='windx.bin';
32     wnd=windx*ones(nx,ny,nt);
33     if kwr > 0,
34     fprintf('write to file: %s\n',namf);
35     fid=fopen(namf,'w','b'); fwrite(fid,wnd,'real*8'); fclose(fid);
36     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='const+40.bin';
61     %u0=0.4;
62     %var=u0*ones(nx,ny);
63     %if kwr > 0,
64     % fprintf('write to file: %s\n',namf);
65     % fid=fopen(namf,'w','b'); fwrite(fid,var,'real*8'); fclose(fid);
66     %end
67    
68     %namf='const-10.bin';
69     %v0=-0.1;
70     %var=v0*ones(nx,ny);
71     %if kwr > 0,
72     % fprintf('write to file: %s\n',namf);
73     % fid=fopen(namf,'w','b'); fwrite(fid,var,'real*8'); fclose(fid);
74     %end
75     %------------------------------------------------------
76    
77 jmc 1.3 namf='ice0_area.bin'; iceC0=1.;
78     iceConc=iceC0*ones(nx,ny); iceConc(:,1)=0;
79     iceConc(:,2)=0.00*iceC0;
80     iceConc(:,3)=0.10*iceC0;
81     iceConc(:,end) =0.00*iceC0;
82     iceConc(:,end-1)=0.01*iceC0;
83     if kwr > 0,
84 gforget 1.2 fprintf('write to file: %s\n',namf);
85 jmc 1.3 fid=fopen(namf,'w','b'); fwrite(fid,iceConc,'real*8'); fclose(fid);
86 gforget 1.2 end
87    
88 jmc 1.3 namf='ice0_heff.bin'; iceH0=0.2;
89     iceVol=iceConc*iceH0;
90     if kwr > 0,
91 gforget 1.2 fprintf('write to file: %s\n',namf);
92 jmc 1.3 fid=fopen(namf,'w','b'); fwrite(fid,iceVol,'real*8'); fclose(fid);
93 gforget 1.2 end
94    
95     %------------------------------------------------------
96    
97 jmc 1.1 dsw0=100;
98     namf=['dsw_',int2str(dsw0),'.bin'];
99     fld=dsw0*ones(nx,ny,nt);
100     if kwr > 0,
101     fprintf('write to file: %s\n',namf);
102     fid=fopen(namf,'w','b'); fwrite(fid,fld,'real*8'); fclose(fid);
103     end
104    
105 jmc 1.3 dlw0=250;
106 jmc 1.1 namf=['dlw_',int2str(dlw0),'.bin'];
107     fld=dlw0*ones(nx,ny,nt);
108     if kwr > 0,
109     fprintf('write to file: %s\n',namf);
110     fid=fopen(namf,'w','b'); fwrite(fid,fld,'real*8'); fclose(fid);
111     end
112    
113     cel2K=273.15; dtx=4; %- dtx = amplitude of air temp variations in X-dir
114     ta_x=cel2K + dtx*sin(pi*(1+2*xc'/nx));
115     ta=repmat(ta_x,[1 ny nt]);
116     namf=['tair_',int2str(dtx),'x.bin'];
117     if kwr > 0,
118     fprintf('write to file: %s\n',namf);
119     fid=fopen(namf,'w','b'); fwrite(fid,ta,'real*8'); fclose(fid);
120     end;
121    
122     cvapor_fac = 640380.000 ;
123     cvapor_exp = 5107.400 ;
124     atmrho = 1.200 ;
125     rh=70; %- specific humid <--> 70.% relative humid
126     tmpbulk = cvapor_fac*exp(-cvapor_exp./ta_x);
127     qa_x = (rh/100.)*tmpbulk/atmrho ;
128     qa=repmat(qa_x,[1 ny nt]);
129     namf=['qa',int2str(rh),'_',int2str(dtx),'x.bin'];
130     if kwr > 0,
131     fprintf('write to file: %s\n',namf);
132     fid=fopen(namf,'w','b'); fwrite(fid,qa,'real*8'); fclose(fid);
133     end;
134    
135     %- salinity
136     sCst=30;
137     so=sCst*ones(nx,ny,nt);
138     namf='socn.bin';
139     %if kwr > 0,
140     % fprintf('write to file: %s\n',namf);
141     % fid=fopen(namf,'w','b'); fwrite(fid,so,'real*8'); fclose(fid);
142     %end;
143    
144     muTf = 5.4e-2;
145 jmc 1.3 tfreeze=-muTf*sCst;
146 jmc 1.1 fprintf('T-freeze = %10.6f\n',tfreeze);
147 jmc 1.3 %- parabolic profile in Y, max @ j=4, min @ j=ny, amplitude=1.K
148     to_y=(yc-3.5)/(ny-4);
149 jmc 1.1 to_y=tfreeze+0.5-to_y.*to_y;
150 jmc 1.3 mnV=min(to_y); MxV=max(to_y); Avr=mean(to_y(2:end));
151     fprintf(' SST* av,mn,Mx: %9.6f , %9.6f , %9.6f , %9.6f\n',Avr,mnV,MxV,MxV-mnV);
152 jmc 1.1 to=repmat(to_y,[nx 1 nt]);
153     namf='tocn.bin';
154     if kwr > 0,
155     fprintf('write to file: %s\n',namf);
156     fid=fopen(namf,'w','b'); fwrite(fid,to,'real*8'); fclose(fid);
157     end;
158    
159     %-- make some plots to check: ----------------
160    
161     hScal=[-1.1 0.1]*abs(H0);
162     figure(1); clf;
163     subplot(211);
164     var=depth;
165     imagesc(xc,yc,var'); set(gca,'YDir','normal');
166     %caxis(hScal);
167     %change_colmap(-1);
168     colorbar;
169     title('Depth [m]');
170    
171     subplot(413);
172     var=depth;
173     j1=2;
174     j2=ny/2;
175     j3=j2+1;
176     plot(xc,var(:,j1),'k-')
177     hold on; j=j+1;
178     plot(xc,var(:,j2),'ro-')
179     plot(xc,var(:,j3),'b-')
180     hold off;
181     axis([-nx/2 nx/2 hScal]);
182     grid
183     legend(int2str(j1),int2str(j2),int2str(j3));
184     title('Depth @ j= cst');
185    
186     subplot(414);
187     i=nx/2;
188     plot(yc,var(i,:),'k-')
189     axis([0 ny H0*1.1 -H0*.1]);
190     grid
191     title(['Depth @ i=',int2str(i)]);
192    
193     %--
194     dewPt=(qa_x*atmrho)/cvapor_fac;
195     dewPt=-cvapor_exp./log(dewPt);
196    
197     figure(2);clf;
198 jmc 1.4 subplot(311)
199     P(1)=plot(xc,ta_x-cel2K,'r-'); hold on;
200     P(2)=plot(xc,dewPt-cel2K,'b-');
201     P(3)=plot(xc,tfreeze*ones(nx,1),'k-');
202     set(P,'LineWidth',1);
203     hold off; AA=axis;
204     axis([-nx/2 nx/2 AA(3:4)]);
205 jmc 1.1 legend('ta','dew');
206     grid
207 jmc 1.4 xlabel('X')
208     title(['Air Temp (^oC): del-Temp-X = ',int2str(dtx),' , RH= ',int2str(rh)]);
209     subplot(312)
210     P(1)=plot(yc,to_y,'b-'); hold on;
211     P(2)=plot(yc,tfreeze*ones(ny,1),'k-');
212     set(P,'LineWidth',1);
213     hold off; AA=axis;
214     L=line([1 1],AA(3:4)); set(L,'LineWidth',2.,'Color',[0 0 0]);
215     axis([0 ny AA(3:4)]);
216 jmc 1.1 grid
217 jmc 1.4 xlabel('Y')
218 jmc 1.1 title('Ocean Temp ^oC');
219 jmc 1.3
220 jmc 1.4 subplot(313)
221     var=iceConc(1,:);
222     P(1)=semilogy(yc,var,'b-x'); hold on;
223     %plot(yc,var,'b-x'); hold on;
224     var=iceVol(1,:);
225     P(2)=semilogy(yc,var,'r-x');
226     %plot(yc,var,'r-x');
227     set(P,'LineWidth',1);
228     hold off; AA=axis;
229     L=line([1 1],AA(3:4)); set(L,'LineWidth',2.,'Color',[0 0 0]);
230     axis([0 ny [0 2]*iceC0]);
231     grid
232     xlabel('Y')
233     legend('iceC','hEff','Location','South');
234     title('Initial ice in Channel : y-section');
235     %-----
236     if kprt == 1, f=2;
237     namfig=sprintf('forcing_%2.2i',f);
238     fprintf([' print fig= %2i to file: ',namfig,' '],f);
239     set(f,'PaperOrientation','portrait')
240     %set(f,'PaperPosition',[0.25 1.5 6. 8.]);
241     set(f,'PaperPosition',[0.25 1.5 5.25 7.]);
242     print(f,'-depsc2',namfig); fprintf('\n');
243     end
244     %-----
245 jmc 1.3
246     figure(3);clf;
247     subplot(311)
248     var=iceConc; ccB=[-1 12]/10;
249     imagesc(xc,yc,var'); set(gca,'YDir','normal');
250     caxis(ccB);
251     %change_colmap(-1);
252     colorbar;
253     title('Ice Concentration in Channel');
254    
255     subplot(312)
256     var=iceVol; ccB=[-1 12]/50;
257     imagesc(xc,yc,var'); set(gca,'YDir','normal');
258     caxis(ccB);
259     %change_colmap(-1);
260     colorbar;
261     title('Effective ice thickness in Channel');
262    
263     subplot(313)
264     var=iceConc(1,:);
265     %plot(yc,var,'b-x'); hold on;
266     semilogy(yc,var,'b-x'); hold on;
267     var=iceVol(1,:);
268     %plot(yc,var,'r-x'); hold off;
269     semilogy(yc,var,'r-x'); hold on;
270     AA=axis; axis([0 ny [0 2]*iceC0]);
271     grid
272     legend('iceC','hEff','Location','South');
273     title('Initial ice in Channel : y-section');
274 jmc 1.4 %-----
275    
276     return

  ViewVC Help
Powered by ViewVC 1.1.22