3.1 Equipment

We used a Geophysical Survey Systems Inc. (GSSI) model SIR 10 B control unit and the GSSI 1.2 W peak power 400 MHz model 5103 antenna transducer unit. The 4.1 ns duration of our mainly 3/2-cycle pulse gives about 0.41 m layer resolution in firn of refractive index n = 1.5 (density, ρ = 600 kg m⁻³). We pushed the 5103 unit ahead of a Tucker Sno-Cat® (Fig. 2); in subsequent years we dragged it in a sled (Reference Arcone, Spikes, Hamilton and MayewskiArcone and others, 2004, in press). We recorded differential GPS position every 30 s, corresponding to about 60 m at typical traverse speeds, and used the data to normalize our profile distances between 1 km segments.

Fig. 2. The 1999 ITASE expedition Tucker Sno-Cat® pushing a tire inner tube within which we placed the antenna unit (inset, lower right). The inset at upper left shows the direct coupling between transmitter and receiver antennas, and a reflection from a metal sheet buried in snow.

We recorded our radar data at 16-bit samples per trace and with time-range gain. Thus our maximum dynamic range in a profile is 96 dB, but the useful range considering noise was probably <90 dB. We recorded at a time range of 600 ns and at 2048 samples per trace. We used a running stack of 32 and a trace acquisition rate of 24 s⁻¹ so that the effective recording rate was about 1 trace per 2.7 m at the traverse speed of approximately 2 ms⁻¹. Profiles recorded in later years extended the time range to 1500 ns (Reference Arcone, Spikes, Hamilton and MayewskiArcone and others 2004, in press), but at a greatly compromised acquisition rate of 1 trace per 15 m, thus sacrificing some data quality.

3.2. Data processing

We alleviated direct-current gain offsets and electronic noise with a wide-band filter, removed constant time-delay nearsurface clutter with a horizontal filter, removed the timevariable gain used during recording, and then applied an inverse range-dependent gain function to compensate for geometric spherical beam spreading losses. The progressive accumulation of two-way electric field transmission losses through each interface in firn is theoretically insignificant.

We used the echo delay formula,

where c = 0.3 m ns⁻¹, to transform echo time delay, t, measured in nanoseconds (ns), into depth, d, in meters. We used an empirically determined relationship (Reference Kovacs, Gow and MoreyKovacs and others, 1996), verified with earlier datasets (Reference CummingCumming, 1952), to transform density, ρ, into dielectric permittivity, e. Equation (1) was applied meter by meter to the densities obtained from our core profiles to produce a slightly nonlinear depth scale that varied insignificantly between core sites. The quantity

Our spatial processing of the profiles consisted of additional stacking (beyond that used during recording), which greatly improved the signal-to-noise quality, and distance normalization using our GPS data to provide an equal number of traces per kilometer. As discussed earlier, we do not correct for elevation, so that the variations in horizon depths relative to the surface are clear. With elevation corrections, it would be difficult to see these depth variations along the stratigraphic dips because they compare with the variations in elevation. The extremely shallow stratigraphic dips precluded any need to migrate our profiles. However, the stacking weakened reflections along steeper dips, as will be discussed.

3.3. Antenna considerations and reflection dynamics

The complexities of the antenna directivity are not important, because of the gentle slope of the layers. We calculate that our reflections originated within an approximately coneshaped volume from beneath the antenna, given the gentle dip of the stratigraphy and the width of the first Fresnel zone. for our situation, the Fresnel zone is the circular area over which distance from the antenna does not vary by more than about λ/8 (0.06 m in firn of ε = 2.6) and has a diameter where d is depth beneath the antenna and λ is the in situ wavelength. At an average ε = 2.6, D is only about 5.2 m at 56 m depth. The recording rate and the Fresnel zone determine the theoretical distance over which the reflections are integrated. This distance varied from 4.3 m at a depth of 5 m to 7.8 m at a depth of 56 m.

4.1. General profile and accumulation rates

We show a 350:1 horizontally compressed version of the entire profile in Figure 3. An important variation feature that we discuss later is the general increase in the dip of the limbs with depth within any one fold. This feature is strongest at 70–130 and 150–177 km, which are the regions of higher ice velocity. A second important feature discussed later is the general failure of the loci of a stratigraphic fold to monotonically migrate down-ice (to the right in Fig. 3) with depth, at the rates expected from the ice speeds. A third feature is the significant change in depth of some horizons over long distances. The horizon marked by vertical up-pointing arrows varies from >60 m depth near 30 and 40 km, where it dips below the record, to 30 m depth at 177 km. This horizon is dated to the year 1734 (personal communication from D. Dixon, 2003). A second horizon, marked by smaller arrows and dated to 1875, shows a similar variable-depth character.

