Constraining model parameters on remotely sensed evaporation : justification for distribution in ungauged basins ?

In this study, land surface related parameter distributions of a conceptual semi-distributed hydrological model are constrained by employing time series of satellite-based evaporation estimates during the dry season as explanatory information. The approach has been applied to the ungauged Luangwa river basin (150 000 (km) 2) in Zambia. The information contained in these evaporation estimates imposes compliance of the model with the largest outgoing water balance term, evaporation, and a spatially and temporally realistic depletion of soil moisture within the dry season. The model results in turn provide a better understanding of the information density of remotely sensed evaporation. Model parameters to which evaporation is sensitive, have been spatially distributed on the basis of dominant land cover characteristics. Consequently, their values were conditioned by means of Monte-Carlo sampling and evaluation on satellite evaporation estimates. The results show that behavioural parameter sets for model units with similar land cover are indeed clustered. The clustering reveals hydrologically meaningful signatures in the parameter response surface: wetlanddominated areas (also called dambos) show optimal parameter ranges that reflect vegetation with a relatively small unsaturated zone (due to the shallow rooting depth of the vegetation) which is easily moisture stressed. The forested areas and highlands show parameter ranges that indicate a much deeper root zone which is more drought resistent. Clustering was consequently used to formulate fuzzy membership functions that can be used to constrain parameter realizations in further calibration. Unrealistic parameter ranges, found for instance in the high unsaturated soil zone values in the highCorrespondence to: H. C. Winsemius (h.c.winsemius@tudelft.nl) lands may indicate either overestimation of satellite-based evaporation or model structural deficiencies. We believe that in these areas, groundwater uptake into the root zone and lateral movement of groundwater should be included in the model structure. Furthermore, a less distinct parameter clustering was found for forested model units. We hypothesize that this is due to the presence of two dominant forest types that differ substantially in their moisture regime. This could indicate that the spatial discretization used in this study is oversimplified.


