IFS Documentation front page

I Observations
II Assimilation
III Dynamics
IV Physics
V Ensemble
VI Technical
VII Waves

5.4 Variational quality control

The variational quality control, VarQC, has been described by Andersson and Järvinen (1999). It is a quality control mechanism which is incorporated within the variational analysis itself. A modification of the observation cost function to take into account the non-Gaussian nature of gross errors, has the effect of reducing the analysis weight given to data with large departures from the current iterand (or preliminary analysis). Data are not irrevocably rejected, but can regain influence on the analysis during later iterations if supported by surrounding data. VarQC is a type of buddy check, in that it rejects those data that have not been fitted by the preliminary analysis (i.e. the current state vector), often because it conflicts with surrounding data.

5.4.1 Description of the method

The method is based on Bayesian formalism. First, an a priori estimate of the probability of gross error is assigned to each datum, based on study of historical data. Then, at each iteration of the variational scheme, an a posteriori estimate of the probability of gross error is calculated (Ingleby and Lorenc, 1993), given the current value of the iterand (the preliminary analysis). VarQC modifies the gradient (of the observation cost function with respect to the observed quantity) by the factor (the QC-weight),which means that data which are almost certainly wrong () are given near-zero weight in the analysis. Data with a are considered `rejected' and are flagged accordingly, for the purpose of diagnostics and feedback statistics, etc.

The normal definition of a cost function is

(5.6)

where is the probability density function. Instead of the normal assumption of Gaussian statistics, we assume that the error distribution can be modelled as a sum of two parts: one Gaussian, representing correct data and one flat distribution, representing data with gross errors. We write:

(5.7)

where subscript refers to observation numer . and are the Gaussian and the flat distributions, respectively:

(5.8)

(5.9)

The flat distribution is defined over an interval which in Eq. (5.9) has been written as a multiple of the observation error standard deviation . Substituting Eqs. (5.7) to (5.9) into Eq. (5.6), we obtain after rearranging the terms, an expression for the QC-modified cost function and its gradient , in terms of the normal cost function

:

(5.10)

(5.11)

(5.12)

(5.13)

where

(5.14)

(5.15)

5.4.2 Implementation

The a priori information i.e. and are set during the screening, in the routine DEPART, and stored in the NCMFGC1 and NCMFGC2-words of the ODB. Default values are set in DEFRUN, and can be modified by the namelist namjo. VarQC can be switched on/off for each observation type and variable individually using LVARQC, or it can be switched off all together by setting the global switch LVARQCG=.F.. Since an as good as possible `preliminary analysis' is needed before VarQC starts, it is necessary to perform part of the minimization without VarQC, and then switch it on.

HJO computes according to Eq. (5.11) and the QC-weight, , according to Eq.(5.14).

The 3D-Var assimilation system (Chapter 3) can handle non-quadratic cost functions because it is using the M1QN3 minimization routine. This allows VarQC to be activated after NITERQC (40 by default) iterations during the minimization process. The 4D-Var assimilation system (Chapter 2) by default use a conjugate gradient minimization method that requires a strictly quadratic cost function. So in 4D-Var we use a quadratic formulation of VarQC, controlled by LQVARQC (default .T.).

In the non-quadratic formulation the QC-weight, , Eq.(5.14) is updated for each simulation based on the normalized departure values calculated from the latest model state. However, to obtain a strictly quadratic cost function the QC-weight is not allowed to change during the minimisation process, because the cost function shape would then change. Also, when LQVARQC = .T. it is not possible to activate VarQC in the middle of the minimisation process, i.e. NITERQC MUST be zero, because this would otherwise introduce a sudden jump in the cost function. Therefore, in 4D-Var assimilations VarQC is not applied during the first outer-loop but only for the second (and possibly subsequent) outer-loop iteration(s).

The quadratic VarQC implementation satisfies the quadratic constraint by calculating the QC-weight, , based on the high resolution trajectory fields and using this constant weight during the next minimisation step. The cost function value is increased/reduced for each simulation during the minimisation by a factor *(norm_dep_LR-norm_dep_HR), linearized in the vicinity of the high resolution cost function at norm_dep_HR. For the minimisation the most important input is the modification of the gradient by the weight . During the trajectory run the high resolution departure is stored in ROBODY(..,MDBIFC1(NUPTRA+1)) in HDEPART. For each simulation during the next minimisation step the high resolution departure is read and normalized by the final observation error in HJO)HJO loop 1.3.2. This is used as input to the VarQC weight calculations (loop 1.5.3-1.5.4 in HJO). It is assured that the same weight is used for each simulation in the minimisation if LLQVARQC_MIN =. T., because ZFJOS_HR(JOBS,JBODY) is constant during the minimisation. The cost function value is calculated like Jo_varqc_LR = Jo_varqc_HR+*(norm_dep_LR-norm_dep_HR)

The a posteriori probability of gross error is stored in the ODB and passed to BUFR feedback. Storing in ODB is done at the final simulation of the last minimisation if LQVARQC =. F., but done during the final trajectory run if LQVARQC=.T. in order to use the updated final trajectory for calculating the QC-weight based on final high-resolution analysis values. This is done in HJO section 1.6.

Variational quality control for time correlated observations: The same method as above is (by default) applied for time correlated observations (LTC=.T.). Here the high resolution departure values are required for all the time correlated observations. This is achieved by storing the normalized high resolution departures in the array RTCNDPHR in HJO) and reading them in COMTC. They are then copied to local array ZTCNDPHR_W and then used instead of the low resolution departures in the calculation (overwriting ZOM1DEP_W by high resolution values).