Fig. 3. The radar profile of the entire transect without elevation corrections. The vertical exaggeration is about 350 : 1. The upper set of small arrows indicates a horizon dated to 1875; the lower set of larger arrows indicates a horizon dated to 1734. Glacier flow is from northeast (left) to southwest (right).

Figure 4 shows the calibration data we used to convert the time delay of any particular horizon in Figure 3 into a cumulative water equivalent (w.e.) value in m a⁻¹. The time delay was then converted to depth using the meter-bymeter procedure discussed in section 3.2, and then to years from the dating. The BSC calibration in Figure 4a was used for the computation of the near-surface accumulation rates in Figure 5a, while all calibrations and an interpolation (Fig. 4b) were used for the deeper rates of Figure 5b.

Fig. 4. The cumulative w.e. accumulation as a function of radar round-trip time delay at the three core sites for the first 60 ns (about 6 m depth) (a) and for 200–580 ns (b). The former was used to convert the 1992–99 horizon time delays into equivalent water accumulation rates in Figure 5a; the latter was used for the accumulation rate calculations in Figure 5b.

Fig. 5. The ice-sheet surface slope and average w.e. accumulation rate per year for 1992–99 (a), and the average w.e. accumulation rate per year for two longer intervals (b). The arrows in (b) indicate sections where the 1734 horizon dipped below the profile record.

The accumulation rates plotted in Figure 5a are 8 year w.e. averages for 1992–99 (start of 1992 to the end of 1999, at the surface). We picked the time delays at every kilometer along the leading edge of a horizon half-cycle that we dated to the start of 1992 using the ice-core data from core site BSC (2.5 m depth; Reference Kreutz, Mayewski, Meeker, Twickler and WhitlowKreutz and others, 2000) and from core site 99-1 at 87 km (4.2 m depth; personal communication from D. Dixon, 2003). We assume that the near-surface density profiles are consistent along the entire transect because cores from 0 km (BSC) and 177 km (99-2) show similar near-surface rates of densification as a function of depth (Fig. 4a). Although density was not obtained for the first 4 m at core site 99-1, the 1992 horizon at BSC does track to about 4.2 m depth at this location.

We chose this shallow horizon so that each accumulation rate would correspond with the topography at that location, whereas deeper horizons include accumulation from some distance up-glacier. The horizon appears to bifurcate at two locations, but where it does, the correct continuation was guided by both the amplitude and general dip of neighboring horizons. Our accuracy in measuring time delay is estimated from the width of a wavelet half-cycle (about 1.3 ns; Fig. 2). Therefore, the error in our accumulation rates varies from a maximum of 6.8% for the minimum time delay of 14.7 ns (at 118 km) to a minimum of 2.5% for the maximum delay of 49.3 ns (77 km). The estimated average error caused by timing inaccuracies is then 4.7% for the entire profile. The accuracy in the dating is about 0.1 year, which increases the 8 year w.e. average error to about 6%.

We calculated along-profile slopes from the GPS elevation profile. The slope at each km represents an average over ±1 km from that point, and we did not smooth the data. The maxima in the plot are not necessarily the absolute maxima because our direction of travel relative to the wind was off by as much as 37°. Slope maxima indicate windward slopes facing up-glacier (to the left) to the northeast, and minima indicate leeward slopes, facing down-glacier (to the right) to the southwest. On average, the transect dipped 0.18° to the southwest. Consequently, leeward slopes are always negative, but some maximum slopes, on windward sides of hills, may still be negative.

The correlation between slope and accumulation rate is strong, regardless of slope value. In some cases, even small perturbations in slope, such as at 106 and 114–123 km, are in phase with corresponding perturbations in accumulation rate. At a few locations, such as at 12 and 136 km, they are out of phase by 180°. At other locations, such as at 40–50 and about 76.3 km, they are out of phase by about 1 km, which may result partially from the 2 km method of computing the slope.

