Transcript ROMS-ESPreSSO-Wilkin-201008

ROMS data assimilation for ESPreSSO
Accomplishments:
* Nested ROMS in larger domain forward simulation
(MABGOM-ROMS) with configuration suitable for IS4DVAR
experimentation. Considerations: boundary conditions,
resolution, computational cost
* IS4DVAR implemented in Slope Sea and MAB shelf waters,
assimilating SST and along-track altimeter sea level anomaly
(SLA). Considerations: tune IS4DVAR horizontal/vertical decorrelation scales, duration of assimilation window, data
preprocessing (error statistics, aliasing, mean dynamic
topography).
* Used withheld data to evaluate how well adjoint propagates
information between variables, and in space and time.
1
ROMS data assimilation for ESPreSSO
Accomplishments:
* Full IS4DVAR reanalysis of NJ inner/mid-shelf for LaTTE
using all data from CODAR, 2 gliders, moored current-meters
and T/S, towed SeaSoar CTD, and satellite SST
* Developed adjoint-based analysis methods for observing
system design and evaluation
* Have an ESPreSSO ROMS system ready for expansion to:
• 2006-2008 reanalysis of ocean physics
• introduction of in situ physical data into reanalysis
• analyze impact of improved physics on ecosystem model
• adjoint/tangent-linear simple optical model, with IS4DVAR
2
Mid-Atlantic Bight ROMS Model for ESPreSSO/IS4DVAR
5 km resolution IS4DVAR
model embedded in …
… ~12 km resolution outer model:
NCOM
global HyCOM/NCODA
ROMS MAB-GoM
3
Mid-Atlantic Bight ROMS
5 km resolution is for IS4DVAR
can use 1 km downscale for
forecast, with forward
ecosystem/optics
• 3-hour forecast meteorology
•
•
•
NCEP/NAM
daily river flow (USGS)
boundary tides (TPX0.7)
nested in ROMS MABGOM
V6 (nested in GlobalHyCOM*) (* which assimilates altimetry)
– nudging in a 30 km boundary
zone
– radiation of barotropic mode
4
Mid-Atlantic Bight ROMS Model for IS4DVAR
5 km resolution IS4DVAR
model embedded in …
… ROMS MAB-GoM V6
which uses global
HyCOM+NCODA
boundary data
5
Sequential assimilation of SLA and SST
Before attempting assimilation of all in situ data for a
full ESPreSSO reanalysis, we are assimilating satellite
SSH and SST to tune for the assimilation parameters
(horizontal and vertical de-correlation scales, duration
of assimilation window, etc.)
Unassimilated hydrographic data are used to evaluate
how well the adjoint model propagates information
between variables, and in space and time.
6
IS4DVAR*
R(xo )
xo
• Given a first guess (the forward trajectory)…
• and given the available data…
*Incremental Strong Constraint 4-Dimensional Variational
data assimilation
7
IS4DVAR
R(xo )
 xo
R(x o   x o )
• Given a first guess (the forward trajectory)…
• and given the available data…
• what change (or increment) to the initial
conditions (IC) produces a new forward trajectory
that better fits the observations?
8
The best fit becomes the analysis
assimilation window
ti = analysis
initial time
tf = analysis
final time
The strong constraint requires the trajectory satisfies the physics in ROMS.
The Adjoint enforces the consistency among state variables.
9
The final analysis state becomes the IC for
the forecast window
assimilation window
tf = analysis
final time
forecast
tf + t = forecast
horizon
10
Forecast verification is with respect to data
not yet assimilated
assimilation window
forecast
verification
tf + t = forecast
horizon
11
Basic IS4DVAR procedure:
 dx

L  J (x)   λ  i  N(xi )  Fi 
 dt

i 1
N
Lagrange function
Lagrange multiplier
J = modeldata misfit
T
i
Fi  F(it )
xi  x(it )
λ i  λ (ti )  λ (it )
N
1
1
T
T
1
J ( x)   x  xb  B  x  xb     H i xi  y i  O 1  H i xi  y i 
2
i 1 2
Jb
Jo
The “best” simulation will minimize L:
model model-data misfit is small and
model physics are satisfied
12
Basic IS4DVAR procedure:
Fi  F(it )
 dx

