EOLDAS¶
Data Assimilation and EOLDAS¶
The purpose of this section is to present an overview of the eoldas coding/configuration approach now that have had the chance in previous sections to use some of its functionality for filtering (smoothing) and for attempting to estimate state variables using a non linear radiative transfer model as an observation oeprator.
Introduction¶
In the previous sections, we have seen that we can phrase the ‘remote sensing problem’ (inference of state variables from remote radiometric measurements) as an optimisation problem with multiple constraints. More general than the remote sensing problem is our desire to provide the most reliable estimate of state vector elements (properties of the land surface) over space/time, given an expression of our current understanding of e.g. biogeochemical processes (encoded as models) and other forms of evidence, including Earth Observation data.
Assuming Gaussian statistics are sufficient to describe state (as is currently done in eoldas, other than allowing state variable transforms), we recall the cost function we want to minimise is:
and that we can simply sum multiple constraints and perform a minimisation of that over the space of the state variables.
We have seen and used this to represent:
An Identity (Prior) Operator
Setting
as our a priori mean expectation of 
, 
as expectation
of how good that expectation is, and the operator itself, 
A Model Operator
Constraining the relationships between state vector elements over space/time. We have seen, and used in eoldas the particular example of a differential operator, setting
as zero,
using a model 
, where 
is the order of the
differential that we expect to be zero, and as the uncertainty in the 
order
derivative (over time, space) of the state variable being zero.An Observation operator
Where
allows us to map between the space of our state vector ( ) and our EO datasets
(
), with being the uncertainty we have in our observations.
Part of the design ideas of eoldas are to allow users to implement their own operators. Many can be easily incorporated into the existing operator classes, but sometimes there are specific issues (e.g. regarding efficiency) that mean that it is better to write a new operator. On the whole, these should be derived from the base Operator class (or classes derived from that) to be consistent with eoldas.
Of particular interest will be to extend the group of model operators beyond the existing differential model.
Operators¶
Within eoldas, we have implementations of these ‘flavours’ of operator within the following modules:
eoldas_Operator
A base class operator, from which others are defined. Functionally, this implements, for an Identity Operator
,
methods to return 
, the differential of with respoect to each element of the state vector, 
(a vector
of the same shape as the state vector), and the second order derivative 
, a matrix with rank the size of the state vector
(though it might be sparse).
This module further provides interfaces for returning and its first and second (only first in this version) derivatives
with respect to the state vector.
It is this core functionality, or a subset of it allowing the calculation of at least and its derivative is required by
all operators. If a method to calculate is not available, a numerical approximation is used ut this can be very slow.
A further function of this base Operator class is to store a list of all of the operators being used and to loop over them to
return the sum of and the sum of and over all operators. Thus, in a sence, a run of eoldas
has only a single operator, but this itself contains further operators.eoldas_DModel_Operator
A module with
The differential operator (to an arbitrary (integer) order) can be be applied to any ‘location’ dimension. It is applied only to a single ‘location’ dimension, so if for example you want to constrain differences in space and time, and your data are described as a function of these ‘locations’, then you would set up three operators, one for each.
A control on the overall uncertainty in the model is provided by the
, which is specified as a state variable (for convenience
and in case we attempt to estimate it more directly in future releases) which must have the name gamma_{location}, e.g.
gamma_time if data are specified as a function of time, or gamma_row or gamma_col for spatial operators.More refined control on the effective ‘degree of smoothness’ can be specified through the
uncertainty
if desired. Note that the default value of is zero (i.e. it is set to zero if you do not
specify it). It is a small step from this class to a generalised linear model (over space and time), but we
have not had the need for that yet.
Various boundary conditions can be specified for the oeprator. The default is None. Also available are ‘periodic’ (with a specified period) or reflexive (though there is currently a bug in that, so don’t use it yet).
The default discrete differiential in this operator is performed between samples in
eoldas_Observation_Operator
The observation operator module, especially tailored to representation of EO-data-like observational data. It has an attempt at a generic interface to radiative transfer models (see the module documentation and code) that should allow users to plug their own observation operators into eoldas.
eoldas_Kernel_Operator
A particular class of observation operators, the family of linear kernels models, Isotropic, Volumetric, and Geometric. These are linear which is a major speed advantage (although we do not as yet bother to use a specific linear solver).
This module makes a call to the kernels class (eoldas_Kernel.py eoldaslib/eoldas_Kernel.py) which impements a large number of linear kernels, but the current interface is only to those defined for the MODIS BRDF/albedo product (also used in the ESA GlobAlbedo product among others).
Although this class contains a form of radiative transfer model, it is termed ‘non spectral’ in eoldas. This means that one must define a state variable (Isotropic, Volumetric, and Geometric) for each waveband in the dataset. This makes formal merging of data with different spectral sampling a little problematic (though not impossible ... see GlobAlbedo).
Observation_Operator_AtCorr
An atmospheric operator on top of a surface operator. This is essentially just an interface to the 6s code. This module works with a Lambertian surface assumption (i.e. does not fully couple the surafce and atmosphere radiation). The reason for this is that it is much faster to compute (but still slow, as this m,odule is making a large number of calls directly to the 6s code, rather than using LUTS or other pragmatic approaches). Really this module should only be used for forward modelling (i.e. simulation of TOA data). It is not recommended to use it for state vector estimation at present.
Observation_Operator_AtCorrCpld
An operator which couples surface atmosphere radiometric interactions.
This is even slower than Observation_Operator_AtCorr because it has to provide integral terms for the coupling. Again, this module should only be used for forward modelling (i.e. simulation of TOA data). It is not recommended to use
it for state vector estimation at present.
Other functions provided by the Operator family of classes (well, provide in the base class) are methods for packing and unpacking data to pass state variables between a ‘master’ representation and the ‘slave’ operators that may function only on a subset of the sull state representation. There are no methods as yet for anything beyond simple subsetting, though one can envisage that spatial modelling of the satellite data footprint would be appropriate in this interface.
Limitations¶
Whilst eoldas is designed to be very flexible and to find the minimum of the function above
the current implemtation has a few limitations worth noting.
- The uncertainty arising from the Operator :math:`mathbf{C_{H}}^{-1} is currently ignored although mechanisms are in place in the code to treat it.
- Whilst full uncertainty matrices can be specified, it is probably safest at present to specifiy only dioagonal uncertainty elements, as not all operators can be guaranteed ‘safe’ for fuller specification of uncertainty.
Using Model and Observation Operators¶
Having outlined something of the design structure and philisophy, and having gained some experience with eoldas functionailty we can now move on to more intersting exanmples merging models and observations.
As mentioned above, the only ‘dynamic’ (well, spatio-temporal) model currently implemented in eoldas is the differential constraint model. That said, we have seen in the various filtering examplesd that it is a powerful tool for constraining DA problems. This is especially true where there is a lack of any current physically0based (space/time) models.
For land surface DA, we are particularly hampered by a lack of any real physical models that we can apply, so such regularisation methods will have to sderve there for some time. One criticism we can put on these models is that they tend to over-smooth sharp changes. This is true of any simple smoothing approach (where the filter ‘influence’ or weight is constant in space/time) but there are many ways of improving on this that can be expplored in future versions of eoldas.
Temporal model and semidiscrete Observation Operator¶
In this section, we demonstrate the applicatiomn of eoldas to the temporal monitoring of vegetaion state variables using the semidiscrete Observation Operator and the differential model.
The dataset we use is MODIS surface reflectance data over a field site in Germany.
We use a generic configuration file semid_default.conf:
# configuration file
# note that some words are protected: values, keys
# so dont use them
# other than that, you can set any level of hierarchy in the
# representation
# this next part serves no purpose other than
# to show how different options can be set
# It is a requirement that parameter.names be defined
[parameter]
location = ['time']
limits = [[1,365,1]]
names=gamma_time,xlai, xhc, rpl, xkab, scen, xkw, xkm, xleafn, xs1,xs2,xs3,xs4,lad
solve = [1]*len($parameter.names)
help_solve='flags for which state vector elements to solve for'
[parameter.result]
filename = 'test/semid/output/test.params'
help_filename="state vector results file"
format = 'PARAMETERS'
[parameter.assoc_solve]
gamma_time = 0
lad = 0
xs3 = 0
xs4 = 0
xhc = 0
rpl = 0
[parameter.x]
datatype = x
names = $parameter.names
default = [0.05]*len($parameter.names)
help_default = "Set the parameter default values"
apply_grid = True
sd = [1.]*len($parameter.names)
bounds = [[0.01,0.99]]*len($parameter.names)
#state = fullState.dat
invtransform=$parameter.names
transform=$parameter.names
[parameter.x.assoc_transform]
xlai=np.exp(-xlai/2.)
xkab=np.exp(-xkab/100.)
xkw=np.exp(-xkw*50.)
xkm=np.exp(-100.*xkm)
[parameter.x.assoc_invtransform]
xlai=-2.*np.log(xlai)
xkab=-100.*np.log(xkab)
xkw=-(1./50.)*np.log(xkw)
xkm=-(1./100.)*np.log(xkm)
[parameter.x.assoc_bounds]
gamma_time = 0.001,1000
xlai = 0.01,0.99
xhc = 0.01,10.0
rpl = 0.01,0.10
xkab = 0.1,0.99
scen = 0.001,1.0
xkw = 0.01,0.99
xkm = 0.3,0.9
xleafn = 0.9,2.5
xs1 = 0.01, 4.
xs2 = 0.01, 5.
xs3 = 0., 0.
xs4 = 0.,0.
lad = 1,5
[parameter.x.assoc_default]
####gamma_time = 0.01
xlai = 0.95
xhc = 5
rpl = 0.01
xkab = 0.95
scen = 0.001
xkw = 0.95
xkm = 0.35
xleafn = 1.5
xs1 = 1.0
xs2 = 0.001
xs3 = 0
xs4 = 0
lad = 5
help_lad = 'lad value'
[general]
datadir = test/semid/input/,.
is_spectral = True
calc_posterior_unc=False
write_results=True
doplot=True
plotmod=10
plotmovie=True
[general.optimisation]
randomise=False
help_randomise="randomise the state vector'
[operator]
modelt.name=DModel_Operator
modelt.datatypes = x
obs.name=Observation_Operator
obs.datatypes = x,y
[operator.obs.rt_model]
model=rtmodel_ad_trans2
use_median=True
help_use_median = "Flag to state whether full bandpass function should be used or not.\nIf True, then the median wavelength of the bandpass function is used"
bounds = [400,2500,1]
help_bounds = "The spectral bounds (min,max,step) for the operator'
ignore_derivative=False
help_ignore_derivative = "Set to True to override loading any defined derivative functions in the library and use numerical approximations instead"
[operator.modelt.x]
names = $parameter.names
sd = [1.0]*len($operator.modelt.x.names)
datatype = x
[operator.modelt.rt_model]
model_order=2
help_model_order='The differential model order'
wraparound=periodic,365
help_wraparound='Wraparound conditions (periodic,reflexive,none). If periodic, then the period can be specified as a second a
rgument (e.g. --operator.modelt.rt_model.wraparound=periodic,365)'
lag=1
help_lag='The grid lag spacing between samples for the differential operator. This can be a list, in which case lag weight is
also required'
lag_weight = 1
help_lag_weight='The weight for the lag. Only of value if multiple lags are specified'
[operator.obs.x]
names = $parameter.names[1:]
sd = [1.0]*len($operator.obs.x.names)
datatype = x
[operator.obs.y]
control = 'mask vza vaa sza saa'.split()
names = "465.6 553.6 645.5 856.5 1241.6 1629.1 2114.1".split()
help_names="set observation band names"
sd = "0.003 0.004 0.004 0.015 0.013 0.01 0.006".split()
help_sd="set observation operator sd"
datatype = y
state = test/data_type/input/test.brf
help_state='set the obs state file'
# angle_transform not yet implemented
#angle_transform = False
#help_angle_transform = "Transform the relative azimuth of the input data by adding 180 degrees. Useful if model and data are not consistently defined."
[operator.obs.y.result]
filename = 'test/semid/output/test-fwd.brf'
help_filename = 'forward modelling results file'
format = 'PARAMETERS'
The configuration file is similar to ther ones we have come across, but we note that by assigning help_filename at various places where we define the input and output files, we can gain control over them from the command line.
We then conduct an experiment where we modify the gamma_time value from the command line:
eoldaslib/eoldas.py –conf=semid_default.conf –logfile=logs/file.log –parameter.result.filename=output/gamma2/state.dat –operator.obs.y.result.filename=output/gamma2/obs.dat –parameter.x.default=2,0.99,5,0.01,0.99,0.001,0.99,0.35,1.5,1,0.001,0,0,5 –operator.obs.y.state=data/modis/brdf_WW_1_A_1.kernelFiltered.dat
This solves for 2920 state variables, being 8 biophysical parameters of the semidiscrete model for 365 days, with a low value of gamma. We expect that this will allow us to arrive at a solution for the state variables for each day, but that the results will be rather noisy.
The resultant state estimate is in output/gamma2/state.dat:
#PARAMETERS time gamma_time xlai xhc rpl xkab scen xkw xkm xleafn xs1 xs2 xs3 xs4 lad sd-gamma_time sd-xlai sd-xhc sd-rpl sd-xkab sd-scen sd-xkw sd-xkm sd-xleafn sd-xs1 sd-xs2 sd-xs3 sd-xs4 sd-lad
62.000000 2.000000 0.878534 5.000000 0.010000 0.984610 0.061855 0.804596 0.432753 1.603578 1.052247 0.128133 0.000000 0.000000 5.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000
63.000000 2.000000 0.937961 5.000000 0.010000 0.990000 0.001000 0.674791 0.434153 1.573094 0.910875 0.010000 0.000000 0.000000 5.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000
64.000000 2.000000 0.965211 5.000000 0.010000 0.989666 0.001000 0.841901 0.388219 1.541949 0.933983 0.010000 0.000000 0.000000 5.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000
65.000000 2.000000 0.982049 5.000000 0.010000 0.989827 0.001000 0.964509 0.353920 1.508224 0.982145 0.010959 0.000000 0.000000 5.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000
66.000000 2.000000 0.990000 5.000000 0.010000 0.989351 0.001000 0.990000 0.342428 1.497282 1.007904 0.010000 0.000000 0.000000 5.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000
67.000000 2.000000 0.954391 5.000000 0.010000 0.983504 0.024033 0.949792 0.345584 1.493383 0.997615 0.049001 0.000000 0.000000 5.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000
68.000000 2.000000 0.795616 5.000000 0.010000 0.972336 0.244144 0.843069 0.416827 1.475381 0.921735 0.200826 0.000000 0.000000 5.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000
69.000000 2.000000 0.569834 5.000000 0.010000 0.964864 0.681852 0.773026 0.583646 1.448639 0.791251 0.416116 0.000000 0.000000 5.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000
70.000000 2.000000 0.711732 5.000000 0.010000 0.975370 0.616745 0.958091 0.604918 1.436592 0.864060 0.391041 0.000000 0.000000 5.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000
71.000000 2.000000 0.880724 5.000000 0.010000 0.982872 0.196240 0.990000 0.439735 1.485515 1.000573 0.227511 0.000000 0.000000 5.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000
72.000000 2.000000 0.973190 5.000000 0.010000 0.988118 0.010014 0.990000 0.355624 1.505714 1.027007 0.073858 0.000000 0.000000 5.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000
and plotted in output/gamma2/state.dat.plot.x.png
Another interesting graphic is the plot of the observations and H(x) at the same locations and wavelengths. The data for this are in output/gamma2/obs.dat and output/gamma2/obs.dat_orig
#PARAMETERS time mask vza vaa sza saa 465.6 553.6 645.5 856.5 1241.6 1629.1 2114.1 sd-465.6 sd-553.6 sd-645.5 sd-856.5 sd-1241.6 sd-1629.1 sd-2114.1
232.000000 1.000000 42.220000 51.240000 42.340000 0.000000 0.045173 0.078350 0.081959 0.235935 0.287922 0.265706 0.172931 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
233.000000 1.000000 31.380000 242.230000 39.510000 0.000000 0.058816 0.088951 0.105734 0.266189 0.324199 0.304242 0.209032 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
233.000000 1.000000 34.710000 -53.000000 42.100000 0.000000 0.056236 0.087703 0.104259 0.260085 0.319270 0.301285 0.210835 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
248.000000 1.000000 42.220000 52.340000 48.040000 0.000000 0.071551 0.103765 0.101696 0.181783 0.276728 0.277846 0.192856 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
254.000000 1.000000 9.980000 -60.890000 47.880000 0.000000 0.044740 0.069006 0.070869 0.148782 0.224932 0.216106 0.136945 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
258.000000 1.000000 37.970000 241.510000 48.590000 0.000000 0.022115 0.034182 0.032888 0.113303 0.166716 0.156624 0.092686 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
261.000000 1.000000 0.590000 -71.270000 50.150000 0.000000 0.048773 0.061643 0.068913 0.115184 0.185096 0.185523 0.126213 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
264.000000 1.000000 39.270000 -246.740000 50.530000 0.000000 0.040916 0.057992 0.069663 0.145711 0.241789 0.256277 0.206893 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
265.000000 1.000000 31.460000 239.490000 51.540000 0.000000 0.040316 0.060021 0.068539 0.141743 0.227941 0.242379 0.201061 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
265.000000 1.000000 34.420000 -60.670000 52.770000 0.000000 0.039633 0.058641 0.068079 0.139110 0.226108 0.241204 0.200986 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
266.000000 1.000000 25.680000 -247.150000 51.480000 0.000000 0.035162 0.049204 0.060163 0.121762 0.203091 0.216431 0.177464 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#PARAMETERS time mask vza vaa sza saa 465.6 553.6 645.5 856.5 1241.6 1629.1 2114.1 sd-465.6 sd-553.6 sd-645.5 sd-856.5 sd-1241.6 sd-1629.1 sd-2114.1
232.000000 1.000000 42.220000 51.240000 42.340000 0.000000 0.053800 0.099900 0.103100 0.303700 0.319300 0.310800 0.204300 0.003000 0.004000 0.004000 0.015000 0.013000 0.010000 0.006000
233.000000 1.000000 31.380000 242.230000 39.510000 0.000000 0.045700 0.085000 0.097700 0.239400 0.291700 0.229000 0.198800 0.003000 0.004000 0.004000 0.015000 0.013000 0.010000 0.006000
233.000000 1.000000 34.710000 -53.000000 42.100000 0.000000 0.054300 0.102100 0.123100 0.266400 0.366400 0.356000 0.226100 0.003000 0.004000 0.004000 0.015000 0.013000 0.010000 0.006000
248.000000 1.000000 42.220000 52.340000 48.040000 0.000000 0.081500 0.101900 0.084800 0.215700 0.298100 0.264200 0.195700 0.003000 0.004000 0.004000 0.015000 0.013000 0.010000 0.006000
254.000000 1.000000 9.980000 -60.890000 47.880000 0.000000 0.040100 0.067000 0.077200 0.155400 0.219600 0.211200 0.138800 0.003000 0.004000 0.004000 0.015000 0.013000 0.010000 0.006000
258.000000 1.000000 37.970000 241.510000 48.590000 0.000000 0.018500 0.032100 0.039500 0.100400 0.161800 0.175700 0.084900 0.003000 0.004000 0.004000 0.015000 0.013000 0.010000 0.006000
261.000000 1.000000 0.590000 -71.270000 50.150000 0.000000 0.050200 0.063800 0.065200 0.122200 0.191000 0.176400 0.128200 0.003000 0.004000 0.004000 0.015000 0.013000 0.010000 0.006000
264.000000 1.000000 39.270000 -246.740000 50.530000 0.000000 0.036800 0.057900 0.073600 0.160900 0.242000 0.249700 0.208600 0.003000 0.004000 0.004000 0.015000 0.013000 0.010000 0.006000
265.000000 1.000000 31.460000 239.490000 51.540000 0.000000 0.033700 0.050700 0.066000 0.125500 0.183700 0.207300 0.183000 0.003000 0.004000 0.004000 0.015000 0.013000 0.010000 0.006000
265.000000 1.000000 34.420000 -60.670000 52.770000 0.000000 0.043800 0.068600 0.081300 0.166300 0.251800 0.269200 0.219300 0.003000 0.004000 0.004000 0.015000 0.013000 0.010000 0.006000
266.000000 1.000000 25.680000 -247.150000 51.480000 0.000000 0.032800 0.049300 0.061100 0.125200 0.200500 0.215000 0.179200 0.003000 0.004000 0.004000 0.015000 0.013000 0.010000 0.006000
And the plots in output/gamma2/obs.dat.plot.y.png and output/gamma2/obs.dat.plot.y2.png, the former showing time sample spectral plots and the latter spectral sample time plots.
The processing for this takes some time, but that is partly because we are using the default convergence thresholds.
We can run this experiment for multiple values of gamma (multiplying by e.g. a factor of 2 each step) e.g. for cross validation exercises. The experiment can be speeded up by passing the estimated state from one value of gamma onto the next step. To do this, we must modify the gamma value in the file output/gamma2/state.dat and replace it with the new gamma. This is not a difficult task using various unix or python approaches, although there is not directly a tool for this in eoldas.
For example we would wanto filter the result of running at gamma 2 and insert gamma = 4 into a new copy of it, called output/gamma4/state-next.dat. Then to run the exercise for gamma 4, we use the command line:
eoldaslib/eoldas.py –conf=semid_default.conf –logfile=logs/file.log –parameter.result.filename=output/gamma4/state.dat –operator.obs.y.result.filename=output/gamma4/obs.dat –operator.obs.y.state=data/modis/brdf_WW_1_A_1.kernelFiltered.dat –conf=confs/set_state.conf –parameter.x.state=output/gamma4/state-next.dat
All that is contained in the new configuration file set_state.conf is:
[parameter.x]
state=fullState.dat
help_state='set the state vector'
which essentially sets up the capability to use –parameter.x.state= on the command line.
Once we have looped over all gamma values that we want to examine, we can visualise the results as a movie (over gamma), using e.g. the unic command convert:
Then we can make an animated gif with e.g:
This makes an animate gif such as:
We can see here the impact of increasing gamma: the high frequency noise is gradually decreased, leaving only rather smooth trejectories of the state variables. There is clearly over smoothing at some sharp trajectories, but otherwise this is a useful result once soime form of cross validation has been conducted.
Whilst we were looping over gamma in our cross validation experiment, we did not calculate uncertainties in the output data to save time. A set of defaault ‘system’ configurations exist in system_confs/eoldas_config.conf, which includes:
calc_posterior_unc=False
help_calc_posterior_unc ="Switch to calculate the posterior uncertainty"
This sets calc_posterior_unc = False by default, but we can override this on the commmand line. If we have generated a result without calculating the posterior uncertainties then we can e.g use:
eoldaslib/eoldas.py –passer –calc_posterior_unc –conf=semid_default.conf –logfile=logs/file.log –parameter.result.filename=output/gamma2/state.dat –operator.obs.y.result.filename=output/gamma2/obs.dat –operator.obs.y.state=data/modis/brdf_WW_1_A_1.kernelFiltered.dat –conf=confs/set_state.conf –parameter.x.state=output/gamma2/state.dat
The flag –passer tells eoldas not to perform an optimisation, but merely to load, check and output the data.
Summary¶
We have seen here how to combine a simple differential model and a physically-based observation operator to provide an estimate of the state of vegetation using eoldas. Once a configuratiuon file is set up for a ‘style’ of experiment, much of the interaction can take place with the command line options, ideally nested in some form of looping structure.
There are many criticisms that could be laid at eoldas ... the only model we use is very simple, the radiative transfer models are slow, but we, the authors believe and hope that it will be a useful tool for the community to gain experience of optimal merging of land surafce EO data and their own models, and that, once the computattional capabilities inevitably increase enough and more is made of faster methods for the DA, this form of variational DA holds a lot of promise for how ‘in the future’ we will be making use of EO (and other data) for monitoring and modelling the Earth land surafce.
There are several directions that eoldas could clearly take from here. There are a few residual bugs that need to be ironed out, but most of all, it needs to be available to users to start to explore these concept further. It should also (with tehse user pages) make a useful teaching tool for those interested in such matter.
Lewis, UCL, September 2, 2001