Urban energy exchanges monitoring from space

Nektarios Chrysoulakis, Sue Grimmond, Christian Feigenwinter, Fredrik Lindberg, Jean-Philippe Gastellu-Etchegorry, Mattia Marconcini, Zina Mitraka, Stavros Stagakis, Ben Crawford, Frans Olofson, et al.

<figure><img src="image_1.png" /></figure>

To cite this version:

Nektarios Chrysoulakis, Sue Grimmond, Christian Feigenwinter, Fredrik Lindberg, Jean-Philippe Gastellu-Etchegorry, et al. Urban energy exchanges monitoring from space. Scientific Reports, 2018, 8, 10.1038/s41598-018-29873-x. insu-03668273

HAL Id: insu-03668273

https://insu.hal.science/insu-03668273

Submitted on 16 May 2022

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>
SCIENTIFIC REPORTS

<figure><img src="image_3.png" /></figure>

OPEN

Urban energy exchanges
monitoring from space

Nektarios Chrysoulakis<sup>1</sup>, Sue Grimmond<sup>2</sup>, Christian Feigenwinter<sup>3</sup>, Fredrik Lindberg<sup>4</sup>, Jean-Philippe Gastellu-Etchegorry<sup>5</sup>, Mattia Marconcini<sup>6</sup>, Zina Mitraka<sup>1</sup>, Stavros Stagakis<sup>1</sup>, Ben Crawford<sup>2</sup>, Frans Olofson<sup>4</sup>, Lucas Landier<sup>5</sup>, William Morrison<sup>1,2</sup> & Eberhard Parlow<sup>1,3</sup>

One important challenge facing the urbanization and global environmental change community is to understand the relation between urban form, energy use and carbon emissions. Missing from the current literature are scientific assessments that evaluate the impacts of different urban spatial units on energy fluxes; yet, this type of analysis is needed by urban planners, who recognize that local scale zoning affects energy consumption and local climate. Satellite-based estimation of urban energy fluxes at neighbourhood scale is still a challenge. Here we show the potential of the current satellite missions to retrieve urban energy budget fluxes, supported by meteorological observations and evaluated by direct flux measurements. We found an agreement within 5% between satellite and in-situ derived net all-wave radiation; and identified that wall facet fraction and urban materials type are the most important parameters for estimating heat storage of the urban canopy. The satellite approaches were found to underestimate measured turbulent heat fluxes, with sensible heat flux being most sensitive to surface temperature variation (–64.1, +69.3 W m<sup>-2</sup> for ±2 K perturbation). They also underestimate anthropogenic heat fluxes. However, reasonable spatial patterns are obtained for the latter allowing hot-spots to be identified, therefore supporting both urban planning and urban climate modelling.

Received: 4 April 2018

Accepted: 17 July 2018

Published online: 31 July 2018

The Urban Energy Balance (UEB) needs to account for the 3D nature of cities, quantifying the fluxes into, out of and the storage change within the control volume (Fig. 1). In the last 15 years, significant advances in understanding urban processes have benefited from enhanced computational capacity, improved resolution of satellite sensors and increased ability to couple advance urban surface parameterization schemes with atmospheric models<fnref n="1" />⁻⁹. Earth Observation (EO) has been widely used to study the Urban Heat Island (UHI) phenomenon, but to a lesser extent to quantify heat fluxes¹⁰⁻¹⁵. Recent studies¹⁶,¹⁷ have investigated the potential of EO to derive turbulent heat fluxes and identify and analyse the associated uncertainties.

Both Earth system science and urban planning communities need spatially disaggregated UEB data at local scale¹⁸⁻²⁰ (neighbourhood, e.g., order (100 m × 100 m) or larger). However, such information is practically impossible to derive for extensive areas by *in-situ* flux measurements and EO-based estimation of spatio-temporal patterns of different UEB components is challenging. Therefore, the question arises whether EO can provide reliable estimates of UEB at the times of satellite acquisitions. Here we address this by investigating the potential of the current satellite missions to retrieve UEB fluxes at local scale, supported by meteorological observations. Recognising the range of city forms, size and settings, we explored three locations: London (UK), Basel (Switzerland) and Heraklion (Greece).

Results and Discussion

EO data from various sources were used to extract urban surface morphology and characteristics. Surface cover and material type were derived using advanced machine learning and fusion techniques and detailed spectral un-mixing approaches²¹ (e.g., Fig. 2a). Using EO-derived high resolution Digital Surface Models (DSM), surface roughness parameters (such as plan area index, frontal area index, roughness length and zero-displacement height) were calculated by morphometric analysis²². Examples of the morphometric analysis results are given in Supplementary Fig. S1. Satellite-derived thermal radiance at 1 km × 1 km was downscaled²³ to 100 m × 100 m. Although uncertainties of the downscaling method exist, as Mitraka et al.²⁴ noted, the error of downscaled