Figure 5b shows plots of average accumulation rates calculated for 1875–1999 and 1734–1999. The rates are based on density values that occur at depths which correspond with time delays of 200–580 ns. We used the density profile at BSC for 0–22 km, an interpolation between the BSC and 99-1 profiles for 23–65 km, the profile at 99-1 for 66–173 km, and the 99-2 profile for 174–177 km. Core site 99-2 appears unique because the density reached 0.82 kg m⁻³ at only 38 m depth. The radar horizons distinctly change in amplitude below this depth (discussed later), and the elevation drops into a basin after about 173 km.

In Figure 5b there is considerable variation for the two horizons. The flat tops within the 1734 horizon represent sections where this horizon dipped below the profile depth. Between 5 and 70 km the 1875–1999 rates vary from about 0.12 to 0.19 m w.e. a⁻¹, and from 1992 to 1999 they have a similar range of about 0.15–0.21 m w.e. a⁻¹. The six prominent maxima for 1892–1999 between 70 and 125 km in Figure 5a also occur in Figure 5b at almost the same locations, but the 1875–1999 range of about 0.085–0.185 m w.e. a⁻¹ is lower than the 1992–99 range of 0.08–0.28 m w.e. a⁻¹. However, from 125 to 177 km the 1875–1999 rates vary from about 0.06 to 0.15, the 1734–1999 rates vary from 0.06 to 0.13 m w.e. a⁻¹, and the 1992–99 rates vary from about 0.09 to 0.135 m w.e. a⁻¹. Therefore, although the tracking of maxima and minima is generally comparable for all horizons to at least 300 years, the larger, long-term ranges of variation between 70 and 125 km appear to diminish with depth. This is reasonable because any rate determined at a single location is an average over an increasingly longer distance upstream as depth increases at that location, and so the variability should be reduced as the depth increases. This smoothing of long-term accumulation rates will be verified in section 5.5 where we compare model and profile rates vs depth.

4.2. Details: 0–45 km

The profile segment presented in Figure 6 contains relatively gentle topography, with slopes varying between about 0.1° and −0.5° (Fig. 5a). We measured the ice speed at BSC at 11 m a⁻¹. According to Reference Price, Bindschadler, Hulbe and BlankenshipPrice and others (2002), the speed at 45 km is about 30 m a⁻¹. Therefore, the maximum down-ice movement at the bottom of this section (about 300 years) should be ~3 km at the left side and ~8 km near the right side of the profile, considering that the transect direction is no more than about 20° to that of the flow.

Fig. 6. The full radar profile (bottom), profile detail (middle) and elevation (solid curve) and slope (dashed curve) profiles (top) for 0–45 km. The arrow indicates the 1992 horizon. The dashed lines track the loci of fold hinges.

The full-depth profile in Figure 6 shows several folds, some of which appear to migrate down-ice. However, any characteristic feature of the folds, such as the loci of hinges (indicated by the dashed lines), shows either no migration or no more than half of what is expected. for example, the hinges starting at 32 and 35–36 km migrate only about 3 km by the bottom of the profile. As shown later by our modeling, this lack of expected movement is caused either by relatively isolated accumulation anomalies or by the effect of ice movement upon the periodicity within the accumulation variation. In comparison with the next 45 km, the periodicity in the accumulation rate is weak, and so the anomalous behavior may dominate in this section.

4.3. Details: 45–90 km

From 45 to 70 km the topography continues to be relatively gentle (Fig. 7). Over this section, slopes vary from about 0.15° to −0.25°. The transect crosses into BIS at about 50 km. After about 67 km the elevation drops steadily and then fluctuates more strongly from 70 to 90 km where the slope varies between about ±0.8°. The hills and valleys continue into the next section and have topographic wavelengths varying from about 8 to 11 km, with an average near 9.6 km.

Fig. 7. Same as Figure 6, but for 45–90 km. The vertical arrows above the lower profile locate vertical sections along which we later model the stratigraphic dips and corresponding accumulation rates.

After 70 km the topography becomes more variable. At 75 km the most negative leeward surface slope of about −0.8° is encountered, along with a sharp accumulation minimum at about 0.12 m w.e.a⁻¹. At 77 km, just slightly beyond a valley bottom, the slope passes through 0°, elevation begins to increase, and the 8 year accumulation rate reaches 0.28 m a⁻¹ (Fig. 5a). At this location the accumulation peak lags behind the peak of the slope by about 2 km, but all subsequent maxima track almost perfectly with slope maxima, as do the rate minima. The accumulation rate remains strong on this complex rise, then drops to a minimum at the minimum in the slope situated at 82 km, and then rises to a maximum at 88 km, just slightly beyond the slope maximum (Fig. 5a).