L  J (x)   λ Ti  i  N(xi )  Fi 
 dt

i 1
N
Lagrange function
xi  x(it )
J ( x) 
Lagrange multiplier
λ i  λ (ti )  λ (it )
J = modeldata misfit
1
T
 x  x b  B 1  x  x b  
2
Jb
N
1
 2 H x
i
i 1
 y i  O 1  H i x i  y i 
T
i
Jo
The “best” simulation minimizes L:
At extrema of L
we require:













L
0
λ i

L
0
xi
dλ  N 
T
1
  i 
λ


H
O
 Hxm  y m 
i
im

dt  x 
dxi
 N(xi )  Fi  0
dt
NLROMS
T
L
 0  B 1  x(0)  xb   λ (0)
x(0)
L
 0  λ (t )  0
x(t )
ADROMS
coupling of NL & AD
i.c. of ADROMS
13
Basic IS4DVAR procedure:
(1)
Choose an x(0)  xb (0)
(2)
Integrate NLROMS t  [0,t ] and save x(t )
(a) Choose a  x(0)
J = modeldata misfit
Inner-loop (3)
Outer-loop (10)
(b) Integrate TLROMS t  [0,t ] and compute J
(c) Integrate ADROMS t  [t , 0] to yield J o  λ (0)
 x(0)
J
(d) Compute
 B 1 x(0)  λ (0)
 x(0)
J ( x) 
1
T
 x  x b  B 1  x  x b  
2
Jb
N
1
 2 H x
i 1
i
 y i  O 1  H i x i  y i 
T
i
Jo
(e) Use a descent algorithm to determine a “down gradient”
correction to  x(0) that will yield a smaller value of J
(f) Back to (b) until converged
(3)
Compute new
x(0)  x(0)   x(0)
and back to (2) until converged
NLROMS = Non-linear forward model; TLROMS = Tangent linear; ADROMS = Adjoint14
xb = model state
(background) at end of
previous cycle, and 1st
guess for the next
forecast
xb
In 4D-Var assimilation
the adjoint gives the
sensitivity of the initial
conditions to mismatch between model
and data
previous
forecast
0
1
2
3
4
time
Observations minus Previous Forecast
A descent algorithm
uses this sensitivity to
iteratively update the
initial conditions, xa,
(analysis) to minimize
Jb+ S(Jo)
x
Adjoint model
integration is forced by
the model-data error
15
Observed information (e.g. SLA, SST) is transferred to
unobserved state variables and
projected from surface to subsurface in 3 ways:
(1) The Adjoint Model
(2) Empirical statistical correlations to generate
“synthetic XBT/CTD”
 In EAC assimilation get T(z),S(z) from
vertical EOFs of historical CTD
observations regressed on SSH and SST
(3) Modeling of the background covariance matrix
 e.g. via the hydrostatic/geostrophic relation
