[Dart-dev] [7985] DART/trunk/models/forced_lorenz_96: updated exercises and docs for DART demos at the STATMOS workshop
nancy at ucar.edu
nancy at ucar.edu
Fri May 15 13:47:25 MDT 2015
Revision: 7985
Author: nancy
Date: 2015-05-15 13:47:24 -0600 (Fri, 15 May 2015)
Log Message:
-----------
updated exercises and docs for DART demos at the STATMOS workshop
Modified Paths:
--------------
DART/trunk/models/forced_lorenz_96/model_mod.html
Added Paths:
-----------
DART/trunk/models/forced_lorenz_96/matlab/make_sequence_input.m
DART/trunk/models/forced_lorenz_96/matlab/set_extended_state.m
DART/trunk/models/forced_lorenz_96/shell_scripts/
DART/trunk/models/forced_lorenz_96/shell_scripts/README
DART/trunk/models/forced_lorenz_96/shell_scripts/advance_model.csh
-------------- next part --------------
Added: DART/trunk/models/forced_lorenz_96/matlab/make_sequence_input.m
===================================================================
--- DART/trunk/models/forced_lorenz_96/matlab/make_sequence_input.m (rev 0)
+++ DART/trunk/models/forced_lorenz_96/matlab/make_sequence_input.m 2015-05-15 19:47:24 UTC (rev 7985)
@@ -0,0 +1,474 @@
+function make_sequence_input(varargin)
+%% make_sequence_input creates input for create_obs_sequence
+%
+% Just as an example, this will create the text input required
+% to create N evenly-spaced observation prototypes. There are no
+% observation values in these. The output file name is 'input_to_create_obs_sequence.txt'
+%
+% num_obs is the number of observations on the unit circle
+% obs_err_var is the observation error variance
+% loc_offset is a static offset for each observation location.
+% loc_offset = 0.0 will make the initial observation
+% coincident with the model grid.
+%
+% Three Examples:
+% make_sequence_input('num_obs',40)
+% make_sequence_input('num_obs',40,'obs_err_var',1.0)
+% make_sequence_input('num_obs',40,'obs_err_var',1.0,'loc_offset',0.02)
+%
+% Example on how to use the output file ...
+%
+% [unixprompt] ./create_obs_sequence < input_to_create_obs_sequence.txt
+%
+% The output file from create_obs_sequence will be called 'set_def.out'
+% and can be used as input to create_fixed_network_sequence.
+
+%% 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
+%
+% DART $Id$
+
+%%------------------------------------------------------------------------------
+% This bit just parses the input and uses default values if need be.
+
+default_num_obs = 40;
+default_obs_err_var = 1.0;
+default_loc_offset = 0.0;
+default_graphic = 0;
+p = inputParser;
+
+addParamValue(p,'num_obs', default_num_obs, @isnumeric);
+addParamValue(p,'obs_err_var',default_obs_err_var, at isnumeric);
+addParamValue(p,'loc_offset', default_loc_offset, @isnumeric);
+addParamValue(p,'graphic', default_graphic, @isnumeric);
+parse(p,varargin{:});
+
+if ~isempty(fieldnames(p.Unmatched))
+ disp('Extra inputs:')
+ disp(p.Unmatched)
+end
+
+Nobs = p.Results.num_obs;
+observation_error_variance = p.Results.obs_err_var;
+loc_offset = p.Results.loc_offset;
+graphic = p.Results.graphic;
+
+locations = (0:Nobs-1)/Nobs + loc_offset; % set the location to be uniform.
+newlocs = mod(locations,1.0);
+filename = 'input_to_create_obs_sequence.txt';
+
+%%------------------------------------------------------------------------------
+% Some stuff gets written to the file once at the beginning.
+
+fid = fopen(filename,'w');
+
+fprintf(fid,'%d \n',Nobs); % upper bound on number of obs in sequence
+fprintf(fid,'0 \n'); % number of copies of data
+fprintf(fid,'0 \n'); % number of quality control values per field
+fprintf(fid,'0 \n'); % 0 is a flag to continue
+
+%%------------------------------------------------------------------------------
+% Some stuff gets written for every observation
+
+for i = 1:Nobs
+
+ fprintf(fid,'%d \n',1); % observing RAW_STATE_VARIABLE
+ fprintf(fid,'%f \n',newlocs(i)); % the unit sphere location
+ fprintf(fid,'0 0 \n'); % days and seconds of observation time
+ fprintf(fid,'%f \n',observation_error_variance);
+
+ if i < Nobs
+ fprintf(fid,'0 \n'); % 0 is a flag to continue
+ end
+
+end
+
+%%------------------------------------------------------------------------------
+% Some stuff gets written to the file at the end.
+
+fprintf(fid,'set_def.out\n'); % name of observation sequence file prototype
+fclose(fid);
+fprintf('Created %s\n',filename)
+
+%%------------------------------------------------------------------------------
+% Optional visualizer of the observation locations.
+
+if graphic ~= 0
+
+ model_size = 40;
+ gridlocs = 2*pi*(0:model_size)/model_size;
+ x = ones(size(gridlocs));
+ h1 = polar(gridlocs,x,'o');
+ set(h1,'MarkerSize',10,'LineWidth',2)
+ relabel(gca)
+
+ hold on;
+ y = ones(size(newlocs));
+ h2 = polar(2*pi*newlocs,y,'r*');
+ set(h2,'MarkerSize',20,'LineWidth',2)
+
+ h = title('forced_lorenz_96');
+ set(h,'Interpreter','none','FontSize',20,'FontWeight','bold');
+ h = xlabel({'unit circle',filename});
+ set(h,'Interpreter','none');
+
+ h = legend('grid locations','observation locations');
+ set(h,'Location','Best','FontSize',14);
+ hold off;
+
+end
+
+end % main function
+
+
+%% Helper function below ...
+
+function relabel(handle)
+
+% remove dotted radii and tickmarks and lables (actuall, all text)
+
+t = findall(handle,'type','text'); %get handles for all text in polar plot
+% set(t,'FontSize',20,'FontWeight','bold')
+set(t,'FontSize',20)
+
+for i = 1:numel(t)
+ switch get(t(i),'String')
+ case {'0'}
+ set(t(i),'String',' 0.0', ...
+ 'HorizontalAlignment','left')
+ case {'90'}
+ set(t(i),'String','0.25')
+ case {'180'}
+ set(t(i),'String','0.50', ...
+ 'HorizontalAlignment','right')
+ case {'270'}
+ set(t(i),'String','0.75')
+ otherwise
+ set(t(i),'String',' ')
+ end
+end
+
+end % function relabel
+
+% <next few lines under version control, do not edit>
+% $URL$
+% $Revision$
+% $Date$
+
+function make_sequence_input(varargin)
+%% make_sequence_input creates input for create_obs_sequence
+%
+% Just as an example, this will create the text input required
+% to create N evenly-spaced observation prototypes. There are no
+% observation values in these. The output file name is 'input_to_create_obs_sequence.txt'
+%
+% num_obs is the number of observations on the unit circle
+% obs_err_var is the observation error variance
+% loc_offset is a static offset for each observation location.
+% loc_offset = 0.0 will make the initial observation
+% coincident with the model grid.
+%
+% Three Examples:
+% make_sequence_input('num_obs',40)
+% make_sequence_input('num_obs',40,'obs_err_var',1.0)
+% make_sequence_input('num_obs',40,'obs_err_var',1.0,'loc_offset',0.02)
+%
+% Example on how to use the output file ...
+%
+% [unixprompt] ./create_obs_sequence < input_to_create_obs_sequence.txt
+%
+% The output file from create_obs_sequence will be called 'set_def.out'
+% and can be used as input to create_fixed_network_sequence.
+
+%% 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
+%
+% DART $Id$
+
+%%------------------------------------------------------------------------------
+% This bit just parses the input and uses default values if need be.
+
+default_num_obs = 40;
+default_obs_err_var = 1.0;
+default_loc_offset = 0.0;
+default_graphic = 0;
+p = inputParser;
+
+addParamValue(p,'num_obs', default_num_obs, @isnumeric);
+addParamValue(p,'obs_err_var',default_obs_err_var, at isnumeric);
+addParamValue(p,'loc_offset', default_loc_offset, @isnumeric);
+addParamValue(p,'graphic', default_graphic, @isnumeric);
+parse(p,varargin{:});
+
+if ~isempty(fieldnames(p.Unmatched))
+ disp('Extra inputs:')
+ disp(p.Unmatched)
+end
+
+Nobs = p.Results.num_obs;
+observation_error_variance = p.Results.obs_err_var;
+loc_offset = p.Results.loc_offset;
+graphic = p.Results.graphic;
+
+locations = (0:Nobs-1)/Nobs + loc_offset; % set the location to be uniform.
+newlocs = mod(locations,1.0);
+filename = 'input_to_create_obs_sequence.txt';
+
+%%------------------------------------------------------------------------------
+% Some stuff gets written to the file once at the beginning.
+
+fid = fopen(filename,'w');
+
+fprintf(fid,'%d \n',Nobs); % upper bound on number of obs in sequence
+fprintf(fid,'0 \n'); % number of copies of data
+fprintf(fid,'0 \n'); % number of quality control values per field
+fprintf(fid,'0 \n'); % 0 is a flag to continue
+
+%%------------------------------------------------------------------------------
+% Some stuff gets written for every observation
+
+for i = 1:Nobs
+
+ fprintf(fid,'%d \n',1); % observing RAW_STATE_VARIABLE
+ fprintf(fid,'%f \n',newlocs(i)); % the unit sphere location
+ fprintf(fid,'0 0 \n'); % days and seconds of observation time
+ fprintf(fid,'%f \n',observation_error_variance);
+
+ if i < Nobs
+ fprintf(fid,'0 \n'); % 0 is a flag to continue
+ end
+
+end
+
+%%------------------------------------------------------------------------------
+% Some stuff gets written to the file at the end.
+
+fprintf(fid,'set_def.out\n'); % name of observation sequence file prototype
+fclose(fid);
+fprintf('Created %s\n',filename)
+
+%%------------------------------------------------------------------------------
+% Optional visualizer of the observation locations.
+
+if graphic ~= 0
+
+ model_size = 40;
+ gridlocs = 2*pi*(0:model_size)/model_size;
+ x = ones(size(gridlocs));
+ h1 = polar(gridlocs,x,'o');
+ set(h1,'MarkerSize',10,'LineWidth',2)
+ relabel(gca)
+
+ hold on;
+ y = ones(size(newlocs));
+ h2 = polar(2*pi*newlocs,y,'r*');
+ set(h2,'MarkerSize',20,'LineWidth',2)
+
+ h = title('forced_lorenz_96');
+ set(h,'Interpreter','none','FontSize',20,'FontWeight','bold');
+ h = xlabel({'unit circle',filename});
+ set(h,'Interpreter','none');
+
+ h = legend('grid locations','observation locations');
+ set(h,'Location','Best','FontSize',14);
+ hold off;
+
+end
+
+end % main function
+
+
+%% Helper function below ...
+
+function relabel(handle)
+
+% remove dotted radii and tickmarks and lables (actuall, all text)
+
+t = findall(handle,'type','text'); %get handles for all text in polar plot
+% set(t,'FontSize',20,'FontWeight','bold')
+set(t,'FontSize',20)
+
+for i = 1:numel(t)
+ switch get(t(i),'String')
+ case {'0'}
+ set(t(i),'String',' 0.0', ...
+ 'HorizontalAlignment','left')
+ case {'90'}
+ set(t(i),'String','0.25')
+ case {'180'}
+ set(t(i),'String','0.50', ...
+ 'HorizontalAlignment','right')
+ case {'270'}
+ set(t(i),'String','0.75')
+ otherwise
+ set(t(i),'String',' ')
+ end
+end
+
+end % function relabel
+
+% <next few lines under version control, do not edit>
+% $URL$
+% $Revision$
+% $Date$
+
+function make_sequence_input(varargin)
+%% make_sequence_input creates input for create_obs_sequence
+%
+% Just as an example, this will create the text input required
+% to create N evenly-spaced observation prototypes. There are no
+% observation values in these. The output file name is 'input_to_create_obs_sequence.txt'
+%
+% num_obs is the number of observations on the unit circle
+% obs_err_var is the observation error variance
+% loc_offset is a static offset for each observation location.
+% loc_offset = 0.0 will make the initial observation
+% coincident with the model grid.
+%
+% Three Examples:
+% make_sequence_input('num_obs',40)
+% make_sequence_input('num_obs',40,'obs_err_var',1.0)
+% make_sequence_input('num_obs',40,'obs_err_var',1.0,'loc_offset',0.02)
+%
+% Example on how to use the output file ...
+%
+% [unixprompt] ./create_obs_sequence < input_to_create_obs_sequence.txt
+%
+% The output file from create_obs_sequence will be called 'set_def.out'
+% and can be used as input to create_fixed_network_sequence.
+
+%% 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
+%
+% DART $Id$
+
+%%------------------------------------------------------------------------------
+% This bit just parses the input and uses default values if need be.
+
+default_num_obs = 40;
+default_obs_err_var = 1.0;
+default_loc_offset = 0.0;
+default_graphic = 0;
+p = inputParser;
+
+addParamValue(p,'num_obs', default_num_obs, @isnumeric);
+addParamValue(p,'obs_err_var',default_obs_err_var, at isnumeric);
+addParamValue(p,'loc_offset', default_loc_offset, @isnumeric);
+addParamValue(p,'graphic', default_graphic, @isnumeric);
+parse(p,varargin{:});
+
+if ~isempty(fieldnames(p.Unmatched))
+ disp('Extra inputs:')
+ disp(p.Unmatched)
+end
+
+Nobs = p.Results.num_obs;
+observation_error_variance = p.Results.obs_err_var;
+loc_offset = p.Results.loc_offset;
+graphic = p.Results.graphic;
+
+locations = (0:Nobs-1)/Nobs + loc_offset; % set the location to be uniform.
+newlocs = mod(locations,1.0);
+filename = 'input_to_create_obs_sequence.txt';
+
+%%------------------------------------------------------------------------------
+% Some stuff gets written to the file once at the beginning.
+
+fid = fopen(filename,'w');
+
+fprintf(fid,'%d \n',Nobs); % upper bound on number of obs in sequence
+fprintf(fid,'0 \n'); % number of copies of data
+fprintf(fid,'0 \n'); % number of quality control values per field
+fprintf(fid,'0 \n'); % 0 is a flag to continue
+
+%%------------------------------------------------------------------------------
+% Some stuff gets written for every observation
+
+for i = 1:Nobs
+
+ fprintf(fid,'%d \n',1); % observing RAW_STATE_VARIABLE
+ fprintf(fid,'%f \n',newlocs(i)); % the unit sphere location
+ fprintf(fid,'0 0 \n'); % days and seconds of observation time
+ fprintf(fid,'%f \n',observation_error_variance);
+
+ if i < Nobs
+ fprintf(fid,'0 \n'); % 0 is a flag to continue
+ end
+
+end
+
+%%------------------------------------------------------------------------------
+% Some stuff gets written to the file at the end.
+
+fprintf(fid,'set_def.out\n'); % name of observation sequence file prototype
+fclose(fid);
+fprintf('Created %s\n',filename)
+
+%%------------------------------------------------------------------------------
+% Optional visualizer of the observation locations.
+
+if graphic ~= 0
+
+ model_size = 40;
+ gridlocs = 2*pi*(0:model_size)/model_size;
+ x = ones(size(gridlocs));
+ h1 = polar(gridlocs,x,'o');
+ set(h1,'MarkerSize',10,'LineWidth',2)
+ relabel(gca)
+
+ hold on;
+ y = ones(size(newlocs));
+ h2 = polar(2*pi*newlocs,y,'r*');
+ set(h2,'MarkerSize',20,'LineWidth',2)
+
+ h = title('forced_lorenz_96');
+ set(h,'Interpreter','none','FontSize',20,'FontWeight','bold');
+ h = xlabel({'unit circle',filename});
+ set(h,'Interpreter','none');
+
+ h = legend('grid locations','observation locations');
+ set(h,'Location','Best','FontSize',14);
+ hold off;
+
+end
+
+end % main function
+
+
+%% Helper function below ...
+
+function relabel(handle)
+
+% remove dotted radii and tickmarks and lables (actuall, all text)
+
+t = findall(handle,'type','text'); %get handles for all text in polar plot
+% set(t,'FontSize',20,'FontWeight','bold')
+set(t,'FontSize',20)
+
+for i = 1:numel(t)
+ switch get(t(i),'String')
+ case {'0'}
+ set(t(i),'String',' 0.0', ...
+ 'HorizontalAlignment','left')
+ case {'90'}
+ set(t(i),'String','0.25')
+ case {'180'}
+ set(t(i),'String','0.50', ...
+ 'HorizontalAlignment','right')
+ case {'270'}
+ set(t(i),'String','0.75')
+ otherwise
+ set(t(i),'String',' ')
+ end
+end
+
+end % function relabel
+
+% <next few lines under version control, do not edit>
+% $URL$
+% $Revision$
+% $Date$
+
Property changes on: DART/trunk/models/forced_lorenz_96/matlab/make_sequence_input.m
___________________________________________________________________
Added: svn:mime-type
+ text/plain
Added: svn:keywords
+ Date Rev Author HeadURL Id
Added: svn:eol-style
+ native
Added: DART/trunk/models/forced_lorenz_96/matlab/set_extended_state.m
===================================================================
--- DART/trunk/models/forced_lorenz_96/matlab/set_extended_state.m (rev 0)
+++ DART/trunk/models/forced_lorenz_96/matlab/set_extended_state.m 2015-05-15 19:47:24 UTC (rev 7985)
@@ -0,0 +1,372 @@
+function set_extended_state(fname, varargin)
+%% set_extended_state modifies the extended state portion of the perfect_ics file.
+%
+% Just as an example, this will create the text input required
+% to create N evenly-spaced observation prototypes. There are no
+% observation values in these. The output file name is 'input_to_create_obs_sequence.txt'
+%
+% fname is the name of the file to read AND WRITE.
+% forcing is the forcing for the extended state - may either be a scalar
+% or an array with exactly 40 elements.
+%
+% (bad) Example (produces a plot):
+% fname = 'perfect_ics';
+% forcing = randn(1,40) * 8.0;
+% set_extended_state(fname,'forcing',forcing)
+%
+% Example - no graphic output:
+% set_extended_state(fname,'forcing',forcing,'graphic',0)
+%
+% Example - unusual forcing:
+% forcing = [ones(1,10)*8 ones(1,30)*7] + randn(1,40);
+% set_extended_state(fname,'forcing',forcing)
+
+%% 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
+%
+% DART $Id$
+
+%%------------------------------------------------------------------------------
+% This bit just parses the input and uses default values if need be.
+
+default_forcing = 40;
+default_graphic = 1;
+p = inputParser;
+
+addRequired(p,'fname', at ischar);
+addParamValue(p,'forcing',default_forcing, at isnumeric);
+addParamValue(p,'graphic',default_graphic, at isnumeric);
+parse(p, fname, varargin{:});
+
+if ~isempty(fieldnames(p.Unmatched))
+ disp('Extra inputs:')
+ disp(p.Unmatched)
+end
+
+forcing = p.Results.forcing;
+graphic = p.Results.graphic;
+
+if (exist(fname,'file') ~= 2), error('file <%s> does not exist',fname), end
+
+%%------------------------------------------------------------------------------
+% Read the existing data
+% Check to make sure the read was successful.
+
+fid = fopen(fname,'r');
+timeline = fgetl(fid); % first record is the valid time of the model state
+[modelstate, nlines] = fscanf(fid,'%f'); % vectorized read of the rest of the file.
+fclose(fid);
+
+if (nlines ~= 80)
+ fprintf('Read %d elements for the model state.\n',nlines)
+ fprintf('Expected to read 80 elements.\n')
+ error('%s not the correct size for forced_lorenz_96.',fname)
+end
+
+%%------------------------------------------------------------------------------
+% Check user input for sanity
+
+if numel(forcing) == 1
+ modelstate(41:80) = forcing;
+elseif numel(forcing) == 40
+ modelstate(41:80) = forcing;
+else
+ fprintf('forcing array can have exactly 1 or exactly 40 values.\n')
+ error('forcing array had %d values.',numel(forcing))
+end
+
+%%------------------------------------------------------------------------------
+% open the file
+% write the timestamp
+% write the model state.
+
+fid = fopen(fname,'w');
+fprintf(fid,'%s\n',timeline);
+fprintf(fid,' %.14g\n',modelstate);
+fclose(fid);
+
+%%------------------------------------------------------------------------------
+% Optional visualizer of the observation locations.
+
+if graphic ~= 0
+
+ xax = 1:40;
+
+ figure; clf
+
+ subplot(3,1,1)
+ h1 = plot(xax,modelstate(1:40),'k*-');
+ set(h1,'MarkerSize',10,'LineWidth',2)
+ title('Original model state','FontSize',14)
+ xlabel('state index')
+
+ subplot(3,1,2)
+ h2 = plot(xax+40,modelstate(41:80),'b*-');
+ set(h2,'MarkerSize',20,'LineWidth',2)
+ title('Forcing','FontSize',14);
+ xlabel('"extended" state index')
+
+ subplot(3,1,3)
+ h3 = plot(modelstate,'b*-');
+ set(h3,'MarkerSize',10,'LineWidth',1.5)
+ title('Entire model state','FontSize',14);
+ xlabel('state index')
+
+end
+
+end % main function
+
+% <next few lines under version control, do not edit>
+% $URL$
+% $Revision$
+% $Date$
+
+function set_extended_state(fname, varargin)
+%% set_extended_state modifies the extended state portion of the perfect_ics file.
+%
+% Just as an example, this will create the text input required
+% to create N evenly-spaced observation prototypes. There are no
+% observation values in these. The output file name is 'input_to_create_obs_sequence.txt'
+%
+% fname is the name of the file to read AND WRITE.
+% forcing is the forcing for the extended state - may either be a scalar
+% or an array with exactly 40 elements.
+%
+% (bad) Example (produces a plot):
+% fname = 'perfect_ics';
+% forcing = randn(1,40) * 8.0;
+% set_extended_state(fname,'forcing',forcing)
+%
+% Example - no graphic output:
+% set_extended_state(fname,'forcing',forcing,'graphic',0)
+%
+% Example - unusual forcing:
+% forcing = [ones(1,10)*8 ones(1,30)*7] + randn(1,40);
+% set_extended_state(fname,'forcing',forcing)
+
+%% 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
+%
+% DART $Id$
+
+%%------------------------------------------------------------------------------
+% This bit just parses the input and uses default values if need be.
+
+default_forcing = 40;
+default_graphic = 1;
+p = inputParser;
+
+addRequired(p,'fname', at ischar);
+addParamValue(p,'forcing',default_forcing, at isnumeric);
+addParamValue(p,'graphic',default_graphic, at isnumeric);
+parse(p, fname, varargin{:});
+
+if ~isempty(fieldnames(p.Unmatched))
+ disp('Extra inputs:')
+ disp(p.Unmatched)
+end
+
+forcing = p.Results.forcing;
+graphic = p.Results.graphic;
+
+if (exist(fname,'file') ~= 2), error('file <%s> does not exist',fname), end
+
+%%------------------------------------------------------------------------------
+% Read the existing data
+% Check to make sure the read was successful.
+
+fid = fopen(fname,'r');
+timeline = fgetl(fid); % first record is the valid time of the model state
+[modelstate, nlines] = fscanf(fid,'%f'); % vectorized read of the rest of the file.
+fclose(fid);
+
+if (nlines ~= 80)
+ fprintf('Read %d elements for the model state.\n',nlines)
+ fprintf('Expected to read 80 elements.\n')
+ error('%s not the correct size for forced_lorenz_96.',fname)
+end
+
+%%------------------------------------------------------------------------------
+% Check user input for sanity
+
+if numel(forcing) == 1
+ modelstate(41:80) = forcing;
+elseif numel(forcing) == 40
+ modelstate(41:80) = forcing;
+else
+ fprintf('forcing array can have exactly 1 or exactly 40 values.\n')
+ error('forcing array had %d values.',numel(forcing))
+end
+
+%%------------------------------------------------------------------------------
+% open the file
+% write the timestamp
+% write the model state.
+
+fid = fopen(fname,'w');
+fprintf(fid,'%s\n',timeline);
+fprintf(fid,' %.14g\n',modelstate);
+fclose(fid);
+
+%%------------------------------------------------------------------------------
+% Optional visualizer of the observation locations.
+
+if graphic ~= 0
+
+ xax = 1:40;
+
+ figure; clf
+
+ subplot(3,1,1)
+ h1 = plot(xax,modelstate(1:40),'k*-');
+ set(h1,'MarkerSize',10,'LineWidth',2)
+ title('Original model state','FontSize',14)
+ xlabel('state index')
+
+ subplot(3,1,2)
+ h2 = plot(xax+40,modelstate(41:80),'b*-');
+ set(h2,'MarkerSize',20,'LineWidth',2)
+ title('Forcing','FontSize',14);
+ xlabel('"extended" state index')
+
+ subplot(3,1,3)
+ h3 = plot(modelstate,'b*-');
+ set(h3,'MarkerSize',10,'LineWidth',1.5)
+ title('Entire model state','FontSize',14);
+ xlabel('state index')
+
+end
+
+end % main function
+
+% <next few lines under version control, do not edit>
+% $URL$
+% $Revision$
+% $Date$
+
+function set_extended_state(fname, varargin)
+%% set_extended_state modifies the extended state portion of the perfect_ics file.
+%
+% Just as an example, this will create the text input required
+% to create N evenly-spaced observation prototypes. There are no
+% observation values in these. The output file name is 'input_to_create_obs_sequence.txt'
+%
+% fname is the name of the file to read AND WRITE.
+% forcing is the forcing for the extended state - may either be a scalar
+% or an array with exactly 40 elements.
+%
+% (bad) Example (produces a plot):
+% fname = 'perfect_ics';
+% forcing = randn(1,40) * 8.0;
+% set_extended_state(fname,'forcing',forcing)
+%
+% Example - no graphic output:
+% set_extended_state(fname,'forcing',forcing,'graphic',0)
+%
+% Example - unusual forcing:
+% forcing = [ones(1,10)*8 ones(1,30)*7] + randn(1,40);
+% set_extended_state(fname,'forcing',forcing)
+
+%% 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
+%
+% DART $Id$
+
+%%------------------------------------------------------------------------------
+% This bit just parses the input and uses default values if need be.
+
+default_forcing = 40;
+default_graphic = 1;
+p = inputParser;
+
+addRequired(p,'fname', at ischar);
+addParamValue(p,'forcing',default_forcing, at isnumeric);
+addParamValue(p,'graphic',default_graphic, at isnumeric);
+parse(p, fname, varargin{:});
+
+if ~isempty(fieldnames(p.Unmatched))
+ disp('Extra inputs:')
+ disp(p.Unmatched)
+end
+
+forcing = p.Results.forcing;
+graphic = p.Results.graphic;
+
+if (exist(fname,'file') ~= 2), error('file <%s> does not exist',fname), end
+
+%%------------------------------------------------------------------------------
+% Read the existing data
+% Check to make sure the read was successful.
+
+fid = fopen(fname,'r');
+timeline = fgetl(fid); % first record is the valid time of the model state
+[modelstate, nlines] = fscanf(fid,'%f'); % vectorized read of the rest of the file.
+fclose(fid);
+
+if (nlines ~= 80)
+ fprintf('Read %d elements for the model state.\n',nlines)
+ fprintf('Expected to read 80 elements.\n')
+ error('%s not the correct size for forced_lorenz_96.',fname)
+end
+
+%%------------------------------------------------------------------------------
+% Check user input for sanity
+
+if numel(forcing) == 1
+ modelstate(41:80) = forcing;
+elseif numel(forcing) == 40
+ modelstate(41:80) = forcing;
+else
+ fprintf('forcing array can have exactly 1 or exactly 40 values.\n')
+ error('forcing array had %d values.',numel(forcing))
+end
+
+%%------------------------------------------------------------------------------
+% open the file
+% write the timestamp
+% write the model state.
+
+fid = fopen(fname,'w');
+fprintf(fid,'%s\n',timeline);
+fprintf(fid,' %.14g\n',modelstate);
+fclose(fid);
+
+%%------------------------------------------------------------------------------
+% Optional visualizer of the observation locations.
+
+if graphic ~= 0
+
+ xax = 1:40;
+
+ figure; clf
+
+ subplot(3,1,1)
+ h1 = plot(xax,modelstate(1:40),'k*-');
+ set(h1,'MarkerSize',10,'LineWidth',2)
+ title('Original model state','FontSize',14)
+ xlabel('state index')
+
+ subplot(3,1,2)
+ h2 = plot(xax+40,modelstate(41:80),'b*-');
+ set(h2,'MarkerSize',20,'LineWidth',2)
+ title('Forcing','FontSize',14);
+ xlabel('"extended" state index')
+
+ subplot(3,1,3)
+ h3 = plot(modelstate,'b*-');
+ set(h3,'MarkerSize',10,'LineWidth',1.5)
+ title('Entire model state','FontSize',14);
+ xlabel('state index')
+
+end
+
+end % main function
+
+% <next few lines under version control, do not edit>
+% $URL$
+% $Revision$
+% $Date$
+
Property changes on: DART/trunk/models/forced_lorenz_96/matlab/set_extended_state.m
___________________________________________________________________
Added: svn:mime-type
+ text/plain
Added: svn:keywords
+ Date Rev Author HeadURL Id
Added: svn:eol-style
+ native
Modified: DART/trunk/models/forced_lorenz_96/model_mod.html
===================================================================
--- DART/trunk/models/forced_lorenz_96/model_mod.html 2015-05-14 22:10:54 UTC (rev 7984)
+++ DART/trunk/models/forced_lorenz_96/model_mod.html 2015-05-15 19:47:24 UTC (rev 7985)
@@ -34,17 +34,18 @@
<A HREF="#PrivateComponents">PRIVATE COMPONENTS</A> /
<A HREF="#Legalese">TERMS OF USE</A>
-<H2>Overview</H2>
-
<P>
-DART interface module for the forced_lorenz-96 model. The 16 public interfaces
+DART interface module for the <em class=program>forced_lorenz_96</em> model. The 16 public interfaces
are standardized for all DART compliant models. These interfaces allow
DART to advance the model, get the model state and metadata describing
this state, find state variables that are close to a given location,
and do spatial interpolation for model state variables.
</P>
+
+<H2>Overview</H2>
+
<P>
-The forced_lorenz_96 model implements the standard L96 equations except
+The <em class=program>forced_lorenz_96</em> model implements the standard L96 equations except
that the forcing term, F, is added to the state vector and is
assigned an independent value at each gridpoint. The result is a
model that is twice as big as the standard L96 model. The forcing
@@ -56,7 +57,137 @@
If the random noise option is turned off (see namelist) the
time tendency of the forcing terms is 0.
</P>
+<center>
+<table width=90% border=0 summary="example of state vector partitioning " cellpadding=8>
+<tr><th colspan = 2>DART state vector composition</th></tr>
+<tr><th width=50% bgcolor=Cyan>state variables<br />indices 1 - 40</th>
+ <th width=50% bgcolor=Gold>forcing terms<br />indices 41 - 80</th></tr>
+<tr><th>traditional Lorenz_96 state</th>
+ <th>"extended" state</th></tr>
+</table>
+</center>
+<P>
+The DART tutorial Section 20 <a href="../../tutorial/section_20.pdf">[pdf]</a>
+explores parameter estimation using the <em class=program>forced_lorenz_96</em> model.
+</P>
+
+<H2>Quick Example</H2>
+<P>To become familiar with the model, try this quick experiment.</P>
+<ol>
+ <li>compile everything in the
+ <em class=file>model/forced_lorenz_96/work</em> directory.<br /><br />
+ <div class=unix>
+ cd models/forced_lorenz_96/work<br />
+ ./quickbuild.csh</div></li>
+<br />
+ <li>make sure the <em class=file>input.nml</em> looks like the following
+ (there is a lot that has been left out for clarity,
+ these are the settings of interest for this example):<br /><br />
+ <div class=routine>
+<pre>&perfect_model_obs_nml
+ start_from_restart = .true.,
+ output_restart = .true.,
+ async = 0,
+ restart_in_file_name = "perfect_ics",
+ obs_seq_in_file_name = "obs_seq.in",
+ obs_seq_out_file_name = "obs_seq.out",
+ ...
+
+&filter_nml
+ async = 0,
+ ens_size = 80,
+ start_from_restart = .true.,
+ output_restart = .true.,
+ obs_sequence_in_name = "obs_seq.out",
+ obs_sequence_out_name = "obs_seq.final",
+ restart_in_file_name = "filter_ics",
+ restart_out_file_name = "filter_restart",
+ num_output_state_members = 80,
+ num_output_obs_members = 80,
+ ...
+
+&model_nml
+ num_state_vars = 40,
+ forcing = 8.00,
+ delta_t = 0.05,
+ time_step_days = 0,
+ time_step_seconds = 3600,
+ reset_forcing = .false.,
+ random_forcing_amplitude = 0.10
+ /</pre> </div> </li>
+<br />
+<li>Run <em class=program>perfect_model_obs</em> to generate
+ <em class=file>True_State.nc</em> and
+ <em class=file>obs_seq.out</em>. The default <em class=file>obs_seq.in</em>
+ will cause the model to advance for 1000 time steps.<br /><br />
+ <div class=unix>
+ ./perfect_model_obs
+ </div>
+</li>
+<br />
+<li>If you have <em class=program>ncview</em>,
+ explore the <em class=file>True_State.nc</em>. Notice that the
+ State Variable indices from 1-40 are the dynamical part of the
+ model and 41-80 are the Forcing variables.<br /><br />
+ <div class=unix>
+ ncview True_State.nc
+ </div>
+</li>
+<br />
+<li>Run <em class=program>filter</em> to generate
+ <em class=file>Prior_Diag.nc</em>,
+ <em class=file>Posterior_Diag.nc</em> and
+ <em class=file>obs_seq.final</em>.<br /><br />
+ <div class=unix>
+ ./filter
+ </div>
+</li>
+<br />
+<li>Launch Matlab and run <em class=program>plot_ens_time_series</em>.<br /><br />
+ <div class=unix>
+ >> plot_ens_time_series<br />
+
+ Input name of prior or posterior diagnostics file:<br />
+ <cr> for Prior_Diag.nc<br />
+ <em class=input>Prior_Diag.nc</em><br />
+ <br />
+ OPTIONAL: if you have the true state and want it superimposed, provide<br />
+ : the name of the input file. If not, enter a dummy filename.<br />
+ : Input name of True State file:<br />
+ <cr> for True_State.nc<br />
+ <em class=input>True_State.nc</em><br />
+ <br />
+ Using state state variable IDs 1 13 27<br />
+ If these are OK, <cr>;<br />
+ If not, please enter array of state variable ID's<br />
+ To choose from entire state enter A 25 50 75 (between 1 and 80)<br />
+ To choose traditional model state enter S 1 23 40 (between 1 and 40)<br />
+ To choose forcing estimates enter F 2 12 22 (between 1 and 40)<br />
+ (no intervening syntax required)<br />
+ <em class=input>A 20 30 40 60 70 80</em>
+ </div>
+ <P>
+ Indices 20, 30, and 40 will be from the dynamical part of the
+ lorenz_96 attractor, indices 60, 70, and 80 will be the
+ corresponding Forcing values. Here are some images for just indices 20
+ and 60. Click on each image for a high-res version.
+ </P>
+</li>
+<br />
+<center>
+<a href="../../doc/images/forced_lorenz_96_evolution_state.png">
+<img src="../../doc/images/forced_lorenz_96_evolution_state_small.png" alt="forced_lorenz_96evolution" /></a>
+<a href="../../doc/images/forced_lorenz_96_evolution_forcing.png">
+<img src="../../doc/images/forced_lorenz_96_evolution_forcing_small.png" alt="forced_lorenz_96evolution" /></a>
+</center>
+</ol>
+
+<P>Repeat the experiment with <em class=code>reset_forcing = .true.</em> when creating the
+true state and <em class=code>reset_forcing = .false.</em> when assimilating. What happens?</P>
+
+<P><!-- make sure the 'top' is aligned correctly --></P>
+
<!--==================================================================-->
<A NAME="OtherModulesUsed"></A>
@@ -620,8 +751,9 @@
</TABLE>
</div>
-<br>
+<P><!-- make sure the 'top' is aligned correctly --></P>
+
<!--===================== DESCRIPTION OF A NAMELIST =====================-->
<A NAME="Namelist"></A>
@@ -634,6 +766,7 @@
Character strings that contain a '/' must be
enclosed in quotes to prevent them from
prematurely terminating the namelist.
+The values shown here are the default values.
</P>
<div class=namelist>
@@ -704,8 +837,7 @@
</TABLE>
</div>
-<br />
-<br />
+<P><!-- make sure the 'top' is aligned correctly --></P>
<!--==================================================================-->
<!-- Describe the Files Used by this module. -->
@@ -715,22 +847,36 @@
<div class="top">[<a href="#">top</a>]</div><hr />
<H2>FILES</H2>
-<TABLE border=0 >
+<TABLE border=0 width=100% summary="list of files used" cellspacing=1 cellpadding=2 width=100%>
<TR><TH align=left>filename</TH>
<TH align=left>purpose</TH></TR>
-<TR><TD>input.nml</TD>
+
+<TR><TD><em class=file>input.nml</em></TD>
<TD>to read the model_mod namelist</TD></TR>
-<TR><TD>Prior_Diag.nc</TD>
+
+<TR><TD><em class=file>True_State.nc</em></TD>
+ <TD>the time-history of the single model state used to generate the synthetic observations</TD></TR>
+
+<TR><TD><em class=file>Prior_Diag.nc</em></TD>
<TD>the time-history of the model state before assimilation</TD></TR>
-<TR><TD>Posterior_Diag.nc </TD>
+
+<TR><TD><em class=file>Posterior_Diag.nc </em></TD>
<TD>the time-history of the model state after assimilation</TD></TR>
-<TR><TD>dart_log.out [default name]</TD>
+
+<TR><TD><em class=file>advance_model.csh</em></TD>
+ <TD>shell script to advance the model as a standalone executable.
+ Each advance will then read the model's namelist.</TD></TR>
+
+<TR><TD><em class=file>dart_log.out</em></TD>
<TD>the run-time diagnostic output</TD></TR>
-<TR><TD>dart_log.nml [default name]</TD>
+
+<TR><TD><em class=file>dart_log.nml</em></TD>
<TD>the record of all the namelists actually USED -
contains the default values</TD></TR>
</TABLE>
+<P><!-- make sure the 'top' is aligned correctly --></P>
+
<!--==================================================================-->
<!-- Cite references, if need be. -->
<!--==================================================================-->
@@ -742,6 +888,8 @@
<li> none </li>
</ul>
+<P><!-- make sure the 'top' is aligned correctly --></P>
+
<!--==================================================================-->
<!-- Describe all the error conditions and codes. -->
<!--==================================================================-->
Added: DART/trunk/models/forced_lorenz_96/shell_scripts/README
===================================================================
--- DART/trunk/models/forced_lorenz_96/shell_scripts/README (rev 0)
+++ DART/trunk/models/forced_lorenz_96/shell_scripts/README 2015-05-15 19:47:24 UTC (rev 7985)
@@ -0,0 +1,18 @@
+The scripts in this directory are:
+
+advance_model.csh -
+Normally small models (including this one) are compiled directly with
+the filter program to make a single executable that does both the data
+assimilation and the model advances. However, to test the mode where the
+model is a separate executable, this script is an example of what is needed.
+See the filter.html documentation for async mode 2 for more details.
+
+The scripts in this directory are:
+
+advance_model.csh -
+Normally small models (including this one) are compiled directly with
+the filter program to make a single executable that does both the data
+assimilation and the model advances. However, to test the mode where the
+model is a separate executable, this script is an example of what is needed.
+See the filter.html documentation for async mode 2 for more details.
+
Property changes on: DART/trunk/models/forced_lorenz_96/shell_scripts/README
___________________________________________________________________
Added: svn:mime-type
+ text/plain
Added: svn:keywords
+ Date Rev Author HeadURL Id
Added: svn:eol-style
+ native
Added: DART/trunk/models/forced_lorenz_96/shell_scripts/advance_model.csh
===================================================================
--- DART/trunk/models/forced_lorenz_96/shell_scripts/advance_model.csh (rev 0)
+++ DART/trunk/models/forced_lorenz_96/shell_scripts/advance_model.csh 2015-05-15 19:47:24 UTC (rev 7985)
@@ -0,0 +1,274 @@
+#!/bin/csh -f
+#
+# 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
+#
+# DART $Id$
+#
+# Standard script for use in assimilation applications
+# where the model advance is executed as a separate process.
+# Can be used as-is with most low-order models and the bgrid model which
+# can be advanced using the integrate_model executable.
+#
+# Arguments are (created by 'filter' or 'perfect_model_obs' and include):
+# 1) the process number of caller,
+# 2) the number of ensemble members/state copies belonging to that process, and
+# 3) the name of the control_file for that process.
+#
+# If this script finishes and the 'control_file' still exists, it is
+# an ERROR CONDITION and means one or more of the ensemble members did
+# not advance properly. Despite our best attempts to trap on this
+# condition, some MPI installations simply hang, some properly terminate.
+#
+# This script loops over all the entries in the control_file to advance
+# any/all of the ensemble members. The number of trips through the
+# loop is the second argument to this script. The control_file contains
+# the information about which ensemble members are to be advanced by THIS TASK.
+# Sometimes it may be just one ensemble member, sometimes all of them.
+# Read DART/doc/html/filter_async_modes.html and the mpi_intro.html
+# for an overview.
+#
+# This script has 4 logical 'blocks':
+# 1) creates a clean, temporary directory in which to run a model instance
@@ Diff output truncated at 40000 characters. @@
More information about the Dart-dev
mailing list