<footnote n="1">¹Foundation for Research and Technology Hellas, Institute of Applied and Computational Mathematics, Remote Sensing Lab, N. Plastira 100, Vassilika Vouton, 70013, Heraklion, Greece. ²University of Reading, Reading, UK. ³University of Basel, Basel, Switzerland. ⁴University of Gothenburg, Gothenburg, Sweden. ⁵Centre d'Etude Spatiale de la Biosphère (CESBIO), Toulouse, France. ⁶German Aerospace Center (DLR), Oberpfaffenhofen, Germany. Correspondence and requests for materials should be addressed to N.C. (email: zedd2@iacm.forth.gr)</footnote>
<figure><img src="image_4.png" /><figcaption>**Figure 1.** Conceptual illustration of the fluxes in the energy balance of an urban building-vegetation-soil-air volume. The equivalent surface energy budget per unit surface area through the top of the volume is: $Q^* + Q_F = Q_H + Q_E + \Delta Q_S + \Delta Q_A$ (W m⁻²), where, $Q^*$ is the net all-wave radiation flux; $Q_F$ is the anthropogenic heat flux (resulting from vehicular emissions, heating and cooling of buildings, industrial processes and metabolic heat release by people and animals); $Q_H$ is the turbulent sensible heat flux; $Q_E$ is the turbulent latent heat flux; $\Delta Q_S$ is the net change in heat storage within the volume; and $\Delta Q_A$ is the net advected flux. Arrows are drawn in the direction to which the corresponding flux is considered positive. $\Delta Q_S$ and $\Delta Q_A$ (= $Q_{A,in}$ – $Q_{A,out}$) are positive if the internal energy of the volume increases.</figcaption></figure>

<figure><img src="image_5.png" /><figcaption>**Figure 2.** Land cover and heat storage for part of London, UK. (a) Land cover (seven classes: buildings, impervious surfaces (i.e., roads, parking lots, sidewalks, etc.), water, bare soil, low vegetation (<2 m) and high vegetation (≥2 m) split into evergreen and deciduous trees). (b) Net change in heat storage ($\Delta Q_s$) (100 m × 100 m) on 19 July 2016 at 11:00 UTC. Highest values are observed in central densely built areas (City of London, Canary Warf (centre right)). River Thames is masked (black). Maps created with QGIS software, version 2.18 (www.qgis.org).</figcaption></figure>
<figure><img src="image_6.png" /><figcaption>**Figure 3.** Spatial patterns of radiation and turbulent heat fluxes. (a–c) EO-derived fluxes for the central area of Basel (on 23 August 2016 at 09:50 UTC, grid reference UTM 32 N (EPSG:32632)); and (d–f) modelled versus observed fluxes at BAES (Basel Aeschenplatz) flux tower for all 2016 cloud-free satellite scenes. (a,d) net all wave radiation (Q*); (b,e) turbulent sensible heat flux (QH); (c,f) turbulent latent heat flux (QE). Comparison with observations are stratified by season (colours, DJF refers to DEC, JAN, FEB, etc.) and wind direction (symbols); 1:1 line (solid), linear regression best fit line (dashed) and extent of 80% of the data points (dotted lines) indicated. Linear regression statistics, a: slope, b: intercept (W m⁻²), R²: coefficient of determination, MAE: mean absolute error (W m⁻²), RMSE: root mean square error (W m⁻²). Maps created with QGIS software, version 2.18 (www.qgis.org). Water surfaces are masked.</figcaption></figure>

temperatures, evaluated using in-situ surface temperature observations from the micrometeorological flux towers, is within 2.4 K (Supplementary Fig. S2).

The study period (January 2016–May 2017) allowed the seasonal cycle to be examined (Supplementary Fig. S3, London). For all cities this was appropriate with net all-wave radiation (Q*), net change in heat storage (ΔQs) and turbulent sensible heat flux (QH) peaking in summer. QH and ΔQs peak in August, whereas turbulent latent heat flux (QE) is relatively minor throughout the year in city centres.

Deriving Q* maps from EO is not straightforward as urban reflectance and thermal emission present anisotropic behaviour, caused by illumination geometry, 3D urban geometry and the distribution of urban material optical and thermal properties². Here, the DART (Direct Anisotropic Radiative Transfer) model⁶ takes into account this anisotropy. It simulates reflected and emitted radiation using information derived from EO satellite images (see more details in Methods). For example, Q* in Basel (Fig. 3a) range between 496 and 633 W m⁻². The highest fluxes occur near the river and in a vegetated region in the upper right of the domain. The mean absolute error (MAE, observed - modelled (EO-derived)) in Q* for 2016 was 22.5 W m⁻² (Fig. 3d), with agreement within 5%. This corresponds to an average root mean square error (RMSE) of 13 W m⁻² and R² > 0.99 (Supplementary Fig. S4). Uncertainty arises from simplifying the complex urban surface structure and the surface elements (e.g., walls) unseen by satellites, except through multiple scattering, although they contribute to Q*.

