Summertime Assessment of an Urban-Scale Numerical Weather Prediction System for Toronto

Abstract

Urban-scale Numerical Weather Prediction (NWP) systems will be important tools for decision-making in and around large cities in a changing climate exposed to more extreme weather events. Such a state-of-the-art real-time system down to 250-m grid spacing was implemented in the context of the Toronto 2015 Panamerican games, Canada (PanAm). Combined with the Global Environmental Multiscale (GEM) model, attention was brought to the representation of the detailed urban landscape, and to the inclusion of sub-daily variation of the Great Lakes surface temperature. Results show a refined representation of the urban coastal environment micro-meteorology with a strong anisotropy of the urban heat island reaching about 2 °C on average for the summer season, coastal upwelling, and mesoscale features such as cumulus clouds and lake-breeze flow. Objective evaluation at the surface with a dense observational network reveals an overall good performance of the system and a clear improvement in comparison to reference forecasts at 2.5-km grid spacing in particular for standard deviation errors in urban areas up to 0.3 °C for temperature and dew point temperature, and up to 0.5 m s${}^{-1}$ for the wind speed, as well as for precipitation with an increased Equitable Threat Score (ETS) by up to 0.3 for the evening accumulation. The study provides confidence in the capacity of the new system to improve weather forecasts to be delivered to urban dwellers although further investigation of the initialization methods in urban areas is needed.

1. Introduction

Crowded cities need early warning systems to face large environmental challenges. Stakeholders will have to adapt to more frequent critical extreme weather events due to climate change [1]. Actions might need to be taken to prevent impacts of heat and cold spells, flooding events associated with heavy precipitation or spring melt, damaging winds, fog, large energy consumption, and unhealthy air quality. Development of urban integrated hydro-meteorological, climate and environmental services for citizens and stakeholders has been recently recognized as an urgent practice for urban metropolitan areas [2,3]. An important piece to achieve such tailored and high-quality services is to refine Numerical Weather Prediction (NWP, see a list of abbreviations at the end of the manuscript) systems to account for the intra urban-scale heterogeneity. Despite the progressive inclusion of urban physical processes in weather, climate and environmental operational activities throughout the world, sufficient resolution for such systems is most of the time missing albeit critical. Due to its predominant impervious nature and its spatial arrangement, the urban landscape interacts with the atmosphere in a relatively different way than other land surfaces. When implemented in high-resolution weather forecasting systems, these effects have been shown to improve local surface energy budgets and near-surface meteorology prediction [4,5,6]. Small-scale surface-induced atmospheric circulations in and around cities are therefore likely to be better represented as shown for example for sea-breeze circulations [7] or precipitation patterns [8]. Based on numerous studies, the general effect of urban areas on meteorological and climatological features can be assessed with the so-called Urban Heat Island (UHI) [1,6]. Thermo-dynamical contrasts between cities and their surroundings can result in meteorological extreme events. For example, enhanced vertical motions lead to enhanced convection with the capacity to produce localized and heavy precipitation [8] or to modify pollutant dispersion [9]. Meteorological impacts of urban areas depend on multiple parameters. For such considerations, operational high-resolution near-real-time NWP remains an efficient way of predicting the meteorology variability at the urban scale.

Recent major international sporting events have provided opportunities to develop high to very-high resolution NWP systems. Advances in wintertime forecasts for mountainous regions with Environment Canada and Climate Change (ECCC)’s NWP systems were provided at the kilometer scale [10] or even with hectometer resolution by use of a uncoupled surface prediction system [11] during the 2010 Vancouver Olympic Games Canada [12]. Grid spacing down to 250 m using the surface-atmosphere coupled system was then used subsequently in real time for the 2014 Sochi Olympic Games, Russia. In particular, visibility was better forecasted with ECCC’s high-resolution system [13]. During other summertime events, sub-kilometer grid spacing as low as 330 m was achieved in support of the 2012 summertime Olympics games in the United Kingdom (UK) for Weymouth Bay to forecast wind conditions for sailing competitions, whereas weather forecasts over London, UK, were provided with 1.5-km grid spacing [14].

To support the 2015 summertime Pan and Parapan American sport Games in Toronto (hereafter referred to as PanAm), a new subkilometer urbanized NWP system was developed and used in real time for a large metropolitan region as part of an ECCC science project. Sports venues and recreation areas during the Games were located in various urban locations. Sailing competitions occurred also in the vicinity of downtown Toronto. Although the primary goal was experimental rather than operational, the main objectives of this effort were to demonstrate the feasibility of implementing such a high-resolution system for real-time forecasting, to assess the quality of the forecasts in an urbanized region, and to initiate collaborations with end-users in order to extract adequate information from such a system.

An overview of the interdisciplinary aspects of this research project is provided in [15]. It involved the deployment of an observational network through the region, a suite of weather, air quality, and environmental forecasts, and products integrated as part of health-related initiatives. For the high-resolution urban forecasting component, the production chain can be seen as a suite of four major requirements: first, the model dependencies (large-scale forcing and advanced methods for surface assimilation), second, the model configuration, then, some post-processing of extreme weather indices or boundary-layer height estimation based on vertical profiles, and finally, sharing with ECCC’s nowcasting platform, and dissemination to the different users. Although the post-processing and distribution were important components to explore the ability of the system to provide valuable operational products and end-users services in the context of a high-attendance event [3], this present work focuses on the scientific modeling components and their dependencies.

Sub-kilometer urban-scale forecasts have long been a research topic. They were shown to capture more details of the sub-mesoscale meteorological features relative to a more detailed surface heterogeneity description, for the city itself, such as the formation of an “urban breeze” similar to a density current [16,17] and for the interaction with the surrounding water bodies [7]. Over London, UK, simulations with 333-m grid spacing could improve the fog forecasting [18]. Analysis of representation of convective rolls in the turbulence grey zone with 50- and 100-m grid spacing has led to the conclusion that, even if turbulent motions are not yet realistically resolved at these scales but rather remain oversized in comparison to reality, traditional forecasts would benefit from hectometre scale NWP [19]. Such systems are known to have currently large computational costs and to the authors’ knowledge there is still no such daily forecasts provided by national weather centres. Due to the growing computational capabilities of national prediction centers, studies with longer periods and demonstration projects are emerging. For example, [5] presented high-resolution forecasts over the Amsterdam metropolitan area during the same summer (2015) showing improvement of 100-m and 500-m grid spacing forecasts at the city scale.

This paper’s objectives are threefold. First, to highlight the advances in the methodology and configuration of this high-resolution urban deterministic system used in the context of the PanAm Games over the Greater Toronto Area (GTA), second, to demonstrate the forecast value for some specific meteorological patterns in and around GTA and finally, to objectively evaluate results during summer 2015. The paper is then organized as follows: Introduction, Materials and methods, Results and discussion, and Conclusions.

2. Materials and Methods

2.1. NWP Systems

