Validating the Mystic Three-dimensional Radiative Transfer Model with Observations from the Complex Topography of Arizona's Meteor Crater

The MYSTIC three-dimensional Monte-Carlo ra-diative transfer model has been extended to simulate solar and thermal irradiances with a rigorous consideration of topography. Forward as well as backward Monte Carlo simulations are possible for arbitrarily oriented surfaces and we demonstrate that the backward Monte Carlo technique is superior to the forward method for applications involving topography , by greatly reducing the computational demands. MYSTIC is used to simulate the short-and longwave radiation fields during a clear day and night in and around Ari-zona's Meteor Crater, a bowl-shaped, 165-m-deep basin with a diameter of 1200 m. The simulations are made over a 4 by 4 km 2 domain using a 10-m horizontal resolution digital elevation model and meteorological input data collected during the METCRAX (Meteor Crater Experiment) field experiment in 2006. Irradiance (or radiative flux) measurements at multiple locations inside the crater are then used to evaluate the simulations. MYSTIC is shown to realistically model the complex interactions between topography and the radiative field, resolving the effects of terrain shading, terrain exposure , and longwave surface emissions. The effects of surface temperature variations and of temperature stratification within the crater atmosphere on the near-surface longwave irradiance are then evaluated with additional simulations.


Introduction
Spatial variations in surface radiation and energy budgets within complex topography often produce a complex mosaic of different microclimates.The radiative energy exchange that plays the important role in producing these microclimates is complicated by the complex earth-atmosphere interface introduced by realistic topography and varying surface properties.Topographic effects such as terrain shading, terrain exposure (slope angle and azimuth of terrain), and terrain reflections greatly affect the shortwave radiation balance, while the longwave energy exchange is affected by emissions from higher surrounding terrain and by the blocking of outgoing radiation.Traditional radiative transfer models, assuming homogeneous flat surfaces, are not capable of reproducing the complex topographic variations of the radiative fluxes observed in three dimensional terrain with varying surface properties (Whiteman et al., 1989;Matzinger et al., 2003;Hoch and Whiteman, 2010).In recent years, 3-D radiative transfer models based on Monte Carlo photon tracing techniques have been developed to model radiative transfer in inhomogeneous cloudy atmospheres.Some of these include also topography, although mainly for the solar spectral range, (e.g., Weihs et al., 2000;Kylling and Mayer, 2001;Müller and Scherer, 2005;Chen et al., 2006;Helbig et al., 2009).For this study we implemented topography for the calculation of solar and thermal irradiance and radiance into MYSTIC (Monte Carlo code for the physically correct tracing of photons in cloudy atmospheres) (Mayer, 2009;Emde and Mayer, 2007).A mathematically rigorous formulation was used in forward and backward Monte Carlo mode.The model is described in Sect. 2. Appendices are provided to document features of the model that are relevant to its use in accurately representing three-dimensional topography.We then applied MYSTIC to simulate radiative transfer within the three-dimensional topography of Arizona's Meteor Crater using measured meteorological inputs and topography data.The simulations are tested against irradiance measurements made within the crater during the METCRAX 2006 meteorological experiment (Whiteman et al., 2008).The tests show that MYSTIC realistically models the complex interactions between topography and the radiative field, resolving the effects of terrain shading, terrain exposure, and longwave surface emissions.Following the successful evaluation of the model, the model is used in a sensitivity study to investigate the role of surface-air temperature discontinuities and temperature stratification on the longwave irradiance in the crater topography.