ΔQs of an urban canopy is approximately 2–6 times larger than for non-urban canopies. Here, ΔQs is estimated using the element surface temperature method (ESTM)³, which reduces the 3D urban structure to one-dimensional (1D) elements (see details in Methods). The sensitivity analysis identified that the fraction of wall facet and materials are the two most important variables in the estimation of ΔQs. Daytime peak (11:00 local time) ΔQs is largest (up to 400 W m⁻²) in densely built areas where tall buildings dominate (e.g., City of London and Canary Wharf in London Fig. 2b). Areas with low building density have, as expected, smaller ΔQs values. Uncertainty is assessed using two alternative approaches implemented within the SUEWS⁷ modelling platform: OHM²⁵, and AnOHM²⁶. All approaches produced similar results, but different spatial patterns are identified, mostly depending on the forcing data. To evaluate ESTM performance, observations from less complex sites with uniform land cover are used, where ΔQs can be directly measured with ground heat flux plates. MAE between 1.1 and 21.7 W m⁻² is found.
<figure><img src="image_7.png" /><figcaption>**Figure 4.** Anthropogenic heat flux ($Q_F$) as a residual from regressing turbulent heat fluxes ($Q_H + Q_E$) with available energy (observed $Q^*$ - modelled $\Delta Q_S$) in central Basel for all 2016 cloud-free satellite scenes stratified by season (colours, as Fig. 3d–f) and wind direction (symbols) for: (a) BAES (Basel Aeschenplatz (394330, 5267379)) and (b) BKL1 (Basel Klingelbergstrasse (393221, 5268560)) flux towers, with 1:1 line (solid), linear regression best fit line (dashed) and extent of 80% of the data points (dotted lines). Linear regression statistics, a: slope, b: intercept (W m⁻²), R², MAE (W m⁻²), RMSE (W m⁻²). (c) Relative spatial patterns of the anthropogenic heat flux ($Q_F$) in central Basel for January and February 2016, scaled from 0 to 1 (unitless). The river is masked (black). Grid reference is UTM 32 N (EPSG:32632). Map created with QGIS software, version 2.18 (www.qgis.org).</figcaption></figure>

The turbulent heat fluxes are strongly modified by large roughness elements and a complex mix of sources/ sinks of heat and water in cities. Here, the Aerodynamic Resistance Method (ARM) is used to derive $Q_H^{16,17}$ and $Q_E^{16}$ at local scale. To implement this, EO-derived downscaled land surface temperature (LST) and surface rough- ness parameters are used combined with air temperature and humidity from meteorological stations (see details in Methods). For example, for 23 August 2016 in Basel (Fig. 3b), the highest $Q_H$ values (up to 270 W m⁻²) are found in the most densely built-up and industrial areas. $Q_H$ values have similar spatial patterns and magnitudes to $\Delta Q_S$, whereas, $Q_E$ (Fig. 3c) is much smaller in the built-up areas (less than 150 W m⁻²) and has largest values in vegetated areas. The eddy covariance flux footprint²⁷ is used to compare measured with modelled turbulent heat fluxes. For the Basel BAES (Aeschenplatz) flux tower, the $Q_H$ MAE for 2016 is 41.2 W m⁻² (Fig. 3e), with modelled $Q_H$ underestimated compared to measured $Q_H$. Similarly, the EO-derived $Q_E$ is underestimated compared to flux tower observations (Fig. 3f). This relates to the uncertainties associated with ARM¹⁶,¹⁷ and the required input variables. Supplementary Figs S5 and S6 compare the EO-based turbulent sensible and latent heat fluxes with the ground truth at the three sites. A sensitivity analysis¹⁷ indicates that, overall, the EO-derived $Q_H$ is most sensitive to LST variation (+69.3, -64.1 W m⁻² for ±2 K perturbation).

Direct measurements of the anthropogenic heat flux $Q_F$ are extremely difficult²⁸. Here, $Q_F$ is obtained from the residual of the UEB, and thus includes the cumulative uncertainties of all the terms. Energy balance closure typically compares available energy ($Q^*$ - $\Delta Q_S$) with turbulent heat fluxes ($Q_H + Q_E$). Basel's two flux towers, have total EC turbulent heat fluxes that tend to be slightly greater than the available energy during winter months (Fig. 4a,b), when $Q_F$ is expected to be highest due to building space heating. The relative spatial pattern of EO-derived $Q_F$ for central Basel (winter 2016, Fig. 4c) appears reasonable with clear correlation with building density. However, the absolute $Q_F$ values for individual pixels include negatives (physically unrealistic). This suggests that individual UEB components are incorrect most likely due to the underestimation of $Q_H$. Highest relative $Q_F$ is found south of the river in the most densely built-up areas (Fig. 4c). Independent $Q_F$ models based on inventories²⁹ also confirm that buildings are the dominant source of $Q_F$. This positive correlation between built-up areas indicates that EO can be used to identify hot-spots with relatively high $Q_F$.

## Conclusions

Based on the analysis of individual UEB components, we conclude that current satellite missions have the potential to provide information about spatial patterns of urban energy exchanges, if supported by suitable meteorological measurements. Synergistic analysis of specific satellite observations is able to identify and map spatial distributions of heat fluxes at local scale, several times per day (if clear weather). The EO-based methods developed are easily transferable to any city and have the potential to support sustainable planning strategies, since
knowledge of UEB patterns at neighbourhood level is needed in urban planning (e.g., to reduce or prevent QH and QF hot spots), health studies (e.g., to estimate impact on thermal comfort) and future proofing (e.g., to plan and implement interventions to reduce heat emissions). Finally, the frequency of heat waves is expected to increase³⁰ with UHI and other urban characteristics exacerbating the respective warming, resulting in increased energy demand for cooling systems in low and mid-latitude cities, which in turn adds to heat emissions and raises temperatures further³¹. Our satellite-based approach, strengthened by uncertainty minimization in future applications, is expected to advance the current knowledge of the impacts of heat fluxes on energy consumption in cities, leading to the development of tools and strategies to mitigate these impacts, improving thermal comfort and energy efficiency.