Daily experimental forecasts were operated at ECCC in addition to operational forecasts provided by the Meteorological Service of Canada (MSC) [15]. Focusing on those used in this study, the High-resolution Deterministic Prediction System (HRDPS) operated at 2.5 km for short-range forecasting in Canada was in the process of becoming operational a few months later [20], and the high-resolution urban deterministic NWP system developed for the Greater Toronto Area (GTA), named TOR-GEM hereafter, was developed to provide a more realistic representation of the surface heterogeneity based on the representation of cities and on refined horizontal and vertical model resolution. Domains are shown in Figure 1 and a summary of physical parametrizations and model grids characteristics are reported in Table 1.

2.1.1. Overview

NWP tools developed at ECCC are based on the Global Environmental Multiscale GEM model [21]. The present systems were built upon one-way nesting limited area model with the fourth GEM generation [22], with a log-hydrostatic pressure as the vertical coordinate following the terrain and with staggered vertical layers (which implies a shift between levels for momentum and thermodynamic variables). HRDPS system configuration was identified as the version v.4.0.0. TOR-GEM uses almost the same GEM model, but with a few other changes in model and configuration. Based on the work of [7] over Vancouver, the number of vertical layers is increased in the lower troposphere. The horizontal resolution is decreased. Current atmospheric parametrization in GEM for the subgrid turbulent transport considers only the vertical component with a scale-aware mixing length. Although the benefit to consider the horizontal components for subkilometer weather prediction has been clearly highlighted [23], recent studies suggest that more research is needed to develop a unified scheme across all scales [24,25,26,27]. Even with that limitation in GEM, previous studies with the model configured with 250-m grid spacing and the same boundary-layer scheme have shown a good performance regarding the near-surface meteorological prediction [7,8,28,29,30].

In urban areas, the atmospheric model vertical layers stand above the rooftop. The canopy morphology (street canyon), a displacement height for the wind profile of 2/3 of the mean building height, and built-up surfaces thermal and radiative properties are considered within the single-layer Town Energy Balance (TEB) scheme [31]. Thermodynamical turbulent fluxes provided to the atmospheric model originate from the street and the roofs as computed with an aerodynamic resistance network and with flow-evolving thermal local roughness lengths [32,33]. Three distinct surface energy budgets are calculated for roofs, roads and walls considering shadowing and radiation trapping effects in the street, leading to prognostic evolution of three-layer temperatures for each of these surfaces. Momentum exchanges and the wind values in the street canyon and above the roofs are based on a canopy-scale roughness length. Additional computations are provided for the evolution of water and snow reservoirs and water runoff.

Figure 1. Geographical location of the forecast domains, (a) over the globe with limits of ECCC’s NWP systems with 10 km (RDPS) and 2.5 km (HRDPS) grid spacing, (b) focus for the two inner grids for TOR-GEM1000 and TOR-GEM250 for which TEB is activated over the gray shaded urban areas, and (c) orography (Above sea level elevation in m) and main landscape highlights represented in TOR-GEM250.

Figure 1. Geographical location of the forecast domains, (a) over the globe with limits of ECCC’s NWP systems with 10 km (RDPS) and 2.5 km (HRDPS) grid spacing, (b) focus for the two inner grids for TOR-GEM1000 and TOR-GEM250 for which TEB is activated over the gray shaded urban areas, and (c) orography (Above sea level elevation in m) and main landscape highlights represented in TOR-GEM250.

Table 1. Summary of physical parametrizations in GEM used for the configurations used.

System	RDPS	HRDPS	TOR-GEM1000	TOR-GEM250
Grid spacing	10 km	2.5 km	1 km	250 m
Vertical levels below 2 km	13	13	28	28
First momentum level	40 m	40 m	10 m	10 m
Time step	450 s	60 s	30 s	12 s
Output Frequency	1 h	30 min	15 min	15 min
Land surface schemes	ISBA 1	ISBA	ISBA + TEB 2	ISBA + TEB
Microphysics	diagnostic scheme 3	double moment 4	double moment 4	double moment 4
Planetary boundary-layer	parameterization 5	parameterization 5	parameterization 5	parameterization 5
Max. turbulence mixing length (neutral)	200 m	200 m	200 m	57.5 m
Boundary-layer clouds	parameterization 5	parameterization 5	parameterization 5	parameterization 5
Shallow convection	parameterization 6	parameterization 6	parameterization 6	parameterization 6
Deep convection	parameterization 7	explicit	explicit	explicit
Radiation	parameterization 8	parameterization 8	parameterization 8	parameterization 8

¹ ISBA [34,35]. ² TEB [7,31]. ³ “Suncon” [36]. ⁴ “MY2” [37]. ⁵ “moistTKE”, vertical diffusion, Boundary Layer clouds [38,39]. ⁶ “Kuo transient” [39]. ⁷ “Kain Fritsch” [40]. ⁸ “Li Barker” [41].

Table 1. Summary of physical parametrizations in GEM used for the configurations used.

System	RDPS	HRDPS	TOR-GEM1000	TOR-GEM250
Grid spacing	10 km	2.5 km	1 km	250 m
Vertical levels below 2 km	13	13	28	28
First momentum level	40 m	40 m	10 m	10 m
Time step	450 s	60 s	30 s	12 s
Output Frequency	1 h	30 min	15 min	15 min
Land surface schemes	ISBA 1	ISBA	ISBA + TEB 2	ISBA + TEB
Microphysics	diagnostic scheme 3	double moment 4	double moment 4	double moment 4
Planetary boundary-layer	parameterization 5	parameterization 5	parameterization 5	parameterization 5
Max. turbulence mixing length (neutral)	200 m	200 m	200 m	57.5 m
Boundary-layer clouds	parameterization 5	parameterization 5	parameterization 5	parameterization 5
Shallow convection	parameterization 6	parameterization 6	parameterization 6	parameterization 6
Deep convection	parameterization 7	explicit	explicit	explicit
Radiation	parameterization 8	parameterization 8	parameterization 8	parameterization 8

¹ ISBA [34,35]. ² TEB [7,31]. ³ “Suncon” [36]. ⁴ “MY2” [37]. ⁵ “moistTKE”, vertical diffusion, Boundary Layer clouds [38,39]. ⁶ “Kuo transient” [39]. ⁷ “Kain Fritsch” [40]. ⁸ “Li Barker” [41].

Soil moisture refinement is critical to representing the Surface Urban Heat Island (SUHI) patterns [4] or the boundary-layer structure as shown for the representation of the low-level jet in Oklahoma city [42]. Initial conditions for soil moisture and temperature for the two soil layers in the vegetation scheme were then prescribed with fields with 2.5-km grid spacing from the Canadian Land Data Assimilation System (CaLDAS) [20,43].

