Sub-kilometer Scale Research Articles

Comparison of Physically Based and Empirical Modeling of Nighttime Spatial Temperature Variability during a Heatwave in and around a City

Abstract The numerical weather prediction model HARMONIE-AROME and a multiple linear regression model (referred to in this article as the TURCLIM model after the local climate observation network) were used to model surface air temperature for 25–31 July 2018 in the City of Turku, Finland, to study their performance in urban areas and surrounding rural areas. The 0200 LT (local standard time) temperatures modeled by the HARMONIE-AROME and TURCLIM models were compared to each other and against the observed temperatures to find the model best suited for modeling the urban heat island effect and other spatial temperature variabilities during heatwaves. Observed temperatures were collected from 74 sites, representing both rural and urban environments. Both models were able to reproduce the spatial nighttime temperature variation. However, HARMONIE-AROME modeled temperatures were systematically warmer than the observed temperatures in stable conditions. Spatial differences between the models were mostly related to the physiographic characteristics: for the urban areas, HARMONIE-AROME modeled on average 1.4°C higher temperatures than the TURCLIM model, while for other land-cover types, the average difference was 0.51°C at maximum. The TURCLIM model performed well when the explanatory variables were able to incorporate enough information on the surrounding physiography. Respectively, systematic cold or warm bias occurred in the areas in which the thermophysically relevant physiography was lacking or was only partly captured by the model. Significance Statement As more and more people are living in an urban environment, the demand for accurate urban climate modeling is growing. This study aims to understand the differences between the numerical weather prediction and multiple linear regression modeling and their limitations in modeling surface air temperature in subkilometer scale. The case study shows that models are capable of predicting the spatial variation of 0200 LT nighttime temperature during a heatwave in a high-latitude coastal city. Both models are therefore valuable assets for city planners who need accurate information about the impacts of the physiography on the urban climate. The results indicate that to improve the performance of the models, more accurate physiographic description and higher spatial resolution of the models are needed.

Wind Speed Downscaling of the WRF Model at Subkilometer Scale in Complex Terrain for Wind Power Applications

Wind Speed Downscaling of the WRF Model at Subkilometer Scale in Complex Terrain for Wind Power Applications

Cold Fog Amongst Complex Terrain

Abstract Cold fog forms via various thermodynamic, dynamic, and microphysical processes when the air temperature is less than 0°C. It occurs frequently during the cold season in the western United States yet is challenging to detect using standard observations and is very difficult to predict. The Cold Fog Amongst Complex Terrain (CFACT) project was conceived to investigate the life cycle of cold fog in mountain valleys. The overarching goals of the CFACT project are to 1) investigate the life cycle of cold-fog events over complex terrain with the latest observation technology, 2) improve microphysical parameterizations and visibility algorithms used in numerical weather prediction (NWP) models, and 3) develop data assimilation and analysis methods for current and next-generation (e.g., subkilometer scale) NWP models. The CFACT field campaign took place in Heber Valley, Utah, during January and February 2022, with support from NSF’s Lower Atmospheric Observing Facilities (managed by NCAR’s Earth Observing Laboratory), the University of Utah, and Ontario Technical University. A network of ground-based and aerial in situ instruments and remote sensing platforms were used to obtain comprehensive measurements of thermodynamic profiles, cloud microphysics, aerosol properties, and environmental dynamics. Nine intensive observation periods (IOPs) explored various mountainous weather and cold-fog conditions. Field observations, NWP forecasts, and large-eddy simulations provided unprecedented data sources to help understand the mechanisms associated with cold-fog weather and to identify and mitigate numerical model deficiencies in simulating winter weather over mountainous terrain. This article summarizes the CFACT field campaign, its observations, and challenges during the field campaign, including real-time fog prediction issues and future analysis.

Detection of Land-Use Change and Rapid Recovery of Vegetation after Deforestation in the Congo Basin

