Aerodynamic roughness variation with vegetation: analysis in a suburban neighbourhood and a city park

Local aerodynamic roughness parameters (zero-plane displacement, zd, and aerodynamic roughness length, z0) are determined for an urban park and a suburban neighbourhood with a new morphometric parameterisation that includes vegetation. Inter-seasonal analysis at the urban park demonstrates zd determined with two anemometric methods is responsive to vegetation state and is 1–4 m greater during leaf-on periods. The seasonal change and directional variability in the magnitude of zd is reproduced by the morphometric methods, which also indicate z0 can be more than halved during leaf-on periods. In the suburban neighbourhood during leaf-on, the anemometric and morphometric methods have similar directional variability for both zd and z0. Wind speeds at approximately 3 times the average roughness-element height are estimated most accurately when using a morphometric method which considers roughness-element height variability. Inclusion of vegetation in the morphometric parameterisation improves wind-speed estimation in all cases. Results indicate that the influence of both vegetation and roughness-element height variability are important for accurate determination of local aerodynamic parameters and the associated wind-speed estimates.


Introduction
The (dis)services of urban vegetation are both context and scale specific, therefore cannot be generalised (Salmond et al. 2016). However, as the socio-environmental and economic benefits of urban 'green spaces' are realised, they are increasingly becoming part of planning agendas to mitigate climate change, improve urban sustainability and improve human well-being (e.g. Gill et al. 2007;Landry and Chakraborty 2009;Roy et al. 2012;Andersson-Sköld et al. 2015;Kremer et al. 2015;Salmond et al. 2016;Ward and Grimmond 2017). Green spaces therefore will continue to be (a greater) part of the urban fabric. Despite this, when modelling the urban environment vegetation is often neglected to simplify the problem (e.g. references within Grimmond et al. 2010Grimmond et al. , 2011. It is imperative that the understanding of the physical implications of urban vegetation is improved across micro-, local-, and regional scales. This extends beyond urban parks and vegetation in street canyonsas the edges of cities are approached vegetation may become the most prominent roughness elements (e.g. Giometto et al. 2017, Kent et al. 2017a.
The presence of urban vegetation has implications for the storage and fluxes of scalar properties (e.g. heat, moisture and pollutants). For example, vegetation can reduce the mean and extreme ambient and indoor temperatures , Schubert et al. 2012, Mavrogianni et al. 2014, Heaviside et al. 2015, whilst also reducing night-time longwave cooling (Coutts et al. 2016). Its presence tends to increase humidity (through increasing evapotranspiration) and is also responsible for precipitation interception, a reduction of run-off and increased soil water storage/ permeability (Stovin et al. 2008, Day et al. 2010, Vico et al. 2014. Vegetation contributes to pollutant absorption and deposition (Tiwary et al. 2009, Tallis et al. 2011, Salmond et al. 2016. Vegetation influences the momentum flux by exerting drag on the mean wind flow (Finnigan 2000, Guan et al. 2003, Krayenhoff et al. 2015, Giometto et al. 2017. At critical aerodynamic porosities (P 3D ) this drag can be as significant as solid structures of the same shape (Hagen and Skidmore 1971, Mayhead 1973, Grant and Nickling 1998, Guan et al. 2000, Rudnicki et al. 2004, Vollsinger et al. 2005, Koizumi et al. 2010, Kent et al. 2017a. Vegetation therefore influences the spatially-averaged mean and turbulent characteristics of the flow in urban areas (Krayenhoff et al. 2015), having implications for in-canopy flow (Salmond et al. 2013), as well as the exchange between in-and above-canopy air masses (Gromke andRuck 2009, Vos et al. 2013).
The influence of a defined surface area upon fluxes of momentum can be indicated using the aerodynamic parameters of the zero-plane displacement (z d ) and aerodynamic roughness length (z 0 ), which are directly related to surface characteristics. Several methods exist to determine these, including algorithms based upon surface form (morphometric methods) or observations (anemometric methods). The presence of all roughness elements is inherently included in anemometric methods, but until recently morphometric methods did not consider both vegetation and buildings in combination. However, Kent et al. (2017a) develop the widely-used Macdonald et al. (1998) (hereafter Mac) morphometric method to include vegetation, which also applies to the Kanda et al. (2013) (Kan) extension of the Mac method.
The objectives of this paper are to use observations at two vegetated urban sites to investigate: (i) the seasonal variability in z d and z 0 with the seasonal change of tree phenology, (ii) Kent et al.'s (2017a) parameterisation of vegetation in the morphometric methods and (iii) the implications of considering vegetation for accurate wind-speed estimation. The interdependence of z d and z 0 means that a single value for each parameter cannot be treated as the 'truth'. Therefore, the analysis provides a comparison between the magnitude and directional variability of roughness parameters determined from the different methods. The wind-speed estimation application provides an independent assessment of the method performance.

Methodology Site description and observations
Measurements from an urban park in Seoul, South Korea (Seoul Forest Park, SFP) and a suburban residential neighbourhood in Swindon, UK (SWD) are used. The obvious contrast of landscape with vegetation phenology means trees and other vegetation are expected to influence the aerodynamic properties of both areas, especially during leaf-on conditions when foliage is at relative maxima. Seoul Forest Park is the third largest park in Seoul (~116 ha), with a dominance of vegetation evident (Fig. 1a-d). The SWD site is typical of UK suburbia, with a slightly larger proportion of buildings than vegetation, but this varies with direction (e.g. Fig. 1e, f). Considerable research at the SWD site means anthropogenic and biogenic controls of energy, water and carbon fluxes and their temporal variability are well understood (Ward et al. 2013(Ward et al. , 2014(Ward et al. , 2015a(Ward et al. , 2015b(Ward et al. , 2015c. In addition, the site has been used during development of the Surface Urban Energy and Water Balance Scheme (SUEWS) (Ward et al. 2016). However, in-depth aerodynamic parameter analysis has not been performed at either the SFP or SWD site.
At each site, fast-response observations of temperature, wind velocity (uhorizontal, vtransverse and wvertical components), CO 2 and H 2 O are processed into 30-min averages (Table 1).

Surface elevation database and differentiation between buildings and vegetation
At both sites, 1-m horizontal resolution digital surface (DSM, ground height + surface features) and digital terrain (DTM, ground height only) models are analysed ( Table 2). The high resolution and accuracy of these data, allow intricacies of surface roughness (e.g. roof pitch) to be resolved. After subtraction of the DTM from the DSM to provide a roughness element surface model (RESM), pixels <2 m high are removed (i.e. street furniture and temporary obstacles, such as vehicles). This retains roughness elements which are most appropriate for application of the morphometric methods. Building and vegetation pixels are differentiated by three techniques.
For the SFP site, initial source area calculations (using the Kormann and Meixner (2001) and Kljun et al. (2015) models) indicate the measurements are consistently influenced by an area within 300 m of the sensor. The area within this radius is classified using a manual and automated technique. The manual technique entails classification of aerial photography (Fig. 2a) into: building, road, impervious, water, forest, grass, bare soil and other (unclassified, but with few roughness elements) (Fig. 2b), with the RESM data overlain to check for inconsistencies. This manual method has some limitations, for example, although buildings (predominantly rectangular with sharp boarders) are mostly captured, those within a waterworks (south of the SFP site) and in a ready mixed concrete (RMC) factory (north-west) are misclassified (Fig. 2b). Additionally, considerable vegetation is missed, especially at land cover interfaces (e.g. along roadsides and bare soil paths, Fig. 2b and c, magenta circles). After re-classification, a surface model of building (BSM) and vegetation canopy (CDSM) heights is created (Fig. 2c).
The automated separation of buildings and vegetation, uses the RGB colour band of aerial imagery, as vegetation tends to be darker (i.e. lower end of the saturation spectrum) for all colour bands. If higher saturation pixels are removed, a binary mask representing pixels which are likely vegetation can be retained (e.g. Fig. 2d) (Crawford et al. 2016). Clouds in the imagery  Raw 10 Hz data processed to 30-min averages with spike detection (Papale et al. 2006;Hong et al. 2009), night-time correction  and double rotation of the wind components, aligning the wind field to the u direction (McMillen 1988, Kaimal andFinnigan 1994 Raw 20 Hz data processed to 30-min averages using EddyPro Advanced (v5-00, LI-COR), which includes de-spiking, double coordinate rotation, humidity correction of sonic temperature and high-and low-frequency spectral corrections (Moncrieff et al. 1997).
* (Stewart and Oke 2012) causes some vegetation to be uncaptured by the mask (cf. Fig. 2a and e, magenta circle). A dark to lighter pixel transition on the edge of vegetation means the mask may be smaller than vegetation's true extent. Therefore, the binary mask and RESM are combined and a filtering algorithm flags pixels as vegetation if they are within ±3 m of another pixel in the binary mask. Pixels not flagged are either buildings or other urban furniture (e.g. cars, street lamps etc.). After removing pixel heights < 2 m, a final CDSM and BSM product is generated (Fig. 2e).
Although the manual (after re-classification) and automated CDSM and BSM products are almost identical (cf. Fig. 2c and e), the latter method is more practical. The remainder of this work uses a combined dataset from both procedures.
At the SWD site, the abundance of vegetation and proximity of built structures makes accurate manual classification difficult. Additionally, the automated technique frequently misclassifies building pixels as vegetation because of the dark roofs and excessive shading (e.g. Fig. 2f). Therefore, a building footprint dataset (OS MasterMap® Topography Layer -Building Height Attribute, Ordnance Survey 2014) (Fig. 2g) was overlain upon the RESM to create the BSM (Fig. 2h, red). The remaining pixels were classed as 'potential' vegetation pixels, with isolated pixels removed if fewer than 6 of the 8 surrounding pixels were not 'potential' vegetation (Goodwin et al. 2009, Lindberg and. The remaining pixels were stored as a CDSM (Fig. 2h, green).

Calculation of aerodynamic roughness parameters
Two anemometric methods are used to determine z d : the temperature variance (TVM, Rotach 1994) and wind variance (WVM, Toda and Sugita 2003) methods. The TVM and WVM are based upon the relation between the nondimensional standard deviation of temperature or vertical wind and stability parameter in the surface layer, during unstable conditions (Wyngaard et al. 1971, Tillman 1972: where σ T and σ w are the standard deviation of temperature and vertical wind velocity respectively, T * is the temperature scale, T * ¼ − w 0 T 0 =u * (with T the temperature, w the vertical wind velocity, u * friction velocity, the overbar representing a mean value and prime indicating deviation from the mean), L is the Obukhov length, L ¼ Tu * 2 κgT * (with g the gravitational acceleration and κ von Karman's constant = 0.4, Högström 1996) and C 1 -C 4 are constants.
The TVM and WVM are amongst the few methods that permit roughness parameters to be derived from single-level turbulence measurements. However, the methods rely on Monin-Obukhov similarity theory and that the resulting flux gradient relations used by the TVM and WVM (Eq. 1 and 2) apply in urban areas (see Roth and Oke 1995). Therefore, the applicability of the similarity relations used by the methods is assessed at both sites during this work. Although the similarity relations are expected to hold where flow is free from roughness-element wakes (i.e. within the inertial sublayer), the TVM is specifically developed to determine z d from measurement locations which may be distorted by local roughness-element wakes (i.e. within the roughness sublayer) (Rotach 1994). Previous analysis indicates results from the WVM are appropriate in similar heterogeneous locations (Toda andSugita 2003, Kent et al. 2017b).
The constants (C 1 to C 4 ) are derived from observations when z d is assumed negligible. Although the constants vary (e.g. Sorbjan 1989, Hsieh et al. 1996, Choi et al. 2004, the z d from the temperature and wind variance methods was found to be relatively insensitive to the range in a dense urban area (Kent et al. 2017b). To assess the effect of constant choice on the final solution to z d the methods are applied with a range of constants (Table 3). Note, if constants are fit to the observations at a site an a priori assumption of z d is required and therefore the z d retrieved is not useful (Kent et al. 2017b).
The right-hand sides of Eq. 1 and 2 are estimated by increasing z d from zero to twice the measurement height (z m ) in 0.1 m increments (producing ϕ est ). The z d is the value which minimises the root-mean-square error (RMSE) between ϕ est and the observed value (ϕ obs ) of σ T /T * or σ w /u * (for the TVM and WVM, respectively). As calculations are undertaken for unstable conditions (0.05 ≤ −z'/L ≤ 6.2, Roth 2000; z' = z mz d ) an initial z d for stability definition is required. Thus, the methods are applied to 10 o wind sectors around the sites with: (i) the z d for stability definition varied from 0 to 10 m in 2m increments (a larger initial z d provides insufficient data to apply the methods); and (ii) different constants (i.e. Table 3).
If measurements are free from roughness-element wakes (i.e. within the inertial sublayer), the 'eddy-covariance (EC) method' can be used to determine z 0 , which is a rearrangement of the logarithmic wind law: where the average wind speed (U z ) and u * are determined from observations at z m . For each 30-min period of observations, z d from both the temperature and wind variance methods are used, providing two z 0 solutions. The EC method, applicable under neutral conditions (|z'/L| ≤ 0.05), requires at least 20 observations to determine z 0 for a directional sector (Beljaars 1987, Grimmond et al. 1998. Additionally, only U z > 1 m s −1 are analysed to ensure sufficient mechanical turbulence (Liu et al. 2009). Stability corrections may be used to apply the EC method outside of neutral conditions. However, these corrections are based upon empirical fits to observed data and vary across studies (Högström 1996). To avoid additional sources of uncertainty only neutral conditions are considered here.
As the SFP site results indicate z d is similar to (or greater than) z m , the EC method to determine z 0 is therefore unusable (and not applied). For both northern-hemisphere sites, leaf-off periods are selected as the (core) winter months of December, January and February; and leaf-on periods are June, July and August. With little solar radiation during winter (leaf-off periods) at the SWD site there are insufficient unstable periods  Table 2; building footprints -Ordnance Survey 2014 to determine z d using the temperature and wind variance methods (and hence z 0 ). Therefore, only leaf-on conditions are analysed at the SWD site.
The Macdonald et al. (1998, Mac) and Kanda et al. (2013, Kan) morphometric methods are used with a new vegetation parameterisation (Kent et al. 2017a). Following the Kent et al. (2017b) methodology, an iterative procedure is applied using the Kormann and Meixner (2001) footprint model with 30-min averaged meteorological observations. Initial rural z d and z 0 values (0.2 and 0.03 m, respectively) are used, as results are independent of these values when applying an iterative procedure (Kent et al. 2017b). Morphometric calculations are only applied to source areas which extend horizontally beyond 50 m from the measurement sensors, as smaller source areas become concentrated upon only a few roughness elements and the morphometric calculations are inappropriate.
For each 30-min observation, the source area weighted geometry is calculated for buildings and vegetation (using the BSM and CDSM). The average, maximum and standard deviation of all roughness-element heights (H av , H max and σ H , respectively) are determined. The plan area index (λ p ) of roughness elements is: where W p,b and W p,v are the sums of weighted pixels in the source area of buildings and vegetation, respectively, W AT is the total sum of weights and P 3D is the aerodynamic porosity of vegetation. The weighted frontal area of buildings and vegetation is determined separately (W f,b and W f,v ), treating vegetation as non-porous. Including vegetation, the Mac method becomes (Kent et al. 2017a): where C Db = 1.2 is the drag coefficient for buildings and α = 4.43 and β = 1.0 are empirical constants for staggered arrays fit to the wind tunnel data of Hall et al. (1996). P v is the ratio between the drag coefficient for vegetation with varying P 3D and buildings (Kent et al. 2017a): derived from experiments with 0 ≤ P 3D ≤ 0.85 (Guan et al. 2000). The Kan method is a development of the Mac method, incorporating roughness-element height variability (Kanda et al. 2013): and where 0 ≤ X ≤ 1, 0 ≤ Y and a 0 , b 0 , c 0 , a 1 , b 1 and c 1 , are regressed constants of 1.29, 0.36, −0.17, 0.71, 20.21 and −0.77. The methods are applicable to any combination of buildings and vegetation, with vegetation phenology and associated drag characteristics being optimisable (through P 3D ). With this information being scarce, and the predominance of deciduous vegetation at both sites it is assumed that all vegetation has a leaf-on porosity of 20% and leaf-off porosity of 60% (i.e. P 3D = 0.2 and P 3D = 0.6, respectively, Heisler 1984; Heisler andDeWalle 1988, Grimmond andOke 1999). During leaf-on and leaf-off transition an intermediate porosity may be used Table 3 Constants (C 1 -C 4 ) for application of the temperature variance (TVM) and wind variance (WVM) anemometric methods (Eq. 1 and 2). For all observations, extensive flat homogeneous terrain is reported. Kaimal and Finnigan (1994) and Toda and Sugita (2003) are after synthesis of coefficients from various studies. At the SWD site, the Choi et al. (2004) constants are not applied, as they predict the scaled σ T and σ w to be much larger and smaller than observations, respectively, meaning z d solutions are consistently zero Tillman ( However, the rapid transition at both sites (< 30 days) means there is insufficient data to investigate the transition periods here. Determination of source-area weighted aerodynamic parameters using the morphometric methods (including vegetation) are implemented into the Urban Multi-scale Environmental Predictor (UMEP, http://www.urban-climate. net/umep/UMEP) climate service plugin for the open source software QGIS (Lindberg et al. 2018).

Impact of roughness elements on observational data
To assess the disturbance to measurements from nearby roughness elements the turbulence data are inspected (Fig. 3). At the SFP site, the data are more variable due to the proximity to roughness elements (measurements are at 1.6H av of all roughness elements in the 300-m radius) ( Fig.  3a-d). In some directions z m is similar to H av (N, SW, W, NW, Table 4b), and H max is always larger than z m . Therefore, the measurements are probably within the roughness sublayer (RSL) and z d is often larger than z m . A peak in the aerodynamic drag coefficient and transverse turbulence intensity between 130 o -180 o is likely caused by the rear sides of the sensor (Fig. 3a, c). In addition, there is a larger proportion of drag between 210 o -330 o where taller roughness elements are located (Fig. 3a).
Although the TVM and WVM methods have been demonstrated to be appropriate in the RSL, the calculation of source areas is more uncertain (e.g. Baldocchi 1997, Sogachev and Lloyd 2004, Vesala et al. 2008, Leclerc and Foken 2014. However, Fig. 3a-d demonstrates there is still some homogeneity to the flow and to characterise the local roughness, the Kormann and Meixner (2001) footprint model is applied at a height likely greater than the RSL (2.5H av = 20 m).
At the SWD site, the measurement height is approximately 2.8H av based on the measurement source area. Combined with a lack of disturbance to turbulence data for most directions (Fig. 3e-h), this indicates measurements are taken within the inertial sublayer (ISL), where it is most appropriate to apply the anemometric methods and source area calculations. Peaks in turbulence data between 100 o and 140 o and at approximately 180 o and 280 o (Fig. 3e-h) are likely caused by houses with maximum heights of up to 8 m -10 m in these directions (within 25 m of the sensor).

Land cover and geometry surrounding the sites
Source areas indicate the likely surface influencing turbulent fluxes measured at a point (Schmid and Oke 1990, Schmid 1997, Leclerc and Foken 2014 and can therefore be used to characterise site surroundings, with varying certainty (Heidbach et al. 2017). Source area characterisation is performed for the SWD site. However, at the SFP site, the    (Table 4b). An exception is to the north and west where the built fraction increases to close to 20% due to the RMC factory (Fig. 2a). Source area size varies with the model and parameters used (e.g. Leclerc andFoken 2014, Heidbach et al. 2017), as demonstrated with aerodynamic parameters (z d and z 0 ) from the Kan and Mac morphometric methods (Fig. 4). Independent of morphometric method, the SWD site source areas include residential housing, back gardens and impervious driveways. The source area climatology is biased towards the predominant south-westerly wind direction, where it also extends further upwind because of the greater wind speeds in this direction.
The Mac method source areas are larger than the Kan method (average upwind extents of 700 m and 400 m, respectively) due to the difference in aerodynamic parameters determined with each method. The peak flux footprint is 50 m upwind for the Mac method and 25 m upwind for Kan (Fig. 4). Although these differences impact the surface geometry and land cover determined within the source area, the consistent fetch at the SWD site means the parameters determined by the Mac and Kan source areas are remarkably similar (Table 5). The source area has 11% built and 2% vegetated roughness elements, with low-level vegetation (i.e. small shrubs), grass, impervious materials (e.g. roads or driveways) and soil forming the remainder. The latter have a comparatively small impact on the aerodynamic roughness parameters at the neighbourhood scale. The average height of vegetation is smaller than buildings (~3.5 m and 4.7 m, respectively), but the maximum tree height (up to 15 m) is slightly larger than buildings (up to 12 m). To the north-east (030 o -090 o ), trees are tallest and most abundant, whilst vegetation is least to the south-west (Fig. 2f, 210 -240 o ).

Seoul Forest Park (SFP)
At the SFP site, the anemometrically determined z d is relatively insensitive to both the 'universal' constants (Table 3) and initial z d used to define stability. Both consistently cause a maximum variability of < 1 m for any 10 o wind sector ( Fig. 5a and b, shading around grey and brown lines) which corresponds to < 10% of the median z d . Despite the proximity of measurements to roughness elements, the maximum RMSE between ϕ est and ϕ obs for the wind variance method is 0.4, which is similar to observations which are higher above roughness elements in other urban areas (Kent et al. 2017a) and provides greater confidence when using the WVM. In contrast, the RMSE for the TVM is much larger (2.0), because of the proximity to roughness elements, the thermal inhomogeneity of the area (i.e., water, grasses, trees, bare-soil, and impervious surfaces, Table 4) and the dissimilarity of roughness length between momentum and heat (e.g. Owen and Thomson 1963, Zilitinkevich 1995, Voogt and Grimmond 2000, Hong et al. 2012. Both the TVM and WVM indicate z d may be larger than z m (12.2 m) for all wind directions (Fig. 5a, b). Both methods have a larger z d during leaf-on (Fig. 5a) than leaf-off (Fig. 5b), which is approximately 1 m larger for the TVM and 2-4 m larger for the WVM. A seasonal increase in z d is also observed by Giometto et al. (2017) for a suburban neighbourhood with the larger amount of leaf-on foliage exerting greater drag upon the flow, acting to raise the centroid of the drag profile (analogous to z d , Jackson 1981).
The anemometric z d is more variable with wind direction during leaf-on, which can be directly related to trees (and Table 5 Roughness-element characteristics in the source areas modelled during leaf-on conditions at SWD using the (a) Kanda et al. (2013, Kan) and (b) Macdonald et al. (1998, Mac) morphometric methods (Fig. 4) for all observations (n = 8787) and in the specified north-easterly (NE) (n = 911) and south-westerly (SW) (n = 2045) directions  (Table 4b). The 6-m variation in z d between 270 o -360 o during leaf-off (Fig. 5b), is attributed to the disturbance to airflow from the RMC factory approximately 250 m upwind (Fig. 2a), and a row of trees (> 20 m) just 60 m upwind. During leaf-on, this variability is not observed because the foliage on the trees dominates airflow disturbance, causing an obvious increase in z d (Fig. 5a). Similar to the anemometric methods, both morphometric methods indicate leaf-on z d is larger than leaf-off z d, (Fig. 5a,  b). The effect is least obvious between 000 o and 120 o due to the lake and open grassed area. However, between 120 o to 280 o the increasing height and proportion of vegetation increases both Mac z d and Kan z d , with leaf-on z d 1 -3 m larger than leaf-off (Table 4b, Fig. 5). Both morphometric methods indicate maximum z d and seasonal signal between 270 o and 330 o , a similar direction to the anemometric methods. Here, Mac z d becomes as large as 7.5 m and Kan z d ; with a more pronounced peak, reaches 14 m. The latter is associated with the increased vegetation plan area, H av , H max and σ H (Table 4b). As the source area rarely extends to the RMC factory (Fig. 2a), the morphometric z d is primarily a function of vegetation in these directions.
Leaf-off z 0 is typically > 0.5 m larger than leaf-on z 0 for both morphometric methods (Fig. 5c, d), as z 0 varies with roughness-element density. In canopies with both vegetated (Shaw and Pereira 1982, Wolfe and Nickling 1993, Raupach 1992, Nakai et al. 2008) and built (Macdonald et al. 1998, Cheng and Castro 2002, Jiang et al. 2008) roughness elements, z 0 has been demonstrated to increase with density until a peak λ f (or leaf area index), beyond which z 0 decreases again. Therefore, the seasonal change of z 0 is expected to be canopy dependent: z 0 will increase with density for sparsely packed canopies, but will decrease with density in dense canopies. The SFP site is an example of the latter, where an already densely packed canopy during leaf-off conditions becomes denser during leaf-on. This effectively closes the canopy creating a smoother surface with a flow more Anemometric methods: TVM (temperature variance, Rotach 1994); WVM (wind variance, Toda and Sugita 2003), with the range of solutions (shading) provided from varying constants used during application of the methods and initial z d used to define stability. Morphometric methods: Kan (Kanda et al. 2013); Mac (Macdonald et al. 1998) applied using the Kormann and Meixner (2001) footprint model for all 30-min observations, with 10th/ 90th percentile for each 10 o sector shaded. z d is for unstable conditions (6.2 ≤ z'/L ≤ 0.05, with z' = z m − z d and L the Obukhov length) and z 0 is for neutral conditions (|z'/L| ≤ 0.05). For morphometric method stability definition, z d in z'/L is determined by the respective morphometric method. Insufficient neutral conditions for southerly winds during leaf-off vegetation state means no morphometric z 0 is reported here characteristic of a skimming regime and reduction in z 0 . The comparatively sparsely packed neighbourhood site analysed by Giometto et al. (2017) is an example of the former, whereby leaf-on transition creates an effectively rougher surface with a flow more characteristic of a wake interference regime and resultant increase in z 0 .
At the SFP site, the seasonal change in Mac z 0 is more obvious than Kan z 0 , as the former is more sensitive to λ f (e.g. Kent et al. 2017b, their Fig. 1) and Kan z 0 considers geometric parameters other than λ f (H max , σ H and λ p ). Kan z 0 is on average 0.2 m less than Mac z 0 because the larger z d determined using the Kan method means physically less frontal area of roughness elements exert drag upon the mean wind flow. For both morphometric methods, leaf-off z 0 is consistent with direction due to the relative lack of foliage (Fig. 5d). However, leaf-on z 0 is much more directionally variable and similarly to z d can be directly related to vegetation geometry and cover. Mac z 0 has greater directional variability because of the aforementioned sensitivity to λ f .
The morphometric methods were applied treating vegetation as buildings (i.e. P 3D = 0) and ignoring it (i.e. P 3D = 1), however, the dominance of vegetation in the area meant the former produced z 0 < 0.25 m and the latter z d < 5 m. Furthermore, applying the methods without the vegetation parameterisation does not produce the seasonal change demonstrated by the observations.

Swindon (SWD)
During leaf-on, median solutions to both the temperature and wind variance methods indicate z d varies between 4 and 10 m surrounding the SWD site (Fig. 6a). However, the range of z d for any 10 o sector is up to 5 m for the temperature variance method and up to 2.5 m for the wind variance method, corresponding to as much as 50% of median z d . This range is larger than previous applications of the methods (e.g. Kent et al. 2017b) and at the SFP site. The z d from the temperature variance method cannot be related to surface characteristics, providing a z d which is consistently close to z m and up to 5 m larger than the wind variance method. The predominantly large z d solutions are likely because of the thermal inhomogeneity of the area, which includes buildings, vegetated, and paved land cover. Therefore, similar to the SFP site, there is considerable variability between ϕ est and ϕ obs for the temperature variance method (RMSE > 0.6) and less confidence in its use.
The wind variance method indicates z d is consistently between 4 and 5 m (i.e. similar to H av ) (Fig. 6a). The directional variability of these results can be directly related to surface Fig. 6 Comparison of anemometric (lines and shading) and morphometric (points) methods to determine the (a) zero-plane displacement (z d ) and (b) aerodynamic roughness length (z 0 ) (note log y axis) surrounding the Swindon site (Fig. 2f). For anemometric methods, z d is the median solution of the temperature variance (TVM) and wind variance (WVM) methods, applied to 30-min observations during unstable conditions (6.2 ≤ z'/L ≤ 0.05, with z' = z mz d and L the Obukhov length) for 10°sectors. The range (shading) represents all possible solutions by varying z d used for stability definition and varying constants used in the methods. z 0 is the median (lines) and upper and lower quartile (shaded) of the eddy-covariance method, during neutral conditions (|z'/L| ≤ 0.05) for each 10 o sector, using z d from the TVM and WVM, respectively. Less than 10 observations in the 110 o sector means no values are reported here. The morphometric methods: Kan (Kanda et al. 2013); Mac (Macdonald et al. 1998) are for each source area during the same conditions as the anemometric methods, applied considering vegetation (subscript bv) and for buildings only (subscript b). For the morphometric method stability definition, z d in z'/L is determined by the respective morphometric method characteristics. Combined with the lower RMSE between ϕ est and ϕ obs (RMSE < 0.2), there is greater confidence in the wind variance than the temperature variance method results. Increases of z d and z 0 of up to 7.5 m and 1 m, respectively between 130 o -180 o and 240 o -280 o are associated with houses within 25 m of the sensor in these directions. Elsewhere, the larger z d from the temperature variance method means that its associated z 0 is consistently 0.25 m -0.5 m less than the wind variance method (Fig.  6b).
Relative minima of both z d and z 0 occur when wind flow is aligned with the smoother road surface to the west of sites between 200 o -210 o and 330 o -360 o (Fig. 6a). A relative increase in z 0 for both methods (to approximately 0.25 m and 0.75 m for the temperature and wind variance methods, respectively) in the 045 o direction is likely because of the taller and more abundant vegetation in the same direction. However, there is not a similar increase of z d .
For both z d and z 0 the morphometric methods have less directional variability than the anemometric methods because of the similarity in geometry surrounding the SWD site. When vegetation is considered Kan z d ranges between 2.5 m and 7 m and is therefore approximately 0.5-1.5H av (Fig. 6a, Kan bv ). In comparison Mac z d is consistently half of this, ranging between 0.25-0.75H av (Fig. 6a, Mac bv ). The z d determined by the morphometric methods is more similar to the wind variance method (than the temperature variance), especially for Kan z d , which has an average difference of 0.2 m. The methods indicate relative minima and maxima of z d in similar directions. A relative reduction between 030 o -090 o occurs where fewer buildings are located (Table 5) values are lower. Morphometrically determined z 0 are within the anemometric range, except when increased friction velocity from nearby roughness elements creates an anemometric z 0 that is double the morphometric results (Fig. 6b, 120 o -180 o ).
When the morphometric methods are applied without considering vegetation, the average reduction of z d is 0.3 m and z 0 is 0.1 m (Fig. 7b, Kan b and Mac b ). However, these values are directionally dependent. For example, not considering the taller and more abundant vegetation between 030 o -090 o means z d and z 0 are reduced with an average of up to 20% and 40%, respectively, for both morphometric methods. Giometto et al. (2017) also demonstrate overlooking vegetation leads to a reduction of up to 50% in both z d and z 0 for a neighbourhood site with a larger plan area of vegetation and taller trees than the SWD site. Both results highlight the importance of considering vegetation during aerodynamic parameter determination.

Similarities in aerodynamic parameter analysis between the sites
The variability in the anemometric methods and the interdependence of z d and z 0 (i.e. the former is used when determining the latter, Eq. 3) means it is difficult to use the anemometric methods as a basis for the most appropriate magnitude of z d and z 0 . However, there are apparent similarities from the aerodynamic parameter analysis performed at both sites. There is greater uncertainty in the application of the temperature variance method to determine z d , than the wind variance method. The RMSE between ϕ est and ϕ obs of the former are consistently twice the latter, which is attributed to the thermal inhomogeneity of both sites and dissimilarity of roughness length between momentum and heat. Greater uncertainty was also found in the temperature variance method in a central urban area (Kent et al. 2017b).
As with previous applications, z d determined using the temperature and wind variance methods at both sites indicates z d is larger than H av (e.g. Grimmond et al. 1998, Feigenwinter et al. 1999, Kanda et al. 2002, Tsuang et al. 2003, Christen 2005, Chang and Huynh 2007, Tanaka et al. 2011, Kent et al. 2017b. Additionally, morphometric z d results are consistently smaller than anemometric results. However, the direct incorporation of height variability in Kan z d means it is more similar to the anemometric methods than Mac z d . Mac z d is less than H av and may be appropriate for homogeneous groups of roughness elements. However, recent literature demonstrates that the disproportionate amount of drag imposed by taller roughness elements in a heterogeneous mix means z d may indeed become larger than H av (Jiang et al. 2008, Xie et al. 2008, Hagishima et al. 2009, Zaki et al. 2011Millward-Hopkins et al. 2011, Tanaka et al. 2011, Kanda et al. 2013. This is particularly evident during leaf-on conditions at the SFP site (Fig. 5a), where z d appraoches 2H av .
The contrast in magnitude of aerodynamic parameters determined using the Kan and Mac morphometric methods (Kan z d is typically twice Mac z d ) relates to the Kan method's direct consideration of roughness-element height variability. However, as the Kan method, is developed from Mac, their directional variability is similar. At both sites, the morphometric methods show similar directional variability to the anemometric methods, indicating sound performance of Kent et al.'s (2017a) vegetation parameterisation.
For both morphometric methods, the range of z d and z 0 for any direction is consistently within ±1 m and 0.2 m of the median, respectively (Fig. 5 shading and Fig. 6 range of points for a direction). The range is attributed to the source area variability with meteorological conditions. For example, a wider range between 000 o -120 o at the SFP site is caused by the proportion of the source area which falls upon the lake, grassed area and surrounding trees (Fig. 5).

Implications for wind-speed estimation
With pre-determined z d and z 0 , the logarithmic wind law can be used to model the neutral vertical profile of wind speed (Tennekes 1973): which theoretically only applies in the inertial sublayer (ISL), where vertical fluxes of momentum can be assumed constant with height (e.g. Tennekes 1973). Closer to a rough surface (i.e. within the RSL) the roughness-element wakes create a highly variable flow which may deviate considerably from Eq. 10 (e.g. Thom et al. 1975, Simpson et al. 1998, Kastner-Klein and Rotach 2004, Christen 2005, Harman and Finnigan 2007, Barlow and Coceal 2009, Giometto et al. 2016. With measurements at the SFP site (1.6H av ) closer to the roughness elements, there is greater confidence to use Eq. 10 at the SWD site (where z m = 2.8H av ).
To assess the vegetation parameterisation within the Kan and Mac morphometric methods, the wind speeds measured at the SWD site during neutral conditions (|z'/L| < 0.05) are estimated using Eq. 10 with the observed u * for each 30-min period and roughness parameters determined with, and without vegetation (Fig. 6, subscript bv with and b without vegetation). The estimated wind speed (U est ) is regressed against the mean observed wind speed (U obs ) for the corresponding time period (Fig. 7). As the RMSE has been demonstrated to disproportionately amplify the error associated with outliers when assessing model performance (Willmott and Matsuura 2005), both the RMSE and mean absolute error (MAE) between U est and U obs are reported.
Wind speeds are overestimated in over 90% of cases, which is more apparent at higher U o b s (Fig. 7). Overestimation could be for several reasons, including uncertainty of the use of the logarithmic wind law closer to roughness elements or the appropriateness of z d and z 0 values obtained from the different methods (e.g. Millward-Hopkins et al. 2012). However, irrespective of the morphometric method U est most resembles U obs when aerodynamic parameters determined considering vegetation and buildings are used. For example, wind speeds estimated using Kan bv and Mac bv have Fig. 7 Observed (U obs ) and estimated (U est , Eq. 10) wind speeds for each 30-min period (at the SWD site) using the Kanda et al. (2013, Kan) and Macdonald et al. (1998, Mac) morphometric methods, considering both buildings and vegetation (subscript bv) and only buildings (subscript b). Data are binned from lowest wind speed in groups of 1250 (30-min) data points. Median (points) and 5th and 95th percentiles shown. The rootmean-square error then mean absolute error (m s −1 ) (between U est and U obs ) are given in the legend MAE from U obs of 0.92 and 1.44 m s −1 , respectively, whereas ignoring vegetation (i.e. Kan b and Mac b ) the MAE is > 0.3 m s −1 larger for both methods (1.31 and 1.73 m s −1 , respectively). The lower errors (both RMSE and MAE) associated with the Kan method indicate that regardless of whether vegetation is considered, incorporating height variability improves windspeed estimates.
Similar comparisons between U est and U obs are performed for wind directions with the least (210 o -240 o ) and greatest (030 o -090 o ) vegetated roughness elements ( Fig. 8a and b, respectively). As the least vegetated directions have similar aerodynamic parameters (Fig. 6) their associated U est are similar irrespective of whether vegetation is considered or not. However, despite the small number of trees, accounting for them still reduces the error in wind speed estimation (i.e. the lower errors of Kan bv and Mac bv ) (Fig. 8a). The importance of considering height variability is apparent again, as the Kan method reduces the errors from U obs by over 0.5 m s −1 , in comparison to the Mac method.
As expected, in directions with maximum vegetation (tree) cover (030 o -090 o ) the impact on estimated wind speeds is greatest. Inclusion of vegetation consistently results in an improvement of wind-speed estimation of over 0.5 m s −1 (Fig. 8b). The smallest differences between the errors associated with Kan bv and Mac bv occur in this direction (0.2 m s −1 ). Combined with the errors for Kan b being larger than Mac bv , the incorporation of vegetation appears more important for accurately estimating the wind speeds than considering height variability (in this case).

Conclusions
Two anemometric and two morphometric methods are used to determine the zero-plane displacement (z d ) and aerodynamic roughness length (z 0 ) for an urban park and a suburban neighbourhood. The anemometric methods use in-situ single-level high frequency observations and therefore inherently include the presence and state of vegetation. The morphometric methods have been developed for bluff bodies only, however a new parameterisation (Kent et al. 2017a) to consider both buildings and vegetation is explored.
At both sites, z d determined using the anemometric methods is larger than the morphometric methods. There is greater uncertainty in an anemometric method based upon scaled temperature variance, as opposed to the vertical wind velocity variance, likely because of the thermal inhomogeneity of the sites. However, the Kanda et al. (2013) morphometric method, which directly considers roughness-element height variability, is consistently most similar to observations, indicating z d is larger than average roughness-element height at the respective sites.
Inter-seasonal analysis is performed at the urban park, which is predominantly vegetation, with few buildings. Both anemometric methods indicate z d during leaf-on vegetation state is up to 1-4 m larger than leaf-off. In addition, leaf-on z d is obviously larger in directions with taller, or a greater proportion of, vegetated roughness elements. The morphometric methods with the vegetation parameterisation have a similar magnitude and directional variability of change, indicating leaf-on z d is 1-3 m larger than leaf-off, which varies with upwind roughness elements. When the anemometric z d is similar to, or larger than, the measurement height z 0 cannot be determined from observations. However, the morphometric methods indicate leaf-on z 0 may be less than half leaf-off z 0 because the additional tree foliage in an already densely packed area creates an effectively smoother canopy.
The suburban neighbourhood has a larger proportion of buildings than trees. Morphometric analyses are undertaken during leaf-on conditions with and without vegetation. Where there is confidence in the anemometric methods, their z d and z 0 can be directly related to surface characteristics surrounding the site. The morphometric methods have similar directional change to the anemometric methods, but with less variability as the geometry of the site surroundings are similar. If vegetation is ignored in the morphometric calculations, z d and z 0 decrease by to 20% and 40%, respectively.
Wind speeds estimated at the suburban site using the logarithmic wind law and aerodynamic parameters from the morphometric methods are compared to observed wind speed at approximately three times the average roughness-element height. Wind-speed estimations most resemble observations when vegetation (in addition to buildings), as well as the height variability of roughness elements are considered. The consideration of vegetation is more important than the roughness-element height variability in directions where vegetation cover is maximal. Kent et al.'s (2017a) extension of the morphometric methods captures the presence and state of vegetation for aerodynamic parameter determination and wind-speed estimation. As green spaces become increasingly part of the urban fabric, understanding the implications of vegetation upon aerodynamic characteristics becomes more important. Further observations with different types, amounts and arrangements of vegetation will allow more thorough assessment of this parameterisation.
The methodology to determine z d and z 0 from surface elevation databases (including vegetation) is freely available in the Urban Multi-scale Environmental Predictor (UMEP, http://www.urban-climate.net/umep/UMEP) for the open source software QGIS.