In contrast, data assimilation techniques for urban surfaces are an emerging research topic but were not available at the time of this study. The urban fabric prognostic variables were then initialized in a simple way. Initialization occurs daily at 0200 Eastern Daylight Time (EDT, or local time corresponding to UTC$-4$ in summer). At this time of the night, the atmosphere is assumed statically stable. Therefore, the built-up surface temperatures are similar to the adjacent air layer (i.e., air canyon for the facades and the roads, and first atmospheric layer for the roof). Internal temperatures below the roads are set to the soil temperature analysis provided by CaLDAS. Finally, considering the absence of a detailed building energy model in this study, the initial indoor temperature is fixed at 23 °C and varies slightly through the day [44]. It is worth mentioning that the release of heat by air conditioning is neglected, but [45] found that it has a small impact on outdoor temperature for the summer of 2015 in Toronto as temperature values were moderate.

A large domain of 250 km${}^2$ was set up for the inner sub-kilometer forecasts for two reasons (Figure 1). First, to match the domain for experimental nowcasting during PanAm dedicated to mesoscale weather pattern tracking (iCAST, [46]). The second reason is to locate the city of Toronto far enough from the lateral boundaries of the model computational domain in order to ensure the quality of the urban forecasts. It was noted that a buffer of several grid points is needed in some conditions for the turbulent motions to stabilize, similarly to features demonstrated by [19] over London, UK.

It has long been recognized that numerical weather and climate prediction systems require even more detailed description of the surface when smaller grid-spacing is used [47] and when physical processes related to the presence of buildings and to the canopy morphology are modeled [48]. An important effort to provide realistic inputs to the surface schemes is presented hereafter.

2.1.2. Surface Description

For domains with 250 m and 1 km grid spacing, orography was determined with the Canadian Digital Elevation Data (CDED-50) and the Shuttle Radar Topography Mission for the United States. Soil texture for the two soil layers used was determined from Agriculture and Agri-Food Canada (AAFC) and the United States Department of Agriculture. Particular attention was brought to TEB parameters and vegetation types at high-resolution. The former method used at ECCC [4], not available for Toronto, was updated by a more general method suitable for ECCC’s systems covering North America.

An overview of the production of high-resolution ancillary data and samples of the results over Downtown Toronto are provided in Figure 2. The method relies on the work of [49] who initially developed ECCC’s GIS tools to automate the production of refined Land Use/Land Cover (LULC) data in urban areas. Microscale (5 m) rasterized maps of detailed urban and natural classes (e.g., paved road, highway, unpaved road, residential area, commercial 24/7 building, city park, school, etc.) were produced by merging pan-Canadian databases (Figure 2a,b). The program was originally designed to facilitate the change or addition of databases. Further developments provided for this study are described hereafter. The national Topographic Data Base (NTDB) was replaced by the vectorial CanVec database provided by Natural Resources Canada to describe detailed LULC and occupancy over Canada. More vegetation attributes are described through the Land Cover Circa 2000 (LCC2000) detailed vectorial database. As the domain used in this study is covering two countries, for the United States the National Land Cover Database (NLCD) 2006 provides 30-m raster information of vegetation types, as well as vegetation and urbanization densities further downscaled at 5 m. In addition, the OpenStreetMap (OSM) database from 2013 is used to extract vectorial objects for roads and more precise land use descriptions where available. Then, the population density (2012 census data from Statistics Canada) has been updated and is used to sub-categorize residential areas. Moreover, detailed individual building heights provided by the City of Toronto for a region covering downtown Toronto with high-rise buildings were used to improve morphological description. Algorithms include a series of consistency verification and priority rules when needed based sometimes on subjective visual inspections.

Final TEB and ISBA input parameters were obtained by assigning 5-m values (following a look-up table, e.g., a 5-m pixel of parking corresponds to 100% of road in TEB, a city park with trees corresponds to the vegetation class of mixed forest), and then by aggregating results as illustrated in Figure 2c–e. Canopy parameters such as the mean building height, the canopy aerodynamic roughness length based on [50] parameterization, and the ratio of wall surfaces to horizontal area following [51] are computed directly on the model grid (Figure 2f,g). Final LULC for the inner grid is presented in Figure 3.

In addition to the presence of the city, the region has heterogeneous terrain and is surrounded by large lakes (including the Great Lakes) with complex water circulations. Therefore, integrated NWP relies on more advanced watershed representation.

2.1.3. Water Surface Temperature

Thermally-driven atmospheric circulations such as lake and land breezes are also likely to be better represented when the water surface heterogeneity is realistically specified. Sea and water surfaces temperature is initialized in HRDPS with a global analysis at 10-km resolution based on buoys measurements and remote sensing products [52] and remains constant through the day. The same method is used for TOR-GEM systems except for the Great Lakes basin which is the dominant watershed in the region. This study took advantage of the development of a coupled ocean-atmosphere system for the Great Lakes [53] that has since become the operational Water Cycle Prediction System [54] for both the Great Lakes (WCPS-GLS) and the Saint Lawrence River watersheds. It uses, respectively, the Nucleus for European Modelling of the Ocean (NEMO) model [55] for the physical processes in the water column and Community Ice CodE (CICE) for the evolution of ice properties above [56] and include river streamflow into the lakes. Lake surface temperature is provided at 2-km grid spacing every hour for the five Great Lakes. Values over the grids with 250-m and 1-km grid spacing were obtained by the use of WCPS specific interpolation-aggregation tools and a linear interpolation was made for time steps in between two values. It is worth mentioning that the system also assimilates in the same way ice extent and thickness throughout the winter.

2.2. Comparison with Measurements

Data from both fixed and mobile measurements were used to evaluate the modelling system from 13 May and 31 August 2015 (corresponding to the period with the largest number of available measurements for PanAm) and their locations are mapped in Figure 4.

2.2.1. Surface Stations

First, most of fixed measurements sites are located in the region of Toronto, and near the Lake Ontario shoreline (Figure 4a). They are a combination of the operational ECCC surface station network that sparsely sample populated regions and the PanAm-specific stations that monitored GTA to a greater extent to capture urban features and track the inland propagation of the lake breeze [15,30]. Stations are classified into two categories, with 63 near-surface (including buoys) and 20 rooftop sites (details of stations are reported in Tables S1 and S2 of Supplementary Materials). Most of the instrumented buildings for rooftop stations were about 10 m tall and they were located in the urbanized regions of Toronto, Hamilton and Oshawa.

In order to evaluate model outputs near the surface, traditional diagnostic heights represented in the NWP systems are respectively about 2 m (10 m) above canopy level for the thermodynamic (wind) variables. A conceptual representation of the location of measurements and their corresponding modelled variables is presented in Figure 5.

Diagnostic wind at 10 m above ground level U${}_{diag}$ was computed directly in the TOR-GEM model depending on the canopy morphology. For grid points with buildings taller than 20 m, the wind speed is the street canyon wind speed [31]. For smaller buildings, the wind velocity value is either the atmospheric wind reduced by the logarithmic law, or a linear interpolation between the reduced wind and the canyon wind speed when 10-m height falls below the top of the canopy. For the ground stations, and due to unavailability of the canyon wind speed output data for the whole period of evaluation, only stations with wind measured at 10-m height are used in the objective evaluation.

