Version 0.0 Monika KRYSTA NCAR/MMM 01-2009

Version 0.1 Yann MICHEL NCAR/MMM 02-2009

Version 0.2 Yann MICHEL NCAR/MMM 03-2009 – update of diagnostics

Stage 0 

gen_be_stage0_wrf.f90 build up "standard perturbation fields"

Purpose : To convert WRF output to "standard perturbation fields" required by the WRF-Var BE covariance generation utility "gen_be".

Input  : WRF forecasts for a specified date (NETCDF format). 

Output : Binary files for use in gen_be_stage1.

[1] Initialize information

• Reads arguments, initializes

[2] Setup up ancillary fields using 1st member values

• Reads fields (lat, map factors). Initializes FFT.

[3] Convert WRF forecast fields to standard fields and output

• Loop over members
• Loop over vertical dimension
• Reads u,v. Compute vorticity/divergence. 
• Calculate streamfunction and potential 
• Read mass fields, and convert to T and r_h.
• Finally, extract surface pressure
• Write out ensemble forecasts for this member

[4] Compute perturbations and output [NMC or Ensemble]

• Take forecast/ensemble difference
• Write out perturbations
• Write out mean/mean square fields

Stage 1  

gen_be_stage1.f90 Produces temporal perturbations

Purpose : to remove the binned mean from the perturbation fields.

Input   : binary files: "pert.ccyymmddhh.e"ce for ENS or "pert.ccyymmddhh.e001" for NMC.

Output : binary files for use of the gen_be_stage2

 [1] Initialize namelist variables and other scalars.

• Reads namelist 'gen_be_stage1_nl.nl' (data dir, bin method...)

[2] Read fields from standard files, and calculate mean fields'

• Loop over dates
• Read perturbation (psi_mean, chi_mean, t_mean, rh_mean, ps_mean...)
• Create bin (if first loop)
• Calculate time/ensemble mean (psi_mean, chi_mean, t_mean, rh_mean, ps_mean)
• Loop over dates (second time =Read fields again, and remove time/ensemble/area) 
• Remove mean
• Write necessary full-fields (ni, nj, nk, latitude, height)
• Write fields to file for further processing
• Write down bin info into 'bin.data'

Stage 2  

gen_be_stage2.f90 Computes regression coefficients: 'covar/var'

Input: Reads 'fullflds/yyymmddhh.fullflds.e+ce'

Reads 'bin.data'

Reads perturbation files 'varyyyymmddhh.var.e+ce'

Output: Write regression coefficients into 'gen_be_stage2.dat '

[1] Initialize namelist variables and other scalars.'

• Read namelist 'gen_be_stage2_nl.nl'

[2] Read fields, and calculate correlations

• Loop over dates
• Loop over members
• Read Full-fields (latitude, height)
• Read bin info (if first time)
• Read perturbations (psi, chi...)
• Calculate psi/chi covariances 3D bins (binning in all three directions: i,j,k)
• (covar1=chi*psi and var1=psi*psi) - depends on b (bin number)
• Calculate psi/ps covariances 2D bins
• (covar2=ps*psi) - depends on b and k (vertical level of psi)
• Calculate  psi/T covariances 2D bins
• (covar3=T*psi) - depends on b and k,k2 (vertical levels of T and psi)
• Calculate  psi/rh covariances 2D bins
• (covar4=rh*T) - depends on b and k,k2 (vertical levels of rh and psi)
• Calculate  psi/psi covariances 2D bins
• (var2=psi*psi) - depends on b and k,k2 (vertical levels of psi and psi)
• Fill in psi/psi covariance by symmetry

[3] Calculate eigenvectors, eigenvalues and inverse for psi/psi covariance

• Loop over bins
• Call da_eof_decomposition = performs eigenvalue/eigenvector decomposition of var2 matrix.
• Truncate eigenvalues to ensure inverse is not dominated by rounding error (eigenvectors are orthogonal)       hard-coded threshold on variance to make it robust
• Inverts 'truncated' var2  <psi,psi>^-1 = E . Lam^-1 . E^T:

[4] Calculate regression coefficients

• psi/chi 3D bins regcoeff1(b) = covar1(b) / var1(b)
• psi/ps 2D bins regcoeff2(k,b) = sum[covar2(k2,b)*var2_inv(k2,k,b)]
• psi/T 2D bins regcoeff3(k,k2,b) =sum [covar3(k,k3,b)*var2_inv(k3,k2,b)] 
• Output regression coefficients for use in 3/4D-Var into gen_be_stage2.dat

Stage 2a  

gen_be_stage2a.f90 Builds unbalanced fields: chi_u, T_u, ps_u

Purpose

Calculates (for every time instant and every member) unbalanced fields (as perturbed field- field explained by perturbed psi via a linear mapping): 

perturbed field - regcoeff*psi (summation over k performed in case of each field)

Input

Reads gen_be stage2.dat [have not checked the filtering yet]

Reads bin.data

Reads perturbation files var/yyyymmddhh.var.e

Output

[1] Initialize namelist variables and other scalars.

• Namelist gen_be_stage2a_nl.nl

[2] Read regression coefficients and bin information

• Read regression coefficients from 'gen_be_stage2.dat'

• Read bin info 'bin.data'