16
MAB Satellite Observations for IS4DVAR
5 km resolution for IS4DVAR
1 km downscale for forecast
SST 5-km daily blended MW+IR
from NOAA PFEG Coastwatch
MAB Sea Level Anomaly (SLA) is
strongly anisotropic with short
length scales due to flowtopography interaction, so use
along-track altimetry (need
coastal altimetry corrections
for shelf data)
• 4DVar uses all data at time of
•
satellite pass
model “grids” data by
simultaneously matching
observations and dynamical
and kinematic constraints
17
Mid-Atlantic Bight ROMS Model for IS4DVAR
Model variance (without
assimilation) is comparable
to along-track in Slope Sea,
but not shelf-break
AVISO gridded SLA differs
from along-track SLA in
Slope Sea (4 cm) and Gulf
Stream (10 cm)
18
All inputs:
NAM
Ocean model based open boundary conditions
River discharge, temperature (USGS)
Altimetry (via RADS; AVISO gridded)
XBT, CTD, Argo
Satellite SST – IR and mWave, passes/blended
HF radar – totals/radials
Cabled observatory time series – MVCO
Glider CTD (and optics)
NDBC buoy time series (T, S, velocity)
tide gauges
waves
Drifters - SLDMB and AOML GDP
Delayed mode
Oleander ADCP
science moorings
19
Assimilation of hydrographic climatology for:
* mean dynamic topography (altimetry)
* removing model bias
• Bias in the background state adversely affects how IS4DVAR
•
•
projects model-data misfit across variables and dimensions
We assimilate a high-resolution (~2-5 km) regional temp/salt
climatology to (i) produce a Mean Dynamic Topography (SSH)
consistent with model physics, and (ii) to remove bias
Climatology computed by weighted least squares (Dunn et al. 2002,
JAOT) from all available T-S data (NODC, NMFS) prior to 2006
(Naomi Fleming)
• Three simulations:
1. ROMS nested in MABGOM V6
2. Free running ROMS initialized with climatology and forced by
climatology at the boundaries and mean surface wind stress
3. ROMS with climatology initial/boundary/forcing and assimilation of
climatology over a 2-day window
20
21
Skill of climatologies and MABGOM-V6 at reproducing all
XBT/CTD from GTS in 2007-2008 in Slope Sea
22
Skill of climatologies and MABGOM-V6 at reproducing all
XBT/CTD from GTS in 2007-2008 in MAB shelf waters
23
24
25
26
27
Mean barotropic velocity from ROMS versus mean alongshelf velocity
from analysis of mooring observations by Lentz (2008)
Blue – mean of ROMS v6
Red – mean of clim ROMS
Black – mean of assim ROMS
Green - observations
28
29
High frequency variability:
model and data issues
ROMS includes high frequency variability typically
removed in altimeter processing (tides, storm surge)
The IS4DVAR cost function, J, samples this high
frequency variability, so it must be either (a) removed
from the model or (b) included in the data
Our approach:
• Run 1-year ROMS (no assimilation) forced by boundary
TPX0.7 tides; compute ROMS tidal harmonics
• de-tide along-track altimetry (developmental in MAB)
• add ROMS tides to de-tided altimeter data
• thus the observations are adjusted to include model tide
• assimilate – high frequency mismatch of model and
altimeter is minimized and cost function is, presumably,
dominated by sub-inertial frequency dynamics
30
High frequency variability:
model and data issues
The IS4DVAR increment is to the initial conditions of
the analysis window, and this itself generates HF
variability (inertial oscillations)
31
High frequency variability:
model and data issues
The IS4DVAR increment is to the initial conditions of
the analysis window, and this itself generates HF
variability (inertial oscillations)
Our approach:
• Apply a short time-domain filter to IS4DVAR initial
conditions
• Reduces inertial oscillations in the Slope Sea
but removes tides
• Tides recover quickly
– approach needs refinement
– possibly using 3-D velocity harmonic analysis of
free running model
32
High frequency variability:
model and data issues
Without a subsurface synthetic-CTD
relationship, the adjoint model can
erroneously accommodate too much of
the SLA model-data misfit in the
barotropic mode
This sends gravity wave at gh along the
model perimeter
Our approach:
• Repeat (duplicate) the altimeter SLA observations at
t = -6 hour, t=0 and t = +6 hour
but with appropriate time lags in the added tide signal
• These data cannot easily be matched by a gh wave
• We are effectively acknowledging the temporal correlation
of the sub-tidal altimeter SLA data
33
High frequency variability:
model and data issues
gh
Our approach:
• Repeat (duplicate) the altimeter SLA observations at
t = -6 hour, t=0 and t = +6 hour
but with appropriate time lags in the added tide signal
• These data cannot easily be matched by a gh wave
• We are effectively acknowledging the temporal correlation
of the sub-tidal altimeter SLA data
34
Sequential assimilation of SLA and SST
Before attempting assimilation of all in situ data for a
full ESPreSSO reanalysis, we are assimilating satellite
SSH and SST to tune for the assimilation parameters
(horizontal and vertical de-correlation scales, duration
of assimilation window, etc.)
Unassimilated hydrographic data are used to evaluate
how well the adjoint model propagates information
between variables, and in space and time.
35
Sequential assimilation of SLA and SST
• Reference time is days
after 01-01-2006
• 3-day assimilation
window (AW)
• Daily MW+IR blended
SST (available real time)
• SSH = Dynamic
topography + ROMS
tides + Jason-1 SLA
(repeated three times)
• For the first AW we just
assimilate SST to allow
the tides to ramp up.
36
37
38
Sequential assimilation of SLA and SST
Assimilation window (3<=t<=6 days)
Observed SST
ROMS SST and currents at 200 m
XBT transect
(NOT assimilated)
Jason-1 data
39
Sequential assimilation of SLA and SST
ROMS solutions along the transect positions [lon,lat,time]
40
Sequential assimilation of SLA and SST
ROMS-IS4DVAR fits the surface observations (SST
and SSH), but how well does it represent
unassimilated subsurface data?
ROMS solutions along the transect positions [lon,lat,time]
41
Forward model
Assimilation of SST and SSH
(no climatology bias correction)
depth (m)
depth (m)
42
ROMS data assimilation for ESPreSSO
Accomplishments:
* Have a system ready for:
1. introduction of in situ physical data into reanalysis
2. 2006-2008 reanalysis of ocean physics
3. analysis of impact of improved physics on
ecosystem (‘fasham’) and optical models
4. construction of adjoint/tangent-linear of optical
model, and subsequent addition of optical data to
cost function and full IS4DVAR
43
IS4DVAR data assimilation
LaTTE: The Lagrangian Transport and Transformation
Experiment
system set-up:
•
•
•
•
•
resolution:
forcing:
rivers:
DA window:
period:
2.5km
NAM model output
USGS Hudson & Delaware gauges
3 days
Apr. 10 – Jun 6, 2006
algorithm:
Incremental Strong-constraint 4DVAR
types and numbers of obs.
(Courtier et al, 1994, QJRMS; Weaver et al, 2003, MWR;
Powell et al, 2008, Ocean Modelling)
1 Nobs
1
J   (Hi Φ  y i )T O1 (Hi Φ  y i )  φT0 B1φ0
2 i 0
2
44
---- reduction of misfit
2006-04-20
06:57:36
IS4DVAR result
model
evolution of cost function
observation
45
IS4DVAR result ---- forecast skills
skill = 1 
RMSafterDA
RMSbeforeDA
RMSafterDA
RMSbeforeDA
1-CC afterDA
1-CC beforeDA
46
Adjoint sensitivity results
J SST
J 
day 0
J
J
J
J
 u (0) 
 T (0) 
 (0) 
 h (0) 
