Localising and quantifying night-time cooling effects from sub-catchments in a mid-European low mountain area

We propose a novel method to analyse and quantify cooling effects provided to a settlement by source areas of nocturnal cold-air drainage. In an interdisciplinary approach, these source areas were defined as hydrological sub-catchments of a complex catchment area in a low mountain range. The cold-air drainage model KLAM_21 was used to exclude the energetic influence of the sub-catchments from the model area by surrounding them with artificial barriers. The outputs of these runs were then compared to a reference run without exclusion to derive the cooling effect of each source area. The results were evaluated at sample points along the main valley and for residential areas of a medium sized city and two smaller settlements. We find that in the complex terrain of the study area, also comparatively remote source areas can have a noticeable cooling effect on the residential neighbourhoods of the target settlements from the middle of the night. The strongest effect however, could be attributed to the sub-catchments in direct vicinity of the target areas. The results at the sample points along the main valley showed that the cooling effect decreased with increasing distance to the sub-catchments and usually gets stronger during the night. The variation in strength of cooling effect between different sub-catchments is likely due to their individual properties such as remoteness, size, terrain, land-cover situation and cold-air exchange with other sub-catchments through overflow effects.


Introduction
With ongoing climate change, the intensity and length of heatwaves are expected to increase in European cities (Guerreiro et al. 2018). There is also evidence of negative synergies between these heat waves and the urban heat island (UHI), with the UHI being more intense during heat waves than during normal periods He et al. 2021;Li and Bou-Zeid 2013). This phenomenon is likely to intensify in case of the night-time UHI of temperate climate cities in the future (Zhao et al. 2018). As the largest temperature differences between a city and the surrounding rural area can be observed at night (Oke 1982), not only the daily temperature extremes but also the lack of nocturnal cooling periods are a major problem for the urban population (Wicki et al. 2018). In cities located in lower parts of river catchments, night-time cold-air drainage flow is known to mitigate the UHI effect (Kuttler et al. 1996;Montvez et al. 2000). Since numerous studies have shown that the UHI increases health risks from heat exposure (Heaviside et al. 2017), cold-air drainage can have a positive impact on the well-being of the urban population.
Cold-air drainage and its effect on urban climate have been a subject of research since the middle of the twentieth century. For example, Balchin and Pye (1947) detected coldair drainage affecting the nocturnal air temperature in the city of Bath, England. Nkemdirim (1980) studied its effect in Calgary, Canada, with the goal to isolate topographical factors from urban climate studies. Here, the problem of reduced turbulence and its consequences for air pollution in the city was addressed but the ecosystem service of UHI mitigation provided by cold-air drainage was not yet of interest. The impact of cold-air drainage on air pollution, however, still remains a research topic, especially in highly polluted areas. It can have both negative and positive effects, depending on the topography of the study area and the sources of air pollution. In Stuttgart for example, an industrial German city located in a narrow basin, city planners intent to preserve open corridors for cold-air drainage from the surrounding slopes in order to ensure ventilation of the urban centre and disposal of pollutants (Reuter and Kapp 2021). The topography of alpine valleys on the other hand can lead to the formation of cold air pools in winter, in which pollutants are trapped (Largeron and Staquet 2016;Quimbayo-Duarte et al. 2021). Kossmann and Sturman (2004) observed a negative effect in Christchurch caused by two converging drainage winds but a positive disposal of pollutants later at night, when the stronger drainage dominated. Impacts of forests in the study area may also play a role for pollutant distribution. Forests are known to remove gaseous pollutants but also reduce wind speed and therefore hinder pollutant dispersion (Nowak 2020). Within the last decades, the reduction of air pollution has been put on the political agenda, for example with the UN sulphur protocols of 1985 (Helsinki) and 1994 (Oslo). While the effectiveness of these protocols cannot be confirmed with significant certainty, SO 2 was especially decreasing in signatory countries (Aakvik and Tjøtta 2011;Barrett et al. 1995). There is also evidence that SO 2 decreased on a city level in temperate climates during at least the last decades (Brimblecombe and Lefèvre 2021;Grøntoft 2021).
In recent years, the positive cooling effect of cold air drainage became a major motivation of studying it in urban environments. Bigg et al. (2014) found that cold-air drainage flows can mitigate the UHI in parts of Sheffield, Son et al. (2022) analysed the cold-air situation in three Korean cities in order to ensure appropriate wind corridor planning, and the effects of future urban land use change on cold-air paths was analysed by Grunwald and Weber (2021). Cold-air flow was also studied in a field experiment by King (1973) where small catchments were blocked by huge curtains and differences to the natural situation were measured. In recent studies, cold-air flow is often calculated with numerical models. One example is KLAM_21, a two-dimensional model developed by the German Weather Service (Deutscher Wetterdienst; DWD) (Sievers and Kossmann 2016). It was used by Grunwald et al. (2019), for instance, to classify cold-air paths in terms of their efficiency for UHI mitigation strategies. This also involved theoretical considerations on the relationship between the city centre and the origin of the cold-air flows as the authors defined cold-air paths as a connection between cold-air reservoir areas and cold-air impact regions. A different numerical model was used to investigate air pollution associated with the formation of cold-air pools in urbanized alpine valleys during wintertime (Quimbayo-Duarte et al. 2019). There, version 3 of the Weather Research and Forecasting (WRF) model (Skamarock et al. 2021) was used for a numerical simulation of the atmospheric conditions in an idealized valley to analyse the transport of passive tracers. Furthermore, the authors discussed theoretical considerations on the influence of orographic variations such as valley width on flow behaviour. Besides model-based approaches, a large number of studies investigated cold-air phenomena with measurements in the field. One example was the Shallow Cold Pool Experiment which was conducted in a small valley in Colorado during October and November 2012 (Mahrt et al. 2014).
Since the UHI effect in combination with global warming is also associated with a higher energy demand (Santamouris et al. 2015), cold-air drainage cooling could also be beneficial to the mitigation of future increased energy demand in the summer. Several studies have shown that the intensity of cold-air drainage is linked to cloud coverage. Clear nights offer ideal conditions for radiative cooling and therefore promoting the formation of cold air, whereas a high degree of overcast weakens the process (Barr and Orgill 1989;Fries et al. 2009;Gudiksen et al. 1992;Novick et al. 2016). It can therefore be assumed for temperate climate cities that cold-air drainage is weaker and less common in winter when cloudy conditions dominate (Meerkötter 2004;Baray et al. 2019). This is supported by a study in Aachen, Germany, where a higher probability of cold-air drainage events in summer than in winter was found (Sachsen 2013). With this in mind, a negative effect of cold air on the energy demand in winter (heating) might be weaker than the positive energetic effect it can have during heat waves in summer (air-conditioning).
Although green spaces within the city play their part in urban ventilation and cold-air production (Eliasson and Upmanis 2000), typical major cold-air flows originate in the rural surroundings. They are therefore an ecosystem service mainly provided by rural catchment areas located upstream of an urban centre. Little research has been done to analyse structures of such cold-air drainage catchments more precisely. Schau-Noppel et al. (2020) calculated linear trajectories between target points and source regions of cold-air drainage in Mainz and Wiesbaden, Germany, in order to visualize the relevance of cold-air production areas for municipal planning processes. Similar to this, a recent study by Steigerwald et al. (2022) identified urban thermal hot spots and cold-air formation areas within a core study area of 12 km × 12 km around Aschaffenburg, Germany. Forward trajectories from formation areas and backward trajectories from hot spots where then used to determine the origin of cooling effects. In a model study, Sachsen et al. (2013) described and quantified the cooling effect of cold-air inflow into an urbanised catchment in Aachen, Germany, by transfluence from other catchments.
However, the question of how much cooling energy at a target point or target area in a settlement originates in a single upstream sub-catchment (hereafter referred to as SC) remains unanswered. It is also not clear how large the cooling effect of SCs are, which are not located in the immediate surrounding of the city but further upstream. The aim of this study is to investigate methods of describing the effective areas of origin and quantify the nocturnal cooling effects in a topographically influenced city on the scale of SC areas in its remote surroundings with a numerical model. We hypothesize that SCs close to a target area have a stronger influence but remote SCs can still play a role for the nocturnal cooling of a city, depending on its individual properties. For this, a river catchment in a typical low mountain range in Central Europe with a city situated at the valley mouth was chosen as the study area. The basic concept to answer the question of how much cooling energy in different settlements of the study area is attributable to different SCs was to divide a given hydrological catchment into SCs and exclude them individually from the cold-air drainage calculations of the numerical model in a predefined sequence, derived from a hydrological SC-hierarchy. Comparisons to a reference run without exclusion would then show the effect of each SC on the study area. The results were analysed for sample points along the main valley and two settlements. The calculations were performed with KLAM_21 which is a two-dimensional cold-air drainage model developed for hilly and mountainous terrain under fair weather conditions (Sievers and Kossmann 2016). A catchment-based analysis of cold-air drainage is possible as it generally follows the hydrological valley structure. It has the additional advantage of clearly defined areas of origin following a typically physical-geographic logic. There are numerous approaches to describe the structure of catchments and their associated streams in a codification system (Verdin and Verdin 1999). In this study, a combination of the German Gewässerkennzahl (waterbody index number, GKZ), first introduced in LAWA (1993), and the stream order after Shreve (1966) was used to sort the SCs prior to the determination of their cooling effects.