[3] Apply ', num_passes, ' pass recursive filter to regression coefficients  if ( num_passes > 0 )

• Call da_filter_regcoeffs

[4] Read standard fields, and compute control variable fields

• Loop over dates
• Loop over members
• Read psi predictor
• Calculate/write unbalanced chi
• Calculate/write unbalanced T
• Calculate/write unbalanced ps

Writes var_u/yyymmddhh.var_u.e+ce (.01 for ps_u)

Stage 3  

gen_be_stage3.f90 Reorganise distribution in vertical direction

Input

Read read fullflds/yyyymmddhh.e+ce

Read var/yyyymmdd.var.e+ce (for all fields; take unbalanced if available; +.01 for ps)

Output

Write var/yyyymmddhh.var.e+ce.ck  

(var=psi,chi_u,t_u,ps_u,rh and ck is vertical level number)

[1] Initialize namelist variables and other scalars

• namelist gen_be_stage3_nl.nl

 if ( .not. data_on_levels ) then

[2] Process fields for variables - loop over variables is in script (gen_be.ksh)

• Loop over dates
• Loop over members
• Create and sort into bins Would be enough to read bin.dat?
• Compute mean field 
• seems mean = average over grid points, not physical cosine weighted mean
• Build 2D perturbations - for each level subtract its spatial mean.
• Fill in upper-right part of BE matrix by symmetry

[3] Calculate vertical eigenvectors and eigenvalues for all variables

• Latitudinally varying BE decomposition  = separately for each 2D bin (local)
• Call da_eof_decomposition 
• Domain-averaged BE decomposition =  mean over all horizontal bins
• Call da_eof_decomposition 
• Output eigenvectors, eigenvalues for use in WRF_Var
• Decide on local or domain-averaged EOFs for horizontal correlations  use_global_eofs

Map binned eigenvectors to x, y grid

Take sqrt of eigenvalues

end (not. data_on_levels)

[4] Transform perturbations (or not), and output.

• Loop over dates
• Loop over members
• Project fields onto vertical modes  call da_transform_vptovv
• Output fields (split into 2D files to allow parallel horizontal treatment)

Stage 4 Regional  

gen_be_stage4_regional.f90 Computes lengthscale for recursive filter

Input: Read rundir/varyyyymmddhh.var.e+ce.ck (ck is vertical level)

Output: Writes ck, coeff, 1/coeff in dir.varck/sl_print.var 

[1] Initialize namelist variables and other scalars.

• reads namelist gen_be_stage4_regional_nl.nl

[2] Input fields and calculate correlation as a function of distance between points.

• Loop over dates
• Loop over members
• if first time; build up new 1D grid
• call get_grid_info distance is in (i,j) space rather than on Earth
• Calculate spatial correlation
• call get_covariance distance is in (i,j) space rather than on Earth

[3] Compute fit of correlation to a straight line

• Calculate fit of input covariance data to Gaussian correlation function

• takes ln and fits a straight line
• call make_scale_length

Stage Diagnostics I  

gen_be_diags.f90 Gather data into file 'be.dat'

Input: 'bin.dat', 'gen_be_stage2.dat', gen_be_stage3.'+var+'.dat', var/'sl_print.'+var

Output: be.dat

[0] Initialize namelist variables and other scalars.

• reads namelist gen_be_diags_nl.nl'

[1] Gather regression coefficients.'

• call da_readwrite_be_stage2

[2] Gather vertical error eigenvectors, eigenvalues

• Loop over unbalanced variables 
• call da_readwrite_be_stage3

[3] Gather horizontal scale length / horizontal error power spectra

• Loop over unbalanced variables 
• call da_readwrite_be_stage4 no bin2d info

Stage Diagnostics II  

gen_be_diags_read.f90 Read file 'be.dat' and print out diagnostics

Input: 'be.dat'

Output: Fortran units

[0] Initialize namelist variables and other scalars.

• read namelist gen_be_diags_nl.nl'

[1] Gather regression coefficients.

• read regcoeff for balanced/physical transform
• call da_print_be_stats_p

[2] Gather vertical error eigenvectors, eigenvalues.

• read vertical eigenvectors/eigenvalues
• call da_print_be_stats_v

[3] Gather horizontal scale length / horizontal error power spectra

• read total power (global) or scale_length (regional)
• call da_print_be_stats_h_global / da_print_be_stats_h_regional

Stage Diagnostics IIIa  

gen_be_cov3d.f90 Covariances between 3D balanced and unbalanced fields

Input: var1/yyyymmddhh.var1.e+ce,  var2/yyyymmddhh.var2.e+ce 

Output: Fortran units

[1] Initialize namelist variables and other scalars.

• Read namelist gen_be_cov3d_nl.nl

[2] Read fields, and calculate correlations

• Loop over date
• Loop over member
• Calculate covariances

Stage Diagnostics IIIb

gen_be_cov2d.f90 Covariances between 2D balanced and unbalanced fields

Input: var1/yyyymmddhh.var1.e+ce,  var2/yyyymmddhh.var2.e+ce 

Output: Fortran units

[1] Initialize namelist variables and other scalars.

• Read namelist gen_be_cov2d_nl.nl

[2] Read fields, and calculate correlations

• Loop over date
• Loop over member
• Calculate covariances