# Presentation: 611184

## Chemcial Data Assimilation for Air Quality ForecastingMethods used in Europe

**H. Elbern**
- with support from
**1****A. Strunk, L. Niradzik, E. Friese, Z. Milbers, **
**2****M. Schrödter-Homscheidt**
*1**Rhenish Institute for Environmental Research *
*at the University of Cologne (RIU)*
- and
*1**Helmholtz virt. Inst. for Inverse Modelling of Atmospheric Chemical Composition (IMACCO)*
*2**German Data Centre for Remote sensing (DLR-DFD)*

## Contents:

- Survey and scope
- Air Quality data assimilation in routine operations today
- Advanced scientific applications
- Kalman filter
- 4dimensional variational data assimilation (4Dvar)
- Example applications
- 4.	Future research directions
- Short term (1-3 years)
- Long term (> 3 years)

## Survey (1) Regional operational data assimilation

- to comply with EU legislation
- core assimilation data of in situ type from observation networks operated by national or regional/state authorities
- Air Quality considered as a high resolution problem  limited area and nested models
- 2 major European efforts in Global Environmental Monitoring and Security (GEMS = European GEOSS)
- ESA funded PROMOTE (_[www.gse-promote.org](http://www.gse-promote.org/)_)
- EU funded GEMS (new!)

**1 Survey an Scope**

## European metropolitan chemograms

- ozone

- NO2

- PM10

- SO2

- CO

- benzene

**1 Survey an Scope**

## Prevair system in PROMOTE

- Method:	Kriging
- Reference: 	Blond et Vautard, JGR, 2003, 2004

- analysis and observations	17.6.05 	forecast and observations

**2 Routine operations**

## EURAD system in PROMOTE

- Method: 	2 and 3Dvar
- Reference:	Reduced version of Elbern and Schmidt, JGR, 2001

- analysis and observations	 28.5.05 	forecast and observations

**2 Routine operations**

**Notes:**
28.5. ozone forecast left ana right

## Survey (2)Global tropospheric partly operational

- MOCAGE (Meteo France) 3DFGAT, (= 3Dvar with “First Guess at Appropriate Time”)
- TM3-5 (KNMI) “Kalman filter” with parameterize covariance update,
- TM3-BETHY (MPI-BGC)  soil column variational (1 spatial D + time D)

**2 Routine operations**

## Special challenges of tropospheric chemistry data assimilation

- The following problems prevail:
- strong influence of manifold processes including *emissions* and *deposition*
*spatially highly variable* “chemical regimes”
- chemical state observability (= “analyseability”) hampered by *manifold hydrocarbon* species
- consistency  with *heterogeneous* data sources: satellite data and in situ observations

**3. Advanced Mthods**
**General remarks**

- This invokes the application space-time data assimilation algorithms
- preserving the **BLUE** property (**B**est **L**inear **U**nbiased **E**stimator)

## Advanced spatio-temporal data assimilation:4D dimensional variational and Kalman filtering

- Model not to  passively accept external 3D analyses (probably engendering spurious relaxations)
- Rather contribute with its dynamics/ chemical kinetics as constraint to **estimate the most probable state** or parameter values

- Best Linear Unbiased Estimate (BLUE)
- while using data ***and*** models consistently
- allows for **Hypothesis testing**

**3. Advanced Mthods**
**General remarks**

## Advanced spatio-temporal methods used in tropospheric chemistry data assimilation 

- 4-dim variational data assimilation (**4D-VAR**, and versions: 4D physical space statistical analysis system (Amodei, 1995; Courtier, 1997); 4D-var incremental (Courtier et al, 1994))
- reduced complexity versions of the Kalman Filter: Reduced Rank Square root KF (**RRSQKF**, Verlaan and Heemink, 1997) , ensemble KF (Evensen, 1994), SEEK, SEIK, (Verron et al., 1999)
**Spacio-temporal BLUEs applied in tropospheric chemistry data assimilation**:
- 4D var:
- with EURAD (Elbern and Schmidt, 1999, 2001),
- with POLAIR (Issartel and Baverel, 2003)
- RRSQKF
- with LOTOS model (van Loon et al, 2000),
- with EUROS model (Hanea et al. 2004)
- _Remark:_* 3D BLUE algorithm analyses like those from Optimal Interpolation, once ingested into a model, do not result in a 4D BLUE analysis*
*	** **For example 3D-FGAT is not!*

**3. Advanced Mthods**
**General remarks**

## Hypothesis: initial state and emission rates are least known

**3. Advanced Mthods**
**General remarks**

##  Kalman filter, basic equations

**3. Advanced Mthods**
**Lakman filter**

## Reduced rank Kalman filter (basic idea)

- Appoximate covariance matrices
- Pb,a (n x n) by a product of suitably low ranked matrix Sb,a (n x p),  and p << n . Same procedure for the system noise matrix Q with T (n x r).
- The forecast step remains unchanged.
- The forecast error covariance matrix rests on 2 x p model integrations only!
- The enlargement of p by r enforces periodic reductions of columns  (with lowest ranked eigenvalues)

**3. Advanced Mthods**
**Lakman filter**

## Reduced rank Kalman filter (calculus)

**Lakman filter**

- In practice, all calculations can be performed without actually calculating matrices P!
- Positive semi-definitenes is maintained!

## Air quality Kalman filter implementationsand comlexity reduction strategy applied

- LOTOS model (van Loon et al., 2000)
- reduced rank square root
- EUROS model (Hanea et al., 2004)
- ensemble
- reduce rank square root

**3. Advanced Mthods**
**Lakman filter**

- Methods used with optimisation of various parameters
- including emissions

## Transport-diffusion-reaction equation and its adjoint

- _Transport-diffusion-reaction equation and its adjoint_

**3. Advanced Mthods**
**4D var**

## In the troposphere for emission rates the product (paucity of knowledge*importance)is high 

**4D var**

**3. Advanced Mthods**
**Lakman filter**

**3. Advanced Mthods**
**Lakman filter**

- In the troposphere for **emission rates** the product (*paucity of knowledge*importance)*is high

## Adjoint integration “backward in time”(slide from lecture 1)

- How to make the
- parameters of resolvents *i*
**M**(ti-1,ti) available in ***reverse***
- order??

## The 4D variational method:Key development: construction of the  adjoint code

- backward model
- (backward differential equation)

- adjoint algorithm
- (adjoint solver)

- adjoint code

**3. Advanced Mthods**
**4D var**

- forward model
- (forward differential equation)

- algorithm
- (solver)

- code

## Source identificationby 4DvarPOLAIR model Issartel and BaverelEcole National de Ponts et Chausees (ENPC)

**3. Advanced Mthods**
**4D var application**

- Calculations of optimal retroplumes
- (Issartel and Baverel, ACP, 2003)

**Notes:**
Fig. 2. Minimum release from a point source a ground or sea level (a)Qmin(x) compatible with the peak measurement of Freiburg
μ1 = 103mBq.m-3 (b)Qmin
4 (x) compatible with the series of measurements in Freiburg μ0 = 6, μ1 = 103, μ2 = 42, μ3 = 6 mBq.m-3;
this value of Qmin
4 (x) would not be altered by considering a fictitious measurement in Stockholm μS = 0mBq.m-3 on 3 February 2000
(c)Qmin
5 (x) compatible with the previous series completed with a fictitious measurement μS = 100mBq.m-3 in Stockholm.

## Mesoscale EURAD 4D-var data assimilation system

**4D var application**

- meteorological
- driver
-  **MM5**

- EURAD
- emission model
**EEM**
**emission 1. guess**

** ****direct CTM**

- emission rates

- Initial
- values

- minimisation

**adjoint**
**CTM**

- observations

- analysis

- gradient

- fore-
- cast

## Computational complexity estimate of the variational algorithm

| Storage strategy | # forward runs/iteration | storage | Complexity
[Tforw] |
| --- | --- | --- | --- |
| total storage | 1 | const*Nx*Ny*Nz*Nc*NT*No
O(1012-1013) | 3 |
| operatorwise
1 level | 2 | Nx*Ny*Nz*Nc*NT*No
O(108-109) | 4 |
| dynamic stepwise
2 levels | 3 | Nx*Ny*Nz*Nc*NT
O(107-108) | 5 |

- Nx*Ny*Nz	 spatial dimensions O(104-105)
- Nc		# constituents O(100)
- NT		# time steps of assimilation window O(10-100)
- No		# operators O(10)
- const		intermediate results O(104)

**3. Advanced Mthods**
**4D var application**

## 2 level forward and backward integration scheme

- time direction

- disk

- forward  integration

- backward integration

**intermediate storage**
**Level 1**

**3. Advanced Mthods**
**4D var application**

## storage sequence: level 2, operator split(each time step)

- storage sequence: level 2, operator split*(each time step)*

**intermediate storage**
**Level 2**

- horiz. Advekt x

- horiz. Advekt y

- Chemie

- vert. Diffusion

- vert.. Advekt

- horiz. Advekt x

- horiz. Advekt y

- vert. Diffusion

- vert.. Advekt

- horiz. Advekt x

- horiz. Advekt y

- vert. Diffusion

- vert.. Advekt

- Adjungierte
- horiz. Advekt x

- Adjungierte
- horiz. Advekt y

- Adjungierte
- vert. Advekt

- Adjungierte
- vert.Diffusion

- horiz. Advekt x

- horiz. Advekt y

- vert. Diffusion

- vert.. Advekt

- Adjungierte
- horiz. Advekt x

- Adjungierte
- horiz. Advekt y

- Adjungierte
- vert. Advekt

- Adjungierte
- vert.Diffusion

- Chemie

- Adjungierte
- Chemie

**intermediate storage**
**Level 1**

**3. Advanced Mthods**
**4D var application**

## Computational resources requested

- with a  14 h assimilation interval about 18 iterations requested
- results in 12  CPU-hours with 121 processors of a T3E

**3. Advanced Mthods**
**4D var application**

## Treatment of the inverse problem to infer emission rates

## Normalised diurnal cycle of anthropogenic surface emissions f(t)emission(t)=f(t;location,species,day) * v(location,species)day in {working day, Saturday, Sunday}	v optimization parameter

**3. Advanced Mthods**
**4D var application**

## Semi-rural measurement site Eggegebirge 

- Semi-rural measurement site **Eggegebirge**

- assimilation interval

- forecast

**+ observations**
**no optimisation**

**initial value opt****.**

**emis. rate opt.**

**joint emis + **
**ini val opt****.**

**3. Advanced Mthods**
**4D var application**

## error statisticsbias (top), root mean square (bottom)

**3. Advanced Mthods**
**4D var application**

## Which is the requested resolution?BERLIOZ grid designs and observational sites  (20.21. 07.1998)

- x=54 km

**Which is the requested resolution?**BERLIOZ grid designs and observational sites  (20.21. 07.1998)

- x=18 km

- x=6 km

- x=2 km

- Control and diagnostics

## Some BERLIOZ examples of NOx assimilation 

**4D var application**

- Some BERLIOZ examples of NOx assimilation

-  Time series for selected NOx stations (upper panel NO, lower panel NO2) on nest 2.
**+** observations,
**- -** no assimilation,**____** N1 assimilation (**18 km**), **____**N2 assimilation.(**6 km**),
- grey shading: assimilated observations, others forecasted.

- NO

- NO2

## Emission source estimates by inverse modellingOptimised emission factors for Nest 3

**4D var application**

- Emission source estimates by inverse modellingOptimised emission factors for Nest 3

- surface

- height layer ~32-~70m

## Status: EURAD assimilation system prepared for 4D-var assimilation production run 1.7.31.12.2003

- Gas phase satellite observations
- KNMI - NO2 tropospheric columns, SCIAMACHY averaging kernel assimilation
- NNORSY neural network based ozone tropospheric columns (revised version based on GOME)
- IMK - MIPAS upper troposphere
- Aerosol phase satellite observations
- DFD - PM10 bulk retrievals SCAMACHY-AATSR

**3. Advanced Mthods**
**4D var application**

## NO2 tropospheric column assimilation IFE:Example: BERLIOZ episode 20.7.1998

**4D var application**

- NO2 tropospheric column assimilation IFE:Example: BERLIOZ episode 20.7.1998

- GOME retieval
- Courtesy: A. Richter, IFE Bremen

- no assimilation

- with assimilation

**NOAA 14**
**11:51 UTC**

**Notes:**
Satellite NOAA 14 Direction Northbound Equator crossing angle 312.898 degrees West Orbit number 18314 Identification 2213/08 (deprecated) Number of lines in pass 4681 Average altitude 856.001 km Time of first line 1998-07-20 11:50:57.812 Time of last line 1998-07-20 12:03:57.767 Duration 13.00278 minutes (theoretical) Equator crossing time 1998-07-20 11:38:49.612 Horizon time 1998-07-20 11:51:00.106 Overhead time 1998-07-20 11:57:43.362 Delay 12.17490 minutes from equator crossing (theoretical) Format HRPT-Packed Number of bad lines 45 CD identification 98 07 20 06 
Satellite NOAA 14 Direction Northbound Equator crossing angle 312.898 degrees West Orbit number 18314 Identification 2213/08 (deprecated) Number of lines in pass 4681 Average altitude 856.001 km Time of first line 1998-07-20 11:50:57.812 Time of last line 1998-07-20 12:03:57.767 Duration 13.00278 minutes (theoretical) Equator crossing time 1998-07-20 11:38:49.612 Horizon time 1998-07-20 11:51:00.106 Overhead time 1998-07-20 11:57:43.362 Delay 12.17490 minutes from equator crossing (theoretical) Format HRPT-Packed Number of bad lines 45 CD identification 98 07 20 06 

## IUP Bremen GOMEGOME forecast validationmodel/retrieval ratio for BERLIOZ 20. (assimilated) + 21.(forecasted)  7.1998

**4D var application**

- IUP Bremen GOMEGOME forecast validation_**model/retrieval ratio**_ for BERLIOZ 20. (assimilated) + 21.(forecasted)  7.1998

- ideal

- no assimilation

- with assimilation

-  20.7.98

- verification
- forecast
- 21.7.98
- (not
- assimilated)

- forecast

- forecast

## KNMI GOMEGOME forecast validationmodel/retrieval ratio for BERLIOZ 20. (assimilated) + 21.(forecasted)  7.1998

- no assimilation

- with assimilation

- ideal

-  20.7.98

- verification
- forecast
- 21.7.98
- (not
- assimilated)

- forecast

- forecast

**3. Advanced Mthods**
**4D var application**

## GOME ozone profile assimilation  BERLIOZ (20.7.1998)Data: Neuronal Network retrieval (Müller et al., 2003)

**4D var application**

- GOME ozone profile assimilation  BERLIOZ (20.7.1998)Data: Neuronal Network retrieval (Müller et al., 2003)

## Ozone observation increment 20.07.1998

**4D var application**

- layer:			200-100 hPa
- radius of influence: 	200 km

- height dependent influence radii

## Future Research directions:short term

- operationalisation of advanced methods: use of complexity constraining algorithms
- aerosol data assimilation
- Model Output Statistics (MOS)

**4. Future Methods**
**short term**

## 4. Future Methods

**short term**

- (NH4)2SO4		watersoluble
- NH4NO3		watersoluble
- H2SO4		watersoluble
- HNO3		watersoluble
- H2O		watersoluble
- prim. organic	insoluble
- Prod. Aromate	watersoluble
- Prod. Alkane	watersoluble
- Prod. Alkene	watersoluble
- Prod. -pinene	watersoluble
- Prod. limonene	watersoluble
- elem. carbon	SOOT
- prim. small aerosols	insoluble
- maritime		SEASALT
- mineral		MINERAL TRANSPORTED
- other anthropogenic	insoluble

**MADE****	****	****SYNAER**

- INSO

- The attribution EURAD-MADE-SORGAM model variables to SYNAER aerosol types

- joint
- SCIAMACHY-AATSR
- retrieval
- (Holzer-Popp et al. JGR, 2003)

## 4. Future Methods

**short term**

- MADE, Level 1, SO4, accumulation mode, mass.-conc. g/m3example 5. August 1997

- Each  MADE-variable modified according to the minimisation result,
-  relative structure of the aerosol vertical profile remains unaffected

- background mass conzentration	      analysis increment           	      new mass conzentration

- 0     10    20   30    40    50

- 0     10    20   30    40    50

- -2    0    2    4    6    8    10

## Future Research directions:long term

- Multi scale assimilation (for error of representativity reduction)  dynamical nesting
- Systems  with non-**Gaußian** error characteristics
- Systems with **phase transitions **(microphysics)
- advanced a posteriori-evaluations

## Necessary a posteriori condition for general consistency test (2 test)

- necessary condition with p observations:
- minimum should be half the number of observations:

- Necessary a posteriori condition for general consistency test (2 test)

- With

**3. Advanced Mthods**
**General remarks**

## Objective function as a sum of different data sources

- _Partition to any data source possible?_
- (Ivanov and Deque, 2001)

- mj number of observations from information source/satellite j
**S**j  source j error covariance matrix
**P**a analysis error covariance matrix

## The EndThank you for attention