Introduction
The Patagonian icefields have been recognized as one of the main temperate mid-latitude ice masses on Earth (Reference Warren and SugdenWarren and Sugden, 1993), with an important contribution to sea- level rise during recent decades, due to the retreat and thinning of their glaciers (Reference Rignot, Rivera and CasassaRignot and others, 2003). This process of deglaciation (Reference Aniya, Sato, Naruse, Skvarca and CasassaAniya and others, 1997) has been associated with the warming trend (Reference Rosenbluth, Fuenzalida and AceitunoRosenbluth and others, 1997) and a reduction in precipitation (Reference Rosenblüth, Casassa and FuenzalidaRosenbluth and others, 1995) observed during the last century in Patagonia. Only a few glaciers of Hielo Patagrónico Sur (HPS; southern Patagonia icefield) (Fig. 1) have been evaluated in relation to ice-elevation changes, and no systematic mass-balance programme is currently available. All of these measurements have been carried out for the ablation areas of the glaciers, indicating high thinning rates in most cases (Reference Naruse, Aniya, Skvarca and CasassaNaruse and others, 1995). Glaciar Chico is selected here because it is one of the few glaciers in Patagonia for which direct global positioning system (GPS) measurements of ice-elevation changes are available in both the accumulation and ablation areas. As well as these direct measurements, several photogrammetric surveys have been carried out on this glacier, allowing investigation of surface elevation and areal changes since 1975. Based upon these topographic changes, the mass balance of the glacier has been inferred using the geodetic method (Reference KrimmelKrimmel, 1999). Glaciar Chico is one of the easternmost glaciers of HPS, and as a consequence of its more continental location it has been considered more sensitive to temperature changes than most of the glaciers located at the maritime margin of HPS (Reference Warren and SugdenWarren and Sugden, 1993). Compared to other HPS glaciers, the glacier front has been retreating at low rates (17-52 m a-1) during the last century. At present, the ice front is only partially calving into one arm of Lago O’Higgins (Brazo Chico of Lago O’Higgins), where the ice is grounded in a small island (‘La Isla’; Fig. 2), which was initially visible in the 1975 aerial photographs. The relatively low calving rates observed on the glacier suggest that Glaciar Chico may be more sensitive to climate-forcing factors than to local factors such as the geometry of the valleys, depth of the fjords or buoyancy at the ice front (Reference Naruse, Aniya, Skvarca and CasassaNaruse and others, 1995). The accumulation area of Glaciar Chico is divided into several basins, the main one located in the so-called ‘Paso de los Cuatro Glaciares’ where Glaciar Chico shares a plateau with Glaciares Viedma and O’Higgins (Fig. 2).
The main aims of this paper are to (i) present new data of surface topography obtained at different dates and with different methods on Glaciar Chico, (ii) analyze and interpret the ice-elevation changes, inferring the mass balance of the glacier, and (iii) discuss possible explanations for the observed changes.
Datasets
Remotely sensed imagery
Several sets of aerial photographs from the period 1945-97 are used here (Table 1). The stereo pairs allow interpretation of frontal and areal changes. Only a few satellite images with reduced cloud coverage have been acquired for HPS (Table 2). Among them, one Advanced Spaceborne Thermal Emission and Reflection Radiometer (ASTER) satellite image of 2001 was obtained from the United States Geological Survey (USGS) EROS Data Center. Due to the gain values used when the ASTER image was acquired (Table 3), the upper accumulation area of HPS is ‘saturated’ (i.e. there is no variation in brightness), preventing discrimination of different snow surfaces or variations in albedo. There is, however, adequate contrast for the ablation area. Five Landsat satellite images were available for the period 1984-2002. They were geo-rectified using the available orbital parameters and 17 ground-control points (GCPs) using ERDAS Imagine 8.2 commercial software, yielding a horizontal error of 43 m. The ASTER scene was ortho- rectified during the digital elevation model (DEM) generation process, yielding a horizontal error <15 m.
Field measurements
Several field campaigns have been carried out to the accumulation area of Glaciar Chico since 1996, when the Chilean Air Force began the so-called ‘Hielo Azul’ annual campaigns which train pilots to perform snow landing and take-off of Twin Otter aeroplanes equipped with skis.
The first campaign was carried out in the accumulation area of Glaciar Chico in January 1996, when a 13 m metal mast was buried to 3 m in snow at the base camp, located at 1440ma.s.l. on a flat area with a 1.4° surface slope, approximately 1 km north of the local ice divide between Glaciares Chico and Viedma (Fig. 2). The second campaign was carried out in October 1997 when topographic quality Trimble Geoexplorer II GPS receivers were used to survey the glacier surface (Fig. 3). The receivers allow singlefrequency data collection in coarse acquisition (C/A) code mode and phase information, and 2-7 m vertical precision after applying differential correction. The third and fourth field campaigns were carried out in September-October 1998 and 2001, when surface topography points (Fig. 3) were measured using geodetic quality GPS receivers (Javed Legacy and Leica model 300 dual frequency with C/A code and phase information, 10-30 cm precision after differential correction applied).
Methods
GPS data
All GPS data were collected from areas covered by snow during the field campaigns, mainly in the accumulation area of the glacier, where snowmobiles transporting the GPS receivers could move easily over snow. As the field campaigns were carried out between September and October (early spring), some measurements took place in the upper ablation area of the glacier.
Most of the GPS data were obtained in kinematic mode, where one receiver was installed on board a snowmobile, collecting data every 5 s. A second receiver was installed simultaneously on a rock outcrop (Nunatak Garcia) located approximately 5 km from the base camp (Figs 2 and 3). The receiver was attached to a 10 cm metal pin glued into a hole drilled in the rock. This second receiver was configured as a base station, collecting data every 5 s. A total of 25 points were measured in 1998 and 2001 using a rapid-static procedure (Fig. 3), and at each point a geodetic quality GPS receiver was installed on a stake for ~15 min, collecting data every 5 s. All GPS data were analyzed using Geogenius commercial software version 1.6, applying a differential correction procedure. The GPS data points were distributed mainly between the temporal snowline and the ice divide with Glaciar Viedma (Figs 2 and 3).
DEM generation from aerial photographs
From vertical aerial photographs (1975, 1981 and 1997), DEMs were generated using Socet Set version 4.4.0, professional digital photogrammetry software, at the University of Zurich, Switzerland. For that purpose, the photographic paper copies were scanned at 800dpi using an EPSON expression 1640 XL scanner, giving a pixel resolution equivalent to 2.5 m on the ground. The calibration information of the camera (coordinates of the fiducial marks and the camera focal length) were used to create the internal orientation of the stereoscopic pair. Tens of tie points were used to connect together all the photographs from each date, generating a photogrammetric mosaic (a so-called image block) with relative orientation. A set of 16 GCPs located on rock available from previous cartography, as well as from some GPS points measured in the surrounding area of the glacier, were used to assign an absolute georeference to the photogrammetric pairs (so-called models). These points were selected from distinguishable landmarks present in the aerial photographs, such as summits, peninsulas and geologic patterns. Triangulation and bundle adjustment was calculated for each model, using GCPs and tie points. With this adjusted model it was possible to calculate the absolute orientations of the entire image block. Once the first image block was solved, the second and third blocks were adjusted to each other, using tie points obtained simultaneously from all the blocks. With the resulting oriented imagery, the photogrammetric software automatically calculated an absolute DEM for each date. A 50 m pixel size was chosen considering possible deformations due to the aerial-photograph scanning process, as well as other horizontal problems derived from geolocation of GCPs. The DEMs covered mainly the ablation area of the glacier, where the stereomatching process was successful. The snow-covered upper parts were not completely covered, especially around the mast installed in 1996, where the low contrast in the aerial photograph prevented any altitude information from being obtained (Fig. 3).
DEM generation from ASTER
Based upon bands 3N and 3B images from one ASTER scene acquired on 14 October 2001 (Table 3), a DEM was generated using PCI Geomatica commercial software at the University of Zurich. The same dataset of GCPs used for the generation of DEMs from aerial photographs was used for the bundle adjustment of the ASTER bands. These points were used in the absolute rectification of the model. Approximately 80 tie points were selected on the screen from all around the image, in order to ensure inter-band registration. By using all the tie points, the epi-polar image was generated by means of rotation and resampling of the bands. An automatic relative DEM was extracted from the epi-polar image, and this was used in combination with the GCPs in the geocoding process to generate an absolute DEM. The resulting DEM was edited in order to eliminate areas with artefacts, especially in the accumulation area of the glacier, where the lack of contrast does not allow any stereo matching. A 50 m pixel size was chosen in order to homogenize all data for further comparisons, and include possible horizontal displacements due to the geolocation of GCPs. The ASTER-derived DEM covered most of the ablation zone of the glacier, but not the snow surface areas, where the saturation of the satellite image prevented any stereo matching (Fig. 3).
Errors
Horizontal errors
The resulting horizontal rms errors for the image rectification, and thus the overall position accuracy of the DEMs generated from the aerial photographs and the ASTER image, were automatically calculated by Socet Set and PCI Geomatica commercial software, yielding 10 and 15m respectively. These horizontal errors are due to the GCPs employed in the DEM generation process, and to possible deformations after scanning the photographs, errors in identifying the GCPs on the imagery and/or deformations of the aerial photograph edges.
The horizontal error of the GPS data was calculated automatically using Geogenius commercial software, yielding on average 0.10 m for points measured with geodetic quality GPS receivers, and 4.5 m for points measured with topographic quality GPS receivers. In order to compare these data with the previous DEMs, the GPS data were transformed to a raster format, where 50 m pixel size DEMs were generated for each dataset.
Vertical error of GPS data
All the GPS measurements were carried out in relatively flat areas away from high mountains, with distances between kinematic and base measurements no longer than 10 km. Under these conditions, and considering only the lowest Precision Orbit Determination (POD) values, the resulting vertical errors were 7 m for measurements based upon topographic quality GPS receivers and better than 0.30 m for geodetic quality GPS receivers.
Vertical error for aerial photograph DEMs
The vertical error obtained for the image rectification from GCPs, and thus the overall vertical accuracy of the digital photogrammetric DEM, was computed by Socet Set, yielding 6.0 m for each of the DEMs. However, when rock outcrops were compared between DEMs derived from aerial photographs from different epochs, an average standard deviation of 12 m was obtained, with biases of 1.4-1.9m. The errors have a normal distribution. The 12 m error is related mainly to the inaccuracies of the GCPs obtained from the regular cartography of the Instituto Geografico Militar of Chile, and is comprised of 8.5 m error for each single DEM (Table 4), which is in good agreement with the result obtained automatically by Socet Set.
Vertical error for the ASTER DEM
The elevations of rock outcrop areas were compared between the ASTER DEM and the aerial photograph DEMs, yielding a mean bias of 22 m and a random error of 24 m. If the random error for the aerial photogrammetric DEMs is taken as 8.5 m, the ASTER-derived DEM random error is 22 m. The differences in rock areas were spatially distributed in relation to the slope aspect derived from the ASTER DEM (Fig. 4), with maximum positive values on western slopes, and maximum negative values on eastern slopes. A best-fit curve was calculated to account for such differences, yielding a cosine relation. We assume that this systematic deviation is due to the ASTER stereo geometry (Table 3). Compared to the ASTER nadir channel (3N), the 3B stereo sensor is directed backwards by 27.6° from nadir. Thus, compared to the nadir image, northern slopes appear significantly shortened and southern slopes significantly stretched in the back-looking image. Such perspective distortion could lead to the systematic errors encountered in the automatic ASTER DEM generation (Reference KääbKaab, 2002).
The cosine relation was applied to the ASTER DEM, giving an adjusted DEM. The rock areas of this adjusted DEM were compared with the aerial photograph DEMs, resulting in a mean error of-3.7 ± 19 m. This total random error is comprised of 8.5 m rms error for the DEMs derived from aerial photographs and a final rms error of 17 m for the ASTER DEM.
Results
Glacier basin
Using the DEMs and the GPS measurements carried out between 1997 and 2001, it was possible to delineate the boundary of the glacier, especially along the ice divide between Glaciares Viedma and Chico (Fig. 2), determining a total ice area of 193 km2 for Glaciar Chico in 2001 (without bedrock exposures).
In order to estimate the equilibrium-line altitude (ELA) of the glacier, we assume that in Patagonia the snowline at the end of the Southern Hemisphere summer (February-March) corresponds to the ELA, which is supported by the fact that temperate ice conditions and high accumulation prevailing at HPS (Reference ShiraiwaShiraiwa and others, 2002) do not allow generation of superimposed ice at the surface. For that purpose, a Landsat Enhanced Thematic Mapper Plus (ETM+) satellite image acquired on 21 March 2001 (Table 2) was digitally analyzed based upon the band ratio 4-5. As a result, the ELA in 2001 was estimated to be at an average altitude of 1320ma.s.l. (Fig. 2), which results in an accumulation-area ratio (AAR) of 0.55.
Areal and frontal changes
Glaciar Chico has been retreating and shrinking since 1945, when the first aerial photograph of the glacier was obtained. In 1945, part of the glacier front was calving into Brazo Chico of Lago O’Higgins (Fig. 5). The rest of the glacier was flowing into a lateral valley, Valle Escondido, where nearby Glaciar GAEA was calving into an ice-dammed lake, Lago Escondido (Fig. 2).
A new set of aerial photographs obtained in 1975 (Table 1) showed slight retreat around the non-calving terminus, but more significant changes at the calving front of the glacier. One of the main changes in the dynamics of the calving front was the appearance of a small island (‘La Isla’) in the aerial photographs of 1975, which was located at the northeast calving glacier front in the middle of Brazo Chico of Lago O’Higgins. This island has increased its surface area since then, reducing the calving area of the tongue to two small sections around ‘La Isla’ (Fig. 6).
The maximum frontal change of the glacier was measured between 1945 and 1975, when the calving frontal tongue retreated at a rate of 52 m a-1. Since then, the glacier has continued to retreat, but at a slower rate (2517m a-1). These frontal retreat rates do not necessarily reflect the real magnitude of changes affecting the glacier. For instance, between 1945 and 2001, Glaciar Chico lost 6.1% of a total surface that was 205 km2 in 1945. The main surface change is taking place in the lower ablation area, where the glacier has been shrinking at all its margins, especially at the calving frontal tongue (Fig. 6). A second significant change is taking place at the ELA of the glacier, where the ice has been receding around Nunatak Garcia (Fig. 2).
Assuming an ice-thickness range of 200-400 m for the wasting area and an ice density of 900 kgm-3 (Reference AniyaAniya, 1999), a volume rate due to area change of -0.064± 0.059 km3 w.e. a-1 was obtained between 1975 and 2001.
Ice-elevation changes since 1975
Table 4 shows the average ice-elevation change measured between different datasets and DEM-1975. In a preliminary analysis of the results, ice thinning is observed in all epochs, with smaller rates between 1975 and 1981 and maximum rates at the end of the study period, in a trend that suggests accelerated thinning. The errors are higher between datasets separated by short periods of time, for example between DEM-1975 and DEM-1981, with thinning rates in agreement with the general trend, but with magnitudes smaller than the errors. The best results were obtained comparing the geodetic quality GPS measurements with the 1975 data.
The differences in altitude between DEMs are more representative for the ablation area, where the coverage is more extensive. Figure 7 shows the differences between DEM-1975 and DEM-1997. The maximum ice thinning is observed at the frontal tongue of the glacier in both the non-calving (250-400 m a.s.l.) and the calving front (250 m a.s.l.), where a maximum thinning rate of 5.4 ± 0.55 m a-1 was measured. At the central flowline of the glacier, between 500 and 1200 m a.s.l., ice thinning rates are smaller, fluctuating between 2.7 and 1.6 ± 0.55 m a-1.
Many hundreds of pixels with GPS altitude values measured between 1997 and 2001 and located between 1100 and 1400 ma.s.l. were compared to DEM-1975, yielding an average thinning rate of 1.7 m a-1 (Table 4).
When the ASTER DEM is compared with the other datasets, the same pattern of thinning rate is obtained for the ablation zone of the glacier, with a mean thinning rate of 1.7 ± 0.7 m a-1 between 1975 and 2001, and 1.6 ± 0.9 m a-1 between 1981 and 2001.
Analyzing 13 rapid-static points measured with geodetic quality GPS receivers in the 1998 and 2001 campaigns in the ablation area of the glacier at 1100-1320ma.s.l., an average thinning rate of 1.5 ± 0.14 m a-1 was obtained. The differences presented a normal distribution, with a minimum thinning rate of 1.3 m a-1 and a maximum of 1.9 m a-1. Considering 12 geodetic quality rapid-static GPS measurements carried out at the accumulation area of the glacier between 1320 and 1450 ma.s.l. in 1998 and 2001, the mean thinning rate was found to be 1.9 ± 0.14 m a-1, with a range of 2.3-1.4 m a-1, and no apparent altitudinal or spatial pattern.
One of the measured points at higher altitude in the glacier was the mast installed at 1440ma.s.l. in January 1996 (Fig. 2). Between 1 October 1998 and 3 October 2001, when the top of the mast was measured with geodetic quality GPS receivers, the structure moved to the north at a horizontal velocity of 53 m a-1 with a total vertical displacement of the mast of -10 m. Considering the surface slope in this area (1.4°), the vertical displacement due to ice advection was estimated to be -3.8 m, resulting in a total ice thinning of 6.2 m (2.0 ± 0.14 m a-1).
Analyzing snow densities and the height of the mast above the snow surface, Reference RiveraRivera (2004) estimated an average snow accumulation of 0.57±0.21ma-1w.e. between 1994/95 and 2000/01. Thinning is therefore three times the annual accumulation, suggesting that in addition to climate, ice-dynamic factors might be responsible for forcing the thinning.
Considering the ice-elevation and areal changes experienced by the whole glacier between 1975 and 2001, a total volume loss of 7.83 km3 w.e. (0.29 ± 0.097 km3 w.e. a-1) was estimated for Glaciar Chico.
Discussion
Ice-elevation changes observed at Glaciar Chico show altitude patterns similar to those of other Patagonian glaciers (Reference Naruse and SkvarcaNaruse and Skvarca, 2000), with maximum thinning at low altitude in the ablation area, and minimum values in the accumulation areas. Nevertheless, the large thinning of 2.0 m a-1 observed in the accumulation area of the glacier is less common. The ice-elevation changes measured for the lower part of the glacier are spatially and dynamically related to calving activity since the maximum thinning rates were measured around the two sections of the ice front where the glacier is losing ice by calving. The changes that have taken place at mid-altitudes of the glacier, with reduction of thinning with altitude, are in agreement with the observed thinning rates of other glaciers (Reference Rignot, Rivera and CasassaRignot and others, 2003). Nevertheless, the estimated warming of 0.02°Ca-1 affecting the region (Punta Arenas’s warming trend;Reference Rosenbluth, Fuenzalida and AceitunoRosenbluth and others, 1997) is not large enough to totally account for such changes.
Ice-elevation changes measured for the accumulation area of the glacier could be partially explained by a warming trend of 0.98°C detected at 850 hPa geopotential height obtained from US National Centers for Environmental Prediction/National Center for Atmospheric Research (NCEP/NCAR)-40 re-analysis data between 1958 and 1985 (Reference Carrasco, Casassa, Rivera, Casassa, Sepulveda and SinclairCarrasco and others, 2002). More probably, ablation due to warming would occur as melting. The sequence of melting and refreezing events, evident from the number of ice layers, ice pipes and ice columns observed at snow pits in the area (Reference RiveraRivera, 2004), is the more likely process affecting the lower part of the accumulation area, increasing the densification rate of the firn. The reduction in the snow accumulation could be a key factor in the glacier response to a warmer climate. A dynamic component such as excess creep (ice thinning from longitudinal stretching) could also be important in explaining ice thinning in Patagonia (Reference Rignot, Rivera and CasassaRignot and others, 2003).
Conclusions
A general warming trend and reduced precipitation have been observed at most stations in Patagonia. This climate- change trend is the background driving factor for the glacier changes in the region, but specific responses are much more complicated, especially given that most glaciers are calving. Frontal and areal responses of Glaciar Chico to the climate change taking place in Patagonia seem to be smaller than the frontal changes of most other glaciers in HPS. These smaller values can be explained mainly by the limited calving activity on the lower tongue of Glaciar Chico. Nevertheless, the rates of ice-elevation changes affecting the ablation area of Glaciar Chico are in agreement with the thinning rates measured for other glaciers of HPS, indicating that thinning of this magnitude expresses regional climatic signals.
The thinning rates measured in the accumulation area of Glaciar Chico could be a response to the warming trend at the 850 hPa geopotential height as detected from the NCEP/NCAR-40 re-analysis data. However, not all the thinning can be explained by warming alone, so it seems likely there is also a dynamic response of the ice mass.
Acknowledgements
This research was done with the support of Centro de Estudios Cientfficos (CECS), the Republic of Chile scholarship of MIDEPLAN, the University of Chile and FONDECYT 1000445. F. Cawkwell edited the manuscript. H. Lange provided valuable acquisition and analysis of GPS data. The Chilean Air Force provided logistic support during the Hielo Azul operations. J. Quinteros, J. Araos, A. Giannini, M. Arevalo and R. Traub (Trauko) gave logistic and mountaineering support during the field campaigns. Global land Ice Measurements from Space (GLIMS) and the USGS gave unrestricted access to ASTER and Landsat ETM+ images. M. Aniya contributed the Landsat TM mosaic of 1986. CECS is a Millennium Science Institute and is funded in part by grants from Empresas CMPC, Fundacion Andes and the Tinker Foundation. The comments of the editor R. Naruse and two anonymous referees are acknowledged.