| 1 |
|
| 2 |
kwr=0; |
| 3 |
%kwr=1; |
| 4 |
|
| 5 |
%gDir='../run_26l/'; |
| 6 |
%gDir='../cs_grid/'; |
| 7 |
gDir='../run/'; |
| 8 |
G=load_grid(gDir,0); |
| 9 |
|
| 10 |
nx=G.dims(1); ny=G.dims(2); nc=ny; |
| 11 |
xc=G.xC; yc=G.yC; xg=G.xG; yg=G.yG; |
| 12 |
|
| 13 |
ccB=[0 0]; shift=-1; cbV=1; AxBx=[-180 180 -90 90]; kEnv=0; |
| 14 |
y1d=[-89.5:1:90]; |
| 15 |
%figure(1);clf; |
| 16 |
%imagesc(msk'); set(gca,'YDir','normal'); |
| 17 |
%colorbar |
| 18 |
|
| 19 |
%-- make SST field (cos shape, no noise) from grid-output file YC: |
| 20 |
yy=yc*pi/90; sst1=9+19*cos(yy); |
| 21 |
%sst=273.15+sst1; |
| 22 |
%fname='SST_cos0.bin'; |
| 23 |
yy=y1d*pi/90; sst1_1d=9+19*cos(yy); %- just for doing a plot |
| 24 |
|
| 25 |
%-- follow APE profile: |
| 26 |
yy=yc*pi/180; |
| 27 |
var=sin(1.5*yy); var=1.-var.*var; |
| 28 |
var(find(abs(yc) > 60.))=0.; sst2=27*var; |
| 29 |
%sst=273.15+sst2; |
| 30 |
%fname='SST_APE_1.bin'; |
| 31 |
%- just for doing a plot |
| 32 |
yy=y1d*pi/180; |
| 33 |
var=sin(1.5*yy); var=1.-var.*var; |
| 34 |
var(find(abs(y1d) > 60.))=0.; sst2_1d=27*var; |
| 35 |
|
| 36 |
%-- close to Symetric zonally-average current SST: |
| 37 |
yyp=abs(yc/90); phi=yyp.^3.; |
| 38 |
sst3=27.8*exp(-7*phi) - 1 ; |
| 39 |
sst=273.15+sst3; |
| 40 |
fname='SST_symEx3.bin'; |
| 41 |
%- just for doing a plot |
| 42 |
yyp=abs(y1d/90); phi=yyp.^3.; |
| 43 |
sst3_1d=27.8*exp(-7*phi) - 1 ; |
| 44 |
|
| 45 |
if kwr > 0, |
| 46 |
fid=fopen(fname,'w','b'); fwrite(fid,sst,'real*8'); fclose(fid); |
| 47 |
fprintf(['write file: ',fname,'\n']); |
| 48 |
end |
| 49 |
|
| 50 |
figure(1);clf; |
| 51 |
subplot(211) |
| 52 |
var=sst; |
| 53 |
grph_CS(var,xc,yc,xg,yg,ccB(1),ccB(2),shift,cbV,AxBx,kEnv); |
| 54 |
title('SST'); |
| 55 |
|
| 56 |
subplot(212) |
| 57 |
plot(y1d,sst1_1d,'g-'); hold on; |
| 58 |
plot(y1d,sst2_1d,'b-'); |
| 59 |
plot(y1d,sst3_1d,'r-'); |
| 60 |
hold off; |
| 61 |
axis([-90 90 -5 29]); grid |
| 62 |
legend('cos','ape','sym') |
| 63 |
title('SST ^oC'); |
| 64 |
%return |
| 65 |
|
| 66 |
%-- make Q-flux file (similar to the one used by Paul, hard coded |
| 67 |
% in original version of mixed_layer.f90) |
| 68 |
qflx_ampl=50. ; |
| 69 |
qflx_width=90.; |
| 70 |
yy=yc/qflx_width; |
| 71 |
yy=yy.*yy; |
| 72 |
qflx=qflx_ampl*(1-2*yy); |
| 73 |
qflx=qflx.*exp(-yy); |
| 74 |
|
| 75 |
figure(2);clf; |
| 76 |
var=qflx; |
| 77 |
grph_CS(var,xc,yc,xg,yg,ccB(1),ccB(2),shift,cbV,AxBx,kEnv); |
| 78 |
title('Q-flux [W/m^2]'); |
| 79 |
|
| 80 |
if kwr > 0, |
| 81 |
fname='Qflux_w90.bin'; |
| 82 |
fid=fopen(fname,'w','b'); fwrite(fid,qflx,'real*8'); fclose(fid); |
| 83 |
fprintf(['write file: ',fname,'\n']); |
| 84 |
end |
| 85 |
%return |
| 86 |
|
| 87 |
%-- make initial pot-temp field by adding noise to T(iter=0) output file: |
| 88 |
|
| 89 |
rDir=gDir; |
| 90 |
namf='T'; it=0; |
| 91 |
tini=rdmds([rDir,namf],it); |
| 92 |
nr=size(tini,3); |
| 93 |
|
| 94 |
var=rand([nx,ny]); var=var-mean(var(:)); |
| 95 |
yy=yc*pi/90; |
| 96 |
var=var.*(2+cos(yy))/3; |
| 97 |
figure(3);clf; |
| 98 |
%var=sst0; |
| 99 |
grph_CS(var,xc,yc,xg,yg,ccB(1),ccB(2),shift,cbV,AxBx,kEnv); |
| 100 |
|
| 101 |
noise=1.e-3; |
| 102 |
tini1=tini+noise*reshape(reshape(var,[nx*ny 1])*ones(1,nr),[nx ny nr]); |
| 103 |
%size(tini1) |
| 104 |
|
| 105 |
if kwr > 0, |
| 106 |
fname='ini_theta.bin'; |
| 107 |
fid=fopen(fname,'w','b'); fwrite(fid,tini1,'real*8'); fclose(fid); |
| 108 |
fprintf(['write file: ',fname,'\n']); |
| 109 |
end |
| 110 |
|
| 111 |
%- spec-humid : put constant Rel-Humid in the lowest troposphere |
| 112 |
relhum=0.8 ; pHum=800.e+2; |
| 113 |
relhum=0.4 ; pHum=700.e+2; |
| 114 |
khum=max(find(G.rC > pHum)); |
| 115 |
|
| 116 |
%- taken from AIM -> qsat in g/kg, pIn = normalised Pressure |
| 117 |
P0=1.e+5; pIn=G.rC/P0; |
| 118 |
qsat=calc_Qsat(1,pIn,tini); |
| 119 |
qsat=reshape(qsat,[nx ny nr]); |
| 120 |
qini=qsat*1.e-3*relhum; |
| 121 |
qini(:,:,khum+1:end)=0; |
| 122 |
|
| 123 |
figure(4);clf; |
| 124 |
pax=G.rC/100; %- in mb |
| 125 |
i1=1; j1=1; |
| 126 |
var=squeeze(qini(i1,j1,:)); |
| 127 |
plot(var,pax,'k-'); hold on; |
| 128 |
i1=nc/2; j1=nc/2; |
| 129 |
var=squeeze(qini(i1,j1,:)); |
| 130 |
plot(var,pax,'r-'); |
| 131 |
i1=nc*2.5; j1=nc*0.5; |
| 132 |
var=squeeze(qini(i1,j1,:)); |
| 133 |
plot(var,pax,'b-'); |
| 134 |
hold off |
| 135 |
set(gca,'YDir','reverse'); |
| 136 |
grid |
| 137 |
legend('mid','eq','pol'); |
| 138 |
title('Q-ini profile'); |
| 139 |
|
| 140 |
if kwr > 0, |
| 141 |
fname=['ini_specQ_',int2str(nr),'l.bin']; |
| 142 |
fid=fopen(fname,'w','b'); fwrite(fid,qini,'real*8'); fclose(fid); |
| 143 |
fprintf(['write file: ',fname,'\n']); |
| 144 |
end |
| 145 |
|
| 146 |
var=reshape(tini,[nx*ny nr]); |
| 147 |
var=mean(var); |
| 148 |
for n=1:ceil(nr/10) |
| 149 |
is=1+(n-1)*10; ie=min(nr,n*10); |
| 150 |
if n == 1, fprintf(' tRef='); else fprintf(' '); end |
| 151 |
fprintf(' %5.1f,',round(10*var(is:ie))/10); |
| 152 |
fprintf('\n') |
| 153 |
end |
| 154 |
|
| 155 |
return |
Parent Directory
| 
Revision Log
| 
Revision Graph