[Dart-dev] [8126] DART/trunk/DART_LAB/matlab: Renamed the files to be in the private directory, nothing more.
nancy at ucar.edu
nancy at ucar.edu
Wed Jun 24 15:21:54 MDT 2015
Revision: 8126
Author: thoar
Date: 2015-06-24 15:21:54 -0600 (Wed, 24 Jun 2015)
Log Message:
-----------
Renamed the files to be in the private directory, nothing more.
Removed Paths:
-------------
DART/trunk/DART_LAB/matlab/advance_oned.m
DART/trunk/DART_LAB/matlab/comp_cov_factor.m
DART/trunk/DART_LAB/matlab/g_prod_plot.m
DART/trunk/DART_LAB/matlab/get_ens_rank.m
DART/trunk/DART_LAB/matlab/get_state_increments.m
DART/trunk/DART_LAB/matlab/kurt.m
DART/trunk/DART_LAB/matlab/lorenz_63_adv_1step.m
DART/trunk/DART_LAB/matlab/lorenz_63_static_init_model.m
DART/trunk/DART_LAB/matlab/lorenz_96_adv_1step.m
DART/trunk/DART_LAB/matlab/lorenz_96_static_init_model.m
DART/trunk/DART_LAB/matlab/norm_inv.m
DART/trunk/DART_LAB/matlab/norm_pdf.m
DART/trunk/DART_LAB/matlab/obs_increment_eakf.m
DART/trunk/DART_LAB/matlab/obs_increment_enkf.m
DART/trunk/DART_LAB/matlab/obs_increment_rhf.m
DART/trunk/DART_LAB/matlab/plot_gaussian.m
DART/trunk/DART_LAB/matlab/plot_polar.m
DART/trunk/DART_LAB/matlab/product_of_gaussians.m
-------------- next part --------------
Deleted: DART/trunk/DART_LAB/matlab/advance_oned.m
===================================================================
--- DART/trunk/DART_LAB/matlab/advance_oned.m 2015-06-24 21:20:35 UTC (rev 8125)
+++ DART/trunk/DART_LAB/matlab/advance_oned.m 2015-06-24 21:21:54 UTC (rev 8126)
@@ -1,28 +0,0 @@
-function x_new = advance_oned(x, alpha, model_bias)
-%% advance_oned(x, alpha, model_bias)
-
-%% 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$
-
-x_new = x + comp_dt(x, alpha, model_bias);
-end
-
-%---------------------------------------------------
-
-% Internal function comp_dt
-function dx = comp_dt(x, alpha, model_bias)
-
-% Compute the time tendency; alpha controls nonlinearity
-% model_bias controls a shift in the model dynamics
-dx = (x + model_bias) + alpha .* x .* abs(x);
-
-end
-
-% <next few lines under version control, do not edit>
-% $URL$
-% $Revision$
-% $Date$
-
Deleted: DART/trunk/DART_LAB/matlab/comp_cov_factor.m
===================================================================
--- DART/trunk/DART_LAB/matlab/comp_cov_factor.m 2015-06-24 21:20:35 UTC (rev 8125)
+++ DART/trunk/DART_LAB/matlab/comp_cov_factor.m 2015-06-24 21:21:54 UTC (rev 8126)
@@ -1,26 +0,0 @@
-function cov_factor = comp_cov_factor(z_in, c)
-%% comp_cov_factor Gaspari Cohn cutoff, z_in is the distance while c is the cutoff
-
-%% 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$
-
-z = abs(z_in);
-
-if( z >= c*2.0)
- cov_factor = 0;
-elseif( z <= c )
- r = z / c;
- cov_factor = ((( -0.25*r +0.5)*r +0.625)*r -5.0/3.0)*r^2 + 1.0;
-else
- r = z / c;
- cov_factor = ((((r/12 -0.5)*r +0.625)*r +5.0/3.0)*r -5.0)*r + 4.0 - 2.0 / (3.0 * r);
-end
-
-% <next few lines under version control, do not edit>
-% $URL$
-% $Revision$
-% $Date$
-
Deleted: DART/trunk/DART_LAB/matlab/g_prod_plot.m
===================================================================
--- DART/trunk/DART_LAB/matlab/g_prod_plot.m 2015-06-24 21:20:35 UTC (rev 8125)
+++ DART/trunk/DART_LAB/matlab/g_prod_plot.m 2015-06-24 21:21:54 UTC (rev 8126)
@@ -1,108 +0,0 @@
-function [prior_mean, prior_sd, obs_mean, obs_err_sd, is_err] = g_prod_plot(h)
-%% g_prod_plot Updates the plot of the prior and observation gaussians
-
-%% 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$
-
-% Successful return as default
-is_err = false;
-
-% Default failed returns for other quantities
-prior_mean = 0; prior_sd = -1;
-obs_mean = 0; obs_err_sd = -1;
-
-h_prior_mean = get(h.edit1);
-prior_mean = str2double(h_prior_mean.String);
-% The mean must be a number
-if(isnan(prior_mean))
- error_banner(h, 'Prior Mean must be a number');
- is_err = true;
- return
-end
-
-h_prior_sd = get(h.edit2);
-prior_sd = str2double(h_prior_sd.String);
-
-% Prior sd must be a number
-if(isnan(prior_sd))
- error_banner(h, 'Prior SD must be a postive number');
- is_err = true;
- return
-end
-
-% Prior sd must also be positive
-if(prior_sd <= 0)
- error_banner(h, 'Prior SD must be positive')
- is_err = true;
- return
-end
-
-hold off
-prior_handle = plot_gaussian(prior_mean, prior_sd, 1);
-set(prior_handle, 'Color', [0 0.73 0], 'LineWidth', 2);
-hold on
-
-h_obs_mean = get(h.edit3);
-obs_mean = str2double(h_obs_mean.String);
-
-% Obs value must be a number
-if(isnan(obs_mean))
- error_banner(h, 'Obs value must be a number');
- is_err = true;
- return
-end
-
-h_obs_err_sd = get(h.edit4);
-obs_err_sd = str2double(h_obs_err_sd.String);
-
-% Obs error sd must be a positive number
-if(isnan(obs_err_sd))
- error_banner(h, 'Obs Error SD must be a positive number');
- is_err = true;
- return
-end
-
-if(obs_err_sd <= 0)
- error_banner(h, 'Obs Error SD must be positive');
- is_err = true;
- return
-end
-
-obs_handle = plot_gaussian(obs_mean, obs_err_sd, 1);
-set(obs_handle, 'Color', 'r', 'LineStyle', '--', 'LineWidth', 2);
-
-% Put on a legend
-legend('Prior', 'Obs. Likelihood');
-
-end
-
-%---------------------------------------------------------
-
-% Internal function to write error banner
-function error_banner(h, message)
-
- hold off
- x= [1 2];
- plot(x, 'Visible', 'off')
- h_fig = get(h.figure1);
- x_limits = get(h_fig.CurrentAxes, 'Xlim');
- y_limits = get(h_fig.CurrentAxes, 'Ylim');
- text(x_limits(1) * 7/8 + x_limits(2) / 8, mean(y_limits), ...
- message, 'fontsize', 16, 'Color', 'r');
-
- % While we're at it, clear out the values for posterior, too
- set(h.text7, 'String', 'Posterior Mean = ');
- set(h.text8, 'String', 'Posterior SD = ');
- set(h.text9, 'String', 'Weight = ');
- return;
-
-end
-
-% <next few lines under version control, do not edit>
-% $URL$
-% $Revision$
-% $Date$
-
Deleted: DART/trunk/DART_LAB/matlab/get_ens_rank.m
===================================================================
--- DART/trunk/DART_LAB/matlab/get_ens_rank.m 2015-06-24 21:20:35 UTC (rev 8125)
+++ DART/trunk/DART_LAB/matlab/get_ens_rank.m 2015-06-24 21:21:54 UTC (rev 8126)
@@ -1,22 +0,0 @@
-function [rank] = get_ens_rank(ens, x)
-%% get_ens_rank
-
-%% 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$
-
-s_ens = sort(ens);
-rank = find(s_ens < squeeze(x),1,'last') + 1;
-if(isempty(rank))
- rank = 1;
-end
-
-end
-
-% <next few lines under version control, do not edit>
-% $URL$
-% $Revision$
-% $Date$
-
Deleted: DART/trunk/DART_LAB/matlab/get_state_increments.m
===================================================================
--- DART/trunk/DART_LAB/matlab/get_state_increments.m 2015-06-24 21:20:35 UTC (rev 8125)
+++ DART/trunk/DART_LAB/matlab/get_state_increments.m 2015-06-24 21:21:54 UTC (rev 8126)
@@ -1,21 +0,0 @@
-function [state_incs] = get_state_increments(state_ens, obs_ens, obs_incs)
-%% get_state_increments Computes state increments given observation increments and
-% the state and obs prior ensembles
-
-%% 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$
-
-% Compute state variance and covariance
-covar = cov(state_ens, obs_ens);
-
-
-state_incs = obs_incs * covar(1, 2) / covar(2, 2);
-
-% <next few lines under version control, do not edit>
-% $URL$
-% $Revision$
-% $Date$
-
Deleted: DART/trunk/DART_LAB/matlab/kurt.m
===================================================================
--- DART/trunk/DART_LAB/matlab/kurt.m 2015-06-24 21:20:35 UTC (rev 8125)
+++ DART/trunk/DART_LAB/matlab/kurt.m 2015-06-24 21:21:54 UTC (rev 8126)
@@ -1,30 +0,0 @@
-function k = kurt(vals)
-% computes the kurtosis of the given input array
-%
-% based on the second formula on this web page:
-% http://www.ats.ucla.edu/stat/mult_pkg/faq/general/kurtosis.htm
-
-%% 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$
-
-% array of diffs from mean
-del = vals - mean(vals);
-
-% compute the square and 4th power of the diffs from mean
-m2 = mean(del .^ 2);
-m4 = mean(del .^ 4);
-
-% compute the kurtosis value. this is the version
-% of the kurtosis formula that is not nonbiased and
-% does not subtract 3.0 from the result.
-k = (m4 ./ (m2 .^ 2));
-
-end
-
-% <next few lines under version control, do not edit>
-% $URL$
-% $Revision$
-% $Date$
Deleted: DART/trunk/DART_LAB/matlab/lorenz_63_adv_1step.m
===================================================================
--- DART/trunk/DART_LAB/matlab/lorenz_63_adv_1step.m 2015-06-24 21:20:35 UTC (rev 8125)
+++ DART/trunk/DART_LAB/matlab/lorenz_63_adv_1step.m 2015-06-24 21:21:54 UTC (rev 8126)
@@ -1,51 +0,0 @@
-function[x_new, time_new] = lorenz_63_adv_1step(x, time)
-%% lorenz_63)adv_1step advances the lorenz convective 3 variable model
-% for a single two step runge-kutta time step
-%
-% x is the 3-vector state, time is the 2-vector days and seconds time
-
-%% 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$
-
-global DELTAT
-
-% Set the fraction for the rk-2 time step
-fract = 1;
-
-% Compute the first intermediate step
-dx = comp_dt(x);
-x1 = x + fract * DELTAT * dx;
-
-% Compute the second intermediate step
-dx = comp_dt(x1);
-x2 = x1 + fract * DELTAT * dx;
-
-% New value for x is average of original value and second intermediate
-x_new = (x + x2) / 2;
-
-% Update the time ; Non-dimensional single unit for development
-time_new = time + 1;
-
-end
-
-
-function[dt] = comp_dt(x)
-
-global SIGMA
-global R
-global B
-
-dt(1) = SIGMA * (x(2) - x(1));
-dt(2) = -1 * x(1)*x(3) + R*x(1) - x(2);
-dt(3) = x(1)*x(2) - B*x(3);
-
-end
-
-% <next few lines under version control, do not edit>
-% $URL$
-% $Revision$
-% $Date$
-
Deleted: DART/trunk/DART_LAB/matlab/lorenz_63_static_init_model.m
===================================================================
--- DART/trunk/DART_LAB/matlab/lorenz_63_static_init_model.m 2015-06-24 21:20:35 UTC (rev 8125)
+++ DART/trunk/DART_LAB/matlab/lorenz_63_static_init_model.m 2015-06-24 21:21:54 UTC (rev 8126)
@@ -1,40 +0,0 @@
-%% lorenz_63_static_init_model Initializes class data for L63, sets up global storage
-% and reads in control data from input file
-
-%% 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$
-
-% Lorenz-63 model parameters
-global SIGMA
-global R
-global B
-global DELTAT
-global TIME_STEP_DAYS
-global TIME_STEP_SECONDS
-
-% Set default values for the model parameters
-SIGMA = 10;
-R = 28;
-B = 8/3;
-DELTAT = 0.01;
-TIME_STEP_DAYS = 0;
-TIME_STEP_SECONDS = 0;
-
-
-% Lorenz-63 fixed model parameters
-global MODEL_SIZE
-
-MODEL_SIZE = 3;
-
-global STATE_LOC
-
-STATE_LOC = (0:2) / 3;
-
-% <next few lines under version control, do not edit>
-% $URL$
-% $Revision$
-% $Date$
-
Deleted: DART/trunk/DART_LAB/matlab/lorenz_96_adv_1step.m
===================================================================
--- DART/trunk/DART_LAB/matlab/lorenz_96_adv_1step.m 2015-06-24 21:20:35 UTC (rev 8125)
+++ DART/trunk/DART_LAB/matlab/lorenz_96_adv_1step.m 2015-06-24 21:21:54 UTC (rev 8126)
@@ -1,70 +0,0 @@
-function[x_new, time_new] = lorenz_96_adv_1step(x, time)
-%% lorenz_96_adv_1step Does a single time step advance for lorenz_96 40-variable model using four step runge-kutta time step
-%
-% x is the 40-vector state, time is the 2-vector days and seconds time
-
-%% 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$
-
-global DELTA_T
-
-% Compute first intermediate step
-dx = comp_dt(x);
-x1 = DELTA_T * dx;
-inter = x + x1 / 2;
-
-% Compute second intermediate step
-dx = comp_dt(inter);
-x2 = DELTA_T * dx;
-inter = x + x2 / 2;
-
-% Compute third intermediate step
-dx = comp_dt(inter);
-x3 = DELTA_T * dx;
-inter = x + x3;
-
-% Compute fourth intermediate step
-dx = comp_dt(inter);
-x4 = DELTA_T * dx;
-
-% Compute new value for x
-x_new = x + x1/6 + x2/3 + x3/3 + x4/6;
-
-% Increment time step
-time_new = time + 1;
-
-
-end
-
-%------------------------------------------------------------------------------
-
-function[dt] = comp_dt(x)
-
-global FORCING
-global MODEL_SIZE
-
-dt = zeros(1,MODEL_SIZE);
-
-for j = 1:MODEL_SIZE
- jp1 = j + 1;
- if(jp1 > MODEL_SIZE), jp1 = 1; end
-
- jm2 = j - 2;
- if(jm2 < 1), jm2 = MODEL_SIZE + jm2; end
-
- jm1 = j - 1;
- if(jm1 < 1), jm1 = MODEL_SIZE; end
-
- dt(j) = (x(jp1) - x(jm2)) * x(jm1) - x(j) + FORCING;
-end
-
-end
-
-% <next few lines under version control, do not edit>
-% $URL$
-% $Revision$
-% $Date$
-
Deleted: DART/trunk/DART_LAB/matlab/lorenz_96_static_init_model.m
===================================================================
--- DART/trunk/DART_LAB/matlab/lorenz_96_static_init_model.m 2015-06-24 21:20:35 UTC (rev 8125)
+++ DART/trunk/DART_LAB/matlab/lorenz_96_static_init_model.m 2015-06-24 21:21:54 UTC (rev 8126)
@@ -1,36 +0,0 @@
-%% lorenz_96_static_init_model Initializes class data for L96, sets up global storage
-% and reads in control data from input file
-
-%% 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$
-
-% Lorenz-96 model parameters
-global FORCING
-global DELTA_T
-global TIME_STEP_DAYS
-global TIME_STEP_SECONDS
-
-% Set default values for the model parameters
-FORCING = 8;
-DELTA_T = 0.05;
-TIME_STEP_DAYS = 0;
-TIME_STEP_SECONDS = 0;
-
-
-% Lorenz-96 fixed model parameters
-global MODEL_SIZE
-
-MODEL_SIZE = 40;
-
-global STATE_LOC
-
-STATE_LOC = (0:MODEL_SIZE - 1) / MODEL_SIZE;
-
-% <next few lines under version control, do not edit>
-% $URL$
-% $Revision$
-% $Date$
-
Deleted: DART/trunk/DART_LAB/matlab/norm_inv.m
===================================================================
--- DART/trunk/DART_LAB/matlab/norm_inv.m 2015-06-24 21:20:35 UTC (rev 8125)
+++ DART/trunk/DART_LAB/matlab/norm_inv.m 2015-06-24 21:21:54 UTC (rev 8126)
@@ -1,83 +0,0 @@
-function x = norm_inv(p)
-% Computes inverse CDF for the normal distribution,
-% evaluated at either the scalar value or array P.
-%
-% based on the code in DART assim_tools_mod.f90, which is
-% in turn based on http://home.online.no/~pjacklam/notes/invnorm
-
-%% 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$
-
-% constants where precision is vital
-a1 = -39.69683028665376;
-a2 = 220.9460984245205;
-a3 = -275.9285104469687;
-a4 = 138.357751867269;
-a5 = -30.66479806614716;
-a6 = 2.506628277459239;
-b1 = -54.4760987982241;
-b2 = 161.5858368580409;
-b3 = -155.6989798598866;
-b4 = 66.80131188771972;
-b5 = -13.28068155288572;
-c1 = -0.007784894002430293;
-c2 = -0.3223964580411365;
-c3 = -2.400758277161838;
-c4 = -2.549732539343734;
-c5 = 4.374664141464968;
-c6 = 2.938163982698783;
-d1 = 0.007784695709041462;
-d2 = 0.3224671290700398;
-d3 = 2.445134137142996;
-d4 = 3.754408661907416;
-
-% Split into an inner and two outer regions which have separate fits
-p_low = 0.02425;
-p_high = 1 - p_low;
-
-% if P is scalar, simplify
-if (numel(p) == 1)
- if (p < p_low)
- q = sqrt(-2.0 * log(p));
- x = (((((c1*q + c2)*q + c3)*q + c4)*q + c5)*q + c6) / ...
- ((((d1*q + d2)*q + d3)*q + d4)*q + 1.0);
- elseif (p >= p_low) && (p <= p_high)
- q = p - 0.5;
- r = q*q;
- x = (((((a1*r + a2)*r + a3)*r + a4)*r + a5)*r + a6)*q / ...
- (((((b1*r + b2)*r + b3)*r + b4)*r + b5)*r + 1.0);
-
- else % (p > p_high)
- q = sqrt(-2.0 * log(1.0 - p));
- x = -(((((c1*q + c2)*q + c3)*q + c4)*q + c5)*q + c6) / ...
- ((((d1*q + d2)*q + d3)*q + d4)*q + 1.0);
- end
-
-else % if P is an array, so do it element by element
-
- inds = (p < p_low);
- q(inds) = sqrt(-2.0 * log(p(inds)));
- x(inds) = (((((c1*q(inds) + c2).*q(inds) + c3).*q(inds) + c4).*q(inds) + c5).*q(inds) + c6) ./ ...
- ((((d1*q(inds) + d2).*q(inds) + d3).*q(inds) + d4).*q(inds) + 1.0);
-
- inds = (p >= p_low) & (p <= p_high);
- q(inds) = p(inds) - 0.5;
- r(inds) = q(inds).*q(inds);
- x(inds) = (((((a1*r(inds) + a2).*r(inds) + a3).*r(inds) + a4).*r(inds) + a5).*r(inds) + a6).*q(inds) ./ ...
- (((((b1*r(inds) + b2).*r(inds) + b3).*r(inds) + b4).*r(inds) + b5).*r(inds) + 1.0);
-
- inds = (p > p_high);
- q(inds) = sqrt(-2.0 * log(1.0 - p(inds)));
- x(inds) = -(((((c1*q(inds) + c2).*q(inds) + c3).*q(inds) + c4).*q(inds) + c5).*q(inds) + c6) ./ ...
- ((((d1*q(inds) + d2).*q(inds) + d3).*q(inds) + d4).*q(inds) + 1.0);
-end
-
-end
-
-% <next few lines under version control, do not edit>
-% $URL$
-% $Revision$
-% $Date$
Deleted: DART/trunk/DART_LAB/matlab/norm_pdf.m
===================================================================
--- DART/trunk/DART_LAB/matlab/norm_pdf.m 2015-06-24 21:20:35 UTC (rev 8125)
+++ DART/trunk/DART_LAB/matlab/norm_pdf.m 2015-06-24 21:21:54 UTC (rev 8126)
@@ -1,32 +0,0 @@
-function [y] = norm_pdf(x, mu, sigma)
-% computes a gaussian (normal) PDF
-% for the points of X with a given mean (mu) and standard deviation (sigma)
-%
-% normal plot, y given x:
-% y = (1 / (sigma * sqrt(2*pi))) * e ^ ((-1/2 * ((x-mu) / sigma)^2)
-% or
-% g(x) = \frac{1}{\sigma\sqrt{2\pi}} e^{ -\frac{1}{2}\left(\frac{x-\mu}{\sigma}\right)^2 }.
-%
-% see: https://en.wikipedia.org/wiki/Probability_density_function
-
-%% 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$
-
-
-e = exp(1);
-
-basen = (1.0 / (sigma * sqrt(2*pi)));
-expon = -0.5 * (((x-mu) / sigma).^2 );
-
-y = basen * (e .^ expon);
-
-end
-
-% <next few lines under version control, do not edit>
-% $URL$
-% $Revision$
-% $Date$
-
Deleted: DART/trunk/DART_LAB/matlab/obs_increment_eakf.m
===================================================================
--- DART/trunk/DART_LAB/matlab/obs_increment_eakf.m 2015-06-24 21:20:35 UTC (rev 8125)
+++ DART/trunk/DART_LAB/matlab/obs_increment_eakf.m 2015-06-24 21:21:54 UTC (rev 8126)
@@ -1,56 +0,0 @@
-function [obs_increments, err] = obs_increment_eakf(ensemble, observation, obs_error_var)
-%% obs_increment_eakf Computes increments for an ensemble adjustment filter
-
-%% 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$
-
-% Set error return to default successful
-err = 0;
-
-% Compute prior ensemble mean and variance
-prior_mean = mean(ensemble);
-prior_var = var(ensemble);
-
-% If both prior and observation error variance are return error
-if(prior_var <= 0 && obs_error_var <= 0)
- err = 1;
- return;
-end
-
-% Compute the posterior mean and variance
-% If prior variance is 0, posterior mean is prior_mean and variance is 0
-if(prior_var == 0)
- post_mean = prior_mean;
- post_var = 0;
-elseif(obs_error_var == 0)
-% If obs_error_var is 0, posterior mean is observation and variance is 0
- post_mean = observation;
- post_var = 0;
-else
-% Use product of gaussians
- % Compute the posterior variance
- post_var = 1 / (1 / prior_var + 1 / obs_error_var);
-
- % Compute posterior mean
- post_mean = post_var * (prior_mean / prior_var + observation / obs_error_var);
-end
-
-% Shift the prior ensemble to have the posterior mean
-updated_ensemble = ensemble - prior_mean + post_mean;
-
-% Contract the ensemble to have the posterior_variance
-var_ratio = post_var / prior_var;
-updated_ensemble = sqrt(var_ratio) * (updated_ensemble - post_mean) + post_mean;
-
-% Compute the increments
-obs_increments = updated_ensemble - ensemble;
-
-
-% <next few lines under version control, do not edit>
-% $URL$
-% $Revision$
-% $Date$
-
Deleted: DART/trunk/DART_LAB/matlab/obs_increment_enkf.m
===================================================================
--- DART/trunk/DART_LAB/matlab/obs_increment_enkf.m 2015-06-24 21:20:35 UTC (rev 8125)
+++ DART/trunk/DART_LAB/matlab/obs_increment_enkf.m 2015-06-24 21:21:54 UTC (rev 8126)
@@ -1,74 +0,0 @@
-function [obs_increments, err] = obs_increment_enkf(ensemble, observation, obs_error_var)
-%% obs_increment_enkf Computes increments for an ensemble Kalman filter with perturbed obs mean correction.
-
-%% 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$
-
-% Set error return to default successful
-err = 0;
-
-% Compute prior ensemble mean and variance
-prior_mean = mean(ensemble);
-prior_var = var(ensemble);
-
-% If both prior and observation error variance are zero return error
-if(prior_var <= 0 && obs_error_var <= 0)
- err = 1;
- return;
-end
-
-% Compute the posterior mean and variance
-% If prior variance is 0, posterior mean is prior_mean and variance is 0
-if(prior_var == 0)
- post_mean = prior_mean;
- post_var = 0;
-elseif(obs_error_var == 0)
-% If obs_error_var is 0, posterior mean is observation and variance is 0
- post_mean = observation;
- post_var = 0;
-else
-% Use product of gaussians
- % Compute the posterior variance
- post_var = 1 / (1 / prior_var + 1 / obs_error_var);
-
- % Compute posterior mean
- post_mean = post_var * (prior_mean / prior_var + observation / obs_error_var);
-end
-
-% Generate the perturbed observations by adding
-% draw from Normal(0, obs_error_sd)
-temp_obs = observation + sqrt(obs_error_var) * randn(size(ensemble));
-
-% Adjust so that perturbed observations have mean = to observation
-% This is technically an enhancement of earliest EnKFs
-temp_obs = temp_obs - mean(temp_obs) + observation;
-
-% Compute new ensemble members by taking product of prior ensemble
-% members and perturbed obs pairs
-updated_ens = post_var * (ensemble / prior_var + temp_obs / obs_error_var);
-
-% Increments are difference between updated and original ensemble
-obs_increments = updated_ens - ensemble;
-
-% Following are enhancements that change characteristics of EnKF
-% Coding is not finalized so do NOT just uncomment
-% Could also adjust the mean and variance of the final sample;
-% This can greatly improve and change the behavior
-% Shift the prior ensemble to have the posterior mean
-%%%updated_ensemble = ensemble - prior_mean + post_mean;
-
-% Contract the ensemble to have the posterior_variance
-%%%var_ratio = post_var / prior_var;
-%%%updated_ensemble = sqrt(var_ratio) * (updated_ensemble - post_mean) + post_mean;
-
-% Compute the increments
-%%%obs_increments = updated_ensemble - ensemble;
-
-% <next few lines under version control, do not edit>
-% $URL$
-% $Revision$
-% $Date$
-
Deleted: DART/trunk/DART_LAB/matlab/obs_increment_rhf.m
===================================================================
--- DART/trunk/DART_LAB/matlab/obs_increment_rhf.m 2015-06-24 21:20:35 UTC (rev 8125)
+++ DART/trunk/DART_LAB/matlab/obs_increment_rhf.m 2015-06-24 21:21:54 UTC (rev 8126)
@@ -1,185 +0,0 @@
-function [obs_increments, err] = obs_increment_rhf(ensemble, observation, obs_error_var)
-%% obs_increment_rhf Computes increments for a rank histogram filter
-% Need to discuss the available options eventually
-% For now this implements the default options
-
-%% 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$
-
-% Set error return to default successful
-err = 0;
-
-% Get the ensemble size
-ens_size = size(ensemble, 2);
-prior_sd = std(ensemble);
-prior_var = prior_sd^2;
-
-% allocate space for some of the density calculation bits
-like_dense = zeros(1,ens_size); % likelihood density
-mass = zeros(1,ens_size);
-height = zeros(1,ens_size);
-
-% Sort the ensemble members and keep the indices
-[x, e_ind] = sort(ensemble);
-
-% Compute the likelihood of each member given the observation
-like = exp(-1 * (x - observation).^2 / (2 * obs_error_var));
-
-% Compute the mean likelihood density in each interior bin
-for i = 2:ens_size
- like_dense(i) = (like(i - 1) + like(i)) / 2;
-end
-
-% For unit normal, find distance from mean to where cdf is 1/(n+1)
-dist_for_unit_sd = -1 * weighted_norm_inv(1, 0, 1, 1/(ens_size + 1));
-
-% Variance of tails is just sample prior variance
-% Mean is adjusted so that 1/(ens_size + 1) is outside
-left_mean = x(1) + dist_for_unit_sd * prior_sd;
-left_var = prior_var;
-left_sd = prior_sd;
-
-% Same for the right tail
-right_mean = x(ens_size) - dist_for_unit_sd * prior_sd;
-right_var = prior_var;
-right_sd = prior_sd;
-
-% Eventually want to support options, for now
-gaussian_likelihood_tails = false;
-
-if(gaussian_likelihood_tails)
- % Need to fill this in
-else
- % Block to do flat tails for likelihood follows
- % This removes assumptions about likelihood and cuts cost
-% new_var_left = left_var;
- new_sd_left = left_sd;
- new_mean_left = left_mean;
- prod_weight_left = like(1);
- mass(1) = like(1) / (ens_size + 1);
-
- % Same for right tail
-% new_var_right = right_var;
- new_sd_right = right_sd;
- new_mean_right = right_mean;
- prod_weight_right = like(ens_size);
- mass(ens_size + 1) = like(ens_size) / (ens_size + 1);
- % End block for flat tail likelihood
-end
-
-% The mass in each interior box is the height times the width
-% The height of the likelihood is like_dense
-% For the prior, mass is 1 / ((n+1) width), and mass = height x width so
-% The height of the prior is 1 / ((n+1) width); Multiplying by width leaves 1/(n+1)
-
-% In prior, have 1/(n+1) mass in each bin, multiply by mean likelihood
-% density to get approximate mass in updated bin
-
-for i = 2:ens_size
- mass(i) = like_dense(i) / (ens_size + 1);
- % Height of prior in this bin is mass/width; Only needed for trapezoidal
- % If two ensemble members are the same, set height to -1 as flag
- if(x(i) == x(i-1))
- height(i) = -1;
- else
- height(i) = 1 / ((ens_size + 1) * (x(i) - x(i-1)));
- end
-end
-
-% Now normalize the mass in the different bins to get a pdf
-mass_sum = sum(mass);
-nmass = mass / mass_sum;
-
-% Get the weight for the final normalized tail gaussians
-% This is the same as left_amp=(ens_size + 1)*nmass(1)
-left_amp = prod_weight_left / mass_sum;
-% This is the same as right_amp=(ens_size + 1)*nmass(ens_size+1)
-right_amp = prod_weight_right / mass_sum;
-
-% Find cumulative mass at each box boundary and middle boundary
-cumul_mass = zeros(1,ens_size+1);
-for i = 1:ens_size + 1
- cumul_mass(i+1) = cumul_mass(i) + nmass(i);
-end
-
-% Begin internal box search at bottom of lowest box
-lowest_box = 1;
-
-new_ens = zeros(1,ens_size);
-obs_increments = zeros(1,ens_size);
-% Find each new ensemble member's location
-for i = 1:ens_size
- % Each update ensemble member has 1/(n+1) mass before it
- umass = i / (ens_size + 1);
-
- % If it is in the inner or outer range have to use normal
- if(umass < cumul_mass(2))
- % It is in the left tail
- % Get position of x in weighted gaussian where the cdf has value umass
- new_ens(i) = weighted_norm_inv(left_amp, new_mean_left, ...
- new_sd_left, umass);
- elseif (umass > cumul_mass(ens_size + 1))
- % It's in the right tail
- % Get the position of x in weighted gaussian where the cdf has value umass
- new_ens(i) = weighted_norm_inv(right_amp, new_mean_right, ...
- new_sd_right, 1 - umass);
- % Coming in from the right, use symmetry after pretending it's on left
- new_ens(i) = new_mean_right + (new_mean_right - new_ens(i));
- else
- % In one of the inner boxes
- for j = lowest_box:ens_size - 1
- % Find the box that this mass is in
- if(umass >= cumul_mass(j+1) && umass <= cumul_mass(j+2))
-
- % Only rectangular is implemented for now
- rectangular_quadrature = true;
-
- if(rectangular_quadrature)
- % Rectangular quadrature block first
- % Linearly interpolate in mass
- new_ens(i) = x(j) + ((umass - cumul_mass(j+1)) / ...
- (cumul_mass(j+2) - cumul_mass(j+1))) * (x(j+1) - x(j));
- else
- % Trapezoidal interpolation block goes here
- end
-
- % Don't need to search lower boxes again
- lowest_box = j;
- break
- end
- end
- end
-end
-
-% Convert to increments for unsorted
-for i = 1:ens_size
- obs_increments(e_ind(i)) = new_ens(i) - x(i);
-end
-
-
-%-----------------------------------------------
-
-
-function [x] = weighted_norm_inv(alpha, mean, sd, p)
-
-% Find the value of x for which the cdf of a N(mean, sd)
-% multiplied times alpha has value p.
-
-% Can search in a standard normal, then multiply by sd at end and add mean
-% Divide p by alpha to get the right place for weighted normal
-np = p / alpha;
-
-% Find spot in standard normal
-x = norm_inv(np);
-
-% Add in the mean and normalize by sd
-x = mean + x * sd;
-
-% <next few lines under version control, do not edit>
-% $URL$
-% $Revision$
-% $Date$
-
Deleted: DART/trunk/DART_LAB/matlab/plot_gaussian.m
===================================================================
--- DART/trunk/DART_LAB/matlab/plot_gaussian.m 2015-06-24 21:20:35 UTC (rev 8125)
+++ DART/trunk/DART_LAB/matlab/plot_gaussian.m 2015-06-24 21:21:54 UTC (rev 8126)
@@ -1,26 +0,0 @@
-function[plot_handle] = plot_gaussian(mymean, sd, weight)
-%% plot_gaussian Plot gaussian over 5 standard deviations
-
-%% 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$
-
-x_min = mymean - 5*sd;
-x_max = mymean + 5*sd;
-x_range = x_max - x_min;
-
-% Number of points is 1001
-num_points = 1001;
-interval = x_range / num_points;
-x = x_min:interval:x_max;
-y = weight * norm_pdf(x, mymean, sd);
-
-plot_handle = plot(x, y);
-
-% <next few lines under version control, do not edit>
-% $URL$
-% $Revision$
-% $Date$
-
Deleted: DART/trunk/DART_LAB/matlab/plot_polar.m
===================================================================
--- DART/trunk/DART_LAB/matlab/plot_polar.m 2015-06-24 21:20:35 UTC (rev 8125)
+++ DART/trunk/DART_LAB/matlab/plot_polar.m 2015-06-24 21:21:54 UTC (rev 8126)
@@ -1,21 +0,0 @@
-function h = plot_polar(y, x, mean_dist, string, model_size)
-%% plot_polar
-
-%% 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$
-
-% Y includes a wraparound point, x does not
-x_t(model_size + 1) = x(1);
-x_t(1:model_size) = x;
-h = polar(y, mean_dist + x_t, string);
-
-end
-
-% <next few lines under version control, do not edit>
-% $URL$
-% $Revision$
-% $Date$
-
Deleted: DART/trunk/DART_LAB/matlab/product_of_gaussians.m
===================================================================
--- DART/trunk/DART_LAB/matlab/product_of_gaussians.m 2015-06-24 21:20:35 UTC (rev 8125)
+++ DART/trunk/DART_LAB/matlab/product_of_gaussians.m 2015-06-24 21:21:54 UTC (rev 8126)
@@ -1,30 +0,0 @@
-function [post_mean, post_sd, weight] = ...
- product_of_gaussians(prior_mean, prior_sd, obs, obs_err_sd)
-%% product_of_gaussians Computes mean, variance and weight of the product of two unit gaussians given the mean and standard deviation of each.
-
-%% 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$
-
-% Get the prior and observational error variance
-prior_var = prior_sd^2;
-obs_err_var = obs_err_sd^2;
-
-% Compute the posterior variance
-post_var = 1. / (1. / prior_var + 1. / obs_err_var);
-post_sd = sqrt(post_var);
-
-% Compute the posterior mean
-post_mean = post_var * (prior_mean / prior_var + obs / obs_err_var);
-
-% Compute the associated weight
-weight = (1. / (sqrt(2. * pi) * sqrt(prior_var + obs_err_var))) *...
- exp(-0.5 * (obs - prior_mean).^2 ./ (prior_var + obs_err_var));
-
-% <next few lines under version control, do not edit>
-% $URL$
-% $Revision$
-% $Date$
-
More information about the Dart-dev
mailing list