Study area
The area of investigation was a horizontal rectangle in the coordinate system ETRS 89 UTM 32 N, which spans between the coordinates 32 N 271,349 5,654,635 and 32 N 346,334 5,579,650 (75 km × 75 km) in the tri-border area Belgium-Netherlands-Germany. The main subject of this study was the upper catchment of the river Rur in the Eifel-Ardenne region in eastern Belgium and western Germany, stretching from the Hautes Fagnes plateau and the low mountain range of the northern Eifel in the southwest to the lowland region of the Jülich-Zülpicher Börde in the northeast (Fig. 1a). This catchment area has a clearly defined valley structure shaped by the River Rur and its tributaries. The stream order after Shreve (1966) based on the main streams of the SCs analysed in this study is given in the lower section of Fig. 1; additionally, the stream order after Strahler (1952). The catchment is roughly funnel shaped with a diameter of approx. 43 km in the southern more elevated areas where the drainage density is relatively high, especially within the catchment of the tributary Urft ( Fig. 1 SCs 6-12). It narrows to the north in stream direction and is approx. 7-9 km wide between Nideggen and Düren. The catchment is influenced by extensive water management with several water reservoirs. The most prominent being the Rur reservoir (SCs 5 and 13) and Urft reservoir (SC 12). The mean slope value of the analysed catchment is 6.2°. Especially high values of over 30° can be observed along the valley of the main river, the Urft reservoir and the right hand tributary Erkensruhr (SC 4), including some parts of natural rock face. The most northern part near Düren, in contrast, is characterized by mostly flat terrain. The landscape is characterized by a combination of forests, mainly at higher altitudes, agricultural areas, scattered settlements and the Rhenish lignite mining region with two large open pit mines in the northeast of the model area but outside the main study area (Fig. 1b). The most populous city directly influenced by the catchment is Düren on the northern edge of the Eifel with 93.660 inhabitants. The height difference between the spring area of the Rur at the Botrange summit (694 m a.s.l.) and the centre of Düren is approximately 592 m.
According to the Köppen-Geiger climate classification system, the study area has a warm-temperate, fully humid climate with warm summers (Cfb) (Kottek et al. 2006). Considerable differences in annual precipitation values can be observed within the area of the upper Rur catchment (Bogena et al. 2005). The western, most elevated areas, which are influenced by the predominant maritime airflow from the west have annual precipitation values of up to 1200 mm/a while a lee effect results in decreasing values in north-eastern direction and less than 700 mm/a around Düren (Bogena et al. 2005).
The study area was chosen for two main reasons. First, the city of Düren was assessed as a climate change precautionary area by the State Office for Nature, Environment and Consumer Protection of North-Rhine Westphalia (LANUV) in its technical report on climate (Johann et al. 2019). The report states that the thermal situation in the city is likely to be unfavourable in the future. It can therefore be considered a relevant target region of nocturnal cold-air drainage. Secondly, the upper Rur catchment, which is the main cold-air source for the Düren area, has a sufficiently large complexity of SCs to assume the applied methodology to be transferable to other complex catchments. Furthermore, the distinct valley structure of this catchment determines clear flow paths of the cold air. This makes it possible to divide the study area in cold-air SCs based on the geomorphological situation, having in mind that according to Sachsen et al. (2013) overflow scenarios from one valley into the other are possible nevertheless.

