01-IMPORTANTE - WindTopo - Downscaling Nearsurface

Download as pdf or txt
Download as pdf or txt
You are on page 1of 21

1477870x, 2022, 744, Downloaded from https://rmets.onlinelibrary.wiley.com/doi/10.1002/qj.4265 by Readcube (Labtiva Inc.), Wiley Online Library on [13/03/2023].

See the Terms and Conditions (https://onlinelibrary.wiley.com/terms-and-conditions) on Wiley Online Library for rules of use; OA articles are governed by the applicable Creative Commons License
Received: 19 August 2021 Revised: 6 March 2022 Accepted: 8 March 2022 Published on: 29 March 2022

DOI: 10.1002/qj.4265

RESEARCH ARTICLE

Wind-Topo: Downscaling near-surface wind fields to


high-resolution topography in highly complex terrain with
deep learning

Jérôme Dujardin1,2 Michael Lehning1,2

1
School of Architecture, Civil and
Environmental Engineering, Swiss Abstract
Federal Institute of Technology in Predicting wind flow in highly complex terrain like the Alps is a challenge for
Lausanne (EPFL), Lausanne, Switzerland
all models. When physical processes need to be resolved in a spatially explicit
2
Institute for Snow and Avalanche
manner, grids with high horizontal resolution of a few hundred meters are
Research (SLF), Swiss Federal Institute for
Forest, Snow and Landscape Research often required and drastically limit, in many cases, the extent and duration of
(WSL), Davos, Switzerland the simulations. Many surface process models, like the simulation of hetero-
Correspondence
geneous snow cover across a season, however, need long time series on large
J. Dujardin, Institute for Snow and domains as inputs. Statistical downscaling can provide the required data, but
Avalanche Research (SLF), Swiss Federal no model can reach the desired resolutions effectively and provide temporally
Institute for Forest, Snow and Landscape
Research (WSL), Davos 7260, Switzerland. resolved wind speed and direction on highly complex topography. The assess-
Email: jerome.dujardin@slf.ch ment of the potential for wind energy in the Alps, a promising player in the
Funding information
energy transition, is an example where the current shortcomings cause strong
Innosuisse through the Swiss Competence limitations. We present “Wind-Topo”, a novel approach based on deep learning
Center for Energy Research, Supply of that discovers some of the interactions between high-resolution topography and
Electricity
coarser-resolution states of the atmosphere to generate near-surface wind fields
with a 50-m resolution. In our test case, we use a large number of measurement
stations in Switzerland to train the model and an operational weather predic-
tion model (COSMO-1) as predictor. Wind-Topo employs a custom architecture
that analyses the state of the atmosphere on various scales and associates it
with high-resolution topography. A dedicated loss function leads to good scoring
metrics as well as accurate wind-speed distributions at 60 independent sta-
tions used for a thorough validation. 50-m resolution wind fields are generated
efficiently and exhibit several expected orographic effects like ridge accelera-
tion, sheltering, and deflection. Furthermore, the bias and mean absolute error
from COSMO-1 at the alpine validation stations, which are 0.72 and 1.77 m⋅s−1
respectively, are reduced to −0.07 and 1.21 m⋅s−1 .

KEYWORDS
complex terrain, convolutional neural network, deep learning, downscaling, orographic effect, wind

This is an open access article under the terms of the Creative Commons Attribution-NonCommercial-NoDerivs License, which permits use and distribution in any medium,
provided the original work is properly cited, the use is non-commercial and no modifications or adaptations are made.
© 2022 The Authors. Quarterly Journal of the Royal Meteorological Society published by John Wiley & Sons Ltd on behalf of the Royal Meteorological Society.

1368 wileyonlinelibrary.com/journal/qj Q J R Meteorol Soc. 2022;148:1368–1388.


1477870x, 2022, 744, Downloaded from https://rmets.onlinelibrary.wiley.com/doi/10.1002/qj.4265 by Readcube (Labtiva Inc.), Wiley Online Library on [13/03/2023]. See the Terms and Conditions (https://onlinelibrary.wiley.com/terms-and-conditions) on Wiley Online Library for rules of use; OA articles are governed by the applicable Creative Commons License
DUJARDIN and LEHNING 1369

1 I N T RO DU CT ION at 171 French stations. These studies calibrated (trained)


point-to-point models that were only capable of predict-
In complex terrain, synoptic wind flows are trans- ing wind speed at specific locations. To obtain spatially
formed by their interaction with the topography and by distributed predictions, Manor and Berkovic (2015) devel-
near-surface processes generating local winds. Numer- oped a Bayesian-aided selection within a catalog of 0.5-km
ical weather prediction (NWP) models are increasingly resolution WRF simulations for northern Israel. The
more capable of reproducing such effects, thanks to the model replicates this high-resolution near-surface wind
constant improvement of their resolution, dynamical speed and direction from a coarser 4.5-km WRF dataset.
cores, and parametrizations. In highly complex terrain For California, Huang et al. (2015) trained a linear model
like the Alps, fine-scale models are required to resolve calibrated with WRF at 3 km to downscale the North
the terrain correctly and obtain flows that are close to American Regional Reanalysis (NARR: 32 km, 3-hourly).
observations. Despite the improvements in NWP models, A deep-learning super-resolution convolutional neural
a trade-off still exists between the extent of the simula- network (CNN) was employed by Höhlein et al. (2020) to
tion’s domain, its resolution, and the modeled period. recreate the ECMWF High-Resolution Forecast (HRES,
In highly complex terrain, given the required resolution, 9 km) from ECMWF ERA5 (30 km) for southern France.
large domains cannot be modeled for long periods. The With the exception of the latter study, which integrates
modeling of many surface processes like snow preferen- a terrain elevation model, these 2D-to-2D methods are
tial deposition and transport, however, requires long time constrained to the domain on which they are calibrated,
series of high-resolution wind fields on large domains because they do not incorporate descriptors of the topog-
(Mott et al., 2010; Vionnet et al., 2021). This is also true for raphy. Furthermore, their targeted resolutions are not
wind-energy potential assessment in such types of terrain, sufficient for highly complex topographies like the Alps.
where spatial variability is high. Other types of models like WindNinja (Wagenbrenner
A solution to the limitation in computational resources et al., 2016; Rios et al., 2018; Hilton and Garg, 2021) and
is to use NWP models at coarser resolution, succeeded TopoSCALE (Fiddes and Gruber, 2014) are based on phys-
by a statistical downscaling model. The former describe ical descriptions of wind–topography interactions and
the state of the atmosphere on a coarse topography and are not bound to the domain of calibration. WindNinja
the latter uses empirical relationships between these data is used frequently for wildfire propagation, but was used
and another source of information that describes the true recently by Vionnet et al. (2021) to model snowdrift. How-
nature of the wind. This ground-truth information can ever, given their generic nature, those methods do not
have various origins, but it may represent the best par- benefit from the valuable information of measurements or
tial description of wind flow in real terrain yet. Statistical high-resolution simulations and thus might not perform
downscaling is employed in many fields and on many well in particular terrain under particular conditions.
types of data. In complex terrain, wind and precipitation An interesting line of research took advantage of the
are the prime examples, because of their high spatial vari- vast network of measurement stations in Switzerland to
ability and importance in surface process modeling. A develop point-to-point methods that can generate spatially
vast repertoire of techniques is also available. Machine distributed predictions in the highly complex Swiss Alps.
learning is gaining importance in environmental sciences All of them use descriptors of the topography to predict
(Hsieh, 2009), and its extension deep learning begins to wind speed at the location of the measurement stations.
offer state-of-the-art results (Reichstein et al., 2019). When calibrated, the models can be used at other locations
For downscaling wind in complex terrain, two dif- and generate maps of wind speed. Etienne et al. (2010)
ferent approaches are employed, depending on whether employed general additive models (GAM) for daily max-
the ground truth comes from measurements or from imum wind speed using landform categories and Fischer
high-resolution NWP models. Philippopoulos and Deli- et al. (2015) employed gradient-boosted regression trees
giorgi (2012) trained an artificial neural network (ANN) to and topographic exposure for the same purpose. Foresti
predict the hourly wind speed at a certain station given the et al. (2011) and Robert et al. (2013), respectively, used
measurements from other stations in Greece. Dupuy et al. support-vector regression on mean wind speed and GAM
(2021) also used an ANN to estimate hourly wind speed at on monthly values, with various terrain descriptors (slope,
a French station, but from values predicted by the Weather difference of Gaussians, and directional derivative). These
Research and Forecasting (WRF) model with 3-km resolu- methods only use topographic descriptors as explanatory
tion. Similarly, Goutham et al. (2021) used analyses from variables and are thus limited to temporally static pre-
the European Center for Medium-Range Weather Fore- dictions. Without time-dependent inputs like the wind
casts (ECMWF, 9 km resolution) as a predictor for random vector from a lower-resolution model, such an approach
forest or gradient-boosting models to predict wind speed can only predict one value per location and not a time
1477870x, 2022, 744, Downloaded from https://rmets.onlinelibrary.wiley.com/doi/10.1002/qj.4265 by Readcube (Labtiva Inc.), Wiley Online Library on [13/03/2023]. See the Terms and Conditions (https://onlinelibrary.wiley.com/terms-and-conditions) on Wiley Online Library for rules of use; OA articles are governed by the applicable Creative Commons License
1370 DUJARDIN and LEHNING

series. Using the Consortium for Small-scale Modeling Intercantonal Measuring and Information System (IMIS)2
at 2- and 7-km resolutions (COSMO-2, COSMO-7)1 and and SwissMetNet3 networks at 7 and 10 m above ground
the directional-dependent terrain parameter Sx (Winstral level (m a.g.l.), respectively (Figure 1). From IMIS, 69 “ex-
et al., 2002), Winstral et al. (2017) trained a linear model to posed” stations are installed on summits or ridges and
generate maps of wind speed at 25-m resolution. As all sta- 86 “sheltered” stations are in locations that are protected
tions were used for calibration, the model’s validation did from the main winds (Lehning et al., 1999). The remaining
not reflect the performance at other locations. The results 166 (“other”) stations are uncategorized and are installed
at the stations, however, were promising and showed the in diverse terrain. As we want the model to predict wind
importance of considering the distribution of wind speed speed and direction, and given the noncontinuity of the
and not only aggregated metrics. From COSMO-2 and a direction, we use u and v, the horizontal components
slope parameter, Helbig et al. (2017) proposed a down- of the wind vector. We use hourly data from October 1,
scaling method that is calibrated from a large catalog of 2015–October 1, 2018. For the same times, the COSMO-1
simulations from the Advanced Regional Prediction Sys- analysis dataset describes the local state of the atmosphere
tem (ARPS) run at 30-m resolution. Using stations for val- surrounding each station, based on the coarse resolution
idation only, the method slightly outperforms COSMO-2 and smoothed topography (slope < 30◦ ) of the model
at exposed sites but is worse at sheltered ones. These (background of Figure 1). Wind-Topo combines this infor-
two methods offer the resolution and efficiency required mation with a high-resolution Digital Elevation Model
for the generation of long time series on large domains (DEM) to predict the station measurements u and v. In
mentioned previously. However, they do not consider the Wind-Topo, “local” means an area of about 21 × 21 km2
change in wind direction when going from coarse to fine (19 × 19 COSMO-1 pixels) surrounding each station, as
resolution and either were not validated spatially or were depicted in Figure 1. The model’s training consists of opti-
ineffective on lee sides. mizing the model’s parameters such that the predictions
We present a novel 2D-to-point statistical downscal- û and v̂ are as close as possible to u and v, using only
ing model: Wind-Topo, first of its kind and based on COSMO-1 and topographic inputs. For training, we use
deep learning, which uses kilometric-resolution NWP the first and last year of the period mentioned and a selec-
model outputs and high-resolution topography to predict tion of 261 stations. The remaining 60 stations and/or the
near-surface wind speed and direction in highly com- middle year are used to test the model, more precisely its
plex terrain. 261 Swiss measurement stations are used for temporal and spatial generalization capabilities.
calibration (training) and 60 different stations are used
for validation. In our test case, the 1.1-km resolution
COSMO-1 is used as predictor and near-surface wind fields 2.2 Model inputs
can be downscaled to 50-m resolution. The custom archi-
The model’s predictors come from either COSMO-1, with
tecture, which relies on multiple CNNs, as well as the
the aforementioned spatial extent and resolution, or the
custom loss function extract meaningful statistical rela-
high-resolution DEM. For the latter, we combined a 2-m
tionships between the low-resolution atmospheric state
DEM4 inside Switzerland and a 1-arcsec DEM5 (about
surrounding the point of interest, wind measurement,
36 m) outside the country. We resampled them such that
and high-resolution topography. The wind–topography
each COSMO-1 pixel contains exactly 21 × 21 DEM pix-
interactions discovered by the model are applied success-
els6 . As COSMO-1 (on a 0.01 × 0.0146 lat/lon grid) has a
fully at the validation stations and to a high-resolution
resolution in our whole Swiss domain of 1,113 m in lat-
grid, and exhibit ridge acceleration, sheltering, and
itude and 1,113 ± 22 m in longitude, the high-resolution
deflection.
topographic inputs have a resolution of 53 × 53 ± 1 m.
Accordingly, our input topographic patches, covering
21.147 × 21.147 ± 0.418 km, have a resolution of 399 × 399
2 M ET H ODS

2.1 Approach overview 2


https://www.slf.ch/en/avalanche-bulletin-and-snow-situation/
measured-values/description-of-automated-stations.html
As ground truth, we use wind speed and direction from 321 3
https://www.meteoswiss.admin.ch/home/measurement-and-
measurement stations spread across Switzerland from the forecasting-systems/land-based-stations/automatisches-messnetz.html
4
https://www.swisstopo.admin.ch/en/geodata/height/alti3d.html#
technische_details
1
https://www.meteoswiss.admin.ch/home/measurement-and- 5
https://www2.jpl.nasa.gov/srtm/
forecasting-systems/warning-and-forecasting-systems/cosmo- 6
This number of pixels (21 × 21) is coincidentally the same as the size of
forecasting-system.html the COSMO-1 input domains. There is no link between them.
1477870x, 2022, 744, Downloaded from https://rmets.onlinelibrary.wiley.com/doi/10.1002/qj.4265 by Readcube (Labtiva Inc.), Wiley Online Library on [13/03/2023]. See the Terms and Conditions (https://onlinelibrary.wiley.com/terms-and-conditions) on Wiley Online Library for rules of use; OA articles are governed by the applicable Creative Commons License
DUJARDIN and LEHNING 1371

F I G U R E 1 Location of the
measurement stations used for
training or validation. The
background shows the topography
used in COSMO-1 and the squares
around selected stations represent the
spatial extent of Wind-Topo’s inputs
to predict the wind at those stations.
The red and grey squares show the
test domain of Section 3.4 and the
domain covered in Figure 2,
respectively [Colour figure can be
viewed at wileyonlinelibrary.com]

pixels. The exact same extent is covered by 19 × 19 ridges and summits. We thus had to use other topographic
COSMO-1 pixels. The two upper left panels of Figure 2 inputs that show the wind–topography interactions more
show an example of such input data from COSMO-1 and directly. Two such descriptors that are well known were
DEM. initially employed: the topographic position index (TPI),
From COSMO-1, we select 2D fields of horizontal indicating the relative height of a location with respect to
components of the wind vector uc and vc (with “c” stand- its neighborhood, and the maximum upwind slope (Sx )
ing for COSMO) for five vertical layers of the model (Winstral et al., 2002), indicating how exposed or sheltered
(terrain-following coordinates) that have average heights a place is, given a specific wind direction. Both of them
of 10, 89, 293, 589, and 1,164 m.a.g.l. We also use the per- share one important parameter: the radius considered for
pendicular deviation of the wind vector from the layer the neighborhood (typically 500–4,000 m). Trained with
surfaces, denoted w′ , which indicates on the five layers them, Wind-Topo obtains scores that are almost as good
where the flow follows the terrain (w′ = 0), separates as the ones presented in Section 3. However, the pre-
upwards (w′ > 0), or converges towards the ground (w′ < dicted wind fields exhibit a strong dependence on the
0). w′ also indicates areas of convection or subsidence. The radius parameter and look like copies of the descriptors,
atmospheric stability is considered via the vertical gradient even if several descriptors with various radii are pro-
of potential temperature between the five layers: Δ𝜃∕Δh. vided. To avoid this problem, we created a set of new
The last inputs are the ground surface sensible heat flux qs topographic descriptors that are free of such parameters
and the elevation of the model’s terrain zCOSMO . Concern- and let Wind-Topo decide how to process them. Given its
ing the high-resolution topographic descriptors, elevation CNN-based nature, Wind-Topo can decide how to inte-
ztopo and the associated slope and aspect are employed. grate them spatially, thus creating an internally controlled
and dynamic radius parameter:
{
2.3 New topographic descriptors E+ = max(sin(𝛼), 0),
E− = min(sin(𝛼), 0),
Ideally, deep learning could allow an end-to-end approach,
where the model learns by itself the operations to apply where
(1)
to the DEM (ztopo ). However, the rather small variability 𝛼 = arctan(tan(slope) cos(𝛿)), with
in topographic inputs (only 321 unique ztopo images, one 𝛿 = arctan2(−vc , −uc ) − aspect.
for each local patch) is far from sufficient for obtaining an
end-to-end model with good spatial generalization. Such
a model could only generate smooth wind fields, lacking Figure 2 shows the variables zCOSMO , ztopo , slope, and
small-scale features, and only increased uc and vc around aspect for the input patch centered around the station
1477870x, 2022, 744, Downloaded from https://rmets.onlinelibrary.wiley.com/doi/10.1002/qj.4265 by Readcube (Labtiva Inc.), Wiley Online Library on [13/03/2023]. See the Terms and Conditions (https://onlinelibrary.wiley.com/terms-and-conditions) on Wiley Online Library for rules of use; OA articles are governed by the applicable Creative Commons License
1372 DUJARDIN and LEHNING

F I G U R E 2 Some topographic inputs for a training station (ELA-1). zCOSMO is the topography in COSMO-1 and ztopo is a 53-m DEM;
slope and aspect are computed from ztopo . (E+ + E− ) reflects places exposed to the wind (red) or sheltered from it (blue). Δutan shows the offset
to apply to u to obtain a flow that is progressively tangential to the terrain as the slope increases, given wind speed and direction [Colour
figure can be viewed at wileyonlinelibrary.com]

ELA-1 (Section 3). From slope, aspect, and the values of uc , descriptors address the challenge of wind deflection: how
vc at the center of the patch we generate six topographic to change uc and vc to obtain a flow that turns when fac-
descriptors. First, we compute E+ and E− (Equation 1) ing a steep slope. We provide the model with the offsets
representing respectively how much a location is exposed for uc and vc , Δutan and Δvtan , needed to obtain such an
to or sheltered from the wind, given its direction. The effect. Equation 2 shows that the steeper the slope, the
best results were obtained when this direction is com- more the flow will have to follow the topographic con-
puted from the second layer above ground in the COSMO-1 tours (be tangential). This is performed via a rotation of
inputs (89 m.a.g.l.). Contrary to Sx , E+ and E− do not incor- the vector (uc , vc ) by an angle 𝛽. This angle is calculated
porate a search distance and are purely local. They are using the slope angle and 𝛿, which is the same angle
computed as the sine of the vertical angle 𝛼 that the wind as in Equation 1. Using the observations as guidance,
vector would have to change to in order to be parallel to Wind-Topo learns how much, where, and when these
the slope while keeping its azimuth. This angle is sim- exposure, sheltering, and theoretical deflections should be
ply computed from the slope angle and from 𝛿, the angle employed:
between the wind direction and the aspect. Splitting this
{
exposure/sheltering into two allows the model to find spe- Δutan = (cos(𝛽) − 1)uc − sin(𝛽)vc ,
cific rules for each of them independently. When summed
Δvtan = sin(𝛽)uc + (cos(𝛽) − 1)vc ,
together (as in Figure 2) and integrated over a certain
radius, the result resembles Sx . As a CNN cannot multi- where
(𝜋 ) (2)
ply some of its inputs together directly, and because our 𝛽= 2
− |𝛿| sign(𝛿) sin(slope), with
model showed difficulties in learning it, we provide E+ uc , 𝛿 ∈] − 𝜋, 𝜋] from Equation 1.
E+ vc , E− uc , and E− vc to the model. The two remaining
1477870x, 2022, 744, Downloaded from https://rmets.onlinelibrary.wiley.com/doi/10.1002/qj.4265 by Readcube (Labtiva Inc.), Wiley Online Library on [13/03/2023]. See the Terms and Conditions (https://onlinelibrary.wiley.com/terms-and-conditions) on Wiley Online Library for rules of use; OA articles are governed by the applicable Creative Commons License
DUJARDIN and LEHNING 1373

F I G U R E 3 Chain of operations to obtain the predicted û and v̂ at the center of the patch of input data [Colour figure can be viewed at
wileyonlinelibrary.com]

2.4 Data flow As presented in Figure 3, we have low-resolution


COSMO-1 data (full and zoom patches), high-resolution
Given the large number of training data (261 stations topographic data (full and zoom patches), and values for all
× 17,520 timesteps × n variables), treating the 399 × 399 those descriptors at the exact location of the station. Time
inputs is an enormous task (3 TB of data). As Wind-Topo of day and day of the year are also provided as pointwise
predicts the wind at the center of the input patch, the information via a cosine–sine transformation for continu-
more relevant information is located close to the cen- ity. Wind-Topo has two branches: one that predicts u and
ter, while the edges of the patch should provide only one that predicts v. Both are independent when making
coarser information. To decrease the amount of informa- a prediction, but they are trained jointly via a common
tion significantly (by a factor of 13) and to help the model loss function (Section 2.6) that optimizes their parame-
focus on high-resolution information towards the center ters simultaneously. This double architecture gives the best
and on progressively coarser information at a distance, results and allows the model to consider u and v jointly.
we apply two consecutive operations to all topographic A widespread technique in machine learning is to aug-
inputs (Figure 3). First, each input is split into two: on ment the amount and diversity of the training data arti-
the one hand the full-extent input is resized to a resolu- ficially. For image processing purposes, for example, with
tion of 77 × 77 pixels, and on the other hand a crop of the CNNs, it is common to rotate, flip, crop, offset, and rescale
central 77 × 77 pixels (a zoom) is extracted. Second, we the images randomly. In our situation, as we only have
apply a foveal blur (see bottom right of Figure 3), which 261 different input topographies, this augmentation is cru-
smooths the information as we move away from the cen- cial. 2D data were randomly rotated by one of 72 angles
ter. For consistency in the model’s architecture, which (0–360◦ every 5◦ ) and all values related to u, v, uc , vc ,
is described below, the same splitting is applied to the and aspect were transformed in a consistent manner that
COSMO-1 inputs. artificially rotates topography and wind together. To stay
1477870x, 2022, 744, Downloaded from https://rmets.onlinelibrary.wiley.com/doi/10.1002/qj.4265 by Readcube (Labtiva Inc.), Wiley Online Library on [13/03/2023]. See the Terms and Conditions (https://onlinelibrary.wiley.com/terms-and-conditions) on Wiley Online Library for rules of use; OA articles are governed by the applicable Creative Commons License
1374 DUJARDIN and LEHNING

F I G U R E 4 Deep-learning architecture of Wind-Topo. Each “CNN low-res” and “CNN high-res” is unique to each input variable. The
two distinct “CNN fusion” blocks treat the concatenated tensors with full or zoom extent, respectively [Colour figure can be viewed at
wileyonlinelibrary.com]

consistent with wind–topography interactions like prefer- Figure 5 details the building blocks of the
ential deflection to the right due to Coriolis, flip, crop, and deep-learning architecture. The CNNs have a classi-
rescale were avoided and proved to be ineffective when cal alternation of convolution layers and pooling layers
tried. Offsets of elevation were investigated but did not in order to extract increasingly larger and higher-level
help. Data augmentation by rotations also allowed us to features from each piece of input data and from the con-
obtain a model with a good rotational invariance (similar catenated data. Special attention was given to the size of
predictions when the inputs are rotated by an angle and the inputs and kernels, as well as to the strides, to avoid
the prediction is then rotated oppositely). padding on the edges and lateral shifts of information. As
Wind-Topo gives a prediction for the exact patch’s center,
it is essential that each operation leads to an odd number
2.5 Deep-learning architecture of rows and columns and that the kernels fall centered on
the central pixel. The sigmoid linear unit (SiLU) function
Figure 4 shows the constitutive blocks of the model, the (silu(x) = x𝜎(x), where 𝜎(x) = 1∕(1 + e−x ) is the logis-
data dimensions, and how the various types of data are tic sigmoid: Hendrycks and Gimpel, 2016) was the best
associated to predict u and v. Each input variable is pre- activation function, outperforming the standard rectified
treated by a dedicated CNN (one CNN for uc , one for vc , linear unit (ReLU) function slightly, most likely because
...) to generate a tensor. Tensors with similar spatial extent of its continuity (important for a regression problem).
(full or zoom) are concatenated and treated by another A particularity of the CNNs treating the high-resolution
type of CNN. As seen in Figures 4 and 5, these latter CNNs topographic data is the addition of a directional mask
and the ones pretreating the high-resolution topographic after the first activation. This mask, which sets to zero the
descriptors generate an additional vector that contains area that is downwind from the patch’s center, has two
values at the exact center of the patch. This procedure advantages: (i) it helps the model understand the notion
retains important high-resolution information at mul- of wind direction and how to deal with windward and
tiple stages of the CNNs. Classically, the outputs of the leeward information, and (ii) it augments the number
CNNs on the concatenated tensors are fed into fully con- of distinct topographic inputs. As mentioned previously,
nected neural networks to generate two 128-entry vectors. central values (red) are extracted before each pooling layer
Finally, all vector data are concatenated and treated by a and before several convolution layers in order to retain
neural network having a unique output, which is û or v̂ low-level high-resolution information about the center of
depending on the model’s branch. the patch.
1477870x, 2022, 744, Downloaded from https://rmets.onlinelibrary.wiley.com/doi/10.1002/qj.4265 by Readcube (Labtiva Inc.), Wiley Online Library on [13/03/2023]. See the Terms and Conditions (https://onlinelibrary.wiley.com/terms-and-conditions) on Wiley Online Library for rules of use; OA articles are governed by the applicable Creative Commons License
DUJARDIN and LEHNING 1375

F I G U R E 5 Constitutive blocks of deep learning in Wind-Topo. The convolution blocks (yellow) indicate kernel size (3 × 3), stride s,
padding size p, and number of filters (e.g., ×16). The average-pooling blocks (purple) depict the same information. Vectors of central
information (red) are extracted at various places. The first convolution layer of “CNN high-res” has an additional directional mask based on
the wind direction from COSMO-1 (e.g., mask in lower right corner). The two types of fully connected neural network are detailed in the
lower left corner [Colour figure can be viewed at wileyonlinelibrary.com]

TA B L E 1 Hyperparameters used to train the model best architecture, input data, and training procedure in
Parameter Value a large design space. In Appendix S1, we present sev-
eral model variants and their performance. This ablation
Optimizer ADAM
study shows the relative importance of various design
𝛽1 , 𝛽2 (0.9, 0.999) decisions by removing one model component at a time.
Initialization Xavier, uniform, gain=1 It appears that the current design is not yet optimized
Learning rate 5 × 10−6 , no weight decay with respect to model complexity and that some building
blocks, despite improving performance slightly, signif-
Early stopping 10th epoch
icantly increase computational requirements. Future
Batch size 128
research will investigate the trade-off between model
ADAM: (Kingma and Ba, 2015). complexity and downscaling performance in more
Xavier: (Glorot and Bengio, 2010) detail.

Training the model consists of optimizing the values of


the model’s parameters (weights and biases for neurons, 2.6 Loss function
filters in a CNN). One needs to adopt a strategy to do so,
in order to reach the best performance on the training and The design of an appropriate loss function was critical to
validation sets. The various elements that compose this obtain the current performance of Wind-Topo. Using any
strategy are generally referred to as the hyperparameters standard loss function like mean-square or absolute error
of the model. For Wind-Topo, an iterative procedure led to (MSE, MAE) led to distributions of û, v̂ and the associated
the hyperparameters presented in Table 1. ̂ that were too narrow. The model
horizontal wind speed vel
The version of Wind-Topo presented here is the result chose to decrease the range of its predictions, favoring val-
of extensive experimentation that aimed at finding the ues near the mean value. This is a well-known problem in
1477870x, 2022, 744, Downloaded from https://rmets.onlinelibrary.wiley.com/doi/10.1002/qj.4265 by Readcube (Labtiva Inc.), Wiley Online Library on [13/03/2023]. See the Terms and Conditions (https://onlinelibrary.wiley.com/terms-and-conditions) on Wiley Online Library for rules of use; OA articles are governed by the applicable Creative Commons License
1376 DUJARDIN and LEHNING

statistical downscaling (Winstral et al., 2017). subpixel resampling is performed efficiently using convo-
lutions with Lanczos kernels (Burger and Burge, 2009). If
1∑
N
the same method is applied to the topographic data as well,
= 𝜏i [(ûi − 𝛽i ui )2 + (̂vi − 𝛽i vi )2 ], it is possible to generate wind fields at a higher resolution
N i=1
{ than the 53-m DEM input data. However, as the model was
𝜖 + veli 𝜏 if ̂ i ≥ veli ,
vel
𝛽i = , 𝜏i = (3) trained on this particular DEM, higher resolutions lead
̂i
𝜖 + vel 1−𝜏 if ̂ i < veli .
vel to smooth outputs with no additional information. In our
case, wind fields generated at a 50- to 100-m resolution
The countermeasure was to compose the loss func- contain the highest level of detail.
tion  in Equation 3, which is essentially a MSE over a The complex pipeline needed to obtain all the required
batch of training data of size N, with a scaling by 𝛽i of input data for every grid point can be performed by the
the target values ui and vi based on the overestimation CPU, in parallel with the predictions performed by the
or underestimation of the predicted velocities vel ̂ i . Addi- GPU. In our case, for Wind-Topo coded with Pytorch
tionally, the term 𝜏i , inspired by the pinball loss function (Paszke et al., 2017),7 it is possible for any domain size
(Takeuchi et al., 2006), allowed the model to remove the to obtain a perfect parallelization of these two tasks,
consistent negative mean bias error (MBE) of vel. ̂ Param- without noticeable overhead and with a GPU running
eters 𝜖, controlling the strength of the scaling, and 𝜏 were constantly at full capacity. As the model is pointwise
tuned iteratively to their respective optimal values of 4 and and time-independent, it is easily parallelizable on mul-
0.425. A lower value of 𝜖 led to an unstable behavior and tiple GPUs and/or machines with, for example, one
incorrect distributions, while larger values led to squeezed GPU computing one part of the domain or a subpe-
distributions. The current value of 𝜏 allowed us to obtain riod and another GPU computing the rest. On a NVIDIA
an unbiased model on the test set. The quality of the dis- RTX2080Ti GPU, one prediction takes 0.6 ms and a
tributions of predicted wind speed is essential for many 300 × 300 domain requires 54 s. This is 30% faster than
applications of high-resolution wind fields. This loss func- on a previous-generation GTX1080Ti, which indicates that
tion is the result of an intense effort to find a differentiable Wind-Topo would be even faster on the latest-generation
function that guides the gradient-descent algorithm of the high-performance GPUs.
model’s optimization efficiently towards distributions that
closely match the measured ones, while obtaining very low
MAE and MBE and high Pearson correlation coefficients 3 RESULTS
between predictions and measurements. These three scor-
ing metrics reflect important qualities of the predictions. This section shows the overall performance of Wind-Topo
A low MBE guaranties no bias between predictions and in highly complex terrain and quantifies its aptitude for
observations, while a simultaneous low MAE and high capturing subgrid-scale wind–topography interactions like
correlation indicate that the model is accurate for a large ridge acceleration, sheltering, and deflection. After train-
range of observations. The development of Wind-Topo was ing Wind-Topo on 261 stations and years 1 and 3, we tested
guided by using these three scoring metrics on the valida- it on separate datasets to assess its performance at new
tion set, and Section 3 focuses on them. locations and/or over a new time period: test set 1 (60 sta-
tions, years 1 and 3), test set 2 (261 stations, year 2), and
test set 3 (60 stations, year 2) being the most informa-
2.7 Making predictions tive because it reflects the spatio-temporal generalization.
In this section, we will evaluate Wind-Topo first quantita-
Wind-Topo was trained to make a point prediction from tively and then qualitatively by discussing some predicted
centered patches of COSMO-1 data and a high-resolution wind fields. From now, when COSMO-1 is mentioned,
DEM covering exactly 21.147 km. Thus, when Wind-Topo it should be understood that a bicubic interpolation was
is used to generate a wind field on, for example, a 100-m applied at the desired locations. Importantly, the 60 test
grid, this is performed in a point-by-point manner. Each stations were selected by an algorithm that ensures an
predicted grid point is an independent output of the model, equitable evaluation. First, all stations were tagged with an
to which we provide the centered patches of data. In this eight-dimensional vector consisting of (1) mean measured
100-m case, COSMO-1 needs to be resampled every 100 m wind speed, (2) MAE (of COSMO-1 with respect to mea-
to avoid 11 consecutive grid points receiving the same surement), (3) MBE, (4) correlation, (5) normalized MAE
COSMO-1 inputs, which would lead to a chessboard-like
pattern in the resulting high-resolution wind field. This 7
https://github.com/pytorch/pytorch/tree/v1.8.1
1477870x, 2022, 744, Downloaded from https://rmets.onlinelibrary.wiley.com/doi/10.1002/qj.4265 by Readcube (Labtiva Inc.), Wiley Online Library on [13/03/2023]. See the Terms and Conditions (https://onlinelibrary.wiley.com/terms-and-conditions) on Wiley Online Library for rules of use; OA articles are governed by the applicable Creative Commons License
DUJARDIN and LEHNING 1377

TA B L E 2 Performance of COSMO-1 (C) and Wind-Topo (WT)

Station Corr. MBE MAE nMAE MAEdir


(vel) Model (DN) (m⋅s−1 ) (m⋅s−1 ) (DN) (deg)
Exposed C 0.51 −0.25 1.83 0.56 41.3
(3.43 m⋅s−1 ) WT 0.66 −0.05 1.49 0.45 32.5
Sheltered C 0.42 1.90 2.27 1.43 67.0
−1
(1.67 m⋅s ) WT 0.45 −0.04 1.03 0.63 58.1
Other C 0.59 0.18 1.34 0.61 50.5
(2.46 m⋅s−1 ) WT 0.64 −0.16 1.16 0.51 47.5
Alps C 0.50 0.72 1.77 0.90 56.4
(2.35 m⋅s−1 ) WT 0.55 −0.07 1.21 0.56 50.1

Note: Each value is the average of the scores at each test station, computed for the test period. The left column gives the average measured wind speed for each
station category.
Abbreviations: Corr., Pearson correlation coefficient; MAE, mean absolute error; MBE, mean bias error; nMAE, normalized MAE (MAE at a station divided by
its average measured wind speed); MAEdir, mean absolute error of wind direction.

(MAE at a station divided by its average measured wind metrics and at all station types. The correlation is improved
speed), (6) elevation, (7) latitude, and (8) longitude. Then, the most at exposed and other sites, while the large bias at
in an iterative manner, the algorithm splits the stations sheltered sites is completely removed, leading to a much
randomly into training and test sets until there is a good lower MAE. The normalized MAE (nMAE) allows us to
representation of this eight-dimensional space in both sets compare all types of station and again shows the lack
and a similar performance of COSMO-1. Concerning the of performance of COSMO-1 at sheltered sites and the
separation of training and test times, our investigation capability of Wind-Topo to correct it. Finally, Wind-Topo
showed that a random split cannot reflect the temporal also corrects wind directions, as expressed by the MAE
generalization: it is easy for a model trained on hours h of wind direction (MAEdir), computed from events with
and h + 2 to predict hour h + 1. Consequently, we kept the a wind speed greater than 1 m⋅s−1 . As Wind-Topo was
complete intermediate year for validation. In this section, designed for complex terrain, the lower line in Table 2 pro-
we employ the same colors as in Figure 1 for stations that vides overall scores for the 44 test stations located in the
are exposed (blue), sheltered (orange), and uncategorized Alps.
(“other”, green). The scores for Wind-Topo were obtained from predic-
tions of the trained model. Figure 6 depicts some of them
during the training. After each epoch (pass through the
3.1 Aggregated scores training set), we evaluated the model on the four datasets.
Correlation, MAE, and MAEdir improved progressively
Table 2 shows various scoring metrics for COSMO-1 and until epoch 9, the model state of which was used for all
Wind-Topo on test set 3. Each score is the average of the predictions in this work. Later on, a typical overfitting
scores at the stations belonging to that category. This is is observed, where the spatial generalization decreases.
more meaningful than computing the scores on all the Interestingly, there is no overfitting with respect to the
data, as, for example, correlations computed on all the data temporal generalization. When training the model further
are significantly higher but are potentially misleading. (not shown), the scores on test set 2 keep improving (MAE
We can observe the importance of the station’s clas- < 1 m⋅s−1 , correlation > 0.8). Temporal generalization was
sification in the left column: the average wind speed easily obtained: many (simpler) models can give accurate
at exposed sites is more than twice that at sheltered predictions on test set 2. In other words, the corrections to
sites. Given its resolution, COSMO-1 cannot discriminate COSMO-1 can be learned easily for a certain period and
between them and has a large positive bias (1.9 m⋅s−1 ) at applied successfully to new periods. Furthermore, a small
sheltered sites. Surprisingly, the highest wind speeds are number of timesteps is required. All results presented
predicted there (Kruyt et al., 2018). Nevertheless, its perfor- are based on Wind-Topo trained with only half the train-
mance is quite high at exposed locations and in flat lands ing times (random selection), and no improvement was
and large valleys (a considerable part of the “other” sta- observed using all times. The strong temporal generaliza-
tions). Wind-Topo improves COSMO-1 significantly for all tion offers a different application: if a short measurement
1477870x, 2022, 744, Downloaded from https://rmets.onlinelibrary.wiley.com/doi/10.1002/qj.4265 by Readcube (Labtiva Inc.), Wiley Online Library on [13/03/2023]. See the Terms and Conditions (https://onlinelibrary.wiley.com/terms-and-conditions) on Wiley Online Library for rules of use; OA articles are governed by the applicable Creative Commons License
1378 DUJARDIN and LEHNING

F I G U R E 6 Quality of the predictions of Wind-Topo compared with COSMO-1 throughout the training phase (one epoch is one pass
through the entire training set). For each quality assessment, the upper left panel corresponds to the behavior of Wind-Topo on the training
set, while the right panels correspond to the test stations and the lower panels to the test period. Scores for wind speed (vel) and direction
(dir) are computed from u and v and are depicted in red. The lower right figure shows, with various colors, the distributions of wind speed for
the three station categories. Various lines styles are used to depict those distributions from measurement, COSMO-1, and Wind-Topo [Colour
figure can be viewed at wileyonlinelibrary.com]

campaign of surface wind is performed somewhere, the 3.2 Disaggregated scores


collected data can be used to train a model like Wind-Topo,
which can generate longer “synthetic” time series cover- We can assess Wind-Topo further by looking at the scores
ing the same period as the NWP model. The lower right at each station. Figure 7 shows such a disaggregated
panel of Figure 6 shows the distributions of measured view for the MAE. Appendix S1 provides similar plots
and predicted wind speeds for the three station types. for correlation, MBE, nMAE, and MAEdir. Each segment
We observe again the lack of discrimination of COSMO-1 represents a station, with its extremities being the scores
between exposed and sheltered sites. It is, however, clear on training and test periods. The test stations are circled
in the measurement, and Wind-Topo can reproduce it in black and are the ones we should focus on. Any point
well. The distributions of predicted wind speed match the located below the black curve indicates that Wind-Topo
observed ones accurately for all station types and for the has a lower MAE than COSMO-1. This is the case for all
four datasets. exposed and sheltered stations. Wind-Topo is, however,
The “early stopping” procedure presented above, in slightly worse at some “other” stations and significantly
association with the chosen learning rate, provides the best worse at one. The clustering of the station categories is
scores and distributions. Other regularization techniques remarkable. Wind-Topo strongly reduces the MAE at the
like weight decay did not show any supplementary advan- sheltered sites with high COSMO-1 MAE. This is also the
tage. On a RTX2080Ti GPU, the 12 epochs of Figure 6 case at exposed sites, albeit less pronounced. The reduced
required 43 hr. performance at the “other” stations is due to the already
1477870x, 2022, 744, Downloaded from https://rmets.onlinelibrary.wiley.com/doi/10.1002/qj.4265 by Readcube (Labtiva Inc.), Wiley Online Library on [13/03/2023]. See the Terms and Conditions (https://onlinelibrary.wiley.com/terms-and-conditions) on Wiley Online Library for rules of use; OA articles are governed by the applicable Creative Commons License
DUJARDIN and LEHNING 1379

FIGURE 7 Mean absolute error of Wind-Topo at all 321 stations compared with COSMO-1. Each segment corresponds to one station.
The color indicates the category of the station. The extremities of each segment show the MAE for the training period (triangles) and the test
period (circles). The 60 test stations are circled in black. Any point below the black line shows an improvement from Wind-Topo [Colour
figure can be viewed at wileyonlinelibrary.com]

good performance of COSMO-1. The orientation of the reflect appropriate corrections. However, it seems that the
segments indicates how the MAE in COSMO-1 changes learned corrections of COSMO-1’s underestimation in this
from one period to another and how Wind-Topo reacts to type of topography are applied to our nearby “outlier”.
it. An orientation parallel to y = x reflects that Wind-Topo Such topographies are not common in our dataset and
follows the variation of score in COSMO-1. An upward Wind-Topo would certainly benefit from more training
deviation from this orientation (triangle to circle) shows stations, located in more diverse terrain.
that Wind-Topo suffers from a lack of temporal general- Figure 8 provides a temporal disaggregation of the
ization and becomes (relatively) worse for the test period. scores for each month of the test period. The 0.25, 0.5,
This is almost never the case and the opposite behavior and 0.75 quantiles, computed for the 60 test stations, are
dominates. represented in bright colors for Wind-Topo and faint col-
The “outlier” mentioned above is an interesting test ors for COSMO-1. The upper panel shows that the win-
station located on the smooth ridge of a hill (700 m ter period is the windiest and that COSMO-1 predicts
above the plains). Its topography is well represented in similar wind speeds at exposed and sheltered sites. This
COSMO-1, which performs well with almost no MBE large positive bias (MBE panel) is always corrected by
(0.2 m⋅s−1 ) and a low MAE (1.2 m⋅s−1 ). Wind-Topo over- Wind-Topo. The negative bias occurring at exposed sites
estimates the wind speed there (MBE of 1.2 m⋅s−1 ), which in winter is also corrected, as well as the slight posi-
leads to a higher MAE (1.5 m⋅s−1 ). Surprisingly, at another tive bias at the “other” stations in summer. The panel
(training) station located 11 km from there and 450 m on correlation shows the significant improvement from
higher on the ridge, COSMO-1 behaves very differently. It Wind-Topo at exposed sites, especially in early sum-
underestimates the wind speed (4.6 m⋅s−1 ) compared with mer, when thermally driven flows become active. At the
measurement (7.7 m⋅s−1 ). This large bias is almost entirely “other” stations, Wind-Topo follows closely, but slightly
corrected by Wind-Topo (MBE of −0.9 m⋅s−1 ). Being in surpasses, COSMO-1. The same happens at sheltered
the training set, this performance does not necessarily sites, except in July–September, when the correlation is
1477870x, 2022, 744, Downloaded from https://rmets.onlinelibrary.wiley.com/doi/10.1002/qj.4265 by Readcube (Labtiva Inc.), Wiley Online Library on [13/03/2023]. See the Terms and Conditions (https://onlinelibrary.wiley.com/terms-and-conditions) on Wiley Online Library for rules of use; OA articles are governed by the applicable Creative Commons License
1380 DUJARDIN and LEHNING

F I G U R E 8 Each panel shows monthly average values for the three types of test station (colors). Those values are the 0.25 quantiles
(lower extremities of the vertical segments), 0.5 quantiles (marks within the segments), and 0.75 quantiles (upper extremities) of the monthly
values for each station. Values for Wind-Topo are depicted with bright colors, while values for COSMO-1 are fainter [Colour figure can be
viewed at wileyonlinelibrary.com]

slightly decreased in Wind-Topo, despite a low MAE. wind speed and direction are systematically lower with
We should note the low average wind speed for this Wind-Topo. This is true for the 0.5 quantile and almost
period (1.2 m⋅s−1 ) and the corresponding large MAEdir always true for the other quantiles, showing that the
from COSMO-1 (68◦ ). Wind-Topo can lower it to 60◦ , downscaling performs well at almost all stations and for all
but at a cost of a lower correlation. Finally, the MAE of seasons.
1477870x, 2022, 744, Downloaded from https://rmets.onlinelibrary.wiley.com/doi/10.1002/qj.4265 by Readcube (Labtiva Inc.), Wiley Online Library on [13/03/2023]. See the Terms and Conditions (https://onlinelibrary.wiley.com/terms-and-conditions) on Wiley Online Library for rules of use; OA articles are governed by the applicable Creative Commons License
DUJARDIN and LEHNING 1381

3.3 Selected stations situation is reversed, with the COSMO-1 direction being
accurate while Wind-Topo rotates it slightly. This can be
For a last quantitative analysis, we selected three test sta- explained by the peculiar location of this station (inter-
tions in the test area of Figure 1: ELA-1, ELA-2, and section of several valleys) and it reflects a weakness of
SAM-0. Indices 1, 2, 0 stand for exposed, sheltered, and Wind-Topo. At ELA-2, Wind-Topo successfully corrects
other, respectively. We can first look at the wind-rose wind speed and direction. These three stations illustrate
diagrams for those stations (Figure 9), generated with the range of potential corrections depending on station
data from measurements, COSMO-1, and Wind-Topo. For characteristics. At exposed sites like ELA-1, COSMO-1 per-
information, their locations are depicted in Figure 10. forms reasonably well but Wind-Topo can still improve the
We can observe that the wind direction in Wind-Topo is flow. At sheltered locations like ELA-2, the coarse reso-
strongly influenced by COSMO-1. This is especially the lution of COSMO-1 cannot reflect the sheltering and this
case at ELA-1, where Wind-Topo can only correct the leads to a strong overestimation of wind speed and wrong
westerly winds from COSMO-1 partially. At SAM-0, the wind direction. Wind-Topo shows a clear ability to lower

F I G U R E 9 Wind-rose diagrams for the three selected test stations during the test period from measurement, COSMO-1, and
Wind-Topo. The color indicates wind speed, while the distance from the center indicates the probability in percent that the wind comes from
one of 16 directions. To ease the comparison, the radial axis has the same extent for the three data sources [Colour figure can be viewed at
wileyonlinelibrary.com]
1477870x, 2022, 744, Downloaded from https://rmets.onlinelibrary.wiley.com/doi/10.1002/qj.4265 by Readcube (Labtiva Inc.), Wiley Online Library on [13/03/2023]. See the Terms and Conditions (https://onlinelibrary.wiley.com/terms-and-conditions) on Wiley Online Library for rules of use; OA articles are governed by the applicable Creative Commons License
1382 DUJARDIN and LEHNING

F I G U R E 10 Wind polar plots for three selected test stations. The upper two rows show the topography used by Wind-Topo (top: full
domain, bottom: central crop). The lower rows are all constructed in the same way: polar plots of a desired quantity by wind speed and
direction. For example, strong northwesterly winds are located in the upper left corner. Events with wind speed smaller than 1 m⋅s−1 are not
depicted. The colors indicate the average value for all events with a given wind speed and direction. Row 3 depicts the probability density, on
which quantile–quantile plots for Wind-Topo and COSMO-1 (of wind speed with respect to measurements) are superimposed. Rows 4 and 5,
respectively, depict the relative error in wind speed from Wind-Topo (blue: underestimation, red: overestimation) and how this error changed
between COSMO-1 and Wind-Topo (green being an improvement). The lower two rows are equivalent for wind direction (azimuth) [Colour
figure can be viewed at wileyonlinelibrary.com]
1477870x, 2022, 744, Downloaded from https://rmets.onlinelibrary.wiley.com/doi/10.1002/qj.4265 by Readcube (Labtiva Inc.), Wiley Online Library on [13/03/2023]. See the Terms and Conditions (https://onlinelibrary.wiley.com/terms-and-conditions) on Wiley Online Library for rules of use; OA articles are governed by the applicable Creative Commons License
DUJARDIN and LEHNING 1383

the wind speed to match observations better and to correct least some of those locally generated flows in the near
the wind direction. Finally, at locations where the topog- future.
raphy is already well described at COSMO-1 resolution, To conclude our analysis for the three stations selected,
Wind-Topo is not able to add relevant information and Figure 11 details the predictions of the two models
might even decrease the quality of the prediction slightly. throughout the test year and for each hour of the day. The
In Figure 10 we present a new type of plot that reveals histograms of wind speed for the whole year reflect the
how a model like COSMO-1 or Wind-Topo performs (with typical behavior of COSMO-1 at two nearby exposed and
respect to observations) across diverse wind conditions. sheltered stations (ELA-1, ELA-2): almost identical, which
The first two rows depict the local topography of the three does not correspond to measurements. As seen previously,
selected stations (full and zoom patches of input data), this is corrected by Wind-Topo. The daily patterns of wind
while the lower rows, coined “wind polar plots”, are polar speed are quite different from site to site and for differ-
plots by wind speed and direction and can be seen as a vari- ent times of the year. SAM-0 has typical, thermally driven
ation of wind-rose diagrams. The color of each point repre- flows in the mid-afternoon, which are stronger in summer.
sents a certain quantity like the bias in wind speed, which COSMO-1 captures them very well and Wind-Topo repli-
is the average value for all events with such wind speed cates them. The measurements at ELA-1 show peculiar
and direction. For example, strong northwesterly wind wind patterns: wind is steadier and stronger in winter, and
events are located in the upper left corner. Additionally, we is significantly reduced in the middle of the day in summer
superimposed on the third row, which shows the probabil- due to complex interactions between the boundary layer
ity density, red curves corresponding to quantile–quantile and the free atmosphere. COSMO-1 predicts the opposite,
(Q–Q) plots for Wind-Topo and COSMO-1 (of wind speed, with an increase of wind speed in the middle of the day in
with respect to measurement). For ELA-1, the Q–Q plot summer. It is remarkable that Wind-Topo can reproduce
shows a small underestimation of COSMO-1 for all but the observed patterns and does not simply copy COSMO-1
high wind speed. Wind-Topo provides an almost perfect at this exposed station, which is located on a type of terrain
correction, except for a small overestimation at low wind where COSMO-1 normally performs quite well.
speed. This quantitative analysis displayed many strengths of
The fourth row confirms the small bias at low wind Wind-Topo, which can (a) distinguish accurately between
speed, but no clear dependence on wind direction. The exposed and sheltered locations and reduce the biases
row below depicts the improvement (green) or worsen- resulting from a simple interpolation of COSMO-1 to
ing (red) of Wind-Topo compared with COSMO-1. Except a higher resolution, (b) generate better wind-speed dis-
for wind speeds below 2 m⋅s−1 , the improvement is consis- tributions in complex terrain, and (c) capture complex
tent across all wind speeds and direction. The bottom two small-scale wind–topography interactions like ridge accel-
rows offer a similar analysis for the error in wind direc- eration, sheltering, and deflection. The model is not per-
tion. Wind-Topo has a small overestimation of the azimuth fect, as it fails to introduce subgrid-scale thermally driven
when the wind is southwesterly to northeasterly, and an flows and might introduce errors in locations where
underestimation otherwise. Compared with COSMO-1, COSMO-1 performs well. Given the current level of per-
the improvement is mixed and presents the same split as formance, based on only 261 training stations, we are
just mentioned. For brevity, the reader is encouraged to confident that more diverse training data (ground truth at
observe the wind polar plots for other stations with various other locations) would benefit Wind-Topo strongly.
topographies. Appendix S1 provides the plots for all 60 val-
idation stations. Looking simultaneously at the probability
density plots and the plots for other metrics, as well as 3.4 Generated wind fields
topography, helps us understand how COSMO-1 performs
in the most frequent wind situations and how Wind-Topo We can finally look at some examples of generated wind
reacts. Wind-Topo seems to perform well across a variety of fields for the test area in Figure 1. COSMO-1 was down-
terrain and wind speed and direction, but it cannot, how- scaled for all hours of the test year to a 100-m grid cov-
ever, introduce small-scale thermally driven flows that are ering this 30 × 30 km2 domain, with a bicubic interpo-
not present in COSMO-1. At some sheltered sites, located lation and with Wind-Topo. The latter required 5.5 days
on slopes exposed to the sun, Wind-Topo fails to recreate on a RTX2080Ti GPU. Figure 12a shows an example of
the observed diurnal slope flows when they are not present an interpolated COSMO-1 wind field and the underlying
in COSMO-1. A decomposition of the wind polar plots for (interpolated) model terrain. Figure 12b is the 100-m DEM,
day and night events in summer (not provided) reveals it and Figure 12c shows the prediction from Wind-Topo
even more clearly. Ongoing effort is addressing this weak- with wind speed enhanced on ridges and exposed slopes
ness and we expect Wind-Topo to be able to introduce at and reduced on sheltered slopes. It is not a pale copy of
1477870x, 2022, 744, Downloaded from https://rmets.onlinelibrary.wiley.com/doi/10.1002/qj.4265 by Readcube (Labtiva Inc.), Wiley Online Library on [13/03/2023]. See the Terms and Conditions (https://onlinelibrary.wiley.com/terms-and-conditions) on Wiley Online Library for rules of use; OA articles are governed by the applicable Creative Commons License
1384 DUJARDIN and LEHNING

F I G U R E 11 Upper panels: distributions of wind speed for the test period from measurement (black), COSMO-1 (red), and Wind-Topo
(green) for the three test stations selected. With the same colors, the lower panels show, for each pair of months (in the test period), the
average wind speed (y-axis) for each hour of the day (x-axis) [Colour figure can be viewed at wileyonlinelibrary.com]

COSMO-1 with a simple scaling by elevation or exposi- and Figure 12f depicts the distributions at every grid point.
tion. A careful observation (full-size panels available in A Weibull distribution is fitted to the histogram of wind
Appendix S1) indicates the presence of flow deflection and speed at each grid point, and the corresponding parame-
recirculation areas on some lee sides. COSMO-1 shows ters are used to determine the color. The scale parameter
such effects on a much larger scale. Wind-Topo can even (resembling the average wind speed) defines the bright-
remove them in some areas, before introducing them in ness, while the shape parameter defines the hue. Blue
others. We should also note the continuity in wind speed areas have a shape parameter close to 1, meaning broad
in the generated field, despite the point-by-point method. distributions, while red areas show narrower distributions,
This is caused by the large overlap of input data between which are more centered around the average. Here, we also
neighboring points. To visualize the flow deflection and observe a good spatial continuity in the fitted distributions.
continuity of wind direction better, Figure 12d depicts the The large valleys all have a low scale parameter, which
wind azimuth with a cyclic color scale and associated increases towards the summits and ridges. However, in
streamlines. many locations the highest values are not found there,
Figure 12e,f provides an insight into the entire down- but slightly below them. Also, the concave areas on slopes
scaled dataset. Figure 12e shows the average wind speed exhibit higher scale parameters than their surroundings.
1477870x, 2022, 744, Downloaded from https://rmets.onlinelibrary.wiley.com/doi/10.1002/qj.4265 by Readcube (Labtiva Inc.), Wiley Online Library on [13/03/2023]. See the Terms and Conditions (https://onlinelibrary.wiley.com/terms-and-conditions) on Wiley Online Library for rules of use; OA articles are governed by the applicable Creative Commons License
DUJARDIN and LEHNING 1385

F I G U R E 12 (a) Wind field from COSMO-1, interpolated (bicubic) on a 100-m grid, for a specific time. The background color indicates
the wind speed and red arrows indicate the wind direction every 500 m. Grey lines show the 100-m contours of the associated topography. (b)
100-m digital elevation model used by Wind-Topo covering the same domain. (c) 100-m wind field generated by Wind-Topo for the same
time. Grey lines show the underlying topography (100-m contours). (d) Associated wind direction (clockwise from North) with a cyclic color
scale. The green lines are the stream lines of the wind field, with arrows indicating the flow direction. (e) Average wind speed from
Wind-Topo over the test period. (f) Weibull parameters of the distribution fitted in each pixel. The brightness indicates the scale parameter
and the color indicates the shape parameter. 100-m contours in grey show the underlying topography [Colour figure can be viewed at
wileyonlinelibrary.com]
1477870x, 2022, 744, Downloaded from https://rmets.onlinelibrary.wiley.com/doi/10.1002/qj.4265 by Readcube (Labtiva Inc.), Wiley Online Library on [13/03/2023]. See the Terms and Conditions (https://onlinelibrary.wiley.com/terms-and-conditions) on Wiley Online Library for rules of use; OA articles are governed by the applicable Creative Commons License
1386 DUJARDIN and LEHNING

A similar plot for COSMO-1 is provided in Appendix S1 are already foreseen, in particular if more diverse loca-
and shows the same increase in concave areas on larger tions with observed or modeled “ground truth” can be used
scales. This local increase can be explained by the chan- for training. Given the rapid evolution of the research on
neling that occurs in such places, which generates many deep learning, we believe that its architecture will also be
situations with relatively high wind speed but few occur- improved in the near future.
rences of a very large one. Wind-Topo is able to reproduce The point-by-point nature of Wind-Topo can be seen
such an effect on smaller scales. as a strength or as a weakness. Using point observations
One might wonder about the physical consistency of as predictands avoids the introduction of inaccuracies or
the wind fields generated. As Wind-Topo only predicts a biases that other (modeled) data sources would have. The
near-surface wind field, we cannot compute a 3D diver- pointwise prediction also simplifies the parallelization of
gence to assess the mass conservation. We can neverthe- the downscaling scheme. However, it prevents the intro-
less estimate the ground-perpendicular motion caused by duction of physical constraints during the training and pre-
the divergence of the 2D field. Appendix S1 describes dictions that would guarantee the generation of physically
this analysis along with figures, which we summarize consistent flows. Ideally, Wind-Topo should be developed
here. The large majority of such motions at 7 m a.g.l. further to incorporate high-resolution 2D or 3D modeled
are between −0.04 and 0.04 m⋅s−1 (negative for motions data as well, which could help identify wind–topography
towards the ground). Rarely, those motions can be 10 times interactions not yet discovered by the model. This new data
larger. When averaged for the entire year, we obtain values source could also provide the physical constraints needed
between −0.06 and 0.05 m⋅s−1 . Furthermore, if we apply a when generating wind fields.
spatial averaging to coarse-grain these maps to the scales Wind-Topo is already quite fast and can be opti-
resolved by COSMO-1, we obtain patterns that are very mized further. In our Swiss case study, it can downscale
similar to those of COSMO-1 in terms of magnitude (−0.01 COSMO-1 (1,110 m) to 6 million grid points in one hour
to 0.006 m⋅s−1 ) and location. on one standard GPU. It is thus possible to downscale the
Swiss Alps (about 29,000 km2 ) operationally (hourly) to a
4 CO N C LU S I O N 50-m grid with a 2-GPU machine. As the computational
requirements scale linearly with the number of grid points
A high spatial resolution in complex terrain is essential desired, Wind-Topo can be applied to much larger domains
to capture small-scale processes shaping the wind flow at this resolution with reasonable resources. In Appendix
that is responsible for phenomena like snow preferential S1, we compare the current version of the model with sim-
deposition or transport. Even state-of-the-art operational pler models and lighter model variants. The 2D nature of
models like COSMO-1 running on supercomputers can- the model provided by its CNNs is of high importance. The
not provide the required resolution for topographies like custom loss function is also crucial in obtaining correct
the Alps. Wind-Topo offers the possibility to extract the distributions of wind speed.
most information from such models and to correct the We hope this work will pave the way for new meth-
near-surface wind fields to match observations closely. ods based on deep learning that will downscale not only
The systematic biases observed in different types of ter- the wind but also other atmospheric variables to higher
rain and at different times of the year and day are strongly resolutions in complex terrain. Such methods will be able
reduced. The downscaled wind fields exhibit some of the to combine state-of-the-art high-resolution NWP models,
expected orographic effects (ridge acceleration, sheltering, observations, and operational models efficiently.
and deflection), and Wind-Topo does not alter the predic-
AU THOR CONTRIBUTIONS
tions from the low-resolution input model in areas where
Jérôme Dujardin: conceptualization; investigation; data
it already performs well. At 44 measurement stations
collection; curating; development; validation; writing -
located in diverse topographies in the Alps and put aside
original draft; writing - review and editing. Michael
for validation, Wind-Topo offers an average correlation of
Lehning: conceptualization; funding acquisition;
0.55, a low bias of −0.07 m⋅s−1 , and a MAE of 1.21 m⋅s−1 .
investigation; supervision; writing – original draft;
This should be compared with the values obtained by a
writing – review and editing.
simple interpolation of COSMO-1: correlation of 0.50, bias
of 0.72 m⋅s−1 , and MAE of 1.77 m⋅s−1 . ACKNOWLEDGEMENTS
With its novel architecture, Wind-Topo not only down- We thank Mathieu Schaer for the first investigations
scales the wind speed to higher resolutions but also incor- on Wind-Topo during his master’s project and Quentin
porates local orographic deflections. Even if subgrid-scale Fisch for his help in parallelizing the code. This work
thermally driven flows are not captured currently, the level was funded by Innosuisse through the Swiss Competence
of performance is remarkable and further improvements Center for Energy Research, Supply of Electricity. Open
1477870x, 2022, 744, Downloaded from https://rmets.onlinelibrary.wiley.com/doi/10.1002/qj.4265 by Readcube (Labtiva Inc.), Wiley Online Library on [13/03/2023]. See the Terms and Conditions (https://onlinelibrary.wiley.com/terms-and-conditions) on Wiley Online Library for rules of use; OA articles are governed by the applicable Creative Commons License
DUJARDIN and LEHNING 1387

Access Funding provided by Ecole Polytechnique Federale Boundary-Layer Meteorology, 179, 133–161. https://doi.org/10.
de Lausanne. [Correction added on 23 May 2022, after 1007/s10546-020-00586-x.
first online publication: CSAL funding statement has been Helbig, N., Mott, R., Van Herwijnen, A., Winstral, A. and Jonas, T.
(2017) Parameterizing surface wind speed over complex topog-
added.]
raphy. Journal of Geophysical Research, 122, 651–667. https://doi.
org/10.1002/2016JD025593.
DATA AVAILABILITY STATEMENT Hendrycks, D. and Gimpel, K. (2016) Gaussian Error Linear Units
The original datasets used to train and validate Wind-Topo (GELUs). https://doi.org/10.48550/arXiv.1606.08415.
can be obtained free of charge from the data own- Hilton, J. and Garg, N. (2021) Rapid wind-terrain correction for
ers MeteoSwiss and the WSL Institute for Snow and wildfire simulations. International Journal of Wildland Fire, 30,
Avalanche Research. Intermediate data products can be 410–427.
obtained from the authors upon request. A complete Höhlein, K., Kern, M., Hewson, T. and Westermann, R. (2020) A
comparative study of convolutional neural network models
use-case dataset is available from the ENVIDAT data plat-
for wind field downscaling. Meteorological Applications, 27, 1–31.
form under the DOI: 10.16904/envidat.301. It includes Hsieh, W.W. (2009) Machine Learning Methods in the Environ-
the architecture of Wind-Topo and its optimized param- mental Sciences: Neural Network and Kernels. Cambridge: Cam-
eters, as well as a python script to downscale uniform bridge University Press. https://doi.org/10.1017/CBO9780511
wind fields with a prescribed vertical profile for any given 627217
53-m DEM. New versions of the code will be available on Huang, H.Y., Capps, S.B., Huang, S.C. and Hall, A. (2015)
the WSL GitLab: https://gitlabext.wsl.ch/dujardin/wind- Downscaling near-surface wind over complex terrain using a
physically-based statistical modeling approach. Climate Dynam-
topo.
ics, 44, 529–542.
Kingma, D.P. and Ba, J.L. (2015). ADAM: a method for stochastic
ORCID optimization. 3rd International Conference on Learning Repre-
Jérôme Dujardin https://orcid.org/0000-0001-5404- sentations, ICLR 2015, pp. 1–15. https://arxiv.org/pdf/1412.6980.
7734 pdf
Kruyt, B., Dujardin, J. and Lehning, M. (2018) Improvement of Wind
REFERENCES Power Assessment in Complex Terrain: The Case of COSMO-1
Burger, W. and Burge, M.J. (2009) Principles of Digital Image Process- in the Swiss Alps. Frontiers in Energy Research, 6, 102. https://10.
ing: Core Algorithms. London: Springer. 3389/fenrg.2018.00102.
Dupuy, F., Duine, G.J., Durand, P., Hedde, T., Pardyjak, E. and Lehning, M., Bartelt, P., Brown, B., Russi, T., Stöckli, U. and
Roubin, P. (2021) Valley winds at the local scale: correcting Zimmerli, M. (1999) SNOWPACK model calculations for
routine weather forecast using artificial neural networks. Atmo- avalanche warning based upon a new network of weather
sphere, 12, 28. and snow stations. Cold Regions Science and Technology, 30,
Etienne, C., Lehmann, A., Goyette, S., Lopez-Moreno, J.-I. and 145–157.
Beniston, M. (2010) Spatial predictions of extreme wind speeds Manor, A. and Berkovic, S. (2015) Bayesian inference aided ana-
over Switzerland using generalized additive models. Journal of log downscaling for near-surface winds in complex terrain.
Applied Meteorology and Climatology, 49, 1956–1970. Atmospheric Research, 164-165, 27–36. https://doi.org/10.1016/j.
Fiddes, J. and Gruber, S. (2014) TopoSCALE v.1.0: downscaling grid- atmosres.2015.04.014.
ded climate data in complex terrain. Geoscientific Model Develop- Mott, R., Schirmer, M., Bavay, M., Grünewald, T. and Lehning,
ment, 7, 387–405. M. (2010) Understanding snow-transport processes shaping the
Fischer, P., Etienne, C., Tian, J. and Krauß, T. (2015) Predic- mountain snow-cover. Cryosphere, 4, 545–559.
tion of wind speeds based on digital elevation models using Paszke, A., Gross, S., Chintala, S., Chanan, G., Yang, E., DeVito, Z.,
boosted regression trees. ISPRS – International Archives of the Lin, Z., Desmaison, A., Antiga, L. and Lerer, A. (2017). Automatic
Photogrammetry, Remote Sensing and Spatial Information differentiation in PyTorch. 31st Conference on Neural Information
Sciences, XL-1-W5, 197–202. https://elib.dlr.de/99861/1/ Processing Systems (NIPS 2017). https://openreview.net/pdf?id=
FischerEtAl.pdf BJJsrmfCZ
Foresti, L., Tuia, D., Kanevski, M. and Pozdnoukhov, A. (2011) Learn- Philippopoulos, K. and Deligiorgi, D. (2012) Application of artificial
ing wind fields with multiple kernels. Stochastic Environmental neural networks for the spatial estimation of wind speed in a
Research and Risk Assessment, 25, 51–66. coastal region with complex topography. Renewable Energy, 38,
Glorot, X. and Bengio, Y. (2010). Understanding the difficulty of 75–82. https://doi.org/10.1016/j.renene.2011.07.007.
training deep feedforward neural networks. Proceedings of the Reichstein, M., Camps-Valls, G., Stevens, B., Jung, M.,
Thirteenth International Conference on Artificial Intelligence and Denzler, J., Carvalhais, N. and Prabhat (2019) Deep learning
Statistics Vol. 9, pp. 249–256. https://proceedings.mlr.press/v9/ and process understanding for data-driven Earth system sci-
glorot10a/glorot10a.pdf ence. Nature, 566, 195–204. https://doi.org/10.1038/s41586-019-
Goutham, N., Alonzo, B., Dupré, A., Plougonven, R., Doctors, R., 0912-1.
Liao, L., Mougeot, M., Fischer, A. and Drobinski, P. (2021) Rios, O., Jahn, W., Pastor, E., Valero, M.M. and Planas, E. (2018)
Using machine-learning methods to improve surface wind speed Interpolation framework to speed up near-surface wind simula-
from the outputs of a numerical weather prediction model. tions for data-driven wildfire applications. International Journal
of Wildland Fire, 27, 257–270.
1477870x, 2022, 744, Downloaded from https://rmets.onlinelibrary.wiley.com/doi/10.1002/qj.4265 by Readcube (Labtiva Inc.), Wiley Online Library on [13/03/2023]. See the Terms and Conditions (https://onlinelibrary.wiley.com/terms-and-conditions) on Wiley Online Library for rules of use; OA articles are governed by the applicable Creative Commons License
1388 DUJARDIN and LEHNING

Robert, S., Foresti, L. and Kanevski, M. (2013) Spatial prediction of Winstral, A., Jonas, T. and Helbig, N. (2017) Statistical downscal-
monthly wind speeds in complex terrain with adaptive general ing of gridded wind speed data using local topography. Journal of
regression neural networks. International Journal of Climatology, Hydrometeorology, 18, 335–348. https://doi.org/10.1175/JHM-D-
33, 1793–1804. https://doi.org/10.1002/joc.3550. 16-0054.1.
Takeuchi, I., Le, V.Q., Sears, T. and Smola, A. (2006) Non-
parametric Quantile Estimation. Journal of Machine Learn- SUPPORTING INFORMATION
ing Research, 7, 1231–1264. https://jmlr.org/papers/volume7/ Additional supporting information may be found online
takeuchi06a/takeuchi06a.pdf. in the Supporting Information section at the end of this
Vionnet, V., Marsh, C.B., Menounos, B., Gascoin, S., Wayand, N.E., article.
Shea, J., Mukherjee, K. and Pomeroy, J.W. (2021) Multi-scale
snowdrift-permitting modelling of mountain snowpack.
Cryosphere, 15, 743–769. How to cite this article: Dujardin, J. & Lehning,
Wagenbrenner, N.S., Forthofer, J.M., Lamb, B.K., Shannon, K.S. M. (2022) Wind-Topo: Downscaling near-surface
and Butler, B.W. (2016) Downscaling surface wind predictions wind fields to high-resolution topography in highly
from numerical weather prediction models in complex ter-
complex terrain with deep learning. Quarterly
rain with WindNinja. Atmospheric Chemistry and Physics, 16,
Journal of the Royal Meteorological Society, 148(744),
5229–5241.
Winstral, A., Elder, K. and Davis, R.E. (2002) Spatial snow model- 1368–1388. Available from: https://doi.org/10.1002/
ing of wind-redistributed snow using terrain-based parameters. qj.4265
Journal of Hydrometeorology, 3, 524–538.

You might also like