How It Works
The PTAA model was developed as an alternative to other methods of measuring glacier mass balances. The data requirements are daily observations of precipitation, maximum and minimum temperatures at a nearby weather station (or an average of two stations), and the area-altitude distribution of the glacier's surface. For some glaciers the weather stations can be as much as 300 km, from the glacier but usually the distance is closer to100 km. or less. The quality of the data appears to be more important than the distance to the glacier. The area-altitude distribution (hypsometery) is developed from DEM models or topographic maps. A finely divided area interval produces the most accurate results, e.g. a 10 meter interval between areas is preferable to 50 meters.
The model relies on the climate history that is embedded in the area-altitude (AA) distribution. The AA distribution is developed by glacier flow and erosion of the underlying bedrock, both of which are dependent on long-term-mass balances and on the climate. The current mass balance is calculated by combining precipitation and temperature observations at a regional weather station with the glacier's current area-altitude distribution to produce daily and annual balances.
To calibrate the model, seven different balance variables are formulated for each day of the period of record and regressed against each other throughout the summer ablation season, for example, the daily balance versus the transient equilibrium line altitude or the snowline altitude versus the AAR. The several hundred regression errors are averaged to produce an optimizing function or calibration error which is minimized. Optimum values are calculated for fifteen coefficients that convert observed precipitation and temperature to daily snow accumulation and snow and ice ablation for each altitude interval. A detailed description of the model is provided in Tangborn, 1999.
The model has been applied to determine annual balances of South Cascade Glacier and compared with manually measured balances for the 1959-1996 period. It has also been applied to several glaciers in Alaska to determine annual balances that are compared with manual balances measured by the USGS (Bhatt, et al, 2007), Korn, (2010), Zhang, et al, (2007a, 2007b) and one glacier in Nepal (Tangborn and Rana, 2000).
Verification of the balances produced by the model is by comparison with annual balances measured the with glaciological method, using snowpits and ablation stakes. To assure an unbiased comparison, the model balances are generated first, before the measured balances are known. A linear regression is then run between measured and modeled balances. Examples of this procedure for the Gulkana and Wolverine Glaciers in Alaska and Vernagtferner in the Austrian Alps are shown on www.ptaagmb.com.
Bhatt U.S., J. Zhang, W.V. Tangborn, and C.S. Lingle, L. Phillips, 2007: Examining Glacier Mass Balances with a Hierarchical Modeling Approach, Computing in Science and Engineering, 9 (2), 61-67. Abstract here.
Korn, D.,2010, "Modeling the mass balance of the Wolverine Glacier Alaska USA using the PTAA model", American Geophysical Union, Fall Meeting 2010http://adsabs.harvard.edu/abs/2010AGUFM.C21C0560K
Tangborn, W. and Rana, B.,2000, Mass Balance and Runoff of the Partially Debris-Covered Langtang Glacier, Debris-Covered Glaciers, Edited by M. Nakawa, C.F. Raymond, & A.Fountain, IAHS Publication 264.
Tangborn, W.V.,, 1997, Using low-altitude meteorological observations to calculate the mass balance of Alaska's Columbia Glacier and relate it to calving and speed. Report of a Workshop, February 28 - March 2, 1997, Byrd Polar Research Center, Report No. 15 PDF of paper here.
Tangborn, W.V., 1999, A Mass Balance Model that Uses Low-altitude Meteorological Observations and the Area-Altitude Distribution of a Glacier , Geografiska Annaler: Series A, Physical Geography Volume 81, Issue 4, December 1999, Pages: 753-765. PDF of this paper is available on www.ptaagmb.com (under How It Works).
Zhang J. , U.S. Bhatt, W. V. Tangborn, and C.S. Lingle, 2007: Response of Glaciers in Northwestern North America to Future Climate Change: an Atmosphere/Glacier Hierarchical Modeling Approach, Annals of Glaciology, Vol. 46, 283 - 290. PDF of paper here.
Zhang, J., U. S. Bhatt, W. V. Tangborn, and C. S. Lingle, 2007: Climate downscaling for estimating glacier mass balances in northwestern North America: Validation with a USGS benchmark glacier, Geophysical Research Letters, 34, L21505, doi:10.1029/2007GL031139. PDF of paper here.
This paper was originally published in Geografiska Annaler: Series A, Physical Geography Volume 81, Issue 4, December 1999, Pages: 753-765
The work done developing this model is the basis of the work presented here and the best explanation of how it works.
A MASS BALANCE
THAT USES LOW-ALTITUDE METEOROLOGICAL
OBSERVATIONS AND THE AREA-ALTITUDE
DISTRIBUTION OF A GLACIER
HyMet Company, Seattle, Washington, USA
Tangborn, Wendell V., 1999: A mass balance model that uses low altitude meteorological observations and the area-altitude distribution of a glacier. Geogr. Ann, 81 A (4): 753-765.
ABSTRACT. A glacier mass balance model that requires only low-altitude precipitation and temperature observations and the glacier's area-altitude distribution is presented as an alternative to direct field measurements. Input to the model for South Cascade Glacier are daily weather observations at stations 30-60 km from the glacier and at altitudes 1300 to l500 m lower than the glacier. The model relies on the internal consistency of mass balance variables that are generated by simulation using the low-altitude weather data. The daily values of such balance variables as snowline altitude, zero balance altitude, glacier balance, balance flux and the accumulation area ratio are correlated throughout the ablation season using two-degree polynomial regressions to obtain the lowest fitting error. When the minimum average error (or maximum R2) is attained, the generated balances and other variables are considered to be real. A simplex optimization technique is used to determine the optimal coefficient values that are used in algorithms to convert meteorological observations to snow accumulation and snow and ice ablation. The independently produced simulation results for the l959-l996 period are compared with balances measured at the glacier. The agreement between annual balances for individual years is fair and between long-term volume changes measured by the geodetic method is excellent.
Key words: glacier, glacier mass balance, climate
Measuring the mass balance of a glacier by the usual ﬁeld techniques of ablation stakes and snow pits is labor-intensive and expensive. A method is needed that will produce balance results with acceptable accuracy but requires only readily-obtained meteorological data. Although most weather stations with reliable records are not located at or near glaciers, there are thousands of stations within 50-100 km of many glaciers worldwide. In addition, the wide availability of long-term (50 to 100 years) weather records would allow reconstruction of mass balances for as long a time period as these weather records exist.
However, weather records from distant stations alone will not determine a glacier's mass balance. Another parameter that is unique to each glacier is required, so that when combined with meteorological observations of daily precipitation and temperature, snow accumulation and snow and ice ablation at speciﬁed altitudes on the glacier can be determined. The distribution of area as a function of altitude (the glacier's shape or AA proﬁle) is a distinct feature of each glacier that is the end result of thousands or millions of years of glacier erosion. Although a link between climate, erosion and mass balance cannot be easily proven, it is reasoned that the combination of Weather observations (representing the climate) and the AA proﬁle (representing erosion) should provide sufﬁcient information about a glacier to calculate its current mass balance.
A computer model that simulates both the current and the historic mass balance using available weather data and the area-altitude distribution has been developed and is presented here. By the application of available meteorological observations, namely precipitation and temperature, from distant and much lower altitude stations, plus the distribution of the glacier's area as a function of altitude, the mass balance is found for each day and for each of the area-altitude intervals taken from topographic maps (designated the PTAA model). To demonstrate model accuracy compared with ﬁeld measurements, the balance is ﬁrst simulated for South Cascade Glacier by using the internal consistency of several generated balance variables to ﬁnd optimum coefﬁcient values. After this task is completed, the simulated balances are then compared with observed balances measured for this glacier since 1959. A similar approach was ﬁrst applied to Columbia Glacier in Alaska (Tangborn 1997). However only a few isolated balance measurements have been made on Columbia Glacier so the model's reliability for this glacier is inconclusive. An earlier balance model for South Cascade Glacier (PT model) also used low-altitude weather data (Tangborn 1980). The difference between the PT and the PTAA model is that the PT model used the measured balances to force-ﬁt balances against meteorological observations; while the PTAA model uses the internal consistency of generated balance variables.
South Cascade Glacier is relatively small (currently about 2 km2) and is losing mass and area rapidly. It is located in the eastern edge of the North Cascade Mountains of Washington State (Figs 1, 2). The glacier has been retreating nearly continuously for the past century after reaching a maximum neo glacial length in the early 1600s, during the Little Ice Age (Miller 1969). The maritime climate of this region is characterized by mild, wet winters and cool, dry summers. Snow depths at higher altitudes often reach 10 m.
Fig. 1. Map of Washington State showing locations of South Cascade Glacier and the two weather stations used in this study to represent the glacier's climate.
Fig. 2. South Cascade Glacier: (a) 29 Septemberl960 (Photo by Austin Post); (b) 9 September 1992 (Photo by Robert Krimmel, US Geological Survey). The glacier lost 30 m (w.e.) of mass during this 32 year period. The glacier was approximately 0.8 km in width at the widest point (at the ELA) and 3.4 km long in 1960.
The glacier has been intensively studied since 1957 and has complete and detailed mass balance, meteorological and hydrological observational records (Meier and Tangborn 1965; Tangborn 1968, 1980; Krimmel 1994, 1995, 1996, 1997). The annual, winter and summer balances are measured each year using stakes, snow pits and by coring. During the 1960s, 10-20 ablation stakes were used each year for calculating the balance. This number has slowly declined and since about 1990, four index stakes located at representative points are used each year. Streamﬂow has been measured continuously since 1960, and meteorological observations have been collected periodically during early years and nearly continuously since about 1993.
Several long-term weather stations in or near the North Cascade Mountains were tested to represent the glacier climate. Only two appear to adequately represent weather conditions at South Cascade Glacier: Concrete and Diablo Dam (Table 1). Stations with long periods of missing observations were omitted from the testing procedure, which is described in more detail below. The correlation of daily precipitation among these stations is poor, suggesting that the correlation between glacier and gage precipitation is unreliable. However, the correlation between observed precipitation and seasonal runoff from nearby basins, some highly glacierized, is high (80 to 90%), signifying that a fair correlation must exist between high-altitude snow accumulation and low-altitude precipitation (Rasmussen and Tangborn 1976; Tangborn and Rasmussen 1977).
Table 1. Precipitation stations.
The temperature records collected at Diablo Dam and Concrete were both tested to represent high-altitude temperatures. Not surprisingly, Diablo Dam, the closest station and with the least elevation difference to the glacier, proved to be the most reliable by producing the highest R2 among mass balance variables. Weather observations at Concrete began in the early 1920s; however, electronic ﬁles for both stations are currently available for only the 1932-1996 period.
The surface area of the glacier between each 20 m altitude increment is determined from the 10 maps that have been made at regular time intervals since 1955 (1955, 1958, 1961, 1964, 1970, 1977, 1980,1985, 1991, 1992) (Krimmel, pers. comm.). Two of the AA distributions, for 195 8 and 1992, calculated from these maps are shown in Fig. 3. The AA distribution that is used for calculating balances for a particular year is simply the one nearest to that year. For example, the 1958 AA distribution is used for 1957, 1958 and 1959 balance calculations. In 1955 the glacier extended from 1610 to 2160 m, a 550 m range, therefore 28 altitude intervals are used.
Fig. 3. South Cascade Glacier area-altitude profiles in 1958 and 1992. The altitude interval is 20 m.
Simulation model structure
Fourteen coefficients are used in algorithms that convert daily meteorological observations to snow accumulation and snow and ice ablation at the glacier. Seven balance variables are then calculated daily from these data for 1955-1996 (this 42-year period was used for calibration because it corresponds to available AA distributions). The balance variables determined each day of the period are:
1. Snow accumulation and/or rain
3. Mass balance
4. Snowline altitude
5. Balance ﬂux
6. Accumulation area ratio
7. Zero balance altitude
The eight main algorithms developed to generate daily balances at each altitude interval are as follows.
I. Precipitation. Four coefficients are required to determine precipitation at the glacier.
P(i,z) glacier : β(z) p(i)gage (1)
where: P(i,z)glacier = precipitation (in m/day) at glacier on day i and altitude z; β(z) = multiplier determined by calibration; p(i)gage = observed precipitation (in m/day) at a low-altitude weather station (weighted average of two stations) on day i.
The multiplier B(z) that converts gage precipitation to high-altitude precipitation occurring over large areas is calculated with three calibration coefficients:
B(z) = (C1 - C2)(E[(z)/C3) + C2 (2)
where: Et (z) = vertical distance above glacier
terminus (m); C1, C2, C3 =
The altitude of maximum precipitation is set by limiting β from the calibrated value of C3. If β is greater than C1, then B = C1 and the value of C3 plus the terminus altitude is the altitude of maximum precipitation. C2 is then the multiplier (times gage precipitation) that determines precipitation at the terminus and C1 is the multiplier for maximum precipitation.
The fourth coefficient determines the relative weights of two potentially useable precipitation stations, and also eliminates unrepresentative or unreliable stations:
P(i)gage = P1(i) (C4) + P2 (i) (1-C4) (3)
where: p1(i), p2 (i) = precipitation at two separate stations on day i; C 4 = coefficient (fraction between 0 and 1).
2. Temperature lapse rate. Four coefﬁcients determine the daily lapse rate, using a simple algorithm that is identical to the one used for Columbia Glacier, Alaska (Tangborn 1997). The lapse rate for each day is calculated on the basis of four general conditions: cloudy/cool (mean temperature is below normal for that day); cloudy/warm (mean temperature above normal); clear skies/cool; clear skies/warm. Cloudiness is estimated from the daily temperature range, which has been shown to be an effective snow and ice melt parameter (Tangborn 1980). The mean temperature deviation is found by subtracting the daily mean (maximum plus minimum/2) from the previously calculated daily temperature averaged for the period of record.
Lr (i) = C5 (Dt(i) + C6) (on days that mean temperature is-above normal) (4)
Lr (i) = C7 (Dt(i) + C8) (on days that mean temperature is below normal) (5)
where: Lr (i) = lapse-rate (in °C per 100 m) on day i; Dt(i) : temperature range (in °C) on day i; C6, C8 = intercept when Dt(i) = 0; C5, C7 = slope of dLr/dDt.
Figure 4 shows the relationship between lapse rate and the diurnal temperature range for above and below normal temperatures.
Fig. 4. The daily lapse rate is calculated from the diurnal temperature range (daily maximum - minimum) and from the deviation of the mean temperature from the normal for that day.
3. Ablation. Four coefficients are required to calculate snow and ice ablation using both the mean daily temperature and the diurnal temperature range. The mean is applied as a turbulent heat transfer index in two ways:
An (i.z) = C9 (Tb(i.z)) (6)
Ap(i,z) = C10 (Tb (i,z)) Pr(i,z) (7)
where: An(i,z) = ablation rate on days without precipitation (m/day) at day i and altitude z; Ap(i,z) = ablation rate on days with precipitation (m/day) at day i and altitude z; Tb (i,z) = mean temperature in °C at day i and altitude z; Pr(i,z) = precipitation as rain (m/day) at day i and altitude z; C9,C10 = coefﬁcients.
The turbulent heat transfer component of ablation is based on heat balance studies conducted at the Blue Glacier (LaChapelle 1961, 1962, 1965).
A short-wave, radiative ablation factor is calculated from the diurnal temperature range, which is a rough measure of cloudiness or solar radiation (Tangborn et al. 1991, Tangborn 1994, p. 223).
Ar(i,z) = C11 (Dt(i)) (CF(z)) (8)
where: A1(i,z) = ablation rate (in m per day) at altitude z; D1(i) = diurnal temperature range (°C); C11 : calibration coefﬁcient; CF (z) = ice ablation factor:
CF(z) = C12 (1 - E(z)/S1) (9)
and: E (z) = altitude at interval z (in m); S1 = snowline altitude (in m); C12 = coefﬁcient.
If CF(z) is less than zero, it is set equal to zero; therefore the ablation factor, CF (z) is applied only below the snowline and decreases with altitude. It is reasoned that radiation is not a signiﬁcant ablation factor for highly reﬂective snow, but is strongly inﬂuential in melting low-albedo ice, and that increased amounts of dust and fewer snow patches in the ablation area reduce the albedo as altitude decreases. Total ablation each day is then the sum of At or An and Ar.
4. Snow accumulation and the snowline altitude. At each 20 m altitude interval, precipitation occurs as snow if the simulated temperature is equal to or less than 0°C, and as rain if the temperature is greater than 0°C. Rain on snow is not treated as accumulation because the snowpack on this glacier is usually not sub-freezing during times when precipitation occurs as rain.
The snowline altitude is controlled by the freezing level altitude and by the ablation rate of snow. If the freezing level falls within the glacier altitude range during the occurrence of precipitation, the snowline is lowered to that altitude. The progression of the snowline up-glacier during the spring and summer exerts considerable inﬂuence over the mass balance. A rapidly rising seasonal snowline exposes darker, low-albedo ice earlier in the season and greatly increases ablation rates. Frequent summer snowstorms do the opposite by covering low-albedo ice with a veneer of highly reﬂective snow. Snow accumulation that occurs during these storms is usually dissipated rapidly. The lowest altitude of these transient snow accumulations is designated the transient snowline. The energy of glacierization (the balance gradient at the equilibrium line altitude (ELA)) has been shown to have a strong effect on the ﬂuctuations of snowline altitude (Shumskiy 1947).
The rate that the snowline rises depends on the ablation rate of snow in the accumulation area, which is determined by two coefﬁcients. One coefﬁcient, C13, controls the rate of rise of the seasonal snowline (Equation 10), and the other, C14, the rate that the transient snowline rises (Equation 11). Visual observations of snowlines on this glacier for several years by the author indicated that there is a signiﬁcant difference in how rapidly each rises under similar weather conditions:
DSs = C13 (Asm) (10)
DST : C14 (Asm) (11)
where: DSS = seasonal snowline altitude increase (in m); DST = transient snowline altitude increase (in m); Asm = ∑ (An + AP), summed from the snowline to the glacier head; C13 = coefﬁcient corresponding to rise of seasonal snowline; C14 : coefﬁcient corresponding to rise of transient snowline.
5. Balance. Based on precipitation and temperature observations from lower elevation stations, the model calculates daily snow accumulation, ablation and the resulting balance for each deﬁned altitude zone. Each of the 28, 20 m altitude intervals is treated as a separate entity, and each interval's accumulation, ablation and resulting balance are independently calculated daily. The total glacier's balance is determined by integrating each altitude interval's balance times its area, from the terminus to the head of the glacier. Thus, each altitude interval develops its own 'summer surface' which occurs at the time its minimum balance is reached each year. The date of minimum balance occurs an average of 21 days earlier at the head of the glacier than at the terminus (the average dates are 20 September and ll October, respectively). The minimum balance for the total glacier each year is then approximately equal to the summer surface as described in standardized mass balance deﬁnitions (Meier 1962; Mayo et al. 1972).
The balance at each altitude zone z and for each day i is:
b(z,i) = c(z,i) + a (z,i) (12)
where: b(z,i) = balance at altitude interval z and on day i; c(z,i) = snow accumulation at altitude interval z and on day i; a(z,i) = ice and snow ablation at altitude interval z on day i (ablation is treated as a negative value).
A summation of each balance is made for the total glacier by
B(i)= ∑ (z,i) aa(z) (13)
where: B (i) = total glacier balance on day i; aa (z) = fraction of total surface area at altitude interval z.
The balance is summed for each altitude interval and for the total glacier, and for both the balance year (1 October - 30 September), and for the period of record.
6. Zero Balance altitude. The altitude at which b (z) = 0 is designated the zero balance altitude (ZBA), and could also be called the transient ELA. When b (Z) at E (z) is positive and b(z-1]) at E (z-I) is negative, the altitude of zero balance is found by interpolation:
ZBA(i) = (E (i,z) - b (i,z)) / (b (i,z) - b (i,z-1) (E (i,z) - E (i,z-1)) (14)
where: ZBA(i) = zero balance altitude (in m) on day i; E (i,z) = altitude at interval z (in m); b (i,z) = balance at interval z (positive); b (i,z-1)) = balance at interval (z-1) (negative); E (i,z-l) = altitude at interval (z-1).
7. Accumulation area ratio. The fraction of the area above the ZBA to total glacier area is deﬁned as the accumulation area ratio (AAR) and is calculated daily from the ﬁrst day that ice becomes exposed at the terminus through September:
AAR(i) = ∑ aa(z) (15)
where: aa (z) = fraction of total glacier area at altitude interval z; n : zero balance altitude on day i;h= glacier head.
8. Balance flux. The movement of glacier mass from higher to lower altitudes produces erosion and ultimately determines the altitude distribution of glacier area. The ﬂux of mass that occurs across the line of zero balance (ZBA) could then be considered an effective means to relate erosion and mass balance. The balance ﬂux is closely related to the energy of glacierization (Shumskiy 1947), or activity index (Meier 1961); a large value for the energy of glacierization translates to a high balance ﬂux.
For a speciﬁc day during the ablation season:
Exm(i) = Vc(i) - Va(i) (16)
where: Va(i) = volume below zero balance altitude (metres water equivalent m w.e.) (negative); VG (i) = volume above zero balance altitude (m w.e.); Exm = balance ﬂux (m w.e.) on day i.
Va(i) = ∑ aa(z)b(i,z) (17)
v,(z) =∑ aa(z)b(i,z) (18)
where: n = zero balance altitude; h = highest elevation interval (glacier head).
The optimum coefﬁcient values are found by minimizing the error (or maximizing the explained variance) that occurs when one balance variable is regressed against another in a non-linear, least squares ﬁtting process. Both linear and non-linear regressions are made: in most cases a two-degree polynomial provided a better least-square ﬁt than a linear regression.
The two ﬁtting functions used for error minimization are:
Linear Y(x) = a2 + a1X (19)
Polynomial Y(x) = a3 + a2X + a1X2 (20)
Fortran 77 program subroutines for polynomial ﬁtting are from Press et al. (1992, p. 665-670).
Nine different (but not independent) regressions are made for each day throughout the summer ablation season, starting on the day that the accumulation area ratio is less than 1.00 and ending on either 30 September or when the balance year ends. Therefore, from 15 June to 30 September, 107 regressions are made for each set of variables, for a total of 963 for each iteration of the simplex. Table 2 gives the results of linear and two-degree polynomial regressions for nine combinations of variables, for just two of the 107 days (30 June and 30 September). It should be noted that the highest R2 shown in Table 2 does not necessarily mean the set producing it is the most valuable contributor to the model's robustness.
An example demonstrating one of the nine regressions that relates daily glacier balance to the zero balance altitude is given in Fig. 5. A similar relationship was found for Vernagtferner relating measured balance to the ELA, which would correspond to the 30 September plot shown in Fig. 5 (Dyurgerov and Rienwarth 1997). The Vernagtferner study demonstrated that there are physically real connections between annual balance and the ELA.
Fig. 5. Balance as a function of zero balance altitude for 30 June and 30 September (1955-1996 period). The R2 for a polynomial regression fit is 0.85 for 30 June and 0.96 for 30 September. Each point represents either 30 June or 30 September balance versus ZBA values for one year during the 1955-1996 period. These are just two of the 963 regressions used to calibrate the model.
Each regression set usually has 42 pairs (corresponding to the 1955-1996 period); however, a year is omitted if the ablation season for that year is delayed beyond 15 June (the AAR is not less than 1.0), consequently some regressions have less than 42 pairs.
The average R2 found by ﬁtting each of the above sets for each day of the summer ablation season is based on an average of the total number of polynomial regressions. The objective error that is minimized in the simplex optimization procedure is simply the mean complement of R2,
where: 0j(r,i) = error for regression r on day i.
The objective function that is minimized in the simplex is determined by:
OBj: ∑ ∑ Oj(r,i)/9f (22)
where: f = number of days for each regression r; OBj = mean regression error used as the objective function.
The complement (1-R2) is used as the simplex's objective function, OBj, rather than the root-mean-square error simply because it is non-dimensional. The simplex ﬁnds the minimum error of a function with several variables (or in this case the maximum R2) by adapting itself to the local landscape formed by multi-dimensional space and contracting to the ﬁnal minimum (Nelder and Mead 1965).
Table 2. The nine regressions used for optimization with the resulting R2
The ﬁnal values for the 14 coefﬁcients, produced when the maximum R2 was attained, are given in Table 3. The ﬁrst four coefﬁcients determine precipitation at the glacier and all have values that produced realistic snow accumulation amounts; the distributions of winter balance with altitude closely resemble ﬁeld measurements determined by depth probing and snow pits (see Fig. 9). The ﬁnal value of C3 indicates that maximum precipitation occurs at 2060 m, or about 100 m below the head of the glacier.
Table 3. Final coefficient values after maximum R2 is attained.
Coefﬁcient C4 determines the relative weights of two precipitation records. For example, if precipitation at station p1 is relatively unrepresentative of glacier precipitation compared to station p2, the value of C4 will be near or equal to zero and station p1 eliminated for further testing. A value of 0.717 was determined for the coefﬁcient C4, when the two most representative stations were used, thus the weighted precipitation used for determining glacier precipitation is:
P(i)w=0.717 (P(i)conc )+0.283 (P(i)diab) (23)
where: P(i)w : weighted precipitation on day i; P(i)conc = observed precipitation at Concrete on day i; P(i)diab = observed precipitation at Diablo Dam on day i.
It is noteworthy that precipitation observations at or near the glacier have not been found to be useful for simulating snow accumulation on the glacier, likely because accurately catching precipitation as snow is nearly impossible at exposed, windy sites.
The next four coefficients (C5 - C8) determine the daily lapse rate, and also have reasonable values that approximately agree with theoretical lapse rates. A cloudy day with below-normal temperatures will have the lowest lapse rate of about 0.52°C per 100 m (corresponding to a saturated adiabatic lapse rate), and a day with clear skies and temperatures above normal will have the highest lapse rates of between 0.80 and 0.90°C per 100 m (corresponding to dry adiabatic conditions). The daily lapse rate between Diablo Dam and South Cascade Glacier averaged O.72°C per 100 m for the 1955-1996 period, which is slightly higher than the average global rate of 0.65°C per 100 m. These average, saturated and dry adiabatic rates approximately agree with measured lapse rates (Hartmann 1994. p. 69-70, 353-355).
Coefﬁcients C9 - C12 determine ablation, which is determined from the mean temperature and temperature range. The validity of the calculated ablation rates can be demonstrated by the relationship between one year of measurements of monthly air temperature and ablation rates in Greenland on Jacobshaven Glacier, which were found to be approximately 7.5 mm of ablation per °C (Braithwaite and Olesen 1989). In comparison, mean monthly ablation (based on Equations 6, 7 and 8) determined for South Cascade Glacier from June through September for the period 1955-1996, varies from 4.2 to 8.3 mm per °C, with an average of 6.2 mm per °C for the summer season - slightly less than that measured in Greenland.
The last two coefficients (C13 - C14) determine the rate that the snowline altitude increases. An example of snowline altitude variations for one summer season is shown in Fig. 6. The effect of summer storms that lower the freezing level altitude and snowline, and the rapid snowline rise following each storm, is apparent. The precipitous drop in the snowline on 1 September 1961 was the result of an unusually severe storm that was observed by the author and others (Springer al. 1969). The value of coefficient C14 is equal to about 11 times C13, indicating that the transient snowline daily increase in altitude is an order of magnitude greater than the seasonal snowline. The reason for this difference is likely because the transient snow cover is thinner, has a lower density and is underlain by darker ice and ﬁrn.
Fig. 6. The simulated daily snowline altitude, 1 May 30 September 1961. Note that the rise in transient snowline altitude after a summer storm is much more rapid than that of the seasonal snowline. An unusual late summer storm on 1 September caused the snowline to drop well below the glacier and essentially ended the summer ablation season that year.
The net balance is deﬁned as the change in glacier volume between two ﬁxed points, either in time or in space. For South Cascade Glacier the minimum balance each year was found to usually occur within a few days of 30 September so that using the ﬁxed date system with 30 September as the ﬁnal day does not produce large errors. Fig. 7 shows the daily variations in the continuous daily balance (Eb (Z, i)) at the terminus (Z = 2, or 1620 m) and in the accumulation area near the head of the glacier (z = 26, or 2100 m) for the 1958-1961 period. Figure 8 is the continuous glacier balance, EB(i) for the same period.
Fig. 7. The cumulative balance at the terminus (1620 m) and on the upper glacier (2100 m) for the 1958 to 1961 period.
Fig. 8. The cumulative glacier balance for the 1958-1961 period. The 1958 balance (-5.7 m) is the lowest calculated for the 1932-96 period.
The distribution of the annual, winter and summer balance with altitude (Eb (z), Ec (z), Ea (z)), averaged for the period of record, is shown in Fig. 9.
Fig. 9. The altitude distribution of the mean winter, summer and annual balances b(z), averaged for the 1955-1996 period. The interval for z is 20 m. The winter balance (3.25 m w.e.) generally agrees with the distribution and amount measured on this glacier; however, the summer balance (-4.50 m w.e.) is approximately 0.35 m (w.e.) more negative than that measured.
Although AA proﬁles before 1955 are not available, the weather records at Concrete and Diablo Dam begin in 1932. Therefore by assuming the 1955 AA proﬁle for 1932 to 1954, the annual balance and other variables can be calculated for the total 1932-1996 period. The cumulative annual balance for the full period is shown in Fig. 10. Figure 11 is winter balances and Fig. 12 is the ELA for each year of this period.
Fig. 10. Cumulative (simulated) balance for the 1932-1996 period.
Fig. 11. The simulated winter balance for the 1932-1996 period. The mean is 3.25 m (w.e.).
Fig. 12. The simulated equilibrium line altitude (ELA) for 1932-1996. The mean ELA is 1935 m for this period.
The zero balance altitude is calculated daily from the time the ﬁrst ice is exposed at the terminus until the glacier is again covered by snow the following winter. The summer balance season ranges from 59 days (beginning on 2 August) to 182 days (beginning on 1 April) and averages 129 days (corresponding to a 24 May start date for the summer season). The mean value of ZBA on 30 September is approximately equal to the ELA. The average for the period of record is 1930 m, which is about 40 m higher than the observed ELA.
The relationship between the zero balance altitude (or ELA) and the snowline altitude is important for understanding the complex interactions between climate ﬂuctuations and the ELA (Kuhn 1989). For unknown reasons, the simulated snowline altitude at South Cascade Glacier is nearly always greater than the ZBA (intuition suggests they should be equal); however, visual observations of glaciers usually place the ELA at the snowline altitude.
The results for determining optimum coefﬁcient values from the maximum average R2 are demonstrated in Fig. 13. The R2 produced after 548 iterations is 0.801, which is an average of 657 polynomial regressions for each iteration. The basis for using a statistical veriﬁcation technique is that an internal consistency is thought to exist among the various parameters that are generated by combining meteorological data with the glacier's area altitude proﬁle. Obviously, for this approach to work, the meteorological records must be from stations that are close enough to sample the same air masses and weather conditions as the glacier.
Fig. 13. Results of the simplex optimization process to obtain the maximum R2. Each point is one iteration and represents from 828 to 963 polynomial regressions (representing 15 June and 1 July starting dates of the summer season). There are 548 iterations and the final average R2 is 0.801.
Comparison with measured balance
A rigorous test of the model's accuracy and reliability is by direct comparison with measured balances produced by South Cascade Glacier ﬁeld surveys. Continuous mass balance data have been collected at this glacier since 1959 (Meier and Tangborn 1965; Krimmel 1994, 1995, 1996, 1997, 1998). These measurements, though labor-intensive and expensive, have been considered the most accurate method for determining a glacier's annual mass balance. A number of analyses have been performed comparing these direct measurements with the geodetic method of determining volume change (by calculating surface elevation changes between mappings). Many of these studies have shown signiﬁcant differences between the direct and geodetic methods (Krimmel 1989).
A year-by-year comparison of simulated and measured balances suggests that there are substantial errors in one or both methods. Figure 14 is a plot of annual measured versus simulated balances for the 1959-1996 period. The explained variance fraction for a linear regression ﬁt is 0.57 with a standard error of 0.75 m w.e).
Fig. 14. The measured versus simulated annual balance for the 1959-1996 period. The R2 for a linear regression is 0.57 and the root mean square error of simulated balance is 0.75 m.
A plot of cumulative balances (measured and simulated) in Fig. 15 indicates that the simulated balance is consistently more negative than measured, chieﬂy due to much greater ice ablation on the lower glacier. The geodetic method used for calculating the volume change between 1975 and 1996 equaled -24 m , compared with a direct measurement of annual balances of -18 m (Krimmel 1999). The model's simulated volume change for this period was also -24 m, indicating that the model and geodetic results are in good agreement. A plot of simulated volume change versus the mean R2 determined by coefﬁcient optimization demonstrates that the -24 m volume change determined by the model is valid and real (Fig. 16).
Fig. 15. Cumulative simulated and measured balance for the 1959-1996 period. The 1975-1996 measured volume change is -18 m (w.e.) and the simulated loss for the same period is -24 m (w.e.), which is equal to the geodetic balance measured for this period.
Fig. 16. The relationship between the R2 for the iterations shown in Fig. 13, and the corresponding 1975-1998 average thickness change calculated for each set of trial coefficients. The thickness change when the maximum R2 was reached is approximately -24 m (w.e.). Note: the measured balance was not used for calibration to obtain the final coefficient values.
Linking climate, mass balance and the AA proﬁle
Intuition suggests that it would be impossible to mathematically relate one physical variable to another and produce realistic results when both are derived from the same data source. A regression of X versus Y if both are based on Z could be manipulated to produce a perfect ﬁt (R2 = 1), but predicting Y from X for values outside the sample range would produce meaningless results. The relationships shown in Table 3 among ﬁve variables generated within the model are all interdependent. However, a likely signiﬁcant ﬁnding is that these regressions tend to counter one another, so that as the R2 for one increases (by altering one or more coefﬁcients), it will tend to decrease for one or more of the others. Therefore, the ﬁnal value of the objective function that is minimized in the simplex is a balance of all regressions that are used. The reason why these regressions generate mass balance results that appear to be physically real is locked in the glacier's area-altitude distribution, which exerts a robust control over the mass balance and related variables. Each glacier's AA proﬁle is unique so that the application of a single set of meteorological observations to several glaciers in a deﬁned region will produce a wide variation in mass balances. The differences between them reﬂect the interaction of the climate and the individual AA proﬁles.
It is proposed that thousands or millions of years of glacial erosion have embedded in the AA proﬁle a key to the climate that formed and now nourishes and sustains the glacier. Glaciated landscapes over much of the earth clearly demonstrate the close link between climate and glacier erosion (Brozovic et al. 1997; Pinter and Brandon 1997). The mass balance (or mass turnover) is an inﬂuential element that relates glacier erosion rates to the climate, as has been demonstrated for glaciated areas in Scotland and Scandinavia (Benn and Evans, 1998, p. 372-375). The dependence of glacial erosion on mass balance has been discussed earlier by others (Andrews 1972; Chemova 1981; Sugden and John 1976).
There can be little doubt that an inextricable climate/erosion/mass balance link exists, and that it is the geometric conﬁguration of the glacier (its AA profile) that connects them. If it were not for the uniqueness of each glacier's area-altitude distribution and its link to the past climate, relating the generated variables to each other to calculate the mass balance (the basis for the PTAA model) would produce erroneous, or even more likely, absurdly nrealistic, results.
The results provided in this study suggest that realistic mass balances can be generated by this approach for any glacier for which reliable meteorological observations from nearby weather stations are available (the maximum permissible distance from the glacier to weather station will vary with the terrain and type of weather system, but is estimated to be about 150 km). The weather data need not be from stations at altitudes similar to that of the glacier - in fact, it is more desirable if the precipitation data were collected at low altitudes because of the notorious inaccuracy of precipitation measurements at exposed, high-altitude sites in mountainous regions. It should be noted that the ﬁnal results for South Cascade Glacier would he identical regardless of the availability of measured balances. Compared to the cost of balance ﬁeld measurements or volume changes by geodetic means, the expense of applying this model is miniscule. In addition, a calibrated PTAA model could be used to predict a g1acier's advance or retreat by increasing or decreasing temperature and/or precipitation. An improved understanding of how climate and glaciers are related is possible if the connection between erosion, mass balance and climate could be further clariﬁed.
My thanks to Sara Zhang and Andrew Tangborn for their invaluable assistance with Numerical Recipes for Fortran 77, and to Bob Krimmel who supplied all the geodetic data and many helpful comments. Funding for this study was provided by HyMet Inc..
Wendell l/. Tangborn,
Andrews, J. T., 1972: Glacier power, mass
balances, velocities and erosion potential.
Zeitschrzflfur Geomorphologie, Suppl, Bd. 13:1-17.
Benn, D.I. and Evans, J.A., 1998: Glaciers and Glaciation. Arnold. London. 734 p.
Braithwaite, R.J. and Olesen, O.B., 1989; Calculation of glacier ablation from air temperature, West Greenland. In: Oerlemans, J. (ed.): Glacier Fluctuations and Climatic Change. Kluwer Academic Publishers. Dordrecht. (417 p.): 219-233.
Brozovic, N.,.Burbank, D.W. and Meigs, A.J., 1997: Climatic limits on landscape development in the Northwestern Himalaya. Science, 276: 571-574.
Chernova, L.P., 1981: Inﬂuence of mass balance and run-off on relief-forming activity of mountain glaciers. Annals of Glaciolagy, 2: 69-70.
Dyurgerov, M. and Reinwarth, 0., 1997: Mass
balance and runoff reconstruction for
Vernagferener Glacier since 1889. Data of Glaciological Studies, Publication 81, Tashkent Symposium '93. 114-116.
Hartman, D.L., 1994: Global Physical Climatology. Academic Press. Sa.n Diego. 411 p.
Krimrnel, R.M. , 1989: Mass balance and volume
of South Cascade Glacier, Washington, 1958-1985. In: Oerlemans, J.
(ed.): Glacier Fluctuations and Climatic Change. Kluwer Academic
Publishers. Dordrecht. (417 p.): 193-206.
- 1994: Runoff, precipitation, mass balance, and ice velocity measurements at South Cascade Glacier, Washington, 1993 balance year. US Geological Survey Water Resources Investigation Report 94-4139. 34 p.
- 1995: Water, ice and meteorological measurements at South Cascade Glacier, Washington, 1994 balance year. US Geological Survey Water Resources Investigation Report 95-4162. 41 p.
- 1996: Water, ice and meteorological measurements at South Cascade Glacier, Washington, 1995 balance year. US Geological Survey Water Resources Investigation Report 96-4174. 37 p.
- 1997: Water, ice and meteorological measurements at South Cascade Glacier, Washington, 1996 balance year. U.S. Geological Survey Water Resources Investigation Report 97-4143. 34 p.
- 1998: Water, ice and meteorological measurements at South Cascade Glacier, Washington, 1997 balance year. US Geological Survey Water Resources Investigation Report 98-4090. 30 p.
- 1999: Mass balance and volume of South Cascade Glacier, Washington. Geografiska Annaler, 81 A: 653-658.
Kuhn, M., 1989: The response of the Equilibrium Line Altitude to climate fluctuations, theory and observations. In: Oerlemans, J. (ed): Glacier Fluctuations and Climatic Change. Kluwer Academic Publishers. Dordrecht. (417 p.): 407-117.
Kuhn, M., Markl,G,, Kaser, G., Nickus, U., Obleitner, F. and Schneider, H., 1985: Fluctuations of climate and mass balance; different responses of two glaciers. Proceedings of Symposium on Climate and Paleoclimate of Lakes, Rivers and Glaciers, Igls, Austria, 4-7 June, 1984. 409-116.
LaChapelle, E., 1961: Annual Mass and Energy
Exchange on the Blue Glacier. Journal of Glaciology, 4:
- 1962: Assessing glacier mass budgets by reconnaissance aerial photography. Journal of Glaciology, 4: 290-297.
- 1965: The mass budget of Blue Glacier, Washington. Journal of Glaciology, 5: 609-623.
Mayo, L.R., Meier, M.F. and Tangborn, W.V., 1972: A system to combine stratigraphic and annual balance systems: a contribution to the International Hydrological Decade. Journal of Glaciology, 11(61): 3-14.
Meier, M.F.,1961: Mass budget of South Cascade
Glacier, 1957-60. US Geological survey Professional Paper 424B.
- 1962: Proposed definitions for glacier mass budget terms. Journal 0fGlc1ciology, 4: 252-261.
Meier, M.F. and Tangborn, W.V., 1965: Net budget and ﬂow of South Cascade Glacier, Washington. Journal of Glaciology, 5: 547-566.
Miller, C.D., 1969: Chronology of Neoglacial Moraines in the Dome Peak Area, North Cascade Range, Washington. Arctic andAlpine Research, 1: 49-66.
Nelder, J./1. and Mead, R., 1965: A simplex
method for function minimization. Computer
Journal, 7: 308-312.
Pinter, N. and Brandon, M ., 1997: How erosion
builds mountains. Scientiﬁc American, 276:
Press, W.H., Teukoisky, S.A., Vetterling, W,T. and Flannery, B.P., 1992: Numerical Recipes in Fortran 77. (2nd edition). Cambridge University Press. 923 p.
Rasmussen, L.A. and Tangborn, W.V., 1976: Hydrology of the North Cascades region, Washington-Part 1. Runoff, precipitation, and storage characteristics. Water Resources Research, 12: 187-202.
Shumskiy, P.A., 1947: The energy of
glacierization and life of glaciers. United States
Publishing House for Geographical Literature. Translated by U.S. Snow, Ice and Permafrost Research Establishment, Translation 7, 1950). 27 p.
Spring, R., Spring, I. and Manning, H., 1969:
The North Cascades National Park. Superior
Publishing Company. Seattle. (143 p.): 35-11.
Sugden, D.E. and John, B.S., 1976: Glaciers and Landscape. Edward Arnold. London. 376 p.
Tangborn, W.V., 1968: Mass balances of some
North Cascade glaciers as determined by hydrologic parameters,
1920-1965. Proceedings of the International Commission of Snow and
Ice, Bern, Switzerland, 1967, International Association of
Scientific Hydrology Publication 79. 267-274.
- 1980: Two models for estimating climate glacier relationships in the North Cascades, Washington, U.S.A. Journal of Glaciology, 25: 3-20.
- 1994: Recent changes in the melt-rate patterns of South Cascade Glacier, Washington (Abstract). American Geophysical Union. Fall Meeting Proceedings. San Francisco, CA.
- 1997: Using low-altitude meteorological observations to calculate the mass balance of A1aska's Columbia Glacier and relate it to calving and speed. In: Van der Veen, C.J. (ed.): Calving Glaciers: RepOI"t of a Workshop, February 28-March 2, 1997. BPRC Report No. 15, Byrd Polar Research Center. The Ohio State University, Columbus, Ohio. (194 p.): 141-161.
Tangborn, W.V. and Rasmussen, L.A., 1976: Hydrology of the North Cascades region, Washington-Part 2. A proposed hydrometeorological streamﬂow prediction method. Water Resources Research, 12: 203-216.
Tangborn, W. V. and Rasmussen, L.A,, 1977:
Application of a hydrometeorological model to
South-Central Sierra Nevada of California. U.S. Geological Survey Journal ofResearch, 5: 33-48.
Tangborn, W.V., Krimrnel, R.M. and Meier, M.F., 1975: A comparison of glacier mass balance measurements by glaciologic, hydrologic and mapping methods, South Cascade Glacier, Washington. Proceedings, International Association of Hydrological Science, Moscow General Assembly, U.S.S.R., 1971. IAHS-AISH Publication 104. 185-196.
Tangborn, W.V., Fountain, A.G. and Sikonia, W.D., 1990: Effect of altitude-area distribution on glacier mass balance - A comparison of North and South Klawatti Glaciers, Washington State, U.S.A. Annals of Glaciol0gy, 14: 278-282.
Tangborn, W.V., Ebbesmeyer, C., Lachapelle, E.R,, 1991: Hidden signals in the Washington State climate record. Presented at Puget Sound Research Conference, January 4-5, 1991, Seattle, WA.
Manuscript received August I998, revised and accepted March 1999.