The ice speed in this section varies from about 30 m a⁻¹ at 45 km, to at least 48 m a⁻¹ from 65 to 87 km. These ice speeds should give a range of displacement for any variation feature of 8–13 km over 300 years, considering that the transect is angled to the ice-flow direction by about 20°. However, although the variation in Figure 6 becomes relatively strong after 70 km, there is no appreciable or consistent migration to the southwest with depth for any of the variation features. for example, the profile detail shows that the synclinal fold-hinge locus starting at about 77 km is vertically stable to about 10 m depth. It then moves about 2 km down-glacier to about 27 m depth, and then is stable to 56 m. Within this fold, the horizon that starts at ~10 m depth at 73 km descends to ~25 m depth by 78 km.

An unusual, lens-like feature occurs at about 43 m depth between 85 and 90 km. In this vicinity the loci of fold hinges noticeably change direction with depth, a feature which occurs often in the rest of the profile. Along any vertical section, the dips of the fold limbs are generally always positive or always negative, with both cases including values of 0° relative to the surface. This feature might appear to represent a single anomaly in accumulation rate or an unconformity. However, our modeling (section 5) will show that it results from the continual superposition of horizons whose spatially varying accumulation rates are progressively out of phase with each other because of ice motion.

4.4. Details: 90–135 km

Within this segment (Fig. 8), surface slopes vary from −0.9° to 0.5° between 90 and 115 km and then vary from about −0.6° to 0.2° between 115 and 135 km. Both the detail in the middle profile and Figure 5a show how well the accumulation rate tracks with the slope until about 128 km. A rare anticorrelation between surface slope and accumulation rate is centered at about 135 km.

Fig. 8. Same as Figure 6, but for 90–135 km.

According to Reference Price, Bindschadler, Hulbe and BlankenshipPrice and others (2002), ice speed at about 92 km is 40 m a⁻¹, and drops to about 25 m a⁻¹ by about 100 km. This latter speed could be expected to give a fold-hinge locus displacement of about 7 km near the bottom of the profile, again considering the 20° misalignment between transect and flow directions. However, the consistent migration with depth of the fold-hinge locus to the southwest at this distance is only 2.5 km. Although the velocity gradient is about −2 m a⁻¹ between 92 and 100 km, any compressional folding caused by this gradient should be an insignificant cause of these direction changes, in view of the strong fluctuations in accumulation rate along this section, and is unlikely because the firn is too young to have developed much strain history. At 125 km, where the ice speed may drop to 20 m a⁻¹, the locus displacement is only 3 km instead of the expected 5–6 km. More common are loci of fold hinges whose directions vary about the vertical near 94, 111, 120 and 129 km.

The radar failed to record stratigraphy in a few places. The blank section below ~35 m depth near 97 km is artificial. At this location the magnitudes of the stratigraphic dips exceed about 0.4° and were sufficient to cause destructive interference between consecutively stacked traces. Reference Arcone, Spikes, Hamilton and MayewskiArcone and others (2004) give quantitative analysis as to why this is the case. However, weak reflections below ~44 m depth at 108–112, 115–121 and 131–135 km are associated with dips of less value. These zones of weak reflections continue near the bottom of the next section (Fig. 9;135–146 km) and are most conspicuous near its end where our core density profile revealed that the firn/ice transition was reached at only 38 m depth (discussed below), and where the accumulation rate was low. Therefore, we speculate that the edges of these zones also represent transitions to ice of density greater than about 820 kg m⁻³, below which density contrasts would be greatly diminished. The arching horizons above them suggest they were also generated in areas of low accumulation rates.

Fig. 9. Same as Figure 6, but for 135–177 km. The boxed area around the varying hinge locus is discussed in section 5.5.

4.5. Details: 135–177 km

The surface slopes in this segment (Fig. 9) vary from about 0.2° to −0.9°. As can be seen in the detail profile or from Figure 5a, the accumulation variations track very well with the slope. We cannot estimate ice speed throughout this section. However, the horizon depth variations between 135 and 155 km are far less severe than between 70 and 135 km. Based on our modeling presented next, this is consistent with the decreased variability in accumulation rates over this section, and with an ice speed that is probably about 20 m a⁻¹. After 155 km the degree of the variation increases, and so does the range of the accumulation rates. At some point the ice speed does also, and we later estimate that by about 160 km it is about 25 m a⁻¹.