Source data and pre-processing
One digital elevation model (DEM), two land cover and two river network datasets were used as input for the cold-air drainage model KLAM_21. KLAM_21 simulates the nocturnal cold-air flow based on the two main inputs land cover and elevation for a standard computational domain of up to 3000 × 3000 grid cells. The model showed good agreement with field observations (Sievers and Kossmann 2016). For geodata processing the programming language R in combination with the integrated development environment RStudio and the geographic information system ArcGIS (version 10.7.1, ESRI) was used. The DEM was the European Digital Elevation Model (EU-DEM, version 1.1) provided by the European Environment Agency (EEA) in the frame of the European Union's Copernicus program (EEA 2016). It is a raster dataset with a spatial resolution of 25 m and is available in the ETRS89 LAEA coordinate system. For the current work, it was projected to ETRS 89 UTM 32 N, which was used for all geodata in this study, and cropped to the area of investigation described above. The result is shown in Fig. 1a.
For land-cover information two sources were used in order to create a contiguous cross-border dataset between the three countries of the study area while having data with high spatial resolution for the majority of the analysed catchment area. For the Belgian and Dutch parts, the Corine Land Cover 2018 (CLC 2018, version 20b2) raster dataset was used (EEA 2019) which has a spatial resolution of 100 m and the coordinate system ETRS89 LAEA. In order to match it to the DEM, the resolution was set to 25 m while preserving the coarser spatial quality of the original data. The 100 m resolution in the Belgian upper catchment part appears to be acceptable as it is mostly a plateau with large scale wood and peat bog. The data set was then projected to ETRS 89 UTM 32 N. Additionally, the Digitales Landschaftsmodell 50 (Digital Landscape Model 50, DLM), a vector dataset provided by the district government of Cologne (Geobasis NRW 2022) was used for the German part of the study area. Among other things, it contains polygon data for different landscape units, has a positional accuracy of 12 m and is already available in the ETRS 89 UTM 32 N coordinate system. A raster dataset that matches the resolution of the DEM was created from the DLM using its OBJART (object type) and BEB (type of building) attributes, which contain the necessary land-cover information for the new raster values. The land-cover values of the CLC 2018 and DLM datasets were then reclassified for their use in the KLAM_21 coldair model, which uses nine predefined land-cover classes. Table 1 shows which CLC 2018 and DLM values were converted into which KLAM_21 classes. Finally, the two datasets were merged and cropped to the study area. Figure 1b shows the resulting land-cover raster for the study area.
A similar approach was used for the river network data. The DLM also contains a set of line vectors that represent watercourses. These were combined with surface water body data from the WISE WFD Reference Spatial Datasets provided by EEA (EEA 2020) for the Belgian and Dutch parts to create a continuous river network for the study area.