Regarding the rooftop stations, velocities values U${}^{ref}$ at the first prognostic level ($z^{ref}$ = 40 m for HRDPS and z${}^{ref}$ = 10 m for TOR-GEM250) were corrected (during post-processing) to provide a fair comparison with the wind measured by the sensor U${}^{sensor}$ at z${}^{sensor}$ = 2.5 m above the roof (or at z${}^{sensor}$ = 5 m for one station as seen in Table S1 of the Supplementary Materials) using the following equation:

$$U^{sensor}=U^{ref}\frac{ln\frac{z^{sensor}+\frac{H}{3}}{z{0}_m}}{ln\frac{z^{ref}+\frac{H}{3}}{z{0}_m}},$$

(1)

with $z{0}_m$ the aerodynamic roughness and H the mean building height at the grid point (in m). Note that H is not defined in HRDPS. For temperature and dew point temperature, 5 m prognostic temperature in TOR-GEM250 is used for the full period comparison.

These features were implemented in the ECCC’s verification package presented in [20] to evaluate pan-Canadian short-range forecasts with the HRDPS. The tool performs several standard skill-scores and uncertainty is assessed with the bootstrapping method. In this study, Bias (BIAS) and standard deviation error (STDE) statistics are analysed for the prediction of near-surface temperature (T${}_{AIR}$), dew point temperature (T${}_D$), wind speed (U) and direction (U${}_{dir}$). Corresponding equations are reported in Appendix A.

Regarding precipitation forecasting skills, the Equitable Threat Score (ETS, [57]) is used to provide an overview of the model performance as detailed in Appendix A. It measures the forecasting performance of a system in its capacity to both represent an existing precipitation event and to avoid prediction of non-existing precipitation. Values fall between −1/3 and 1 and the closer that ETS is to 1, the better the precipitation forecast.

2.2.2. Mobile Measurements

The second database is the ECCC Automated Mobile Meteorological Observation System (AMMOS, [15]) which provides near-surface meteorological data for 21 days (mostly daytime) during the summer of 2015 for three instrumented cars. Routes were large arterial roads but not highways, sampling the central part of GTA (Figure 4b, and see Section S1 of Supplementary Materials). They were following three parallel routes from the shoreline and the downtown to about 19 km north, sometimes up to 35 km at King City (at the location of ECCC’s radar), and more rarely, up to Lake Simcoe. A few itineraries were covering the area parallel to the Lake Ontario shore (at about 19 km away). Data used in this study are from quality-controlled meteorological data averaged at 1-min interval (original measurements were at 1-s interval). Temperature and dew-point temperature from AMMOS are compared with model values at the nearest grid points. Given the fact that there were four times fewer quality-controlled measurements for wind than for temperature and that more analyses would need to be performed, wind data are not retained in this study. In addition to bias and STDE, the linear correlation coefficient r is computed for the scatter plots evaluation of these data (Appendix A).

3. Results and Discussion

3.1. High-Resolution Weather Prediction

A selection of results demonstrating the capacity of TOR-GEM250 to represent more realistic meteorological patterns is presented in this section. Due to the numerous outputs and physical processes to be considered with this system, the reader is invited to consult other specific companion studies where model outputs from this system were used. For example, Ref. [30] have evaluated thermal stress indices newly implemented in the system. Ref. [58] have analyzed the modeled lake-breeze system for two cases with light and strong inland penetration with both Lidar and surface measurements.Prediction of lightning occurrences has been evaluated at different resolutions [59].

3.1.1. Coastal Upwelling

The impact of the integration of improved lake surface temperature on the NWP system is investigated here. In particular, the northern Lake Ontario shore is subject to frequent upwelling resulting in cooler near-shore water temperature values in summer [60], and water flows between Lakes Erie and Ontario at different background temperatures. Figure 6a shows a map of water skin surface temperature on the grid with 250-m grid spacing on a day with offshore north-west synoptic wind (23 July 2015), which exhibits colder temperatures near the shore as a signature of an upwelling process in the watershed. Buoy measurement close to the Toronto islands is in agreement with the assimilated temperature. An important gradient of 9 °C can be observed between the north and the south shores of Lake Ontario for this event as the opposite downwelling process occurs on the southern shore. Lake Erie is warmer but still exhibits lower temperatures on the north shore. As can be seen in Figure 6b, the operational global analysis of water surface temperature (Figure 6b, see Section 2.1.3) exhibits less spatial variability over Lake Ontario and does not consider colder temperatures on the north shore. Similar differences can be seen for Lake Erie. This analysis is, however, used to prescribe daily Lake Simcoe (located in Figure 1c) and other water bodies surface temperatures. The hourly assimilation of the lake water temperature then allows for a more detailed representation of gradients occurring in the region.

3.1.2. Urban-Scale Heterogeneity

Urban-scale summertime 2015 general features are presented in this section. Averaged 15-min results were computed for canopy-top turbulent fluxes and near-surface meteorological variables for the period between 13 May and 31 August 2015. Time evolution of the different variables reflects the energy budget partitioning in and out urban areas [4,33].

Detailed results over the GTA are presented in Figure 7 at 2000 EDT corresponding to 18 h after initialization. At this time of the day slightly before sunset, spatial distribution of values of turbulent fluxes, surface temperature at different heights and dew point temperature varies. Figure 7d exhibits larger and positive sensible heat flux values over the GTA with values up to 60–70 W m${}^{-2}$ in the downtown region with high-rise buildings near the shore. More energy is released into the atmosphere after being stored during the day in the tall canopy where radiation fluxes were trapped previously. Over natural surfaces, sensible heat flux returns to a slightly negative value which is an expected behavior at that time [4]. Map of the latent heat flux (Figure 7e) shows values below 40 W m${}^{-2}$ in GTA as compared to values in the surroundings reaching 100 W m${}^{-2}$ over agricultural land outside the GTA, with minimum values in downtown and around the metropolitan area ring, where more impervious surfaces are present (such as industrial and commercial areas and large highways, as seen in Figure 3d). Over Lake Ontario, sensible and latent heat fluxes have small absolute values due to small vertical gradient of temperature and humidity at that time.

Skin radiative surface temperature (Figure 7a) exhibits a clear distinction of urban areas similarly to the sensible heat flux. Values reach 25 °C downtown but no more than 13 °C above Lake Ontario. These large temperature gradients through the warm season contribute to the enhancement of lake breezes processes in the region.