The righthand side of this segment ends at core site 99-2, which was located in a basin where we measured ice speed at 48 m a⁻¹. The ice core shows the firn/ice transition at a relatively shallow 38 m depth (unpublished data), and the location itself has a very low accumulation rate (Fig. 5a). We traced the age of the horizon at this depth back to core site 99-1 and found it is about 300 years old. Therefore, approximately 14 km of displacement after 300 years might be expected for the locus of any fold feature here. However, as with the previous sections, there is no consistent migration down-ice of any variation feature, and the few loci that we traced for the fold hinges are primarily vertical.

5.1. Model parameters

We seek to explain the unusual profile features noted above with a simple model. We take the glacier radar surface to be flat, as viewed in a radar profile, and determine the thickness, y, of an annual deposition with a cosinusoidal accumulation rate such that

The coefficient A is the amplitude of the variable part of the accumulation rate and C is a constant, background rate. We chose A = 0.125 m a⁻¹ and C = 0.225 m a⁻¹ for snow at density ρ = 400 kg m⁻³. These parameters assure there is always horizon separation and that the rate will always be positive because it varies between 0.1 and 0.35 m a⁻¹. The quantity x is distance, v is ice speed (such as 10 m a⁻¹), a is time in years, and λ is the wavelength of topographic variation, which we assume to be 10 km (peak to peak) to simplify comparisons. The actual average hill spacing near 99-1 is about 9.6 km.

We make no assumptions regarding the stability or shape of the surface topography. We seek only the effects of a systematically displaced series of sinusoidally varying thickness functions upon a radar record. We do not correct our model profile for elevation, because this would make the effects we seek less clear, and difficult to compare with the radar profile. Consequently the introduction of any distance phase lag between the thickness function and the surface slope would be meaningless because we only seek the effects of a variable accumulation rate and ice movement at one point in time. Migration of the surface topography (Reference Black and BuddBlack and Budd, 1964; Reference WhillansWhillans, 1975a) would complicate a time evolution picture of the stratigraphy, but whether this truly happens is debatable (Reference WhillansWhillans, 1975a).

The first accumulation surface is measured relative to the instantaneous topographic surface and is followed by succeeding surfaces at depths that are characterized by greater ages, a, and greater displacement down-ice, va. The depth, y _n, to the nth surface at any distance x is the addition of all the layer thicknesses above it, and is expressed mathematically as

This summation produces some interesting features that will be seen in the examples to follow.

The initial density of 400 kgm⁻³ for all layers is a good approximation for that of the first 10 m of snow, which is the thickness of all our modeled layers. We then multiplied the depths, y _n, by a densification (or compaction) function, which simulates the densification at site 99-1, to find a new depth , and applied this function to all of our models. The mathematical form of this empirical function is

where y _n and are negative numbers. The expression in braces was derived by integration of an empirical exponential density function,

to derive the total mass at a depth y. Although this compaction consideration provides the approximately correct depths of 40–50 m for the 300 year surface, it marginally affects the qualitative results discussed below. Our densification function is continuous and provides no abrupt density contrasts. Therefore, the individual horizons we plot could represent very thin density discontinuities, such as those caused by the millimeter-thick ice crusts and associated hoar layers that occur so frequently in ice cores (Reference GowGow, 1968; Reference Arcone, Spikes, Hamilton and MayewskiArcone and others, 2004).

5.2. Changes in uniform ice speed

Figure 10a–c show cases where the ice speed is steady. In Figure 10a the speed is relatively low at 10 m a⁻¹. The accumulation rate differences appear to propagate down-ice uniformly with depth, the dips of the folds become steeper with depth, and the loci of the fold-hinge apices appear to form straight lines. However, after 300 years the deepest hinge axis has moved only 1.5 km down-ice instead of the 3 km expected from the velocity, and reproduces the unexpectedly reduced migration of fold hinges down-ice seen in the profiles. This reduced migration is explained mathematically by either a simple trigonometric argument or an integration of Equation (3). Trigonometrically, the addition of two or more sinusoids, each with a progressively increasing incremental change in phase, ϕ, produces a new function with half the phase shift. By trigonometric identities