Calculation of SC areas
The SC areas were calculated following a standard DEM based watershed delineation approach using the hydrology toolset of ArcGIS with the implemented D8 algorithm introduced by Jenson and Dominique (1988). Prior to the watershed delineation the combined river network dataset was used to prepare the DEM for input, using a so called "stream burning" procedure, whose method was introduced by Hutchinson (1989) and discussed by Saunders (1999). The idea is to artificially lower the DEM at the location of the streams prior to watershed delineation procedure in order to force the calculated flow directions to the correct watercourses. A comparison between the results of the watershed delineation with and without stream burning in the study area has shown that more realistic catchments were calculated using the stream burning approach, especially in its flat parts. Based on flow direction and flow accumulation, a raster cell was defined as part of a stream if it accumulates the flow of at least 50,000 raster cells which is an area of 31.25 km 2 . The watersheds were then calculated using flow direction and the intersection points of the derived stream network. The relevant 16 SCs of the upper Rur were selected and saved as a separate vector dataset (Fig. 1a).

Preparation of barriers around the SCs and implementation of KLAM_21 model runs
Different technical approaches for the exclusion of areas within the model were tested. The most promising results were archived by surrounding an SC with a two-dimensional barrier, similar to the experimental approach of King (1973). These barriers are normally used to represent objects like dams or walls. They have a blocking effect on the cold-air flow and can only be overflowed if the height of the cold-air layer is higher that the barrier itself (Sievers 2008). A barrier of sufficient height that encloses an SC should trap the cold air that develops within it. The geometries of the calculated SCs for the upper Rur were used to create the two-dimensional barriers that will exclude them from the numerical model. KLAM_21 uses a format called PL-file to implement vector data (Sievers 2008). A wall-like barrier can be constructed in this format specifying the xy-coordinates and height of each point that makes up the line representing the barrier. KLAM_21 then aligns the given shape to the model grid. Because the barriers are artificial objects, which interfere with the cold-air flow, their influence on the model environment should be minimised, except for their purpose of trapping the cold air in the excluded area. It would, e.g., not make sense to exclude a single SC in the middle of the study area, because cold-air flowing down from higher altitudes would be blocked by the barrier and alter its course arbitrarily. That is why every SC had to be excluded together with all upstream SCs at the same time. To ensure this, the SCs were excluded in a distinctive order and the shape of a newly excluded catchment was merged with the shape of the ones excluded before. To determine this order, the  (Verdin and Verdin 1999), using a reversed numbering in accordance to the non-systematic but rather cold-air runoff related downhill oriented research question. In order to realise this automatically, based on available datasets, an algorithm was developed in R. It takes two variables, the GKZ, a German watercourse identifier, and the hierarchy number after Shreve (1966) of the main watercourse in each SC. The GKZ was obtained from river network features in the DLM dataset.
The GKZ for the Rur is 282 while the GKZ of all its tributaries also begins with 282. The order of these tributaries is determined by the fourth digit. A tributary having a GKZ beginning with 2821 flows into the Rur before a tributary beginning with 2822. Tributaries of the tributaries are ordered in the same way. In the Shreve hierarchy, all source streams start with the value 1. If two streams merge, their values are added.
The SCs were first sorted by GKZ, interpreted as character string, and then by their Shreve values. The first SC of this pre-order list (POL) was therefore the upmost SC of the Rur (GKZ = 282, Shreve = 1). This SC was put in first position of the final-order list (FOL). Before adding the next SC of the Rur main stream to the FOL, the algorithm looks for a tributary, which would be the next SC in the POL that is one level higher in the GKZ system than the one last added to the FOL. If a tributary was found, the algorithm added it to the FOL and looked for a tributary of that tributary and, if it existed, added it to the FOL as well. This was done until no higher GKZ level was found. In that case, the algorithm looked for the next catchment in the POL with a GKZ value one level lower than the one of the SC last added to the FOL, which was the highest possible. In other words, we proceeded down the main flow. This was done as long as no further catchment was found by the algorithm. In this case, the lower end of the study area was reached. The final order of the 16 Rur SCs can be seen in Fig. 1a.
This order was now used to create the barriers for the KLAM_21 model runs. A total of 17 runs were made, one for the exclusion of each SC (run 1-16) and one reference run without barriers (run 0). The area excluded in run 1 is equal to the area of the first SC. For run 2-16, the area of the respective SC was merged with the area that was excluded in the run before. In run 3, e.g., the areas of SC1, 2 and 3 were merged and excluded jointly. An R-Script was written that automatically generated the required vector shapes for run 1-16 from ESRI Shapefiles of the ordered SCs and saved them in KLAM_21 readable PL-file format. Because the PL-files have a maximum number of possible points, the polygons of the SCs have been simplified in advance. The height of the barriers was set to 1000 m to ensure that no cold-air flow could pass them.
Finally, three data inputs were provided for all 17 model runs: DEM, land cover and the PL-file representing the excluded area of the respective run. The DEM and landcover datasets were the same for all runs and have been transformed to KLAM_21 readable ASCII files. The model area was defined by the study area described above and had an extent of 3000 × 3000 cells, each with a size of 25 × 25 m.