Introduction
Hydrological models in data sparse areas are often oversimplified.This is partly due to the lack of observational data to justify more complexity.As a result, parsimony in model parameters is often advocated to prevent the undesirable occurrence of equifinality (e.g.Beven and Binley, 1992;Beven and Freer, 2001;Savenije, 2001).Although parsimony results in simple and to a certain extent identifiable models, their predictive capacity, for instance of land cover changes, is rather small, because parameters usually have little physical meaning and cannot represent the variability inherent in the landscape.Even with simple models, parameters are often poorly identifiable (e.g.Uhlenbrook et al., 1999) and cannot be justifiably distributed in space in view of the problem of equifinality.A related issue is model structural uncertainty, which is probably even more difficult to define and quantify (Wagener and Gupta, 2005;Young, 2001;Fenicia et al., 2008).Furthermore, in many remote river basins, especially in developing countries, measurement networks are collapsing.Sometimes only old records (often Published by Copernicus Publications on behalf of the European Geosciences Union.1404 H. C. Winsemius et.al.: RR parameter constraints on evaporation from colonial periods) exist which cannot be confronted with new (remotely sensed) data sources.This further reduces our capability of hydrological understanding and compromises prediction in areas with emerging water scarcity, where these tools are needed most (Sivapalan, 2003).The low understanding we have of "ungauged basins" forces modellers and experimentalists to look beyond the classical approach of rainfall-runoff "curve-fitting" and to come up with either robust alternative strategies to find behavioural parameter distributions or to use other data sources than classical streamflow series.Kuczera (1983) showed that additional information (not necessarily streamflow) can be used to reduce uncertainty through updating of prior likelihoods by means of Bayes' law. with where the left-hand side represents the posterior probability of parameter set i of a given model M and the right handside represents the process of Bayesian updating with observation Y and joint prior distribution p(Y ).More and more modellers are applying this or similar methods, some for the purpose of taking into account new stream flow data as time proceeds (e.g.Freer et al., 1996), some for the purpose of presenting uncertainty based on a joint posterior parameter distribution (Kuczera, 1983) and some for the purpose of learning and consequently detecting model structural deficiencies and henceforth improving the model structure (Vaché and McDonnell, 2006;Son and Sivapalan, 2007;Fenicia et al., 2008).
In the era of remote sensing, new potentially interesting data sources emerge that may allow us to step-wise infer constraints on parameter distributions based on satellite measurements.Even though these data sources are often subject to a great deal of noise, resulting in a substantial and ill-quantifiable uncertainty, they can be employed as "soft data" (Seibert and McDonnell, 2002) to update prior likelihoods of parameters.Although this cannot be done by formal Bayesian updating, which would require the use of a formal likelihood measure, including knowledge about the nature of the model residuals, at least it can be done by imposing a certain degree of acceptance of a parameter set.Franks et al. (1998) made an attempt to constrain a model's parameter space by sequentially updating the parameter distribution using the Generalized Likelihood Uncertainty Estimation (GLUE, Beven and Binley, 1992;Beven and Freer, 2001), first using streamflow and consequently Synthetic Aperture Radar (SAR) estimates of saturated areas.In their study, a fuzzy estimate (rather than the estimate itself) of the total saturated area was used as a way to deal with uncertain-ties related to subjective choices within the estimation procedure and the observations used.Franks and Beven (1997) included fuzzy estimates of Landsat TM derived evaporation estimates in the uncertainty reduction of a land surface model.
Other examples where remotely sensed information is used for calibration (rather than updating of parameter probability distributions) are described by Campo et al. (2006), Johrar (2002) and Immerzeel and Droogers (2008).In the former, a correlation of SAR backscatter with the top-soil moisture content of a distributed hydrological model of the Arno (∼8230 (km) 2 ) was assumed.The authors recognised the problems of SAR soil moisture inference over vegetated areas.Moreover, SAR observations are difficult to apply on the large scale.In the latter two, first attempts of calibration with remotely sensed evaporation were performed.Johrar (2002) calibrated an agro-hydrological model on remotely sensed evaporation rates and validated with in-situ measurements of groundwater levels.Immerzeel and Droogers (2008) calibrated a SWAT model for the Bhima catchment (∼45 000 (km) 2 ) on remotely sensed evaporation estimates by using a global optimisation algorithm.Here, streamflow is regulated by reservoirs and thus unsuitable for calibration of natural hydrological processes.A marked difference between the first and the last two studies is that actual evaporation is a flux, completely equivalent to evaporation from a hydrological model, while SAR soil moisture estimates represent a state, which is not equivalent to the soil moisture state in a hydrological model (i.e. the problem of representativeness of the measurement, e.g.Liu and Gupta, 2007).
In this paper, an attempt is made to transfer a prior distribution of parameters that relate transpiration to soil moisture states, into a justifiable fuzzy posterior distribution, by constraining the priors on satellite-based evaporation estimates during a dry season in the Luangwa basin in Zambia.Evaporation in our terminology means the sum of all processes where water is tranferred from the liquid to the gas phase (Brutsaert, 1982).The satellite-based evaporation is based on thermal-infrared imagery, which compromises application in the cloud-covered wet season, which is the reason why only dry-season estimates have been used.Although formally not applied in this study, we could say that Y in Eq. ( 1) is the remotely sensed evaporation and M is the given model structure.We emphasize that the parameters that determine the depletion of the soil moisture zone in hydrological models also determine to a large extent the separation of rainfall into soil moisture and streamflow.If the approach followed is successful, then this opens up new opportunities for constraining rainfall-runoff models in ungauged basins.
The analysis is performed on a semi-distributed conceptual model of the Luangwa river basin, a semi-arid area in Zambia, where recent information on stream flow is not available.At the moment, the basin is clearly ungauged, having no reliable streamflow records after 1980 and poorly concomitant available time series of streamflow and rainfall (see  Huffman et al., 2007).
Fig. 1 for an overview of available data).The hypothesis is that the evaporation estimates impose a two-fold constraint on model parameters.First, the model is forced to obey the water balance and in particular, to follow a realistic depletion of soil moisture within the dry season.According to old monthly records of rainfall and streamflow, evaporation accounts for about 85% of the annually averaged water balance in this river basin so it is a strong prescriptor for the water balance.Second, the modeller can attempt to regionalize evaporation sensitive parameters making use of observed land cover.An additional benefit of this approach is that the accuracy of the rainfall estimates that is used to force these models, is not of direct importance for this step in the reduction of parameter uncertainty, as long as the moisture status in the end of the wet season is more or less accurate and the evaporation estimates do not exhibit significant bias.Evaporation estimation based on thermal-infrared satellite imagery can only be performed during cloud-free moments and is based on indirect estimation procedures, which may introduce a great deal of ill-quantifiable uncertainties (Wagener and Gupta, 2005).These are caused for instance by transferring of radiometric surface temperatures to land surface temperatures; undetected low clouds or aerosols; roughness and emissivity parameterisation and; in the case of the Surface Energy Balance Algorithm for Land (SEBAL, used in this study and described in a later section) the somewhat subjective choice of a "wet" and "dry" pixel as extremes in the surface energy balance.Therefore, the evaporation estimates obtained, are in this study used as a proxy for the response of the land surface, assuming that errors are uncorrelated in time and space and are of a random nature.

