The Cryosphere Multi-channel ground-penetrating radar to explore spatial variations in thaw depth and moisture content in the active layer of a permafrost site

Multi-channel ground-penetrating radar (GPR) was applied at a permafrost site on the Tibetan Plateau to investigate the influence of surface properties and soil texture on the late-summer thaw depth and average soil moisture content of the active layer. Measurements were conducted on an approximately 85 × 60 m2 sized area with surface and soil textural properties that ranged from medium to coarse textured bare soil to finer textured, sparsely vegetated areas covered with fine, wind blown sand, and it included the bed of a gravel road. The survey allowed a clear differentiation of the various units. It showed (i) a shallow thaw depth and low average soil moisture content below the sandcovered, vegetated area, (ii) an intermediate thaw depth and high average soil moisture content along the gravel road, and (iii) an intermediate to deep thaw depth and low to intermediate average soil moisture content in the bare soil terrain. From our measurements, we found hypotheses for the permafrost processes at this site leading to the observed latesummer thaw depth and soil moisture conditions. The study clearly indicates the complicated interactions between surface and subsurface state variables and processes in this environment. Multi-channel GPR is an operational technology to efficiently study such a system at scales varying from a few meters to a few kilometers. Correspondence to: U. Wollschl̈ager (ute.wollschlaeger@ufz.de)


Introduction
Permafrost on the Tibetan Plateau reacts very sensitively to climate change (Liu and Chen, 2000;Zhao et al., 2010).Current observations show, for instance, an increase in ground temperatures and active layer thickness (Cheng and Wu, 2007;Zhao et al., 2010), and a related lowering of ground water and lake water levels and desertification that is caused by the drier ground surface (Cheng and Wu, 2007).
In cold regions, the active layer is "the layer of ground that is subject to annual thawing and freezing in areas underlain by permafrost" (van Everdingen, 1998).In this environment, it plays a significant role since almost all biogeochemical, hydrological, ecological, and pedogenic processes take place in this uppermost part of the soil (Kane et al., 1991;Hinzman et al., 2005).The thaw depth of the active layer determines, for instance, the partitioning of surface and groundwater runoff, it influences the rate of water and gas exchange between the soil and the atmosphere (Duguay et al., 2005) and the occurence of specific vegetation communities (Kane et al., 1991;Walker et al., 2003).
From a large-scale perspective, the thaw depth of the active layer primarily depends on temperature and the length of the thaw season (Hinzman et al., 2005).The seasonal thaw depth, however, is locally influenced by various additional factors that affect the local microclimate and surface energy balance.Important factors are surface temperature, physical and thermal properties of the surface cover and the substrate, vegetation cover, soil moisture content, albedo, and thickness and duration of snow cover (Duguay et al., 2005;Hinzman et al., 2005).
Published by Copernicus Publications on behalf of the European Geosciences Union.
Monitoring of active layer freeze and thaw processes is most often based on point measurements like boreholes (e.g., Brown et al., 2000;Wang et al., 2000).The standard protocol of the Circumpolar Active Layer Monitoring (CALM) program, for instance, recommends mechanical probing of the end-of-season thaw depth of the active layer with a frost probe as standard monitoring procedure (Brown et al., 2000;Smith and Brown, 2009).However, at point locations, yearto-year and microtopographic variations in surface temperature and soil moisture content can induce large spatial variations in thaw depth (Lemke et al., 2007).Examples for active layer thaw depth variability are given, e.g., by Wright et al. (2009) for small-scale variations and by Nelson et al. (1998) and Hinkel and Nelson (2003) for larger areas.Hence, in order to obtain a deeper understanding of spatial active layer processes, efficient methods which allow a high-resolution mapping of active layer conditions are desired.
In order to derive large-scale spatial information about permafrost conditions, various satellite remote sensing techniques are currently being applied (Zhang et al., 2004;Duguay et al., 2005).So far, the challenge with mapping near-surface permafrost conditions from space is that only information about surface properties is available which then has to be related to complicated subsurface phenomena (Zhang et al., 2004;Duguay et al., 2005).
At an intermediate scale, ranging from a few meters to a few kilometers, ground-based geophysical measurements yield spatially highly resolved information on material properties and state variables like soil moisture content.In permafrost regions, ground-penetrating radar (GPR) is well suited for studying near-surface phenomena which are related to soil moisture and ice content, respectively.GPR by now has been applied in various permafrost studies, e.g., for mapping the depth of the permafrost table (Arcone et al., 1998;Hinkel et al., 2001;Moorman et al., 2003;Schwamborn et al., 2008) or layer boundaries in frozen ground (Hinkel et al., 2001), exploring massive ground ice (Annan and Davis, 1976) and pingo ice (Yoshikawa et al., 2006), ice wedges (Hinkel et al., 2001;Munroe et al., 2007), or thaw depth beneath streams (Bradford et al., 2005;Brosten et al., 2006Brosten et al., , 2009)).In addition, GPR has been employed to infer the near-surface moisture content of the active layer from the analysis of ground wave measurements (Moorman et al., 2003).
The measurement principle of GPR is based primarily on the analysis of the propagation velocity of electromagnetic waves travelling through the ground.The propagation velocity predominantly depends on the soil's dielectric permittivity which itself is directly related to the water and ice content, respectively (Davis and Annan, 1989).Reflections occur at interfaces that are related to abrupt changes in dielectric permittivity, e.g., at layer boundaries, at a groundwater table in a coarse-textured soil or at the frozen-unfrozen ground interface.Due to the strong dielectric contrast between the frozen (dielectric permittivity of permafrost: 4...5 (Sharma, 1997)) and the overlying unfrozen active layer (dielectric permittivity of wet sand/gravel: 10...20 (Sharma, 1997)), the frost table usually appears as a well detectable, continuous reflector in the radargram.
The major drawback of the typical common offset, single channel GPR measurements is, that the travel time measured to estimate the depth of a reflector strongly depends on the dielectric permittiviy of the ground as well.Hence, there are two unknown quantities, typically variable both in space and time, that determine the single measurement quantity.A common work-around is to calibrate dielectric permittivity and the related reflector depths with additional measurements, e.g., in boreholes (e.g., Arcone et al., 1998;Lunt et al., 2005), at excavates of soil profiles (e.g., Wollschläger and Roth, 2005) or common-midpoint (CMP) GPR measurements (e.g., Schwamborn et al., 2008).All this is practically feasible at only a few locations.However, due to the heterogeneous texture of natural soils, soil moisture content and consequently dielectric permittivity typically vary horizontally over rather short distances (e.g., Greaves et al., 1996).This results in distortions of reflector depth as displayed in a radargram if a constant signal propagation velocity for the complete transect is assumed (Neal, 2004).Hence, estimates of reflector depth which are calibrated at only a few points along the measurement transect are only able to provide a rough estimate of the true reflector structure.
An efficient way to obtain more accurate information about the true shape of the reflector and additionally infer average soil moisture content are continuous multi-offset measurements obtained using multi-channel GPR systems (Bradford, 2008;Gerhards et al., 2008).Here, the survey is done in the profiling mode with an array of several coupled GPR antennas which allows to measure the travel times from a number of different antenna separations at once.The acquired data can be evaluated in analogy to a standard CMP survey while the distance between the CMP data sets is determined by the pre-selected trace interval (e.g., each 0.1 m).Given suitable, continuous reflectors, this new measurement technique allows to infer high-resolution spatial information about the reflector structure and average soil moisture content without requiring much additional data for calibration.In permafrost regions, the technique allows to non-invasively map both, the thaw depth and the vertically integrated soil moisture content of the active layer which is important information, for instance for the interpretation and modelling of thawing and freezing processes.
The objective of our study is to apply multi-channel GPR at a continuous permafrost site in order to efficiently infer spatial variations in thaw depth and average volumetric soil moisture content of the active layer.In addition, we use this advanced geophysical measurement technique to obtain deeper understanding of field-scale active layer processes.Therefore, we first present an extended measurement technique with respect to the one described in Gerhards et al. (2008) by using an eight-channel GPR system and adding a detailed topographical survey which allows for an accurate three-dimensional representation of both, the ground surface and the frost table topography.Furthermore, we stabilize the inversion procedure by applying a small modification in the determination of the air wave travel times.Finally, we briefly discuss the observed patterns in thaw depth and soil moisture content with respect to the surface and soil properties of the investigated site.