Fig. 10. (a–c) Three models of stratigraphic variability caused by ice movement and a sinusoidally varying snow accumulation rate. The ice speeds are labeled. The straight line in (b) does not track the locus of the fold hinges. The scale of the elevation profile is not relevant to the model and is presented only to show where the accumulation maxima and minima would align with the slope. The parameters in the labeled accumulation function are defined in the text, and negative signs are used to plot depths below the surface.

Figure 10b shows the case where ice speed is 30 m a⁻¹. The down-ice migration of the last surface is 4.5 km instead of the expected 9 km. The straight line shows that a locus connecting the fold hinges would deviate from linearity, so that the fold feature is beginning to deform. In addition, the variation in depth of the lowest horizon is far less than in Figure 10a, so that the glacial speed relative to the accumulation wavelength smooths the variation.

Figure 10c shows the case of v = 50 m a⁻¹. The variable direction of the fold-hinge loci is now significant. This ice speed gives one wavelength (10 km) of ice travel after 200 years. The 200 year surface is the horizontal (relative to the surface) horizon at about 29 m depth. Along any vertical section near and across this horizon, the hinge orientation periodically changes from anticlinal above to a down-ice displaced, synclinal below. As shown previously, changes in direction of the fold-hinge loci occur at 80–90, 103–135 and 157–167 km. We model the one between 80 and 90 km below.

5.3. Anomalous changes in ice speed

Figure 11a presents another possible way in which the fold-hinge loci could change direction between relatively up- and relatively down-glacier. In this example, the speed varies with age, starting out at 50 m a⁻¹ for the most recent 100 years, then 10 m a⁻¹ for the second, and finally 100 m a⁻¹ for the third most recent 100 years. These changes translate into additional phase lags in Equations (2) and (3). for example, at year 210, the total translation results from 100 years of motion at 50 m a⁻¹, plus 100 years of motion at 10 m a⁻¹, and 10 years of motion at 100 m a⁻¹. Rather than produce more down-ice translation at the bottom, the pattern appears to produce up-ice translation. for example, the loci of the anticlinal fold hinges that start at 5 and 15 km form a curved line that first slopes down-glacier, then tends toward vertical at about 25 m depth, then starts to slope back up-glacier at about 29 m depth, and then tends back towards vertical at 32 m depth.

Fig. 11. Three models of stratigraphic variability caused by anomalous changes in ice speed (a) and accumulation rate (b, c). The ice speeds are labeled. The parameters in the labeled accumulation function are defined in the text, and negative signs are used to plot depths below the surface. Bold contours are where the accumulation rates change their functional dependence. In (b) the rate is localized in space but not in time. In (c) it is localized in space for 10 years.

The 10 m a⁻¹ speed model of Figure 10a predicts the largest depth fluctuations of the deeper horizons in Figure 3, but the model is unrealistic because such slow speeds occurred only in the first quarter of the profile where topographic and accumulation variations were relatively small. The more realistic model of Figure 10c predicts depth variations of about 7.5 m for any horizon older than 120 years. However, between 70 and 130 km the horizon dated to 1875 varies as much as 14 m, and that to 1734 much more. Although the variable-speed model of Figure 11a produces up to 10 m change in horizon depth after 300 years along with variable directions in fold-hinge loci, and many speed variations are possible, any significant change in speed should cause significant distortion throughout the entire profile in Figure 3, which is not the case. In addition, there is no other evidence that such changes have ever occurred in such a short time-span. Therefore, we think that the situation depicted in Figure 10c provides a more reasonable explanation for loci direction changes because it assumes a single glacial speed.

5.4. Anomalous change in accumulation rate

We have argued that a periodically varying accumulation rate precludes tracing the burial history of an ice particle along the loci of any particular feature, such as a fold hinge. Figure 11b shows that the correct burial trajectory can be traced from an accumulation anomaly localized only in space and superimposed on a uniform deposition rate, as might occur year after year on an isolated windward, monocline surface. A nearly identical model was developed and used by Reference Leonard, Bell, Studinger and TremblayLeonard and others (2004) to explain an accumulation anomaly seen in the Vostok (East Antarctica) ice core. As with the previous models, the surface is normalized to horizontal and we assume each successively deeper layer had the same anomalous accumulation distribution when it was at the surface. The anomalous peak amplitude of 0.25 m a⁻¹ is superimposed on a constant rate of 0.35 m a⁻¹, all for snow at 400 kgm⁻¹ density. The anomaly is represented by a displacement y _a, and has a Gaussian distribution expressed as