As felt by standing persons at about 2-m above ground level, air temperature intra-urban variability has clearly a different spatial distribution (Figure 7b) than the skin surface temperature. Larger values are still found in correlation with highly urbanized pixels, but the maxima are found over a region slightly north of the downtown and extends into surrounding urban areas due to advective air movements. In particular, maximum temperature values are seen in an industrial region with few vegetation and more absorbing materials at the south-west of the highways intersection (see Figure 3d, and identified in Figure S4 of the Supplementary Materials by the business improvement area in the city yard 7). A stable air layer is found over Lake Ontario with an air temperature minimum value of 15 °C, 2 °C larger than surface temperature. In between the lake and the GTA, a sharp temperature gradient is found near the shore. It is likely that colder air from the lake, and the air layer over the densely built areas warmed earlier in the daytime, were advected to the North by the often occurring lake-breeze flow.

As a result of more complex mesoscale circulations through the day, less spatial correlation between the urban area and the warmer air temperature is seen at the first prognostic atmospheric level at 5-m above the canopy level (Figure 7c). A region with diluted contours with larger temperature remains over the urban areas and illustrates the regular presence of a canopy urban heat island (CUHI) over Toronto. One can note, however, that the maximum values are also reached in the vicinity of the two green valleys respectively along the Don (East side of downtown) and Humber (west side of downtown) rivers (Figure 1c). A wider region than in Figure 7b of colder air is found near the shore, as cool air from the lake is advected aloft an internal urban boundary-layer.

Dew point temperature (Figure 7f) at the first atmospheric level exhibits the lowest values over the central areas of Toronto and north of it, even lower than over the lake. On the contrary to the lowest values, the largest dew point values are correlated with regions of the largest latent heat fluxes over the vegetation north and west of Toronto.

Daily time series for the corresponding domain and period are represented in Figure 8. Values are averaged based on the LULC clusters following Figure 3b urban fraction (URBF) thresholds. As seen in Figure 8a, near-surface air temperature reaches maximum values of 23.5 °C around 1500 and 1600 EDT for the most urbanized grid points. The effect of using initialization based on a lower-scale model is seen as air temperature is relatively homogeneous during the 7–8 h after the initialization over land. Later and until the end of simulations, values decrease as the urban fraction decreases. Regarding the quantitative description of the daytime CUHI, it is respectively 8 °C (2.5 °C) between the densest area of the city and the watershed (and the rural grid points). Early in the night, air temperature above the lake and in the surroundings of the city (with URBF < 30%, see Figure 3b) becomes similar but the regions forming the core of the city remains warmer, by about 1.8 °C. Standard deviation values among clusters are small, about 0.2–0.3 °C but are slightly larger for the densest cluster (URBF > 80%). This might be due to the effect of other urban descriptors or to the effect of the distance to the lake.