Processing of the model results
Six model output parameters were returned per raster field in hourly temporal resolution for 1-8 h after sunset for each run: heat deficit E(x, y) in J/m 2 , which is the total negative energy balance and main subject of the evaluation, effective cold-air layer depth H eff (x, y) in dm, vertically averaged horizontal wind components u(x, y) and v(x, y) in cm/s as well as the horizontal wind components at 2 m above ground uz(x, y) and vz(x, y) in cm/s (Sievers and Kossmann 2016). The outputs were raster datasets with the same extent and resolution as the inputs.
The hourly outputs for E(x, y) of reference run 0 were subtracted from the corresponding outputs of run 1-16 in order to get the respective change in heat deficit (CHD) attributable to the excluded area of each run. The CHD is interpreted as the cooling effect of the excluded catchment in the study area at a given hour. CHD it would be the effect of the excluded area of run i at t hours after sunset. Because the excluded area of run i was the combined area of SC1 to SC i , another step was necessary to determine the CHD for a single SC defined as CHDS it , with i being the number of the respective SC. For SC1 CHDS 1t = CHD 1t applies. For all other SCs CHDS it was calculated by formula (1).
Here, the difference between run i − 1 and run 0 was subtracted from the difference between run i and run 0 for a given hour after sunset and therefore determining the effect of the area that was additionally excluded in run i , which was the area of SC i . The result was 256 raster datasets representing the cooling effect of all 16 SCs for 1-8 h after sunset in absolute and relative values.

Evaluation
To quantify the cooling effect of each SC, the CHDS raster datasets were analysed in two different ways. First, the changes were extracted for sample points along the Rur valley. The sample points were created by setting one point at each river kilometre starting at the Perlenbach spring (SC2) and continuing down the main Rur valley. Sample point 1 (1) CHDS it = CHD it − CHD (i−1)t is the first sample point outside the borders of SCs 1 and 2. Additionally, the changes within the residential areas of the target settlements were analysed. Here the mean CHDS was calculated. Düren was chosen as an example of a city which is potentially vulnerable for heatwaves and therefore dependent on mitigation effects such as cold-air drainage (see Sect. 2.1). Furthermore, the effects in residential areas in the low mountain area itself were compared by data extracted for two districts of the city of Nideggen. Here the situation within the Rur valley, represented by Nideggen-Brück (referred to as Brück), was compared to the situation in the main district of Nideggen (referred to as Nideggen), above the valley.