Study area description
The Luangwa basin (Fig. 2) is a relatively pristine and remote area with only a small amount of old hydrometric station data (at the time of writing only one was fully operational, installed in November 2007).Dominant land cover characteristics are given in Fig. 3.The approximately 100 km wide river valley (here referred to as riverine) consists of sandy/loam soils (among which black cotton soils), covered by typical tropical savanna vegetation such as low shrubs, interspersed with Mopane (Frost, 1996) and Acacia.This part of the catchment dries out completely during the dry season.Some low lying areas are interspersed with wetlands, locally called dambos, which are also used for land cultivation.In the Western part of the catchment, some areas with mixed forests were identified, mostly consisting of Miombo and Mopane forest.The North-Western boundary (the Muchinga escarpment) is a moist highland covered by densely forested pristine wetland areas.Here vegetation consists of a mixture of Miombo forest and other, partly evergreen species.It has a completely different hydro-climatology than the riverine savannas.Temperatures in the highlands are much lower and given the type of vegetation present, these areas have a higher capability of retaining moisture during the dry season than the lower savannah regions.In the last 10 years, the annual rainfall in the catchment was around 1000 mm per year (Fig. 2).Rainfall is concentrated in one wet season from November until April.Although the basin is largely ungauged, the heterogeneity of this area and the moisture limited evaporation, typical for semi-arid areas, makes it an excellent site for research on the applicability of spatially distributed remotely sensed data in hydrological models.
The hydrological response of this area is crucial to the operation of the Cahora Bassa reservoir in Mozambique, the downstream riparian country.The Luangwa can generate critical and unexpected peak flows during Zambezi floods.The Luangwa joins the Zambezi, closely upstream lake Cahora Bassa.Operators of Cahora Bassa are sometimes forced to release large amounts of water from the reservoir, not knowing the exact magnitude of the Luangwa floods (e.g. in February andMarch 2001, February 2007).The lack of knowledge of the Luangwa compromises optimal operation of the reservoir during these floods.

Evaporation estimates
To derive spatially and temporally variable estimates of the evaporation, the Surface Energy Balance Algorithm for Land (SEBAL) has been applied.An elaborate description of the current state of SEBAL is given by Allen et al. (2007).SE-BAL was validated both on field and catchment scale e.g. by Bastiaanssen et al. (1998) and Bastiaanssen et al. (2005).It was applied in several studies with varying applications.Some of them are described by Bastiaanssen et al. (2002); Schuurmans et al. (2003); Mohamed et al. (2004); Immerzeel and Droogers (2008); Gragne et al. (2008).15 MODIS TERRA images, ranging from May 2006 until October 2006 (dry season) have been processed.Unfortunately, evaporation cannot be estimated for a complete hydrological year, because during the wet season, no cloud-free images can be found for this region.SEBAL is a residual based energy bal-ance approach in which instantaneous estimates of the energy balance are made based on the energy balance: where ρ w is the density of water and "dry" extremes in the satellite image).Then, the following assumption is made: where a and b are calibration coefficients that can be found through calibration on the two anchor points.With this linear equation, T s −T a,ref is found for the whole satellite image.Finally, λE is found as the residual of Eq. (3).24-h evaporation is found by assuming that the evaporative fraction, ρ w λE/(R n −G) is constant over the day as given below (with time (t) in hours): Here, the subscript "inst" stands for "instantaneous".To yield period-averaged evaporation, daily surface conductance g s [L T −1 ] estimates were derived by inverting the Penman-Monteith equation (Monteith, 1981;Penman, 1948)).
s and γ are the slope of the vapour pressure curve at given air temperature and the psychrometric constant  physically based approach was followed to downscale the coarse grids of ECMWF to near-surface 1×1 (km) 2 variables (Voogt, 2006).A Jarvis model (Jarvis, 1976) was used to correct g s for temporal variability of the meteorological conditions within the period of time (Bastiaanssen et al., 1994;Farah, 2001).This estimate of g s was inserted in Penman-Monteith with the concurrent period-averaged meteorological conditions to yield period-averaged (typically 10-15 days) estimates of E over the dry season of 2006.