u (0)
T (0)
 (0)
 h (0)
Upstream
temperature
Density
Surface
current
J X
104
105
2  10 4
X
2
1
J
 X (C 2 )
X
2  10 4
105
X
SSH
Viscosity
Diffusion
3 105
1
0.3
101
102
105
106
2  10 5
3 107
105
3 107
47
Ensemble measure of the influence of glider MURI track at
the end of the glider mission
t
2
1
2
2

 dLdt
Cost function: J 
(
T

T
)

(
S

S
)



L(t2  t1 ) t1 L
Covariance between
J and temperature,
cov( J , T ( x, y, z, t )),
reflects the
influence of glider
observation, as
plotted in the right.
t: the finish time of a
glider mission.
48
Ensemble measure of the influence of glider MURI
track 5 days after the glider mission
t: 5 days after the
mission is finished.
49
Observation evaluation
 (T  T )2 ( S  S )2 
1
J


 dtdV
V t V t  OT
OS 
Assuming: model error ~
ocean state anomaly
glider
 (T  T )2 ( S  S )2 
1
J


 dtdV
V t V t  OT
OS 
Mooring
observation window
forecast window
50
Observation evaluation (cont’d)
northerly wind
southerly wind
 (T  T )2 ( S  S )2 
1
J


 dtdV


V t V t  OT
OS 
observation window
forecast window
51