where v =50 m a⁻¹. This distribution is seen most clearly in the first accumulation horizon. Each horizon represents 10 years of accumulation. The depths of the successive 10 year surfaces are summed as in Equation (3). The anomalous deposition is centered at 4.5 km, as expressed by the delay of 4.5 km in the exponential. The velocity then centers the anomaly at 5 km distance after 10 years.

The maximum displacement of any feature within the 300 year surface might be expected to be 15 km after 300 years. However, only the leading edge of the anomalies (delineated by the dashed line in Fig. 11b) propagates at the ice speed. In the lowest surface it has moved 15 km down-ice. Other features have spread with depth but are displaced far less. Therefore, trajectories that indicate the correct ice speed could only be recovered from the leading edge of local anomalies that are well isolated from any other variations in accumulation. The spread of the variation at depth is reasonable because the displacements of the accumulation disturbances in the recent, upper layers translate through to the lower ones.

Radar profiles of down-glacier spreading of stratigraphic variability with depth caused by a single anomaly are given by Reference Vaughan, Corr, Doake and WaddingtonVaughan and others (1999, fig. 2c) and by Reference Leonard, Bell, Studinger and TremblayLeonard and others (2004). Vaughan and others show that an anomaly of about 0.5 km width at the surface becomes almost 1.3 km wide at about 96 m depth. Our modeling suggests their interpretation of ice speed from the locus connecting the crests of the synclinal hinges is approximately correct since it does not differ significantly from a locus joining the down-glacier, leading edge of the disturbance, while the up-glacier edge appears vertical.

Figure 11c shows a similar, Gaussian-shaped anomaly that is localized in time (120–130 years BP), as well as in space (centered at 12 km distance), and is superimposed upon a periodic accumulation rate and a steady ice flow at 30 m a⁻¹. Mathematically, this anomalous thickness, y _a, is expressed as

The anomaly doubles the accumulation rate to a peak of 0.7 m a⁻¹ during a 10 year period. As in Figure 11b, the extra accumulation delays the reflections from all subsequent layers. Consequently the anomaly distorts the stratigraphy beneath it and creates a vertical series of secondary fold hinges. A similar effect would result if the anomaly were expressed as an increase in density, which would also supply extra time delay.

There is no indication in our data or modeling of any feature that consistently tilts down-glacier, starting from the bottom of the glacier, as suggested by Reference Whillans and JohnsenWhillans and Johnsen (1983) and Reference Jacobel, Gades, Gottschling, Hodge and WrightJacobel and others (1993). In these cases, they surmised that the cause of the folding occurred far upstream, where basal boundary conditions changed during the transition from inland ice to ice-stream flow.

5.5. Comparisons with profile features

We recreate two features that appear to be central expressions of the stratigraphic variation process. They are (1) steepening fold dips associated with asymptotically vertical hinge loci; and (2) periodically changing fold dips associated with fold hinges whose loci directions waver about vertical. We first attempt to reproduce the slopes and the accumulation rates associated with the fold limbs located near 77 km (within BIS) in Figures 3 and 7, with the isolated accumulation anomaly model of Figure 11b. We chose this feature because the change from a dipping to a vertical direction of the fold-hinge locus at about 25 m depth appears similar to that of our model, and because the ice speed is nearly 50 m a⁻¹ across this section. The synclinal hinges begin near the start of the hilly area crossing BIS and there is little topographic variation up-glacier from this section. The section is not yet in the midst of where the strong, periodically varying accumulation rate occurs. Consequently, this feature should have some aspects of an anomaly localized in space but not in time, because the hillside is a likely permanent feature in response to sub-ice topography.

Figure 12a and b compare the model with the feature. The vertical dashed lines are where the comparisons between accumulation rates (Fig. 12c) and slopes (Fig. 12d) were measured. In Figure 12b the vertical line is placed at 76.5 km where the anomaly appears to originate, and about 260 years of accumulation are considered. Similarly, the vertical line in the model is placed where the anomalous accumulation is centered at the surface. The model uses a densification function that simulates the density profile at 99-1 (at 87 km), and the maximum age is 300 years. As mentioned earlier, in the real profile the horizon at 10 m depth at about 73 km drops to about 25 m depth at 78 km. In the model, the 12 m horizon on the left side drops to about 23 m depth.