Site description
The study site is located near the Qitedaban mountain pass (35 • 45 N, 79 • 26 E, 4950 m a.s.l.), in the northern Aksai Chin Region, Xinjiang Uigur Autonomous Region, W-China (Fig. 1).This western part of the Tibetan Plateau is mainly characterised by high-cold desert landscapes (Jin et al., 2007).Being situated in the rain shadow of the Karakorum and Kunlun mountain ranges, annual convective precipitation in the area is less than 50 mm which is occasionally enhanced by monsoon rainfall (Gasse et al., 1991).The ratio of evaporation to precipitation ranges between 20 and 50 (Gasse et al., 1991) indicating the extremely arid climate of the study area.Zhou et al. (1998) provide a lower permafrost limit of 4400 m to 4500 m a.s.l. for the W-Kunlun region.The site is situated within an area of continuous permafrost (Jin et al., 2000).
Measurements were conducted at the foot of an alluvial fan (Fig. 2).The soil surface at the upper part of the slope is almost bare.It is covered with a number of small flow channels which were dry during the survey.During an earlier measurement campaign in early October 2006 (Gerhards et al., 2008) some of the channels were water bearing but surface water infiltrated completely before reaching the foot of the slope.In a soil profile which was excavated during the 2007 measurement campaign a few hundred meters upslope of the measurement plot a groundwater table was observed at a depth of 0.76 m below the ground surface.The soil of the alluvial fan primarily consists of sand and gravel.Several finer textured areas can be identified by salt precipitations appearing at the soil surface.They indicate the higher water retention capacity of the fine grained soil which makes the soil water available for evaporation close to the ground surface.At the foot of the alluvial fan, on both sides of the Xinjiang-Tibet Highway, small vegetated areas occur.In the valley bottom, a dry river bed extends parallel to the gravel road.In some sections it is located upslope of the gravel road, then crosses the highway and continues in the downslope area of the road.
In contrast to the former study of Gerhards et al. ( 2008) who applied data from a single GPR transect crossing a small creek at this permafrost site to demonstrate the new multichannel GPR technique, this study aims at the spatial exploration of active layer thaw depth and soil moisture content and the investigation of the relationship of these param- eters to soil surface properties and soil texture.Therefore, we selected an approximately 85 × 60 m 2 sized plot at the foot of the slope where the soil surface was partially bare and partially covered with sparse grass (Fig. 2).In addition, the area is traversed by the Xingjiang-Tibet Highway, which represents another textural unit of our study site.The soil beneath the bare surface basically consists of medium to coarse grained sand and gravel.Below the vegetated area it is characterised by rather homogeneous fine to medium grained sand which is covered by a 10 to 15 cm thick layer of even finer, presumably wind-blown sand that has been accumulated by the sparse grass cover.Within the vegetated area, precipitations of salt were observed on the soil surface.The surface salt indicates the high amount of evaporation from the fine textured soil in this area.The sand accumulations within the vegetated patch produce a rough, slightly elevated surface topography with a maximum height of about 0.5 m (Fig. 2).