Conceptual model
A semi-distributed conceptual model has been set up for this study.Current state of the art land cover maps do not sufficiently describe the wetland and Miombo land cover in Eastern Zambia.Therefore the spatial distribution of the model was based on a one-year time series of decadal SPOT NDVI images.An unsupervised classification was performed to differentiate between different land covers.The characteristics of the main dominant land covers was roughly defined through a short field investigation, elevation differences and, where available through investigation of high resolution, google earth overflight information http://earth.google.com/.
Based on this information, 4 dominant land covers were defined: riverine, dambos (or wetlands), forested and highlands as described in Sect. 2. The regions were manually delineated into polygon-shaped model units to decrease the computation time (Fig. 3).The model structure applied, is a simplified version of the 1-dimensional box-model HBV (Fig. 4, Lindström et al., 1997).The goal of this study is in first principle to investigate what the information density of the evaporation data is with respect to model identification (i.e. to what extent can the evaporation data explain the spatial and temporal variability of the land surface response to moisture availability), not to find the optimal model structure or parameter set for a given catchment.Therefore the amount of parameters was kept as low as possible in order to gain parameter identifiability and to be able to physically interpret the results in terms of land surface response to rainfall, drought and meteorological forcing.The model now consists of an interception store with a parameter D [L T −1 ], an un- Non-linearity of the runoff generation is represented by the power function in Eq. ( 9) (given as a curve in Fig. 4).The runoff component is transferred to an upper zone S q [L].
Streamflow is generated from this zone, assuming it behaves as a linear reservoir with residence time K q [T].S q represents the fast flow generated from water bodies or dambos.Finally, a lower zone S s [L] (conceptualising groundwater) receives a maximum amount of percolation R [L T −1 ] from S q , determined by the parameter F perc [L T −1 ].
R= min(S q , F perc ) This zone also behaves as a linear reservoir, contributing to the base flow with one parameter K s [T] representing the average residence time of reservoir S s .An overview of all parameters, their physical meaning and unit is given in Table 1.Note that in this study, the focus is only on the soil reservoir of the model, not on the runoff generating reservoirs.An upward flux from the runoff generating reservoirs to the soil reservoir was deliberately excluded.This results in a model structure in which there are only 2 parameters, S max and l p and one state, S u , that influence the transpiration when there is no rainfall.Within the dry season, there is no significant sensitivity for parameter D and B, since they do not influence the depletion.D was fixed on 2 mm day −1 .The prior value of the parameter B was constrained, by making it dependent on soil texture.A normalised soil texture map was derived from the WISE-ISRIC dataset (Batjes, 2006)  types with their respective texture class -coarse (0), medium (0.5) or fine (1).B was roughly estimated by assuming that it has a value ranging between 1 and 4.5, linearly depending on the normalised soil texture between 0 and 1.Also the routing parameters F perc , K s and K q do not influence our results, because there is no feedback from the discharge generating reservoir towards the soil moisture.Therefore these parameters are not mentioned in the remainder of this article.
Daily rainfall input was taken from the Tropical Rainfall Measuring Mission (TRMM, Huffman et al., 2007).Product 3B42 was used, which is a gauge-corrected product based on merged rainfall estimates from several sensors.For potential evaporation, local daily measurements of wind speed, minimum and maximum temperature and relative humidity were obtained from the Zambian Meteorological Office.These  New et al., 2002) following the approach presented by Reynolds (1988).All meteorological fields were used to compute potential evaporation, following Allen et al. (1998).

Parameter distribution
The next step involves estimation of the two remaining evaporation-sensitive parameters S max and l p using the SE-BAL estimates as evaluation data set.Parameter sampling and model evaluation was done separately for each model unit.To this end, SEBAL evaporation estimates were lumped per model unit.SEBAL evaporation maps correspond to the average evaporation over a period of 10 to 15 days.When more than 25% of the pixels of a model unit within an evaporation map appeared to be cloud-covered (and thus were excluded from the SEBAL computation), the value was discarded for the evaluation.The evaporation, generated from each model unit was averaged over the same periods as the SEBAL estimates.