## Methods

**Case studies and datasets.** Three cities were selected: highly urbanized mega-city (London, UK); typical central European medium size city, that requires substantial winter heating (Basel, Switzerland); and a smaller, low latitude Mediterranean city that requires substantial summer cooling (Heraklion, Greece). Meteorological data from wireless automatic weather station networks (AWSN) and flux towers (net all-wave radiation and turbulent heat fluxes) in all cities were supplemented by Large Aperture Scintillometers²⁸ in London. High spatial resolution satellite imagery (between 10 m × 10 m and 30 m × 30 m) from Copernicus Sentinel 2 and Landsat 5/7/8, as well as very high resolution (VHR) data (<10 m × 10 m) from Copernicus Contributing Missions, such as SPOT (Satellite Pour l'Observation de la Terre) and WorldView-2, RapidEye, TerraSAR-X and TanDEM-X are used. The EO thermal infrared data are from MODIS (Moderate Resolution Imaging Spectroradiometer), Sentinel 3 (1 km × 1 km) and ASTER (Advanced Spaceborne Thermal Emission and Reflection Radiometer) observations (at 90 m × 90 m)²¹. Synergies between Sentinels 2 and 3 proved useful to identify fractional surface cover needed in local scale surface emissivity and LST²³ estimates. All UEB fluxes are estimated at 100 m × 100 m resolution.

**EO data analysis.** Surface characteristics include morphology, cover and biophysical parameters. Morphology requires the location and dimensions of urban objects. For this purpose, DSMs are derived from VHR optical stereo imagery (Heraklion) or airborne LiDAR observations (Basel²¹ and London³²). Surface structure parameters (i.e., building volume, sky view factor, plan area index, frontal area index, roughness length and zero-plane displacement height) are derived by morphometric analysis²², using the Urban Multi-scale Environmental Predictor (UMEP)³³. Surface cover (plus its spatial and temporal changes) are obtained from EO using advanced machine learning techniques and detailed spectral mixture models³⁴. Multi-temporal acquisitions are used to update the surface cover fraction maps for the classes of interest. By fusing VHR data with Sentinel-2 observations it is possible to categorize different roof materials³⁵ and build spectral libraries³⁶. Other required parameters (e.g., emissivity, albedo and vegetation indices) are derived as outlined in Marconcini et al.²¹. As daily high spatial resolution thermal imagery is unavailable, the thermal radiance retrieved at 1 km × 1 km resolution (MODIS, Sentinel 3) are downscaled to 100 m × 100 m LST and emissivity using Mitraka et al.'s²³ novel technique. Uncertainties in this approach are quantified by Mitraka et al.²⁴. The downscaled LST values are compared to in-situ surface temperature observations from the micrometeorological flux towers.

**Estimation of net all-wave radiation.** Q* is the difference of the incoming solar (shortwave) radiation (Kl) and atmosphere thermal (longwave) emission (Ll), minus the outgoing reflected solar radiation (Kt) and urban surface thermal emission (Lt). Satellites observe in few spectral bands from restricted viewing directions, whereas Q* is an integral over a wide spectral domain and the whole hemisphere. Here, the DART model⁹,³⁷ is used to estimate local scale Q*. DART simulates radiative transfer of the urban surface-atmosphere system. It operates on scenes simulated as 3D arrays of rectangular cells. Here a processing chain is developed that calibrates DART with satellite images: the iterative comparison of the model output with images allows the material optical properties (OP) per urban element at the spatial resolution of the satellite data to be obtained. The iterative procedure accounts for the multiple scattering mechanisms and has six steps: 1) Urban morphology and material type are given as input to DART. 2) DART simulates at high spatial resolution (e.g., 2.5 m) all spectral bands of the satellite sensor being used for comparison. The simulation accounts for the specific atmosphere and illumination conditions of each input satellite scene. In the first iteration, spatially constant OP are set per urban element type. 3) DART-simulated spectral image output are georeferenced and spatially resampled to the input satellite scene. 4) DART-simulated and satellite images are compared on a pixel basis to improve the OP map per urban element, using the area of the urban elements within each satellite pixel (derived from the urban morphology and material database). If the DART-simulated and satellite images do not match, the procedure re-iterates from step 2. 5) DART simulates the angular spectral radiative flux along all directions sampling the upper hemisphere. 6) Shortwave K* (Kl - Kt) and longwave L* (Ll - Lt) exitance maps are computed as a double integral over the spectrum and the upper hemisphere. Finally, Q* (K* + L*) is resampled to 100 m × 100 m spatial resolution and compared against observations from the micrometeorological towers.

**Estimation of the net change in heat storage.** ΔQs is the net flow of heat stored in the urban volume (i.e., air, trees, buildings, ground, etc.). In urban areas, the net heat stored in the canopy is a relatively large fraction of Q* and directly evaluating ΔQs in the urban canopy is difficult³⁸. To determine ΔQs we use the Element Surface Temperature Method (ESTM), which reduces the 3D urban volume to four elements, i.e., building roofs, walls, and internal mass and ground (road, vegetation, etc.)⁶:

$$ \Delta Q_s = \sum_i \frac{\Delta T_i}{\Delta t} \omega_{ci} \Delta x_i f_i, \qquad (1) $$
where $\Delta T_i/\Delta t$ is the rate of temperature change over the period for each element $i$, $\rho c_i$ is the volumetric heat capacity, $\Delta x_i$ is the element thickness and $f_i$ is the fraction of each surface type. Each element type has sublayers (e.g., a wall can be built up by brick, insulation and wood). Without internal element temperature data the average is determined:

$$ \rho C \frac{\partial T}{\partial t} = -\frac{\partial Q}{\partial x} = -\frac{\partial}{\partial x} \left( -k \frac{\partial T}{\partial x} \right), \qquad (2) $$

where $Q$ is the conductive heat flux through the surface and $k$ the thermal conductivity. For the internal surfaces (i.e., roof, exterior and interior walls, and floors) the surface temperature of element $i$ is determined by setting the conductive heat transfer out of (in to) the surface equal to the radiative and convective heat losses (gains), as described by Offerle et al.⁶.

The morphology of the urban surface was derived from high resolution DSM, including 3D information of vegetation³²,³³. As land cover thermal properties differ, detailed land cover information was required. The seven land cover classes (Supplementary Fig. S1) have three surface materials types for impervious covers and five building categories. These are derived from Urban Atlas³⁹ land use, street view photographs (ground level provided by Google Maps) and local knowledge. The meteorological forcing from the AWSN with the EO down-scaled LST are used.

**Estimation of the turbulent heat fluxes.** The methodology uses the ARM approach⁴⁰ to estimate $Q_H$:

$$ Q_H = \rho c_p \frac{LST - T_a}{r_a}, \qquad (3) $$

where $\rho$ is the density of air, $c_p$ the specific heat of air at constant pressure (1005 J kg⁻¹ K⁻¹), $T_a$ is the air temperature provided by the AWSN and $r_a$ is the aerodynamic resistance (s m⁻¹). Analogously, $Q_E$ is expressed as:

$$ Q_E = \frac{\rho c_p e_s^* - e_a}{\gamma r_a + r_s}, \qquad (4) $$

where $e_s^*$ is the saturation water vapour pressure (hPa) at surface air temperature, $e_a$ is the atmospheric water vapour pressure (hPa), $\gamma$ is the psychrometric constant (0.67 hPa K⁻¹) and $r_s$ is the stomatal resistance (s m⁻¹). Stomatal resistance is calculated following Kato et al.⁴¹ using a simplified equation from Nishida et al.⁴²:

$$ \frac{1}{r_s} = \frac{f_1(T_a)f_2(PAR)}{r_{sMIN}} + \frac{1}{r_{cuticle}}, \qquad (5) $$

where PAR is the photosynthetic active radiation, $r_{sMIN}$ is the minimum stomatal resistance and $r_{cuticle}$ is the canopy resistance related to the diffusion through the cuticle layer of leaves (10⁵ s m⁻¹). Functions $f_1$ and $f_2$ are from Nishida et al.⁴² and $r_{sMIN}$ is for each vegetation type⁴¹. $Q_E$ is calculated by land cover type and weighted by the previous $r_{sMIN}$ type in every pixel. The aerodynamic resistance $r_a$ is⁴⁰:

$$ r_a = \frac{1}{u_* \kappa} \left[ \ln\left(\frac{z_{ref} - z_d}{z_{0m}}\right) - \psi_h\left(\frac{z_{ref} - z_d}{L}\right) + \ln\left(\frac{z_{0m}}{z_{0h}}\right) \right], \qquad (6) $$

with $u_*$ the friction velocity (m s⁻¹):

$$ u_* = U\kappa \left[ \ln\left(\frac{z_{ref} - z_d}{z_{0m}}\right) - \psi_m\left(\frac{z_{ref} - z_d}{L}\right) - \psi_m\left(\frac{z_{0m}}{L}\right) \right]^{-1}, \qquad (7) $$

where $\kappa$ is von Kármán's constant (0.4), $z_{ref}$ is height of wind measurements, $z_d$ is the zero-plane displacement height, L is the Obukhov length, $z_{0m}$ ($z_{0h}$) is the roughness length and $\psi_m$ ($\psi_h$) is the stability function for momentum (heat). $z_{0h}$ is often reported within the dimensionless number $k\beta^{-1}$, defined as

$$ \kappa\beta^{-1} = \ln\left(\frac{z_{0m}}{z_{0h}}\right). \qquad (8) $$

EC and scintillometry²⁸ measurements can be used to determine $k\beta^{-1}$ for their flux footprint. Here, flux tower measurements are used as reference values for the magnitude of the fluxes of momentum, sensible and latent heat during the satellite overpass. $z_{0h}$ is calculated from the roughness Reynold's number ($Re_*$) as⁴³:

$$ z_{0h} = z_{0m}[7.4 \exp(-\alpha Re_*^{0.25})], \qquad (9) $$

where $\alpha = 1.29$ and $Re_* (=z_{0m}u_*/\nu)$ uses the kinematic molecular viscosity $\nu$ (1.461 × 10⁻⁵ m s⁻¹).