Abstract The Congo Basin is severely understudied compared to other tropical regions; this is partly due to the lack of meteorological stations and the ubiquitous cloudiness hampering the use of remote sensing products. Clustering of small-scale agricultural deforestation events within the basin may result in deforestation on scales that are atmospherically important. This study uses 500-m MODIS data and the Global Forest Change dataset (GFC) to detect deforestation at a monthly and subkilometer scale and to quantify how deforestation impacts vegetation proxies (VPs) within the basin, the time scales over which these changes persist, and how they are affected by the deforestation driver. Missing MODIS data meant that a new method, based on two-date image differencing, was developed to detect deforestation on a monthly scale. Evaluation against the yearly GFC data shows that the highest detection rate was 79% for clearing sizes larger than 500 m2. Recovery to predeforestation levels occurred faster than expected; analysis of postdeforestation evolution of the VPs found 66% of locations recovered within a year. Separation by land-cover type also showed unexpected regrowth, as over 50% of rural complex and plantation land recovered within a year. The fallow period in the study region was typically short; by the sixth year after the initial deforestation event, ∼88% of the locations underwent a further considerable drop. These results show the importance of fine spatial and temporal information to assess Congo Basin deforestation and highlight the large differences in the impacts of land-use change compared to other rain forests.

Typhoon Track, Intensity, and Structure: From Theory to Prediction

To improve understanding of essential aspects that influence forecasting of tropical cyclones (TCs), the National Key Research and Development Program, Ministry of Science and Technology of the People’s Republic of China conducted a five-year project titled “Key Dynamic and Thermodynamic Processes and Prediction for the Evolution of Typhoon Intensity and Structure” (KPPT). Through this project, new understandings of TC intensification, including outer rainband-driven secondary eyewall formation and the roles of boundary layer dynamics and vertical wind shear, and improvements to TC data assimilation with integrated algorithms and adaptive localizations are achieved. To promote a breakthrough in TC intensity and structure forecasting, a new paradigm for TC evolution dynamics (i.e., the correlations, interactions, and error propagation among the triangle of TC track, intensity, and structure) is proposed; and an era of dynamic-constrained, big-data driven, and strongly coupled data assimilation at the subkilometer scale and seamless prediction is expected.

PyNAPLE: Lunar Surface Impact Crater Detection

ABSTRACT In the last 20 yr, over 600 impact flashes have been documented on the lunar surface. This wealth of data presents a unique opportunity to study the meteoroid flux of the Earth–Moon environment, and in recent years the physical properties of the impactors. However, other than through serendipitous events, there has not been yet a systematic search and discovery of the craters associated to these events. Such a meteoroid-crater link would allow us to get insight into the crater formation via these live observations of collisions. Here, we present the pynaple (python NAC Automated Pair Lunar Evaluator) software pipeline for locating newly formed craters using the location and epoch of an observed impact flash. We present the first results from pynaple, having been implemented on the 2017-09-27 impact flash. A rudimentary analysis on the impact flash and linked impact crater is also performed, finding that the crater’s ejecta pattern indicates an impact angle between 10–30°, and although the rim-to-rim diameter of the crater is not resolvable in current LRO NAC images, using crater scaling laws we predict this diameter to be 24.1–55.3 m, and using ejecta scaling predict a diameter of 27.3–37.7 m. We discuss how pynaple will enable large scale analyses of sub-kilometer scale cratering rates and refinement of both scaling laws, and the luminous efficiency.

An Application of Sea Ice Tracking Algorithm for Fast Ice and Stamukhas Detection in the Arctic

For regional environmental studies it is important to know the location of the fast ice edge which affects the coastal processes in the Arctic. The aim of this study is to develop a new automated method for fast ice delineation from SAR imagery. The method is based on a fine resolution hybrid sea ice tracking algorithm utilizing advantages of feature tracking and cross-correlation approaches. The developed method consists of three main steps: drift field retrieval at sub-kilometer scale, selection of motionless features and edge delineation. The method was tested on a time series of C-band co-polarized (HH) ENVISAT ASAR and Sentinel-1 imagery in the Laptev and East Siberian Seas. The comparison of the retrieved edges with the operational ice charts produced by the Arctic and Antarctic Research Institute (Russia) showed a good agreement between the data sets with a mean distance between the edges of <15 km. Thanks to the high density of the ice drift product, the method allows for detailed fast ice edge delineation. In addition, large stamukhas with horizontal size of tens of kilometers can be detected. The proposed method can be applied for regional fast ice mapping and large stamukhas detection to aid coastal research. Additionally, the method can serve as a tool for operational sea ice mapping.

Spatially distributed simulations of the effect of snow on mass balance and flooding of Antarctic sea ice