Results
From the large amount of calculated results (effects of 16 SCs for 1-8 h after sunset for sample points and residential areas in absolute and relative values), a selection is presented here. In a spatial perspective, Fig. 2 shows the CHDS caused by the exclusion of SCs 2 (Fig. 2a) and 14 (Fig. 2b) compared to the reference run. The situation at the end of a typical mid-latitude summer night, 8 h after sunset, is displayed which is expected to show the maximum effect of the two SCs on the study area. Negative values indicate less heat deficit at the location of the pixel when the SC was excluded. Therefore, the opposite number of these values can be interpreted as cooling effect originating in the specific SC. SC2 is an example for an SC in the most elevated area of the analysed upper Rur catchment. SC14 represents an SC which lies in close proximity of Düren. The sample points along the Rur valley, on which the x-axis of Figs. 3 and 4 is based, are marked as black dots.
The effect is most intense in the parts that are located directly adjacent downstream of the SCs. This is also the case for all other SCs. Along the Rur valley, the effect got weaker until it faded out in the northern parts of the study area. While 8 h after sunset, the cooling effect caused by SC2 reached the city of Düren, these effects are between 0.1 and 1% and are therefore marginal. A different picture can be seen for SC14, where at the end of the night, the CHDS caused by its exclusions ranged between − 10 and − 15% around Düren. There were also effects on areas around the SCs where the cold-air drainage did not follow the hydrological valley structure. This can be seen at the southern and eastern border of SC2 and northwest of SC14. In longitudinal valley profiles, Fig. 3 compares the CHDS for three selected SCs at the locations of the sample points along the Rur valley at 4 h, 6 h, and 8 h after sunset in absolute and relative values. All of these SCs represent different parts of the study area. SC10 is part of the Urft catchment (SC6-12) which makes up the whole north-eastern part of the upper Rur catchment. SC15 is one of the SCs in the lowest part of the study area, already enclosing parts of Düren's residential areas. The spatial effect of SC2 at 8 h was already depictured in Fig. 2b which can now be compared to the value curve at the sample points and different time steps. The effects were stronger in close proximity to the SCs and generally stronger later at night, except in the case of SC15 where CHDS near the SC was stronger in the middle of the night. Also noticeable is the difference between the SCs for absolute and relative values. In case of SC15, for example, low absolute heat-deficit values in the reference run, downstream of the SC, led to high relative change, while the absolute effects were far weaker compared to SCs 2 and 10. SC10 is worth noticing because despite its remoteness, it was the SC of the Urft catchment with the largest CHDS within the main Rur valley. It also played an important role for heat deficit in Nideggen which is shown at the end of the chapter.  Figure 4 is the equivalent to Fig. 3 for SCs 13 and 14. These neighbouring SCs represent a similar level in catchment order and were those with the largest effect on the residential areas of Düren besides SC15. However, the effects at the sample points differed between the two SCs. The exclusion of SC14 generally produced larger CHDS compared to the reference run than the exclusion of SC13. For SC14, the relative effect on sample points number 43-46 was larger earlier in the night, which could not be seen in the absolute values.
In addition to the effects along the Rur valley, the CHDS was analysed for residential areas. Figure 5 shows the CHDS accountable to SCs 13 (a-c) and 14 (d-f) in spatial perspective for 4, 6 and 8 h after sunset with the residential areas of Düren shown as black polygons. The larger effect of SC14 compared to SC13 which could already be seen for the sample points, could also be observed here. At 4 h after sunset the cold air from both SCs led to similar effects in the range of 1-3%. At 6 h, the effect from SC14 was already stronger, a trend that persisted until the end of the night. Also noticeable is that a large part of heat deficit originating in SC14 did not even affect the main Rur valley but rather areas northwest of the SC. This phenomenon is weaker for SC13. At the north-western border of SC13, values were above 0, which means that after the exclusion of SC13, the heat deficit in this area increased. Figure 6 shows the mean CHDS within residential areas of Düren for different SCs in chronological sequence for all time steps in absolute and relative values. SC15, which is directly adjacent to Düren's residential areas and even contains some parts of it, has by far the strongest effect. For the first 4 h after sunset, it is also the only one with a relevant cooling effect on the target area. Among the other, more remote source areas SCs 13 and 14 had the strongest effect. Also, SC10, which has a remote location from the city, produced a mean effect of 4% at the end of the night. Leaving SC15 aside, a noteworthy effect could only be seen 4 h after sunset even for SCs 13 and 14. However, the mean effect of the more remote SCs summed up to 11.62% after 6 h, 19.93% after 7 h and 28.84% after 8 h.
In addition to Düren, the heat-deficit effects were analysed for Nideggen and Brück, where the situation above  Figure 7 shows spatial characteristics of cold-air drainage for SCs 4 and 10. Despite the fact that SC10 is further away from Nideggen (mind the scale difference between a-c and d-f), the effect caused by SC10 in the surrounding of the town was generally stronger than that of SC4, especially later at night. The CHDS values on the areas directly downstream of SC10 were the strongest of all analysed SCs. In case of Fig. 6 Mean change in heat deficit (mean CHDS) caused by exclusion of different SCs for residential areas in Düren. Note that only 72% of Düren's residential area was used for the calculation of the mean CHDS for SC15 as the other 28% lay within the SCs boundaries and where therefore left out Fig. 7 Change in heat deficit (CHDS) 4, 6 and 8 h after sunset after exclusion of SCs 4 and 10 compared to the reference run SC4, the effect of the directly adjacent downstream parts weakened until the end of the night but got stronger for the area around Nideggen. This could be either due to less cold-air drainage from SC4 (which is supported by sample point data not shown here) or because of more cold-air drainage from the Urft catchment, weakening the relative share of the cold-air drainage from SC4 or due to a combination of both. Figure 8a and b shows the mean absolute values in CHDS for the areas of Brück and Nideggen, attributable to the given SCs for 1-8 h after sunset. Brück, which lies at the bottom of the valley, was much more affected by coldair drainage than Nideggen centre. However, those SCs which are accountable for most of the heat deficit were the same in both areas. Also, the absolute values for Nideggen were more evenly distributed between the SCs than for the relative values ( Fig. 8c and d). For the mean relative values of CHDS, it is noticeable that the heat-deficit effect in Brück was more evenly distributed between the different SCs. SC10 had the highest relative effect in the last 3 h, but the difference to SCs 4 and 3 was only a few percent. A different picture can be seen for Nideggen where SC4 had the highest relative mean CHDS after 4 h and 5 h, SC3 after 6 h and 7 h and SC10 only after 8 h.