Model diagnostics and parameter identifiability
Of each model unit for which evaporation estimates of satisfactory quality were found, the best 2% of the parameter realizations (consisting of parameters S max and l p ) was rescaled to a posterior likelihood L s by where n is the number of realizations, belonging to the best 2%.Equation ( 12) is a transformation function that tranfers the likelihood measure L( i ) such, that the models with highest values for L s are most likely.The sum of all values of L s is equal to 1.The results from all model units that have the same land cover class were plotted together in Fig. 5.The number of model units belonging to the same class is indicated on the right-side of each sub-figure.For each land cover class, an example of a well-performing model is given in Fig. 6.The corresponding parameter sets and performances are given in Table 2.
S max is clearly clustered for all land cover classes.High optimal values are found in forested and highland regions, while riverine areas and dambos show relatively low values for the optimum of S max .Forested model units show least clustering, although the parameter ranges suggest that S max should at least be higher than 1000 mm.For 2 land cover classes, dambos and highlands, the parameter response surface reveals a clear clustering of optimal parameter ranges for both parameters.Forested and riverine areas show less clustering for l p .

Physical interpretation of parameter posteriors and validity of the model structure
Although the parameters represent area-averaged and thus effective values, the relation of the parameter response surface with the identified land cover follows the knowledge we have on the land cover regime.For instance, the dambo dominated areas do not have deep rooting vegetation, since water levels are too shallow in the wet season for deep rooting trees to survive.Shallow rooting grasses, reeds and shrubs dominate these areas and it is well known that the grass wilts very soon after the rains have passed.This explains the rather persistent optimal value of l p for all 3 dambo model units being close to 1 (Fig. 5d).It implies that transpiration rates decline immediately when S u <S max .A relatively low optimal value for S max (Fig. 5c) is found which also concurs with the small root depth of the dambo-type vegetation.Riverine areas show more or less a similar pattern (Fig. 5a and b), although the amount of model units used to sample from the parameter space was limited to only 2 and the parameter clustering is less pronounced.This is clearly visible in Fig. 5b, where it is clear that l p still shows variability.An obvious reason for this is that not all heterogeneity in soil type and vegetation coverage within the chosen dominant land cover classes is accounted for.
The highlands on the other hand, are covered with (partly evergreen) forests that apparently have a large reservoir of water available for transpiration.l p is rather low (Fig. 5h), indicating drought resistance, and S max (Fig. 5g) is high in these areas (2000 mm or even beyond the prior range).There is no optimum found for S max and we believe that a greater prior parameter range would not yield any better results given that the annual rainfall is in general much less than the maximum prior value for S max of 2000 mm.This means that as S max becomes higher, modelled evaporation becomes less sensitive for variability in S max .The highlands are located on the rather isolated Muchinga escarpment, where the evergreen forests on the edge of the escarpment may act as a sink to which groundwater converges, perhaps even from outside the Luangwa catchment itself.It is therefore likely that the model structure applied, may not be suitable for these areas: first, the Miombo woodlands may tap from groundwater, and second, there may be a lateral influx of groundwater, which cannot be modelled with a 1-dimensional box-type model such as the one we present here.A perception could be that the model should be replaced by a 2-box configuration where soil moisture is replenished in the dry season by uptake from the groundwater and discharge is generated from a groundwater reservoir.There are however no measurements done in this area to support this hypothesis.
The posterior parameter distribution for forested regions (Fig. 5e and f), shows a bi-modal distribution in the combined response surface of the four forested model units.This dual mode may be due to the two main forest types in the basin.In the field, we have observed large areas covered with multiple species of Miombo (Brachystegia), some species suited to the hotter lower areas, some to the colder, higher elevated areas, exhibiting far less seasonality (Fuller, 1999).Some of the lower lying and hot areas are covered by Colophospermum Mopane, known to intersperse dambos (Chidumayo, 2005).Mopane is known to dominate areas with relatively shallow and poor soils that are not well drained (Lewis, 1991) whereas in areas with more favourable growing conditions, other (for instance Miombo) species will dominate.Miombo is known to root deeply and use deeper soil moisture or groundwater reserves.The reason for their dry-season dormancy may well be temperature related rather than soil moisture related (Chidumayo, 2005), which could mean that in these woodlands, we should include temperature as a transpiration constraint, which may lead to further model structure improvement.If for instance average temperature will rise in the future, these trees may keep on transpiring at the cost of deeper groundwater.Our misconceptualized model would actually stop transpiring in this future climate due to moisture constraints.It would require a far more detailed land cover map to identify what type of forest we are dealing with and what its coverage is.