AbstractSouthern Ocean sea ice can exhibit widespread flooding and subsequent snow-ice formation, due to relatively thick snow covers compared to the total ice thickness. Considerable subkilometer scale variability in snow and ice thickness causes poorly constrained uncertainties in determining the amount of flooding that occurs. Using datasets of snow depth and ice thickness acquired in the Weddell Sea during austral winter 2013 (AWECS campaign) from three floes, we demonstrate large spatial variability of a factor 10 and 5 for snow and combined snow and ice thickness, respectively. The temporal evolution after the floe visit was recorded by automatic weather station and ice mass balance buoys. Using a physics-based, multi-layer snow/sea ice model in a one-dimensional and distributed mode to simulate the thermodynamic processes, we show that the distributed simulations, modeling flooding across the entire heterogeneous floe, produced vastly different amounts of flooding than one-dimensional single point simulations. Three times the flooding is produced in the one-dimensional simulation for the buoy location than distributed (floe-averaged) simulations. The latter is in close agreement with buoy observations. The results suggest that using point observations or one-dimensional simulations to extrapolate processes on the floe-scale can overestimate the amount of flooding and snow-ice formation.

Satellite-derived quantification of the diurnal and annual dynamics of land surface temperature

The diurnal and annual cycles of Earth drive the land surface temperature (LST) variations in space and time but their combined effects are conveyed differently from place to place as a function of climate, local weather conditions, and surface characteristics. Fitting a model of the diurnal or annual LST cycle (DTC and ATC, respectively) to satellite data yields several parameters that summarize the thermal dynamics of each pixel. These parameters are distinct for each place and reveal how the ATC and DTC differ between pixels. Currently, due to the lack of LST with high spatial and temporal resolution, it is impossible to characterize both cycles at kilometer or sub-kilometer scale. This work addresses this problem and presents a method for simultaneously quantifying the annual and diurnal LST dynamics at kilometer scale. The proposed method is based on Annual Cycle Parameters (ACP) retrieved from geostationary satellite data for every 30 min between 00:00 and 23:30 local time and uses statistical downscaling to enhance their coarse spatial resolution. We apply our method to continental Europe for the period 2009–2013 and validate it using independent ACP retrieved from satellite and in-situ LST. Our results reveal how the climate and surface conditions influence the shape of the ACP diurnal profiles over Europe and illustrate for the first time the excellent agreement between satellite-derived ACP and ACP derived from in-situ LST. Even though our method is limited by the spatial coverage of the geostationary data, it provides a novel way for describing Earth's thermal response to solar heating.

Estimating root zone soil moisture across diverse land cover types by integrating in-situ and remotely sensed data

Many soil moisture networks monitor only one land cover type, typically grassland, and the availability of in-situ soil moisture data in other land cover types is severely limited. Satellite-based radiometers lack adequate resolution to match the spatial variability in land cover, which often occurs at the sub-kilometer scale. Thus, spatial and temporal dynamics of root zone soil moisture in regions with heterogeneous land cover types remain poorly understood. Our objective was to determine how effectively root-zone soil moisture for diverse land cover types can be estimated using a water balance model driven by normalized high-resolution, remotely sensed vegetation indices (VI) data and in-situ meteorological data. Root zone soil moisture dynamics under four different land cover types were estimated using normalized VI data as a proxy for the basal crop coefficient. Correlation coefficients (r) between measured and modeled soil moisture ranged from 0.50–0.92, mean absolute error (MAE) ranged from 0.03–0.06 m3 m−3, and mean bias error (MBE) ranged from -0.05–0.02 m3 m−3 across tallgrass prairie, cropland, mixed hardwood forest, and loblolly pine plantation sites. Model-estimated soil moisture under each land cover type was more accurate than both measured data from the nearest long-term grassland monitoring site and data from the NASA-USDA Enhanced Soil Moisture Active-Passive (SMAP) soil moisture product, providing evidence that in-situ meteorological data and remotely sensed VI data may be integrated into a simple water balance model to better estimate root zone soil moisture across diverse land cover types.

An Experiment on the Prediction of the Surface Wind Speed in Chongli Based on the WRF Model: Evaluation and Calibration

In this study, the ability of the Weather Research and Forecasting (WRF) model to generate accurate near-surface wind speed forecasts at kilometer- to subkilometer-scale resolution along race tracks (RTs) in Chongli during the wintertime is evaluated. The performance of two postprocessing methods, including the decaying-averaging (DA) and analogy-based (AN) methods, is tested to calibrate the near-surface wind speed forecasts. It is found that great uncertainties exist in the model’s raw forecasts of the near-surface wind speed in Chongli. Improvement of the forecast accuracy due to refinement of the horizontal resolution from kilometer to subkilometer scale is limited and not systematic. The RT sites tend to have large bias and centered root mean square error (CRMSE) values and also exhibit notable underestimation of high-wind speeds, notable overestimation or underestimation of the near-surface wind speed at high altitudes, and notable underestimation during daytime. These problems are not resolved by increasing the horizontal resolution and are even exacerbated, which leads to great challenges in the accurate forecasting of the near-surface wind speed in the competition areas in Chongli. The application of postprocessing methods can greatly improve the forecast accuracy of near-surface wind speed. Both methods used in this study have comparable abilities in reducing the (positive or negative) bias, while the AN method is also capable of decreasing the random error reflected by CRMSE. In particular, the large biases for high-wind speeds, wind speeds at high-altitude stations, and wind speeds during the daytime at RT stations can be evidently reduced.