The displacement height ($z_d$) and $z_{0m}$ are calculated with morphometric parameterization of Kanda et al.⁴⁴, using urban characteristics (mean and maximum building height, standard deviation of building height, plan area index and frontal area index) derived from 1 m × 1 m DSM, using UMEP²²,³³.

To evaluate the EO-derived turbulent heat fluxes, EC observations are compared to modeled fluxes in the flux tower source area. Here, the Kormann and Meixner²⁷ footprint model is used.
**Estimation of the anthropogenic heat flux.** With the other UEB components known, Q<sub>F</sub> is estimated as a residual and hence contains the net errors of all the components. The spatial Energy Balance Closure (EBC) can be determined from regression between (Q<sub>H</sub> + Q<sub>E</sub>) and (Q<sup>*</sup> - ΔQ<sub>S</sub>) to give an estimate of Q<sub>F</sub> and uncertainty. This requires three assumptions: 1) All sources of anthropogenic heat are released into the environment. 2) Advection is negligible at the scale of interest or cancels between scales: at the microscale (e.g., Q<sub>H</sub> between shadowed and sunlit patches, Q<sub>E</sub> between wet and dry patches); at the local scale (e.g., between parks, water bodies, and built-up areas of different density); and at the meso-scale (e.g., between city and surrounding rural environment: urban breeze; coastal cities sea breezes; topographic induced anabatic/katabatic effects). Moisture advection enhances latent heat flux and can be of similar size to Q<sub>H</sub>, but opposite in sign, thereby essentially offsetting each other<sup>45</sup>. For this local scale study the advection error is within the Q<sub>F</sub> estimate from energy balance closure. 3) Any unmeasured terms are incorporated in the error of the Q<sub>F</sub> estimate. The resulting Q<sub>F</sub> is evaluated by comparison to Q<sub>F</sub> estimates from alternative approaches, based on inventories<sup>29,33,46</sup>.

**Data availability.** This research was conducted under the framework of the URBANFLUXES project that received funding from the European Union's Horizon 2020 Research and Innovation Programme and joined the H2020 Pilot on Open Research Data. Therefore, all datasets related to the present study are included in the URBANFLUXES Data Repository as Open Data and can be accessed at: http://urbanfluxes.eu/data/.

**References**

1. Arnfield, A. J. Two decades of urban climate research: a review of turbulence, exchanges of energy and water, and the urban heat island. *Int. J. Climatol.* **23**, 1–26 (2003).

2. Voogt, J. A. & Oke, T. R. Thermal remote sensing of urban climates. *Remote Sens. Environ.* **86**, 370–384 (2003).

3. Offerle, B., Grimmond, C. S. B. & Fortuniak, K. Heat Storage and Anthropogenic Heat Flux in Relation to the Energy Balance of a Central European City Centre. *Int. J. Climatol.* **25**, 1405–1419 (2005).

4. Weng, Q. Thermal infrared remote sensing for urban climate and environmental studies: Methods, applications, and trends. *ISPRS J. Photogramm. Remote Sens.* **64**, 335–344 (2009).

5. Grimmond, C. S. B., Roth, M. & Oke, T. R. et al. Climate and More Sustainable Cities: Climate Information for Improved Planning and Management of Cities (Producers/Capabilities Perspective). *Procedia Environ. Sci.* **1**, 247–274 (2010).

6. Gastellu-Etchegorry, J. P., Yin, T. & Lauret, N. et al. Discrete Anisotropic Radiative Transfer (DART 5) for Modeling Airborne and Satellite Spectroradiometer and LIDAR Acquisitions of Natural and Urban Landscapes. *Remote Sens.* **7**, 1667–1701 (2015).

7. Ward, H. C., Kotthaus, S., Järvi, L. & Grimmond, C. S. B. Surface Urban Energy and Water Balance (SUEWS): development and evaluation at two UK sites. *Urban Climate* **18**, 1–32 (2016).

8. Ando, T. & Ueyama, M. Surface energy exchange in a dense urban built-up area based on two-year eddy covariance measurements in Sakai, Japan. *Urban Climate* **19**, 155–169 (2017).

9. Oke, T. R., Mills, G., Christen, A. & Voogt, J. A. Urban Climates. Cambridge University Press (2017).

10. Chrysoulakis, N. Estimation of the all-wave Net Radiation Balance in Urban Environment with the combined use of Terra/ASTER multispectral imagery and in-situ spatial data. *J. Geophys. Res.* **108**(D18), 4582 (2003).

11. Kato, S. & Yamaguchi, Y. Analysis of urban heat-island effect using ASTER and ETM+ Data: Separation of anthropogenic heat discharge and natural heat from sensible heat flux. *Remote Sens. Environ.* **99**, 44–54 (2005).

12. Kato, S. & Yamaguchi, Y. Estimation of storage heat flux in an urban area using ASTER data. *Remote Sens. Environ.* **110**, 1–17 (2007).

13. Xu, W., Wooster, M. J. & Grimmond, C. S. B. Modelling of urban sensible heat flux at multiple spatial scales: A demonstration using airborne hyperspectral imagery of Shanghai and a temperature-emissivity separation approach. *Remote Sens. Environ.* **112**, 3493–3510 (2008).