Effect of bias on posterior parameter distribution
All previous analyses have been carried out under the assumption of no bias.Nevertheless, there are realistic possibilities that evaporation is sometimes biased.In particular over areas where thin cloud-cover may be expected, SEBAL can result in overestimation of evaporation due to the contamination of land surface temperature estimates by cloudiness.Immerzeel et al. (2008) for instance suspected SEBAL to provide overestimated evaporation over forested areas in India.Figure 7 shows the effect of such overestimation on the parameter inference of 1 forest-covered model unit, by multiplying SEBAL evaporation estimates with a multiplication actor β and consequently perform the Monte-Carlo analysis described in Sects.3.3 and 4.1.Although sets of behavioural parameters remain within the same part of the parameter space, it is clear that this bias has a significant effect on parameter estimation.A generally lower evaporation rate will reduce S max and increase l p , which can physically be interpreted as a lower capacity of retaining moisture and higher sensitivity to moisture stress.Bias may be accounted for if the annual water balance at the catchment scale is known (i.e. if annual precipitation and streamflow are known with enough accuracy).In case annual evaporation can be derived by SEBAL, this water balance allows for a bias correction by assuming that total annual evaporation should be equal to precipitation minus discharge.In this case, β can be imposed on evaporation estimates before taking them into account in the parameter inference under the assumption that the bias is equal over the catchment.However, in Southern Africa, SEBAL cannot be applied over the entire year due to persistent cloud cover during the rainy season.In this case, a better parameter estimation can only be found if additional information to constrain the water balance, in particular discharge, is used jointly with evaporation in the parameter estimation process.Then, β may be introduced as a latent variable, as for instance proposed by Kavetski et al. (2006).In this way bias can be estimated as part of the parameter inference process, while maintaining the temporal and spatial moisture depleting behaviour.

Outlook
The results of this study show that, even with limited groundtruth knowledge, remotely sensed evaporation can both reveal model structural deficiencies and condition model parameters.When observing a natural river basin from above, patterns can be observed that are the result of evolution, which has resulted in a co-existence of ecosystems and hydrological behaviour.As a result, there is interdependence between vegetation, evaporation and runoff.This interdependence can clearly be identified in our model structures, which represent a simplified perception of nature, for instance in the parameters S max and l p in our conceptual representation of the soil moisture zone.These parameters influence both runoff and evaporation a great deal.It means that, although we condition the range of possible parameter realizations in a period without any forcing in terms of rainfall, it will have a great impact on parameter realizations in terms of discharge as well.The wide range of likely parameter realizations, found for forested regions, also reveals that we have not yet learned enough about the land surface to condition our models reliably on these data.Besides uncertainty in the evaporation data, there are still a lot of hypotheses on which we may further condition our model structure and parameters (e.g.difficulties in identification of the type of trees within forested model units, groundwater convergence in the highlands, missing knowledge about the geology, redistribution of runoff in wetlands).Nonetheless, it is far better to use these posterior distributions as a soft constraint (Seibert and McDonnell, 2002) than not to use them at all.Although evaporation estimates from space may be noisy and result in discontinuous response surfaces (as seen in Fig. 5), they offer one of the few, maybe even the only, opportunities to justifiably distribute parameters in ungauged basins.This would improve the decision making that can be done based on such hydrological models.Even as a soft constraint, for instance in the form of a fuzzy measure that may be employed as a first measure of degree of acceptance, these posteriors will impose a constraint for distributed land-surface related parameters in basins with little gauging.In this case, the red trapezoidal functions, displayed in Fig. 5 could be used as fuzzy measures of acceptability of priors in a subsequent calibration step on for instance (old) streamflow records.The form of such fuzzy measures is somewhat arbitrary, but the trapezoidal shape makes sure that models on the left or right slope are not immediately rejected.Instead, the trapezoid can be used to give a penalization score to models that are on the slopes.In addition, effects of bias in the evaporation estimates may be compensated for if multiple (water balance) constraints are jointly introduced in this parameter inference process.It is the view of the authors that when more and more of such soft constraints are included and combined with hard constraints in the modelling process, we will eventually be able to make better predictions in ungauged basins, including the necessary uncertainty assessments.

