Challenges of modelling soaring flight in humid landscapes A sensitivity analysis of vegetation parameters in the modelling of sensible heat flux

(1)

Challenges of modelling soaring flight

in humid landscapes

A sensitivity analysis of vegetation parameters in the modelling of sensible heat flux

(Source: http://klimaatadaptatie.nl/wp-content/uploads/2017/05/succesfactoren.jpeg)

BSc Future Planet Studies Earth sciences major

Thesis by Rens van der Veldt – 10766162 Supervisor: Willem Bouten

Secondary Supervisor: Emiel van Loon July 4th_{, 2017 Amsterdam}

(2)

Some species of birds use thermal convection to soar. Modelling this could help better understand bird behaviour. However, modelling sensible heat flux, a proxy for thermal convection, has not been done on high enough resolutions that provide insight into bird behaviour. By means of solving the energy balance (Formula 1) SHF can be calculated as a residual. For this purpose two models have been created on a suitable resolution of 100x100m and hourly resolution. These models for evapotranspiration and ground flow created by Sweerts (2016) and Borgman (2016) respectively, are combined to provide a model that calculates SHF. Unfortunately, the parameters concerning the land cover, calculated by Kooij (2016), contain uncertainties due to the dynamics of the landscape. Thus the aim of this thesis is to analyse how the energy balance behaves in various situations and how sensitive the integrated model is to changes in crop factors (Kc) and leaf area indexes (LAI) under these various conditions. It was found that Kc and LAI have an impact of at least 30% in the reference cells, a short period of drought was found to cause an increase in all sensitivity. This short period of drought did not cause noticeable amounts of water stress, but did decrease water fraction of the soil and thus thermal conductivity. In Veluwe forests the sensitivity is highest (>100%) due to the minimal role of ground flow (G). Well-drained soils are expected to have higher sensitivity due to a lower role of G, this effect was not found in any of the reference cells. This indicates that other mechanisms overshadow the effect of soils on the tested parameters as sensitivity is larger in polders. The sensitivities indicate the importance of using accurate land use data for parameterisation.

Keywords: energy partitioning, sensible heat flux, ground flow, evapotranspiration, modelling, crop factors, leaf area index, landscape dynamics, vegetation, sensitivity analysis.

(3)

Page | 3

Introduction

Soaring flight and modelling background

Many bird species use natural processes that cause air to “flow” upward in order to increase their flight efficiency. One of these processes is thermal convection, which is used by birds to gain height and to soar almost effortlessly. This so called thermal convection can occur in a landscape when there is energy left over after evapotranspiration and the heating of the soil, which is then instead used to warm the air (Delworth & Manabe, 1988). The energy used for warming the air is the sensible heat flux (SHF) which serves as a proxy for thermal convection (Bohrer et al, 2011). Although modelling SHF has been done separately on a daily basis (Bastiaanssen, 1999) or on a 1x1m resolution for an urban area of 0.650 km2_{(Voogt & Grimmond, 1999), these models remain unfit to assist in predicting bird behaviour. As}

bird behaviour is relevant at high resolutions, it is more meaningful to calculate SHF on higher temporal resolution and for larger and more diverse areas combined. In order to accomplish this, two models were created to model SHF based on the energy balance products of evapotranspiration (ET) modelled by Sweerts (2016), ground flow (G) modelled by Borgman (2016) and the measured net radiation (Rn). The evapotranspiration is the energy flow product of evaporation and transpiration combined. Ground flow is the energy flow that heats the ground. By calculating the energy fluxes of both the overlying atmospheric system and, the hydrological underlying system the energy balance can be solved (Formula 1).

𝑅𝑛 = 𝐺 + 𝐸𝑇 + 𝑆𝐻𝐹 [1]

By deducting the ET and G fluxes from the input radiation data the SHF is calculated. The models for evapotranspiration and ground flow were created on an hourly basis with a 100x100m spatial scale which are far more suitable for analysing soaring behaviour than the limited previous models.

The landscape dynamics

The focus of this paper lies in the fact that in humid landscapes soaring flight can be challenging. In these landscapes a high fraction of the available energy is spent on evapotranspiration, leaving little energy for the sensible heat fluxes required for soaring flight. The landscape which will provide context for this research is The Netherlands. In this country water is abundant and with large areas close to and even below sea level, groundwater is often in proximity to the surface. The humid, temperate climate also provides ample water in the form of precipitation. The large water availability means that most of the time vegetation may transpire unhindered by water stress. Furthermore, wet soils are better thermal conductors, channelling energy away from the surface during daytime. The result is that many different landscape properties play an important role in the modelling of SHF. Some elements of this landscape change constantly. If the model is to be used to accurately predict SHF it is useful to know how sensitive the model is to certain changes in this dynamic landscape.

Various changes may constitute the dynamics of the landscape. Vegetation changes throughout the year as in forests leaves fall in autumn and regrow during spring. In agricultural fields this change may be more abrupt when crops are planted and harvested after reaching adulthood. With the landscape becoming instantly barren, and with no transpiration happening, the SHF may increase sharply. These oscillations are captured by the leaf area index (LAI), a measure of how much of the ground is shaded by vegetation and influencing the SHF through co-determining how much plants can transpire. The amount of transpiration a plant can do is captured by the values of crop factors (Kc), this factor determines how much a certain vegetation cover can transpire. The result is that different crops may show various patterns in transpiration dependent on Kc and LAI and due to crop cycles the vegetation in agriculture can change much from year to year. To illustrate this dynamic, two landscapes captured by Google Earth, are compared. The first at the end of June, 2015 (Figure 1). The landscape is

(5)

Page | 5 mostly shown as a pallet of greens with only a few barren plots. At the exact same location by the end of August, 2012 (Figure 2), many of the fields have been harvested and the landscape turns into a greyish brown. It also becomes visible, that the division of many areas changed. For example, in the red highlighted area the border between two fields shifted between these snapshots.

Figure 1. Agricultural landscape at the end of June, 2015 in Flevoland.

Figure 2. Agricultural landscape at the end of August, 2012 in Flevoland.

Larger changes of the land use may also play a role. Cities that expand, agriculture that is turned into recreation or forests that get chopped down. This is emphasized looking at data from the Centraal Bureau voor de Statistiek (2017), which states that over 1.214 square kilometres of agricultural land

(6)

Page | 6 disappeared between 1996 and 2012. Of this area, 642 square kilometres was transformed into building or industrial terrain, housing areas, or roads. These areas constitute the lowest range of values in Kc and LAI, due to absence of vegetation, therefore only little can be transpired which positively affects the SHF.

The described dynamics would not be problematic if parameterisation was dependant on frequently updated data of the landscape. However, the state of the art land use data that the vegetation parameters are dependent on is currently updated every 6 years. This data is provided by the EEA (European Environment Agency) in the CORINE land cover (CLC) datasets, which exist in the required resolution of 100x100m. Using this data, the landscape was parameterized by Kooij (2016). This poses a problem in relation to the described dynamics of the landscape and vegetation. As LAI and Kc may both change with changing vegetation cover they become uncertain with dynamic land uses.

The model is dependent on many other factors that together calculate the SHF. The goal of this paper is to find out how important the inclusion of accurate vegetation parameters are in the model. This is done through analysing the sensitivity of the landscape energy balance to the LAI and Kc parameters in a spatial and temporal high resolution model. Because other factors also play a role this is done for varying landscapes which will be elaborated on in the methods. In doing this it is desired to find that the model is not sensitive in all situations. This would indicate that in order to gain viable results from the model, the landscape data to be put in would not have to be revised all the time.

Research questions and structure

The research questions aim to explain the behaviour and sensitivity in different situations. The first research question is how the agricultural, forest and urban areas’ energy balance is influenced by varying conditions. Secondly, this research will focus on how sensitive SHF in these areas is under changing vegetation parameters. These question will help answer how important it is to include accurate vegetation parameters for the calculation of SHF in this humid climate.

(7)

Page | 7

Methods

Theoretical framework

For the theoretical framework the theories behind the two models created by Sweerts (2016) and Borgman (2016) will be discussed. Furthermore, the theory on which parameterization was performed by Kooij (2016) will be introduced. Furthermore, the abbreviations found in this chapter are in line with variable names in the model, thus also serving as a documentation for the model. If available, formulas are accompanied by references to the function names in the supplementary material. Finally, it must be noted that most of the theoretical background is retrieved from the theses by Sweerts (2016), Borgman (2016) and Kooij (2016) and this chapter serves as a note on their theoretical descriptions.

Evapotranspiration

The evapotranspiration model uses the potential or reference evapotranspiration (ET0) for a cell in combination with reducing factors to determine the actual evapotranspiration (ETact). ET0 is calculated based on an adaptation of the Penman-Monteith method, done by the American Society of civil engineers (Allen et al.,1998)(Formula 2).

𝐸𝑇0 =0.408𝛥(𝑅𝑛−𝐺)+ 𝛾 𝐶𝑛

𝑇+273.3 𝑈2 (𝑒𝑠−𝑒𝑎)

𝛥+ 𝛾(1+𝐶𝑑 𝑈2) [2]

In this formula Δ is slope vapour pressure (d) and is calculated using “calc_d.m” (Formula 3). 𝛥 = 4098 (0.6108𝑒

17.27𝑇 𝑇+237.3)

(𝑇+237.3)2 [3]

Rn (Rn) is the incoming solar radiation, G is the soil heat flux (GroundFlow), which will be discussed in a separate paragraph about its theoretical background. Furthermore, (γ) is the psychometric constant, T is the temperature at 2 meters (T2m) converted to Celsius, Cn (Cdn and Cn) are denominators for different vegetation heights and day/night cycles. U2 Is the wind speed at 2 meters (U2m) which was

converted from 10 meter wind speeds using the FAO conversion curve. Es and Ea (Es, Ea) are the saturation vapour pressure and actual vapour pressure (Formula 4 and 5) in “calc_es_ea.m”. These formulas use the same temperature at 2 meters in combination with dew point temperature at 2 meters (D2m).

𝑒𝑠 = 0.6108 𝑒 17.27 𝑇

𝑇 + 237.3 [4]

𝑒𝑎 = 0.6108 𝑒 17.27 𝑇𝑑𝑒𝑤

𝑇𝑑𝑒𝑤 + 237.3 [5]

To get from ET0 to ETact, it is multiplied by reducing factors (Formula 6). The first being a crop coefficient (Kcb) according to the FAO dual crop method. There are Kcb values for midseason (KcbMid) and full grown crops (KcbFull). The calculation of these was done by Kooij (2016) and will be explained in the landscape parameter section.

𝐸𝑇𝑐 𝑎𝑐𝑡 = (𝐾𝑠𝐾𝑐𝑏 + 𝐾𝑟𝐾𝑒)𝐸𝑇0 = 𝐾𝑐𝑎𝑐𝑡 𝐸𝑇0 [6] The Kcb values are reduced by multiplying with a water stress coefficient (Ks) (Formula 7). TAW is the total available water in the root zone and Dr the water shortage relative to field capacity (Dr). The 𝑝-values determine the volumetric wilting point for a type of vegetation (Fpwp). Ks is calculated in the model through the function “calc_ks.m”.

(8)

Page | 8 Furthermore, the soil evaporation factor (Ke) is calculated (Formula 8). Ke is determined through subtracting the current vegetation coefficient which is KcbMid from the value at maximum growing conditions, KcbFull.

𝐾𝑒 = (𝐾𝑐𝑚𝑎𝑥 − 𝐾𝑐𝑏) [8]

This Ke is reduced by the corresponding soil evaporation reduction factor (Kr) (Formula 9). The total evaporable water in the soil at field capacity (TEW) and the water that can be evaporated before water stress occurs (REW) are used in combination with the soil moisture depleted the last day (De,j) in the function “calc_kr.m”.

𝐾𝑟 =𝑇𝐸𝑊−𝐷𝑒,𝑗−1

𝑇𝐸𝑊−𝑅𝐸𝑊 [9]

The final step is combining the calculated factors with ET0 (Formula 6). This is done through the function “calc_etact_evap.m”, providing us with the required evapotranspiration (ETact).

Ground flow

The ground flow model created by Borgman is based on the method proposed by Bouten (2016), seen in formula 10 and functions somewhat differently from the evapotranspiration model. This model works with an extra “dimension” where layers are introduced, using a factor (Ks) multiplied by the rate of ground flow (GroundFlow) in every layer.. From this, the exchange between surface and ground is the net flux of energy, or ground flow.

𝐶 ∗ 𝜕𝑇

𝜕𝑡= 𝐾𝑠∗ 𝜕2𝑇

𝜕𝑧2 [10]

C, Volumetric heat capacity (VHC) is determined through the soil characteristics such as the VHC of mineral content and their fractions, which are provided in the model and calculated through the function “calc_vhc.m”. 𝜕𝑇 is the change in temperature over 𝜕𝑡, the change in time. These are the rates for each soil layer with height z (ThickL), as found in the model configuration. The actual heat content (HeatCont) in each layer is calculated (Formula 11) as proposed by Bouten (2016).

𝐻𝑐 = 𝑇 ∗ 𝐶 ∗ 𝑧 [11]

The gradient in heat content drives the flows between layers. The calculation in the model is done through the function “calc_flows.m”. Fluctuation in the heat content is co-determined by the boundary conditions present in the model. These conditions allow energy to leave and enter the system, thus causing flow. For the bottom boundary condition, the lowest layers temperature is set to 0. The initial top layer is set to match the temperature at 2 meters (T2m).

The calculation of flows requires some preparation, one of the main drivers here is the net radiation (Rn) which is multiplied by a cover fraction (CoverFraction) (Formula 12).

𝑐𝑓 = exp (−0.5 ∗ 𝐿𝐴𝐼

𝑍𝑒𝑛𝑖𝑡ℎ) [12]

In the model this is done using the weekly values for the leaf area index (LAI) and the constant for Zenith (Solar_Zenith) (Yang et al., 1999) in the function “calc_cover_frac.m”. Flow in heat content is also influenced by aerodynamic resistance of a landscape, such can occur due to the vegetation layer, which is aerodynamic resistance (Rh) is introduced as another boundary condition. The aerodynamic resistance values from Goudriaan (1977) were assigned in the model to the corresponding classes of

(9)

Page | 9 the Corine Land Cover (CLC) map. This is done in the model in the function “convert_aero_resistance.m”.

𝐻 =𝜌𝐶𝑝 (𝑇𝑎𝑒𝑟𝑜− 𝑇𝑎)

𝑅𝑎ℎ [13]

The resulting Rh can be used to finally create the aforementioned boundary condition (Formula 13) from Zuang et al. (2015). However, in the model itself the aerodynamic dynamic is simply used as a boundary condition in reduction of the top layers’ temperature flow. This is conducted during the calculation for flows.

Landscape parameters

The landscape parameters concerning vegetation are the Leaf area index (LAI) and Kcb values (KcbMid and KcbFull). LAI was calculated by Kooij (2016) using the method proposed by Hagemann et al. (1999) (Formula 14).

𝐿𝐴𝐼𝑖 = 𝐿𝐴𝐼𝑚𝑖𝑛+ 𝐹𝑖∗ (𝐿𝐴𝐼𝑚𝑖𝑛− 𝐿𝐴𝐼max) [14]

Where LAI in week i is determined by its maximum and minimum value and a factor Fi which is the

growth factor at week i. This factor is determined by temperature data. The LAI can be used to calculate the crop factors KcbMid and KcbFull (Formula 15-16), according to the FAO method from Allen et al. (1998).

𝐾𝑐𝑏𝑚𝑖𝑑 = 𝐾𝑐𝑚𝑖𝑛+ (𝐾𝑐𝑏𝑓𝑢𝑙𝑙∗ 𝐾𝑐𝑚𝑖𝑛) ∗ (1 − 𝑒−0.7 𝐿𝐴𝐼) [15]

𝐾𝑐𝑏𝑓𝑢𝑙𝑙 = 𝐾𝑐𝑏ℎ+ (0.04 ∗ (𝑈2− 2) − 0.004 ∗ (𝑅𝐻𝑚𝑖𝑛− 45)) ∗ ( ℎ3)

0.3 _[16]

KcbMid contains mid-season crop coefficients, which are due to the implementation of LAI, applicable through all seasons. Kcbh is a proxy for KcbMid derived from maximum plant height and U2 the wind

speed at 2 meters (U2m). RHmin is the 30-year average minimum daily humidity during the mid-season.

(10)

Page | 10

Data and tools

The model was written in Matlab version R2016b on a 64 bits Windows version using 12 GB RAM. The integrated model was constructed according to the models created by Sweerts (2016) and Borgman (2016), also in Matlab. The model was used with data spanning May and June of 2012. This data was originally retrieved from multiple sources (Table 1). Most of its content is retrieved from the theses of Sweerts (2016) and Borgman (2016). As some data resolutions, both temporal and resolution wise, did not match the provided data in for the models they were edited here to match the models’ input.

Source Data Spatial resolution Temporal resolution

ERA‐interim

(reanalysis) by ECMWF

-Temp at 2m

- Dew point temp at 2m -Wind speed at 10m -Total cloud cover

0.0125° (Degrees) 1-Hourly

NHI by Deltares -Soil moisture content -Soil moisture deficit -Root zone depth -Soil properties

250 m Daily

PCR-GLOBWB (model) by University of Utrecht

-Soil moisture content 5’ (geographical minutes)

Daily

CORINE land cover programme by EEA -Land use 100 m - Satellite Application Facility on Climate Monitoring (CM SAF) -Surface incoming shortwave radiation 0.05° (Degrees) 1-Hourly

GLCCD by USGS -Global ecosystems map 1 km -

European Soil Data Centre

-Soil texture 1 km

(11)

Page | 11

Workflow

The workflow diagram (Figure 3) helps to explain the steps taken from data loading and manipulation, to the modelled output. In this section of the methods the workflow diagram will be explained step by step.

Figure 3. Workflow diagram.

Adjustments

The first step is loading the data, control constants and control points. The following primary data adjustments are applied first. Wind speed data was converted from 10 to 2 meter in the script according to the FAO conversion curve (Allen et al., 1998). Incoming radiation was adjusted from mega joule to joule. Furthermore, European scale data is clipped to match Dutch spatial dimensions. This converts the model input to the proper units and size.

(12)

Page | 12

Interpolation and area selection

The next step is interpolation, as resolutions do not yet match (Table 1). The interpolation was done so that all data matches the CORINE land cover data resolution of 100x100m. This is done using two separate functions that simplify readability of the script. Interp_3d_data and interp_2d_data, for data with that changes over time (3d) and two dimensional data respectively. The resulting interpolated maps have a 100 by 100 meter resolution and are now resolution wise uniform.

The last part of the preparation is the selection of an area. This step is optional and included in the main script in the supplementary material. It allows the user to select an area from the available data for which the model will be ran. This can drastically improve runtime as all data is clipped to this size. The model functions by initializing the data and variables used for dynamic calculation, listed in the workflow diagram.

Integration

The combining of models for ET and G by Sweerts (2016) and Borgman (2016) respectively was realized in order to calculate SHF based on the energy balance equation (Formula 1). The models were restructured and split into the functions mentioned in the theoretical framework. This can be found in “evapotranspiration_functions” and “groundflow_functions” in the supplementary material. Furthermore, all hardcoding was taken out of the model input.

After re-writing, the integration itself became fairly simple. The dynamic calculations part was split into a hierarchical time loop, meaning the main daily loop contains 24 hourly loops to calculate temporally more detailed data. This is different from the original model structures as these calculate many variables outside of such loops. The hierarchical structure of time loops was introduced bearing in mind that the interception module as originally created by Mol (2016) and edited by Zijlstra (2017), with a 5 minute time step, could be easily included. This is illustrated in the time hierarchy structure (Figure 4). Furthermore, in regards to the calculation sequence, G was calculated before the ET0. This is because G is used in the adapted Penman-Monteith method. This provides the model as seen in the model workflow (Figure 3).

(13)

Page | 13

Period of analysis

Analysis of the model is conducted by running the model for multiple cells of interest. This is done for both dry and wet periods to assess the increase or decrease of sensitivity with changing conditions. During May, 2012 De Bilt experienced a period of relative drought between the 24th_{and 31}st_{(Figure 5)}

after intense rainfall. This will be the period for which the model is ran.

Figure 5. Rainfall in De Bilt (mm, left) and duration (minutes, right) (KNMI, 2012).

To see whether this drought is expressed in the data the decrease in root zone saturation between the 25th_{and 30}th_{of May is generated (Figure 6). Unfortunately, the expected decrease is not}

represented everywhere. This is because many areas will still be provided with water by either natural conditions, such as high groundwater levels or retention rate, but also in the Dutch polders, where water levels are often regulated by humans. Other areas do show a decrease in soil moisture, such as the Veluwe and the Utrechtse Heuvelrug. These moraine areas were formed by glacial activity and contain well drained sandy soils (van Buuren, 2004). This explains how soils and location may also influence moisture availability. However, some vegetation cover also possess larger root zones, reaching into deeper water and allowing transpiration to be maintained longer during drought.

(14)

Page | 14 Figure 6. Moisture decrease between 25th_{until 30}th_{of May in mm per square meter of the soil.}

Analysis

The model is ran for day 24 until 30, allowing one calibration day. The model analysis will be focussing on cells at the end and start of the discussed dry period. For the duration of this period of drought, agricultural fields and forests in two contrasting areas are analysed: the Veluwe and Flevoland, these are contrasting areas of soil properties and moisture availability. Finally, the situation for cities is discussed. This amounts to a total of 5 reference cells to be analysed. A sensitivity analysis can be used to quantify the influence of landscape related parameters. Global sensitivity analysis in Matlab requires data to be uniformly distributed. However, the parameters of interest were not (Chi-Squared GoF, alpha = 0.05). Another attempt at setting realistic ranges was to use a confidence interval. However, Kolmogorov Smirnov tests for normal distribution (alpha = 0.05) also proved insignificant. Without any statistic grounds to set the range from, this forces the range to be the maximum and minimum of parameter values. These values influence the modelled G and ET to a certain degree. By using this G and ET to calculate SHF the sensitivity of the model in this thesis can be defined as: the percentage difference in SHF between the modelled SHF and the SHF with the maximum and minimum parameter as input. This will allow for drawing conclusions about the importance of incorporating accurate landscape factors in a humid landscape in different settings, thus answering the final research question.

(15)

Page | 15

Misprojection

After looking at the data, a misprojection problem was identified. The data that was retrieved from PCR-GLOBWB was in a different projection (Figure 7). More specifically, most data was provided in the European Terrestrial Reference System 1989 (ETRS89) also known as EPSG:3035. The root zone saturation and saturation deficiency, total evaporable water (TEW) and total available water (TAW), are projected in the World Geodetic System 1984 (WGS84) also known as EPSG:4326. Without spatial reference in the PCR-GLOBWB dataset, this data had to be reprojected for the model results to make sense. This was done manually using the Geometric transformation toolbox. Geo-referencing was done using the Control Point selection toolbox. The total available water map and ρ-values (Fpwp/Wilting point values) were used with this toolbox to select the points required.

Figure 7. Original misprojection shown as difference between map cover in both blue and yellow. Because the difference in projection is “curved”, and not simply rotated, a 3rd_degree

polynomial transformation was applied using 30 control points set at distinct landscape features. After the transformation the images have to be clipped in order to match the other data. The data boundaries on the all sides but north can be used directly to clip the maps. The northern boundary however does not contain the Wadden islands in the transformed data. The Northern border clip size was eventually estimated visually.

(16)

Page | 16 Figure 8. Misprojection after reprojection shown as difference between map cover in both blue and yellow.

After reprojection the match between the data has improved (Figure 8). Most of the difference in projection that seems to remain can be explained by the difference in resolution and the imperfection of the applied technique. This partly answers the first research question if the models can be

(17)

Page | 17

Results

The first result is the average energy balance for a large section of the Netherlands (Figure 9). The different areas were distinguished based on the CLC classes. The urban environment is distinguished by CLC classes 1-11. The forest environment is limited to natural environment classes at CLC classes 23-25. And in the middle, agricultural areas, classes 12-22.

Figure 9. The energy balances for different areas on the 25th of May.

For two different cells with agricultural landscape classes, the cells are located in Flevoland and the Veluwe. The energy balances were generated both before (Figure 10) and after drought (Figure 11). The radiation data for the day after drought contains a dip due to clouds moving over the area.

(18)

Page | 18 Figure 11. Energy balance after drought in agricultural fields on different soils, the dip in radiation is caused by clouds moving over the landscape.

This result (Table 2) shows percentage change in SHF when Kc and LAI are the minimum value, similar to the urban situation, on the left-hand. SHF increases because ET decreases with lower parameter values. The right-hand number is the change that happens when Kc and LAI are at a maximum, similar to a forest at peak transpiration.

SENSITIVITY AGRICULTURE BEFORE DROUGHT AFTER DROUGHT

VELUWE SHF % CHANGE +30.52 -31.21 +35.49 -42.04 POLDER SHF % CHANGE +39.98 -37.69 +39.16 -56.39 Table 2. Sensitivity of the agricultural cells before and after drought on different soils.

(19)

Page | 19 The energy balance for two cells with forest land cover classes in Flevoland and the Veluwe were generated for the same period as the agricultural cells. Both before (Figure 12) and after (Figure 13) drought are shown together with the relevant sensitivity analysis (Table 3).

Figure 12. Energy balance before drought in forests on different soils.

Figure 13. Energy balance after drought in forest areas on different soils, the dip in radiation is caused by clouds moving over the landscape.

(20)

Page | 20

SENSITIVITY FOREST BEFORE DROUGHT AFTER DROUGHT

VELUWE SHF % CHANGE +68.41 -4.92 +104.77 -8.36 POLDER SHF % CHANGE +39.90 -40.28 +43.53 -54.93 Table 3. Sensitivity of the forest cells before and after drought on different soils.

The final energy balance that was generated (Figure 14) from the model shows the values for a cell in Utrecht, an urban area, which serves as a reference to what a situation with naturally low Kc and LAI looks like.

(21)

Page | 21

Discussion

Energy balance in various situations

The energy balances for the different areas (Figure 9), urban, agricultural and forest, show how different land use can lead to varying energy partitioning. The urban area constitutes the lower range of LAI and Kc values, which leads to the expected result of low ET and high SHF. The higher range values are found in the forests, as expected the ET is thus much higher on average. In agricultural areas divergent values of both high and low LAI and Kc are found, leading to an ET that seems average between the urban and forest situation. However, it also becomes visible that G is much higher. This could be due to higher moisture content leading to increased thermal conductivity. Although this may be temporary and due to the intense precipitation event that preceded this period, it is not visible in the other areas. Because SHF is reduced by ET and G, the forest area, with high ET and low G, has lower SHF than the agricultural area, with average ET and G.

In the compared agricultural areas (Figure 10 and 11) an important difference is noted. Before drought occurs ground flow is higher in the polder which causes the difference in SHF between the two areas. This could be explained by the fact that the Veluwe, as mentioned before, contains sandy soils which are well-drained, this causes saturation in the soil to be lower in the Veluwe. Therefore the thermal conductivity and thus the ground flow is lower in the Veluwe. As this is the situation before drought, both areas have ample water available and water stress levels (Ks) are very low (<0.002). After the drought Ks does not increase noticeably (Figure 11), this means that a week of drought in the model is not enough time to cause sufficient water stress to be noticeable in the agricultural ET. Moreover, ET has even increased. This was expected in the polder, due to increasing soil moisture (Figure 6), however it is also seen in the Veluwe. This increase can be due to any number of reasons but is possibly caused by the natural increase of LAI and Kc as after each week passes both increase. The impact of drought on G, in contrary to ET is more evident. As the period before drought was preceded by intense precipitation, soil saturation was high. This caused thermal conductivity to be high, which decreases a lot after drought. In contrast, soil moisture increased in the polder, but also shows the same decrease in G. With G decreasing in both areas the relative influence of ET on SHF grows over this period of drought.

In the forested Veluwe area (Figure 12) the ET is much higher than in the agricultural part. Although G is decreased by naturally higher LAI, preventing large fractions of the radiation from reaching the surface. Due to the well-drained soils and less radiation reaching the soil the role of G is relatively very low. Looking at the energy balance of the forest in Flevoland, it is visible that ET is not as high and very similar to the agricultural area. This is unexpected as in forests LAI and Kc should be relatively high. The high thermal conductivity in pre-drought conditions might overshadow the LAI, or the area might not be as densely forested as its Veluwe counterpart. The first explanation would mean that the fact that the soil in the Veluwe is well-drained plays an important role in the sensitivity to changes in ET. After the period of drought (Figure 13) the same pattern occurs as for the agricultural areas. ET increases and G is reduced, increasing the relative influence of ET on SHF. The reduction for the Veluwe forest, however, is barely visible due to the originally low values.

The urban energy balance (Figure 14) has much higher SHF than any other of the environments. This corresponds to the assumption made in the first paragraph that in cities the values of Kc and LAI are at their lowest. However the difference between this energy balance (Figure 14) and the first energy balance in the result (Figure 9) is that G behaves much more expectedly. The pattern seen in all other areas is that G increases as the radiation rises and decrease over the afternoon. This is inverted in the first urban energy balance. The primary explanation would be the inversion of the flux, as directions of

(22)

Page | 22 these can be confused. However, the agricultural and forest balances (Figure 9) were generated from the same data. This behaviour could still be explained by the fact that it was generated as an average of all the urban classes in the CLC map. The second urban energy balance (Figure 14) shows more expected behaviour, generated from the city centre of Utrecht. Although it is out of the scope of this research it would be interesting for the understanding of the model dynamics to know which classes are responsible for this pattern.

Sensitivity

The analysis proposed in the methods was applied to the compared agricultural and forest areas, both before and after drought. This resulted in the sensitivity of SHF, expressed as the percentage that SHF goes up or down due to change of LAI and Kc to both their lowest and highest value. Broadly speaking, high values of LAI and Kc cause increased ET, and thus decreased SHF. The sensitivity rises when the relative influence of ET increases. The influence increases when G is relatively small to ET.

Contrary to the expectations, the SHF in the agricultural Veluwe is less sensitive to the vegetation parameters despite the smaller influence of G than the polder (Table 2). This indicates that in the Veluwe other conditions limit the influence of Kc and LAI, even overshadowing the fact that the relative influence of ET is lower to begin with. With further declining influence of G over drought, this sensitivity does increase however as is expected with declining conditions for thermal conductivity in the soil. A possible cause for this behaviour is that in Flevoland the sensitivity is larger due to higher moisture availability, which allows vegetation with high values of Kc to transpire more than in the Veluwe. This does not account for the larger sensitivity in the lower range of Kc and LAI however. The lower values cause an increase in SHF over drought however, the sensitivity for these values remains around 39%. This indicates that the moisture conditions in the polder do not influence the SHF while Kc and LAI are at their lowest. What does become clear is that with a two-sided influence range between 37% and 56%, even as drought occurs, the polder is more sensitive to changes in Kc and LAI. Between the Veluwe and Flevoland, the latter is also the area in which agriculture is the most dominant form of land use. Because of the dynamics and abrupt changes in agricultural areas it becomes that much more important to incorporate accurate vegetation parameters into the model.

The sensitivity in the Veluwe’s forest environment is much more in line with expectations. The sensitivity analysis (Table 3) shows that with the highest parameters the Veluwe SHF changes only slightly. This is due to the fact that these parameters are already naturally high in these areas. With low Kc and LAI however, the sensitivity increases dramatically. Especially after drought, when ET increases and the influence of ET is the highest, SHF can be 104% higher. This can occur when the area is deforested for example. The sensitivity in the polder before drought is more evenly spread and similar to the agricultural polder. Moreover, sensitivity is much more two sided as Kc and LAI values are naturally lower than in the Veluwe. This even spread causes the cumulative sensitivity in the polder to be slightly higher before drought. In the after situation, the difference is relatively smaller. The relative role of G has decreased in the polder, resulting in two sided increase in sensitivity. The unexpected behaviour of the polder forest could be caused by the method of selecting reference cells. Although the analysed cell might be categorized as a forest, its actual parameters might not always be in line with such environments. If the role of other parameters that might cause this effect could be ruled out, this would emphasize the importance of accurate parameterisation.

Integration

Much more time than expected went into the integration process due to some unforeseen problems. Primarily, the integration would have been done in cooperation with other members of the research group, which was eventually confined to this research. Furthermore, the problem with misprojection

(23)

Page | 23 (Methods, Misprojection) occurred. For a more broadly orientated study and more accurate results it would be interesting to develop integration further. This should be done by integrating the model developed by Zijlstra (2017). Creating an even more realistic model as it incorporates interception evaporation. To make a step further towards predicting bird behaviour, the model for orographic lift developed by Hoekstra (2017) should also be integrated.

Better metadata

Some concerns were raised by the fact that the data provided for the model lacks metadata. Unfortunately, the current data can hardly be provided with sufficient metadata within the timespan of this thesis. The data collected from PCR-GLOBWB lacks referencing for instance. Without geo-referencing this had to be done manually, with its inherent imperfections. For further research new data collection would be advised.

Higher temporal resolution

Furthermore, in the current model the ground flow often jumps at the start of the day. Inside the scope of this research it is impossible to actually find the parameter or mechanism that causes this. However, this is probably due to the fact that many soil properties are calculated on a daily basis, although the changes that happen over a day are gradually. Although the early day is not as interesting for SHF, the model would be more reliant if more variables were calculated.

Conclusion

It was concluded that the analysis period of one week was insufficiently long to cause water stress and reduce ET. Due to drought the soil thermal conductivity, and therefore G, decreased. In both agricultural areas this causes the relative influence of ET on SHF to grow over this period of drought, this drought thus increases sensitivity. Although G decreases more in the Veluwe, the sensitivity in Flevoland is higher. Higher moisture availability for transpiration in the polder might cause this, though other mechanisms might be at play. In the forest areas of the Veluwe sensitivity is at its highest (>100% in the reference cell), this is due to the small role G plays and the naturally high LAI and Kc. G might be low because leaves shade most of the ground. In the urban area, ET and G are naturally low and SHF is higher than in the other areas. The low values of Kc and LAI reflect the situation that occurs due to urbanisation, deforestation, but also upon harvesting crops. This emphasizes the importance of using accurate land use data for parameterisation.

From the analysis it has become clear that the modelled SHF is sensitive to changes in these parameters. The sensitivity increases over drought and from well-drained to polder soils. Forests are more sensitive to lower range parameters, while agricultural areas are more two-sided. For better understanding of the mechanisms behind the sensitivity more research has to be done with more parameters. Moreover, the research would have to be conducted over longer periods of drought and potentially with higher resolution data and calculations. It has become clear that the sensitivity of SHF to Kc and LAI in the analysed areas can be influenced by various conditions. Nevertheless, in each of the areas sensitivity larger than 30% remains, indicating the importance of accurate vegetation parameters.

(24)

Page | 24

References

Allen, R.G., et al. (1998). Crop evapotranspiration – guidelines for crop water requirements – FAO irrigation and drainage paper 56. FAO.

Bastiaanssen, W. G. M. (1999). SEBAL-based sensible and latent heat fluxes in the irrigated Gediz

Basin, Turkey. Journal of Hydrology 229 (2000) 87–100.

Bohrer, G. et al. (2011). Estimating updraft velocity components over large spatial scales: contrasting

migration strategies of golden eagles and turkey vultures. Ecology Letters, 15(2), pp. 96-103.

Borgman, C. (2016). High Resolution modelling of the Soil Heat Flux with the inclusion of a Vegetation

Layer.

Bouten, W. (2016) Modelling and Simulating Spatial Processes in Landscape Ecosystems. Course manual for Modelling and Simulation, BSc major Earth Sciences, Universiteit van Amsterdam. Buuren, van, S. (2004). Grondsoorten in Nederland. Deltawerken online: Geology of the Netherlands. Retrieved from: http://www.deltawerken.com/Geology-of-the-Netherlands/112.html

Centraal Bureau voor de Statistiek (2017). Minder landschap, meer bebouwing, 1950-2015. Retrieved from: https://www.cbs.nl/nl-nl/nieuws/2017/10/minder-landschap-meer-bebouwing-1950-2015 Delworth, T.L., Manabe, S. (1988). The Influence of Potential Evaporation on the Variabilities of

Simulated Soil Wetness and Climate. American meteorological Society.

Goudriaan, J. (1977). Crop micrometeorology: a simulation study (No. 551.5 G688). Centre for Agricultural Publishing and Documentation, Wageningen.

Hagemann, S., Botzet, M., Dümenil, L. & Machenhauer, B. (1999). Derivation of global GCM boundary

conditions from 1 km land use satellite data. Max-Planck-Institute for Meteorology, Report 289,

Hamburg, Germany.

Hoekstra (2017). T.B.A. Submitted to UvA.

Klimaatadaptatie (2017). Klimaatadaptatie. Retrieved from: http://klimaatadaptatie.nl/wp-content/uploads/2017/05/succesfactoren.jpeg

KNMI (2012). Archief maand/seizoen/jaaroverzichten. Retrieved from:

https://www.knmi.nl/nederland-nu/klimatologie/maand-en-seizoensoverzichten/

Kooij, D. (2016). Parameterization of the Corin Land Cover dataset for the calculation of the energy

balance on a high spatial resolution.

Mol, M. (2016). Interception and the sensible heat flux.

Sweerts, B. (2016). High resolution modelling of evapotranspiration and sensible heat flux in the

Netherlands: Explorative research for a global thermal convection model for predicting bird behaviour.

Voogt, J.A., Grimmond, C. S. B. (1999). Modeling Surface Sensible Heat Flux Using Surface Radiative

Temperatures in a Simple Urban Area. American meteorological Society, Journal of Climate volume 1

issue 5.

Yang, Z-L., Dai, Y., Dickinson, R.E., Shuttleworth, W. J. (1999). Sensitivity of ground heat flux to

(25)

Page | 25 Zijlstra (2017). T.B.A. Submitted to UvA.

Zuang, Q., Wu, B., Yan, N., Zhu, W., Xing, Q. (2015). A method for sensible heat flux model

parameterization based on radiametric surface temperature and environmental factors without involving the parameter KB-1. International Journal of Applied Earth Observation and Geoinformation,

42, 50-69.

Acknowledgements

This research is done in cooperation with the UvA-BiTS (Bird Tracking System) program, which tracks birds in order to understand and study their migration, foraging and navigation behaviour. This cooperation was achieved by means of supervision by Willem Bouten. Furthermore, the research group in which cooperation this thesis was developed consists of 3 Future Planet BSc students. Zijlstra looked at the importance of incorporating high resolution precipitation data and interception evaporation into the calculation of SHF. Hoekstra examined the behaviour of 2 types of gulls in relation to orographic lift in a relatively flat area.

Appendix

Data management

All data for the model was retrieved from an external hard drive provided by Bouten. This data is copied and stored in combination with all other relevant material(Figure 15). All supplementary material referenced to in the text is found here as well. Copies will be stored on the hard drive the data was retrieved from as well as the following GitHub repository:

https://github.com/rensvanderveldt/rensvanderveldt_Bsc

Challenges of modelling soaring flight in humid landscapes A sensitivity analysis of vegetation parameters in the modelling of sensible heat flux