Figure 8b exhibits similar trends for the radiative surface temperature values with average maximum values of 34 °C. Diurnal variations are however much larger (up to 20 °C as compared to about 9.5 °C except for the lake water. The largest heterogeneity among clusters is seen for densest areas and for the lake water with the standard deviation reaching 3 °C. During daytime, resulting skin–air surface temperature difference presented in Figure 8c remains negative over water and increases progressively with the urban fraction categories with large positive values up to 10 °C. All clusters identify the air layer to return to a nighttime stable state but with a time lag of 4 h from the rural regions to the densest grid points. This is in agreement with the sensible heat flux evolution in the evening with remaining positive values until respectively 1900 LST (2300 LST) for the rural (most urbanized) category (Figure 8e). Sensible heat flux also exhibits the maximum 2 h later for the most urbanized pixels than for the rural grid points, with about double values (260 W m${}^{-2}$ and 120 W m${}^{-2}$). Heat trapped into the urban canopy during daytime is then released to the atmosphere in the evening and in the early night through sensible heat in addition to surface radiative cooling. By contrast, the latent heat flux becomes homogeneous through the urban regions at 2100 EDT (Figure 8d). Daytime values are maximum over rural regions (350 W m${}^{-2}$), with large difference among clusters.

The Bowen ratio (sensible to latent heat flux ratio) in Figure 8f highlights the evolution of the partitioning of the thermodynamic turbulent fluxes. During daytime when both fluxes are positive above land, it varies from a plateau not exceeding 0.5 for the fewer urbanized pixels (less than 30%) but reaching 1 and above for denser grid points. The more the urban fraction, the less regular is the curve. It is obvious for pixels with URBF > 80%, where two maximum values of about 3.5 are reached at 1230 EDT and 1800 EDT. In between, values slightly smaller highlight the time for the maximum of the latent heat flux at about 1530 EDT whereas maximum for the sensible heat flux occurs earlier at 1400 EDT.

Impact of lake-breeze flows on GTA’s intra-urban CUHI is further investigated in the following paragraph. During the Games period (from 8 July to 15 August), 95% of the days featured lake-breeze development [15]. Figure 9 presents maps averaged for days with low deformation lake-breeze over the lake Ontario (10, 11, 27, 28 July and 15 August 2015). Low deformation lake breeze means that the influence of the synoptic-scale wind is low and that the lake breeze propagates inland parallel to the coastline [61]. Maps of air temperature are provided during daytime (1500 EDT) and night-time (0100 EDT). During daytime (Figure 9a), air temperature maximum values are found over the city but at about 10 km north of the shoreline. During nighttime (Figure 9b), maxima are found close to the shore and concentrated around Toronto downtown, due to the interaction of the land breeze with the persistent heat released in urban areas (release of energy stored during the day). These results are in agreement with the expected climatology for coastal city under lake-(or sea-) breeze regime (see Figure 2 in [1]). In those situations, the maximum CUHI over land reaches respectively 5 °C (7 °C) for daytime (nighttime), but with a strong anisotropy for daytime with almost 14 °C as compared to the water side.

As shown in this section, the system is able to capture the heterogeneity of the surface-atmosphere interactions due to detailed representation of the urban, rural and water surfaces at very high resolution. The amplitude of the CUHI varies through the time of the day and the level of urbanization, but also to the meteorological conditions and the proximity to water bodies. Comparison with previous studies is therefore exploratory. For example, Ref. [62] found a CUHI of 0.7 °C (2 °C) respectively for daytime (night-time) for Toronto during a heat-wave with the Weather Research and Forecasting (WRF) model. However, intensities were based on a few grid points with locations that were not especially representative of the maximum/minimum in the region and the urban fabric heterogeneity for Toronto and its surrounding remained unrefined even with 333-m grid spacing. For Chicago which is another coastal urban agglomeration under the influence of the Great Lakes, and with the same model with 1-km grid spacing, similar nocturnal CUHI values of about 2 were found in [63]. Maximum temperature values during daytime were also found away from the coast, at about 15–30 km, with results sensitive to the description of the urban canopy. GTA’s summertime meteorology is influenced by the lake-breeze flow. Characteristics of such mesoscale circulations are further analyzed in the following section.

3.1.3. Mesoscale Meteorology

The TOR-GEM250 grid covers the region of Southern Ontario between Lakes Ontario, Erie, St. Clair, Huron (including Georgian Bay), and Simcoe (starting south of GTA and going clockwise) which has been intensively studied for its frequent occurrence of thermally-induced pressure gradients at the shores responsible for mesoscale circulations such as lake-breeze (and land-breeze) often responsible for convective patterns [61,64,65] and air chemistry features [64].

An illustration of the development of multiple lake-breeze fronts in the study region on 18 July 2014 is shown in Figure 10. Due to weak synoptic flow on that day, low-deformation lake-breeze circulations developed allowing the lake breeze to penetrate well inland around each lake, as seen on a visible satellite imagery at 1450 EDT from the MODIS Aqua radiometer (Figure 10a). Cumulus cloud development over the region in between the lake-breeze fronts can be identified with the spatial organization of significant surface-based convection, with more intense convection to the south. At 10 min apart, Figure 10b shows an analysis of the lake-breeze fronts as retrieved with the iCAST nowcasting system and using both GOES satellite, radar observations and in-situ measurements. These two products provide coherent information on the location of the fronts, including the merger of two fronts in one thin frontal zone in the Niagara peninsula (between Lakes Ontario and Erie). It should be noted that the location of lake breeze fronts is not identical to the visible cloud edges.

Forecasting of this situation is shown at the same time as the satellite with grid spacing of 250 m (Figure 10c,d). Predicted cloud cover exhibits the same organization as the visible satellite image. Multiple lake breezes with inland propagation can be identified with a South-East wind blowing over the GTA (Lake Ontario), a westerly wind blowing over the West of the domain (Lake Huron) and a southerly wind blowing over Lake Erie’s north shore. In agreement with MODIS, a region with less clouds along the southern Lake Huron’s coast can be seen (with the cloud edge forming a curve). This region corresponds to convective rolls in the region with moderate wind (due to the lake-breeze flow) upwind of the formation of inland cellular convective cells (as seen in Figure S2 in Supplementary Materials). Enhanced cloud development over the Oak Ridges Moraines can be seen as well, although slightly more than in MODIS. Clouds forming at the edge of Lake Erie’s front are slightly larger than in the north in both MODIS and TOR-GEM250. As was observed that day, the onset of the lake breeze occurred earlier on that shore (1200 EDT), allowing sufficient time for the convection to intensify.

Predicted surface winds with manually identified lake-breeze fronts annotated are presented in Figure 10d. Such evaluation of subkilometer-scale weather predictions can be difficult in a highly turbulent atmosphere and has some degrees of uncertainties [58,66]. While the synoptic wind was weak that day, indicators such as horizontal wind velocity, wind confluence, as well as associated updrafts and downdrafts, were found to be the most relevant (see Figures S2 and S3 in Supplementary Materials). Results of TOR-GEM250 suggest that lake-breeze fronts between Erie and Ontario did not fully merge, but the frontal zone is more clearly defined for the Lake Ontario shore (Figure S5c). The mesoscale analysis suggests that lake-breeze fronts have merged (Figure 10d), but this was not yet the case one hour before which suggests that the situation was evolving fast. Front collision in the Niagara peninsula is regularly occurring as described for two case studies with 200-m grid spacing simulations with WRF and mesoscale analysis [66]. Interestingly, as the grid spacing is increasing in the predictions for the nesting grids up to 10 km, cumulus clouds are also found in the same region but they are less numerous and tend to be larger, and the representation of the frontal zone becomes less detailed (Figure S4 in Supplementary Materials). TOR-GEM250 provides the best comparison with the satellite image. It is worth mentioning that the clouds in the satellite image are even smaller in some areas, which suggests that the model with 250-m grid spacing still does not have sufficient resolution to fully represent their organizations. Characteristics of the modelled lake-breeze front, more likely the turbulent structures than the location, might be impacted by the current one-dimensional planetary boundary-layer scheme used, as shown recently in [25] for idealized simulations using several turbulent schemes.

In summary, despite the complexity of the local scale atmospheric daytime evolution in summer, this section has highlighted the capacity of the TOR-GEM250 model to reproduce the impact of the lake-breeze currents propagation on the convective organization in the region. These mesoscale features are important in this region in summertime and can be associated with diverse effects, such as cooling of the coastal urban regions, air pollutant transport, or the occurrence of convective precipitation. Preliminary analysis of the precipitation patterns are examined hereafter.

With seasonal averaging (from 13 May to 31 August 2015) similar to that in the previous Section 3.1.2, Figure 11 shows afternoon and evening 6 h accumulated precipitation amounts between 1400 and 2000 EDT (Figure 11a), and between 2000 and 0200 EDT (Figure 11b). In the afternoon (Figure 11a), larger values are found generally inland, at some distance from the lakes coasts. Maxima are found in the south-west region and in the centre of the land areas, with largest values in the region of the city of London, ON. Lower values are found over Toronto and north of it, with slightly more accumulation around the Oak Ridges Moraine outcrop north of Toronto. Minimum is found over the middle of the Great Lakes. These patterns can be explained by the prevalence of lake breezes inland and air mass advection towards slightly more elevated regions (Figure 3c).

In the evening and later (Figure 11b), the accumulation of precipitation exhibits different spatial patterns. Large values are still found on the west part. A west–east oriented band with large values can be seen from Burlington to a region near Kitchener and the south of the small Conestogo Lake. It is possible that this pattern corresponds to the area where the lake breezes from different lakes sometimes interact and the convergence intensification is likely to produce more precipitation, but it also likely corresponds to advection of leftover convection from storms during the day. Over the Great Lakes, larger values were found than earlier in the day, likely due to end of the lake breeze flow and possibly the onset of land breeze flow. The Toronto area did not receive much precipitation that season as compared to the rest of the domain, as the lake breeze is often occurring and propagates beyond the agglomeration borders, and as latent heat flux is smaller over the GTA. To assess the realism of these maps, results were compared with Canadian Precipitation Analysis (CaPA) available for that period with 10-km grid spacing [67]. Unfortunately, the scale resolved by the background field is too coarse to provide a meaningful evaluation of TOR-GEM250 precipitation. Evaluation at surface stations will be therefore performed in the next section.

This section has provided an analysis of the capacity of the TOR-GEM250 system to provide detailed and realistic representation of the meteorology in and around a city for the GTA area. Further evaluation with in-situ measurements is provided in the following section.

3.2. Objective Evaluation

In this section, following the methodology in Section 2.2, evaluations are provided over the inner modelling domain with 250-m grid spacing (Figure 1) with both the surface monitoring network and the mobile meteorological measurements, and TOR-GEM250 is compared with the reference ECCC’s HRDPS system for short-range forecasting in Canada.

Time series evaluations are presented for the 24 h of the simulation for TOR-GEM250 and for the operational HRDPS system initialized at the same time, for both surface and rooftop stations. Corresponding local time starts at 0200 EDT.

3.2.1. Temperature and Dew Point Temperature

Statistics for temperature are shown in Figure 12. For surface stations, STDE values for 2.5-m above ground level $T_{2.5m}$ are consistently lower for TOR-GEM250 than in HRDPS (Figure 12a). Improvement is small but existing for the first 12 h. It is likely that the simplified initialization of TEB’s prognostic variables limits the performance of TOR-GEM250. STDE improvement increases however with time and reaches 0.15 °C by the end of the simulation, highlighting the better representation of spatial variability in temperature. Diurnal cycle of the errors is similar between the two experiments, with larger errors during daytime. Corresponding bias (Figure 12b) shows a slight underestimation of air temperature through the day for TOR-GEM250 of less than 0.5 °C and decreases for the second night to about 0.15 °C. In comparison, bias in HRDPS is either negative during the night with the same value than TOR-GEM250 after sunset and it features a positive maximum of 0.6 °C occurring after sunrise. Ground stations are located in mixed environments with significant vegetation fractions. Therefore, the role of evapotranspiration processes after sunrise in ISBA plays an important role in air temperature warming and might be overestimated in HRDPS.

With overall similar trends for the rooftop stations, statistical improvement of TOR-GEM250 over HRDPS is even larger than for ground stations, up to 0.3 °C in the STDE (Figure 12c) in late afternoon. There is less amplitude in the error diurnal cycle for TOR-GEM250, in particular for hours following sunrise, probably due to the impact of a better representation of physical processes over built-up surfaces. Bias is improved during nighttime (less negative). During daytime HRDPS shows substantial temperature overestimation (up to 0.2 °C). TOR-GEM250 exhibits an underestimation in a larger extent (up to 1.3 °C). As was mentioned earlier, the modeled variable from HRDPS and TOR-GEM250 above the roof are slightly different (see Section 2.2.1). As buildings are represented in TOR-GEM250 the variable compared with measurements is featured at 5-m above roof level.

The impact of such a difference between the temperature at the first prognostic level (T${}_{5m}^{prog}$) and a temperature diagnosed directly above the roofs at the sensor level (T${}_{2.5m}^{diag}$) is shown in Figure 13 on averaged for 5 days when the variable was stored from the system for further investigations. Unstable conditions are found with (T${}_2$${}^{diag}$) up to 1.35 °C warmer than (T${}_5$${}^{prog}$) for the maxima, which is likely to explain the daytime negative bias in Figure 12d. Although this is shown for a few days only, it is argued that it might have an impact for the four and a half months period.

Air water vapor content is traditionally evaluated in NWP with the dew point temperature, as it impacts numerous surface and boundary-layer processes. On Figure 14 is shown STDE and bias of T${}_D$ at 2.5-m above ground level and 5-m above roof level. Near the ground (Figure 14a), STDE are lower in TOR-GEM250 than in HRDPS during the first 6 h of the simulations, slightly larger during the next 6 h, and lower again until the end of the simulation by 0.15 °C. This result indicates that the spatial variability is better represented in the new system. Bias is relatively small (Figure 14b) for the two systems between −0.6 and 0.5 °C, but it is even smaller after 6 h of simulation in TOR-GEM250, almost nonexistent until the sunset after which the bias is improved by up to 0.4 °C in comparison to HRDPS. The bias has fewer diurnal variations as the fraction of impervious area in the surrounding of the station is better represented. The difference of representation of the buildings and roads between the two systems is obvious over the rooftop stations (Figure 14c,d). STDE errors remain improved during all the simulation and up to 0.3 °C during the afternoon, which is statistically significant. Due to the overestimation of the evapotranspiration of the soil and vegetation in HRDPS the bias for rooftop station is large and reaches 1.4 °C in the morning. In TOR-GEM250 the bias is improved by 1.3 °C during daytime and still by 0.4 °C during the second night. The diurnal cycle is less present in TOR-GEM250 as the fraction of vegetation is likely to be better represented.

To complement the objective evaluation, comparison of model outputs with mobile measurements is performed and provides sampling of different areas than the fixed mesonet. They covered the denser regions of the city from the downtown to city suburban rings (Section 2.2.2). Predicted versus measured scatter plots for all measurements are shown in Figure 15 for the mobile stations. Good correlations are found for the temperature (Figure 15a) and dew point temperature (Figure 15b) with relatively similar statistics for the three cars. The bias (STDE) for temperature is −0.52 °C (1.76 °C) which is rather consistent with results in Figure 12; and for the dew point temperature, respectively 0.13 °C (1.58 °C), which demonstrates a better STDE than for Figure 14. Good linear correlation coefficient is found for air temperature (0.88) and dew-point temperature (0.91).

This section has demonstrated the good performance of TOR-GEM250 for temperature and dew point temperature forecasting during the summer of 2015 for various locations. Heterogeneity in suburban areas is also improved as compared with HRDPS in particular for the second part of the runs. The wind speed is the next variable to be evaluated.

3.2.2. Wind Field

Wind vectors exhibit more details when forecasted at high resolution. The representation of built-up surfaces also tends to reduce the mesoscale wind [7]. These changes are highly dependant on the location of the measurement sites. This section is however exploring if objective evaluation can be informative on the new system performance.

As seen in Figure 16a for the ground stations, STDE errors are relatively similar for the two systems and are comprised between 1.3 m s${}^{-1}$ during nighttime and 1.75 m s${}^{-1}$ during daytime. The error increases earlier in TOR-GEM250 and remains above the one in HRDPS between 8 to 16 h of simulation (10 to 18 EDT). This time corresponds however to an improved bias as seen in Figure 16b with increased wind speed by up to 0.4 m s${}^{-1}$. Bias tends to remain smaller in TOR-GEM250 during the night, and slightly larger at sunset, despite a similar STDE than HRDPS. Differences in these figures might be due to the larger spatial variability of the wind field resolved in TOR-GEM250 and in particular due to the detailed representation of the urban canopy.

For the rooftop stations, STDE errors for TOR-GEM250 are significantly improved as compared with HRDPS during all the time series by up to 0.5 m s${}^{-1}$ in the middle of the day and of about 0.3–0.4 m s${}^{-1}$ during the nights (Figure 16c). Values are relatively similar in HRDPS for ground and rooftop stations. By contrast, values drop to 0.9–1 m s${}^{-1}$ during the night for TOR-GEM250, and 1.25 m s${}^{-1}$ during the day. This result indicates that TOR-GEM250 provides a better performance in the forecasting of the wind field heterogeneity in the coastal urbanized areas where rooftop stations are located. Corresponding bias values in Figure 16d suggest however a wind speed underestimation in both systems, by 1.2 m s${}^{-1}$ in HRDPS and 1 m s${}^{-1}$ in TOR-GEM250 in the middle of the day and respectively by 0.25 m s${}^{-1}$ to 0.4 m s${}^{-1}$ during the night. It is possible that the method used to adapt the forecasted wind field to the roof level decelerates the wind too much. It is worth mentioning that a revision of Equation (1) taking into account the static stability was tested for several cases resulting in relatively few differences. Evaluation of the wind direction was relatively similar between the two systems (Figure S5 of the Supplementary Materials).

Overall, results suggest objectively that the wind speed is generally improved in TOR-GEM250 for rooftop stations with a remarkable improvement of the STDE all day in comparison to HRDPS. STDE is increased for the ground stations in particular during daytime. It is possible that the very large wind heterogeneity in TOR-GEM250 in the complex summertime situations with cellular convection organization augment the statistical errors found at point stations. Summertime simulation with WRF with 333-m grid spacing over Toronto also revealed more difficulties to compare the wind with surface measurements than the temperature [62].

3.2.3. Precipitation

The spatial distribution of surface stations reported in Figure 4 is relatively dense in the centre of the domain, in and around the Toronto Greater area and along the lake Ontario’s coast. There are unfortunately less stations available in the western agricultural landscape that would help evaluating precipitation patterns in Figure 11. ETS (Section 2.2.1) is presented for results exceeding an accumulation of 5 mm during 6 h for the forecasting lead times 6, 12, 18 and 24 h, for the ground (Figure 17a) and rooftop (Figure 17b) stations. During the first part of the simulation (6 and 12 lead time), TOR-GEM250 exhibits a similar score than HRDPS or slightly below (between 0.22 and 0.29 for ground stations and between 0.25 and 0.3 for the rooftop stations). Later, ETS is significantly improved in TOR-GEM250. ETS is respectively improved for ground (rooftop) stations from 0.23 (0.23) in HRDPS to 0.3 (0.35) in TOR-GEM250 for precipitation accumulation between 12 and 18 h lead time and from 0.18 (0.19) in HRDPS to 0.38 (0.55) in TOR-GEM250 for precipitation accumulation between 18 and 24 h lead time.

The results indicate that TOR-GEM250 clearly outperforms HRDPS in precipitation forecasting after 12 h of simulations. This increase of performance with time indicates that TOR-GEM250 requires sufficient spin-up period or it needs improved initial conditions. Other results considering various precipitation amount thresholds confirm the good performance of the system for the second part of the simulations, with the exception of an overestimation of events with very small amount (see Figure S6 in Supplementary Material), which could be due to measurement uncertainties. The results also underline that TOR-GEM250 is more capable of representing the summertime convective precipitation occurring in late afternoon and evening.

4. Summary and Conclusions

Real-time NWP at high-resolution for urban areas and their surroundings is an important step towards the delivery of multiple services to city stakeholders, with model outputs directly available, or as a component of comprehensive decision-making tools. This paper presents results from the experimental TOR-GEM250 system set-up over the GTA and run in real-time in the context of the PanAm Games in Toronto prior to and during the summer of 2015. A detailed procedure to obtain geographical surface descriptors to be used by the surface schemes was followed. This method is now further used in all-scale NWP and air-quality systems with TEB at ECCC (currently up to 2.5 km). Forecasts were produced at grid-spacing as low as 250 m with optimized configuration. Their quality has been demonstrated both objectively, with results at surface stations and mobile measurements, and subjectively with analysis of small-scale meteorological patterns. The main results and perspectives are summarized as follows:

• The intra-urban scale micro-meteorology is predicted thanks to the representation of the urban fabric heterogeneity. Analysis of thermodynamic patterns over the summer season provides information on the characteristics of the GTA UHI which is important for both day, with a strong anisotropy, and night. On average for the season, CUHI is about 2 °C considering the rural surroundings. In details, it is about 1.8 °C during the night and 2.5 °C during the day as compared to the northern side of the city, and up to 8 °C during the day on the Lake Ontario side of the city. Contrasts can even reach 14 °C for some days in July when a low-deformation lake breeze penetrates inland. Urban-water surface contrasts are particularly large for this region which increase the impact of the city in the surface differential heating responsible for air advection and convergence zones. In the evening, heat stored in the urban canopy is restituted to the atmosphere through sensible heat flux in addition to surface radiative cooling;
• A further step towards an integrated system is reached in this study with the prescription of time-evolving Great-Lakes surface temperature based on realistic physical processes in the watershed. This additional factor in the representation of small-scale heterogeneity is influencing atmospheric circulations. Further development could include inline coupling of atmosphere and ocean processes at the same high-resolution (as is the case at ECCC for the global scale), or assimilation of in-situ buoys temperature. It is worth mentioning that surface temperature evolution over the Great Lakes have been successfully evaluated at a larger spatio-temporal scale [53], but more evaluation at the microscale would be needed as well. In particular there still can exists situation with some bias could be responsible for a lack of representation of the lake breeze system with potential impacts on the simulation of air pollutant dispersion in the GTA [9];
• Lake-breeze systems are important patterns in the GTA region and impact the local meteorology. Subjective comparison of TOR-GEM250 cloud patterns with satellite imagery on a day with lake breezes on multiple Great Lakes highlights the capacity of the system to provide realistic meteorological results which are not achieved by lower resolution systems. A complete objective analysis of lake breezes for the whole summertime season would be of interest;
• The real-time experimental system, TOR-GEM250, outperformed MSC’s HRDPS for all variables and in particular for the second part of the simulation, as evidenced by objective evaluation with 83 in-situ surface stations and mobile measurements. The system improves particularly the standard deviation errors through the day, which indicate that such a system would be valuable to represent small-scale variability in and around the city. This improvement reaches 0.3 °C for temperature and dew-point temperature and 0.5 m s${}^{-1}$ for wind speed. A method to evaluate model output with measurements located both at the ground and at rooftop level has been suggested in this study. The authors advise to carefully consider the sensor height to evaluate the model output with the variable that is intended to represent the measured value. Additional evaluation of the wind field’s intra-urban variability will also be required. Further investigation would be needed to assess the representativeness of the wind measured very close to surface (at about 2 m) as is the case for the AMMOS mobile stations and for some of the fixed stations.
  It also appears that the system is not optimal in the first few hours of the simulation. In the near-future, emphasis on the more accurate initialization of TEB prognostic variables and soil moisture in urban areas will be deployed. Method currently under investigation include the recycling of variables, and/or the implementation of the urban areas in the CaLDAS assimilation system [68];
• Precipitation patterns include many more details in TOR-GEM250 than in lower resolution systems and in particular for the late afternoon and evening events. Improvement of objective evaluation at stations with ETS reaching 0.55 for the evening period in urban areas indicates that the refined system better represents the meteorology. The city of Toronto and its surroundings do not receive a lot of precipitation as compared to the rest of the domain. This kind of system could be used to develop a very high-resolution precipitation analysis that would provide important information for many sectors of activities. For example, recent development of an ECCC’s 2.5 km ensemble-based precipitation analysis system still fail to fully represent summertime convective precipitation [69];
  Finally, the ongoing replacement of the microphysics scheme and, more generally, the modernization of the physics in GEM [70] for ECCC’s NWP systems is also likely to improve the forecasting of clouds and precipitation properties. Enhancement in several aspects of the predictions can be expected with an adaptation of the turbulent subgrid scheme to the grey zone of turbulence [24,71].