14. Weng, Q., Hu, X., Quattrochi, D. & Liu, H. Assessing intra-urban surface energy fluxes using remotely sensed aster imagery and routine meteorological data: A case study in Indianapolis, USA. *IEEE J. Sel. Top. Appl. Earth Obs. Remote Sens.* **7**, 4046–4057 (2014).

15. Chen, S. & Hu, D. Parameterizing Anthropogenic Heat Flux with an Energy-Consumption Inventory and Multi-Source Remote Sensing Data. *Remote Sens.* **9**, 1165 (2017).

16. Feigenwinter, C., Parlow, E. & Vogt, R. et al. Spatial Distribution of Sensible and Latent Heat Flux in the City of Basel (Switzerland). *IEEE J. Sel. Top. Appl. Earth Obs. Remote Sens.* https://doi.org/10.1109/JSTARS.2018.2807815 (2018).

17. Crawford, B. et al. Variability of urban surface temperatures and implications for aerodynamic energy exchange in cities. *Q. J. R. Meteorol. Soc.* in press (2018).

18. Seto, K. C., Dhakal, S. & Bigio, A. et al. Human Settlements, Infrastructure and Spatial Planning, Client-Mate Change 2014: Mitigation of Climate Change. *Contribution of Working Group III to the Fifth Assessment Report of the Intergovernmental Panel on Climate Change* (2014).

19. Creutzig, F., Baiocchi, G. & Bierkandt, R. et al. Global typology of urban energy use and potentials for an urbanization mitigation wedge. *Proc. Natl. Acad. Sci. USA* **112**, 6283–6288 (2015).

20. Seto, K. C. & Christensen, P. Remote sensing science to inform urban climate change mitigation strategies. *Urban Climate* **3**, 1–6 (2013).

21. Marconcini, M., Heldens, W. & Del Frate, F. et al. EO-based Products in support of urban heat fluxes estimation. Proceeding of Joint Urban Remote Sensing Event JURSE 2017, Dubai, UAE, March 30 – April 1, https://doi.org/10.1109/JURSE.2017.7924592 (2017).

22. Kent, C. W. et al. Evaluation of urban local-scale aerodynamic parameters: Implications for the vertical profile of wind speed and for source areas. *Boundary-Layer Meteorol.* **164**, 183–213 (2017).

23. Mitraka, Z., Chrysoulakis, N. & Doxani, G. et al. Urban surface temperature time series estimation at the local scale by spatial-spectral unmixing of satellite observations. *Remote Sens.* **7**, 4139–4156 (2015).

24. Mitraka, Z., Doxani, G., Del Frate, F. & Chrysoulakis, N. Uncertainty Estimation of local-scale land surface temperature products over urban areas using Monte Carlo simulations. *IEEE Geosci. Remote Sens.* Lett. **13**, 917–921 (2016).

25. Grimmond, C. S. B. & Oke, T. R. Heat storage in urban areas: Local-scale observations and evaluation of a simple model. *J. Appl. Meteorol.* **38**, 922–940 (1999).

26. Sun, T., Wang, Z. H., Oechel, W. C. & Grimmond, S. The Analytical Objective Hysteresis Model (AnOHMv1.0): methodology to determine bulk storage heat flux coefficients. *Geosci. Model Dev.* **10**, 2875–2890 (2017).

27. Kormann, R. & Meixner, F. X. An analytical footprint model for non-neutral stratification. *Bound.-Layer Meteorol.* **99**, 207–224 (2001).

28. Crawford, B., Grimmond, C. S. B. & Ward, H. C. et al. Spatial and temporal patterns of surface-atmosphere energy exchange in a dense urban environment using scintillometry. *Q. J. R. Meteorol. Soc.* **143**, 817–833 (2017).

29. Gabey, A. M., Grimmond, C. S. B. & Capel-Timms, I. Anthropogenic Heat Flux: advisable spatial resolutions when input data are scarce. *Theor. Appl. Climatol.* https://doi.org/10.1007/s00704-018-2367-y (2018).

30. Russo, S., Sillmann, J. & Fischer, E. M. Top ten European heatwaves since 1950 and their occurrence in the coming decades. *Environ. Res. Lett.* **10**, 124003 (2015).

31. Founda, D. & Santamouris, M. Synergies between Urban Heat Island and Heat Waves in Athens (Greece), during an extremely hot summer (2012). *Sci. Rep.* **7**, 10973 (2017).
32. Lindberg, F. & Grimmond, C. S. B. Nature of vegetation and building morphology characteristics across a city: influence on shadow patterns and mean radiant temperatures in London. *Urban Ecosyst.* **14**, 617–634 (2011).

33. Lindberg, F., Grimmond, C. S. B. & Gabey, A. *et al*. Urban Multi-scale Environmental Predictor (UMEP): an integrated tool for city-based climate services. *Environ. Modell. Softw.* **9**, 70–87 (2018).

34. Mitraka, Z., Del Frate, F. & Carbone, F. Nonlinear Spectral Unmixing of Landsat Imagery for Urban Surface Cover Mapping. *IEEE J. Sel. Top. Appl. Earth Obs. Remote Sens.* **9**, 3340–3350 (2016b).