Thickness of orthopyroxene-rich materials of ejecta deposits from the south pole-Aitken basin

The South Pole-Aitken (SPA) basin is the largest impact structure on the Moon and is believed to have excavated the orthopyroxene (Opx)-rich lower crust and/or upper mantle materials. For the complex craters outside of the SPA transient cavity, the origin of the Opx-rich central peaks is possibly from either the Opx-rich materials of the SPA ejecta deposits or the unexcavated lower crust and/or upper mantle. To estimate the thickness of the Opx-rich materials of the SPA ejecta deposits, this study investigated large complex craters (dozens of kilometers) that have penetrated the Opx-rich materials and exposed deeper mafic-poor crustal material based on spectra extracted from small fresh craters (sub-kilometer scale). The amount of foreign material introduced to these large complex craters by other lunar impact events was estimated to guarantee the least influence on the compositional analysis. The study results suggest that a 56 km-diameter crater at ~640 km northwest of the SPA center is large enough to penetrate the Opx-rich materials of the SPA ejecta deposits, which are thinner than ~4.7 km at this location. This result also indicates that the intense bombardment history of the large craters and basins outside of the SPA transient cavity excavated and redistributed a large amount of crustal material across the basin, possibly resulting in the heterogeneous distribution of mafic-rich and mafic-poor materials on the SPA surface.

Evaluating Climate Change Impacts on Soil Moisture and Groundwater Resources Within a Lake‐Affected Region

AbstractThe impact of climate change on surface water resources is reasonably well studied. However, the impact on groundwater resources has only been considered by a few studies worldwide. Here we present an analysis of climate change impacts on groundwater resources in a well‐instrumented 6,800‐km2 watershed in the Laurentian Great Lakes Basin. We employ a physics‐based modeling pipeline consisting of an ensemble of high‐resolution regional climate model projections based on the Weather Research and Forecasting model and the fully integrated three‐dimensional hydrologic model HydroGeoSphere. The Weather Research and Forecasting model is run at a resolution as fine as 10 km using two different physics configurations, while HydroGeoSphere simulates the terrestrial hydrosphere at subkilometer scale, from deep groundwater to surface water, including surface water‐groundwater interactions. The two Weather Research and Forecasting model physics configurations exhibit opposite climate change responses in summer precipitation. The hydrologic simulations follow the climate forcing, but due to the memory of the subsurface, differences in summer affect the entire seasonal cycle. In the drier climate scenario groundwater levels and recharge decline, while in the wetter scenario groundwater levels rise (recharge remains unchanged). Soil moisture changes accordingly, but primarily in late summer. It is also shown that the magnitude of climate change impacts on groundwater is strongly modulated by local physiographic features. In particular, regions where the groundwater table is deep (below 2 m; 15% of the area) show a high sensitivity to changes in climate forcing. Furthermore, changes in groundwater levels, recharge, and soil moisture typically occur in the same regions, suggesting potentially compounding impacts.

Finite-Element barotropic model for the Indian and Western Pacific Oceans: Tidal model-data comparisons and sensitivities

In this study, a 9.6 million node large-scale unstructured grid finite-element forward barotropic model is developed and applied to understand the tidal dynamics and dissipation mechanisms of the Indian and western Pacific Oceans down to sub-kilometer scale at the coast. Tidal model-data comparisons are presented to assess the capabilities and limitations of our large-scale barotropic model. The average root-mean-square (RMS) discrepancies of tidal elevations at coastal tide gauges is 14 cm, which is ∼ 3 cm smaller than those of a state-of-the-art global data assimilated barotropic tidal model. Sensitivities to lateral boundary conditions, bathymetry, and dissipative processes are explored to guide future endeavors related to large-scale barotropic modeling in the region and other regions throughout the world. Lateral boundary conditions are found to induce adverse resonant effects on the lunar semi-diurnal modes when poorly placed elevation specified boundary conditions are used. This problem is largely resolved by using an absorption-generation layer at the boundary. Parameterization of internal tide energy conversion is identified as the most important aspect to control deep water solutions, and help reduce the RMS discrepancies of the entire system. Two forms of this parameterization are presented and their spatial distributions of dissipation are compared. Bathymetry has a negligible effect on the tidal solutions in deep water, but local high resolution bathymetry results in significant reductions to the average RMS discrepancies on the continental shelf (26%) and at the coast (30%). Implementing a spatially varying bottom friction coefficient based on sediment types decreases the average RMS discrepancy at the coast by 9% predominantly due to its positive effects in the Yellow Sea. The model is shown to capture a large amount of the tidal physics and has the potential for application to a range of barotropic problems such as wind-driven surge and tidal processes.