Sample points
The detailed analysis of the cooling effects at the sample points along the Rur valley suggests two general rules. First, a large influence was noticeable at points close to the respective SCs during the whole night. This influence decreased with increasing distance. This can also be seen on the maps of Figs. 2, 5 and 7, most clearly in Fig. 2b where the cooling effect originating in SC14 fanned out widely into the flat areas at the mountain foreland. Second, the cooling effect of a single SC (absolute and relative) usually gets stronger during the night, which was true in all cases for sample points with a distance > 10 km from an SC along the valley. These two basic rules show that the combination of distance and time plays an important role in the question of cold-air origin. A target area closer to an up-valley SC can benefit to a higher degree from cooling effects originating from it, conversely, an SC more remote to a target point generally has a lower effect. Also, the total cooling effect produced within an SC is not immediately available downvalley, because cold air from areas located far away from the mouth of the respective SCs main valley simply needs time, according to its velocity, to stream through the SC before it can affect other parts of the study area. Therefore, a point closer to an SC can benefit from a larger amount of cold-air produced within it. However, within this general framework, there were differences visible between the SCs which can presumably be explained by variations in their characteristics like shape, size, terrain and land cover as well as their spatial arrangement to each other. Also, the location of a target point within the catchment structure is important. Sample point 1 for example is most likely only affected by SCs 1 and 2 while sample point 43 is potentially affected by 14 different SCs.
In terms of terrain characteristics, SC15 was different from the other SCs as it is situated in the mountain foreland and, thus, relatively flat with few elevated parts in the south. In this setting, the cold air flow from the higher areas is likely to disperse more evenly throughout the flat terrain and towards the border of the SC and is not confined to narrow valleys like in SCs 1-14. This might have resulted in more stable availability of cooling energy at the sample points over time. This could partly explain the differences in Fig. 3 where the absolute cooling energy from SC15 stayed relatively constant throughout the night compared to SCs 2 and 10. The decreasing relative effect on near sample points during the night, which was noticeable for SCs 15 (Fig. 3) and 14 (Fig. 4), indicates that the absolute heat deficit originating from other SCs in this part of the valley grew faster than the effect of SCs 15 and 14.
Minor irregularities in Figs. 3 and 4 can be attributed to an uneven distribution of the sample points in the valley cross section as they are located along the main river, which is meandering with higher amplitude than the main valley. This could also explain a similar uneven distribution in Fig. 7. In general, the cooling effect was lower at the centre of the valley and higher at its margins. A possible explanation might be a lower wind speed at the edges, resulting in a more stable cold air layer and less mixing with the surrounding air.

Düren
The most relevant target areas to benefit from cold air ecosystem services in the form of nocturnal cooling effects are settlements. The analysis of settlements in Düren and Nideggen showed that the origin of these is often a question of size and distance of the source area, but not in all cases. For Düren, the most relevant of all considered source areas was SC15 (Fig. 6). It lies in direct vicinity to the target area and therefore has an immediate effect which is also stronger than that of any other SCs. This observation is consistent with the findings from the sample point analysis, where closer points where affected stronger than distant ones. The decreasing relative effect from 5 h after sunset can be explained by the incoming cold air flows from other SCs, which increase their share in total cooling effect at the same time.
Besides SC15, SCs 13 and 14 had the largest effect on residential areas in Düren (Fig. 6). Both are large enough to produce relevant amounts of cold air and lie in close proximity to the city's residential areas. SC14 has a larger share of unsealed open areas, which might be the reason for its stronger effect. However, a noteworthy effect was only present 4 h after sunset. For the first 3 h after sunset, when parts of the urban population are still likely to be outdoors, these source areas did not play an important role. Nevertheless, they can contribute to the cooling of the city at hours when people are dependent on a comfortable temperature to sleep, which is important as heat exposure has a negative effect sleeping quality (Okamoto-Mizuno and Mizuno 2012).
Looking at the other SCs, it could be seen that SC12, which is the next closest to Düren, and SC6, which is the largest SC, had a weaker effect than SCs 3, 4 and 10 which are more remote and smaller, respectively. For SC12, its large share of forest compared to unsealed open areas which indicates the production of less cold air (Kiefer and Zhong 2015) could provide an explanation, as well as its smaller size. In case of SC6 (upper valley of side river Urft), large parts of its cooling energy affected areas east of the Rur catchment (not shown here). A reason could be the terrain situation of SC6, with large parts of potential cold-air source areas lying relatively far west of the Urft valley, resulting in a strong cold-air drainage in eastern directions that is able to overflow the eastern watershed. Furthermore, the exclusion of SC6 counterintuitively led to more heat deficit from 3 h after sunset in the valley of the lower Urft (SC12) and its tributary Olef (SC11). The wind field of the reference run (not shown here) explains this phenomenon. After 3 h, large overflow effects from the Olef into the Urft valley and further over the eastern watershed of SC6 take place. This results in upstream flow of cold air within SC6 and drains off large amounts of cooling energy from the Olef valley. When SC6 was excluded and this overflow effect therefore prevented, the cooling energy from the Olef valley would affect the main Rur valley instead. The loss of cooling energy east of SC6 should therefore be lower by the amount of the blocked overflow effects in its west. SC10 is a medium sized SC and relatively remote to Düren. Despite its remoteness, it has a substantial cooling effect 7 and 8 h after sunset. It has a high share of unsealed open areas and only fragmented forest and settlements which could promote the production of large amounts of cold air and their undisturbed outflow.

