[Dart-dev] [7732] DART/trunk: This is the long-awaited merge of branches/ cam back onto the trunk.
nancy at ucar.edu
nancy at ucar.edu
Fri Mar 13 16:53:45 MDT 2015
Revision: 7732
Author: nancy
Date: 2015-03-13 16:53:44 -0600 (Fri, 13 Mar 2015)
Log Message:
-----------
This is the long-awaited merge of branches/cam back onto the trunk.
This was developed within the CESM1_2_1 framework.
The major upgrades and changes include the following.
More detailed comments are interspersed in the log since late 2013.
1) Addition of WACCM to the list of CAM variants which can be used in assimilations.
This required the addition of a new vertical coordinate; scale height ('log_invP'),
which has a small but widespread footprint.
2) Addition of CAM-SE to the list of CAM variants.
This required new interpolation routines and grid file management to handle
the logically non-rectangular grid (cubed sphere).
This consists mostly of a large new section of the module, but some scattered
changes as well.
3) Fixing the misinterpretation of DART KINDs vs TYPEs.
This was fixed independently on the trunk.
4) Code review leading to some structural and extensive cosmetic changes to bring
it into line with DART practices.
5) Rewriting of get_close_obs to handle the DART QC status better and
improve the clarity and logic, as well as accomodate items 1)-4).
6) Implementation of an improved damping algorithm for the top layers of the state,
where CAM applies enhanced diffusion and the assimilation should be weaker.
This algorithm has the right units, works for vert_coord = 'pressure' and
'log_invP', and allows more of the state to feel the observations.
7) Inclusion of scripting which will enable running multiple assimilation cycles
in a single job, which eliminates much of the waiting in the queue, and saves
~30% in computing time by running st_archive.sh offline as a single processor job.
This is accomplished through the creation of scripts ${CASE}.submit_cycles and
${CASE}.run_cycles by the new CESM1_2_1_setup_advanced. CESM1_2_1_setup_hybrid
now exists as the streamlined, 'beginners' script, which trades some of the
capabilities in CESM1_2_1_setup_advanced in favor of ease of use.
The CAM-FV version from branches/cam was tested against the trunk.
When the improvement in the choice of PS location in convert_vert is taken
into account, the results are identical.
No direct comparisons against the trunk are possible for WACCM or CAM-SE.
Those were compared qualitatively against the branches' CAM-FV assimilations.
A notable limitation discovered during testing is that filter can fail
for task counts > ~4200, due to an MPI imbalance. Filter completes all of
its work, but never returns from MPI_finalize. Work-arounds are available,
but are not included in the trunk.
At this time the scripts may only work for setting up the CAM5 version of WACCM.
The CAM4 version forecasts fail, probably due to input/forcing file incompatibility.
This is at the top of the 'fix-it' list.
Several additions and upgrades are planned for the near future on the trunk:
> many cosmetic changes to bring it into compliance with the DART style manual,
> inclusion of CAM-Chem, once the selective covariance mechanism is in place
elsewhere in DART.
In the medium term, cam/model_mod.f90 will be split into versions specialized
for various selected flavors of CAM, characterized by:
> dy-core (FV, SE, SE+refined_grid, eulerian, ...?),
> physics at vertical structure (CAM4 at 26, CAM5 at 30, WACCM#@(66,70), CAM#-Chem(s),...?).
Modified Paths:
--------------
DART/trunk/adaptive_inflate/adaptive_inflate_mod.f90
DART/trunk/diagnostics/matlab/linked_observations.m
DART/trunk/diagnostics/matlab/plot_rank_histogram.m
DART/trunk/diagnostics/threed_sphere/obs_diag.f90
DART/trunk/ensemble_manager/ensemble_manager_mod.f90
DART/trunk/filter/filter.f90
DART/trunk/location/threed_sphere/location_mod.f90
DART/trunk/models/cam/dart_to_cam.f90
DART/trunk/models/cam/model_mod.f90
DART/trunk/models/cam/model_mod.html
DART/trunk/models/cam/model_mod.nml
DART/trunk/models/cam/work/input.nml
DART/trunk/models/cam/work/path_names_cam_to_dart
DART/trunk/models/cam/work/path_names_closest_member_tool
DART/trunk/models/cam/work/path_names_create_fixed_network_seq
DART/trunk/models/cam/work/path_names_create_obs_sequence
DART/trunk/models/cam/work/path_names_dart_to_cam
DART/trunk/models/cam/work/path_names_fill_inflation_restart
DART/trunk/models/cam/work/path_names_filter
DART/trunk/models/cam/work/path_names_obs_common_subset
DART/trunk/models/cam/work/path_names_obs_diag
DART/trunk/models/cam/work/path_names_obs_seq_to_netcdf
DART/trunk/models/cam/work/path_names_obs_sequence_tool
DART/trunk/models/cam/work/path_names_perfect_model_obs
DART/trunk/models/cam/work/path_names_restart_file_tool
DART/trunk/mpi_utilities/mpi_utilities_mod.f90
DART/trunk/obs_kind/DEFAULT_obs_kind_mod.F90
DART/trunk/utilities/closest_member_tool.f90
Added Paths:
-----------
DART/trunk/models/cam/deprecated/advance_model.csh
DART/trunk/models/cam/deprecated/advance_model.retry.csh
DART/trunk/models/cam/deprecated/check_model.csh
DART/trunk/models/cam/deprecated/job.simple.csh
DART/trunk/models/cam/deprecated/run-cam.csh
DART/trunk/models/cam/doc/CESM_2013_Homme_DA.pptx
DART/trunk/models/cam/doc/DReaD_2013_cs_algorithms.pptx
DART/trunk/models/cam/doc/cubed_sphere_algorithm.pptx
DART/trunk/models/cam/doc/highest_state_p_Pa.pptx
DART/trunk/models/cam/doc/initial_vs_restart_files_for_CAM-FV.pdf
DART/trunk/models/cam/doc/restart_contents_SE
DART/trunk/models/cam/doc/restart_file_contents.txt
DART/trunk/models/cam/doc/roots_of_m_equation.pptx
DART/trunk/models/cam/matlab/BottomAnnotation.m
DART/trunk/models/cam/matlab/PlotCubedSpherePatches.m
DART/trunk/models/cam/matlab/cesm.stdout.273578.array
DART/trunk/models/cam/matlab/compare_hop_tests.m
DART/trunk/models/cam/matlab/compare_hop_tests_SE.m
DART/trunk/models/cam/matlab/err_cam-se.csh
DART/trunk/models/cam/matlab/get_info_index.m
DART/trunk/models/cam/matlab/my_getVar.m
DART/trunk/models/cam/matlab/project_cells.m
DART/trunk/models/cam/matlab/project_cells_make_array.csh
DART/trunk/models/cam/matlab/quadratic_eq_coeffs.m
DART/trunk/models/cam/matlab/se_plot_prep.csh
DART/trunk/models/cam/matlab/two_experiments_evol_setup.m
DART/trunk/models/cam/work/path_names_model_mod_check
DART/trunk/obs_def/obs_def_CO_Nadir_mod.f90
DART/trunk/observations/SABER/progs/
Removed Paths:
-------------
DART/trunk/models/cam/clm_ens_avg.f90
DART/trunk/models/cam/work/mkmf_clm_ens_avg
DART/trunk/models/cam/work/path_names_clm_ens_avg
Property Changed:
----------------
DART/trunk/adaptive_inflate/
DART/trunk/assim_tools/
DART/trunk/assim_tools/assim_tools_mod.f90
DART/trunk/matlab/GetClmInfo.m
DART/trunk/models/bgrid_solo/model_mod.f90
DART/trunk/models/wrf/model_mod.f90
DART/trunk/obs_kind/DEFAULT_obs_kind_mod.F90
DART/trunk/observations/MADIS/
DART/trunk/observations/NCEP/ascii_to_obs/real_obs_mod.f90
DART/trunk/utilities/
DART/trunk/utilities/closest_member_tool.f90
-------------- next part --------------
Property changes on: DART/trunk/adaptive_inflate
___________________________________________________________________
Modified: svn:mergeinfo
- /DART/branches/development/adaptive_inflate:4680-6255
/DART/branches/gitm/adaptive_inflate:5143-6215
/DART/branches/gitm_lanai/adaptive_inflate:6571-6652
/DART/branches/helen/adaptive_inflate:5995-6161
/DART/branches/inf_restart:4784-4812
+ /DART/branches/cam/adaptive_inflate:6639-7729
/DART/branches/development/adaptive_inflate:4680-6255
/DART/branches/gitm/adaptive_inflate:5143-6215
/DART/branches/gitm_lanai/adaptive_inflate:6571-6652
/DART/branches/helen/adaptive_inflate:5995-6161
/DART/branches/inf_restart:4784-4812
Modified: DART/trunk/adaptive_inflate/adaptive_inflate_mod.f90
===================================================================
--- DART/trunk/adaptive_inflate/adaptive_inflate_mod.f90 2015-03-13 22:44:50 UTC (rev 7731)
+++ DART/trunk/adaptive_inflate/adaptive_inflate_mod.f90 2015-03-13 22:53:44 UTC (rev 7732)
@@ -360,28 +360,45 @@
trim(det), trim(tadapt), trim(sadapt), trim(akind)
call error_handler(E_MSG, trim(label) // ' inflation:', msgstring, source, revision, revdate)
-! if inflation was set from a namelist, that message has already been
-! printed (further up in this routine). for values set from a restart
-! file, if the inflation flavor is 2, the values printed are the min and
-! max from the entire array. for flavors 1 and 3 there is only a single
-! value to print out.
+! if inflation was set from a namelist a message with the values has been
+! printed (further up in this routine). for values set from a restart file
+! we now know the min/max values, so print them out here.
+! if the inflation flavor is 2 there are 2 values: the min and max from the
+! state-vector-sized array. for flavors 1 and 3 there is only a single value.
+! also check for bad values (missing_r8s should not be found in the inflation files).
if (inf_flavor > 0) then
- if (mean_from_restart) then
+ ! if my task owns the mean/sd, test for any missing_r8 and error out if found
+ call get_copy_owner_index(ss_inflate_index, owner, owners_index)
+ if (owner == ens_handle%my_pe) then
+ if (any(ens_handle%vars(:, owners_index) == MISSING_R8)) then
+ call error_handler(E_ERR, 'adaptive_inflate_init', 'illegal missing values found in inflation mean file', &
+ source, revision, revdate)
+ endif
+ endif
+ call get_copy_owner_index(ss_inflate_sd_index, owner, owners_index)
+ if (owner == ens_handle%my_pe) then
+ if (any(ens_handle%vars(:, owners_index) == MISSING_R8)) then
+ call error_handler(E_ERR, 'adaptive_inflate_init', 'illegal missing values found in inflation sd file', &
+ source, revision, revdate)
+ endif
+ endif
+ ! task 0 knows the min and maxes and needs to print them for the log
+ if (my_task_id() == 0 .and. mean_from_restart) then
if (inf_flavor == 2) then
- write(msgstring, '(A, F8.3, A, F8.3)') &
+ write(msgstring, '(A, F12.3, A, F12.3)') &
'inf mean from restart file: min value: ', minmax_mean(1), ' max value: ', minmax_mean(2)
else
- write(msgstring, '(A, F8.3, A, F8.3)') &
+ write(msgstring, '(A, F12.3, A, F12.3)') &
'inf mean from restart file: value: ', minmax_mean(1)
endif
call error_handler(E_MSG, trim(label) // ' inflation:', msgstring, source, revision, revdate)
endif
- if (sd_from_restart) then
+ if (my_task_id() == 0 .and. sd_from_restart) then
if (inf_flavor == 2) then
- write(msgstring, '(A, F8.3, A, F8.3)') &
+ write(msgstring, '(A, F12.3, A, F12.3)') &
'inf stddev from restart file: min value: ', minmax_sd(1), ' max value: ', minmax_sd(2)
else
- write(msgstring, '(A, F8.3, A, F8.3)') &
+ write(msgstring, '(A, F12.3, A, F12.3)') &
'inf stddev from restart file: value: ', minmax_sd(1)
endif
call error_handler(E_MSG, trim(label) // ' inflation:', msgstring, source, revision, revdate)
Property changes on: DART/trunk/assim_tools
___________________________________________________________________
Modified: svn:mergeinfo
- /DART/branches/development/assim_tools:4680-6255
/DART/branches/gitm/assim_tools:5143-6215
/DART/branches/gitm_lanai/assim_tools:6571-6652
/DART/releases/Kodiak/assim_tools:5292-5577
+ /DART/branches/cam/assim_tools:6639-7729
/DART/branches/development/assim_tools:4680-6255
/DART/branches/gitm/assim_tools:5143-6215
/DART/branches/gitm_lanai/assim_tools:6571-6652
/DART/releases/Kodiak/assim_tools:5292-5577
Property changes on: DART/trunk/assim_tools/assim_tools_mod.f90
___________________________________________________________________
Modified: svn:mergeinfo
- /DART/branches/development/assim_tools/assim_tools_mod.f90:4680-6255
/DART/branches/gitm_lanai/assim_tools/assim_tools_mod.f90:6571-6652
/DART/branches/helen/assim_tools/assim_tools_mod.f90:5995-6161
/DART/releases/Kodiak/assim_tools/assim_tools_mod.f90:5020-5583
+ /DART/branches/cam/assim_tools/assim_tools_mod.f90:6639-7729
/DART/branches/development/assim_tools/assim_tools_mod.f90:4680-6255
/DART/branches/gitm_lanai/assim_tools/assim_tools_mod.f90:6571-6652
/DART/branches/helen/assim_tools/assim_tools_mod.f90:5995-6161
/DART/releases/Kodiak/assim_tools/assim_tools_mod.f90:5020-5583
Modified: DART/trunk/diagnostics/matlab/linked_observations.m
===================================================================
--- DART/trunk/diagnostics/matlab/linked_observations.m 2015-03-13 22:44:50 UTC (rev 7731)
+++ DART/trunk/diagnostics/matlab/linked_observations.m 2015-03-13 22:53:44 UTC (rev 7732)
@@ -81,6 +81,7 @@
figure2 = figure(2); clf(figure2); orient tall; wysiwyg
%% Create axes for time VS. QC
+% FIXME ... choose a format for datetick based on the timespan
axes4 = axes('Parent',figure2,'OuterPosition',[0 0.80 1 0.175],'FontSize',14);
set(axes4,'XAxisLocation','top')
box(axes4,'on');
@@ -96,6 +97,7 @@
ylabel(obs.colnames{obs.qcindex});
%% Create axes for observation index VS. time
+% FIXME ... choose a format for datetick based on the timespan
%axes3 = axes('Parent',figure2,'OuterPosition',[0 0.575 1 0.175]);
axes3 = axes('Parent',figure2,'OuterPosition',[0 0.600 1 0.2],'FontSize',14);
box(axes3,'on');
Modified: DART/trunk/diagnostics/matlab/plot_rank_histogram.m
===================================================================
--- DART/trunk/diagnostics/matlab/plot_rank_histogram.m 2015-03-13 22:44:50 UTC (rev 7731)
+++ DART/trunk/diagnostics/matlab/plot_rank_histogram.m 2015-03-13 22:53:44 UTC (rev 7732)
@@ -3,7 +3,7 @@
% Part of the observation-space diagnostics routines.
%
% 'obs_diag' produces a netcdf file containing the diagnostics.
-% 'obs_diag' condenses the 'obs_seq.final' information into summaries
+% 'obs_diag' condenses the 'obs_seq.final' information into summaries
% for a few specified regions - on a level-by-level basis.
%
% USAGE: plotdat = plot_rank_histogram(fname, timeindex, [,obstypestring]);
@@ -14,10 +14,10 @@
% Must match something in the 'ObservationTypes' variable.
% (ncdump -v ObservationTypes obs_diag_output.nc)
% If specified, only this observation type will be plotted.
-% If not specified, ALL observation types incluced in the
+% If not specified, ALL observation types incluced in the
% netCDF file will be plotted.
%
-% OUTPUT: two files will result for each observation type plotted. One is a
+% OUTPUT: two files will result for each observation type plotted. One is a
% postscript file containing a page for each level - all regions.
% The other file is a simple text file containing summary information
% about how many obs were assimilated, how many were available, etc.
@@ -26,7 +26,7 @@
% EXAMPLE 1 - plot the rank histogram for ALL observation types, ALL levels
%
% fname = 'obs_diag_output.nc'; % netcdf file produced by 'obs_diag'
-% timeindex = 1; % plot the histogram for the first timestep
+% timeindex = 1; % plot the histogram for the first timestep
% plotdat = plot_rank_histogram(fname, timeindex);
%
%
@@ -35,7 +35,7 @@
% observation type in the netCDF file.
%
% fname = 'obs_diag_output.nc'; % netcdf file produced by 'obs_diag'
-% timeindex = 3; % plot the histogram for the third timestep
+% timeindex = 3; % plot the histogram for the third timestep
% plotdat = plot_rank_histogram(fname, timeindex, 'RADIOSONDE_TEMPERATURE');
%
%
@@ -123,6 +123,8 @@
skip_seconds = time_to_skip(4)*3600 + time_to_skip(5)*60 + time_to_skip(6);
iskip = time_to_skip(3) + skip_seconds/86400;
+% set up a structure with all static plotting components
+
plotdat.timecenters = plotdat.timecenters + timeorigin;
plotdat.timeedges = plotdat.timeedges + timeorigin;
plotdat.Ntimes = length(plotdat.timecenters);
@@ -138,10 +140,6 @@
datestr(plotdat.timeedges(timeindex,2),21));
end
-% set up a structure with all static plotting components
-
-plotdat.linewidth = 2.0;
-
if (nvars == 0)
[plotdat.allvarnames, plotdat.allvardims] = get_varsNdims(fname);
[plotdat.varnames, plotdat.vardims] = FindTemporalVars(plotdat);
@@ -158,29 +156,39 @@
plotdat.NQC6index = get_copy_index(fname,'N_DARTqc_6');
plotdat.NQC7index = get_copy_index(fname,'N_DARTqc_7');
+figuredata = setfigure();
+
%----------------------------------------------------------------------
% Loop around (time-copy-level-region) observation types
%----------------------------------------------------------------------
for ivar = 1:plotdat.nvars
-
+
% create the variable names of interest.
% netCDF can only support variable names less than 41 chars.
-
- plotdat.myvarname = plotdat.varnames{ivar};
+
+ plotdat.myvarname = plotdat.varnames{ivar};
plotdat.guessvar = sprintf('%s_guess',plotdat.varnames{ivar});
plotdat.rhistvar = BuildFullVarname(plotdat.varnames{ivar});
+ [present, ~] = nc_var_exists(fname, plotdat.rhistvar);
+ if ( ~ present )
+ fprintf('Could not find %s in %s ... skipping\n',plotdat.rhistvar, fname)
+ continue
+ end
+
% remove any existing postscript file - will simply append each
% level as another 'page' in the .ps file.
-
- psfname = sprintf('%s_rank_hist.ps',plotdat.varnames{ivar});
- fprintf('Removing %s from the current directory.\n',psfname)
- system(sprintf('rm %s',psfname));
- % remove any existing log file -
-
+ for iregion = 1:plotdat.nregions
+ psfname{iregion} = sprintf('%s_rank_hist_region%d.ps',plotdat.varnames{ivar},iregion);
+ fprintf('Removing %s from the current directory.\n',psfname{iregion})
+ system(sprintf('rm %s',psfname{iregion}));
+ end
+
+ % remove any existing log file -
+
lgfname = sprintf('%s_rank_hist_obscount.txt',plotdat.varnames{ivar});
fprintf('Removing %s from the current directory.\n',lgfname)
system(sprintf('rm %s',lgfname));
@@ -207,36 +215,38 @@
end
plotdat.nlevels = length(plotdat.level);
- % Here is the tricky part. Singleton dimensions are auto-squeezed ...
+ % Here is the tricky part. Singleton dimensions are auto-squeezed ...
% single levels, single regions ...
- guess_raw = nc_varget(fname, plotdat.guessvar);
+ guess_raw = nc_varget(fname, plotdat.guessvar);
guess = reshape(guess_raw, plotdat.Ntimes, plotdat.ncopies, ...
plotdat.nlevels, plotdat.nregions);
- rhist_raw = nc_varget(fname, plotdat.rhistvar);
+ rhist_raw = nc_varget(fname, plotdat.rhistvar);
rhist = reshape(rhist_raw, plotdat.Ntimes, plotdat.Nrhbins, ...
plotdat.nlevels, plotdat.nregions);
-
- % Collapse the time dimension if need be.
- if ( timeindex < 0 )
- guess = sum(guess,1);
- rhist = sum(rhist,1);
+
+ % Collapse the time dimension if need be.
+ if ( timeindex < 0 )
+ guess = sum(guess,1);
+ rhist = sum(rhist,1);
plotdat.timeindex = 1;
end
% check to see if there is anything to plot
- nposs = sum(guess(plotdat.timeindex,plotdat.Npossindex,:,:));
+ nposs = sum(guess(plotdat.timeindex,plotdat.Npossindex,:,:)) - ...
+ sum(guess(plotdat.timeindex,plotdat.NQC5index ,:,:)) - ...
+ sum(guess(plotdat.timeindex,plotdat.NQC6index ,:,:));
if ( sum(nposs(:)) < 1 )
fprintf('%s no obs for %s... skipping\n', plotdat.varnames{ivar})
continue
end
-
+
for ilevel = 1:plotdat.nlevels
fprintf(logfid,'\nlevel %d %f %s\n',ilevel,plotdat.level(ilevel),plotdat.level_units);
-
+
plotdat.ges_Nqc4 = squeeze(guess(plotdat.timeindex,plotdat.NQC4index ,ilevel,:));
fprintf(logfid,'DART QC == 4, prior %d\n',sum(plotdat.ges_Nqc4(:)));
@@ -249,7 +259,8 @@
plotdat.ges_Nqc7 = squeeze(guess(plotdat.timeindex,plotdat.NQC7index ,ilevel,:));
fprintf(logfid,'DART QC == 7, prior %d\n',sum(plotdat.ges_Nqc7(:)));
- plotdat.ges_Nposs = squeeze(guess(plotdat.timeindex,plotdat.Npossindex, ilevel,:));
+ plotdat.ges_Nposs = squeeze(guess(plotdat.timeindex,plotdat.Npossindex, ilevel,:)) ...
+ - plotdat.ges_Nqc5 - plotdat.ges_Nqc6;
fprintf(logfid,'# obs poss, prior %d\n',sum(plotdat.ges_Nposs(:)));
plotdat.ges_Nused = squeeze(guess(plotdat.timeindex,plotdat.Nusedindex, ilevel,:));
@@ -257,15 +268,10 @@
% plot by region
- if (plotdat.nregions > 1)
- clf; orient tall; wysiwyg
- else
- clf; orient landscape; wysiwyg
- end
-
for iregion = 1:plotdat.nregions
+ figure(iregion); clf; orient(figuredata.orientation); wysiwyg
- plotdat.region = iregion;
+ plotdat.region = iregion;
plotdat.myregion = deblank(plotdat.region_names(iregion,:));
if ( isempty(plotdat.level_units) )
plotdat.title = plotdat.myvarname;
@@ -275,19 +281,18 @@
plotdat.level(ilevel), ...
plotdat.level_units);
end
-
+
plotdat.rank_hist = squeeze(rhist(plotdat.timeindex, :, ilevel,iregion));
-
- myplot(plotdat);
- end
- % create a postscript file
- print(gcf,'-dpsc','-append',psfname);
+ myplot(plotdat,figuredata);
- % block to go slow and look at each one ...
- disp('Pausing, hit any key to continue ...')
- pause
+ % create a postscript file
+ print(gcf,'-dpsc','-append',psfname{iregion});
+ % block to go slow and look at each one ...
+ % disp('Pausing, hit any key to continue ...')
+ % pause
+ end
end
end
@@ -295,57 +300,50 @@
% 'Helper' functions
%----------------------------------------------------------------------
-function myplot(plotdat)
+function myplot(plotdat,figdata)
- % 6 Rejected because of incoming data qc.
- % 5 Not used because of namelist control.
+ nobs_poss = plotdat.ges_Nposs(plotdat.region);
- nobs_poss = plotdat.ges_Nposs(plotdat.region) - ...
- plotdat.ges_Nqc5( plotdat.region) - ...
- plotdat.ges_Nqc6( plotdat.region);
-
% nobs_used = plotdat.ges_Nused(plotdat.region);
nobs_used = sum(plotdat.rank_hist(:));
obsstring = sprintf('%d obs possible, %d obs binned',nobs_poss, nobs_used);
% Determine some quantities for the legend
- % Plot
-
- subplot(plotdat.nregions,1,plotdat.region);
+ % Plot
+
+ ax1 = subplot('position',figdata.position);
bar(plotdat.rank_hist, 1.0);
- set(gca,'TickDir','out')
+ set(ax1,'TickDir','out','FontSize',figdata.fontsize,'YGrid','on')
axlims = axis;
axlims = [0 plotdat.Nrhbins+1 0 axlims(4)];
axis(axlims)
h = text(plotdat.Nrhbins/2, 0.9*axlims(4),obsstring);
- set(h,'FontSize',12,'FontWeight','Bold')
+ set(h,'FontSize',figdata.fontsize,'FontWeight','Bold')
set(h,'HorizontalAlignment','center')
- xlabel('Observation Rank (among ensemble members)')
-
if ( strcmp(plotdat.outlierstring, 'TRUE'))
- ylabel('count - including outlier observations')
+ ylabelstring = 'count - including outlier observations';
else
- ylabel('count')
+ ylabelstring = 'count';
end
- % more annotation ...
+ set(get(gca,'Ylabel'),'String', ylabelstring, ...
+ 'Interpreter','none','FontSize',figdata.fontsize)
+ set(get(gca,'Xlabel'),'String', {'Observation Rank (among ensemble members)', plotdat.timespan}, ...
+ 'Interpreter','none','FontSize',figdata.fontsize)
- if (plotdat.region == 1)
- title({plotdat.title, plotdat.myregion}, ...
- 'Interpreter', 'none', 'Fontsize', 14, 'FontWeight', 'bold')
- BottomAnnotation(plotdat.timespan, plotdat.fname)
- else
- title(plotdat.myregion, 'Interpreter', 'none', ...
- 'Fontsize', 14, 'FontWeight', 'bold')
- end
-
-
+ title({plotdat.myregion, plotdat.title}, ...
+ 'Interpreter', 'none', 'Fontsize', figdata.fontsize, 'FontWeight', 'bold')
+ BottomAnnotation(plotdat.fname)
-function BottomAnnotation(timespan, main)
+
+%=====================================================================
+
+
+function BottomAnnotation(main)
% annotates the directory containing the data being plotted
subplot('position',[0.10 0.01 0.8 0.04])
axis off
@@ -362,18 +360,14 @@
string1 = sprintf('data file: %s',fullname);
end
-h = text(0.5, 0.67, timespan);
+h = text(0.5, 0.5, string1);
set(h,'HorizontalAlignment','center', ...
'VerticalAlignment','middle',...
'Interpreter','none',...
'FontSize',8)
-h = text(0.5, 0.33, string1);
-set(h,'HorizontalAlignment','center', ...
- 'VerticalAlignment','middle',...
- 'Interpreter','none',...
- 'FontSize',8)
+%=====================================================================
function [y,ydims] = FindTemporalVars(x)
@@ -386,7 +380,7 @@
for i = 1:length(x.allvarnames)
indx = strfind(x.allvardims{i},'time');
- if (indx > 0)
+ if (indx > 0)
j = j + 1;
basenames{j} = ReturnBase(x.allvarnames{i});
@@ -404,7 +398,9 @@
end
+%=====================================================================
+
function s = ReturnBase(string1)
inds = strfind(string1,'_guess');
if (inds > 0 )
@@ -427,20 +423,45 @@
end
+%=====================================================================
+
function s = BuildFullVarname(varname)
% netCDF restricts the length of a variable name to
-% be 40 characters.
+% be 40 characters.
rankvar = sprintf('%s_guess_RankHist',varname);
-if (length(rankvar) > 40 )
+if (length(rankvar) > 40 )
s = rankvar(1:40);
else
s = rankvar;
end
+%=====================================================================
+
+function figdata = setfigure()
+%%
+% figure out a page layout
+% extra space at the bottom for the date/file annotation
+% extra space at the top because the titles have multiple lines
+
+orientation = 'portrait';
+fontsize = 16;
+position = [0.15 0.25 0.7 0.6];
+linewidth = 2.0;
+
+figdata = struct('expcolors', {{'k','r','g','m','b','c','y'}}, ...
+ 'expsymbols', {{'o','s','d','p','h','s','*'}}, ...
+ 'prpolines', {{'-',':'}}, 'position', position, ...
+ 'fontsize',fontsize, 'orientation',orientation, ...
+ 'linewidth',linewidth);
+
+
+%=====================================================================
+
+
function value = local_nc_varget(fname,varname)
%% If the variable exists in the file, return the contents of the variable.
% if the variable does not exist, return empty value instead of error-ing
@@ -456,7 +477,9 @@
end
+%=====================================================================
+
function value = local_nc_attget(fname,varid,varname)
%% If the (global) attribute exists, return the value.
% If it does not, do not throw a hissy-fit.
Modified: DART/trunk/diagnostics/threed_sphere/obs_diag.f90
===================================================================
--- DART/trunk/diagnostics/threed_sphere/obs_diag.f90 2015-03-13 22:44:50 UTC (rev 7731)
+++ DART/trunk/diagnostics/threed_sphere/obs_diag.f90 2015-03-13 22:53:44 UTC (rev 7732)
@@ -147,7 +147,8 @@
! 3 Evaluated only, but the posterior forward operator failed
! --- everything above this means only the prior is OK
! 4 prior forward operator failed
-! 5 was not included to be assimilated or evaluated in the namelist
+! 5 not assimilated or evaluated because of namelist specification of
+! input.nml:obs_kind_nml:[assimilate,evaluate]_these_obs_types
! 6 prior QC rejected
! 7 outlier rejected
! 8+ reserved for future use
@@ -880,7 +881,7 @@
obs_used_in_epoch(iepoch) = obs_used_in_epoch(iepoch) + 1
!--------------------------------------------------------------
- ! If it is a U wind component, all we need to do is save it.
+ ! If it is a U wind component, we need to save it.
! It will be matched up with the subsequent V component.
! At some point we have to remove the dependency that the
! U component MUST preceed the V component.
@@ -2654,6 +2655,14 @@
! If the prior mean cannot be calculated (i.e. is missing) we put it in the
! last 'bin' of the crude histogram. This is pretty much a 'z' score in the
! statistical sense.
+!
+! In Jan of 2014 I ran into a special case. We are performing a perfect model
+! experiment - perturbing a single state and then assimilating. We wanted
+! to compare against the true observation value (errvar = 0.0) and for the
+! first time step, there is no ensemble spread either. In this case the denom
+! was really zero and the calculation blew up. Since we only use this to track
+! how far apart these things are, we can restrict the distance to the worst-case
+! scenario ... the last bin.
real(r8) :: InnovZscore
real(r8), intent(in) :: obsval, prmean, prspred, errvar
@@ -2661,12 +2670,22 @@
real(r8) :: numer, denom
+InnovZscore = real(MaxSigmaBins,r8) ! worst-case ... really far apart
+
if ( qcval <= qcmaxval ) then ! QC indicates a valid obs
numer = abs(prmean - obsval)
denom = sqrt( prspred**2 + errvar )
- InnovZscore = numer / denom
-else
- InnovZscore = real(MaxSigmaBins,r8)
+
+ ! TJH FIXME ... test this before putting on the trunk!!!!
+ ! At worst, the InnovZscore can be 'MaxSigmaBins' i.e. 100
+ ! protect against dividing by a really small number
+ ! if numer/denom < 100 then go ahead and calculate
+ ! if numer < 100 * denom ...
+
+ if ( numer < InnovZscore*denom ) then
+ InnovZscore = numer / denom
+ endif
+
endif
end Function InnovZscore
Modified: DART/trunk/ensemble_manager/ensemble_manager_mod.f90
===================================================================
--- DART/trunk/ensemble_manager/ensemble_manager_mod.f90 2015-03-13 22:44:50 UTC (rev 7731)
+++ DART/trunk/ensemble_manager/ensemble_manager_mod.f90 2015-03-13 22:53:44 UTC (rev 7732)
@@ -307,7 +307,7 @@
endif
! Store this copy in the appropriate place on the appropriate process
- ! map from my_pe to physical task number all done in send and recieves only
+ ! map from my_pe to physical task number all done in send and receives only
call put_copy(map_task_to_pe(ens_handle,0), ens_handle, i, ens, ens_time)
end do
@@ -1337,7 +1337,7 @@
if (use_copy2var_send_loop .eqv. .true. ) then
! Switched loop index from receiving_pe to sending_pe
-! Aim: to make the commication scale better on Yellowstone, as num_pes >> ens_size
+! Aim: to make the communication scale better on Yellowstone, as num_pes >> ens_size
! For small numbers of tasks (32 or less) the receiving_pe loop may be faster.
! Namelist option use_copy2var_send_loop can be used to select which
! communication pattern to use
Modified: DART/trunk/filter/filter.f90
===================================================================
--- DART/trunk/filter/filter.f90 2015-03-13 22:44:50 UTC (rev 7731)
+++ DART/trunk/filter/filter.f90 2015-03-13 22:53:44 UTC (rev 7732)
@@ -1224,7 +1224,7 @@
'Reading in initial condition/restart data for all ensemble members from file(s)')
else
call error_handler(E_MSG,'filter_read_restart:', &
- 'Reading in a single ensemble and perturbing data for the other ensemble members')
+ 'Reading in a single member and perturbing data for the other ensemble members')
endif
endif
Modified: DART/trunk/location/threed_sphere/location_mod.f90
===================================================================
--- DART/trunk/location/threed_sphere/location_mod.f90 2015-03-13 22:44:50 UTC (rev 7731)
+++ DART/trunk/location/threed_sphere/location_mod.f90 2015-03-13 22:53:44 UTC (rev 7732)
@@ -1836,7 +1836,6 @@
if (howmuch > 1) then
! DEBUG
write(errstring,"(A,I8,A,36(I8,1X))") ' obs_box(',i,') =', gc%obs_box(1:min(i,36)) ! (nobs)
- !write(errstring,*) ' obs_box(',i,') =', gc%obs_box ! (nobs)
call error_handler(E_MSG, 'locations_mod', errstring)
else if(howmuch > 0) then
write(errstring,*) ' obs_box(',i,') =', gc%obs_box(1:min(i,sample+1))
Property changes on: DART/trunk/matlab/GetClmInfo.m
___________________________________________________________________
Modified: svn:mergeinfo
- /DART/branches/development/matlab/GetClmInfo.m:4680-5731
/DART/branches/gitm_lanai/matlab/GetClmInfo.m:6571-6652
/DART/branches/helen/matlab/GetClmInfo.m:5995-6161
+ /DART/branches/cam/matlab/GetClmInfo.m:6639-7729
/DART/branches/development/matlab/GetClmInfo.m:4680-5731
/DART/branches/gitm_lanai/matlab/GetClmInfo.m:6571-6652
/DART/branches/helen/matlab/GetClmInfo.m:5995-6161
Property changes on: DART/trunk/models/bgrid_solo/model_mod.f90
___________________________________________________________________
Modified: svn:mergeinfo
- /DART/branches/development/models/bgrid_solo/model_mod.f90:4680-6255
/DART/branches/gitm_lanai/models/bgrid_solo/model_mod.f90:6571-6652
/DART/branches/helen/models/bgrid_solo/model_mod.f90:5995-6161
+ /DART/branches/cam/models/bgrid_solo/model_mod.f90:6639-7729
/DART/branches/development/models/bgrid_solo/model_mod.f90:4680-6255
/DART/branches/gitm_lanai/models/bgrid_solo/model_mod.f90:6571-6652
/DART/branches/helen/models/bgrid_solo/model_mod.f90:5995-6161
Deleted: DART/trunk/models/cam/clm_ens_avg.f90
===================================================================
--- DART/trunk/models/cam/clm_ens_avg.f90 2015-03-13 22:44:50 UTC (rev 7731)
+++ DART/trunk/models/cam/clm_ens_avg.f90 2015-03-13 22:53:44 UTC (rev 7732)
@@ -1,722 +0,0 @@
-! DART software - Copyright 2004 - 2013 UCAR. This open source software is
-! provided by UCAR, "as is", without charge, subject to all terms of use at
-! http://www.image.ucar.edu/DAReS/DART/DART_download
-!
-! $Id$
-
-program clm_ens_avg
-
-!----------------------------------------------------------------------
-! purpose: Average an ensemble of CLM initial/restart NetCDF files
-! for use with the ensemble mean CAM files to be used as ICs for forecasts.
-!
-! method: Read new ensemble file of CLM, 'initial', snow and water fields.
-! Shift snow layers so that top layers line up. Average them. Shift them back down.
-! Preserve special values in some columns.
-! (Other fields will be averaged using NCO ncea). Write the averaged snow and water
-! fields into the pre-existing (from NCO ncea) ensemble averaged file.
-! Not generically coded to handle arbitrary lists of fields in arbitrary order.
-! Little abstraction/object-oriented programming.
-! But does distinguish CLM3 from 3.6 by reading :source attribute from initial files.
-! Updated for CAM3.6; mss_XXX,flx_XXX, etc. fields
-! Note that at least 1 dimension name changed (levsoi -> levgrnd),
-! and that the the dimension names in 3.5 which were not used (levsoi, ...?),
-! hence not passed to the new file, may be different in 3.6.
-! Flx and mss fields don't have special values anywhere, but snw_rds should be
-! averaged only over members with snow in the relevant layer(s).
-! ? flx_ are dimensioned (levsno+1,column), the extra being the first ground layer.
-! NetCDF will allow just (levsno,column) of them to be read,
-! and write those back out again. (The ground layer will be averaged by ncea.)
-!
-! author: Kevin Raeder 4/22/09
-! upgrade to CLM 3.6 Kevin Raeder 8/20/09
-! based on model_mod:read_cam_init
-!
-!----------------------------------------------------------------------
-
-use types_mod, only : r8
-use utilities_mod, only : get_unit, file_exist, nc_check, &
- initialize_utilities, finalize_utilities
-use netcdf
-use typeSizes
-
-! eventually use model_mod, only : read_clm_init
-! use assim_model_mod, only : &
-! will want eventually for CLM; write it now? get_model_size , &
-
-implicit none
-
-! version controlled file description for error handling, do not edit
-character(len=256), parameter :: source = &
- "$URL$"
-character(len=32 ), parameter :: revision = "$Revision$"
-character(len=128), parameter :: revdate = "$Date$"
-
-character (len = 128) :: file_in = 'snow_water_ens.nc', file_out = 'clm_ens_avg.nc' &
- ,glob_att_name = 'source', glob_att_val, model_ver
-
-! lo_top is the index of the top layer of snow when all the layers on "on the ground"
-! hi_bot is the index of the bottom layer of snow when all the layers have been "raised up"
-! to allign all the top layers of snow into index 1.
-! zero_ are the complements to lo_top and hi_bot, used for zeroing the vacated layers.
-integer :: file_unit, ncfldid, ncfileid, e,c,k, glob_att_len
-integer :: num_col,num_ens, num_levsno, t_soisno_num, lo_top, hi_bot, zero_bot, zero_top
-
-logical :: do_output = .false. &
- ,lcam4 = .true.
-
-! Temporary allocatable storage to read in a native format for cam state.
-
-! 2D fields (col,lvl) with ens dimension
-real(r8), allocatable, dimension(:,:,:) :: &
- dzsnoE,zsnoE,zisnoE,t_soisnoE,h2osoi_iceE,h2osoi_liqE &
- ,snw_rdsE,qflx_snofrz_lyrE &
- ,mss_bcphoE,mss_bcphiE,mss_ocphoE,mss_ocphiE &
- ,mss_dst1E,mss_dst2E,mss_dst3E,mss_dst4E &
- ,flx_absdvE,flx_absdnE,flx_absivE,flx_absinE
-! 2D fields
-real(r8), allocatable, dimension(:,:) :: &
- dzsno,zsno,zisno,t_soisno,h2osoi_ice,h2osoi_liq &
- ,snw_rds,qflx_snofrz_lyr &
- ,mss_bcpho,mss_bcphi,mss_ocpho,mss_ocphi &
- ,mss_dst1,mss_dst2,mss_dst3,mss_dst4 &
- ,flx_absdv,flx_absdn,flx_absiv,flx_absin
-! 1D fields (col) with ens dimension
-real(r8), allocatable, dimension(:,:) :: snowdpE,h2osnoE
-! 1D fields (col)
-real(r8), allocatable, dimension(:) :: snowdp,h2osno
-! 1D field with ens dimension for storing # of snow layers
-integer, allocatable, dimension(:,:) :: snlsnoE
-! 1D field
-integer, allocatable, dimension(:) :: snlsno
-real(r8) :: snowdp_sum, ICE_sum, LIQ_sum, DZ_sum, SNOWmax, TSOI_sum, SNODIFF, &
- snw_rds_sum,qflx_snofrz_lyr_sum, &
- mss_bcpho_sum,mss_bcphi_sum,mss_ocpho_sum,mss_ocphi_sum, &
- mss_dst1_sum,mss_dst2_sum,mss_dst3_sum,mss_dst4_sum, &
- flx_absdv_sum,flx_absdn_sum,flx_absiv_sum,flx_absin_sum
-
-print*,'==========================================================================='
-print*,'clm_ens_avg.f90:'
-
-call initialize_utilities('clm_ens_avg')
-
-! Static init assim model calls static_init_model which now (merge/MPI) calls read_cam_init)
-! call static_init_assim_model()
-
-! Read in the array dimensions from the ens file
-call nc_check(nf90_open(path = trim(file_in), mode = nf90_write, ncid = ncfileid), &
- 'clm_ens_avg', 'opening '//trim(file_in))
-
-! Get model version from the CLM initial file global attribute :source
-! :source = "Community Land Model: CLM3" ; or CLM3.6
-call nc_check(nf90_inquire_attribute(ncfileid, nf90_global, trim(glob_att_name), len = glob_att_len), &
- 'clm_ens_avg','inquire; get CAM version')
-call nc_check(nf90_get_att(ncfileid, nf90_global, glob_att_name, glob_att_val), &
- 'clm_ens_avg','get_att; get CAM version')
-! assim_model example:
-! call nc_check(nf90_put_att(ncFileID%ncid, NF90_GLOBAL, "assim_model_source", source ), &
-
-model_ver = glob_att_val(22:glob_att_len)
-print*,'glob_att_name, Model version are ',trim(glob_att_name),' ', trim(model_ver)
-
-call read_clm_init_size(ncfileid, num_ens, num_levsno, num_col)
-
-! ens fields
-allocate (dzsnoE (num_levsno,num_col,num_ens), zsnoE (num_levsno,num_col,num_ens) &
- ,zisnoE (num_levsno,num_col,num_ens), t_soisnoE (num_levsno,num_col,num_ens) &
- ,h2osoi_iceE(num_levsno,num_col,num_ens), h2osoi_liqE (num_levsno,num_col,num_ens) &
- ,snw_rdsE (num_levsno,num_col,num_ens), qflx_snofrz_lyrE(num_levsno,num_col,num_ens) &
- ,mss_bcphoE (num_levsno,num_col,num_ens), mss_bcphiE (num_levsno,num_col,num_ens) &
- ,mss_ocphoE (num_levsno,num_col,num_ens), mss_ocphiE (num_levsno,num_col,num_ens) &
- ,mss_dst1E (num_levsno,num_col,num_ens), mss_dst2E (num_levsno,num_col,num_ens) &
- ,mss_dst3E (num_levsno,num_col,num_ens), mss_dst4E (num_levsno,num_col,num_ens) &
- ,flx_absdvE (num_levsno,num_col,num_ens), flx_absdnE (num_levsno,num_col,num_ens) &
- ,flx_absivE (num_levsno,num_col,num_ens), flx_absinE (num_levsno,num_col,num_ens) )
-
-allocate (snowdpE(num_col,num_ens), h2osnoE(num_col,num_ens), snlsnoE(num_col,num_ens))
-! ens avg fields
-allocate (dzsno (num_levsno,num_col), zsno (num_levsno,num_col) &
- ,zisno (num_levsno,num_col), t_soisno (num_levsno,num_col) &
- ,h2osoi_ice(num_levsno,num_col), h2osoi_liq (num_levsno,num_col) &
- ,snw_rds (num_levsno,num_col), qflx_snofrz_lyr (num_levsno,num_col) &
- ,mss_bcpho (num_levsno,num_col), mss_bcphi (num_levsno,num_col) &
- ,mss_ocpho (num_levsno,num_col), mss_ocphi (num_levsno,num_col) &
- ,mss_dst1 (num_levsno,num_col), mss_dst2 (num_levsno,num_col) &
- ,mss_dst3 (num_levsno,num_col), mss_dst4 (num_levsno,num_col) &
- ,flx_absdv (num_levsno,num_col), flx_absdn (num_levsno,num_col) &
- ,flx_absiv (num_levsno,num_col), flx_absin (num_levsno,num_col) )
-
-allocate (snowdp( num_col), h2osno( num_col), snlsno( num_col))
-
-! Initialize output fields to 0.
-dzsno = 0.0_r8
-zsno = 0.0_r8
-zisno = 0.0_r8
-t_soisno = 0.0_r8
-h2osoi_ice = 0.0_r8
-h2osoi_liq = 0.0_r8
-snowdp = 0.0_r8
-h2osno = 0.0_r8
-snlsno = 0
-snw_rds = 0.0_r8
-qflx_snofrz_lyr = 0.0_r8
-mss_bcpho = 0.0_r8
-mss_bcphi = 0.0_r8
-mss_ocpho = 0.0_r8
-mss_ocphi = 0.0_r8
-mss_dst1 = 0.0_r8
-mss_dst2 = 0.0_r8
-mss_dst3 = 0.0_r8
-mss_dst4 = 0.0_r8
-flx_absdv = 0.0_r8
-flx_absdn = 0.0_r8
-flx_absiv = 0.0_r8
-flx_absin = 0.0_r8
-
-
-! Read the ensembles of fields from the file of concatenated clminput_#.nc files
-! (member,levels,field)
-! [?generalize program with: call read_clm_init(ens(1,1,fld),field_name(fld), num_ens, levsno, fld)]
-
-! Multi-level fields
-call read_real(ncfileid, 'DZSNO', num_levsno,num_col,num_ens, dzsnoE )
-call read_real(ncfileid, 'T_SOISNO', num_levsno,num_col,num_ens, t_soisnoE )
-call read_real(ncfileid, 'H2OSOI_ICE', num_levsno,num_col,num_ens, h2osoi_iceE )
-call read_real(ncfileid, 'H2OSOI_LIQ', num_levsno,num_col,num_ens, h2osoi_liqE )
-if (trim(model_ver) == 'CLM3.6') then
- call read_real(ncfileid, 'snw_rds', num_levsno,num_col,num_ens, snw_rdsE )
- call read_real(ncfileid, 'qflx_snofrz_lyr', num_levsno,num_col,num_ens, qflx_snofrz_lyrE )
- call read_real(ncfileid, 'mss_bcpho', num_levsno,num_col,num_ens, mss_bcphoE )
- call read_real(ncfileid, 'mss_bcphi', num_levsno,num_col,num_ens, mss_bcphiE )
- call read_real(ncfileid, 'mss_ocpho', num_levsno,num_col,num_ens, mss_ocphoE )
- call read_real(ncfileid, 'mss_ocphi', num_levsno,num_col,num_ens, mss_ocphiE )
- call read_real(ncfileid, 'mss_dst1', num_levsno,num_col,num_ens, mss_dst1E )
- call read_real(ncfileid, 'mss_dst2', num_levsno,num_col,num_ens, mss_dst2E )
- call read_real(ncfileid, 'mss_dst3', num_levsno,num_col,num_ens, mss_dst3E )
- call read_real(ncfileid, 'mss_dst4', num_levsno,num_col,num_ens, mss_dst4E )
- call read_real(ncfileid, 'flx_absdv', num_levsno,num_col,num_ens, flx_absdvE )
- call read_real(ncfileid, 'flx_absdn', num_levsno,num_col,num_ens, flx_absdnE )
- call read_real(ncfileid, 'flx_absiv', num_levsno,num_col,num_ens, flx_absivE )
- call read_real(ncfileid, 'flx_absin', num_levsno,num_col,num_ens, flx_absinE )
-endif
-
-! Single level fields
-call read_2Dreal(ncfileid, 'SNOWDP', num_col,num_ens, snowdpE)
-call read_2Dreal(ncfileid, 'H2OSNO', num_col,num_ens, h2osnoE)
-
-! Integer single level fields
-call read_2Dint(ncfileid, 'SNLSNO', num_col,num_ens, snlsnoE)
-
-call nc_check(nf90_close(ncfileid), 'clm_ens_avg', 'close input file')
-
-! find maximum number of snow layers amoung the ensemble members for each column.
-snlsno = minval(snlsnoE,2)
-
-! Shift snowy columns up to top levels, so that all top layers are in index 1.
-do c=1,num_col
- ! average the snow depth for additional filtering of calculations
-! Led to wrong values? snlsno(c) = minval(snlsnoE(c,:))
-! if (c == 1) write(*,'(A,5I4)') 'snlsno, snlsnoE = ',snlsno(c),(snlsnoE (c,e),e=1,num_ens)
-! write(*,'(A,(5I5))') 'col,snlsno, snlsnoE = ',c,snlsno(c),(snlsnoE (c,e),e=1,num_ens)
- SNOWmax = maxval(snowdpE(c,:))
- if (SNOWmax > 0.0_r8) then
- ! only need to test first member because these spvals come from a mask, not physical variables
- if (h2osoi_liqE(1,c,1) < 1.e+35) then
- snowdp_sum = 0.0_r8
- do e=1,num_ens
- snowdp_sum = snowdp_sum + snowdpE(c,e)
-! if (c == 1) print*,'snowdpE, snowdp_sum = ',snowdpE(c,e), snowdp_sum
-! lo_top, hi_bot, ... refer to positions in array as found on the NetCDF file.
-! Those positions are counted in the opposite direction from how CLM numbers them.
- lo_top = num_levsno + snlsnoE(c,e) + 1
- if (lo_top > 1 .and. snlsnoE(c,e) < 0) then
- hi_bot = -snlsnoE(c,e)
-! write(*,'(A,5I4)') 'snlsnoE,lo_top,hi_bot = ',snlsnoE (c,e),lo_top,hi_bot
-! write(*,'(A,5F8.5)') 'dzsnoE before shift up = ',(dzsnoE(k,c,e),k=1,num_levsno)
- dzsnoE (1:hi_bot,c,e) = dzsnoE (lo_top:num_levsno,c,e)
- h2osoi_liqE (1:hi_bot,c,e) = h2osoi_liqE (lo_top:num_levsno,c,e)
- h2osoi_iceE (1:hi_bot,c,e) = h2osoi_iceE (lo_top:num_levsno,c,e)
- t_soisnoE (1:hi_bot,c,e) = t_soisnoE (lo_top:num_levsno,c,e)
- if (trim(model_ver) == 'CLM3.6') then
- snw_rdsE (1:hi_bot,c,e) = snw_rdsE (lo_top:num_levsno,c,e)
- qflx_snofrz_lyrE(1:hi_bot,c,e) = qflx_snofrz_lyrE(lo_top:num_levsno,c,e)
- mss_bcphoE (1:hi_bot,c,e) = mss_bcphoE (lo_top:num_levsno,c,e)
- mss_bcphiE (1:hi_bot,c,e) = mss_bcphiE (lo_top:num_levsno,c,e)
- mss_ocphoE (1:hi_bot,c,e) = mss_ocphoE (lo_top:num_levsno,c,e)
- mss_ocphiE (1:hi_bot,c,e) = mss_ocphiE (lo_top:num_levsno,c,e)
- mss_dst1E (1:hi_bot,c,e) = mss_dst1E (lo_top:num_levsno,c,e)
- mss_dst2E (1:hi_bot,c,e) = mss_dst2E (lo_top:num_levsno,c,e)
- mss_dst3E (1:hi_bot,c,e) = mss_dst3E (lo_top:num_levsno,c,e)
- mss_dst4E (1:hi_bot,c,e) = mss_dst4E (lo_top:num_levsno,c,e)
@@ Diff output truncated at 40000 characters. @@
More information about the Dart-dev
mailing list