The prevalence of kilometer-scale heterogeneity in the source region of MORB upper mantle.

The source regions of mid-ocean ridge basalts (MORB) are heterogeneous, consisting of chemically and lithologically distinct domains of variable size. Partial melting of such heterogeneous mantle sources gives rise to diverse isotopic compositions of MORB and abyssal peridotites. Variations in radiogenic isotope ratios in MORB are attributed to mixing of melts derived from enriched and depleted mantle components. However, melt mixing alone cannot fully account for the difference between the average 143Nd/144Nd in abyssal peridotites and their spatially associated MORB. We show that the more depleted Nd isotope composition in abyssal peridotites is a natural consequence of melt migration-induced mixing or smearing in the melting column. Sub-kilometer scale enriched mantle components or heterogeneities are significantly damped or homogenized in both the residue and erupted melt during their transit through the melting region. Heterogeneities with larger size and higher incompatible trace element abundance are more resistive to the mixing processes. The size-sensitive mixing depends on a parameter called the enrichment strength, which is the product of the heterogeneity size and the ratio between incompatible trace element abundance in the enriched and depleted mantle sources. Observed Nd-Hf isotope variations in MORB and abyssal peridotites can be reproduced if the enrichment strength is 20 to 60 km. These heterogeneities could be on the kilometer scale and have similar isotope ratios to but less incompatible trace element abundances than recycled oceanic crust.

Distributed sensing of ionospheric irregularities with a GNSS receiver array

AbstractWe present analysis methods for studying the structuring and motion of ionospheric irregularities at the subkilometer scale sizes that produce L band scintillations. Spaced‐receiver methods are used for Global Navigation Satellite System (GNSS) receivers' phase measurements over approximately subkilometer to kilometer length baselines for the first time. The quantities estimated by these techniques are plasma drift velocity, diffraction anisotropy magnitude and orientation, and characteristic velocity. Uncertainties are quantified by ensemble simulation of noise on the phase signals carried through to the observations of the spaced‐receiver linear system. These covariances are then propagated through to uncertainties on drifts through linearization about the estimated values of the state. Five receivers of SAGA, the Scintillation Auroral Global Positioning System (GPS) Array, provide 100 Hz power and phase data for each channel at L1 frequency. The array is sited in the auroral zone at Poker Flat Research Range, Alaska. A case study of a single scintillating satellite observed by the array is used to demonstrate the spaced‐receiver and uncertainty estimation process. A second case study estimates drifts as measured by multiple scintillating channels. These scintillations are correlated with auroral activity, based on all‐sky camera images. Measurements and uncertainty estimates made over a 30 min period are compared to a collocated incoherent scatter radar and show good agreement in horizontal drift speed and direction during periods of scintillation for which the characteristic velocity is less than the drift velocity.

Stalled Improvement in a Numerical Weather Prediction Model as Horizontal Resolution Increases to the Sub-Kilometer Scale

Stalled Improvement in a Numerical Weather Prediction Model as Horizontal Resolution Increases to the Sub-Kilometer Scale

ENCELADUS’ GEYSERS: RELATION TO GEOLOGICAL FEATURES