Materials and methods
Multi-channel GPR measurements were performed as a modification of the method presented by Gerhards et al. (2008).In the following, we give an overview about the new equipment together with the central equations of the multi-channel evaluation of Gerhards et al. (2008) and the modification of the evaluation procedure as adapted for this study.

Materials
Measurements were conducted using a multi-channel GPR array (Fig. 3a, b) that consisted of three standard 200 MHz shielded antenna systems (TR 200 K2) in combination with a DAD K2-MCH control unit and the distribution box (4CH), all manufactured by IDS (Ingegneria dei Sistemi S.p.A., Italy).The antenna boxes were coupled in a row with a www.the-cryosphere.net/4/269/2010/The Cryosphere, 4, 269-283, 2010 fixed distance of 1.01 m between the two front antennas and a distance of 2.015 m separating both rear antennas.As in a standard common-offset survey, the whole antenna array was pulled along the transect while signals were triggered by the survey wheel which was mounted at the back of the rear antenna.
The employed multi-channel GPR system allows to acquire data from a total of nine transmitter-receiver combinations, in the following referred to as "channels" (Fig. 3b), while the applied measurement software is capable to simultaneously record signals arriving from eight channels.For our survey we chose to measure radargrams with the following antenna separations: 2 × 0.19 m (T2-R2; T3-R3), 0.82 m (T2-R1), 1.2 m (T1-R2), 1.83 m (T3-R2), 2.21 m (T2-R3), 2.84 m (T3-R1), and 3.22 m (T1-R3).Compared to the study of Gerhards et al. ( 2008) which was performed with a 4-channel setup, the availability of eight channels is expected to stabilise the multi-channel evaluation since it is less sensitive to outliers.
As a new part of the measurement technique, the position of the GPR array and surface topography were recorded using a TCRA1102 tachymeter (Leica GeoSystems, Germany).The survey lines were mapped by the tachymeter by automatically tracking a 360 • prism mounted on top of the middle antenna (Fig. 3a, b).These additional data allow us to effi-ciently create maps of the depth and absolute topography of the frost table and the spatial distribution of active layer soil moisture content throughout the measurement area.

Methods
Multi-channel GPR measurements were acquired on 31 August 2007 when the annual thaw depth of the active layer was presumed to be close to its deepest position.The site was explored by measuring eight parallel lines (Fig. 2) starting upslope of the vegetated area with one profile covering a complete transect of bare soil.The following three lines started and finished on bare soil while crossing the vegetated area in-between.The next profile was measured in the roadside ditch of the Xinjiang-Tibet Highway and was followed by a line located directly on the gravel road.Finally, two more transects on the opposite side of the highway were explored.Here, in the valley bottom, the soil surface was again partially covered by sparse vegetation.
Radargrams were recorded using a time window of 100 ns, with 10 samples per ns and 24 stacks per trace.For the multi-channel evaluation, measured data from the different channels first have to be resorted such that they share the same common measurement point and reflections occur from a similar area of the reflector (Fig. 3c).Hence, in case of good data quality, ideally data recorded with seven different antenna separations (by design, two of the eight available channels have the same antenna separation, cf.Sect.3.1) are available for the inversion.In this investigation, we used the position of the 360 • prism of the laser tachymeter (Fig. 3b) as reference position for the relocation procedure.Since the survey is run in the normal profiling mode, the separation between the CMP data sets is equal to the spatial trace increment which was set to 0.05 m for this survey.By this moving CMP technique, we are able to efficiently estimate relative dielectric permittivity, reflector depth and average soil moisture content for each individual position along the radargram.
Data processing was minimal and consisted of (i) application of a dewow filter in order to remove low frequency noise, (ii) clipping of all radargrams to a total time window of 80 ns, and (iii) adjustment of the color density for each radargram individually for visualisation purposes.No amplification was applied to the data set.Next, the travel times of the frozen/unfrozen boundary reflection were picked in each radargram and air wave travel times were determined as a reference for the calculation of the absolute travel times of the reflected waves (Fig. 4).
Based on a ray approach, the multi-channel evaluation method (Gerhards et al., 2008) assumes that the two-way travel time t of the electromagnetic signal to a reflector at a position x in the surroundings of the reference position x 0 and antenna separation a is given as where ε c is the relative dielectric permittivity of the soil, c 0 is the speed of light in vaccuum (0.3 mns −1 ), α is the reflector's inclination angle, and d its depth at x 0 .Relative dielectric permittivity of the soil, depth and inclination angle of the reflector for every position along the radargram are estimated from the absolute travel times of the signals measured with all available transmitter-receiver combinations.This is done by minimizing the objective function where b={ε c ,d,α} is the parameter vector, t refl and t model are the measured and modelled reflected wave travel times for N measurements around x 0 obtained from K antenna separations, x n (n = 1, ..., N) are measurement points around x 0 , and a k (k = 1, ..., K) are the antenna separations.If, due to low data quality, picking of reflections along certain sections of single radargrams is not possible, the evaluation program automatically conducts the inversion with all data which are available for the pre-defined distance interval around x 0 which was set to 0.5 m for this survey.With respect to the trace interval of 0.05 m and ideally data from eight channels, for this setup up to 72 data points are available to minimise the objective function for one single position along the radargram.Air wave travel times are required since, in most GPR surveys, the absolute start-of-trace signal (time zero) is not available.Instead, the travel time of the air wave, with its known velocity, is used as a reference.Since, with the system used here, amplitudes of air wave signals measured with long antenna separations are often very low, picking of the accurate air wave travel time is difficult or even impossible.To overcome this difficulty we used the air wave adaptation procedure of Gerhards et al. (2008) which is based on the multichannel evaluation.In order to confine the inverse adaptation of air wave travel times to a reasonable range, as a modification to the work of Gerhards et al. (2008), air wave travel times of the channels with 0.19 m antenna separation (T2-R2 and T3-R3), which were well identifiable in the radargrams, were first picked and then fixed to a constant average travel www.the-cryosphere.net/4/269/2010/The Cryosphere, 4, 269-283, 2010 time which was calculated for both channels individually.Subsequently, air wave travel times of all other channels were adjusted by using the air wave adaptation algorithm.This procedure was done for each of the eight transects individually in order to account for a drift of the antenna electronics.
For data evaluation, we use constant air wave travel times for each individual radargram.This time-zero determination procedure may cause a certain bias on the estimated absolute relative dielectric permittivities and the related reflector depths and average soil moisture contents.However, relative changes in these parameters within the data set will still be recognisable after this correction (Gerhards et al., 2008).The average volumetric water content θ of the unfrozen active layer was calculated from the estimated relative dielectric permittivities according to Roth et al. (1990) for each position along the radargram by using the CRIM (Complex Refractive Index Method) formula: where ε c is the composite dielectric permittivity of the soil, composed of the volume fractions of solid matrix ε s , water ε w , and air ε a , and φ is porosity.Re-arrangement of Eq. ( 3) for calculating volumetric soil moisture content leads to The dielectric permittivity of water was set to a constant value of 86.1 corresponding to a soil temperature of 5 • C (Kaatze, 1989).For the dielectric permittivity of the solid matrix we initially assumed a constant value of ε s = 5, and set the porosity to an estimated fixed value of φ = 0.4.In addition, the uncertainty of these estimates was further evaluated by running the analysis with porosities and dielectric permittivities of the soild matrix ranging from 0.3 to 0.5 and 4 to 6, respectively.These values are assumed to be reasonable estimates for the outer limits of φ and ε s for the soil at the investigated site.
For plotting of the radargrams we considered topography by shifting each trace vertically by 0.1 ns per meter height difference.For spatial visualisation of the data acquired from all eight transects, contour plots of surface and reflector topography, depth of the reflector below the ground surface and average volumetric soil moisture content were generated by a simple bi-linear interpolation of the values measured along the different transects.