Nideggen
The municipality of Nideggen with its districts Nideggen (uphill) and Brück (valley bottom) has completely different preconditions compared to Düren. Nideggen has less builtup area and is situated within the Eifel low mountains, in direct influence of the major cold-air flow in the Rur valley. The influence can be seen in the mean absolute values of CHDS in Fig. 8a and b, which were significantly stronger compared to Düren (Fig. 6a). The comparison of the two districts is also interesting. The same SCs had much greater absolute effects in Brück, which lies at the bottom of the valley and therefore within the major cold-air flow. Nideggen on the other hand is located above the valley and is therefore less influenced by cold air. The more evenly distributed relative values in Brück might indicate a more stable flow of cold air, while the higher and frequently changing relative values in Nideggen is a sign of a more unstable and unsystematic availability. An explanation for that might be that the Nideggen residential area is much larger and heterogeneous than that of Brück. In combination with overall less heat-deficit values it might increase the chance of outliers affecting the mean relative values. Nideggen and Brück lie at the northern margin of SC13. This SC, representing the directly adjacent upstream surroundings of the municipality, is likely to have a large effect in terms of heat deficit, similar to SC15 and Düren's residential areas.

Method limitations
A considerable part of the cooling energy originating in the SCs of the Rur catchment also affected areas in surrounding valleys, which might be attributed to overflow effects. Such effects have been observed in the area around Aachen and could partly be attributed to blocking effects of vegetation (Sachsen 2013). They may also occur if more cold air is produced within a catchment than it can contain  or might be due to interactions of cold-air streams during the confluence with other valleys. In some situations, the approach of using barriers for SC exclusion is able to visualize these overflow effects from the excluded SCs into adjacent valleys (Figs. 2, 5 and 7), which is a positive side effect of the approach. However, overflow effects from SCs into areas which were already excluded in the previous model run, cannot be seen. The possible overflow from SC2 into SC1 for example is not visible in Fig. 2. This is because the overflow from SC2 into SC1 is blocked in the model run in which only SC1 was excluded. Therefore, this flow cannot take its natural path through the study area. This may lead to counterintuitive results, especially in the areas with complex overflows patterns as described for SC6 in Sect. 4.2.1. Also, minor interactions of the barrier with the modelled cold air flow due to inaccuracies in the computation of the watersheds are to be expected. This has to be kept in mind when interpreting the results.
An adjustment of the method should be aimed at by which a comparison of the effect from the direct vicinity of the target areas can be compared with the effects from the more distant source areas mainly presented here. For example, the effects produced within SC13 are likely to be important for the first hours after sunset in the target areas of Nideggen and Brück because larger effects from the considered SCs could only be seen 4 h after sunset (Fig. 8). This effect however cannot be shown with the presented method, because heat deficit pattern within an excluded area can hardly be interpreted.

Conclusion
The results show that cooling effects of cold-air drainage at a specific target point or area neither originate equally from all parts of the catchment nor only from the direct target vicinity, but there is a complex relation of nearby and remote source areas several kilometres away. Generally, all source areas within the complex terrain of the study area, analysed by a system of SCs, prove to contribute to the major coldair flow towards and the cooling effect in the target areas. Usually, the cold-air flow followed the hydrological valley structure. In detail, a relationship between distance of the SC from the target area and the intensity and delay of the effect was visible.
The method described here can be used to analyse the cold-air heat deficit in an urban target area, attributable to a single SC in its remote surroundings. The connection of cold-air source areas and their effect on urban target areas was also addresses by previous studies. Grunwald et al. (2019) quantified the efficiency of cold-air paths that connect cold-air reservoir areas (CARA) with cold-air impact regions (CAIR) within the city borders of Braunschweig, Germany. Schau-Noppel et al. (2020) and Steigerwald et al. (2022) identified the connection between the most relevant source and target areas with linear trajectories. Both studies examined source areas within or near the city borders. These areas close to the target settlement seem to play the major role in providing cooling energy, especially in the early night. This could be seen when looking at the effects of SC15 on the residential areas of Düren in Fig. 6. By including a greater catchment area and showing the effect of more remote source areas, this study presents how the question of cold-air related ecosystem services could be analysed on a larger spatial scale within complex terrain. The presented method also indicates how cold-air source areas can be understood as whole parts of a greater cold-air catchment. This makes it possible to also include the effects of land cover types into the analysis which are less relevant for cold-air formation.
The case of Nideggen shows that settlements which are exposed to strong cold-air related cooling effects have a lower vulnerability to heat stress. This is a fact that might become a valuable quality for local recreation areas in the future. Consequently, a strong influence of cold-air drainage could have a positive impact on tourism.
Despite the limitations mentioned in Sect. 4.3., it has been shown that the method is suitable to show the amounts of cooling energy induced by remote cold-air drainage source areas on the scale of SCs. It can locate the origin of the cooling effects in a target area and enables to compare the impacts of different source areas with each other.