The 3-D radiative transfer model MYSTIC
MYSTIC is capable of radiative transfer calculations in three-dimensional atmospheres in plane-parallel and spherical geometry (Mayer, 2009;Emde and Mayer, 2007;Kylling and Mayer, 2001).The model, developed as a solver for the freely available libRadtran (Mayer and Kylling, 2005) radiative transfer package handles three-dimensional clouds, inhomogeneous surface albedo, and topography.The accuracy of MYSTIC was demonstrated previously by comparison with the one-dimensional DISORT code (Stamnes et al., 1988) and with other three-dimensional solvers during the I3RC (intercomparison of 3-D radiation codes) campaign, where deviations of less than 1% were found for well-defined conditions (Cahalan et al., 2005).Experimental validation of 3-D radiative transfer codes is very difficult in general: pronounced 3-D atmospheric effects on radiation are usually found only in cloudy cases, and the three-dimensional characterization of clouds is a complex task.Usually, the uncertainty in the cloud structure used as input to the calculation is known only with large uncertainty and the comparison between the simulated and measured radiation does not allow firm conclusions to be drawn concerning the accuracy of the model.Nevertheless, a comparison of observations and calculations during a total solar eclipse has shown that MYSTIC compares very well with observations also for three-dimensional setups (Emde and Mayer, 2007;Kazantzidis et al., 2007).
MYSTIC has recently been described in detail by Mayer (2009).Here we concentrate on the specific requirements for use of topography as a lower boundary condition in the model.In our case, elevation is specified on a rectangular grid with constant grid spacing dx and dy in x-and ydirections.Between the grid points, elevation is interpolated bi-linearly (see Fig. where the coefficients a, b, c, and d for each grid cell are unambiguously defined by the altitudes of the four enclosing grid points.To determine the coefficients, ( 1) is evaluated at the four vertices of each individual grid cell using the coordinates of each grid point (x ij ,y ij ,z ij ) which yields a linear equation system with four equations which is solved for the four unknowns a ij , b ij , c ij , and d ij .Thus, elevation in the model domain is completely defined by the set of coefficients a ij , b ij , c ij , and d ij of each grid cell (alternatively one could have split each grid cell into two triangles and thus describe the surface by a set of 2N triangular planes instead of N bent surfaces).With this description, reflection at slopes can be accurately simulated.As in any regular MYSTIC run, photons are traced from scattering to scattering where the jump width is sampled from Lambert-Beer's law.At the end of each path element the photon altitude is checked to see if it is below the highest surface altitude.If it is, the last photon path element is re-traced to determine whether the photon intersects the surface in any of the elevation grid cells along the photon path (Appendix A describes a numerically stable calculation of the intersection point).If an intersection is found, the photon is reflected according to the surface properties at the intersection point.Currently Lambertian reflection is assumed.Non-Lambertian bidirectional reflectance distribution functions (BRDFs) would be possible as well, but we are not aware of any BRDF data extending to zenith angles larger than 90 • which are needed when a photon is reflected down-slope.Reflection is determined by a random process: if a random number between 0 and 1 is smaller than the surface albedo, the photon is reflected and a new random photon direction is chosen; otherwise the photon path ends.Alternatively, the photon may be reflected always while multiplying the photon weight with the surface albedo.Which of those methods is more efficient depends on the application.At the surface we have to distinguish between slopeparallel and horizontal irradiance, E s and E h .The slopeparallel irradiance is the radiant flux on the sloped area divided by the surface area and describes the actual energy flux received by the surface.While this is the correct quantity for energy budget considerations, instruments often observe the horizontal irradiance, that is, the radiant flux on a horizontal surface divided by the horizontal surface area.
By counting of all photons received by the sloped surface we obtain the slope-parallel irradiance E s .In particular, where E 0 is the extraterrestrial irradiance, cosθ 0 is the cosine of the solar zenith angle which corrects for the slant incidence at top-of-atmosphere, A d is the area of the total model domain, and N 0 is the total number of photons traced.
is thus the radiant flux associated to the individual photon incident at TOA.By counting the weights w i of all photons reaching the respective surface pixel and dividing by the surface area A s of the sample pixel we finally obtain the slope-parallel irradiance E s .For our bi-linear description of the surface, the calculation of the surface area A s turns out to be rather complex and is described in Appendix B.
The calculation of the horizontal irradiance E h is a bit more complicated: for each photon we need to consider the probability that the photon had actually not hit the sloped surface but the horizontal surface.Figure 2 illustrates what needs to be done: to consider that the photon would have hit surface element dA h instead of dA s we need to multiply the photon weight with the ratio of the projected areas (using dA h =dA s cosδ).We finally obtain for the horizontal irradiance Please note the index i at θ h,i , θ s,i , and δ i which indicates that all three angles differ from photon to photon because the slope varies over each pixel.The θ h,i <90 • under the sum indicates that only photons with incident angles <90 • ("coming from above") are counted into the downward irradiance.While these angles are easily calculated by simple geometry, such a pre-factor may cause problems in terms of numerical noise.In particular the factor 1/cosθ s,i may become arbitrarily large.As a consequence, the result is affected by spikes which slow down the convergence or may even completely prevent convergence of the result.This is a strong motivation for using the backward Monte Carlo technique where the photons are started from the detector and traced backwards through the atmosphere towards the sun or the point where they are emitted thermally, see e.g., Emde and Mayer (2007).Please note that in the following we discuss only the diffuse irradiance component.The direct irradiance is easily evaluated by Lambert-Beer's law.
Backward horizontal irradiance E h is calculated exactly as without topography: photons are either emitted from a point location or from a random location within a given sample pixel, if the average over a certain area is desired as in the forward calculation.The photons are emitted at the surface elevation, into a random upward direction with a probability proportional to the cosine of the zenith angle θ h .Please note that some of the thus emitted photons immediately hit the slope and are reflected or absorbed accordingly.Slopeparallel irradiance needs some special attention: After selection of the start location (distributed randomly in the horizontal) each photon is assigned a start weight of cosδ i where δ i is the local slope angle, see Fig. 2.This (a) accounts for the fact that actually more photons start at steeper slopes which actually have a larger surface area; and (b) makes sure that surface-parallel irradiance relates to the sloped surface rather than to the horizontal surface.Photons are started in a random outward direction with a probability proportional to the cosine of the zenith angle θ s .Forward and backward irradiances calculated by MYSTIC generally agreed to much better than 1% and in all cases within the Monte Carlo photon noise.
Figure 3 illustrates that backward Monte Carlo is the method of choice for irradiance calculations with topography.The upper plot shows the results of 10 000 individual MYSTIC calculations with 1000 photons each.For this purpose, the cloudless US standard atmosphere (Anderson et al., 1986) was chosen.Spectral irradiance was calculated at 320 nm for solar zenith angle 60 • for a 2-D mountain with triangular shape (height 1 km, 45  that actually an array of mountains is calculated).The plot shows the diffuse downward irradiance averaged over the entire domain.The scatter of the blue dots (forward) is considerably larger than the scatter of the red dots (backward) and a number of large outliers is clearly visible in the forward case.These are caused by the 1/cosθ s,i in (4).Actually, the standard deviation is more than a factor of 2 smaller in the backward simulation compared to the forward calculation.This, together with the absence of spikes in the backward calculation and the fact that the computational times of forward and backward were nearly identical is a clear motivation for using backward rather than forward Monte Carlo for calculating horizontal irradiance in structured terrain.
Another important advantage of backward Monte Carlo which is relevant for slope-parallel as well as for horizontal irradiance calculations is the fact that with backward Monte Carlo the irradiance has to be computed only for those pixels which are really needed.With forward Monte Carlo, photons are started evenly distributed over the whole model domain; unless an average over the whole model domain is desired, only a small fraction actually hits the sample pixel which may be very small if e.g.irradiance at a specific location such as a mountain top is to be calculated.In contrast, all photons in backward Monte Carlo are started at the sample pixel and each photon contributes to the result.Also, the sample pixel in backward Monte Carlo may be arbitrarily small (even a point location) while in forward mode the noise increases approximately with 1/ sample pixel area.
In this paper, MYSTIC is used to simulate radiative transfer in the three-dimensional terrain of Arizona's Barringer Meteor Crater.The backward Monte Carlo technique has been used for all simulations to allow faster and computationally less expensive calculations of radiance and irradiance for individual sites where radiation measurements were made in the Meteor Crater topography.Instead of calculating irradiances for an entire large domain, which would be computationally expensive, backward calculations were made only for locations of interest.All calculations used a topographic domain of 4×4 km 2 , with a horizontal grid resolution of 10 m. Figure 4 shows a centered cut-out of the topography.For every irradiance calculation, 10 6 photons were traced.The isotropic reflectivity factors used for the different sites were determined from the measured monthly mean surface albedos at these individual sites.

Observations during METCRAX
During METCRAX 2006, extensive radiation data were collected in a unique, idealized topographic basin formed by the impact of a meteorite 50 000 years ago on the Colorado Plateau, 40 km east of Flagstaff, Arizona (Whiteman et al., 2008).Strong temperature inversions form in this basin during synoptically undisturbed nights in this basin, making it an ideal venue to study the different processes involved in stable boundary layer development.The idealized bowl shape facilitates studying the interactions between topography and the radiative field.Detailed observations of the shortwave and longwave components of the surface radiation budget were made at six different sites within the crater topography (Figure 4) as described in detail by Hoch and Whiteman (2010).Two of the sites were located over quasi-horizontal surfaces, one on the crater floor (FLR) and one on the crater rim (RIM).The remaining four sites were located on the sloping crater sidewalls, two on the west slope (West Upper, WU, and West Lower, WL) and two on the east slope (East Upper, EU, and East Lower, EL).At the slope sites, the instruments were oriented parallel to the underlying slope as estimated by eye over a slope area of 10s of meters around the radiometer.The site characteristics, including the inclination and azimuth angles of the radiometers are listed in Table 1.The four main components of the radiation balance (shortwave incoming, shortwave reflected, longwave incoming and longwave outgoing radiation) were measured individually, quasi-parallel to the slope of the underlying terrain, using Eppley PSP or Kipp & Zonen CM21 pyranometers and Eppley PIR or Kipp & Zonen CG4 pyrgeometers.Diffuse radiation was measured at the crater rim and floor sites with LiCor pyranometers and shadowbands.The uncertainties of the radiation instruments used during METCRAX were discussed in detail by Hoch and Whiteman (2010) and are based on a detailed study by Kohsiek et al. (2007) who reported the following accuracies: incoming shortwave: maximum of 5 W/m 2 or 1% of value; incoming longwave: 10 W/m 2 (daytime), 5 W/m 2 (nighttime); outgoing shortwave: maximum of 5 W/m 2 or 6% of value; outgoing longwave: 10 W/m 2 (daytime), 5 W/m 2 (nighttime).All instrumentation was supplied and installed by the National Center for Atmospheric Research (NCAR).
The detailed observations of the atmospheric temperature and humidity structure from a 10-m meteorological tower on the bottom of the crater, frequent tethersonde ascents throughout the crater atmosphere, and radiosoundings outside of the crater as described in Whiteman et al. (2008) are well-suited to construct the model input necessary for MYS-TIC.

Comparisons
Calculations with MYSTIC in both the shortwave and longwave spectrum were compared with observations made during METCRAX.The clear sky day of 21 October 2006 was selected for the comparison of solar irradiance using the Air Force Geophysics Laboratory (AFGL) mid-latitude standard atmosphere (Anderson et al., 1986) as input to the model.For the comparison of longwave irradiance, the night of 22-23 October 2006 was chosen.During this night, atmospheric profiles of temperature and humidity were available from three-hourly radiosonde ascents between 15:00 and 09:00 Mountain Standard Time (MST).These profiles, combined with those collected from tethersonde flights and a 10 m meteorological mast at the crater floor site were used as input data sets for the model calculations.

Shortwave irradiance
The modeled downward shortwave components of the radiation balance include direct and diffuse solar radiation.The modeled components were compared to observations at the crater floor (Fig. 5).There, prominent shadows cast from the surrounding crater rim strongly affect the radiation field in the morning and afternoon, a feature simulated by the model.The sum of the downward direct and diffuse components (global radiation), shortwave reflected radiation, and albedo (the ratio between shortwave reflected and global radiation) can be compared at all sites.
Observed and modeled direct and global solar radiation agree within 7±8 W/m 2 (mean bias±RMS difference) and 7±10 W/m 2 , respectively, except at the East Upper site where the RMS difference of the global irradiance reaches 19 W/m 2 (which is still less than 5% of the average irradiance).Possible reasons for this larger discrepancy are discussed below.Otherwise, the curves overlap at most times.In particular, the topographic effect of terrain shading is very well reproduced by MYSTIC: The morning shading by the east rim and the afternoon shading by the west rim both match well with the observations.In all cases except the East Upper site, the difference between observed and calculated direct and global irradiance is within or close to the uncertainty of the measurements (5 W/m 2 or 1% of the value, according to the manufacturer which does not include possible misalignment errors).
Global radiation (Fig. 6) is symmetric about solar noon at the crater floor and rim, where the instruments were mounted horizontally.The asymmetry of the diurnal variation of global radiation at the sloping sites -the maximum global radiation is received prior to solar noon on the west sidewall and after solar noon on the east sidewall -is well represented by the MYSTIC calculations.The topographic influence on the radiation balance by terrain exposure (slope and azimuth angle) is very well reproduced by MYSTIC.In addition, the timing of the shadows cast from surrounding topography in the morning and evening (terrain shading) matches well with the observations.Figure 7 compares the diurnal variation of observed and modeled shortwave reflected radiation.The model calculation is influenced by the choice of the isotropic reflectance factor that is prescribed for each model run.A value of 0.2 was used for all calculations, with the exception of the crater rim site, where 0.3 was chosen.These values represent the observed mean surface albedos reported by Hoch and Whiteman (2010).The use of these time-independent values with MYSTIC yielded shortwave outgoing irradiances that are slightly higher than those observed at most sites.With the exception of the East Upper site, however, MYSTIC realistically reproduced the site-specific diurnal variation seen in the observations.At East Upper, the influence of the site exposure was slightly underestimated, which may partly be due to a misrepresentation of the instrument measurement plane or due to a misrepresentation of the true local topography by the relatively coarse digital elevation model.
Figure 8 compares the diurnal variations of observed and modeled albedo at the horizontally oriented sensors 2 m above the underlying terrain at two selected sites, Rim and Floor.Although a constant surface reflectivity is prescribed in MYSTIC for each terrain pixel regardless of its orientation, both the modeled and simulated albedos are seen to vary with time.This diurnal variation in the observations is produced by the fact that the horizontal radiometers are exposed to reflections that come not only from the underlying near-horizontal surface but also from the surrounding complex terrain within the field of view of the sensor surfaces.The similarity in shape between calculated and observed diurnal variation is a strong test of MYSTIC and shows that the model successfully accounts for the reflections from all terrain pixels in the field of view of the 2-m sensor.MYS-TIC reproduced the smooth diurnal albedo cycle at the Floor site where reflections from the surrounding terrain affect both the incoming and outgoing irradiances.MYSTIC matched the distinctive shape of the diurnal albedo curve at the Rim where, because the site is on a high elevation ridgeline, the terrain reflections affect primarily the upwelling irradiance.At the Rim site, the terrain slopes downward to the west.Afternoon reflections from this lower terrain may enhance the upwelling irradiance, leading to higher albedos at the horizontally oriented instrument plane.At East Upper, the previously mentioned mismatch in reflected shortwave radiation led to a mismatch of the diurnal albedo pattern (not shown).Despite the relatively crude terrain resolution, MYSTIC is shown to successfully simulate the effects of topography on the diurnal cycle of albedo.

Longwave irradiance
Atmospheric profiles used as input for the MYSTIC longwave calculations were obtained from data collected during METCRAX.Half-hourly means of the temperature and humidity profiles from 4 levels of the 10 m meteorological tower on the crater floor were combined with temporal interpolations of the frequent sounding data of the central tethersonde inside the crater and the 3-hourly radiosonde profiles obtained just outside the crater.The modeled longwave irradiances were then compared with 30-min mean observations.For the MYSTIC calculations, the atmospheric profiles of temperature and humidity were assumed to be invariant across the crater from sidewall to sidewall, with the ground temperature being equal to the air temperature at that elevation.In future work, it will be possible to input a grid of surface temperatures.The present simulations thus do not take into account the shallow stable and unstable atmospheric layers that were observed over the crater sidewalls, something that we hope to improve in future work.
Figures 9 and 10 show the diurnal variation of the incoming and outgoing longwave irradiances as observed and modeled with MYSTIC.Similar to the shortwave irradiance, we find a close agreement between the observed and modeled incoming longwave irradiances: the mean difference in the diurnal variation between 15:00 MST and 09;00 MST does not exceed 8.5 W/m 2 at any of the 6 sites.The differences are mainly due to offsets, as indicated by standard deviations of the differences of less than 2 W/m 2 .Again, the differences are within or close to the instrument uncertainty: 10 W/m 2 (day) or 5 W/m 2 (night).
Differences between the modeled and observed outgoing irradiances are more apparent, but they are site specific.At Floor, the difference is only 3 W/m 2 ±3.6 W/m 2 (mean bias ±RMS difference) -better agreement than at the other sites is to be expected as the temperature input to the model is compiled from data taken not far away from the observation.At all other sites the model tends to overestimate the outgoing radiation during the night.This is attributed to the aforementioned near-surface inversions over the slopes of the crater topography that are not resolved with the 1D temperature profile observed over the crater center.Accordingly, shallow super-adiabatic layers that exist over the sunlit slopes in the afternoon and during post-sunrise periods are not resolved by the model.An underestimation of the outgoing longwave irradiance of between 20 and 60 W/m 2 is seen for these short time intervals, for example at West Upper and West Lower between 07:00 and 09:00 MST and between 15:00 and 17:00 MST at East Upper (no observation data at East Lower at this time).

Topographic effects on the radiation balance
Because MYSTIC performed well in simulating radiation observations in the basin, a second step was taken to run sensitivity simulations with MYSTIC focused on gaining improved understanding of the influences of atmospheric stability and air-ground temperature differences on the radiation field within the crater.For this purpose, the 18:00 MST sounding was modified so that an isothermal and constant mixing ratio atmosphere extended down into the crater from an elevation 500 m above the crater floor.The temperature of the lowest 500 m layer was 10 • C and the water vapour mixing ratio was 2.1 g/kg.Sensitivity simulations were then  Figure 11 shows results for the isothermal atmosphere with the three air-ground temperature discontinuities on a south-north transect across the crater basin.All curves show that longwave incoming radiation peaks when the terrain is steeply inclined.In steep terrain, more radiation is received from surrounding terrain and less from the sky.This topographic effect of terrain emissions increases (decreases) when the surrounding terrain is warmer (colder), as seen by larger (lower) variations in the longwave incoming radiation under 20 • C warmer (colder) terrain relative to the air.The upward longwave irradiance shows no widespread topographic effect.Only the increase (de-crease) with higher (lower) surface temperature according to the Stefan-Boltzmann Law is seen.
Figure 12 illustrates the effects of differing atmospheric stability on the longwave irradiance transect.Firstly, the mean longwave downward irradiance decreases under stable stratification.This reduction is partially due to the colder air overlying the terrain, and this effect is seen at the flat sections of the transect as well.Secondly, the effect of terrainemitted radiation and its growing influence with steepening terrain is seen.The total variation in longwave incoming irradiance along the transect is reduced under stable conditions.When moving along the transect toward the crater floor a compensating effect comes into play -the increase in irradiance caused by the larger contribution from terrain-emitted radiation is reduced when this contribution originates from terrain with colder surface temperatures.Colder surface temperatures are found at the lower elevations in a stable atmosphere.Thus, the incoming irradiance at the bottom of the crater can be lower than the incoming irradiance on the plain outside the crater.The With an assumption of a linear distribution of temperature with terrain height, the transect mirrors the crater topography.

Conclusions
The MYSTIC 3-D radiative transfer model has been extended to allow forward and backward calculations of horizontal and slope-parallel irradiance in structured terrain.It was shown that for allmost all applications involving topography the backward Monte Carlo method is to be preferred over the forward technique because (a) horizontal irradiance is affected by spikes in the forward mode; and (b) the backward technique simulates only those sample pixels or even point locations which are really needed while in the forward technique all sample pixels of the whole model domain need to be calculated.The backward Monte Carlo technique thus represents an enormous reduction in computational costs when compared with the more traditional forward approach.MYSTIC was validated by comparing simulations for the complex terrain environment of Arizona's Meteor Crater with data collected during METCRAX.MYSTIC provides accurate simulations of shortwave and longwave irradiances, reproducing the feedbacks between topography and the radiation field.Terrain shading, effects of varying terrain exposure, terrain reflection and longwave counter-radiation from surrounding terrain are represented by the model.
The comparison also pointed out areas for improvement.For the present study we used a constant temperature profile which means that inversion layers or superadiabatic nearsurface layers can not be prescribed in detail.As an improvement, a two-dimensional surface temperature map could be used in these situations, prescribing the surface temperature across the entire model domain and thus improving the representation of the outgoing longwave irradiances.
In two simple case studies we demonstrate how MYSTIC can be utilized to evaluate topographic feedback on the radiation field.Model calculations under different atmospheric stratifications, and with varying surface temperature boundary conditions, showed how the slope-parallel longwave irradiances are influenced by terrain and atmospheric emissions.
Future research with MYSTIC will include the calculation of radiative cooling and heating rates within the crater basin under varying atmospheric conditions.Comparisons with the observed temperature tendencies during METCRAX will help to understand the role of radiative cooling in stable boundary layer development.Parametric studies will further address the influences of basin size and shape, and influences of near-surface temperature gradients due to nighttime inversions and daytime superadiabatic sublayers.

Fig. 1 .
Fig. 1.Implementation of topography in MYSTIC.The illustration to the right shows one surface element where the elevation is interpolated bi-linearly between the four grid points.

Fig. 2 .
Fig. 2. Translation from slope-parallel to horizontal irradiance.dA S is the sloped surface element, dA h is the corresponding horizontal surface element.dA S and dA p are the corresponding projections normal to the photon direction.n S and n h are the normals to the sloped and horizontal surface elements.

Fig. 3 .
Fig. 3. (Top) Domain average diffuse downward horizontal irradiance calculated in forward and backward mode; for each of the 10 000 data points, 1000 photons were traced.(Bottom) Histogram of the diffuse irradiances shown in the top plot.

Fig. 4 .
Fig. 4. Topographic map showing the radiation measurement sites within the crater.The dashed line shows the location of the transect discussed later.

Fig. 6 .
Fig. 6.Observed (black) and modeled (red) slope-parallel global radiation at the six sites in the crater for 21 October 2006.

Fig. 8 .
Fig. 8. Observed (black) and modeled (red) albedo at the crater floor and rim sites for 21 October 2006.

Fig. 9 .
Fig. 9. Observed (black) and modeled (red) slope-parallel longwave incoming radiation at the six crater sites for the night of 22-23 October 2006.

Fig. 10 .
Fig. 10.Observed (black) and modeled (red) slope-parallel longwave outgoing radiation at the six crater sites for the night of 22-23 October 2006.

Fig. 11 .
Fig. 11.Modeled longwave downward (left) and upward (right) irradiances at 2 m above the topography under isothermal atmospheric stratification.Black: Surface temperature same as air temperature.Blue: Surface temperature 20 • C colder than air temperature.Red: Surface temperature 20 • C warmer than air temperature.Irradiances are represented by the scales on the left.Elevation cross sections are shown as thin dotted lines with the scale given on the right.

Fig. 12 .
Fig. 12. Modeled longwave downward (left) and upward (right) irradiances 2 m above the topography under different atmospheric stratifications, isothermal (black) and stable (5 • C/100 m, blue).The surface temperature was set to the temperature of the adjacent air.Irradiance values are shown on the left axes, while elevation values for the terrain cross section (thin dotted line) are shown on the axes.

Table 1 .
Radiation measurement sites and their characteristics.