We apply histogram analysis, photogeological methods, and tidal stress modeling to Porco et al.'s survey of 101 Enceladus South Polar Basin geysers and their three-dimensional orientations to test if the jet azimuths are influenced by their placement relative to surface morphology and tectonic structures. Geysers emplaced along the three most active tiger stripe fractures (Damascus Sulcus, Baghdad Sulcus, and Cairo Sulcus) occur in local groupings with relatively uniform nearest-neighbor separation distances (~5 km). Their placement may be controlled by uniformly spaced en echelon Riedel-type shear cracks originating from left-lateral strike-slip fault motion inferred to occur along tiger stripes. The spacing would imply a lithosphere thickness of ~5 km in the vicinity of the tiger stripes. The orientations of tilted geyser jets are not randomly distributed; rather their azimuths correlate with the directions either of tiger stripes, cross-cutting fractures, or else fine-scale local tectonic fabrics. Diurnal tidal stress modeling suggests that periodic changes of plume activity are significantly affected by cross-cutting fractures that open and close at different times than the tiger stripes that they intersect. We find evidence of sub-kilometer scale morphological modification of surface geological features surrounding geysers from sublimation-aided erosion, and ablation, and scouring. We propose that the simultaneous crushing and shearing action of periodic transpressional tidal stress on ice condensing on the inside walls of geyser conduits is the mechanism that extrudes the peculiar, paired narrow ridges known as shark fins that flank the medial tiger stripe fissures. We present a gallery of high-resolution image mosaics showing the placement of all the jets in their source region and consequently their geological context.

Steep spatial gradients of volcanic and marine sulfur in Hawaiian rainfall and ecosystems

Sulfur, a nutrient required by terrestrial ecosystems, is likely to be regulated by atmospheric processes in well-drained, upland settings because of its low concentration in most bedrock and generally poor retention by inorganic reactions within soils. Environmental controls on sulfur sources in unpolluted ecosystems have seldom been investigated in detail, even though the possibility of sulfur limiting primary production is much greater where atmospheric deposition of anthropogenic sulfur is low. Here we measure sulfur isotopic compositions of soils, vegetation and bulk atmospheric deposition from the Hawaiian Islands for the purpose of tracing sources of ecosystem sulfur. Hawaiian lava has a mantle-derived sulfur isotopic composition (δ34S VCDT) of −0.8‰. Bulk deposition on the island of Maui had a δ34S VCDT that varied temporally, spanned a range from +8.2 to +19.7‰, and reflected isotopic mixing from three sources: sea-salt (+21.1‰), marine biogenic emissions (+15.6‰), and volcanic emissions from active vents on Kilauea Volcano (+0.8‰). A straightforward, weathering-driven transition in ecosystem sulfur sources could be interpreted in the shift from relatively low (0.0 to +2.7‰) to relatively high (+17.8 to +19.3‰) soil δ34S values along a 0.3 to 4100ka soil age-gradient, and similar patterns in associated vegetation. However, sub-kilometer scale spatial variation in soil sulfur isotopic composition was found along soil transects assumed by age and mass balance to be dominated by atmospheric sulfur inputs. Soil sulfur isotopic compositions ranged from +8.1 to +20.3‰ and generally decreased with increasing elevation (0–2000m), distance from the coast (0–12km), and annual rainfall (180–5000mm). Such trends reflect the spatial variation in marine versus volcanic inputs from atmospheric deposition. Broadly, these results illustrate how the sources and magnitude of atmospheric deposition can exert controls over ecosystem sulfur biogeochemistry across relatively small spatial scales.

Wintertime Subkilometer Numerical Forecasts of Near-Surface Variables in the Canadian Rocky Mountains

Abstract Numerical weather prediction (NWP) systems operational at many national centers are nowadays used at the kilometer scale. The next generation of NWP models will provide forecasts at the subkilometer scale. Large impacts are expected in mountainous terrain characterized by highly variable orography. This study investigates the ability of the Canadian NWP system to provide an accurate forecast of near-surface variables at the subkilometer scale in the Canadian Rocky Mountains in wintertime when the region is fully covered by snow. Observations collected at valley and high-altitude stations are used to evaluate forecast accuracy at three different grid spacing (2.5, 1, and 0.25 km) over a period of 15 days. Decreasing grid spacing was found to improve temperature forecasts at high-altitude stations because of better orography representation. In contrast, no improvement is obtained at valley stations due to an inability of the model to fully capture at all resolutions the intensity of valley cold pools forming during nighttime. Errors in relative humidity reveal that the model tends to overestimate relative humidity at all resolutions, without improvement with decreasing grid spacing. Wind speed forecasts show large improvements with decreasing grid spacing for high-altitude stations exposed to or sheltered from wind. However, no systematic improvement with decreasing grid spacing is found for all stations, which is similar to previous studies. In addition, the model’s sensitivity at subkilometer grid spacing is investigated by evaluating the effects of (i) accounting for additional drag generated by subgrid orographic features, (ii) considering slope angle and aspect on surface radiation, and (iii) using high-resolution initialization for the surface fields.