Conclusions
In this study, we presented a method with which remotely sensed evaporation can be employed to construct posterior parameter distributions of land-surface related parameters in hydrological models.A simple conceptual representation of the soil moisture dynamics has been used in order to present clearly identifiable parameter constraints.The results show that there is a clear consistency in the posterior likelihoods of parameters, belonging to the same land cover class.The consistent modes of the response surface are hydrologically meaningful in the sense that these modes concur with what one can expect from the land surface: landscapes covered by deep rooting vegetation reveal high optimal values for S max , the unsaturated zone capacity, and areas with shallow rooting vegetation, dambos and riverine areas, show much smaller values for S max .Furthermore, the dambo-dominated areas, typically covered with seasonal grasslands, are easily stressed for moisture, which corresponds with high values for l p .Because the evaporation estimates used are typically noisy, we feel that the posterior likelihoods found, can only be used in a "soft" mode for instance as a fuzzy measure of acceptability, combined with likelihood measures constructed from other information.
We argue that S max may reach unrealistically high values, especially for the forest covered highlands, which may point in the direction of deficiencies in the model structure: the deep-rooting forests are probably capable of tapping from groundwater sources during the dry season, which may be replenished laterally.This is a process that is not yet included in the model structure and which could lead to better process understanding if included in the model structure.
Furthermore, the response surface for forested area exhibits a bi-modal distribution.This teaches us that we probably have oversimplified the variability in forest types and their coverage.Multiple species of Miombo and Mopane forests co-exist in these areas where in general Mopane favours shallow and poor soils and Miombo has a deeper rooting system, dominating richer soils.This underlines the need for collaboration between hydrologists and ecologists.Collaboration may improve the understanding of the synergy between hydrology and ecology in large eco-systems and hence improve the spatial discretisation of our models.
The approach followed, enabled us to spatially distribute and constrain some crucial parameters that determine rainfall-runoff behaviour in an ungauged basin of considerable size without any direct calibration on streamflow and, interestingly, outside periods of direct forcing.Moreover, in semi-arid areas, where evaporation is a much larger water balance term than streamflow, this step should be combined with calibration on discharge in order to justifiably include spatial distribution of parameters that determine the water balance.The derived constraints will prove useful in ungauged catchments, especially if used in a Bayesian updating framework, combining multiple constraints.

Fig. 5 .
Fig.5.Likelihoods of S max and l p given the SEBAL evaporation estimates.The results from different model units with similar land cover (the amount of model units is given on the right side between brackets) are combined in one figure to reveal resemblances in the response surface.The red lines indicate possible trapezoidal fuzzy measures that could be applied in a later calibration step as parameter constraint.

Fig. 6 .
Fig. 6.Examples of the simulation performance of transpiration for well-performing parameter sets per land cover class.

Fig. 7 .
Fig. 7. Effect of bias on posterior parameter distributions of S max and l p

Table 1 .
Lindström et al. (1997)ls, their physical meaning and units.[L](completelyequal to the HBV soil zone), consisting of 3 parameters, in this paper referred to as S max [L], B [-] and l p [−], equivalent to the abbreviations "FC", "BETA" and "LP" in the publication byLindström et al. (1997).Outgoing fluxes are transpiration T a [L T −1 ] and runoff generation r c [−], computed as

Table 2 .
Parameter values of S max and l p and corresponding performance, belonging to simulations displayed in Fig.6.
Subsequently a Monte-Carlo framework was applied to estimate posterior parameter distributions for each model unit, given the model structure and given that the quality of the SEBAL data is reasonable and at least unbiased.Uniform prior distributions were imposed for both S max (varying between 200 and 2000 mm) and l p (between 0.3 and 0.95).A model run time from September 2005 (driest moment in the year) until October 2006 was used.This period covers one rainy season, which allows a spin-up of the soil zone for the model evaluation period, May 2006 until October 2006.The following objective function L was used: 2 (11)where t p is a time period[T]for which a SEBAL map was generated, m is the number of observation time periods available and E o and E s are the observed and simulated evaporation respectively [L T −1 ].L becomes lower as the model performs better.