Transect crossing bare soil and vegetation
Radargrams for all eight channels derived from a multichannel measurement along one selected transect are shown in Fig. 4. The data were recorded along a measurement line upslope of the roadside ditch partly traversing bare soil and partly the vegetated area (cf.Fig. 2).In the middle of the transect a non-vegetated section within the vegetated part of the site was crossed.
All radargrams show one dominant reflection from the frost table which is well detectable on each channel.In general, the shape of the reflector indicated by its two-way travel time follows the surface topography, with travel times slightly longer below the bare soil sections than below the vegetated sections.This initially indicates a deeper frost table and/or a higher soil moisture content in this area compared to the vegetated patch.Below the vegetated area, at long antenna separations, the reflection becomes weak or even disappears.We presume that this is mainly due to the weaker antenna-ground-coupling at the rough ground surface in combination with high geometrical spreading at long antenna separations.Fig. 5. Surface topography, reflector depth (top graph), average volumetric soil moisture content and relative dielectric permittivity (bottom graph) inferred for the example transect shown in Fig. 4; red: raw data, dark blue: data averaged using a linearly weighted running mean filter over a distance of 2 m, cyan: uncertainty limits for the averaged soil moisture contents (dark blue) using the CRIM formula and assuming variations in porosity and dielectric permittivity of the solid matrix of 0.3...0.5 and 4...6, respectively, green: reflector depth obtained from a hypothetical single-channel GPR measurement assuming a constant dielectric permittivity of ε c = 15.
The relative permittivity axis refers to the averaged moisture contents printed in dark blue assuming φ = 0.4 and ε s = 5.
Figure 5 shows the topographically corrected reflector depth, relative dielectric permittivity and average volumetric soil moisture content of the unfrozen active layer as derived from the multi-channel analysis for the survey line displayed in Fig. 4. Both, the topography of the reflector below the ground surface and also the dielectric permittivity estimated from the multi-channel inversion are rather noisy (Fig. 5, red lines).We explain the origin of these small-scale variations by small changes in antenna separations, instrument noise, and the picking accuracy which is on the order of 0.2 ns.All these factors induce small-scale statistical variations in the measured travel times of the various channels which are expressed as high-frequency fluctuations in the raw data and are carried further into the multi-channel evaluation.For the interpretation of the variations of reflector depth and dielectric permittivity, we smoothed both data sets with a linearly weighted running mean filter calculated with 2 m support (Fig. 5, dark blue lines) and concentrate in the following on the discussion of major variations in reflector depth, dielectric permittivity and average soil moisture content, respectively, which are large compared to the statistical noise.In addition, a systematic error may be induced by the air wave adaptation method which may lead to a small shift in absolute travel times of single channels but will not affect significantly the inferred relative changes of reflector depth and soil moisture content along single transects (see also Sect.Fig. 6.Number of measurement points available for each position along the transect (top) and errors in reflector depth (middle) and relative dielectric permittivity (bottom) resulting from the inverse parameter estimation.The number of data points is determined by the number of channels where picks of the frost table reflection are available (Fig. 4).
wave determination and adaptation is given by Gerhards et al. (2008).
As it was already observed from the travel times displayed in the radargrams (Fig. 4), the inverted frost table depth (dark blue line) generally follows the surface topography but differences between bare soil and vegetated sections along the transect appear not as obvious as indicated by the signal travel times.Nevertheless, it is still more shallow below the vegetated sections compared to the bare soil area which becomes more obvious when looking at the spatial data (Fig. 7) which are described in Sect.4.2.This thaw depth variability can be related to strong changes in dielectric permittivity along the transect which can primarily be referred to lateral changes in average soil moisture content (Fig. 5).On the first approximately 10 m of the radargram, below a bare soil section, moisture contents are rather high and then decrease sharply when entering the vegetated area.Below the section of the profile where vegetation and accumulations of fine sand at the soil surface are missing (∼ 30...40 m), the dielectric permittivity, hence average soil moisture content, rises abruptly and drops again as the following vegetated patch is crossed.Within the adjacent bare soil section, average volumetric soil moisture content and dielectric permittivity are again high and remain at these values for the final section of the radargram.
To visualise the impact of the varying dielectric permittivity on the estimations of the depth of the frost table, re-flector depth was calculated for a constant dielectric permittivity of ε c = 15 as inferred from the multi-channel inversion for most sections of the bare soil areas.It corresponds to a constant signal velocity of 0.08 m ns −1 (Fig. 5, green line).Compared to the estimates with variable dielectric permittivity, this leads to a frost table that is 0.5 m shallower beneath the vegetated patch.
The analysis clearly shows, that differences in observed travel times (Fig. 4) are not exclusively related to changes in reflector depth.They strongly depend on lateral variations in average soil moisture content which differs significantly between the bare and vegetated sections.Short travel times beneath the vegetated area are strongly related to low dielectric permittivities which can be referred to low volumetric soil moisture contents in this area.Longer travel times in the bare soil sections are strongly influenced by a higher soil moisture content in this region.Assuming one single propagation velocity for the GPR signal for the complete measurement line (Fig. 5, green line) leads to severe deviations of the estimated reflector depth from the real situation.We comment that similar observations with high relative dielectric permittivities, hence average volumetric soil moisture contents below the bare soil areas and low values below the vegetated area were made for the other radargrams crossing the vegetated area.They are visualised in the spatial representation of the data which is discussed in Sect.4.2.

Uncertainty estimation
In the following, we will use the results from the transect presented in Sect.4.1 for discussing the uncertainty of the different steps of the multi-channel evaluation procedure.Since picking of the reflection from the frost table was not possible for all channels along the complete measurement line (Fig. 4) a variable number of measurement points along the transect is available for the inverse parameter estimation.We will first discuss the uncertainty of the inversion procedure with respect to the estimated reflector depths and dielectric permittivities taking also into account the available number of measurement points.Afterwards, we will evaluate the uncertainty in the soil moisture content estimates using the CRIM formula.

Errors resulting from the inverse parameter estimation
Errors (Fig. 6) were calculated as confidence limits on the estimated model parameters (Press et al., 2002) assuming that the problem around the parameter optimum (b opt ) can be linearised and that the measurement errors are normally distributed.The error b i alias the confidence interval for parameter b i was calculated by where OF(b opt ) is the objective function at its minimum value at the end of the inverse parameter estimation, and J is the Jacobian matrix.The value S − 3 represents the number of the used data points reduced by the number of the parameters which is 3 in our case.The value below the square root is a diagonal entry of the covariance matrix.The factor √ 11.3 sets the confidence level to 99% with respect to the number of parameters.
The errors for reflector depth and dielectric permittivity are low remaining below 0.08 m and 0.3, respectively along most parts of the transect.In addition, the uncertainty reduces with more data points being available for the inverse parameter estimation.In cases where all eight channels provide data (high number of measurement points), the errors are lowest.In contrast, for sections, where only few data points are available, as in the vegetated sections at around 30 m and between 40 m and 45 m, errors are largest.This www.the-cryosphere.net/4/269/2010/The Cryosphere, 4, 269-283, 2010 U. Wollschläger et al.: Multi-channel GPR for exploring thaw depth and moisture content of the active layer result initially supports the application of a multi-channel GPR system with a large number of different antenna separations since the large number of available measurement points stabilises the parameter estimation.However, in our example, the large error along sections with a small number of available data points can additionally be referred to the fact that, due to site-specific conditions, only data from short antenna separations are available resulting in travel times which are not too much different from each other and which increases the uncertainty of the parameter estimation.
In cases where only a system with a small number of channels is available and subsurface conditions allow large antenna separations, the multi-channel setup should be chosen such that the travel times measured on the various channels cover a broad range.Finally we have to state that the errors given in Fig. 6 only represent the uncertainty resulting from the inverse parameter estimation.Further uncertainty which e.g. may result from changes in antenna separations during the measurement is not considered here.

CRIM formula
While the estimation of reflector depth and dielectric permittivity only depends on the measured travel times, the calculation of volumetric soil moisture contents using the CRIM formula introduces further uncertainties.These result from the constant estimates for porosity, dielectric permittivity of the solid matrix and soil temperature which influences the value of the dielectric permittivity of water.A detailed error propagation analysis of the single components of the CRIM formula is provided by Roth et al. (1990).In order to visualise the range of uncertainty for the measurements of this study, average soil moisture contents were calculated using extreme values for porosity (0.3...0.5) and dielectric permittivity of the solid matrix (4...6).The resulting outer limits for soil moisture content are included in Fig. 5 (cyan lines).Variations of porosity and dielectric permittivity of the solid matrix along the survey track may lead to a reduction of the differences in the inferred soil moisture contents between bare soil sections and the vegetated area.However, even when considering extreme values, soil moisture content beneath the vegetated plot will remain lower than in the bare soil sections.A complete error analysis for the application of the CRIM formula would demand information on correlations between the various parameters as well as the temperature field which enters through ε w (T ).With the available data, this is out of reach.However, since some of the correlations lead to a reduction of the total error, the given bounds are considered as upper limits.We further state that the multichannel GPR method as it has been applied in this study is a non-invasive large-scale method which, by nature, is not as accurate as point measurements like gravimetric sampling or time-domain reflectometry.It includes a certain level of uncertainty where detailed spatial information of specific values such as porosity, dielectric permittivity of the solid ma-trix and the vertical soil temperature distribution are usually unknown and have to be estimated.These estimates could be specified if additional measured values from small-scale investigations such as soil samples or temperature profiles were available.

Spatial data
Contour plots of surface and reflector topography, thaw depth below ground surface, average volumetric soil moisture content and, in addition, total soil moisture content of the unfrozen active layer are shown in Fig. 7.The total soil moisture content of the active layer can be calculated directly from the multi-channel GPR measurements by multiplying the measured thaw depth and the average soil moisture content.It is an important measure in permafrost research as it provides direct information about how much water is stored in the active layer which can be applied to improve hydrological calculations.Furthermore, it gives a rough estimate of the latent heat content of the active layer or, in other words, the amount of energy which has to be released for freezing the active layer if no water is being added or removed (e.g., by rain or groundwater inflow) during freezing (Westermann et al., 2010).Like for the transect described in Sect.4.1, estimated frost table depths and water contents of all transects were smoothed with a linearly weighted running mean filter calculated with 2 m support.Volumetric water contents were calculated assuming a porosity of 0.4 and a dielectric permittivity of the solid matrix of 5.
Generally, the measured thaw depth and average soil moisture content of the active layer show a high spatial variability throughout the measurement plot which can be related to the observed differences in soil textural and surface properties.The large-scale topography of the reflecting frost table generally follows the surface topography (Fig. 7).Neglecting the slope and considering only the thickness of the thawed active layer, one can clearly observe a shallower frost table beneath the vegetated area.Here, where the ground surface is also higher than in the surrounding region, the frost table builds a shallow mound below the vegetation.A more shallow frost table was also observed along the roadside ditch which results from the lower topographic height compared to the adjacent measurement lines.Here, the frost table apparently is not able to penetrate deeper into the ground.Upslope of the road, in the more coarse textured, bare soil sections surrounding the vegetated area, the depth of the frost table is about 0.3 m to 0.7 m deeper than underneath the vegetated patch.Here, the estimated absolute late summer thaw depth of the active layer ranges between 1.5 m and 2.0 m which is in good agreement with the depths reported by Gerhards et al. ( 2008) one year earlier further upslope of the alluvial fan.For the transect along the gravel road and the two profiles further downslope, an intermediate thaw depth of about 1.5 m was detected.A similar small-scale variablility in active layer The Cryosphere, 4, 269-283, 2010 www.the-cryosphere.net/4/269/2010/thaw depth as observed for this site which could be related to microclimate, hydrology, soil moisture content, substrate properties, vegetation and topography has e.g.been observed by Nelson et al. (1998) and Hinkel and Nelson (2003) at several sites in Alaska.Measured average soil moisture contents underneath the vegetated patch are low with values of less than 0.2.They rise beneath the directly surrounding bare soil area up to values between 0.25 and 0.3.With greater distance to the vegetated patch, average soil moisture contents decrease again in the upslope area of the bare soil sections down to values of less than 0.15.They agree well with the soil moisture contents measured by Gerhards et al. (2008) in the distant areas of the investigated flow channel.Average soil moisture contents larger than 0.25 were detected below the gravel road, the adjacent roadside ditch, and along the profiles measured in the valley bottom downslope of the road.
The calculated total soil moisture contents of the thawed active layer show a similar spatial pattern as thaw depth and average soil moisture content.The total amount of moisture stored in the active layer is high in the valley bottom along the gravel road with values of larger than 0.4 m.High total moisture contents were also found below the bare soil sections upslope of the gravel road in the direct surroundings of the vegetated patch bordering the shallow mound in the frost table.With greater distance to the vegetated area, the total amount of moisture stored in the profile decreases down to values of less than 0.3 m.The lowest total amount of water was calculated for the area beneath the vegetation with values of mostly less than 0.25 m.

Discussion of the observed frost table depths and average soil moisture contents
The observed late-summer patterns of thaw depth and soil moisture content are the result of various factors and short and long-term processes which determine the thawing rates of the active layer throughout the site.Since our data only represent one point in time and additional continuous observations like meteorological variables, soil temperature or soil moisture content from this site are not available, we can only hypothesise about the processes that led to the observed thaw pattern.One patricular observation demands explanation, however: The measured distribution of soil moisture content throughout the site is somewhat unexpected with low moisture contents in fine textured soils and high moisture contents in coarse textured areas which usually should be reverse (e.g.Jury et al. 1991, Ch. 2;Hillel, 1998, Ch. 18).Hence, in the following, we will attempt to find some possible explanations for the observed patterns of thaw depth and soil moisture content.
From the factors controlling the thaw depth of the active layer (cf.Sect. 1) we presume soil texture, soil moisture content and topography to be most important for our site.Additional subordinate factors may be vegetation and albedo.Snow cover which can have an important effect on the areal thawing process and consequently on thaw depth (e.g., Smith, 1975) should not be important at the Qitedaban site since, first, the site is located in an extremely arid environment with very low precipitation (Gasse et al., 1991), hence, the total amount of snowfall is low.Second, due to its low latitude and high atmospheric transparency at high elevations, the Tibetan Plateau receives a large amount of solar radiation (e.g., Wang and French, 1994;Ma et al., 2009) which hampers the formation of a persistant snow cover and which, at the same time, supplies a large amount of energy for thawing of the active layer.
For the following discussion, the site is divided into three major units which exhibit significant differences in thaw depth and soil moisture content (Fig. 7): (i) the bare soil area upslope of the Xinjiang-Tibet Highway which is characterised by a relatively deep frost table and low to intermediate average soil moisture contents, (ii) the vegetated, fine sand area with shallowest thaw depths and low average soil moisture contents, and (iii) the gravel road and its downslope area with high average soil moisture contents and intermediate thaw depths.We will first discuss the observed soil moisture content pattern since soil moisture content influences the thermal properties of the soil and consequently affects thawing rates.Afterwards we will use this information for the discussion of the observed frost table depths.

Soil moisture content
We first recall that multi-channel GPR yields a vertically averaged moisture content for the complete active layer, ranging from the ground surface down to the frost table.Typically, water is not distributed uniformly in this layer, and the measured value includes contributions from groundwater (fully saturated) as well as from the unsaturated zone.The presence of groundwater at the site has been identified in the soil profile further upslope (Sect.2) and is also expected for our measurement plot in the valley bottom.We expect groundwater to be the principal factor determining the observed soil moisture content pattern: By nature, a groundwater table is rather flat while the ice table can be rough with mounds and troughs.Hence, areas with a deeper frost table will contribute a larger amount of groundwater to the measured average soil moisture content which will thus be larger than at a location with a shallower frost table.This is illustrated schematically in Fig. 8 where the groundwater body at the site can be expected to be considerably thicker in the bare soil sections and beneath the gravel road compared to the vegetated areas with shallower frost table.Similarly, we explain the enhanced average and total soil moisture contents in the bare soil regions upslope of the vegetated plot by pondig of groundwater in front of the shallow mound in the frost table.
With respect to the vadose zone soil moisture content, we generally can expect a higher soil moisture content and a  Evaporation is expected to be higher in the fine grained and compacted soils at the site due to the higher capillary rise of water from the groundwater table towards the soil surface.
higher capillary rise of water towards the soil surface to occur within the fine textured areas (e.g., Jury et al. 1991, Ch. 2;Hillel, 1998, Ch. 18).In general, this initially leads to a wetter soil in the vegetated area compared to the bare regions.There, water it is available (i) for evaporation and (ii) for root water uptake.This situation is confirmed by the presence of plants and the appearence of salt precipitations at the soil surface which indicate the higher amount of evaporation in this part of the site compared to the bare soil areas.Both, evaporation and root water uptake finally lead to a reduction of soil moisture contents in the vegetated area andin combination with the thin groundwater body on top of the shallow frost table -may cause the low average soil moisture contents measured with multi-channel GPR in this region.Small average pore sizes and a high capillary rise of soil moisture from the groundwater table towards the soil surface can also be expected for the compacted sediments of the gravel road and its downslope area.This, together with the location in the valley bottom, where groundwater will accumulate, may explain the high average soil moisture contents in this section of the experimental site.

Thaw depth
The observed thaw depth can be related to various surface and subsurfeace factors and processes occuring within the major units of the site which influence the energy transfer from the atmosphere to the frost table to different extents.
In the vegetated area, there are many factors that may lead to the observed shallow thaw depths.Grassy surfaces, as found as rather sparse cover at the experimental site, only have a very small shading effect and influence mean ground temperatures insignificantly (Yershov, 1998, pp. 362-369).We consider the effect of enhanced evaporation of soil moisture being present close to the soil surface in this fine textured area (cf.Sect.4.3.1)and, to a small extent, also transpiration from the grass to be more important for the measured shallow thaw depths.During times when water is available close to the ground surface, evaporation leads to a cooling of the soil surface and thus to lower subsurface temperatures (Williams and Smith, 2008, pp. 68-70) and, consequently promotes a shallow frost table.Another factor which may favour shallow frost table depths in the vegetated area may be the enhanced albedo caused by the preciptations of white surface salts (cf.Fig. 2).
Evaporative cooling is obviously much lower in the bare soil units of the site located upslope of the gravel road since salt precipitations at the soil surface are missing.Here, capillary rise is lower due to the coarser soil texture which limits the evaporation from the ground surface.Consequently, higher surface and active layer temperatures can be expected, leading to deeper frost tables in these regions.
The intermediate thaw depth together with the high water content observed along the gravel road and its downslope region may result from different processes whose contributions cannot be quantified from a one-time multi-channel GPR measurement alone.These processes include evaporative cooling which is again indicated by salt precipitations on the ground surface (cf.Fig. 2) which would lead to lower active layer temperatures and shallower frost table depths compared to the bare soil regions further uplsope.The observation that thaw depth is deeper here compared to the vegetated patch may be explained by the runoff of ground water and surface water which is concentrated in the valley bottom.Deep active layers in places of flowing water are reported from other permafrost regions, e.g. on Svalbard (Brown et al., 2000).

Summary and conclusions
Multi-channel ground-penetrating radar was applied at a continuous permafrost site on the Tibetan Plateau in order to spatially explore late-summer thaw depth, average and total moisture content of the active layer.For this, an improved version of the method presented in Gerhards et al. (2008) was used.The additional application of a laser tachymeter provided detailed spatial information about the topography of the soil surface and the frost table and assisted in the interpretation of the observed water content pattern.Even if conducted at a rather small example site with strongly varying surface and soil textural conditions, this investigation demonstrates the potential of multi-channel ground penetrating radar to efficiently and non-invasively map thaw depth, average and total soil moisture content of the active layer with high spatial resolution.This method is reasonably accurate and orders of magnitude more efficient, in terms of measuring time, than traditional methods like time domain reflectometry or gravimetric sampling.
From our measurements, we deduced hypotheses for the observed late-season thaw depths and soil moisture contents at the site.Surface and subsurface factors (soil texture and soil moisture content, topography and potentially and to a lesser extent vegetation and albedo) influence the shape of the frost table.This influence is far from trivial, however, since surface and subsurface factors interact with processes in a complicated manner.
The multi-channel GPR technique efficiently covers the intermediate scales between traditional point measurements and space-based remote sensing.The method yields spatial information of subsurface active layer parameters and mesoscale heterogeneity that is hard to obtain otherwise.Hence, it can e.g.be used as ground-truth information to calibrate and interpret remote sensing data and permafrost models.In a modelling context, especially the average soil moisture content for the complete active layer is of relevance since it significantly influences the soil's thermal properties and cannot be inferred in a comparable fashion from any other field method at the scale covered by multi-channel GPR.More importantly, the technique opens a way for spatially resolved, detailed process studies, of ecosystem dynamics or hydrology in permafrost environments.In electromagnetically lowloss media such as sand or peat soils with reasonably deep active layers which allow for a clear identification of the frost table in the radargram, multi-channel GPR may be considered as an alternative, non-invasive method for active layer monitoring to the standard frost-probe measurements suggested by the CALM protocol (Brown et al., 2000;Smith and Brown, 2009).
In order to quantify hillslope hydrological processes such as runoff or subsurface storage, knowledge of the temporal development of the frost table depth throughout the thawing season is required (Hinzman et al., 1998;Wright et al., 2009).The spatial maps of active layer thaw depth and soil moisture content (Fig. 7) e.g.indicate that the major groundwater discharge on top of the frost table at the investigated site predominantly occurs in the bare soil regions while the mound in the frost table below the vegetated part will more of less be passed by by the groundwater.Due to its fast and efficient applicability, time series of multi-channel GPR measurements may in future provide important information for a more detailed process understanding.

Fig. 2 .
Fig. 2. Photographs of the study area: GPR measurements were acquired at the foot of an alluvial fan, partly across bare soil, a small vegetated area and the roadbed of the Xinjiang-Tibet Highway (red box; arrows indicate the direction of GPR lines as shown in Fig. 4; thin dashed line: transect discussed in Sect.4.1, Figs. 4 to 6, transitions between vegetated and non-vegetated areas are marked separately).A detailed photograph of the vegetated area and the adjacent bare soil is provided on the left photograph, the black arrow indicates the location and viewing direction of the photograph at the left.
Fig. 3. (a) Photograph of an 8 channel setup.(b) Sketch of the 9 possible transmitter-receiver combinations.For analysis, radargrams are shifted to a common reflection point, here to the position where the 360 • prism of the tachymeter is mounted on the middle antenna.(c) CMP pathways after lateral shift of radargrams.

Fig. 4 .
Fig. 4.Radargrams from all 8 channels showing a measurement which traversed both, bare soil and the vegetated area (Fig.2).Data are plotted including topographical correction.Vegetation and bare soil sections are indicated in the topmost radargrams.Air wave and ground wave travel times are drawn in light blue, antenna separations (AS) are given at the lower right of each radargram.Below the vegetated areas the reflection from the frost table is not clearly visible at long antenna separations.
Figure5shows the topographically corrected reflector depth, relative dielectric permittivity and average volumetric soil moisture content of the unfrozen active layer as derived from the multi-channel analysis for the survey line displayed in Fig.4.Both, the topography of the reflector below the ground surface and also the dielectric permittivity estimated from the multi-channel inversion are rather noisy (Fig.5, red lines).We explain the origin of these small-scale variations by small changes in antenna separations, instrument noise, and the picking accuracy which is on the order of 0.2 ns.All these factors induce small-scale statistical variations in the measured travel times of the various channels which are expressed as high-frequency fluctuations in the raw data and are carried further into the multi-channel evaluation.For the interpretation of the variations of reflector depth and dielectric permittivity, we smoothed both data sets with a linearly weighted running mean filter calculated with 2 m support (Fig.5, dark blue lines) and concentrate in the following on the discussion of major variations in reflector depth, dielectric permittivity and average soil moisture content, respectively, which are large compared to the statistical noise.In addition, a systematic error may be induced by the air wave adaptation method which may lead to a small shift in absolute travel times of single channels but will not affect significantly the inferred relative changes of reflector depth and soil moisture content along single transects (see alsoSect.3.2).A detailed discussion on errors for the air

TheFig. 7 .
Fig. 7. Surface and reflector topography, thaw depth, average soil moisture content and total soil moisture content calculated for the thawed active layer within the measurement area.Black lines indicate positions of GPR profiles.For better orientation, the positions of the gravel road, vegetated and bare soil area, and the line discussed in Figs. 4 and 5 are indicated separately.The dashed line roughly marks the vegetated area upslope of the gravel road.

Fig. 8 .
Fig.8.Schematic illustration of the interpretation of observed thaw depth and groundwater level beneath the major units of the investigated site.The measured average soil moisture content is composed of both, groundwater and water being present in the unsaturated zone.Evaporation is expected to be higher in the fine grained and compacted soils at the site due to the higher capillary rise of water from the groundwater table towards the soil surface.