1 |
% |
% |
2 |
% [] = A_compute_potential_density(SNAPSHOT) |
% [ST,LON,LAT,DPT] = A_compute_potential_density(SNAPSHOT) |
3 |
% |
% |
4 |
% For a time snapshot, this program computes the |
% For a time snapshot, this program computes the |
5 |
% 3D potential density from potential temperature and salinity. |
% 3D potential density from potential temperature and salinity. |
6 |
% THETA and SALTanom are supposed to be defined on the same |
% THETA and SALTanom are supposed to be defined on the same |
7 |
% domain and grid. |
% domain and grid. |
8 |
% |
% SALTanom is by default a salinity anomaly vs 35. |
9 |
|
% If not, (is absolute value) set the global variable is_SALTanom to 0 |
10 |
|
% |
11 |
% Files names are: |
% Files names are: |
12 |
% INPUT: |
% INPUT: |
13 |
% ./netcdf-files/<SNAPSHOT>/<netcdf_THETA>.<netcdf_domain>.<netcdf_suff> |
% ./netcdf-files/<SNAPSHOT>/<netcdf_THETA>.<netcdf_domain>.<netcdf_suff> |
47 |
ncSALTa = netcdf(ferfile,'nowrite'); |
ncSALTa = netcdf(ferfile,'nowrite'); |
48 |
SALTavariables = var(ncSALTa); |
SALTavariables = var(ncSALTa); |
49 |
|
|
50 |
|
global is_SALTanom |
51 |
|
if exist('is_SALTanom') |
52 |
|
if is_SALTanom == 1 |
53 |
|
bS = 35; |
54 |
|
else |
55 |
|
bS = 0; |
56 |
|
end |
57 |
|
end |
58 |
|
|
59 |
%% Gridding: |
%% Gridding: |
60 |
% Don't care about the grid here ! |
% Don't care about the grid here ! |
61 |
% SALTanom and THETA are normaly defined on the same grid |
% SALTanom and THETA are normaly defined on the same grid |
85 |
for iz = 1 : nz |
for iz = 1 : nz |
86 |
if toshow,disp(strcat('Compute potential density at level:',num2str(iz),'/',num2str(nz)));end |
if toshow,disp(strcat('Compute potential density at level:',num2str(iz),'/',num2str(nz)));end |
87 |
|
|
88 |
S = SALTavariables{4}(iz,:,:) + 35; % Move the anom to an absolute field |
S = SALTavariables{4}(iz,:,:) + bS; % Evantualy move the anom to an absolute field |
89 |
T = THETAvariables{4}(iz,:,:); |
T = THETAvariables{4}(iz,:,:); |
90 |
SIGMATHETA(iz,:,:) = densjmd95(S,T,zeros(ny,nx)) - 1000; |
SIGMATHETA(iz,:,:) = densjmd95(S,T,zeros(ny,nx)) - 1000; |
91 |
|
|
92 |
% Eventualy make a plot of the field: |
% Eventualy make a plot of the field: |
93 |
if 0 |
if 0 & iz==1 |
94 |
clf;pcolor(squeeze(SIGMATHETA(iz,:,:))); |
clf;pcolor(squeeze(SIGMATHETA(iz,:,:))); |
95 |
shading interp;caxis([10 40]);colorbar |
shading interp;caxis([20 30]);colorbar |
96 |
drawnow |
drawnow |
97 |
M(iz)=getframe; % To make a video |
%M(iz)=getframe; % To make a video |
98 |
end %if1 |
end %if1 |
99 |
end %for iz |
end %for iz |
100 |
|
|
156 |
|
|
157 |
nc=close(nc); |
nc=close(nc); |
158 |
|
|
159 |
|
|
160 |
|
% Output: |
161 |
|
switch nargout |
162 |
|
case 1 |
163 |
|
varargout(1) = SIGMATHETA; |
164 |
|
case 2 |
165 |
|
varargout(1) = SIGMATHETA; |
166 |
|
varargout(2) = lon; |
167 |
|
case 3 |
168 |
|
varargout(1) = SIGMATHETA; |
169 |
|
varargout(2) = lon; |
170 |
|
varargout(3) = lat; |
171 |
|
case 4 |
172 |
|
varargout(1) = SIGMATHETA; |
173 |
|
varargout(2) = lon; |
174 |
|
varargout(3) = lat; |
175 |
|
varargout(4) = dpt; |
176 |
|
end |