35. Wang, Q., Shi, W., Li, Z. & Atkinson, P. M. Fusion of Sentinel-2 images. *Remote Sens. Environ.* **187**, 241–252 (2016).

36. Kotthaus, S., Smith, T. E., Wooster, M. J. & Grimmond, C. S. B. Derivation of an urban materials spectral library through emittance and reflectance spectroscopy. *ISPRS J. Photogramm. Remore Sens.* **94**, 194–212 (2014).

37. Gastellu-Etchegorry, J.-P., Demarez, V., Pinel, V. & Zagolski, F. Modeling radiative transfer in heterogeneous 3-D vegetation canopies. *Remote Sens. Environ.* **58**, 131–156 (1996).

38. Meyn, S. K. & Oke, T. R. Heat fluxes through roofs and their relevance to estimates of urban heat storage. *Energy and Build.* **41**, 745–752 (2009).

39. European Environment Agency. Urban Atlas, https://www.eea.europa.eu/data-and-maps/data/urban-atlas (2012).

40. Voogt, J. A. & Grimmond, C. S. B. Modeling surface sensible heat flux using surface radiative temperatures in a simple urban area. *J. Appl. Meteor.* **39**, 1679–1699 (2000).

41. Kato, S., Yamaguchi, Y., Liu, C. C. & Sun, C. Y. Surface Heat Balance Analysis of Tainan City on March 6, 2001 Using ASTER and Formosat-2 Data. *Sensors* **8**, 6026–6044 (2008).

42. Nishida, K., Nemani, R. R., Running, S. W. & Glassy, J. M. An operational remote sensing algorithm of land surface evaporation. *J. Geophys. Res.* **108**, 4720 (2003).

43. Kanda, M. & Moriizumi, T. Momentum and heat transfer over urban-like surfaces. *Boundary- Layer Meteorol.* **131**, 385–401 (2009).

44. Kanda, M., Inagaki, A., Miyamoto, T., Gryschka, M. & Raasch, S. A new aerodynamic parametrization for real urban surfaces. *Boundary- Layer Meteorol.* **148**, 357–377 (2013).

45. Pigeon, G., Lemonsu, A., Grimmond, C. S. B., Durand, P. & Masson, V. Divergence of turbulent fluxes in the surface layer: case of a coastal city. *Bound.-Layer Meteorol.* **124**, 269–290 (2007).

46. Lindberg, F., Grimmond, C. S. B., Yogeswaran, N., Kotthaus, S. & Allen, L. Impact of city changes and weather on anthropogenic heat flux in Europe 1995-2015. *Urban Climate* **4**, 1–15 (2013).

## Acknowledgements

The research was conducted in the framework of the project URBANFLUXES (http://urbanfluxes.eu), that has received funding from the European Union's Horizon 2020 Research and Innovation Programme, under Grant Agreement No 637519. We would like to thank all URBANFLUXES team members for their contributions to the project and especially Judith Klostermann from ALTERRA; Fabio Del Frate and Danile Latini from GEOK; Andy Gabey, Ting Sun and Christoph Kent from University of Reading; Roland Vogt and Andreas Wicki from University of Basel; Thomas Esch and Wieke Heldens from DLR; Ahmad Al Bitar and Nicolas Lauret from CESBIO; and Giannis Lantzanakis, Nektarios Spyridakis, Giorgos Kogxylakis and Nikos Manioudakis from FORTH. We would also like to thank Michael Abrams from Jet Propulsion Laboratory, California Institute of Technology, for the provision of time series of ASTER images. Other support was received by NERC (TRUC, AirPro).

## Author Contributions

Individual contribution of each co-author to the reported research: N.C. led the Heraklion case study, contributed to satellite data analysis, EC observations analysis, interpretation of results and wrote the first version of the paper. S.G. was responsible for the anthropogenic heat flux estimation, contributed to interpretation of results and led the London case study. C.F. was responsible for the turbulent heat fluxes estimation, contributed to interpretation of results and led the Basel case study. F.L. was responsible for the net change in heat storage estimation and contributed to interpretation of results. J.P.G.E. was responsible for the net all-wave radiation estimation and contributed to interpretation of results. M.M. was responsible for satellite data analysis and contributed to interpretation of results. Z.M. developed the surface temperature downscaling method, contributed to the satellite data analysis, to the interpretation of results and to the preparation of the Figures. S.S. contributed to the EC observations analysis, to the satellite data analysis, to the interpretation of results and to the preparation of the Figures. B.C. performed the Scinitillometry and ARM analyses in London and error analyses and Q~F~ analyses. F.O. implemented the ESTM approach for net change in heat storage estimation. L.L. performed the DART simulations for net all-wave radiation estimation. W.M. performed the surface temperature observations in London and contributed to the interpretation of results. E.P. supported the turbulent heat fluxes estimation and contributed to the interpretation of results.

## Additional Information

Supplementary information accompanies this paper at https://doi.org/10.1038/s41598-018-29873-x.

**Competing Interests:** The authors declare no competing interests.

**Publisher's note:** Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

<figure><img src="image_8.png" /></figure>

**Open Access** This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons license, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons license and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this license, visit http://creativecommons.org/licenses/by/4.0/.

© The Author(s) 2018