Fig. 12. (a–d) A qualitative comparison between the model of Figure 11b (a) and the profile section located near 77 km (b), and quantitative comparisons between yearly average accumulation rates relative to the surface (c) and the slopes of the fold limbs (d). The compared calculations were made along the vertical dashed lines in (a) and (b). The vertical streak near 78 km is a noise band. The curves with data points are for the profile.

The yearly-averaged w.e. accumulation rates in Figure 12c are referenced to the surface (end of year 1999). We obtained horizon dates by tracking the dated horizons from site 99-1 at 87 km. The accumulation rate function (given in Fig. 12a) was chosen to generate the measured rate of 0.28 m a⁻¹ for the 8 year old surface at 77 km (Fig. 5a). Qualitatively, the model variations track those of the profile. Quantitatively, the model gives similar peak values, as it should, and a similar asymptotic rate of about 0.12 that is consistent with the calculations at this distance in Figure 5b. The rates decrease from a peak value because lower-accumulation-rate snow has moved beneath the peak of the accumulation. Figure 12d shows that the model closely simulates the relatively rapid decrease in slope, from 0° (the normalized surface) to about 0.3° at 15–16 m depth, the knee in the slope variation, and then the near-steady value of about −0.22° to the bottom of the record.

Figure 13a and b compare the model in Figure 10c with a segment of the profile located near core site 99-1 at 87 km. We chose this feature because the change from anticlinal folding relatively up-glacier, to synclinal relatively down-glacier at about 4 km before the topographic peak, occurs in both model and profile, and because the ice speed across this section is close to 50 m a⁻¹. The model uses a sinusoidal accumulation rate amplitude of 0.125 m a⁻¹ superimposed upon a steady rate of 0.30 m a⁻¹. A horizontal (relative to the surface) horizon occurs at about 38.3 m depth (time delay of 386 ns) in the profile. This depth dates to about 1813 (personal communication from D. Dixon, 2003), or an age of 186 years. As discussed above, this horizon should occur at about 200 years BP (1999) given the glacial speed of about 50 m a⁻¹ and the accumulation wavelength variation of about 10 km.

Fig. 13. (a–d) A qualitative comparison between a model slightly different than the one shown in Figure 10c (a) and the profile section located near 87 km (b), and quantitative comparisons between yearly average accumulation rates relative to the surface (c) and the slopes of the fold limbs (d). The compared calculations were made along the vertical dashed lines in (a) and (b). The curves with data points are for the profile.

The slight discrepancy in years between the model and profile horizontal stratum is explained by the slightly less hill spacing of the profile data. The age, a _n, in years of the horizontal line is found from the simple formula

At a wavelength of 10 km and speed of 50 m a⁻¹, a _n = 200 years. If the wavelength is reduced to the actual value of 9.6 km, then a _n = 192 years, which improves the agreement with the data, but at a speed of 48 m a⁻¹ it is still 200 years. We do not think that this computation should consider the vector component of the speed in the direction of the transect (46 m a⁻¹), because the depth of this horizon would then steadily increase to infinity as the transect changed direction.

Figure 13c shows good agreement between the model and the profile yearly-averaged w.e. accumulation rates, and peak and asymptotic values at 86.5 km. The asymptotic value of about 0.11–0.12 agrees with the calculations for this distance in Figure 5b. The w.e. values decrease from the peak because lower-accumulation-rate snow has moved down-glacier and beneath the peak of the accumulation. Figure 13d shows good agreement between fold slopes, which vary sinusoidally and are always negative.

A second possible example of this phenomenon appears in the boxed area in Figure 9 between about 158 and 165 km. It occurs at about 35 m depth at the midpoint of the fold-hinge locus that changes direction. This segment was recorded within a series of four hills whose average peak-to-peak spacing is 4.9 km. From 160.5 to 162 km the folds change from anticlinal to synclinal going down-glacier, with an intervening horizon at about 35 m depth. This horizon can be traced back to core site 99-1 where it coincides with the same 200 year horizon discussed in Figure 13. Using Equation (9), the ice speed at this location should then be about 25 m a⁻¹. Although ice speed increases to 48 m a⁻¹ by another 16 km, this is not an unreasonable value given the decrease in ice speed seen in the data of Reference Price, Bindschadler, Hulbe and BlankenshipPrice and others (2002) after 110 km.