ESA.enkf.etkf_cor
index
/home/lfeng/otool/ESA/enkf/etkf_cor.py

Solve Ensemble transform matrix by using  approach LU solver
 
Authors: L. Feng, Edinburgh University
History: v0.9, 2012.09.11
History: v0.95, 2013.02.15
 
Classes:
===============================================
1. etkf_cl:class for solving ETKF data assimilation equation

 
Modules
       
scipy.sparse.linalg.dsolve.linsolve
numpy.linalg
numpy
sys

 
Classes
       
etkf_cor_cl

 
class etkf_cor_cl
    class for Ensemble Transform Kalman Filter solver for correlated observations
 
See Feng et al., ACP 2009; Palmer et al., 2011.
 
 
Members:
-----------------------------------------------------
 
# control variables 
1. nx:<int>: size of control variable 
2. ne:<int>: size of perturbation ensemble 
3. x:<array, (nx, ne)>: perturbation ensemble 
4. x_mean:<array, (nx)>: mean value of control variable 
5. xnorm:<float>: normalized factor 
 
# model observation 
6.  y:<array, (nobs, ne)>: ensemble for model observation perturbation 
7.  ny:<int>: ny=nobs
8.  y_mean:<array, (nobs)>: model observations for self.mean_x
9.  sp_y: <array, (nobs, ne)>: sparse matrix of self.y
 
# observation errors 
 
9.  r:<array, (nobs, nobs)>: observation err covariance 
10. sp_r:<array, (nobs,)>: sparse matrix of r
 
11. sp_b:<array, (nobs, nobs)>: background error covariance in obsevation space 
    sp_b=sp_y*sp_y^T
12. sp_r:<array, (nobs, nobs)>: sp_r=sp_r+
 
 
 
Functions
---------------------------------------------------------------
1.  __init__: initialization
2.  reset: re-do initialization
3.  get_posterior: get posterior 
4.  get_tm: get transform matrix
6.  get_inc_m: get increment matrix 
7.  get_x_inc: get analysis increment
8.  get_y_inc: get analysis increment in observation space
9.  do_y_transform: transform maxtrix in y space
 
  Methods defined here:
__init__(self, x_mean, y_mean, x, y, r, xref=1, xnorm=1.0)
Initialization
 
Inputs:
---------------------------------------------------------
1. x_mean:<array, (nx)>: mean values of control variables x
2. Y_mean:<array, (ny)>: mean model observation values.  
3. x: <array, (nx, ne)>:  ensemble of x perturbations 
4. y: <array, (nobs, ne)>:  ensemble of model forecast perturbations. ny=nobs 
5. r: <array, (nobs, nobs)>: observation error
6. xref:<float>: value used to scale mean_x 
7. xnorm:<float>: value used to scale ensemble purbations
get_inc_m(self, yobs)
get the analysis increment matrix 
 
 
Inputs:
-------------------------------------------------
1. yobs:<array, (nobs)>: the observations
 
Returns:     
------------------------------------------------
1. inc_m :<array, (ne,)>:increase matrix 
 
 
Notes:
--------------------------------------------------------
1. inc_m=y^T*(R+Y*Y^t)^-0.5*dy
get_posterior(self, yobs)
calculate posterior from assimilation of observations 
Inputs:
-------------------------------------------------
1. yobs:<array, nobs>: observations
 
Returns:
---------------------------------------
1. postprior:<array, nx>: 
 
Notes:
-----------------------------------------
1. xa=xf+k*(yobs-f(x))=xf+dx*inc_m
get_tm(self)
the analysis increment on the variables  
 
Returns:
----------------------------------------------
1. tm :<array, (ne, ne)>: matrix for posterior ensemble: x(f)=x(a)*tm 
 
Notes: 
-----------------------------------
1. tm*tm=[1.0+(y^t*R^-1*y) ]
get_x_inc(self, inc_m, x_in=None, mean_x_in=None)
the analysis increment of control variables  
Inputs:
------------------------------------------
1. inc_m:<array, (ne)>: increment matrix 
2. x_in:<array, (nx, ne)>: ensemble perturbations
3. mean_x_in:<array, (nx)>: mean value
 
Returns:
---------------------------
1. xinc: <array, nx>: the x increments
get_y_inc(self, y_in, inc_m, mean_y=None)
the analysis increment on the variables  
 
Inputs:
--------------------------------------------
1. y_in:<array, (new_ny, ne)>: the 'new' observation perturbation 
2. inc_m :<array, (ne)>: increment matrix 
 
Returns:
-------------------------------------------------------
1. yinc: <array, new_ny>: y increments from xinc
get_y_transform(self, y_in, tm, mean_y=None)
multiply y_in by tm 
 
 
Inputs:
--------------------------------------------------------------
1. y_in:<array>, (new_ny, ne)>: 'new' observation perturbations 
2. tm  :<array, (ne, ne)>: Transformation matrix
 
Returns:
---------------------------------------------------------
1. yt:<array, (new_ny, ne)>: yt=y_in*tm
reset(self, x_mean, y_mean, x, y, r, xref=1.0, xnorm=1.0, do_debug=False)
Re-Initialization
 
Inputs:
---------------------------------------------------------
1. x_mean:<array, (nx)>: mean values of control variables x
2. Y_mean:<array, (ny)>: mean model observation values.  
3. x: <array, (nx, ne)>:  ensemble of x perturbations 
4. y: <array, (nobs, ne)>:  ensemble of model forecast perturbations. ny=nobs 
5. r: <array, (nobs, nobs)>: observation error
6. xref:<float>: value used to scale mean_x 
7. xnorm:<float>: value used to scale ensemble purbations