Open Collections

UBC Theses and Dissertations

UBC Theses Logo

UBC Theses and Dissertations

The relative importance of calving and surface ablation at a lacustrine terminating glacier : a detailed… Chernos, Matthew 2014

Your browser doesn't seem to have a PDF viewer, please download the PDF to view this item.

Item Metadata


24-ubc_2014_september_chernos_matthew.pdf [ 24.1MB ]
JSON: 24-1.0166946.json
JSON-LD: 24-1.0166946-ld.json
RDF/XML (Pretty): 24-1.0166946-rdf.xml
RDF/JSON: 24-1.0166946-rdf.json
Turtle: 24-1.0166946-turtle.txt
N-Triples: 24-1.0166946-rdf-ntriples.txt
Original Record: 24-1.0166946-source.json
Full Text

Full Text

The Relative Importance of Calvingand Surface Ablation at a LacustrineTerminating GlacierA Detailed Assessment of Ice Loss at Bridge Glacier,British ColumbiabyMatthew ChernosB.Sc., The University of British Columbia, 2010A THESIS SUBMITTED IN PARTIAL FULFILLMENT OFTHE REQUIREMENTS FOR THE DEGREE OFMASTER OF SCIENCEinThe Faculty of Graduate and Postdoctoral Studies(Geography)THE UNIVERSITY OF BRITISH COLUMBIA(Vancouver)August 2014c© Matthew Chernos 2014AbstractBridge Glacier is a lacustrine calving glacier located in the southern Coast Mountainsand terminates in a 6.2 km2 proglacial lake. The glacier has retreated more than 3.55km up-valley since 1984, the majority of the retreat having occurred since 2003. Whilesurface melt may have contributed to the retreat, calving allowed for an additional annualvolume of ice loss. The relative contributions from surface melt and calving to the totalvolume of ice loss is examined for the 2013 melt season. Surface melt is quantified usingon-glacier meteorological data to drive a distributed energy balance model. The calvingflux is quantified using field measurements of lake bathymetry, terminus area change, andice thickness. Calving flux estimates are completed by daily measurements of terminussurface velocity derived from manual feature tracking using oblique time lapse cameraimagery. Calving accounts for 23% of the total ice loss in the 2013 melt season, suggestingthat surface melt is the main driver of mass loss at Bridge Glacier.Data from the 2013 field season is used to inform historical calving flux and surface meltestimates from 1984 to 2013. The calving flux is minor until 1991, at which point the glacierterminus achieves flotation, and begins to discharge large tabular icebergs. Calving wascharacterized by large, multi-annual retreats, alternating with periods of relative stability.The calving flux peaked from 2005 to 2010, when it was roughly equal to the mass lossdue to surface melt. Calving was a much smaller contributor of mass loss from BridgeGlacier, except for a transient high-calving period in the late 2000s. Looking forward,Bridge Glacier will retreat into shallower water where the terminus will no longer float,and calving losses should decrease substantially. Although calving losses will become anincreasingly minor portion of the mass balance, future retreat is expected at Bridge Glacierdue to a legacy of dynamic thinning brought about by its transient calving phase.iiPrefaceThis research was made possible with the help of many people. Technical and academicadvice from Dr. Michele Koppes and Dr. R. Dan Moore significantly improved fieldworkplanning, focus, and the writing and content of this thesis. Alistair Davis and MelanieEbsworth assisted with field work. Data collection, field work, and data analyses wereimproved by advice and support from Lawrence Bird.A paper is currently being written that encompasses Chapters 3, 5, and 6. It is cur-rently titled “The Relative Contributions of Calving and Surface Ablation to Ice Loss at aLacustrine Terminating Glacier” and is co-authored by M. Koppes and R.D. Moore.iiiTable of ContentsAbstract . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . iiPreface . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . iiiTable of Contents . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . ivList of Tables . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . viiiList of Figures . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . ixAcknowledgements . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xiDedication . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xii1 Introduction: Calving Glaciers and Their Relationship to Climate . . 11.1 Motivation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11.2 Glacier Mass Balance . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 21.3 Glacier Flow . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 41.4 Glacier Calving . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 51.4.1 Calving Mechanisms . . . . . . . . . . . . . . . . . . . . . . . . . . . 61.5 An Intricately Coupled System . . . . . . . . . . . . . . . . . . . . . . . . . 91.6 Research Gaps and Study Goals . . . . . . . . . . . . . . . . . . . . . . . . 111.7 Study Location . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 121.8 Thesis Structure . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 142 A First Order Estimate of the Influence of Climate on Retreat . . . . . 162.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 162.2 Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 172.3 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 182.3.1 Retreat and Climate . . . . . . . . . . . . . . . . . . . . . . . . . . 182.3.2 Climate . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 192.3.3 Inverse Linear Model . . . . . . . . . . . . . . . . . . . . . . . . . . 222.4 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 242.5 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 26ivTable of Contents3 A Distributed Energy Balance Model for Ice Ablation . . . . . . . . . . 283.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 283.1.1 Surface Energy Balance Modelling . . . . . . . . . . . . . . . . . . . 283.2 Instrumentation and Field Methods . . . . . . . . . . . . . . . . . . . . . . 313.3 Energy Balance Modelling Methods . . . . . . . . . . . . . . . . . . . . . . 333.3.1 Net Radiation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 343.3.2 Turbulent Heat Fluxes . . . . . . . . . . . . . . . . . . . . . . . . . 373.3.3 Snowline Retreat . . . . . . . . . . . . . . . . . . . . . . . . . . . . 443.3.4 Complementary Melt Rates and Meteorology . . . . . . . . . . . . . 453.4 Results and Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 463.4.1 Volumetric Ice Melt . . . . . . . . . . . . . . . . . . . . . . . . . . . 463.4.2 Relative Energy Contributions . . . . . . . . . . . . . . . . . . . . . 463.4.3 Modelled Glacier Temperature Distribution . . . . . . . . . . . . . . 483.4.4 Comparative Model Performance . . . . . . . . . . . . . . . . . . . 483.4.5 Comparison of FPL Model Uncertainty . . . . . . . . . . . . . . . . 523.4.6 Climatic Sensitivity . . . . . . . . . . . . . . . . . . . . . . . . . . . 553.4.7 Model Uncertainty . . . . . . . . . . . . . . . . . . . . . . . . . . . 583.5 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 584 Measuring Glacier Velocity with Time Lapse Imagery . . . . . . . . . . 604.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 604.2 Study Site and Instrumental Setup . . . . . . . . . . . . . . . . . . . . . . 624.3 Measuring Glacier Velocity from Time Lapse Images . . . . . . . . . . . . . 654.3.1 Tracking Features . . . . . . . . . . . . . . . . . . . . . . . . . . . . 654.3.2 Converting pixel displacement into distance displacement . . . . . . 674.3.3 Calculating the horizontal distance from camera . . . . . . . . . . . 674.3.4 Correcting from non-perpendicular flow . . . . . . . . . . . . . . . . 704.3.5 Data Processing . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 704.4 Results and Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 714.4.1 Ground Control Points . . . . . . . . . . . . . . . . . . . . . . . . . 714.4.2 Terminus Velocity . . . . . . . . . . . . . . . . . . . . . . . . . . . . 724.4.3 Nunatak Velocity . . . . . . . . . . . . . . . . . . . . . . . . . . . . 734.5 Error Analysis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 774.5.1 Human Measurement Error and Camera Shift . . . . . . . . . . . . 784.5.2 Photographic and Feature Quality . . . . . . . . . . . . . . . . . . . 794.5.3 The Geometry of Control Points . . . . . . . . . . . . . . . . . . . . 794.5.4 Comparison with Surveyed Results . . . . . . . . . . . . . . . . . . 814.6 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 81vTable of Contents5 The Magnitude and Timing of Calving . . . . . . . . . . . . . . . . . . . . 835.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 835.2 Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 855.2.1 Terminus Area Change . . . . . . . . . . . . . . . . . . . . . . . . . 855.2.2 Terminus Flow Velocity . . . . . . . . . . . . . . . . . . . . . . . . . 865.2.3 Bathymetry and Ice Thickness Estimate . . . . . . . . . . . . . . . 865.2.4 Limnology and Calving Events . . . . . . . . . . . . . . . . . . . . . 885.3 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 895.3.1 Calving Events . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 895.3.2 Calving Correlations . . . . . . . . . . . . . . . . . . . . . . . . . . 935.3.3 Calving Flux . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 955.3.4 Uncertainty Analysis . . . . . . . . . . . . . . . . . . . . . . . . . . 975.4 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 975.4.1 Potential Calving Mechanisms . . . . . . . . . . . . . . . . . . . . . 975.4.2 Calving Forcings . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 995.4.3 The Relative Importance of Calving . . . . . . . . . . . . . . . . . . 1025.5 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1036 Calving, Retreat, and Surface Melt: Present and Historical Conditions 1056.1 The Relationship Between Calving and Climate . . . . . . . . . . . . . . . 1056.2 Conditions at Bridge Glacier - 2013 . . . . . . . . . . . . . . . . . . . . . . 1066.2.1 Ablation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1066.2.2 Glacier Flow . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1066.2.3 Calving at Bridge Glacier . . . . . . . . . . . . . . . . . . . . . . . . 1076.3 Relative Importance and Long Term Implications . . . . . . . . . . . . . . 1096.3.1 Relative Importance of Calving at Bridge Glacier 2013 . . . . . . . 1096.3.2 Long-Term Mass Loss Fluctuations . . . . . . . . . . . . . . . . . . 1106.3.3 Bridge Glacier Relative to Other Lacustrine Calving Systems . . . . 1136.3.4 Future Implications . . . . . . . . . . . . . . . . . . . . . . . . . . . 1157 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1177.1 Major Findings . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1177.2 Future Work . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 118Bibliography . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 119AppendicesA Distributed Energy Balance Modelling Data Processing . . . . . . . . . 127A.1 Pressure Transducer Post-Processing . . . . . . . . . . . . . . . . . . . . . 127A.2 Solar Geometry Calculations . . . . . . . . . . . . . . . . . . . . . . . . . . 128viTable of ContentsA.3 Modelled Data Gaps . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 130A.3.1 Incoming Longwave Radiation Gaps . . . . . . . . . . . . . . . . . . 130A.3.2 Incoming Wind Speed and Direction Gaps . . . . . . . . . . . . . . 131B Time Lapse Camera Methodology and Walkthrough . . . . . . . . . . . 132B.1 Getting Started . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 132B.2 Importing Pictures . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 132viiList of Tables1.1 Some notable statistical relations for mass balance in western North America 31.2 Characteristics of selected major lacustrine calving glacier studies . . . . . . 102.1 Linear inverse model parameters . . . . . . . . . . . . . . . . . . . . . . . . 223.1 Measured variables from on-site automatic weather stations . . . . . . . . . 323.2 Variables used in this study . . . . . . . . . . . . . . . . . . . . . . . . . . . 343.3 Statistical coefficients used in FPL modelling . . . . . . . . . . . . . . . . . 414.1 Time lapse camera variables . . . . . . . . . . . . . . . . . . . . . . . . . . . 654.2 Terminus tracked points . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 724.3 Nunatak centreline velocity time series . . . . . . . . . . . . . . . . . . . . . 754.4 Nunatak Camera distributed points . . . . . . . . . . . . . . . . . . . . . . . 754.5 Comparison of GPS surveyed stakes and time lapse camera derived studyperiod displacements . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 816.1 Characteristics of selected major lacustrine calving glacier studies with addedBridge Glacier variables . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 114viiiList of Figures1.1 Regional variations in the calving rate and water depth . . . . . . . . . . . 61.2 Calving along a ‘preferential line of weakness’ . . . . . . . . . . . . . . . . . 71.3 The force imbalance at the calving front . . . . . . . . . . . . . . . . . . . . 81.4 The relative importance of calving controls related to glacier flow speed . . 91.5 Study site location . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 132.1 Bridge Glacier retreat from Landsat imagery . . . . . . . . . . . . . . . . . 182.2 Winter accumulation predictors . . . . . . . . . . . . . . . . . . . . . . . . . 192.3 Summer melt predictors . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 202.4 Bridge Glacier climate and retreat . . . . . . . . . . . . . . . . . . . . . . . 212.5 Modelled and observed retreat . . . . . . . . . . . . . . . . . . . . . . . . . 242.6 Calving-climate feedback loop . . . . . . . . . . . . . . . . . . . . . . . . . . 263.1 Bridge Glacier 2013 study area and instrumentation . . . . . . . . . . . . . 313.2 On-glacier instrumentation photos . . . . . . . . . . . . . . . . . . . . . . . 333.3 Flow Path Lengths for Bridge Glacier . . . . . . . . . . . . . . . . . . . . . 383.4 Temperature depression as a function of ambient off-glacier temperature . . 393.5 Wind speed as a function of ambient air temperature . . . . . . . . . . . . . 423.6 Modelled snowline retreat from Landsat data . . . . . . . . . . . . . . . . . 453.7 Modelled ice melt for the study period . . . . . . . . . . . . . . . . . . . . . 473.8 Average daily energy fluxes contributing to melt within the ablation area . 483.9 Air temperature distributions for 13:00, July 31, 2013 . . . . . . . . . . . . 493.10 Cumulative melt from DEBM, point EBM and measured sources . . . . . . 503.11 Comparison between Point EBM and DEBM predictions at Glacier AWS . 513.12 Cumulative ice melt for DEBM and ablation stakes . . . . . . . . . . . . . . 523.13 Comparison of DEBM turbulent heat fluxes with more frequently employeddistributions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 533.14 Changes in turbulent heat fluxes using different models for spatially dis-tributing temperature and wind speed . . . . . . . . . . . . . . . . . . . . . 553.15 Calving loses relative to climatic variability scenarios affecting the melt rate 563.16 Predicted additional melt from a 1◦C increase in air temperature . . . . . . 574.1 Time Lapse camera locations . . . . . . . . . . . . . . . . . . . . . . . . . . 62ixList of Figures4.2 Photos of time lapse camera and velocity stake . . . . . . . . . . . . . . . . 634.3 Time lapse camera views and tracked points . . . . . . . . . . . . . . . . . . 644.4 Time lapse image processing workflow . . . . . . . . . . . . . . . . . . . . . 664.5 Camera geometry for converting pixel displacement into pixel values . . . . 684.6 Camera depth geometry . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 694.7 Distance from camera as a function of the height above glacier surface andzenith angle . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 704.8 Displacements for stationary Ground Control Points . . . . . . . . . . . . . 714.9 Terminus velocity time series . . . . . . . . . . . . . . . . . . . . . . . . . . 734.10 Nunatak centerline velocity time series . . . . . . . . . . . . . . . . . . . . . 744.11 Nunatak Camera lateral velocity time series . . . . . . . . . . . . . . . . . . 764.12 Tracked points from the Nunatak Camera . . . . . . . . . . . . . . . . . . . 774.13 Lateral profile from Nunatak Camera . . . . . . . . . . . . . . . . . . . . . . 784.14 Ground Control Point standard deviations . . . . . . . . . . . . . . . . . . . 794.15 Standard deviations for 8 tracked points from each time lapse camera . . . 805.1 Bird’s eye view of a conceptual model of calving flux . . . . . . . . . . . . . 855.2 Photograph showing terminus inflection point indicating flotation . . . . . . 875.3 Processed and unprocessed water level data . . . . . . . . . . . . . . . . . . 895.4 Large tabular calving event, June 23, 2013 . . . . . . . . . . . . . . . . . . . 905.5 Moderate calving event, August 2, 2013 . . . . . . . . . . . . . . . . . . . . 915.6 Large calving event, August 21, 2013 . . . . . . . . . . . . . . . . . . . . . . 925.7 Measured limnological, meteorological and glaciological variables at BridgeGlacier . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 945.8 Bathymetry for Bridge Lake . . . . . . . . . . . . . . . . . . . . . . . . . . . 965.9 Total volumetric mass loss at Bridge Glacier during the study period . . . . 975.10 Tilting ice front and exposed meltwater notches . . . . . . . . . . . . . . . . 985.11 Photograph of a rising ice front . . . . . . . . . . . . . . . . . . . . . . . . . 995.12 Photograph of terminus melt notches . . . . . . . . . . . . . . . . . . . . . . 1006.1 Photograph of a terminus force imbalance . . . . . . . . . . . . . . . . . . . 1086.2 Modelled mass balance gradients . . . . . . . . . . . . . . . . . . . . . . . . 1116.3 Historical ice loss from calving and surface melt . . . . . . . . . . . . . . . . 1126.4 Bridge Glacier from above . . . . . . . . . . . . . . . . . . . . . . . . . . . . 115A.1 Difference between maximum and minimum depth readings for each day forthe pressure transducer. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 128A.2 Pressure transducer observed melt rate plotted against Distributed EnergyBalance Model melt rate . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 129B.1 Memory in use . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 133B.2 Manually increase the memory size . . . . . . . . . . . . . . . . . . . . . . . 133xAcknowledgementsAlthough only my name appears on this thesis, what is written would not have beenpossible without the help of many people. I would like to thank Michele Koppes fortaking me on as a student, and for guiding me through two exciting years of science anddiscovery, from field work, to computer work, to actually trying to make sense of, andcommunicate, the findings. I would like to also thank Dan Moore for providing anotherguiding voice in the process, giving a broader context to the work, new ideas, and manyhours of guitar playing. This work was supported by a Discovery Grant from NaturalSciences and Engineering Research Council (NSERC) to Dr. Koppes and Dr. Moore aswell as a Graduate Support Initiative scholarship to M. Chernos.I would like to thank the Geography Department and its faculty for not only two greatyears as a master student, but four more as an undergraduate. The courses and peoplehere have taught me to think critically, and have given me a new appreciation and passionfor the science that governs our natural world. I would like to thank Dr. Brett Eaton,who initially convinced me after a Geography field course, that I should be pursuing aBSc. (not a BA. for which I was originally registered), and helped navigate the innerworkings of the university to make it happen. I would also like to thank the professors forwhom I was a teaching assistant (Drs. Khaled Hamdan, Simon Donner, Greg Henry, MikeBovis, Andreas Christen, and Dan Moore) who shared with me their passion for teaching,and their field. Finally, I would like to thank Mike Church his guidance, passion, andperspective on science and research through work on the Fraser River and beyond.I would also like to acknowledge Lawrence Bird, for his support navigating BridgeGlacier field work, Geography Building construction, R tutorials, field gear construction,and the spatial distribution of UBC pubs. I would also like to thank Alistair Davis for hisfield work acumen, great field stories, and pretty clever Glacier AWS tripod design.I would like to thank my family for their support and encouragement, and for givingme perspective on my journey to become the Master of Glaciers (or something like that).Finally, I would like to thank Melanie Ebsworth, who has graciously put up with manyhours of typing and incessant conversations about glaciers. Her support and patience hasnever been more evident than when our hikes get side-tracked to get “just another view”of a glacier.To my family and friends, colleagues, and supervisory committee:Thank you, and happy trails!xiTo Bridge Glacier,and everyone else’s Bridge Glacier,wherever that may be.xiiChapter 1Introduction: Calving Glaciers andTheir Relationship to Climate1.1 MotivationSince the 19th century Little Ice Age, glaciers across the globe have been shrinking at anaccelerated rate. Although this retreat has been irregular, a general trend of 20th centuryretreat is pervasive, and well correlated with an increase in global mean temperatures (Oer-lemans, 2005). The reduction in glacial area has raised considerable concern about potentialchanges in the timing, volume, and duration of summer streamflow (Marshall et al., 2011;Stahl et al., 2008). These changes have major implications for hydroelectric projects, agri-culture, and water quantity and quality for general consumption. Similarly, glacier retreat,and subsequent hydrological changes, also has significant ramifications for eustatic sea levelrise (Dyurgerov et al., 2002; Barry, 2006). While recent glacier retreat is well documented(Kaser et al., 2006), the projection of future glacier retreat is critical to the managementof water resources, projection of future channel morphology (Cowie et al., 2014), and theestablishment of ecological and aquatic habitats (Milner and Bailey, 1989).Due to their direct relationship with air temperatures and precipitation, glaciers serveas important high altitude weather stations (Kaser et al., 2006). As barometers of achanging climate, fluctuations in glacier mass balance, area, or length can be related to airtemperature (Oerlemans, 2005), supplementing a global network of temperature records,much of which is biased towards lower elevation sites (Christy et al., 2000; Oerlemans, 2005).Finally, glacier fluctuations provide an additional record of a changing global climate,independent of other climatic indicators such as tree rings, sediment cores, and isotoperecords (Oerlemans, 2005; Dyurgerov et al., 2002).Glacier fluctuations are predominantly a product of changing climatic conditions, al-though glaciers that terminate in bodies of water respond at least partially independently11.2. Glacier Mass Balanceof climate on decadal timescales (Warren and Kirkbride, 2003; Post et al., 2011). This de-coupling of the climate-glacier signal is due to the fact that calving, the process by whichicebergs are mechanically discharged from the glacier terminus, can be an important addi-tional source of ice loss (Benn et al., 2007a). While calving glaciers have a more complexclimatic signal than one from glaciers that terminate on land (Van der Veen, 2002; Motykaet al., 2003), their inherent instability suggests that they have the potential to contributedisproportionately to eustatic sea level rise (Meier and Post, 1987; Dyurgerov and Meier,2005), making them an important ‘grey area’ in our understanding of how glaciers respondto climate.While previous work has focused on the mechanisms and rate of calving in both marineand lacustrine glacier systems (Benn et al., 2007b; Post et al., 2011), more work is neededto assess the role and importance of calving on glacier mass balance. The goal of this studyis to quantify the contributions of ice loss from calving and surface melt over a single meltseason. This study then aims to use field data from the melt season to inform estimatesof how the relative importance of calving and surface melt has changed over the transientcalving phase of a lacustrine calving glacier.1.2 Glacier Mass BalanceBefore 1950, there were few in situ observations of glacier mass change, while before 1900,direct glacier observations were extremely sparse (Rasmussen and Conway, 2004; Oerle-mans, 1994), limited to anecdotal evidence and artwork (Nussbaumer et al., 2007). It wasonly with the inception of the International Hydrological Decade in 1965 that a standardmethodology was instituted to quantify glacial change in Canada, and consistent observa-tions began in the southwestern Canadian Cordillera. Today there are continuous recordsfrom 1965 to present of Peyto Glacier in the Rocky Mountains, as well as similar lengthrecords from Place Glacier and Helm Glacier in the southern Coast Mountains (Demuthet al., 2006; Moore and Demuth, 2001). These records complement a long-term datasetat South Cascade Glacier in Washington state, USA, from 1959 to present (Rasmussenand Conway, 2004; Beedle et al., 2009), among other North American datasets in Al-berta, northern British Columbia, and Alaska (Pelto, 1987; Rasmussen and Conway, 2004;Dyurgerov and Meier, 2005; Braithwaite, 2002).Glacier change is quantified by measuring the net mass balance, a glacier’s annualchange in volume, as measured in units of water equivalent melt. Mass balance can be21.2. Glacier Mass Balancepartitioned into seasonal balances, where the net balance for the year (bn) is equal to thesum of the winter (bw) and summer (bs) balances:bn = bw + bs. (1.1)In many regions, including North America, the winter balance is a period of net accumu-lation, while the summer is a period of net loss (ablation).Land-terminating glaciers provide an relatively unambiguous measure of climatic vari-ability due to their strong response to atmospheric indicators. Glacier mass balance ismeasured using ablation stakes (Østrem and Brugman, 1991; Braithwaite, 2002), in con-junction with snow depth and density measurements. Although the most direct way tomeasure glacier mass balance involves manual measurements, simple predictive models haveshown that changes can be well predicted by winter precipitation, summer temperatures(see Table 1.1), and glacier geometry (Moore and Demuth, 2001; Letre´guilly, 1988). Airtemperatures have consistently been shown to produce accurate estimates of summer abla-tion; however, they lack the explanatory power to describe the main drivers of melt. Whilethe energy pathways to generate melt are more complex (termed the energy balance), airtemperatures are still an effective proxy for melt primarily because they are significantlycorrelated with many of the meteorological variables used in energy balance modelling.Table 1.1: Some notable statistical predictors of mass balance in western North America usingprecipitation (P ), temperature (T ), synoptic indicators, and glacier co-variability. a: Letre´guilly(1988), b: Moore and Demuth (2001), c: Beedle et al. (2009), d: Shea and Marshall (2007), e:Rasmussen and Conway (2004)Predictand Predictors Location R2 Sourcebn Pw, Ts (Vancouver, BC ) Place Glacier, BC 0.82 abn min.Ts (Jasper, AB) Peyto Glacier, AB 0.86 abn Ts, bw (Agassiz, BC ) Place Glacier, BC 0.63 bbn Pmarch (Prince George, BC ) Castle Creek Glacier, BC 0.51 cbw NDJF 500 mb heights Peyto Glacier, AB 0.65 dbs JJA 500 mb heights Peyto Glacier, AB 0.59 dbn, bs bn, bs co-variability Gulkana, Wolverine, AK 0.50 ebn, bs, bw bn, bs, bw co-variability Peyto, South Cascade, WA 0.60 e31.3. Glacier Flow1.3 Glacier FlowIn glacial environments, flow is regulated by the driving force of gravity, controlled by theup-glacier mass, and resisted by friction due to drag at the base and sidewalls of the glacier(Nye, 1952). The driving force responsible for down-valley flow is made possible by glaciersgaining mass in the upper reaches of the glacier (accumulation area), while losing mass inthe lowest reaches (ablation area). This allows for a mass imbalance causing the glacier tore-adjust its geometry and flow downslope.For land-terminating glaciers, flow is achieved through a combination of deformationand sliding, resisted partially by lateral and basal drag. The driving and resisting forces arecontrolled by topography, subglacial hydrology and mass balance (Benn and Evans, 2010).Seasonal variations in glacier velocity are related to melt rates, where higher melt ratesincrease the amount of water at the glacier bed, in turn reducing basal drag (Bartholomewet al., 2010; Zwally et al., 2002). Sub-annual variations in glacier flow speed are com-monly observed for land-terminating glaciers; however, over annual timescales, flow viadeformation and sliding is roughly equal (Nye, 1952; Benn and Evans, 2010).In calving glacier systems flow is complicated by regions where the glacier is in contactwith a water body. Water at the bed of the glacier leads to a significant reduction in basaldrag, or the complete removal, if the terminus is floating. This reduction in friction resultsin a substantial increase in velocity in the floating terminus (Van der Veen, 2002; Rignotand Kanagaratnam, 2006; Meier and Post, 1987).Calving itself has been shown to increase velocity as well, due to the excess ice lossfrom the terminus enhancing the longitudinal mass imbalance. Similar to a large ablationyear, significant ice loss from the terminus increases the average glacier slope, in turnenhancing the driving force and accelerating glacier flow (Van der Veen, 1996; Kirkbrideand Warren, 1999; Benn et al., 2007a). Variations in velocity following calving events havebeen documented in Greenland (Joughin et al., 2004), as well as Alaska, South Americaand New Zealand (Pelto and Warren, 1991; Brown et al., 1982). The connection betweencalving and velocity suggests that calving has the potential to enhance velocity, in turnenhancing calving, leading to a runaway effect.41.4. Glacier Calving1.4 Glacier CalvingFor glaciers that terminate in lacustrine or marine waters, calving can be an importantadditional component of the mass balance (bn, Equation 1.1), as it leads to greater ice lossthan possible through surface ablation alone (Van der Veen, 2002). The calving rate, Uc,is calculated asUc = UT −dLdt(1.2)where UT is the terminus-averaged surface velocity, and dLdt is the change in glacier lengthwith time.The calving rate can be calculated from field measurements (or estimates) of the surfacevelocity near the terminus, and the change in terminus position (Benn et al., 2007a).Calving glaciers are observed in both marine and fresh water settings. Calving ratesare roughly an order of magnitude greater in marine (tidewater) settings (Venteris, 1999;Haresign, 2004), suggesting that while land-terminating glaciers respond most readily toclimatic changes, and marine systems are relatively insensitive, lacustrine glaciers occupya ‘middle ground’. Differences between marine and lacustrine glaciers can be attributedto some processes only observed in marine environments, such as tides, greater watercirculation, and water density differences (Rignot et al., 2010; Benn et al., 2007b), whileother factors, such as flow speed and longitudinal strain rates, vary only in magnitude(Van der Veen, 2002; Warren and Kirkbride, 2003).Multiple studies have found significant empirical correlations between terminus waterdepth and calving rates (Skvarca et al., 2002; Warren and Kirkbride, 2003; Haresign, 2004).The relationship between calving rate and water depth applies in a variety of settings,from grounded termini to buoyant or near-floating termini across the globe. However,these ‘calving laws’ do not provide a physically based explanation for calving dynamics,nor does it account for observed seasonal fluctuations (Van der Veen, 2002). Furthermore,regional variations in the calving rate (Figure 1.1), as well as greater calving rates inmarine environments, indicate that water depth does not provide a full explanation forcalving variability.One of the reasons water depth is important for glacial stability is the potential for deepwater to promote terminus buoyancy due to the relative density difference between ice andwater. In many cases, flotation of the terminus has coincided with large-scale retreat (Meierand Post, 1987; Warren and Sugden, 1993). Flotation can produce immediate calving due51.4. Glacier CalvingWater Depth (m)Calving Rate (km/yr)Alaska (marine)Svalbard (marine)Alps (lacustrine)West Greenland (marine)New Zealand (lacustrine)Patagonia (lacustrine)0123450 100 150200250 300Figure 1.1: Regional variations in the calving rate with water depth. Adapted from Haresign(2004).to torque, while under low velocities a floating terminus may develop. If calving is notimmediate, the terminus becomes more sensitive to small perturbations in water level andstrain, eventually leading to terminus disintegration and substantial retreat (Boyce et al.,2007).1.4.1 Calving MechanismsCalving events are driven by an interplay of forces acting on the terminus. Longitudinalstrain, terminal force imbalance, subaqueous melt, and buoyant torque have been identifiedas the primary mechanisms responsible for individual calving events (Benn et al., 2007b).Longitudinal strain near the terminus, the product of large velocity gradients, has beenidentified as a primary driver of calving events. If velocity gradients are high enough, strainon the glacier leads to crevasse formation. If a crevasse is located relatively close to theterminus, it becomes a ‘preferential line of weakness’. This weakness provides the fracturemargin for calving events (see Figure 1.2), producing large tabular icebergs in floatingterminus environments.Calving can also result from an outward force imbalance at the terminus. Above the61.4. Glacier CalvingFigure 1.2: Calving along a ‘preferential line of weakness’, Bridge Glacier, June 20, 2013.waterline, outward-acting cryostatic pressure is opposed only nominally by the atmosphere.Below the waterline, cryostatic pressure is greater than the resisting hydrostatic pressure,except at the base of the ice front (Paterson, 2000). This imbalance is greatest at thewaterline, where the resisting water pressure is zero (Figure 1.3), leading to the generationof small volume, high frequency calving events, where ice ‘topples over’ into the waterbody.In cases where melt of the ice face is greater below the waterline than above, under-cutting will occur. This uneven melt progressively increases the force imbalance at theterminus. Conversely, in environments where the water temperature is low or stronglystratified, an ‘ice foot’ can develop (Dykes et al., 2011; Rohl, 2006). The ice foot is lessdense than the surrounding water, and buoyant forces cause it to fail, shooting to thesurface in dramatic fashion (Skvarca et al., 2002; Motyka et al., 2003).Although the calving mechanisms have been presented as separate entities, it is rarefor a situation to exist without multiple mechanisms acting on the terminus. However, the71.4. Glacier CalvingPice  PwaterPwaterPiceFigure 1.3: The force imbalance at the calving front, adapted from Benn et al. (2007b).dominant mechanism differs depending on the current glacier, water body and climatolog-ical conditions. The mechanisms acting on the terminus have significant implications forthe magnitude and character of calving events.In floating terminus environments, a hierarchy in driving calving mechanisms exists.Longitudinal strain is the only mechanism that creates a physical margin along whichfailure can propagate, making it a first order control on the calving rate. The magnitude ofevents triggered due to strain is often orders of magnitude greater than other mechanisms.These large, relatively infrequent events supersede the low volume, high frequency eventsof terminus force imbalances, even often occurring on overlapping areas.Conversely, for grounded termini, longitudinal strain is less important in determiningthe calving flux, since glacier flow patterns do not allow crevasses to penetrate the fulldepth of the terminus (Benn et al., 2007a). Instead, the terminus force imbalance becomesthe main driver of individual calving events. This change in driving mechanisms affectsthe magnitude and character of calving. While floating termini in some environmentsexperience low frequency, large tabular calving events, calving driven by a force imbalanceis defined by low magnitude, high frequency events.Water temperature also plays a significant role in determining the character and magni-tude of calving events due to its influence on terminus geometry. High water temperaturespromote undercutting, which in turn enhances the terminus force imbalance, promoting81.5. An Intricately Coupled Systemfurther high frequency, low magnitude events. Alternatively, low water temperatures in-hibit subacqueous melt, and if above-water melt is sufficient, can allow for the developmentof an ice-foot. Ice feet fail due to buoyant torque, and supply low frequency, moderate tolarge magnitude events in floating terminus environments.Although water temperature has an effect on terminus geometry, its influence is oftendampened by terminus ice velocity. Low velocities allow for longer exposure of the ice to thewater body, allowing greater melt rates, and a quicker establishment of instability (Rohl,2006). High glacier terminus velocities dampen this effect by limiting the amount of timethe ice face is exposed to warmer waters before being transported to deeper water, wherethe glacier calves more readily. The relative importance of individual mechanisms havebeen related to terminus velocities (Figure 1.4), where fast glaciers respond more stronglyto topographic controls, while slower glaciers are more affected by lacustrine processes.Glacial )low 2versteepening7hermal 8ndercutting/ongitudinal Stretching%uoyancyRelative ,mportanceGlacier )low SpeedFigure 1.4: The relative importance of calving controls related to glacier flow speed at the terminus.Adapted from Haresign (2004).1.5 An Intricately Coupled SystemAttempts to formulate a portable ‘calving law’ to predict calving rates have proven elusivedue to significant variability in calving environments (see Table 1.2). Lacustrine Patagonian91.5. An Intricately Coupled SystemGlaciers Upsala and Leon can be characterized by high flow speeds, high lake temperatures,and deep waters, promoting high calving rates. Conversely, small lacustrine glaciers Maud,Grey, Godley, Hooker, and Ruth in the Southern Alps (New Zealand) have small waterdepths, low flow speeds, and low lake temperatures, all of which restrict calving losses.While some studies have found topography to be the major control on calving rates andterminus position (Koppes et al., 2011), others have attributed changes to flow gradients(Tsutaki et al., 2011), water temperature (Dykes et al., 2011; Rohl, 2006; Robertson et al.,2012; Motyka et al., 2003), and buoyancy (Boyce et al., 2007).Table 1.2: Characteristics of selected major lacustrine calving glacier studies. Dw is the meanwater depth, Tw is the mean water (depth averaged or range) temperature, Uc is the calving rate,and UT is the terminus averaged flow speed. Citations: a: Boyce et al. (2007), b: Motyka et al.(2003), c: Warren and Kirkbride (2003), d: Dykes et al. (2011), e: Warren and Aniya (1999), f :Stuefer et al. (2007), g: Haresign (2004)Location Year Dw (m) Tw (◦C) Uc (ma−1) UT (ma−1) SourceAlaskaMendenhall 1997 - 2004 45 - 52 1 - 3 12 - 431 45 - 55 a, bNew ZealandMaud, Grey,Godley, Ruth,Hooker1995 4 - 20 1.7 - 4.3 14 - 88 5 - 151 cTasman 1995 10 0.5 28 11 c2000 - 2006 50 1 - 10 78 69 d2006 - 2008 153 1 - 10 227 218 dPatagoniaUpsala West,Ameghino,Grey, Nef1995 153 - 325 2 - 7 355 - 2020 370 - 1620 ePerito Mereno 1995 - 2006 175 5.5 - 7.6 510 535 e, fLeon 2001 65 4.5 - 7.0 520 - 1770 520 - 1810 gIcelandFjallsjokull 2003 75 1.5 - 3.0 582 258 gAlong with regional heterogeneity, there is substantial variability in the temporal con-trols on calving rates. Over long timescales (decadal or longer), climatic warming has beensuggested as the dominant control on glacial mass balance and calving rates (Warren andKirkbride, 2003; Motyka et al., 2003). Sustained negative mass balance for multiple yearsresults in thinning, most pronounced near the terminus, which can initiate buoyancy-drivenretreat from deep to shallower water. Over annual timescales there is a seasonal increase in101.6. Research Gaps and Study Goalsthe calving rate in the late summer, attributed to surface ablation and thinning (Haresign,2004; Van der Veen, 2002). Thinning acts in concert with an observed increase in lake level(and most likely subglacial meltwater enhanced velocity), increasing the probability of theterminus achieving flotation.For lacustrine glaciers, temporal variability is further accentuated by the transient ‘life-cycle’ of calving retreat over an over-deepened lake or freshwater fjord. Controls, mech-anisms, and calving rates vary as the glacier retreats from deep waters, where geometryand submarine melt play dominant roles, and calving rates are high (Meier and Post, 1987;Rignot et al., 2010), to shallower waters, where flow gradients and water level play a largerrole and calving losses are significantly smaller (Tsutaki et al., 2011; Boyce et al., 2007).Eventually, the terminus retreats to dry land, where calving is no longer a significant sourceof ice loss, allowing the glacier to thicken and possibly re-advance (Van der Veen, 2002), orcontinue to retreat due to a legacy of dynamic thinning and continued climatic warming(Oerlemans, 2010).1.6 Research Gaps and Study GoalsWhile advances have been made in predicting the magnitude and timing of calving inwater-terminating glaciers, there remains many questions about their nature, timing, andimportance. There is an insufficient theoretical understanding of the drivers of calvingto make reliable predictions and quantify ’calving laws’. Similarly, a robust empiricalframework has so far proven elusive due to a complex set of interconnected mechanisms,diverse range of environments, and relative paucity of study sites worldwide. In orderto better understand each unique calving system, a better understanding of the broadcommonalities in lake-terminating glacier response is required.Although understanding lacustrine glaciers is of critical importance for better watershedmanagement and for isolating the climatic signal in calving glaciers, few lacustrine glaciersare studied worldwide. Most notably, work in North America has focused on MendenhallGlacier in Alaska (Motyka et al., 2003; Boyce et al., 2007), while international efforts havefocused on Tasman Glacier in New Zealand (Warren and Kirkbride, 2003; Dykes et al.,2011; Dykes and Brook, 2010) and a few Patagonian glaciers, most notably Perito MerenoGlacier (Warren and Sugden, 1993; Warren and Aniya, 1999; Stuefer et al., 2007). BridgeGlacier presents another valuable study site to supplement the worldwide database, andprovides information from a mountain range where no other calving glacier is currently111.7. Study Locationbeing studied.Few studies have focused on relating the calving flux to mass balance in order to quantifythe relative importance of calving on the mass balance of a lake-terminating glacier. Abetter understanding of the glaciological, lacustrine, and climatological conditions, and howthey relate to calving, will help outline driving conditions promoting ice loss. Furthermore,this data will help elucidate the broad commonalities between calving glaciers worldwide,allowing for a more universal understanding of lacustrine glacier calving.This thesis explores the relative importance of calving and surface ablation on thesummer mass loss at Bridge Glacier, BC. Specifically, it investigates the following questions:• What is the relationship between calving and climate in a lacustrine calving glacier?How much of the observed retreat at Bridge Glacier can be attributed to calving?• How much ice loss at Bridge Glacier is due to surface melt? How sensitive is BridgeGlacier to changing climatic conditions? What are the main meteorological variablesdriving melt?• Can surface terminus velocities, integral to an accurate estimate of the calving flux,be adequately measured using oblique time lapse imagery? How accurate is thismethod?• How does the calving flux compare to the volume of ice lost through surface melt?What are the main variables controlling the calving rate over the 2013 melt season?What is the relationship between glacier dynamics, surface melt, and calving events?• How representative is 2013 relative to the calving history of Bridge Glacier? Howhas the calving flux varied over the last 30 years. How do changes in the historicalcalving flux relate to changing environmental conditions?• What does the future hold? Will the relative importance of calving and surface meltchange over the coming decade(s)?1.7 Study LocationBridge Glacier (50◦48’11”N, 123◦38’40”W) is located in the Pacific Ranges of the CoastMountains of southwestern British Columbia, Canada (Figure 1.5). The glacier is an outlet121.7. Study Locationof the Lillooet Icefield, one of the largest icefields in British Columbia, and is locatedroughly 175 km north of Vancouver, Canada. Bridge Glacier has an area of 81 km2,extending from an elevation of over 2900 m at Bridge Peak to 1390 m, where it terminatesin a large proglacial lake. The lake has grown from under 2 km2 in 1972 to over 6 km2in 2013. The glacier has experienced large tabular calving events since the early 1990s,indicative of a floating terminus. The far (east) end of the lake houses numerous large(several hundred m2) icebergs, which have been present for multiple (in some cases morethan 10) years.VancouverPembertonBridge Glacier20 kmComoxUSACOAST MOUNTAINSLillooetNFigure 1.5: Study site location taken from Google Map satellite imagery of southwestern BritishColumbia.The region is characterized by a large precipitation gradient, and Bridge Glacier lieson the divide between the heavy precipitation coastal Pacific Ranges, and the arid interiorChilcotin Ranges. Synoptic air flow is predominantly from the west, leading to the highestelevation, most westerly reaches experiencing heavy snowfall, while the eastern flank of theglacier is more arid. Manual snow surveys at the eastern edge of Bridge Lake show mean131.8. Thesis StructureMay 1 SWE of 600 mm (, while theregion is characterized by cool moist winters, and warm dry summers.Bridge Glacier is also a particularly valuable site to study due to its importance forindustry, and as a site of previous glaciological work. Bridge River, which begins as anoutlet of the proglacial lake in which Bridge Glacier terminates, supplies water for a hy-droelectric complex consisting of 3 dams. The project has been operational since 1960,and currently supplies 6 - 8% of British Columbia’s electricity ( Bridge Glacier was also a siteof Environment Canada annual mass balance surveys from 1977 until the program wasdiscontinued in 1985 (Mokievsky-Zubok et al., 1985). More recently, the glacier has beenthe location of several glacio-meterological, hydrological and paleo-environmental studies(Ryder, 1991; Stahl et al., 2008; Shea et al., 2009; Allen and Smith, 2007).1.8 Thesis StructureThis thesis aims to shed light on the research goals (Section 1.6) by examining the pertinentconcepts thematically. Chapter 2 presents an introduction to Bridge Glacier and outlinesits recent history from 1972 to 2012. The chapter reconstructs Bridge Glacier’s retreatrate, and uses a linear inverse model to estimate the climatic-driven retreat, giving a firstorder estimate of the relative importance of calving.Chapter 3 estimates surface melt during the 2013 melt season. Surface melt is modelledusing a distributed energy balance model that spatially distributes meteorological variables.The model follows, and evaluates, recent work on modelling katabatic flow to account forthe relatively pronounced effect of katabatic winds on temperature, humidity and windspeed across the glacier. Finally, this chapter outlines the main meteorological variablesdriving surface melt, and estimates how sensitive Bridge Glacier is to changing climaticconditions.Chapter 4 presents a method for producing glacier velocity measurements from a timelapse camera. The method requires only a single camera, open source software, and basicgeometry to manually measure features on-glacier. The chapter then evaluates the potentialand observed errors in the method and compares results to in situ measurements. Finallythis chapter presents daily and seasonal displacements from 2 sites near the terminus ofBridge Glacier.Chapter 5 quantifies the calving flux for the 2013 season. The volume of ice discharged141.8. Thesis Structurefrom the terminus is reconstructed from satellite imagery, time lapse camera derived ter-minus flow speeds, in situ ice front height observations, and lake bathymetry. Individualcalving events, captured using time lapse cameras, are also related to meteorological andlimnological variables, and glacier velocities. The character and driving forces of calvingevents are examined and related to the current terminus environment. Finally, the relativeimportance of calving during the 2013 melt season is explored.Chapter 6 synthesizes the main findings and implications from the 2013 melt season.This chapter explores how these findings relate to historical rates of calving and surfacemelt, what these findings mean for the future of Bridge Glacier, and how the relativeimportance of calving is be expected to change over the coming decades. Finally, thischapter examines how findings from Bridge Glacier relate to other lacustrine calving glaciersin order to identify broad commonalities and unique qualities in each system.Chapter 7 summarizes the main findings from this study, and provides potential avenuesof future work.15Chapter 2A First Order Estimate of theInfluence of Climate on Retreat2.1 IntroductionGlaciers serve as valuable barometers of a changing climate due to the fact that variationsin their length, area and mass can be tied directly to regional and global air temperatures.While this strong correlation provides valuable information for remote locations whereweather station records do not exist, the climatic signal in glacier systems is obscured whencalving is present. This is due to the fact that calving allows for an additional mechanismof ice loss, leading to larger volumes of ice loss than would be otherwise possible throughsurface melt alone. In order to isolate the climatic signal in calving glacier systems, andaccurately predict future retreat, a better understanding of the relative importance andinfluence of calving is needed.In other lacustrine calving systems, dramatic retreat has been found to coincide withterminus floatation (Warren and Kirkbride, 2003; Boyce et al., 2007; Dykes et al., 2011).In many cases, terminus flotation may be achieved through thinning near the terminus dueto successive years of high melt rates. Enhanced thinning due to a run of hot summershas been observed at Mendenhall Glacier (Motyka et al., 2003). This climate inducedthinning led to increased instability and propensity to calve, and eventually to the collapseand widespread retreat into shallower waters (Boyce et al., 2007). Similar findings havebeen made at Tasman Glacier, New Zealand (Warren and Kirkbride, 2003; Dykes andBrook, 2010), and in Patagonia (Skvarca et al., 2002), suggesting that retreat from climaticwarming may enhance calving rates over decadal time scales.While a few studies have examined lake-terminating glacier retreat (Haresign, 2004),the data are limited both by the number of glaciers observed and temporal extent, wheremost studies have focused on sub-decadal responses (Benn et al., 2007b). In order to162.2. Methodsbetter understand the long-term impact of climate on glacier retreat, the full range ofcalving behaviour needs to be taken into account, from calving in lake-contact termini, togrounded termini, and to buoyant termini. The goal of this chapter is to broadly outline therelationship between calving and climate in a lacustrine terminating calving glacier. First,this research aims to reconstruct and quantify retreat at Bridge Glacier, British Columbia,between 1972 and 2012. Second, this chapter aims to quantify the role of calving on retreatby comparing observed rates to the expected climatic response. Finally this chapter willmake a first order estimate on the relationship between climate and calving, and how thatvaries through the ‘life-cycle’ of a lacustrine calving glacier.2.2 MethodsTerminus positions were obtained from Landsat imagery available using LandsatLookViewer ( through MSS, TM, ETM+ and ETM+ SLCOFF sensors. Bridge Glacier imagery was obtained between 1972 and 2012; however, be-fore 1984 usable imagery was sparse, and terminus positions are missing from 1978 to 1981,as well as 1983. The terminus position for each year was recorded using imagery taken atthe end of the melt season (August 26 to October 16), using the last available image beforethe glacier became snow-covered. Notable exceptions were 1972, 1973, 1975 and 1977,when imagery was only available for mid-August dates, and 1982, when a December imagewas the only image that season that was not cloud-covered. Three points were taken foreach season’s terminus and were connected to delineate the shape of the terminus. Retreatwas calculated by finding the average down-glacier change between the two annual terminipositions.Climate data were obtained from Environment Canada for Whistler, BC (50◦08’ N122◦57’ W, elevation = 659 m, ID #1048898), and Vancouver International Airport, BC(49◦12’ N 123◦11’ W, elevation = 4 m, ID #1108447), the two closest stations with rela-tively complete records for the period of study ( had monthly records of precipitation and temperature since 1977; however, datasince 2007 had not been verified or calibrated, and contained gaps and physically improb-able values in the record. The Vancouver International Airport data spaned the periodfrom 1937 to present, but contained incomplete entries for September 2010 as well as forJanuary to May 2012. Annual winter precipitation was computed by summing precipita-tion between November and April (inclusive). Mean summer temperatures were obtained172.3. Resultsby averaging mean monthly temperatures between May and October (inclusive).2.3 Results2.3.1 Retreat and Climate(a) September 22, 1985 (b) September 11, 1995.(c) October 23, 2005 (d) October 10, 2012.Figure 2.1: Landsat scenes from late fall 1985, 1995, 2005, and 2012. Note the substantial thinningand disintegration of the terminus once it reaches flotation.Landsat imagery (see, Figure 2.1) showed a relatively stable terminus position between 1972 and1984. During this period the proglacial lake was relatively small and did not contain largeicebergs. In 1989, icebergs were consistently observed in the lake, suggesting that calvingwas occurring. Larger icebergs, however, only began to appear in the lake in the summerof 1991. In 1994, and again in 1997, two large calving events led to a more substantialretreat. The glacier entered a period of relative stability from 1997 until 2004, at whichpoint large crevasses formed near the terminus and became water-filled. In the summer of2005, the terminus disintegrated and retreated 469 m, the largest calving rate during the182.3. Resultsstudy period. After 2005, the glacier remained at a relatively stable position until 2010,when it underwent another substantial retreat toward its present day margin.2.3.2 ClimateFor Whistler and Vancouver, annual winter precipitation is significantly correlated withthe May 1 snow water equivalent for Bridge Glacier terminus, obtained between 1995 and2012 by the Water Stewardship Division of the Government of British Columbia Ministry ofEnvironment (Station #1C39, Therelation with Whistler (n = 13) has a larger r2 value than that for Vancouver (Figure 2.2);however, caution is required since the sample size is larger for the Vancouver dataset(n = 18).500 700 900 11004006008001000Whistler Winter Precipitation (mm)Bridge May 1 SWE (mm)1995-2007R2 = 0.84500 700 900 11004006008001000Vancouver Winter Precipitation (mm)Bridge May 1 SWE (mm)1995-2011R2 = 0.43Figure 2.2: Relationship between winter precipitation at Whistler and Vancouver, BC and BridgeGlacier May 1 SWE.Values for both Vancouver and Whistler are significant predictors of mean annual flow(1979-2010) recorded at the outlet of Bridge Glacier (Figure 2.3), indicative of summer meltat Bridge Glacier. Hydrometric data was obtained from Environment Canada gauging sta-tion 08ME023 (50◦51’22” N, 123◦27’01” W, ). Vancouver is a better predictor of Bridge River mean annual192.3. Resultsflow than Whistler, suggesting it is a strong predictor of summer melt. Due to the factthat Vancouver has a longer climate record than Whistler, and a strong correlation withBridge Glacier climatology, it has been used for all further analyses in this study.11 12 13 14 15 165101520Whistler Summer Mean Temperature ( ° C)Mean Annual Flow (cms)R2 = 0.591977-200711 12 13 14 15 165101520Vancouver Summer Mean Temperature ( ° C)Mean Annual Flow (cms)R2 = 0.651977-2007Figure 2.3: Mean Summer Temperature (May - October) for Whistler and Vancouver related toMean Annual Flow at Bridge Glacier’s outlet.Winter precipitation (Figure 2.4a) during the period oscillated in short, 3 to 5 year cy-cles of high and low snowfalls, with relatively large variations. Average winter precipitationfor the period was 819 mm, but varied up to 400 mm in individual years. Four consecutiveheavy snowfall years were observed between 1980 and 1984, while from 2000 to 2009 onlytwo seasons experienced above average precipitation. In general, the relative lack of a longterm trend or shift in precipitation patterns suggest that its impact on retreat would havebeen minimal, a conclusion echoed in other studies (Oerlemans, 2005; Post et al., 2011).Summer temperatures (Figure 2.4b) show a general trend of warming during the period.A long run of warmer than average summers was observed from 1987 to 1998, with onlytwo years falling below normals. From 1999 to 2001, three consecutive summers were coolerthan normal, followed by 3 of the 5 warmest summers in the record from 2003 to 2005.Both winter precipitation and summer temperature patterns are significant predictors ofBridge Glacier equilibrium line altitudes.The equilibrium line altitude (ELA, Figure 2.4c) showed two strongly below averageseasons in 1974 and 1976. The ELA was anomalously high for 9 out of 11 seasons, between202.3. ResultsPrecipitation Anomaly (mm)-2000200 a.Temperature Anomaly ( ° C)- Anomaly (m)-1000100 c.MAF Anomaly (m3 /s)-10-50510d.Retreat Rate (m/yr)-2000200600e.1980 1990 2000 2010Figure 2.4: a. Vancouver winter precipitation anomaly (x¯ = 819 mm) b. Vancouver summertemperature anomaly (x¯ = 14.8◦C), c. Bridge Glacier equilibrium line altitude (x¯ = 2089 m), d.Bridge River mean annual flow anomaly (x¯ = 10.7 m3s−1), e. Annual retreat rate of Bridge Glacier,dashed line is loess-smoothed retreat (span = 0.5).1985 and 1999. Save one anomalously low ELA season in 2000, only minor deviations fromthe mean were observed in the last 20 years. The ELA trend is well correlated with bothtemperature and mean annual flow anomalies.212.3. ResultsMean annual flow from the gauging station at the outlet of Bridge Lake (Figure 2.4d.)follows the general pattern of summer temperatures and is significantly correlated with bothWhistler and Vancouver, however Vancouver is a stronger predictor. Four consecutive lowflow years were observed between 1983 and 1985, while between 1989 and 1998, 9 of the10 years were above average. The period from 2003 to 2006 was characterized by anotherrun of 4 consecutive high flow years.Bridge Glacier’s retreat (Figure 2.4e.) was relatively consistent until 1991 when calvingbegan to have a larger impact on terminus position. Bridge Glacier’s post-1991 retreat ischaracterized by large retreat years followed by periods of relative stability. The rate ofretreat between 1972 and 1991 was relatively low (21 ma−1), but accelerated substantiallyafter 1991, with an average rate of 144 ma−1. The rate of retreat accelerated again since2009, retreating over 1200 m in 3 years.2.3.3 Inverse Linear ModelIn order to estimate the relative roles of calving and climate, Bridge Glacier’s retreat wasmodelled using an inverse linear method (Oerlemans, 2005, 2007). The model treats theglacier as a one-dimensional (longitudinal) flowline that responds only to changes in climaticconditions. By modelling retreat solely as a function of climatic forcing, we can projecthow the glacier would have responded had it not been affected by calving. Through the useof this model, the influence of climate can be quantified, while deviations from modelledresults can be hypothesized to be due to calving influence.Table 2.1: Parameters used in the linear inverse model. Both c and τ fit well within the rangefound in Oerlemans (2005).Model Parameters Symbol ValueSlope θ 0.0818 m/mGlacier Reference Length L 18 500 mMean Annual Winter Precipitation P 819.31 mm/yearMass Balance Gradient β 0.00543 K/mClimatic Sensitivity c 2366.36 m/KResponse Time τ 139.03 yearsIn the linear inverse model, the retreat rate (dLdt ) is calculated asdLdt= −1τ[cT ′(t) + L′(t)] (2.1a)222.3. Resultswhere T ′ is the annual summer temperature anomaly (◦C), c is the glacier’s individualclimate sensitivity (m/K), and τ is response time in years, and L′ is the glacier’s length,relative to a reference length (see Table 2.1). The relationship can be rewritten to solvefor L′(t) asL′(t) = (1− e− tτ )cT ′(t). (2.1b)which implies that when the time (t) is equal to the response time (τ), roughly two-thirdsof the climatic forcing has manifested itself as a change in glacier length.In order to reconstruct retreat rates, a response time and climatic sensitivity are cal-culated asτ =13.6βθ(1 + 20θ)−0.5L−0.5 (2.2a)c =2.3P 0.6θ(2.2b)where, following Oerlemans (2005), θ is the average glacier slope and L is the glacierreference length in metres. The mass balance gradient (K/m), β, is calculated asβ = 0.006P 0.5 (2.2c)and is a function of P , the average annual winter precipitation (mm/year).Values of response time and climate sensitivity were derived using 1972 as the referenceyear (see Table 2.1). The derived response time for Bridge Glacier is 139 years, whichis relatively long, but consistent with other low slope, large glaciers (Oerlemans, 2005).Climate sensitivity was increased by 50% because the accumulation area is roughly fourtimes the ablation area, a method that is consistent with other studies, and has beenconfirmed in numerical simulations (Oerlemans, 2007, 2011). Therefore, the value of c iswell within the range of average values derived for glaciers around the world (Oerlemans,2005).The model is generally robust to changes in parameters (Figure 2.5). Sensitivity anal-yses suggest that changes in climate sensitivity and response time have a minor impact onthe amplitude of response, and cannot account for the far larger observed retreat rates.Since modelled results from Equation 2.1b are noisy, and retreat rates fluctuate byas much as 571 m, length changes are smoothed (Oerlemans, 2007) using local regressiontechniques (loess in the R statistical package). The loess-smoothed modelled annualretreat predicts a cumulative retreat of 613 m between 1972 and 2012. Modelled retreat232.4. Discussion1980 1990 2000 20100100020003000YearCumulative Retreat since 1972 (m)Observed RetreatModeled RetreatLoess Smoothed Modelled Retreatτ = 75 years, c = 2366 m/Kτ = 200 years, c = 2366 m/Kτ = 139 years, c = 3000 m/Kτ = 139 years, c = 1500 m/KFigure 2.5: Modelled and observed retreat at Bridge Glacier. The model is run with variedresponse time (τ) and climatic sensitivity (c) parameters.agrees quite well with observed rates until 1991, at which point observed retreat ratesincrease substantially while modelled retreat increases only gradually. This departurebetween observed and modelled retreat coincides with the onset of large calving eventsin Bridge Lake. In the summer of 2005, large calving events further increases the rate ofterminus retreat, while the modelled rates do not vary significantly, leading to a furtherdeparture between the two retreat rates.2.4 DiscussionThe general agreement between the observed and modelled results up to 1991 suggests that,between 1972 and 1991, Bridge Glacier’s retreat was directly related to climatic conditions.242.4. DiscussionBetween 1972 and 1991, anomalously hot summers produced retreat consistent with thederived response times and climate sensitivity. Since 1991, however, observed retreat was atleast partially independent of climatic forcing. The acceleration of retreat, asynchronousto climate trends, suggests that Bridge Glacier’s retreat was less directly influenced byclimatic conditions once calving began. After the onset of calving, the glacier exhibited astaggered retreat, characterized by larger annual variations than could be expected fromclimate-driven retreat alone, followed by periods without any significant retreat.Although calving modifies the relationship between climate and glacier retreat, climatemay still play an important role in the long term retreat of calving glaciers. Many studieshave found that flotation is a key variable in determining the stability and retreat of water-terminating glaciers (Meier and Post, 1987; Van der Veen, 2002); however, the mechanismsof thinning are diverse.For Bridge Glacier, calving only had a measurable effect on the terminus position onceit retreated into deeper water and achieved floatation, a finding is supported by the large,infrequent tabular icebergs produced from the terminus. Similar observations have beenmade at Mendenhall Glacier in Alaska (Motyka et al., 2003; Boyce et al., 2007), whereseries of high melt summers preceded flotation and large scale retreat. Climatic indicatorsfor Bridge Glacier suggest strongly negative mass balance years from 1986 to 1991 (Figure2.4), which would have increased the rate of glacier thinning and created the conditionsnecessary for terminus flotation, and enhanced the rate of calving.While it is possible that other factors could have affected the terminus thinning duringthe period leading up to the floatation of Bridge Glacier’s terminus, these scenarios aremuch less probable. It is possible that increased basal water due to enhanced melt rates(or an earlier, more prolonged melt season) would have allowed for greater flow speeds,enhancing strain rates and dynamic thinning. It is also possible that thinning could havebeen produced by a decrease in flow speeds in the lowest reaches of the glacier, leading to alower supply of ice to the terminus, and ice stagnation. Both of these conditions, however,require a systematic shift in the dynamics of the glacier system before the onset of calving.While they cannot be ruled out completely, it is unlikely that these these scenarios can fullyexplain the significant increase in Bridge Glacier’s retreat rate since they do not correspondwith the observed climatic trend of warming preceding the increased rate of retreat.Once calving begins to have a significant impact on the retreat rate of the glacier italso plays a role in driving for glacier dynamics. Once Bridge Glacier reaches flotation,calving becomes a more important source of ice loss. Increased ice loss, leading to terminus252.5. Conclusionretreat, increases near-terminus velocities and longitudinal strain rates. This process, inturn, leads to dynamic thinning, crevassing, increased instability, and more calving events,instigating a positive feedback loop independent of climate (Figure 2.6). While climateis responsible for the initial thinning, once the terminus achieves flotation, retreat can beenhanced without any further climatic forcing. In this mechanism, climate is responsiblefor creating the conditions necessary for a calving-driven retreat, but once floatation-drivencalving begins, the glacier becomes more vulnerable to calving events, making it inherentlyunstable and unable to regain stability until it becomes grounded.Figure 2.6: Schematic representation of the processes that drive the calving-climate feedback loop.2.5 ConclusionIn this chapter, the retreat of Bridge Glacier has been reconstructed and has been related toclimate through a simulation of its linear response to climatic forcing. From 1972 to 2012,Bridge Glacier has retreated 3.55 km. The retreat rate was modest up until 1991, with lessthan 600 m of cumulative retreat. Since 1991, however, retreat accelerated, culminating inretreat rates of over 400 ma−1 between 2009 and 2012. The post-1991 retreat is irregular,262.5. Conclusioncharacterized by large tabular calving events and high retreat years, and periods of relativeterminus stability. The high magnitude, low frequency tabular character of observed calvingevents suggest that calving is primarily due to terminus floatation (Benn et al., 2007b).Observed retreat follows modelled retreat projected from climatic conditions until 1991;at which point retreat rates far exceed modelled results. Bridge Glacier was retreating insync with climate until 1991, however, since then the retreat rate is driven primarily bycalving. Modelled results suggest that calving is responsible for an additional 2.5 - 3.0km of retreat, suggesting that since 1991, calving has been the dominant control on theterminus position.Although these findings emphasize the relative importance of calving in Bridge Glacier’srecent retreat, climatic records and modelling suggest that an initial climatic warming isnecessary to create the conditions necessary for calving to become an important componentof retreat. A series of consecutive high ablation summers is necessary to thin the lowestreaches of the glacier enough to allow the the terminus to become buoyant. Once theterminus becomes buoyant, however, retreat can continue independent of the larger climatictrend due to the vulnerability of the terminus to calving. Terminus flotation and increasedcalving rates in turn promote higher terminus velocity, accelerating terminus thinning andinstigating a positive feedback loop where a calving-driven retreat can continue independentof climatic conditions.This first order estimate of calving at Bridge Glacier suggests that while the inverselinear model employed in this chapter attributes the vast majority of observed retreat inthe last 40 years to calving, calving itself acts within the initial bounds of climate. Withinthe ‘life-cycle’ of a lacustrine calving glacier, climatic forcing is necessary to create theconditions necessary to promote further calving retreat.This work suggests that calving is responsible for the majority of the observed retreatduring this period; however, the same implication does not necessarily hold true for BridgeGlacier’s mass balance during the period. Beyond the implication that terminus buoyancyis critical to determining the calving rate, this work does not elucidate any of the controlsof calving itself on Bridge Glacier. In order to better constrain the role and importance ofcalving to the overall ‘health’ of Bridge Glacier, a quantitative account of the volumetric iceloss from both mechanisms is required. These questions, among others, are the overarchingconcepts that drive research in subsequent chapters.27Chapter 3A Distributed Energy BalanceModel for Ice Ablation3.1 IntroductionSince the Little Ice Age, the cryosphere, and mountain glaciers in particular, have beenreceding at an accelerated pace (Barry, 2006). Although irregular, this retreat is pervasiveacross regions within Canada and across the globe, and is well correlated to rising meanglobal temperatures (Oerlemans, 2005). This retreat has significant ramifications for futurevolume and timing of water resources, as well as future eustatic sea level rise (Dyurgerovand Meier, 2005), and requires a detailed predictive system to ensure proper managementover short and long timescales.Numerical melt models allow for a quantification of meltwater generation from glacierand snow surfaces by relating surface melt to meteorological or climatological variables.Melt models vary substantially in complexity and data requirements and have been shownto have high predictive value over a range of temporal and spatial scales (Hock, 2005).Empirical frameworks, such as Positive Degree Day models (Braithwaite and Olesen, 1989)and Temperature Index models (Hock, 2003), which relate melt to air temperature, canprovide very good daily to seasonal predictions. Empirical temperature models are desir-able in many cases due to few data requirements and their ability to out-perform physicallybased models in locations where detailed reliable meteorological data cannot be obtained(Shea, 2010).3.1.1 Surface Energy Balance ModellingSurface energy balance models offer an alternative to empirical models by explicitly ac-counting for the fluxes of energy available for melt at the glacier surface. While energybalance models have higher predictive power over short timescales, they also require sub-283.1. Introductionstantially more input variables, increasing the complexity of the model (Hock, 2005). Fur-thermore, the transferability of energy balance model input variables is limited, necessi-tating detailed in situ meteorological observations to accurately capture melt processes(Anslow et al., 2008; MacDougall and Flowers, 2011). Notwithstanding, energy balancemodels are preferable in many cases because they also offer the added benefit of quantifyingthe relative importance of individual melt processes, allowing for future projections of thesensitivity of glacier melt to a changing climate.Surface energy balance models are most commonly made at a single point on the glacier(Munro, 2006; Brock et al., 2000). Typically, a weather station is installed at a locationwithin the ablation area that measures all relevant energy balance variables. While thismethod is logistically superior, it does not account for spatial variability in relevant me-teorological variables, and relies on the representativeness of a single point to accuratelycapture the main drivers of melt over the glacier.As an extension of point energy balance models, distributed energy balance modelspresent a more complex model that can account for the spatial variation in melt acrossthe glacier. The energy available for melt is calculated for each grid point on the glacierby interpolating or distributing meteorological parameters spatially. With this increase inspatial complexity, distributed models have larger data requirements, including a digitalelevation model, and a routine to spatially distribute turbulent flux components.Turbulent Heat Parameter DistributionsSensible and latent heat fluxes are important components of melt energy supply, and arecontrolled by temperature and vapour pressure gradients, as well as wind speed. Theseturbulent fluxes distribute predominantly as a function of elevation; however, in glacier-ized basins this spatial and altitudinal distribution is disrupted by the development andpresence of katabatic winds (Shea and Moore, 2010). Air directly above the glacier sur-face is cooled, creating denser, less buoyant air, and initiating down-valley flow. As airflows downslope it is further cooled, creating a persistent down-glacier wind (Munro, 2006).Katabatic flows create a sharp, shallow thermocline (between 5 and 20 m above the sur-face) where temperatures have been observed to be significantly cooler than the ambientair temperature (Oerlemans, 2010; Shea and Moore, 2010).During the melt season, winds over the glacier surface are predominantly characterizedby a down-glacier katabatic flow (Oerlemans, 2010; Oerlemans and Reichert, 2000), while293.2. Instrumentation and Field Methodsobservations in the Coast Mountains found that the presence of katabatic flow is dependenton a threshold temperature (Shea and Moore, 2010). Using piecewise regression, katabaticeffects were found only when the ambient air temperature was greater than 4◦C. Thisobservation was corroborated on Pasterze Glacier, Austria, where higher downslope windspeeds correlate well with higher ambient air temperatures (Greuell and Smeets, 2001).While many studies have used a constant wind speed measured at on-glacier weatherstations to distribute wind speed across the glacier (Hock and Holmgren, 2005; MacDougalland Flowers, 2011), this distribution scheme does not take into account the inherent com-plexity of katabatic flow. Expected spatial variability in katabatic wind speed magnitudes,as well as terrain sheltering and exposure across the glacier, suggest that a constant windspeed across the glacier is unrealistic, and will ultimately reduce the accuracy of turbulentheat calculations. Similarly, temperature is commonly distributed using a −6 ◦C km−1lapse rate (Hock, 2005), or an average measured seasonal lapse rate (Anslow et al., 2008;MacDougall and Flowers, 2011). However, the cooling of on-glacier air temperatures bykatabatic flow ultimately suggests that standard atmospheric lapse rates are not suitableto distribute on-glacier air temperatures (Munro, 2006; Shea and Moore, 2010).As an alternative, recent work has focused on modelling the magnitude and extent ofkatabatic winds across a range of glacier sizes (Shea, 2010). Using terrain morphometry,katabatic wind is shown to be stronger in the lowest reaches of the glacier, and can signif-icantly depress temperature and vapour pressure (Shea and Moore, 2010), and dominatethe character of valley winds (Petersen and Pellicciotti, 2011). Since katabatic wind cansignificantly enhance turbulent heat fluxes (Munro, 2006), a full evaluation of katabaticforcing is required in order to accurately spatially and temporally model melt for glacierswith a strong katabatic influence.This chapter addresses three main themes. First a distributed energy balance modelis used to calculate the volume of ice lost at Bridge Glacier during the 2013 melt season.Secondly, this model elucidates the main drivers of melt during the 2013 melt season, andexamines how those drivers vary both spatially and temporally. Finally, the work of Sheaand Moore (2010) is used to account for the strong influence of katabatic wind on thedistribution of turbulent heat flux parameters, and compare the melt projections of thismodel against more commonly used turbulent parameter distribution schemes.303.2. Instrumentation and Field MethodsEasting (m)Northing (m)Elevation (m asl)562600056280005630000563200056340005636000450000 455000 460000 465000Glacier AWSLake AWS   Ridge AWS     Ablation Stake1300 m1500 m1700 m1900 m2100 m2300 m2500 m2700 m2900 mFigure 3.1: Bridge Glacier study area, instrumentation, and 2013 terminus position. Contourintervals are 100 m.3.2 Instrumentation and Field MethodsBetween June 20 and September 12, 2013, data were collected from 3 automatic weatherstations (AWS) (see Figure 4.1). The majority of the data were taken from an AWS in-stalled on-glacier (Glacier AWS, 50◦49’02” N 123◦33’35” W, 1607 m a.s.l.). A second AWS(Ridge AWS, 50◦50’33” N 123◦31’23” W, 1639 m a.s.l.) was used for ambient temperaturedata, and was shielded from strong, persistent katabatic flow observed throughout the meltseason, measuring ambient meteorological conditions. A third weather station was locatedat the shore of the lake (Lake AWS, 50◦50’25” N 123◦30’07” W, 1393 m a.s.l.), along aterminal moraine, and was used to measure incoming longwave radiation. All meteorolog-313.2. Instrumentation and Field MethodsTable 3.1: Variables measured at on-site automatic weather stations. All variables were measuredat 1.75 m above the surface.Variable Symbol Units AccuracyGlacier AWSReflected Radiation K ↑ Wm−2 <1%Temperature Tg ◦C 0.21◦CRelative Humidity RH% % 2.5%Wind Speed u ms−1 4%Wind Direction ud ◦ 5◦Atmospheric Pressure P Pa 3.0 mbarRainfall Rate R mmhr−1 1.0%Melt Rate M Pa 0.05 md−1*Ridge AWSTemperature Tamb. ◦C 0.2◦CIncoming Shortwave Radiation K ↓ Wm−2 <1%Lake AWSIncoming Longwave Radiation L ↓ Wm−2 <1%*Approximate value, see Appendix A.1 for further information.ical variables were measured 1.75 m above the surface (see Table 3.1), and computed into10 minute means. For energy balance modelling, data was further aggregated into 1 hourmeans (sums for rainfall data) to reduce computation time.Elevation data was obtained using a 25 m resolution lidar DEM1. The elevation modelwas rescaled to 50 m to reduce computation time and reduce some digital artifacts in thedata. In order to ground-truth velocity measurements (see Chapter 4), and modelled meltestimates, 3 m ablation stakes were installed at 6 locations. Due to logistical challenges, andobtaining results that could also be usable for velocity estimates, the stakes were limitedto a relatively small spatial extent within the ablation area (see Figure 4.1). The stakeswere surveyed June 18, July 19, and September 13; however, in September, all the stakeshad melted out, and were unusable to estimate melt. As well, stake F, the lowest elevationstake, had moved into a heavily crevassed area where it was unable to be retrieved. Apressure transducer was installed at Glacier AWS in a 3 m water-filled borehole, whichmeasured 10 minute water depths. The borehole was re-drilled July 19. Post-processingfound that there was a mean diurnal variation of 10 cm water depth, but variations up to80 cm were observed. The data were further processed by selecting the maximum water1The authors gratefully acknowledge M. Demuth, C. Hopkinson, and B. Menounos for the LiDAR dataused in this study, which was collected as part of a C-CLEAR effort to develop LiDAR environmentalapplications, and funded in part by the Western Canadian Cryospheric Network (WCCN).323.3. Energy Balance Modelling Methods(a) Glacier AWS in June. Gimbel joint used to measure reflectedshortwave is visible on the left end of the tripod.(b) Ablation and velocitystake C, September 13, 2013.Figure 3.2: On-glacier instrumentation photosdepth for each day, allowing daily melt rates to be obtained (see Appendix A.1 for furtherinformation). Using the pressure transducer data, the date of the ablation stake (A) atGlacier AWS melting out was determined to be August 28 (the date at which the pressuretransducer also reached the surface). This date is assumed constant for all other stakes,although it is likely that there is an uncertainty of up to two to three days for individualstakes.3.3 Energy Balance Modelling MethodsGlacier surface melt (M) is expressed in metres of water equivalent melt per day and iscalculated asM =QMLfρi(3.1a)where QM is the sum of available energy at the surface (Wm−2), Lf is the latent heatof fusion, and ρi is the density of ice. Energy supplied to the glacier surface is taken aspositive, while energy flux away from the surface is negative. The sum of available melt333.3. Energy Balance Modelling Methodsenergy can be further divided asQM = Q∗ +QH +QE +QR (3.1b)where the net radiation (Q∗), sensible heat flux (QH), latent heat flux (QE), and rain heatflux (QR) can all contribute melt energy to the glacier surface. It is assumed that all energyfluxes occur at the glacier surface (Oerlemans, 2010; Munro, 2006), and subglacial melt isneglected. All calculations were completed using R Statistical Software unless otherwisenoted.Table 3.2: Variables used in this study.Symbol Variable ValueIo Solar Constant 1366.5 W m−2ϕ Glacier Latitude 50.81◦ψa Clear Sky Transmissivity 0.75Po Mean Sea Level Pressure 1013.15 hPaΓP Pressure Lapse Rate -0.12 hPa m−1ΓT Temperature Lapse Rate -0.006 ◦C m−1εterrain Terrain emissivity 0.95εice Ice emissivity 0.98αterrain Terrain albedo 0.17g Gravitational constant 9.81 m s−1ca Specific heat capacity of air 1006 J kg−1K−1cw Specific heat capacity of water 4200 J kg−1K−1k von Karman constant 0.41zo Roughness length of momentum 2.5× 10−3 mz Measurement height 1.75 mLv Latent heat of vapourization 2.26× 106 J kg−1Lf Latent heat of fusion 3.34× 106 J kg−1ρi Density of ice 917 kg m−3ρw Density of water 1000 kg m−33.3.1 Net RadiationNet radiation (Q∗) is calculated as the sum of incoming (↓) and outgoing (↑) shortwave(K) and longwave (L) radiation as follows:Q∗ = (K ↓ −K ↑) + (L ↓ −L ↑). (3.2)343.3. Energy Balance Modelling MethodsTo spatially distribute net radiation, it is rewritten asQ∗ = (S ↓ +D ↓)(1− αice) + (L ↓ −L ↑) (3.3)where shortwave radiation is separated into direct (S) and diffuse (D) components.Shortwave RadiationDirect radiation (Wm−2) for each point, S ↓, is calculated asS ↓= (1−D ↓K ↓)K ↓Kexi,jKex(3.4)by relating the ratio of direct shortwave radiation at Glacier AWS and observed incomingshortwave radiation to the ratio of potential direct radiation at the grid point (Kexi,j )and at Glacier AWS (Kex). The potential solar radiation at the top of the atmosphere iscalculated asKex = Io(RmR)2 cos θ (3.5a)where Io is the solar constant, Rm and R are the mean and instantaneous Sun-Earthdistances, and θ is the angle of incidence (see Appendix A.2 for solar geometry calculations).Potential direct radiation was corrected for slope geometry, and the angle of incidence (θ)wherecos θ = cosZ cosβ + sinZ sinβ cos(Ω−A) (3.5b)and Z is the elevation angle, β is the grid cell slope, Ω the solar azimuth, and A is theaspect.Diffuse shortwave radiation for all cells when Kex is greater than 0 is calculated asD = K ↓ (D ↓K ↓)φ+ αterrainK ↓ (1− φ) (3.6)where φ is sky view factor for each grid cell (Hock and Holmgren, 2005; MacDougall andFlowers, 2011).Due to the complications and heterogeneity involved in measuring the albedo (α) forthe surrounding terrain, a constant value of 0.17 was assumed, which is well within therange for dark, rocky surfaces (Oke, 1988). The errors associated with this assumptionshould be minor, since the sky view factor for most grid cells is greater than 0.9.353.3. Energy Balance Modelling MethodsThe ratio of incoming diffuse to shortwave radiation is calculated asD ↓K ↓=0.15 K↓Kex ≥ 0.800.929 + 1.134 K↓Kex − 5.111(K↓Kex)2 + 3.106( K↓Kex )3 0.15 < K↓Kex > 0.801.0 K↓Kex ≤ 0.15(3.7)where K↓Kex is the ratio of observed to potential shortwave radiation at the weather station.The method follows empirical relationships derived by Collares-Pereira and Rabl (1979),and re-confirmed with a different dataset by Hock and Holmgren (2005).In order to spatially distribute incoming shortwave radiation, each grid point is modelledas either shaded or sunny. A shading algorithm was implemented that calculates themaximum horizon angle for each grid point within a 10 km2 window, using 10◦ azimuthbins. At each time step, if the horizon angle is greater than the elevation angle (Z), thegrid point is shaded, and only receives diffuse radiation. For times when the horizon angleis smaller than elevation angle, the grid point receives both direct and diffuse radiation.Sky view factor was calculated using SAGA GIS software and a 25 m lidar DEM. Thealgorithm,φ =12pi∫ 2pi0cos2HdΩ, (3.8)integrates the maximum horizon angles (H) for each grid cell, for each azimuth angle (1◦interval). A maximum search window of 10 km2 was implemented to reduce computationtime.Longwave RadiationIncoming longwave radiation (L ↓) was measured directly at Lake AWS. In order to dis-tribute longwave radiation across the glacier, it is scaled by the sky view factor (φ) asL ↓= Linmeasuredφi,jφaws+ Lterrain(1− φi,j) (3.9)where additional longwave input is supplied by the surrounding terrain (Lterrain). Terraintemperature is assumed to equal air temperature, and the Stefan-Boltzmann equation isused with an emissivity (εterrain) of 0.95 (Oke, 1988).363.3. Energy Balance Modelling MethodsOutgoing longwave can be calculated asL ↑= εiceσT4 + (1− εice)L ↓ (3.10)where an ice emissivity (ε) of 0.98 is used, σ is the Stefan-Botlzmann constant, and T isthe ice temperature. It is assumed that during the melt season the glacier surface wasconstantly at the melting point (273.15 K).Albedo and Net RadiationReflected shortwave radiation was measured on-glacier, in the ablation area, throughoutthe melt season using a gimbal joint designed to keep the pyranometer pointed directlydownwards. Due to difficulties ensuring a pyranometer attached to the glacier weatherstation would remain level, incoming shortwave measurements were taken from the off-glacier Ridge AWS. The differences in shading were found to be minor. To eliminate smalldiscrepancies in shading, uneven cloud patterns, and low solar angle errors (Oerlemans,2010), ice albedo (α) is calculated asαdaily =∫K ↑ dt/∫K ↓ dt (3.11)for all daylight hours, and is assumed constant throughout the day. The weather stationwas located on bare ice for the entire study, and since this study only concerns ice loss,snow albedo was not considered.3.3.2 Turbulent Heat FluxesSensible heat (QH) and latent heat (QE) are important components of melt energy avail-ability, and can be expressed according to the bulk transfer approach asQH = ρaircaCu(Ta − Tg) (3.12)QE = ρairLvCu(0.622(ea − eg)P) (3.13)where cair is the specific heat capacity of air, C is a turbulent transfer coefficient, u isthe windspeed, T is air temperature, Lv is the latent heat of vaporization, e is the vapourpressure, and P is the atmospheric pressure at the surface (see Table 3.2). The gradients for373.3. Energy Balance Modelling Methodsair temperature and vapour pressure are denoted where the subscript g refers to on-glacier,and a refers to off-glacier ‘ambient’ conditions.In order to reflect the significant forcing of katabatic winds on turbulent flux variables,glacier morphometry is used to model the extent and magnitude of katabatic forcing. Usingflow path lengths, defined as the average flow distance to a given point from an upslopesummit or ridge, statistical relationships are formulated to spatially distribute temperatureand wind speed across the glacier. Flow path lengths (FPL) for the glacier were calculatedusing a 25 m lidar DEM and a routine in the Terrain Analysis - Hydrology module ofSAGA GIS (Quinn et al., 1991), and is shown in Figure 3.3. Since downslope flow driveskatabatic forcing, flow path lengths allow for a quantification of the magnitude of relativeforcing across the glacier (Shea and Moore, 2010).Easting (m)Northing (m)450000 455000 460000562600056280005630000563200056340005636000020004000600080001000012000Flow Path Length (m)Flow Path Length (m)Figure 3.3: Flow Path Lengths for Bridge Glacier.383.3. Energy Balance Modelling MethodsDistributing TemperatureThe magnitude of katabatic forcing was modelled as a function of the temperature dif-ference (∆T ) between on-glacier (Glacier AWS) and off glacier (Ridge AWS) (outside thekatabatic boundary layer) weather stations. Temperature differences were separated intoupslope (northeasterly) and downslope katabatic (southwesterly) periods, based on thewind directions of Glacier AWS. Linear regressions of both groups against ambient (off-glacier) temperature (Figure 3.4, Ta, from Ridge AWS) show a strong linear increase in∆T with increasing ambient temperatures for katabatic flows, indicating the magnitudeof katabatic forcing increases with increasing ambient temperatures. Conversely, ∆T doesnot significantly vary as a function of ambient temperatures during upslope flow, althoughambient temperature above 10 ◦C during these periods was rare. The elevation of bothweather stations are within 100 m, and small corrections to potential temperature, usinga −6◦C km−1 lapse rate, did not produce a meaningful difference in the linear fit.5 10 15 20051015Ambient Temperature (°C)Ambient − Glacier Temperature ( °C)Downslope WindUpslope WindFigure 3.4: Katabatic on-glacier cooling, as a function of ambient temperature from Ridge AWS.393.3. Energy Balance Modelling MethodsOn-glacier temperature for each grid point is modelled as a function of the katabatictemperature depression whereTg = Ta − (k1Ta + ∆T∗) (3.14a)and k1 is the magnitude of katabatic forcing for each point on the glacier, Ta is the off-glacier (ambient) air temperature, and ∆T ∗ is the threshold temperature differential atwhich katabatic flow is observed.The magnitude of katabatic forcing for each point on glacier is calculated by usingstatistical coefficients and glacier flow path lengths (FPL) wherek1 = β1 exp(β2FPL) (3.14b)and the empirical coefficient β2 is taken from Shea (2010). The coefficient is derivedusing data from multiple elevations and multiple glacier basins in the southwestern CoastMountains (see Table 3.3). The coefficient β1 is calculated asβ1 =mTexp(β2FPLaws)(3.14c)where FPL and the degree of katabatic forcing (mT ) values are extracted from the locationof Glacier AWS.The coefficients mT and ∆T ∗ can be derived from linear regression (Figure 3.4), where∆T = mTTa + ∆T∗ (3.14d)and the difference between on-glacier and off-glacier air temperatures (∆T ) is regressedagainst ambient air temperatures (from Ridge AWS).In the rare case that wind direction is upslope, temperatures are distributed usingon-glacier temperature, Tg, and a standard temperature lapse rate of −6◦C km−1.403.3. Energy Balance Modelling MethodsTable 3.3: Statistical coefficients used in FPL modelling.Variable Valueβ1 3.90× 10−1β2 4.43× 10−5*β3 6.70× 10−2β4 -3.39× 10−1*∆T ∗ -1.67 ◦Cu∗ 1.08 m s−1*From Shea and Moore (2010)Distributing Vapour PressureVapour pressure is calculated from measured relative humidity and saturation vapour pres-sure (es). Saturated vapour pressue is calculated using Teten’s formula, wherees ={6.11× 10(7.5T237.7+T ) T ≥ 06.11× 10(9.5T265.5+T ) T < 0(3.15a)and T is the air temperature. Vapour pressure (e) can be related to relative humidity(RH%) as follows:e = esRH100. (3.15b)Since only a weak relationship was found between on-glacier vapour pressures and off-glaciertemperatures, relative humidity measured at Glacier AWS is held spatially constant acrossthe glacier for each timestep, and saturation vapour pressure is calculated from distributedon-glacier air temperatures.Distributing Wind SpeedGlacier wind speed was distributed in a similar fashion to air temperatures, as a functionof katabatic forcing and ambient temperatures, following Shea (2010). On-g;acier windspee is regressed against ambioent air temperature. For downslope winds, wind speedincreases linearly with increasing ambient air temperature, while upslope wind speeds showno significant change with ambient temperature (see Figure 3.5).On-glacier wind speed for each grid point is modelled as a function of off-glacier tem-413.3. Energy Balance Modelling Methods5 10 15 2002468Ambient Temperature (°C)Glacier Windspeed (m/s)Downslope WindUpslope WindFigure 3.5: Wind speed as a function of ambient (off-glacier) air temperature.perature whereug = u1Ta + u∗ (3.16a)and u1 is the magnitude of katabatic forcing for each point and u∗ is the katabatic windspeed when the off-glacier temperature is 0◦C.The magnitude of katabatic forcing on wind speed for any point on the glacier calculatedusing statistical coefficients and glacier flow path length asu1 = β4 + β3 log(FPL) (3.16b)where the empirical coefficient β4 is taken from Shea and Moore (2010). The coefficient β3is calculated asβ3 =mu − β4log (FPLaws)(3.16c)423.3. Energy Balance Modelling Methodswhere FPL and the degree of katabatic forcing on wind speed (mu) values are extractedfrom the location of Glacier AWS.The coefficients mu and u∗ can be derived from linear regression (Figure 3.5), whereug = muTa + u∗ (3.16d)and the on-glacier wind speed is regressed against off-glacier air temperatures. When thewind direction is upslope, on-glacier wind speed is held constant, using measured windspeeds from Glacier AWS.Physical Limitations of Flow Path Length ModellingThis model has the major drawback in that FPL < 1 m results in a negative wind speed. Inorder to mitigate this drawback, points with no katabatic forcing (u1 ≤ 0) are assumed tohave the same wind speed as the off-glacier Ridge AWS. The vast majority of these pointsoccur near ridges and peaks, and are unlikely to be affected by katabatic flow. Secondly, thephysical meaning of u∗ breaks down at low temperatures, theoretically producing katabaticwinds when the ambient air temperature is colder than the glacier surface (Tamb. < 0◦C).Since sub-zero temperatures were not observed during the study period, the potential of anon-linear fit at low temperatures cannot be tested. It is likely, however, that the surfaceof the glacier would be less than 0◦C during events with low ambient air temperatures,and would not be melting substantially.Similar limitations exist in distributing on-glacier temperatures. The calculated value of∆T ∗ is < 0 ◦C, suggests that a katabatic temperature depression occurs even when glacierair temperatures are higher than ambient values. This phenomenon is not observed duringthe study period, and seemingly contradicts expected katabatic behavior. However, forcingthe linear fit through the expected value of ∆T ∗ = 0◦C substantially reduces the linear fit.Similar to wind speed limitations, it is likely that katabatic forcing is uncommon at lowertemperatures during the summer melt period, an assumption supported by temperaturesof roughly and below 5◦C, exhibiting predominantly upslope wind directions.Finally, while β1 and β3 were recalculated, and could be compared with values derivedfrom Shea and Moore (2010), variables β2 and β4, which scale FPL magnitude spatially,could not be tested with the current 2013 field dataset since it requires data from multipleelevations, and were therefore assumed to be the same.433.3. Energy Balance Modelling MethodsTurbulent Transfer CoefficientUsing the bulk aerodynamic method (Moore, 1983), the turbulent transfer coefficient, C(dimensionless), is calculated using a stability correction (Θ) and surface roughness lengths(zo and zx) as follows:C = Θk2[log (z/zo) log (z/zx)]−1. (3.17a)The magnitude of the stability correction is calculated asΘ ={(1− 5.2Rb)2 Rb > 0(1− 16.0Rb)0.75 Rb < 0(3.17b)where positive bulk Richardson Numbers (Rb) indicates a stable atmosphere. The bulkRichardson Number is calculated asRb =gTg + 273.15(Ta − Tg)zu2(3.17c)where g is the gravitational constant, and z is the height of measurement (m).The roughness length for momentum (zo) is taken as 2.5 mm for ice (Munro, 1989),while the roughness length for temperature and vapour pressure (zx) is calculated as zo/300following Hock (1998). Both roughness lengths are assumed constant in space and timewithin the ablation area.Rain Heat FluxAlthough relatively minor, energy supplied to the surface due to rain is calculated asQR = ρwcwRTR (3.18)where R is the rainfall rate (ms−1), TR is the temperature of rain, which is assumed thesame as the ambient (off-glacier) air temperature, and ρw and cw are the density andspecific heat of water (see Table 3.2).3.3.3 Snowline RetreatSince an adequate spring snow survey was not logistically possible, and proxy snow sur-veys were considered unreliable due to large precipitation gradients in the region (Stahl443.3. Energy Balance Modelling MethodsJun Jul Aug Sep14001600180020002200Snowline (m asl)Figure 3.6: Modelled snowline retreat from Landsat data. Dashed line is loess smoothed al., 2008; Shea and Moore, 2010), snowline retreat was reconstructed using Landsatimages obtained from the LandsatLook Viewer ( Mul-tiple measurements were taken for each scene (9 unique dates throughout the study period),and averaged to produce a basin-wide snowline elevation. Interpolation between points wasachieved using the loess smoothing function in R (Figure 3.6). Substantial variability ex-ists in the snow coverage in the north and south channels; however, the errors associatedwith this interpolation should be relatively minor in the overall summer ablation. In orderto calculate the volume of ice melt, any melt above the snowline, at each time step, is setto Complementary Melt Rates and MeteorologyAs an alternative method of deriving and constraining 2013 melt, melt rates are alsomodelled using a daily temperature index model (TIM) (Hock, 2003; Shea et al., 2009) asfollows:M = DDF × Ta. (3.19a)453.4. Results and DiscussionThe degree day factor (DDF) is calculated asDDF =∑M∑T+a(3.19b)where T+a is the positive off-glacier air temperature (from Ridge AWS), and M is themeasured seasonal ablation from ablation stakes near Glacier AWS.3.4 Results and Discussion3.4.1 Volumetric Ice MeltThroughout the study period of June 20 to September 12, 2013, the model predicted icemelt of 1.0 m (w.e) near the ELA, to 5.9 m (w.e) near the terminus. Melt is greatestalong the main tongue, displaying the signature of enhanced sensible heat flux due topersistent katabatic flow (Figure 3.7). The southernmost tributary glacier shows relativelylow melt rates relative to its elevation, due predominantly to the fact that it is shelteredfrom high winds, and its north-facing aspect allows for substantial shading throughout themelt season.The volume of ice loss during the study period was 0.124 km3. By the end of the meltseason, the ablation area was 27.7 km2, accounting for 34% of the total glacier surface area(81.3 km2).3.4.2 Relative Energy ContributionsDuring the study period, net radiation (Q∗) was the primary driver of ice melt (Figure 3.8).Due to the relatively low measured albedo of ice (values between 0.4 and 0.2), net radiationis strongly controlled by incoming shortwave radiation (K ↓). The relative importance ofQ∗ decreases throughout the melt season from roughly 80% of melt energy to less than60% by September. This decrease is driven primarily by the decrease in sunlight hourslater in the season and a decrease in the incidence angle.Sensible heat flux represents between 20% and 40% of the melt energy supplied tothe ice surface during the melt season. Driven primarily by peaks in wind speed andtemperature, the relative importance of QH coincides with the warmest days in late Julyand early August. While a significant decrease in Q∗ is observed during storm events, bothsensible and latent heat show only very minor decreases in supplied energy, suggesting that463.4. Results and DiscussionNorthing (m)449000 451000 453000 455000 457000 459000 4610005625000562700056290005631000563300056350005637000Easting (m)0 1 2 3 4 5 5.9Summer Ice Melt (m w.e.)Figure 3.7: Modelled ice melt (water equivalent) for the study period from June 20 to September12, 2013.inclement weather does not significantly change the supply of turbulent heat to exposedice on the glacier.Latent heat flux was a relatively minor source of energy throughout the melt season.The relatively dry air on the lee side of the Coast Mountains made for little condensationon the glacier surface in July and early August. Evaporation was predicted for several daysin mid-July, leading to latent heat transfer away from the glacier surface, indicating a lowair moisture (roughly 5.98 hPa). Late August rain storms supplied high moisture contentin the air, leading to slightly elevated latent heat flux and notable rain events supplyingadditional energy to the glacier surface.473.4. Results and DiscussionJul Aug Sep0100200300400Melt Energy (Qm)Net Radiation (Q*)Sensible Heat Flux (Qh)Latent Heat Flux (Qe)Rain Heat Flux (Qr)Average Daily Melt Energy (Wm−2 )Figure 3.8: Average daily energy fluxes contributing to melt within the ablation area.3.4.3 Modelled Glacier Temperature DistributionDuring hot, sunny days (see Figure 3.9), modelled temperature shows a large depressionof up to 14 ◦C on-glacier, relative to off-glacier values. The magnitude of katabatic forcingis substantially higher further down-glacier, while the uppermost reaches of the glacier arerelatively insensitive to katabatic forcing, and reflect temperatures that would be expectedfrom a standard environmental lapse rate. The magnitude of this forcing is scaled byoff-glacier temperatures (Tamb.), where warmer ambient temperatures create a strongerkatabatic forcing, and a larger temperature depression on-glacier (see Equation 3.14a).3.4.4 Comparative Model PerformanceAlthough there was insufficient data to fully ground-truth the distributed energy balancemodelled (DEBM) melt, the model was tested against sparse data from several sourcesto constrain the range of uncertainty in predicted melt. Modelled results were comparedagainst a point energy balance model using the same modelling methodology, and derivedsolely from observations at Glacier AWS (i.e. no FPL distributed T , e, and u). The DEBM483.4. Results and DiscussionKatabatic Temperature Depression 2013−07−31 13:00:00Easting (m)Northing (m)Temperature (C)562600056280005630000563200056340005636000450000 455000 4600008101214161820Figure 3.9: Air temperature distributions for 13:00, July 31, 2013.was also evaluated against daily melt measured at Glacier AWS using a pressure transducerin a water-filled borehole. Summer cumulative melt was further constrained by 6 ablationstakes.A comparison of cumulative melt (Figure 3.10) shows that the Point EBM predicts thehighest melt (5.10 m w.e.), while both measured (4.39 m w.e.) and DEBM (4.67 m w.e.)agree relatively well after 67 days (June 20 - August 28). All three methods reproducea sharp decline in the melt rate for several days in late June. However, measured valuesshow a substantial decrease in the melt rate in mid-July that is not reproduced in eitherenergy balance model. This mid-July decrease in the melt rate corresponds with very largeobserved diurnal variations in the measured borehole depth (see Appendix Figure A.1),suggesting that the borehole may not have been refilling completely during the period,obscuring the melt rate signal. All three methods also show a mid-August decrease in themelt rate, although the measured value shows significantly more scatter, most likely alsodue to changes in the water level within the borehole.The Point EBM produces slightly higher cumulative melt and a consistently higher493.4. Results and DiscussionJul Aug Sep0123456Cumulative Melt (m w.e.)DEBMPoint EBMMeasuredFigure 3.10: Cumulative melt at Glacier AWS from the Distributed Energy Balance Model(DEBM), a point Energy Balance Model, and from measured daily melt taken from a pressuretransducer in a borehole.melt rate throughout the summer study period (Figure 3.11). Net radiation is consistentlyslightly higher in the Point model than the DEBM, due to the fact that the point modeldoes not take into account the surface slope at Glacier AWS, leading to a slightly higherestimate of incoming shortwave radiation.Evaluation of the sensible and latent heat fluxes from the DEBM and Point EBM givesa sense of how well the FPL reproduces turbulent heat parameters at Glacier AWS. Whilesensible heat is well reproduced by the DEBM model, latent heat is less well reproduced,showing a consistently lower estimate than from the Point EBM. This suggests that holdingrelative humidity constant across the glacier, and calculating vapour pressure using FPLdistributed temperatures, consistently underestimates the vapour pressure at Glacier AWS.Although further in situ data are necessary to confirm this finding, it is likely that thisunderestimate is greater near the terminus, where the FPL model predicts the strongestkatabatic winds and lowest air temperatures, in turn producing the lowest vapour pressures(Shea and Moore, 2010).Summer melt measured at the ablation stakes is slightly less than that predicted fromthe DEBM; however, the discrepancy is less than 0.20 m w.e. (3.5%) for all stakes except503.4. Results and Discussion0.04 0.08 Rate (m/day)PointDistributed50 100 150 200 250 30050100150200250300Net Radiation (Wm−2)PointDistributed20 40 60 80 100 12020406080100120Sensible Heat Flux (Wm−2)PointDistributed−10 10 20 30 40 50 60−100102030405060Latent Heat Flux (Wm−2)PointDistributedFigure 3.11: Comparison between Point EBM and DEBM predictions at Glacier AWS.D (Figure 3.12). Although all ablation stakes are assumed to have melted out at thesame time, it is possible that there was up to several days between each stake meltingout, leading to underestimates of the total melt during the period. Furthermore, surfaceheterogeneity, such as small meltwater pools forming at the surface, or patches of dirty ice(which would lower the albedo) were observed across the glacier, and could have affectedthe melt rate at individual stakes, including stake D, which shows substantially less meltthan nearby stakes (12%). Overall, total melt at the stakes was well reproduced by theDEBM, suggesting that the model accurately captures the main drivers of melt in this partof the ablation area.513.4. Results and DiscussionA B C D E FAblation Stake0100200300400500Summer Melt (cm w.e.)ModeledMeasuredFigure 3.12: Cumulative ice melt between June 20 and August 28 from DEBM and measuredablation stakes.3.4.5 Comparison of FPL Model UncertaintyIn order to test the sensitivity of the distributed energy balance model to the choice ofturbulent parameter distribution, the flowpath length model (see Sections 3.3.2 and 3.3.2)was compared with other frequently employed temperature and wind speed distributionmodels (e.g. MacDougall and Flowers (2011); Hock and Holmgren (2005)). The FPL modelis tested against turbulent heat flux calculations obtained by using a spatially constant windspeed (measured at Glacier AWS) distributed across the glacier, and by using a standardtemperature lapse rate (see Table 3.2) to distribute measured on-glacier temperatures (fromGlacier AWS).Comparing wind speed distribution models, the FPL model predicts measurably lesstotal melt than using a constant wind speed (Figure 3.13). While the main (southern)channel of Bridge Glacier seems relatively insensitive to the choice of wind speed distribu-tion model, the model predicts over 1.5 m of additional melt for the north-flowing channelfor a constant wind speed. This additional ice melt prediction is due to the fact that theFPL model predicts significant sheltering in this side channel, shielding it from the muchgreater wind speeds experienced near the terminus and Glacier AWS. Although not tested,523.4. Results and Discussion450000 452000 454000 456000 458000 460000 4620005624000562600056280005630000563200056340005636000Constant WindspeedNorthing (m)450000 452000 454000 456000 458000 460000 462000Constant Temperature Lapse Rate562600056280005630000563200056340005636000Easting (m)-2-1012Model−FPL Model Predicted Melt (m w.e.)Figure 3.13: Plots of the difference in ice melt (in m w.e.) between spatially constant windspeedand temperature lapse rate distributions, and FPL modelled temperature and is highly likely that wind speeds are much higher at Glacier AWS than further up-glacierdue to a combination of being located in a constrained area with high sidewalls (allowingwinds to be funnelled), and close to the terminus, where katabatic winds are expected tobe strongest. This FPL model also suggests that, if using a constant wind speed across theglacier, the model will be very sensitive to the location of the wind speed measurements.Temporally, both sensible and latent heat fluxes increase across the glacier using aconsistent wind speed throughout the melt season (Figure 3.14). This is most observableduring late July through early August, coinciding with the highest off-glacier air temper-atures and highest observed on-glacier wind speeds. Since a constant wind speed modeldoes not take into account any potential for sheltering or spatial variation in wind speeds,it likely overestimates observed wind speeds in sheltered areas.The FPL model for katabatic temperature distribution predicts more melt overall, com-pared to a constant temperature lapse rate distribution. Near the terminus there is littlechange in modelled melt, except for immediately at the terminus, where a constant temper-ature lapse rate predicts and additional 0.3 m w.e of melt throughout the summer (roughly533.4. Results and Discussion5% m). Above Glacier AWS, higher melt is predicted using the FPL model, with a pre-dicted difference of up to 1 m (18%). This discrepancy is due to the fact that the FPLmodel predicts that the magnitude of katabatic forcing generally decreases with increasingelevation, leading to slightly warmer temperatures further up-glacier. Distributing tem-peratures by using an on-glacier weather station and standard temperature lapse rate willpredict much lower temperatures than would be expected from distributing temperatureswith off-glacier weather stations, even in areas such as peaks and ridges, where katabaticcooling would not be expected. By using on-glacier temperatures from an AWS in an areastrongly affected by katabatic cooling, a standard temperature lapse rate will lead to asignificant underestimate of higher elevation temperatures, and subsequently sensible heatfluxes. It also suggests that due to the potential for differential katabatic forcing, the loca-tion of an on-glacier weather station will ultimately have a significant effect on modelledtemperature distributions using a standard lapse rate.Sensible heat flux derived from a constant temperature lapse rate is smaller than pre-dicted from a FPL model throughout the melt season (Figure 3.14). This trend is mostapparent when katabatic forcing is greatest, due to high ambient temperatures and on-glacier wind speeds. During cool or rainy weather, the difference in sensible heat flux fromboth models is minor, suggesting that katabatic forcing is negligible during poor weather.This assumption is supported by less frequent observations of downslope winds duringthese periods.Since the relative humidity was held constant across the glacier, changes in the latentheat flux are themselves a product of air temperatures. Since the constant temperaturelapse rate model predicts very low temperatures at higher elevations during strong katabaticflow conditions, as a by-product, it also predicts very dry conditions, most apparent duringthe warm, high katabatic flow period of late July through mid-August. Since the constanttemperature lapse rate most likely underestimates air temperatures at higher elevations,it is likely that it also underestimates vapour pressure at those elevations. This suggeststhat using a constant temperature lapse rate coupled with a constant relative humidityacross the glacier is sensitive to measurement location, and has the potential to significantlyunderestimate the latent heat flux, specifically at higher elevations, where katabatic forcingis expected to be less significant.543.4. Results and DiscussionLatent Heat Flux (Wm−2 )-20020406080100Flowpath Length Constant Windspeed Constant Temperature Lapse RateJul Aug SepSensible Heat Flux (Wm−2 )-20020406080100Figure 3.14: Changes in turbulent heat fluxes using different models for spatially distributingtemperature and wind speed.3.4.6 Climatic SensitivityResearch has shown the length of time ice is exposed is by far the most critical factor indetermining cumulative ice melt during the summer (Moore and Demuth, 2001). In orderto further position the 2013 study period within a longer term context, the DEBM was re-run with climatic perturbations to simulate the relative importance of climatic variability.The snowline position throughout the summer is altered, by delaying the retreat by oneand two weeks, as well as accelerating it by two weeks. Climatic variability is also testedby increasing ambient temperatures for each timestep by 1◦C.In all cases, the absolute changes in total volumetric ice loss are smaller than thecalving flux (see Chapter 5). The absolute change in volume of ice loss due to a 2 weekdelay snowline retreat; however, is comparable (72.5%) to the calving flux (Figure 3.15).While snowline retreat (and subsequent ice exposure) has long been demonstrated as oneof the most important controls on glacier mass balance, the magnitude of this variabilityis well within the expected annual observed variability in spring snowpack and snowlineretreat.553.4. Results and Discussion0. ° CSnowline RetreatCalvingFlux1 wk. later 2 wks. later 2 wks. earlierVolume of Ice (km3 )Additional GainAdditional LossFigure 3.15: Calving loses relative to climatic variability scenarios affecting the melt rate. TheDEBM is modified by increasing ambient air temperatures by 1◦C, as well as delaying the snowlineretreat by 1 and 2 weeks, as well as accelerating it by 2 weeks.In order to assess the relative sensitivity of Bridge Glacier to future changing climaticconditions, the DEBM also run with increased air temperature (Figure 3.16. Throughoutthe summer, the ambient temperature (at Ridge AWS) was increased by 1◦C for eachtimestep. Results show that this increase in off-glacier air temperature yields an additionalmelt of over 0.25 m w.e. in the lower bench of the glacier, and an additional 5.73 ×10−3km3 of ice loss.Higher ambient temperatures increase melt directly due to higher on-glacier tempera-tures, while it also indirectly enhances melt by intensifying katabatic flow. This feedbackis borne out in the model by enhanced melt along the largest flow path lengths along themain southern tongue of the glacier, while sheltered areas show a minimal increase in melt.This feeback mechanism also suggests that ambient temperatures have a stronger effect onmelt through their ability to enhance wind speeds rather than air temperatures themselves.While an increase in ambient air temperature does produce a meaningful increase inice melt, the impact is comparatively modest when taking into account the magnitude ofchanges in other melt variables. This increase in ambient air temperature only results in563.4. Results and DiscussionNorthing (m)449000 451000 453000 455000 457000 459000 4610005625000562700056290005631000563300056350005637000Easting (m)0.05 0.15 0.25Additional Ice Melt (m w.e.)Figure 3.16: Predicted additional melt from a 1◦C increase in air temperature during the 2013melt season.a 5% increase in cumulative melt. While this increase is significant, it could also be easilyachieved by a small increase in incoming radiation, due to fewer cloudy days, or a fewadditional melt days.Additionally, since it has been shown that sensible heat flux has a larger impact onsnow rather than ice (Munro, 2006), an increase in air temperatures would acceleratethe snowline retreat and provide additional ice melt on top of what has been estimatedhere. Ultimately, however, this test shows that the Bridge Glacier’s climatic sensitivity isrelatively large, and small changes in summer air temperatures, winter precipitation, ordecreased cloudiness can all provide significant additional ice loss.573.5. Conclusions3.4.7 Model UncertaintyTo quantify the level of uncertainty in melt and volume calculations presented above,the DEBM uncertainty is calculated based on changes in model parameters and in situmeasurements of ablation. The DEBM was re-run using both a steady temperature lapserate, and constant wind speed (see Section 3.4.5). The volumetric difference is computed,relative to the FPL model used, and the difference is found to be 0.4% for the standardtemperature lapse rate, and 6.5% for the constant wind speed. These potential errors arecombined with the observed error at ablation stakes near Glacier AWS. The observed errorat the ablation stakes is calculated at 5.5% (0.29 m), giving a net uncertainty in melt of12.4%.3.5 ConclusionsBetween June 20 and September 12, 2013, up to 5.9 m w.e. of ice was lost in the ablationarea, leading to a total melt volume of 0.124 km3. The melt rate was highest in Julythrough mid-August, corresponding with warmest temperatures, and predominantly clear-sky conditions. Modelled melt during the study period corresponds well with measuredmelt from both ablation stakes, and cumulative daily melt from a pressure transducer in awater-filled borehole at Glacier AWS.Climatic sensitivity analyses suggest that given the current geometry of Bridge Glacier,the volumetric ice loss variability is in the range of ±0.03 km3. Ice loss is most stronglycontrolled by the length and timing of ice exposure, while it is only moderately sensitiveto changes in mean air temperature.During the study period, net radiation, driven primarily by incoming shortwave radi-ation, is the main source of melt energy across the ablation area. Net radiation accountsfor almost 80% of the ice melt in June and July, becoming less important throughout thesummer as the number of sunlight hours decreases. Net radiation is very variable, showinglarge decreases during summer storms, a trend that is not mirrored in latent and sensibleheat fluxes.Contrasting with many other DEBM temperature and wind speed distribution mod-els (MacDougall and Flowers, 2011; Hock and Holmgren, 2005; Anslow et al., 2008), thisformulation considers observed katabatic forcing (Petersen and Pellicciotti, 2011; Jiskootand Mueller, 2012). Katabatic flow exerts a significant decrease in on-glacier air temper-583.5. Conclusionsature and an increase in on-glacier wind speed, while both trends are strongly enhancedby high off-glacier air temperatures. Using flow path lengths, the magnitude of katabaticforcing is modelled, allowing for an account of katabatic influence on air temperature. Themodel produces results that allow for significant cooling in the lowest on-glacier reaches,while producing high elevation results that are in line with what would be expected fromoff-glacier standard temperature lapse rates.59Chapter 4Measuring Glacier Velocity withTime Lapse Imagery4.1 IntroductionMeasuring glacier velocity has long been an important part of glacier study, and is anintegral variable in calculating calving fluxes in tidewater and lacustrine systems (Meier andPost, 1987), strain rates and estimating basal water pressures. Original measurements ofglacier velocity involved repeat manual surveys of on-glacier targets. While this techniquebenefits from a simple design and low instrumentation requirements, it is less desirable froma glaciological perspective because it provides a lower time resolution since the temporalresolution depends on the frequency of visits to the glacier.With the increased reliability in Differential Global Positioning Systems (DGPS), pro-viding accuracy of up to 10 cm, the temporal resolution of a velocity dataset can be sub-stantially augmented by installation of a DGPS on-glacier (Danielson and Sharp, 2013).While this technique can yield high accuracy measurements, translating into hourly tosub-daily resolutions, it only provides displacements for one location, and is limited by thehigh costs of instrumentation and maintenance.The use of remotely sensed images, such as Landsat or SPOT, to track glacier flowovercomes the spatial constraints of point measurements. However, it too suffers from lowtemporal resolution, dictated by the frequency of repeat images. In environments whereheavy clouds and consistent precipitation patterns are expected, the number of usableimages is further reduced due to cloud-cover obscuring the glacier surface.As a relatively low cost alternative, time lapse photography offers a remote option thatcan provide sub-daily velocity measurements over an array of targets distributed acrossthe field of view of the camera. The use of time lapse photography to calculate glaciervelocity has been in practice for several decades, and has yielded valuable data, albeit604.1. Introductionwith significant challenges, such as camera displacement, inclement weather, and poorbattery performance (Krimmel and Vaughn, 1987; Harrison et al., 1992). However, recenttechnological advances increasing image quality and digital storage, as well as the design ofautomatic camera triggers and robust weather enclosures, have made repeat photographya convenient and reliable method of glaciological observation.Monoscopic (single) cameras are the most basic time lapse set-up that can yield dis-placement values, although without additional surveyed points, displacement values areonly relative to the camera. To overcome this limitation, ground control points must besurveyed, or estimated, in order to calculate true displacements. Another method that hasbeen implemented to avoid this complication involves using a DEM to orthorectify eachimage (Corripio, 2004), yielding a distributed velocity time series from a monoscopic timelapse dataset (Rivera et al., 2012). This method is reliant on having an accurate DEM,which is unfeasible in many instances due to substantial glacial thinning. Since errors inglacier surface topography will propagate into orthorectification, and subsequently, mea-sured displacements, values derived using a DEM are, in some cases, no more accuratethan assuming a planar glacier surface.To overcome this limitation, stereoscopic (two overlapping images) photography presentsan alternative that allows for the derivation of real world coordinates directly. This methodis preferable for some applications (Sund et al., 2011; Rivera et al., 2012) and has the addedbenefit of estimating glacier thinning rates. However, this method requires additional geo-metric and geographic inputs, a full account of lens distortion to accurately produce directdisplacements, as well the increased cost of deploying two systems.More recent advances have employed computer vision algorithms instead of manualimage-pixel tracking (Ahn and Howat, 2011). Automated algorithms can provide greaterspatial coverage, at a fraction of the time required for manual pixel tracking, and in somecases can provide more accurate solutions (Eiken and Sund, 2012). However, the dynamicnature of image lighting, changing glacial features throughout the season, and inclementweather end up reducing their effectiveness and reliability unless substantial image pre-processing is applied (Ahn and Howat, 2011).The purpose of this chapter is to present a method for producing glacier velocity mea-surements from a monoscopic time lapse camera using manual feature tracking, outlininga low-cost method with minimal external input requirements. Secondly, this chapter willpresent results from two monoscopic cameras and quantify the errors in this method. Fi-nally, this chapter aims to produce daily velocities at two locations on Bridge Glacier that614.2. Study Site and Instrumental Setupcan be used to quantify the calving flux and strain rates (see Chapter 5).4.2 Study Site and Instrumental SetupThree cameras were installed at Bridge Glacier between June 17 and September 14, 2013.Two cameras were installed roughly perpendicular to glacier flow, while a third was locatedwith a full view of the terminus, roughly 2 km from the terminus, and was not used tomeasure glacier velocities. One camera (referred to as ‘Terminus Camera’) was installedperpendicular to the ice front along Bridge Lake, and captured the floating terminus of theglacier. A second camera (referred to as ‘Nunatak Camera’) was installed perpendicular toflow on top of a prominent (former) nunatak, with a view along the south (main) branchof the glacier (see Figure 4.1 and 4.3). Both cameras were programmed to shoot every 30minutes between 6 am and 9 pm, producing 34 images per day.459000 460000 461000 462000 4630005628000562900056300005631000563200014001500160017001800190020002100             Elevation          (m asl)Bridge Glacier Time Lapse Camera LocationsEasting (m)Northing (m)Bridge LakeBridge GlacierNunatakTerminusFigure 4.1: Bridge Glacier extent at the end of the 2013 study period, and the location of the twotracking time lapse cameras, Terminus Camera and Nunatak Camera. Contour intervals are 40 m624.2. Study Site and Instrumental SetupThe time lapse camera setup was designed ready-made by Harbortronics. A CanonEOS Rebel T3i camera with a 10-20 mm wide-angle lens was installed in a custom-madeenclosure powered by two lithium batteries, which are recharged by a solar panel. Thecamera is controlled using custom command line software designed to program a DigiSnap2000 intervalometer.(a) Nunatak Camera, September 12, 2013. (b) Velocity Stake survey, June20, 2013.Figure 4.2: Photos of time lapse camera and velocity stakeCameras were installed on June 17, 2013. On July 18, 2013, the cameras were checked,and the SD card was changed. On September 14, 2013, the cameras were uninstalled andremoved from the field. During the June and July visits, the location and elevation of eachcamera was checked using a handheld Garmin GPS. Other geographic parameters (Table4.1) were derived from BC TRIM and Google Earth SPOT data.Each camera produced over 3000 images during the study period. Time lapse imageswere culled before feature tracking was applied. Daily images from Nunatak Camera at18:00 hrs were selected. Terminus Camera images were selected bi-daily, at 06:00 and21:00 hrs. Images from both cameras were manually checked, and several images werereplaced by successive images in the time series due to fog, poor lighting, water droplets,or even in one case, a marmot peering into the inclosure. Images were selected to minimize634.2. Study Site and Instrumental Setup(a) Nunatak Camera view. The blue line represents 595 m at the centrelinepoints (1030m from the camera).(b) Terminus Camera view. The blue line represents a scale of 690 m at the icefront (970 m from the camera).Figure 4.3: Camera views and tracked points from images taken on June 17, 2013. Points in theforeground correspond to ground control points.644.3. Measuring Glacier Velocity from Time Lapse ImagesTable 4.1: Time lapse camera variables.Terminus Camera Nunatak CameraGlacier Slope 0.060 0.075Flow Direction (δ) 10.6◦ 13.0◦Focal Length 15 mm 16 mmHorizontal Angle of View 100.4◦ 96.7◦Height Above Glacier (h) 65 m 130 mthe effects of changing solar angles, and all replaced images were within 2 hours of the‘ideal’ time. Typical camera views are shown in Figure 4.3, along with targets trackedthroughout the study period.Six stakes were also installed on the glacier, two of which were within the viewshed ofthe Nunatak Camera. These stakes were surveyed at visits in June, July and Septemberwith a Garmin eTrex Legend Cx commercial grade GPS. The GPS is WAAS enabled andhas been tested to be accurate up to 1 m.4.3 Measuring Glacier Velocity from Time Lapse ImagesThe overall prodecure for measuring glacier velocity from oblique time lapse images isoutlined in Figure Tracking FeaturesFeatures were tracked manually using Tracker (,an Open Source program that produces X and Y pixel values for a tracked point in suc-cessive frames. While Tracker has an automatic tracking feature, it was found to beunreliable due to the relatively large changes in target lighting and shape throughout theseason. Changing shadows (due to time of day, weather, seasonality), changing featureshape (due to melt or movement), and disappearance of the target altogether (calving ormelt) made tracking a single feature through the whole time series difficult, and sometimesimpossible. In the case that a feature was lost, a nearby feature was selected for the re-mainder of the dataset (the time-step of the change was flagged so as not to produce anartificial displacement).In order to pick the best possible feature to track, a few simple guidelines were followed:• Tracked features were selected to have a substantially different colour than the back-654.3. Measuring Glacier Velocity from Time Lapse ImagesTrack on-glacier featuresCalculate pixel displacementCorrect for non-perpendicular flowCorrect for camera motionFilter erroneous valuesTrimmed mean of target displacementsCalculate pixel displacement Track Ground Control Points (GCP)Create baseline from mean GCP displacement Glacier velocity time-seriesCalculate time between framesCamera/Glacier geometryThreshold and standard error estimatesConvert to distance using pixel-distance factorFigure 4.4: Time lapse image processing workflowground ice. Rocks or crevasses worked best since the black features stuck out againstthe grey/white glacial ice in a wide range of lighting conditions.• The middle of a feature was tracked. Shadows tended to distort the edge of features,making it harder to discern shadow from feature edge.• Features were selected near the centre of the image to minimize the potential forimage distortion, which was not accounted for in this study.• Late afternoon images had the least amount of light angle changes throughout theseason, while early morning and evening pictures went from direct to diffuse sunlightas the days got shorter, while lengthening shadows also made features more difficultto discern. The lighting for midday images was found to be most affected by changes664.3. Measuring Glacier Velocity from Time Lapse Imagesin cloud cover.In order to account for camera motion due to winds, temperature changes, and otherfactors shifting the camera orientation (Ahn and Howat, 2011; Eiken and Sund, 2012; Ahnand Box, 2010), several ground control points were tracked. Ground control points wereselected as distinguishable features that would not be expected to move, such as large rocksor lichen on a rock in the foreground. The pixel displacement for each ground control pointwas checked manually, and a mean background pixel displacement was computed. Thisbaseline value was then subtracted from each tracked feature pixel displacement.4.3.2 Converting pixel displacement into distance displacementConverting a pixel displacement into distance displacement was achieved by creating apixel to distance scale (m/pixel). By knowing the location of two control points (CP inFigure 4.5), the distance between them perpendicular to the camera (x) can be computed(Harrison et al., 1992) or measured as follows:x′ = xy′y. (4.1)In this study, initial values of x and y were measured using Google Earth. Using thedistance from the line x to the camera (y), and from the tracked feature point to thecamera (y′), a scale can be created where x′ is proportional to x, which can be related tothe number of pixels in that distance to calculate a distance per pixel translation factor.Assuming that the glacier flow is perpendicular to the camera (parallel to and at linex′), we can reference points (in real world distance) with respect to x′, where 0 correspondswith the centre of the frame (the intersection of lines x′ and y′). In this convention, positivedisplacement is attributed to movement from left to right, and was corrected since glacierflow was in the opposite direction.4.3.3 Calculating the horizontal distance from cameraIn order to accurately calculate the down-glacier displacement in distance units (Equation4.1), an adequate distance per pixel ratio must be calculated. This calculation relies onan accurate measurement of the distance between the camera and the tracked point (yi).Estimations of this distance can be aided and constrained by glacier features that can beseen from georeferenced remotely sensed data (such as medial moraines, or large crevasses).674.3. Measuring Glacier Velocity from Time Lapse ImagesCx’xCP CPy’yθδvFigure 4.5: Camera geometry for converting pixel displacement into pixel values.However, errors in the distance between the camera and tracked feature will scale linearlywith measured displacement (Eiken and Sund, 2012), making an accurate calculation of yimore important.By assuming a planar glacier surface, geometric relationships can be used to calculatethe distance of any point on the glacier along the line y (Figure 4.6). By knowing theelevation of the camera above the glacier plane (h), and the angle of the point of interestfrom horizontal (α), referred to as the zenith angle, a target’s distance from the camera(yi) can be calculated as follows:yi =htanαi. (4.2)In order to solve Equation 4.2, the zenith angle (α) must be calculated for the point ofinterest. This can be accomplished by having two known points, and their distances fromthe camera (y1, y2), as well as their corresponding Y pixel values (Y1, Y2). The zenith angle(α) is first calculated for each known point as follows:α = tan−1hy. (4.3a)684.3. Measuring Glacier Velocity from Time Lapse Imagesyy y KCα2y α1Figure 4.6: Camera depth geometry, where yo is horizontal to the camera.Since the zenith angle scales linearly with Y pixel values (Y ), the zenith angle for any Ypixel location (Yi) can be calculated asm =Y2 − Y1α2 − α1(4.3b)where m is the slope of the regression line. Individual pixel points can be converted intozenith angles asαi =Yi − bm(4.3c)using the intercept of the regression line, b.Finally, zenith angles for tracked points can be converted into distances from the camerausing Equation 4.2. This relationship suggests that the degree of horizontal foreshorteningis strongly dependent on the height of the camera above the glacier plane (h). For relativelylow heights, the degree of foreshortening is great, and a small decrease in the zenith anglewill result in a large increase in distance from the camera. Since picture quality is confinedto a finite number of pixels, the accuracy of any pixel measurement in the y directiondecreases substantially with low camera heights and low zenith angles (Figure 4.7).694.3. Measuring Glacier Velocity from Time Lapse Images0 20 40 60 8001000200030004000α (°)Distance from Camera (m)10000500020001000500100 Height Above Surface (m)Figure 4.7: Distance from camera as a function of the height above glacier surface (h) and zenithangle (α).4.3.4 Correcting from non-perpendicular flowIn the likely case that flow is not perfectly perpendicular to the camera (represented asv in Figure 4.5), knowledge of the target’s X pixel location and direction of flow yieldsa corrected displacement value (d). The raw displacement value (d′) can be corrected,following Eiken and Sund (2012), whered =d′cos δ + sin δ tan θ(4.4)and θ is the angle of the target relative to the centre of the camera frame and δ the directionof glacier flow relative to perpendicular to the camera.4.3.5 Data ProcessingDisplacement values were filtered, where negative (up-glacier) displacements and valuesgreater than a specified threshold (1.5 md−1 in this study) were classified as erroneous.Velocities were computed by dividing the displacement by the individual time interval704.4. Results and Discussionbetween images. Finally, if multiple points were tracked in the same approximate location,measurements for each time step were compared. A final velocity was calculated by usinga trimmed mean. In each case, if the number of velocity measurements for a time-stepwas greater than half the total possible measurements, and the spread was larger than3 standard deviations, the value furthest from the mean was removed. This process wasrepeated a second time, and second outlier was removed if the spread was still larger than3 standard deviations.Outliers were predominantly substantially larger than the mean tracked velocity. Theseoutliers are most likely due to the fact that if a feature is not captured in the previousframe, the displacement captures 2 (or more) days, artificially overestimating the dailyvelocity. The majority of these values were removed by applying a threshold (1.5 md−1 isroughly 3 times the average flow speed). However, some outliers remained, and were mostoften visibly divergent from the majority of tracked points.4.4 Results and Discussion4.4.1 Ground Control PointsJul Aug Sep−15−10−5051015 Noise from Stationary Nunatak Control PointsPixel Displacement(a) Nunatak Camera GCP noiseJul Aug Sep−15−10−5051015 Noise from Stationary Terminus Control PointsPixel Displacement(b) Terminus Camera GCP noiseFigure 4.8: Displacements for stationary Ground Control Points for the two cameras. Note: theanomalously large displacement in late July is due to the camera being moved when the SD cardswere changed.The displacements for the ground control points show significantly larger displacementsin the nunatak dataset (Figure 4.8). The average displacement in the terminus data was inthe order of ±1-2 pixels, while the nunatak data showed an average displacement of ±4-5714.4. Results and Discussionpixels. The nunatak camera was substantially more exposed, and had a larger wind dragdue to a wood backing on the tripod, making it more susceptible to shaking. The terminuscamera, on the other hand, was sheltered by a large boulder and was located in a smalldepression.Ground control point displacements on Terminus Camera were found to have a diurnalpattern, with movements downglacier during the afternoon, and upglacier in the morning.This pattern was most strongly observed during periods of high pressure, sunny weather,coiciding with the ideal conditions that promote strong katabatic winds, suggesting thatwind was the main factor responsible for camera stability.4.4.2 Terminus VelocityTracked points from Terminus Camera showed relatively consistent velocities of roughly 0.4md−1 throughout the study period, corresponding with total displacements of just under34 m. Between 6 and 14 frames were discarded from each target due to anomalously high,or negative daily displacement values (see Table 4.2).Table 4.2: Terminus tracked points. All points are located along the terminus ice front, roughlyperpendicular to glacier flow, 972 m from the camera.Point Discarded Frames Displacement (m) Flow Speed (ma−1)a 13 34.28 140.60b 14 31.97 131.13c 13 33.56 137.62d 6 33.70 138.20e 14 32.52 133.38f 12 34.30 140.65g 9 36.40 149.27h 9 33.91 139.08Mean 11 33.83 138.74Raw velocity data were very noisy, and were smoothed using a running mean (spanof 5) in order to attempt to isolate seasonal variations. Average velocity measurementsshow little variability throughout the summer (Figure 4.9). A small increase in velocity isapparent throughout July and mid-August. In early September, terminus velocity decreasesroughly 0.1 md−1.724.4. Results and DiscussionJul Aug Sep0. (m/day)Figure 4.9: Terminus velocity time series. The grey lines represent each tracked target, while thebold line is the trimmed mean. Data are smoothed into a running mean (5 time steps).This relatively steady summer velocity profile is consistent with what would be physi-cally expected from a floating terminus. Since there is no basal friction, seasonal changesin meltwater channels, and subsequently basal water pressure, do not affect the terminus.Absent changes in basal friction (or in this case lateral friction as well, since the terminusdoes not abut valley side-walls), we would expect ‘plug flow’, where ice velocity is deter-mined primarily by the channel gradient, which did not measurably change throughout thestudy period.4.4.3 Nunatak VelocityTracked points along the glacier centreline from Nunatak Camera showed very similaraverage velocities during the study period, with substantially higher peaks, and morevariability. Similar to the terminus dataset, raw daily velocity measurements were noisy,and had to be smoothed using a running mean to isolate weekly and seasonal trends (Figure4.10). A peak is located at the beginning of the dataset (June 18-20) and then a shortperiod of low velocities at the end of June. July 1-12 is marked by well above averagevelocities, reaching running 5 day averages of over 0.85 md−1. From July 20 onwards,734.4. Results and Discussionthe velocity time series settled into a variable pattern of relatively low values, oscillatingbetween 0.3 and 0.5 md−1. Between 6 and 15 frames were discarded from each trackedpoint due to anomalously high, or negative displacement values (see Table 4.3).Jul Aug Sep0. (m/day)Figure 4.10: Nunatak centerline velocity time series. The grey lines represent each tracked target,while the bold line is the trimmed mean. Data are smoothed into a running mean (5 time steps).Tracked points along the glacier centreline from Nunatak Camera showed very similaraverage velocities during the study period, with substantially higher peaks, and morevariability. Similar to the terminus dataset, raw daily velocity measurements were noisy,and had to be smoothed, using a running mean to isolate weekly and seasonal trends(Figure 4.10). A peak is located at the beginning of the dataset (June 18-20) and then ashort period of low velocities at the end of June. July 1-12 is marked by well above averagevelocities, reaching running 5 day averages of over 0.85 md−1. From July 20 onwards,the velocity time series settled into a variable pattern of relatively low values, oscillatingbetween 0.3 and 0.5 md−1. Between 6 and 15 frames were discarded from each trackedpoint due to anomalously high, or negative displacement values (Table 4.3).The increase in velocity in late June is consistent with classical glacier flow theorywhere an increase in basal water pressure builds up early in the melt season. High basalwater pressure reduces friction at the bed, allowing for faster flow. As the summer pro-gresses, and more melt water is stored in the glacier, sub-glacial channels are scoured,744.4. Results and DiscussionTable 4.3: Nunatak centreline velocity time series. All points are from near the centerline, 897 mto 980 m from the camera, roughly perpendicular to glacier flow.Point Target Distance (m) Discarded Frames Displacement (m) Flow Speed (ma−1)a 897 10 33.61 137.85h 911 6 34.04 139.60d 920 12 32.67 133.98g 946 14 37.81 155.06b 951 13 32.26 132.31c 958 12 33.84 138.79f 976 13 30.97 127.01e 980 7 33.74 138.36Mean 942 11 33.62 137.87eventually allowing the stored basal water to drain, decreasing the basal water pressure,and subsequently, flow speeds.Since Nunatak Camera was located high enough above the glacier surface, it was pos-sible to track points distributed across the width of the glacier. Tracked points show verylow velocities very close to the north edge of the glacier (Table 4.4). Velocities increasedwith distance from the camera (southward across the glacier), except for the two furthestpoints (b and l), which were both within 150 m of the south margin of the glacier, whichshowed slightly lower flow speeds than those close to the centreline.Table 4.4: Nunatak Camera distributed points. All points are within 200 m of the centre of theimage. The glacier extends from 337 m to 1698 m from Nunatak Camera.Point Target Distance (m) Discarded Frames Displacement (m) Flow Speed (ma−1)c 422 22 13.00 53.30a 494 21 13.02 53.38e 576 10 21.88 89.75g 722 6 27.07 111.00f 1080 12 32.56 133.52h 1169 18 37.74 154.79k 1262 12 39.02 160.03d 1425 14 39.49 161.94b 1568 11 37.20 152.56l 1642 12 37.99 155.81754.4. Results and DiscussionSince Nunatak Camera was located high enough above the glacier surface, it was pos-sible to track points distributed across the width of the glacier. Tracked points show verylow velocities very close to the north edge of the glacier (Table 4.4). Velocities increasedwith distance from the camera (southward across the glacier), except for the two furthestpoints (b and l), which were both within 150 m of the south margin of the glacier, whichshowed slightly lower flow speeds than those close to the centreline. The velocity timeJul Aug Sep0. (m/day)Distance from Camera (m)422494576722108011691262142515681642Figure 4.11: Nunatak Camera velocity time series (running mean, span = 5) for points across thewidth of the glacier.series for all tracked points showed a very similar pattern to the centreline points (Figure4.11). All points show a decrease in flow velocity at the end of June, followed by a largepeak in early July. Throughout the rest of the study period, velocity changes were mini-mal, with some points further from the camera (closer to the centreline and south margin)showing small increases in flow speed in late August. More generally, flow speeds weregreater further southward, and faster points showed larger amplitude changes throughoutthe study period, while low velocity points were less variable.The observed lateral velocity pattern (Figure 4.12,4.13) is consistent with glacial flowmodels, where lateral drag reduces flow speeds close to the margins of the glacier, andthe highest flows are found along the thalweg. In this case, the highest flow speeds arefound further towards the south margin, most likely due to the fact that the glacier is764.5. Error Analysiscoming out of a southeasterly bend, although it is also possible that this is a reflection ofthe longitudinal stress coupling between the floating terminus, and grounded up-valley ice.Extremely low flow speeds observed at the north (close) margin reflect the shallow, almoststagnant state of the ice, and the eddy created by the nunatak and flow bend.460000 460500 461000 461500 462000562900056295005630000563050056310001400160018002000             Elevation          (m asl)Bridge Glacier Nunatak Target Velocity (ma-1)Easting (m)Northing (m)Bridge LakeNunatakFigure 4.12: Tracked points from the Nunatak Camera (40 m contour intervals).4.5 Error AnalysisAlthough the measurements presented for Terminus and Nunatak camera conform rela-tively well with expected results, there are several sources of error that affect both thetotal and relative value of measured feature displacement, and are discussed below.774.5. Error Analysis400 600 800 1000 1200 1400 160001020304050Distance From Camera (m)Summer Displacement (m)Figure 4.13: Lateral profile from Nunatak Camera. In reality, points can be up to 200 m fromcentre of the frame in the (x) direction, but have been lined up to estimate the lateral velocitydistribution. Glacier margins are represented by the solid black lines at 337 m and 1698 m.4.5.1 Human Measurement Error and Camera ShiftManual tracking has benefits over automated matching solutions in that it can produceresults when a statistical match is unlikely, benefitting from a human ability to discernchanging feature shapes and shadows in each frame. While this suggests a superior method,this rational selection is at least partially offset by a loss in precision. Although Trackersoftware enables very high (sub-pixel) precision in manual measurements, in practice thosemeasurements will be limited by an ability to select the ideal location. This selection errorwill be compounded by the fact that, given the finite nature of pixels, sub-pixel precisionwill be limited by having to manually interpolate within the pixel to select the centroid ofthe tracked feature.The standard deviation was calculated for each time step from each camera, compar-ing the relative displacement (in pixels) of measured stationary Ground Control Points(Figure 4.14). The mean standard deviation for Nunatak Camera was 0.52 pixels/day,while Terminus Camera was 0.41 pixels/day. These values suggest that human selectionerror is roughly half a pixel per displacement step, and represents the relative uncertaintyin correcting raw displacement values for camera shift. This uncertainty can be reduced784.5. Error Analysissomewhat by filtering data using a threshold, and removing negative displacements in orderto mitigate erroneous corrections.Jul Aug Sep0. Ground Control PointsStandard Deviation (pixel/day)=Jul Aug Sep0. Ground Control PointsStandard Deviation (pixels/day)Figure 4.14: Standard deviations for 3 tracked ground control points from each time lapse camera.4.5.2 Photographic and Feature QualityBy assuming that the 8 points tracked within 100 m from Nunatak Camera and within200 m on the floating margin of Terminus Camera have the same displacement, comparingtheir relative pixel displacements yields an estimate of the uncertainty due to human errorsas well as error induced by photographic quality and the changing shape and orientationof tracked features. The quality of time lapse data is strongly dependent on the qualityof images and camera used, both in absolute pixel density, and in sensor quality. Thecameras used in this study are 12.1 megapixels, and are moderate quality and relativelyinexpensive. Due to limited image storage capacity, there is a tradeoff between imagequality and shooting frequency that must be balanced. In this study, standard deviationsfor both cameras are between 0.44 pixels/day (Terminus) and 0.46 pixels/day (Nunatak),suggesting this quality of image produces errors of roughly 0.12 md−1 at 1 km from thecamera (Figure 4.15).4.5.3 The Geometry of Control PointsAlthough the relative uncertainty of individual velocity time steps is relatively minor, overthe season, error is introduced when attempting to derive absolute distances from timelapse data. In order to produce accurate absolute results, the distance to pixel ratio mustbe computed, which is itself dependent on the geometry of the tracked point relative to the794.5. Error AnalysisJul Aug Sep0. Velocity Error AnalysisStandard Deviation (pixel/day)Jul Aug Sep0. Velocity Error AnalysisStandard Deviation (m/day)Figure 4.15: Standard deviations for 8 tracked points from each time lapse camera. Points removedfrom the outlier mean calculated average are also omitted from this calculation. Distance estimatesare derived from the mean distances of centreline tracked While the relationships for computing this geometry are elementary, substantialdifficulties remain in selecting image control points and accurately surveying them.In order to derive an initial distance to pixel ratio, two control points must be surveyed,and the distance between them calculated perpendicular to the camera (see Figure 4.5).In this study, Google Earth was used to calculate the distance between points, and theaccuracy of the measurements are estimated to be ±5 m. This estimate is compounded bythe accuracy of the surveyed location of the camera, which is given as 1 m. Uncertainty insurveying the elevation of the camera is assumed at ±10 m, producing an additional errorof 0.7%. The total estimate from daily displacement measurements translates into an errorof 2.7% on absolute displacement.Another geometric uncertainty is in calculating the distance of a tracked feature fromthe camera. The first component in calculating this value is deriving accurate distancesalong the centre of the frame for two known control points. Since features used had tobe along the glacier plane, both the close and far edge of the glacier were selected. Theuncertainty involved in delineating the margin of the glacier and accurartely surveyingits position is estimated at ±10 m, translating to an error of 6% of the final distance topixel translation factor. This error is likely a slight underestimate in the total error in thisprocedure, due to the fact that it is difficult to quantify the effect of a non-planar glaciersurface, which will introduce further uncertainty in the calculated distance of a featurefrom the camera. The total error associated with geometric simplifications and surveyingis estimated at 9% of the absolute measured displacement.804.6. Conclusions4.5.4 Comparison with Surveyed ResultsComparison between surveyed stake displacement, and time lapse camera measured dis-placement over then entire study period suggests that the time lapse camera method pro-duces greater displacement values (Table 4.5). A comparison of surveyed stake ‘B’, whichis within 100 m of Nunatak Camera centreline measurements, suggests absolute displace-ments are within the estimated error. Comparison between surveyed stake ‘A’ and Nunatakdistributed point ‘h’, however, shows a significantly higher camera-derived displacement.Table 4.5: Comparison of GPS surveyed stakes and time lapse camera derived study period dis-placements.Distance From Camera (m) Displacement (m)Nunatak PointsNunatak Averaged Centreline 942 33.6 ±3.1Nunatak Distributed point ‘h’ 1169 37.7 ±3.5Surveyed PointsB 830 27 ±5A 1222 23 ±2In this study, the GPS values are relatively unreliable for a full evaluation of the timelapse camera method due to the high errors involved. The quality of the GPS surveyingwas substantially reduced by high valley side walls, short occupation times, and the rela-tive logistical challenges in bringing higher accuracy (D)GPS systems into the field. Themeasured stakes suggest; however, that tracked displacement values agreed with surveyedvalues within a relatively narrow range. Differences within a few metres suggest that thediscrepancies may just as likely be due to the uncertainties associated with the GPS systemas values measured from time lapse imagery manual feature tracking.4.6 ConclusionsIn this chapter, a method has been presented for producing usable glacier velocity measure-ments from a monoscopic time lapse camera using manual feature tracking. The methodhas proven reliable at daily timescales, although it has potential to provide sub-daily dis-placement measurements under ideal conditions.The estimated error for this method is 0.45 pixels per displacement measurement dueto human error, and 9% of the measured absolute displacement due to uncertainty in814.6. Conclusionsthe measurement of control points. This error, already acceptable for many glaciologicalapplications, could be substantially reduced by using high quality GPS equipment to surveyground control points and the use of high resolution, up-to-date remotely sensed imageryto estimate glacier margins. This error could also be substantially reduced with the use ofa high quality, up-to-date DEM, although this method is not discussed here (see Corripio(2004)).Between June 17 and September 14, 2013, Bridge Glacier flowed at a speed of 139 ma−1in the terminus region, and up to 162 ma−1 in the fastest section of the nunatak region.Averaged centreline speeds were 138 ma−1. Flow speeds were relatively constant alongthe floating terminus, however, flow increased dramatically in early July all along theablation area, before settling into significantly lower flow speeds from late July throughthe remainder of the study period.82Chapter 5The Magnitude and Timing ofCalving5.1 IntroductionFor glaciers that terminate in lacustrine or marine waters, calving can be an importantcomponent of mass balance, as it leads to greater ice loss than what would be possiblethrough surface ablation alone (Van der Veen, 2002). The calving rate, Uc (ma−1), can beexpressed asUc = UT −dLdt− Uwm (5.1)and is a function of the average surface velocity at the terminus, UT (ma−1), the changein glacier length with time, dLdt (ma−1), and the terminus melt below the waterline, Uwm(ma−1).Field measurements (or estimates) of the surface velocity near the terminus, and thechange in terminus position, can be used to calculate the calving rate (Benn et al., 2007b).For continental lacustrine systems, where lakes freeze in the winter, it is assumed thatcalving events only occur in the summer months, while in low elevation or marine systemscalving may be a source of ice loss throughout the year (Warren and Aniya, 1999; Skvarcaet al., 2002; Nick et al., 2010; Van der Veen, 2002).Multiple studies have emphasized the importance of water depth to the rate and mag-nitude of calving (Skvarca et al., 2002; Warren and Kirkbride, 2003; Post et al., 2011).Significant positive linear relationships have been found between water depth and calv-ing rates (Pelto and Warren, 1991; Haresign, 2004). However, this relationship does notfully explain the complexity of calving responses in different environments, or the seasonalvariation in calving rates (Van der Veen, 2002).Water depth remains an important variable in determining calving rates due to thepotential for deep water to promote terminus buoyancy. Terminus buoyancy is found to835.1. Introductionhave, in many cases, coincided with disintegration and rapid retreat of the ice terminus(Boyce et al., 2007). Other calving systems have been found to withstand terminus flotationfor several seasons, although they have become more sensitive to small perturbations interminus conditions, and more likely to experience large calving losses (Benn et al., 2007b;Post et al., 2011).Water depth is recognized as an important control on the multi-annual calving flux;however, calvng loses are ultimately determined by a series of discrete, stochastic events.Current literature has made great strides in explaining the overarching controls on annualto decadal calving rates, but current literature lacks predictive systems to capture individ-ual events. Attempts have been made to relate calving events to changes in water level(Haresign, 2004), water temperatures and subaqueous melt (Motyka et al., 2003; Rohl,2006; Robertson et al., 2012), and strain rates (Boyce et al., 2007; Benn et al., 2007a).While these studies have elucidated several important controls on calving rates in lacus-trine and marine environments, they have also served to emphasize the uniqueness andvariability within and between each system. This uniqueness also extends into the relativeimportance of calving in each system. While some glaciers have shown massive recent icelosses due to calving (Van der Veen, 1996; Warren and Kirkbride, 2003), others have beenshown to lose relatively little mass to calving processes (Boyce et al., 2007; Tsutaki et al.,2011).This chapter aims to contribute to a better understanding of the roles and importanceof calving at Bridge Glacier during the 2013 melt season, and has three main goals. First,it aims to identify the timing, nature, and relative size of individual calving events observedduring the study period. Secondly, this chapter examines whether individual calving eventscan be related to meteorological, limnological or glaciological observations on a sub-seasonaltimescale. Finally this chapter sets out to quantify the volume of ice lost through calvingduring the 2013 melt season and relate it to the volume of ice lost through surficial meltin order to assess the relative importance of calving on the glacier’s mass balance.845.2. Methods5.2 MethodsThe volume of ice discharged through calving from the glacier terminus, Qcalving (m3a−1),termed the calving flux (Figure 5.1), can be quantified asQcalving =(dATdt+ Ux)HI (5.2)where dATdt is the rate change in glacier area at the terminus (m2a−1), U is the terminus flowvelocity (ma−1), HI is the height of ice at the terminus, and x is the cross sectional areaat the terminus (m). Ice front melt is assumed negligible compared with the magnitudesof flow velocity and terminus retreat, as well as the errors associated with estimating eachvariable.8G$GWCXU UHQ W*OD FLHU PRVL WLRQxFigure 5.1: Bird’s eye view of a conceptual model of calving flux.5.2.1 Terminus Area ChangeThe rate of change in glacier area (dATdt ) during the study period was found by comparingLandsat images from June 23 and September 11, 2013. Historical terminus positions (1984- 2012) were derived from Landsat images taken from mid-September to late October ofeach year (see Section 2.2). Landsat images were acquired using the LandsatLook Viewer( Shapefiles for end-of-season terminus positions weregenerated by manual deliniation from Google Earth imagery. The change in area was thencalculated using routines developed in R, using the rgeos package.855.2. Methods5.2.2 Terminus Flow VelocityThe terminus flow velocity was measured by tracking features from two time lapse cam-eras pointed at the terminus and roughly 1 km above the floating terminus (Terminus andNunatak TLC). Points were manually tracked using Tracker video analysis and modellingtool ( Raw pixel displacement was con-verted into real-world displacement values using known camera and control point geome-tries (see Chapter 4).Eight on-glacier points are tracked from each camera throughout the study periodusing daily noon-time images. All tracked points from the Terminus TLC are located onthe floating terminus, and within 200 m of each other, while Nunatak TLC points wereall selected near the centre flowline, within 300 m. Filtering routines discarded roughly10% tracked data points due to negative displacement, or loss of target. Final terminusand nunatak velocity time series are generated by averaging the daily displacement foreach tracked point, and an averaged summer velocity is found by averaging the totaldisplacement for each tracked point during the study period.5.2.3 Bathymetry and Ice Thickness EstimateIn order to convert changes in area into volumetric ice loss, a suitable estimate of glacierthickness must be derived. The terminus of Bridge Glacier sits in a large lake, and has anotable inflection point (see Figure 5.2) where it is assumed that the terminus transitionsfrom grounded to floating. Large tabular calving events during the melt season have shownlimited mobility immediately after calving, suggesting that the ice is right at the boundarycriterion for flotation. Given these observations, the thickness of ice at the terminus canbe approximated by assuming the glacier is at the threshold of flotation. Using the heightabove buoyancy criterion (Van der Veen, 1996; Benn et al., 2007b) whereHI = Hb +ρwρiDW (5.3)and Hb is the height of ice above the waterline (m), DW is the depth of water (m), and ρis the density of ice (i, 917 kgm−3) and lake water (w, 1000 kgm−3) respectively, the icethickness can be calculated.865.2. MethodsFigure 5.2: Photograph of the terminus inflection point, indicating flotation noted by red arrows.Yellow arrow indicates large crevasse that eventually leads to calving of the large tabular icebergto the left of arrow.Assuming the terminus is floating, the minimum ice thickness can be calculated asHI = Hb(ρwρi− 1)−1 (5.4)which is preferable in many cases because it allows for supplemental in situ observations ofthe height of ice above the waterline, which can constrain our estimation of ice thickness.To acquire water depth near the terminus, a bathymetric survey was conducted using aLowrance HDS Gen2 fishfinder. The vertical range is given at 3000 ft while the horizontalGPS accuracy is estimated at 5 m. A total of 893 points were taken in an irregularpattern, while some areas in immediate proximity to the terminus and near Lake AWS,were choked full of icebergs and were unable to be surveyed. An additional 47 manualpoints were taken by linear interpolation between surveyed points to improve coverage inlow spatial coverage areas. Bathymetric points were processed using R including routinesfrom the gstat package (Pebesma, 2004). Irregular point data was interpolated onto aregular 10 m grid using inverse distance weighting (maximum 4 neighbours). The points875.2. Methodswere overlain on a shapefile of Bridge Lake generated through Google Earth and Landsatimagery.5.2.4 Limnology and Calving EventsWater level and temperature were recorded using a HOBO U20 Water Level Logger. Thepressure transducer measured 10 minute mean pressure (which is converted into a depth)and water temperature, and is accurate to 0.5 cm and 0.44◦C. The pressure transducer washoused inside a radiation shield, and attached to an anchor. The transducer was deployedJune 18, 2013 into 1 m of water, roughly 5 m offshore, in the north end of the lake nearthe terminus.Water level, WD (m), was calculated asWD =PL − Pgρw(5.5)where PL is the measured transducer pressure, P is the atmospheric pressure obtainedfrom the Lake AWS barometer (hPa), g is the gravitational constant (9.81 ms−2), and ρwis the density of water (1000 kgm−3).Water levels were found to fluctuate due to waves, and also to rise improbably oververy short time steps, most likely due to waves displacing the transducer or forcing it toslip and sink deeper into the lake. In order to remove these suspected artificial signals, thewater level data was filtered. A routine was written to separate artificial from true waterlevel fluctuations based on a 3 cm threshold. If the differences in water level between atime step (ti) and the subsequent datapoint (ti+1) were greater than the threshold, andtwo time steps past the data point of interest (ti+2) and the time step preceding the datapoint of interest (ti−1) were also greater than the threshold, the time step is flagged as an‘artificial shift’. In order to correct for this shift, the magnitude of the shift is subtractedfrom all data points after the event. Several different thresholds were tested, and smallerthresholds were found to over-prescribe ‘artificial shifts’, while larger thresholds did notcapture clearly visible shifts in the data. Given the current dataset and lake size, a 10minute mean shift of over 3 cm is highly improbable. In total, 7 time steps were flaggedas ‘artificial shifts’ (see Figure 5.3).In order to capture the timing and relative size of individual calving events, the terminuswas tracked by a time lapse camera located at Ridge AWS. The camera acquired pictures885.3. ResultsJul Aug Sep1. Level (m)ProcessedRawArtificial ShiftFigure 5.3: Processed and unprocessed water level data. Darker vertical lines indicate that two‘artificial shifts’ were detected in near proximity, suggesting the shift was a large wave affectingmeasurement, rather than a shift in the position of the transducer.every hour of the day between June 18 and July 20, and between 06:00 and 22:00 from July21 to September 14. Pictures were visually scanned to find the timing of calving eventswithin the hour.5.3 Results5.3.1 Calving EventsDuring the 2013 summer study period, the vast majority of the calving flux can be explainedby three major calving events. On June 23 at 01:00 hrs the first major calving event isobserved (Figure 5.4a and 5.4b). This event is characterized by a large tabular iceberg(roughly 400 m × 500 m) failing along a series of interconnected crevasses, removing theV-shaped protrusion from the centre of the terminus. The failure was preceded by severaldays of the slow propagation of the crevasses towards one another.The second major calving event was observed August 2 at 13:00 hrs (Figure 5.5a and5.5b). This calving event (roughly 100 m × 150 m) failed on the south side of the terminus.Similar to the previous major event, this event was preceded by up to a week of the slowpropagation of a crevasse. During that period, the prow of the terminus in that regionraised substantially above the surrounding ice, tilting backwards towards the terminusbefore ultimately failing.895.3. Results(a) 18:00 hrs, June 22, 2013.(b) 06:00 hrs, June 23, 2013.Figure 5.4: Large tabular calving event, June 23, 2013. Note large crack propagating along thefull extent of the terminus protrusion.905.3. Results(a) 12:00hrs, August 2, 2013.(b) 13:00 hrs, August 2, 2013.Figure 5.5: Moderate calving event, August 2, 2013. Note the uplifted prow of the terminus at farlooker’s left (south).915.3. Results(a) 06:00 hrs, August 21, 2013.(b) 07:00hrs, August 21, 2013.Figure 5.6: Large calving event on looker’s right sidewall, August 21, 2013. Note large crackpropagating up-glacier along terminus (looker’s right).925.3. ResultsThe third and final major calving event during the summer was observed on August 21at 07:00 hrs. This event was roughly 800 m × 100 m, and involved a failure along the fullnorth facing (longitudinal) length of the terminus, starting near the identified inflectionpoint between floating and grounded terminus (Figure 5.2). The ice appears to have failedalong a series of interconnected crevasses. In the days leading up to the event, the terminusice front rose several metres, and fell backwards as it calved. The event was characterizedby multiple smaller events that occurred within several hours of each other (Figure 5.6aand 5.6b).During the summer, there were also numerous small calving events observed wheremasses of ice less than 10 m in width failed along the ice-front. In most cases, this calvingevent exhibited different behaviour than the high magnitude, low frequency nature of themajor calving events. These events appear to calve almost exclusively above the waterline,failing along a crevasse. Rough estimates of the size of calving events yields a volume of0.032 km3 assuming an average 110 m ice thickness.5.3.2 Calving CorrelationsThroughout the summer of 2013, there does not appear to be any clear relationship be-tween major calving events, and pertinent glaciological, limnological or meteorologicalvariables. Throughout the study period, the melt rate obtained from both a distributedenergy balance model (DEBM) and temperature index model (TIM) is relatively consis-tent, oscillating between 5 and 10 cm/day (Figure 5.7a.). Melt peaks in early July, andagain in early August. Several large storm events depress the melt rate for several days inlate June, mid-August and early September. No large calving event is preceded by higherthan normal melt days. In fact, the highest summer melt rates coincide with the longestnon-calving period (July) during the study period.The lake level has both diurnal and seasonal oscillations (Figure 5.7b.). Daily variationsshow a peak in lake level in the mid-to-late afternoon, roughly 17:00 hrs, with minimumlevels occurring near sunrise (06:00 hrs). Seasonally, lake level rises through June andearly July, before dropping almost to pre-summer levels in mid-July. The lake level risesagain in August, before another drop in late August. Rises in lake level follow sustainedperiods of high melt rates, while drops follow storm events, which have low melt rates (evenwhen significant rain occurs). Daily mean lake level is significantly positively correlatedwith mean daily air temperature and modelled glacier melt rates. Both the June and935.3. Results4681012dat$timeMelt Rate (cm/day)DEBMTIMa.01020304050dat$timeWater Level (cm)b.10 min MeanDaily Mean0. Temperature ( ° C)c.10 min MeanDaily MeanJul Aug Sep0. Velocity (m/day)TerminusNunatakd.Figure 5.7: Measured melt rate, water level, water temperature and running mean (3 days) ter-minus velocities for two sites at Bridge Glacier. Calving events are shown as thick grey lines.early August calving events coincide with an increase in daily mean water level; however,all three events also coincide with near-minimum diurnal (10-minute mean) water levelvariations.During the study period, water temperature is found to decrease from 1.5◦C to 0.9◦C;however, this decrease is far less than the observed daily variability (±1.5◦C). More in-945.3. Resultsdepth lake temperature observations (L. Bird, pers. comm.) suggest that the lake is wellmixed, and that temperatures are relatively uniform both spatially and with depth. It islikely that the variability observed here is slightly elevated due to the sensor’s proximity tothe shore. Daily maximum water temperatures are observed in the late afternoon (roughly16:00 hrs), while minima are observed close to sunrise (05:00 hrs) (Figure 5.7c.).Overall, there does not appear to be any discernible short-term relationship betweencalving events and water temperature, except for the final (late-August) calving event,where water temperature drops substantially following the event. This event is closest tothe location of the water temperature probe, making it likely that newly discharged ice fromthe calving event may have flowed quite close to the sensor, depressing water temperaturesin the vicinity. Seasonally, water temperature decrease coincides with higher melt ratesand more calving. Both of these processes allow for water near the freezing point to beadded to the lake. Minimal stratification and winds blowing icebergs several 100 m wouldenhance mixing and further depress temperatures.Finally, there does not appear to be any relationship between calving events and trackeddaily glacier flow velocities at either the terminus or up-glacier near Nunatak TLC (Figure5.7d.). Terminus velocity remained relatively constant throughout the study period, oscil-lating between 0.35 and 0.50 m/day. There is no discernible seasonal trend, and observedvelocities do not appear to respond to either changing melt rates or water level. However,nunatak measured velocities show a doubling of flow speed in early July, following a largespike in melt rate. The elevated flow speeds remain for over a week, before settling to a rel-atively consistent velocity, slightly lower on average than those observed along the floatingterminus. Observed velocities show a slight decrease immediately following calving events;although this trend is minor, and the change in velocity is well within the measurementerror for the method.5.3.3 Calving FluxThe change in terminus area in 2013 is 2.97×10−1 km2 over the 89 day period. The averageterminus velocity during the study period was measured as 139 ma−1, while the terminuscross-sectional distance was 1055 m, giving an additional area change of 3.42×10−2 km2.Due to the paucity of water depth measurements close to the terminus, water depthwas estimated from a cross-section adjacent to the June 2013 terminus position (see Figure5.8). The median depth was 109 m, corresponding to a height above buoyancy of 9.9 m,955.3. Results0 50 100 150 186Water Depth (m) 20  20  20  40  60  80  100  120  140  140  160 Bridge Glacier461000 462000 463000 464000 465000 466000Easting (m)562900056300005631000563200056330005634000Northing (m)Figure 5.8: Bathymetry for Bridge Lake (20 m contour intervals).and an estimated ice thickness of 109 m. Combining these terms in Equation 5.2 yields anestimated calving flux of 3.62×10−2 km3 for the 85 day study period.Comparing the volume of mass lost through calving with the volume of surface icemelt during the same period (Chapter 3) yields a total mass loss of 1.60×10−1 km3 of ice(Figure 5.9). For the study period, calving accounts for 23% of the mass loss, equivalentto an additional 1.3 m of surficial ice melt in the ablation area.965.4. DiscussionNet Ice LossSurface MeltCalving0.00 0.02 0.04 0.06 0.08 0.10 0.12 0.14 0.16 0.18 0.20Volumetric Ice Loss (km3)Figure 5.9: Total volumetric mass loss at Bridge Glacier during the study period, June 20 toSeptember 12, 2013.5.3.4 Uncertainty AnalysisThe rate of change in terminus position is measured as a product of Landsat 7 resolution,which is given as 30 m. The error associated with measuring features from Landsat imageryis estimated as 2 pixels, which gives a total error of 7200 m2a−1 or 2.4%. The cross-sectionalarea, which is also estimated from Landsat imagery is calculated at 5.7%, while the terminusdisplacement is calculated in Chapter 4 as 9%. Finally, the uncertainty in ice thickness iscalculated as a product of water depth. Errors in water depth are calculated by summingthe absolute value of the residuals of a loess smoothing algorith used to estimate waterdepth at the June 2013 terminus. The absolute error in water depth is estimated ±6.7 mat the terminus, giving an ice thickness uncertainty of 5.6%. The total uncertainty in thecalving flux calculation is estimated to be 9.3%.5.4 Discussion5.4.1 Potential Calving MechanismsCalving events observed during the 2013 melt season share several commonalities in char-acter and magnitude, and are consistent with observations from other lacustrine glacierswith floating termini (Dykes et al., 2011; Boyce et al., 2007). First, the observed calvingevents are relatively long processes, taking upwards of several days for crevasses to fullypropagate to allow the ice to be discharged. Secondly, in the days and weeks leading upto the calving event, the ice margin tilted backwards against the terminus, exposing small975.4. Discussionmeltwater notches as the ice separates from the terminus (see Figure 5.10). Both of theseobservations suggest that buoyancy induced torque is the primary mechanism responsiblefor calving.Figure 5.10: Ice front tilting backwards against the terminus prior to calving, exposing smallmeltwater notches, June 17, 2013, 5 days prior to calving event.A backwards tilt of the ice margin, most clearly observed in the August 2 and August21 events (see Figures 5.10, 5.11), suggest the presence of an ice-foot, similar to what hasbeen observed at Tasman Glacier (Rohl, 2006; Dykes et al., 2011). Cold water near theterminus likely depresses melt rates below the waterline, save for a shallow near-surfacelayer that experiences daily warming. Given minimal expected subaqueous melt, it is likelythat the ice front extends further into the lake below the waterline. This greater volumeof ice below the waterline destabilizes the terminus due to an upwards buoyant force. Thisprocess is in contrast to warmer tidewater settings (Rignot et al., 2010; Motyka et al.,2003), where subaqueous melt is a significant source of ice loss, and can lead to terminusundercutting.The observed tilt of the ice margin can also be explained by basal crevassing, createdas the glacier terminus transitions from grounded to floating (Benn et al., 2007a). A clearinflection point (see Figure 5.2) along the terminus delineates the point at which the glacierbegins to float in the lake, allowing for substantial torque and possibly basal crevassing.This mechanism also explains why all major calving events are observed in the buoyantregion of the terminus.985.4. DiscussionFigure 5.11: Photograph of the ice front tilting backwards as the prow rises, July 20, 2013, roughly2 weeks before this section of the terminus calved.Ultimately, timing and the position of individual calving events appears to be moststrongly controlled by the location, orientation, and size of crevasses in the buoyant ter-minus. One commonality with all three major calving events in 2013 was the slow calvingof major tabular icebergs. This suggests that the gradual expansion, propagation, andjoining of individual crevasses, due to shear and surface melt, form a ‘preferential line ofweakness’ which allows a margin for calving (Benn et al., 2007a). This suggests that thenature of buoyant calving at Bridge Glacier is progressive and gradual, rather than a rapid,catastrophic failure brought about by abrupt perturbations in conditions.5.4.2 Calving ForcingsAn examination of daily and weekly (sub-seasonal) meteorological, glaciological and limno-logical variables yields no significant relationship with the timing or magnitude of calvingevents. This finding reflects the significant control of geometry and crevassing on the lo-cation and extent of calving. The relative insensitivity of calving events to changes inwater level contrasts with observations in Patagonia and Iceland (Haresign, 2004), but isconsistent with observations from Mendenhall Glacier, AK (Boyce et al., 2007). Whilechanges in lake level should affect the balance of hydrostatic forces against the terminusand thus promote buoyant torque, it appears that these factors are muted when the ter-995.4. Discussionminus has already achieved flotation. A combination of the relatively low magnitude ofdiurnal and seasonal changes (< 5 cm, < 50 cm, respectively), combined with a floatingterminus, creates conditions where lake level changes have a relatively small influence onthe terminus force imbalance. This influence is likely further muted because the glacier isnot close enough to being grounded to allow changes in lake level to affect the buoyancyof the terminus.Figure 5.12: Photograph of terminus melt notches, July 19, 2014, which are exposed as theterminus flexes upwards due to buoyancy.The relatively low water temperatures observed during the summer suggest that un-dercutting due to melt below the waterline is not a significant factor affecting calving atBridge Glacier. Meltwater notches are observed in close proximity to the waterline; how-ever, these notches are in the order of several centimetres (see Figure 5.12), and do notappear to measurably affect the terminus force imbalance. This hypothesis correspondswith very few observations of oversteepened ice-front failure, indicative of high subaqueousmelt rates. Conversely, cool water promotes the growth of an ice foot along the terminus,due to lower subaqueous melt, relative to that above the waterline. The presence of anice foot should serve to destabilize the terminus, and is consistent with observations ofpre-calving terminus tilting (Figure 5.11). While the low lake temperatures should desta-bilize the terminus, and alter the character of calving events, its overall effect on the rate1005.4. Discussionof calving is muted since ice foot development is a slow process and also relies on velocityand surface melt to create the necessary force imbalance to promote failure.Changes in glacier flow speed near the terminus did not concide with observed ma-jor calving events. The observed differences in variability and average magnitude of thetwo velocity time series suggest that there is substantial strain between the floating andgrounded terminal zone, and is supported by the presence of heavy crevassing as the glacierspills over a final rock bench into the lake. However, the majority of the glacier terminus incontact with the lake is not significantly crevassed, suggesting that strain is not significantat the surface in the floating portion of the terminus.In general, crevassing patterns are complicated, although it is expected that strainshould be highest at the transition from the grounded to floating ice margin, due to thereduction of basal drag. This signal is delayed, however, by the need for crevasses topropagate and connect before calving can occur. Furthermore, crevassing is also likelydue to buoyant torque, something that is not adequately captured from surface velocitygradients. While glacier velocity, and more generally strain in the terminus region, willhave an impact on the frequency, location and magnitude of individual calving events, theprocess has a much longer time scale than can be adequately captured by discrete changesin velocity patterns.Finally, while the intensity and cumulative effects of surface melt certainly play aphysical role in calving dynamics at the terminus, they do not show any relationshipto individual calving events. Surface melt should reduce the volume of ice above thewaterline, in turn increasing buoyant torque and the potential for calving induced failure,however, this forcing is cumulative rather than intensity based. Furthermore, melt is alsodirectly related to other variables affecting forces acting on the terminus, such as lakelevel, and indirectly related with others, such as grounded ice velocity (due to the factthat more basal water supply reduces basal drag). Both of these factors complicate anymeteorological control on major calving events.While calving in the floating margins of Bridge Glacier is the product of multipleinterrelated sub-seasonal variables, calving events are ultimately triggered by changes thatare cumulative, rather than stochastic. Furthermore, the interrelation of relevant calvingvariables further complicates matters, and makes isolating individual factors and assigningcausality difficult. Ultimately, although there is a physical basis for certain meteorologicalvariables triggering individual calving events, these factors, at least at Bridge Glacier, actwithin the imposed bounds of the geometry and bathymetry of the terminus region. These1015.4. Discussionfindings suggest that calving in floating lacustrine glacier termini does not significantlyrespond to sub-seasonal variability, and instead is the product of cumulative, multi-seasonalchanges in crevassing, ice thickness, and buoyancy.5.4.3 The Relative Importance of CalvingWhile calving still contributes significant ice loss that would not otherwise be possible,it is a comparatively small contributor of mass loss relative to climatic melt at BridgeGlacier during the 2013 melt season. Calving losses in this system are controlled by severalglaciological, and topographical controls that ultimately limit the magnitude of the calvingflux. While both high frequency, low magnitude and high magnitude, low frequency eventswere observed, rough calculations of the size of the three main calving events (Section5.3.1), and a detailed calculation of the total summer calving flux (Section 5.3.3), suggestthat roughly 90% of the calving losses can be explained by the three major events.The calving flux is limited due to the fact that it only occurs over a very small percentageof the total glacier area. The cross-sectional area at the calving margin is just over 1 kmwide, which strongly limits the volume of ice that can reach the floating terminus, in turnlimiting the size of calving events. Conversely, the end-of-season ablation area was over27.6 km2 in 2013. The large contributing area for melt, relative to calving, allows for meltprocesses to contribute a substantially larger volume of ice loss than from calving processes.The calving flux is also limited by relatively modest glacier flow speeds at the terminus.Velocity is limited at Bridge Glacier by a gentle gradient in the lower reaches of the glacier,and relatively narrow valley side-walls, both of which inhibit the speed at which ice canbe delivered to deep lake waters at the terminus where it becomes buoyant and calves.Near-terminus flow speeds are 1 to 2 orders of magnitude smaller than observed at largetidewater calving glacier systems in Patagonia and Greenland (Rivera et al., 2012; Koppeset al., 2011; Meier and Post, 1987).The calving flux is also limited by the bathymetry of Bridge Lake. While the terminusof Bridge Glacier is, at least partially, buoyant, it is dependent on the ice thickness ofthe terminus region. Any significant thickening, which would theoretically increase thepotential volume of ice lost due to calving, would also serve to reduce buoyancy of theterminus. If the terminus were grounded, it would serve to stabilize, and significantlyreduce calving losses. Therefore, any major increase in terminus thickness is more likelyto reduce, rather than enhance, the calving flux.1025.5. ConclusionsThese limitations on the magnitude of calving flux suggest that climate is the maindriver of ice loss at Bridge Glacier. Although calving supplies a significant portion of iceloss at Bridge Glacier, three times as much ice is removed due to surface ablation. Thishas significant implications for the future health of Bridge Glacier as it retreats further up-valley, its terminus eventually grounds, and it no longer experiences major calving events.While it would be tempting to suggest that the glacier will begin to gain mass again, thesefindings show that calving losses, at least for the 2013 melt season, are within the rangeof expected variability in seasonal melt. While the future health of Bridge Glacier maydepend in part on the changing magnitude and frequency of calving, ultimately it is moredependent on future climatic conditions.5.5 ConclusionsDuring the 2013 melt season, Bridge Glacier experienced three major calving events, whichwere responsible for the vast majority of the seasonal calving flux. The events, in late June,early August, and late August, are characterized by the discharge of large tabular icebergs.The events were preceded by up to several weeks of crevasses propagating and joining,ultimately causing to the ice front to tilt backwards against the terminus and calve. Thecharacter and magnitude of the three major calving events suggest that buoyancy is theprimary driver of calving, while crevassing provides the margin where ice calves.A high temporal resolution examination of meteorological, limnological and glaciologicalvariables related to calving found no significant relationships with individual calving events.The surface melt rate, lake water level and temperature and terminus velocity ultimatelyhave an effect on the nature and timing of calving events; however, their contributionsare cumulative rather than intensity-based. The current environment, characterized by afloating terminus in cold, well mixed water, dampens any influence of individual variables.Calving events are more strongly controlled by seasonal and annual changes, rather thanrapid perturbations in conditions.During the 2013 melt season, Bridge Glacier lost 3.62×10−2 km3 of ice through calving,accounting for 23% of total ice loss between June 20 and September 12, 2013. While calvingis an important component of mass loss at Bridge Glacier, it contributes much less massloss than surface melt in the much larger ablation area. Calving is limited by a modestterminus flow velocity, a relatively narrow cross-sectional area, and lake depth, all of whichlimit the amount of ice that can be supplied to the floating terminus.1035.5. ConclusionsThis research suggests that calving is a significant source of ice loss at Bridge Glacier;however, the glacier is still substantially more influenced by climatic melt. While annualvariations in calving can affect the mass balance of Bridge Glacier, the glacier is moresensitive to changes in the melt rate. As Bridge Glacier retreats to shallower waters (andeventually dry land), this research suggests that although future calving losses may bereduced, and ultimately cease completely, the mid to long-term mass budget of BridgeGlacier is at the mercy of the much larger influence of a changing climate.104Chapter 6Calving, Retreat, and SurfaceMelt: Present and HistoricalConditions6.1 The Relationship Between Calving and ClimateIn Chapter 2, the retreat of Bridge Glacier between 1971 and 2012 was reconstructedusing satellite imagery. The climatic influence on this retreat was examined using proxyrecords from nearby stations, and a first order estimate of the impact of calving was madeby modelling the expected retreat, had Bridge Glacier not been terminating in a lake.This study suggests that climatic warming was responsible for roughly 700 m (20%) ofthe observed 3.55 km retreat between 1971 and 2012; however, it does not consider howchanges in retreat relate to volumetric changes.Observed retreat was found to correspond quite well with the modelled climatic responseuntil 1991. This transition coincided with the beginning of observations of tabular icebergsin Bridge Lake, suggesting the terminus had achieved flotation. From 1991 to present, theretreat rate diverges from the expected climatic response, suggesting that once a glacierachieves flotation, its retreat rate is driven primarily by calving rates. This finding is inagreement with studies in other lacustrine and marine environments (Post et al., 2011;Boyce et al., 2007; Warren and Kirkbride, 2003), which shows that once a glacier achievesflotation, its climatic signal is modified by the overarching control of calving dynamics.That chapter also suggests that while the response of calving glaciers that achieveflotation can become desensitized to climatic forcing, the initial conditions to achieve thisstate are controlled by long-term climatic trends. In order to achieve flotation, the glaciermust thin substantially in order to become buoyant. Given no likely increase in flow speedsand extensional flow to drive dynamic thinning before 1991, the most probable source of1056.2. Conditions at Bridge Glacier - 2013thinning is due to a run of years with strongly negative mass balances.Climate is responsible for the initial thinning required for the terminus to achieveflotation; however, once the terminus is buoyant, glacier velocity is expected to increase,due to a decrease in basal drag. This increase in velocity is enhanced by calving eventsat the terminus, creating a mass deficit near the terminus, and ‘drawing down’ more icefrom the upper reaches of the glacier. In turn, the ‘draw-down’ promotes further thinningof the terminus, and further calving. This ‘runaway effect’ emphasizes the dominant roleof buoyancy on terminus stability, and the degree to which buoyant calving glaciers canretreat more or less independently of climate.6.2 Conditions at Bridge Glacier - 20136.2.1 AblationIn Chapter 3, ice melt for the 2013 study period is reconstructed through a combination ofin situ observations and modelling. Melt was reconstructed to be up to 5.9 m (w.e.) of iceloss in the lowest reaches of the glacier, translating to a volumetric ice loss of 0.124 km3.Melt is shown to be strongly dependent on solar radiation, particularly over low albedo ice.The spatial distribution of temperature, humidity and wind speed is found to be stronglycontrolled by high magnitude katabatic winds that persist through 95% of the summerstudy period.The geometry of Bridge Glacier plays a substantial role in both the magnitude of iceloss, and the spatial distribution of meteorological parameters. A large, low-elevation basin,combined with deep valley walls in the lowest reaches promote strong, persistent katabaticwinds. The winds are funnelled by glacier side-walls, originating in a large up-glacier basin.Flow path lengths are an effective predictor with which to model the spatial variation ofkatabatic winds, and thus improve turbulent heat flux estimates in the ablation area.6.2.2 Glacier FlowIn Chapter 4, a method is presented to measure glacier velocity using time lapse cam-eras and manual feature tracking. The method is computationally simple and relativelycost-effective, relying only on consumer grade hardware. The method produces velocitymeasurements that have an uncertainty of 0.45 pixels/day, corresponding to 0.12 md−1at 1 km from the target. This uncertainty allows for a reliable velocity time series at a1066.2. Conditions at Bridge Glacier - 2013temporal resolution of one measurement per day, with the potential for higher frequencymeasurements for closer targets or faster moving glaciers.Bridge Glacier flow speed was measured at two sites: along the north flank of thefloating terminus, and 1 km up-glacier, along a horizontal transect. Velocities were foundto average 139 ma−1 over the study period, while faster flow was observed closer to thecentreline of the glacier, reaching average speeds of over 168 ma−1. The time series showsa doubling of flow speed in early July, with highest magnitude changes closest to theglacier centerline, coinciding with high melt rates. This observation follows classical glacierflow theory (Nye, 1952), which suggests that the drainage of basal water, collected fromearly summer melt, temporarily reduces basal drag and allows for faster flow. Up-glaciervelocities became more consistent later in the summer, and oscillated around 0.4 m/day.Measured velocities in the terminus region were found to be, on average, 138 ma−1.While the net displacement is found to be similar between both grounded and floatingsites, the timing of those flow speeds varied. In contrast to the up-glacier measurements,velocities along the terminus were relatively consistent, reflecting a lack of drag. In allperiods except early July, flow speeds were roughly 0.1 m/day greater at the terminus. Thedifference in flow speed is enough to promote the strain necessary to instigate crevassing inthe lower reaches of the glacier. These surface crevasses were found primarily along a benchbetween the two camera sites, where the glacier spills into the lake. The relatively smallvelocity gradient between the two sites suggests that there is only minor spatial variabilitybetween in the lowest reaches of the glacier.6.2.3 Calving at Bridge GlacierIn Chapter 5, the calving flux is reconstructed for the 2013 melt season. The volume ofice loss due to calving is found to be 0.036 km3, the majority of this loss occurring inthree major calving events (June 20, August 2, and August 22, 2013). All three eventsfailed along a series of interconnected crevasses, and occurred in the buoyant region of theterminus.Observed Calving MechanismsOver the study period, calving was characterized by two distinct types of events. Lowmagnitude, moderate-high frequency events occurred throughout the study period, failingalong the ice front. This event was found to involve a failure along fractures several metres1076.2. Conditions at Bridge Glacier - 2013from the ice front, which was most likely due to a force imbalance along the terminal icecliff (see Figure 6.1). Slightly warmer surface water undercut the ice front due to thermalerosion, and ultimately created an overhanging ice cliff that failed along a point of weakness,most often a near-terminus crevasse.Figure 6.1: Photograph of a terminus force imbalance. The potential locations of futurefailure are suggested by red arrows. For scale, the ice cliff is approximately 10 m tall.The second type of calving observed at Bridge Glacier in 2013 was characterized bya high magnitude, low frequency event that occurred along a series of interconnectedcrevasses. The size of icebergs discharged from this mechanism, relative to the calcu-lated calving flux (roughly 90%, see Section 5.3.1 and 5.3.3), suggest that this mechanismwas responsible for the majority of the calving flux during the 2013 study period. Allthree events had similar character: slow (several days) processes where the ice cliff failedalong a series of crevasses. In all three cases, the terminal ice cliff tilted backwards in thedays and weeks preceding the event, exposing thermal notches below the waterline. Theseobservations suggest that buoyant torque was the primary force initiating calving events.It is likely that in order for buoyant torque to trigger calving events, crevasses arerequired along the base of the glacier as well as at the surface, similar to what has beenproposed elsewhere Benn et al. (2007a). During the glacier’s transition from grounded tofloating, torque at the base of the glacier is most likely great enough to create crevassesalong the inflection point. As the glacier advances further into the deep lake waters, basalcrevasses can eventually connect with surface crevasses, providing a preferential line ofweakness on which calving can occur.1086.3. Relative Importance and Long Term ImplicationsPrimary Controls on CalvingIn Chapter 5, the timing of the three major calving events is compared with meteorological,limnological and glaciological variables. Over the study period, short, rapid perturbationsin water level and temperature, surface melt and glacier flow speed have only a minorassociation with individual calving events. Although there is a physical basis for thresholdconditions to promote calving, those conditions are achieved cumulatively rather thanthrough daily (or sub-daily) perturbations.The cumulative conditions required to trigger calving events are likely dampened byterminus flotation. Since the glacier is not grounded, changes in water level and ice thick-ness, within the range observed during the study period, do not appear to alter the forcebalance at the terminus. These findings suggest that buoyancy is the main control on ter-minus stability in lacustrine environments. In environments where the terminus is buoyant,calving events, and ultimately the flux of ice discharged, is determined by the cumulativeeffects of glacier flow, water depth, and surface melt creating the conditions required toprovide the necessary buoyant torque for individual calving events.6.3 Relative Importance and Long Term Implications6.3.1 Relative Importance of Calving at Bridge Glacier 2013By combining surface melt calculations from Chapter 3 with the calving flux calculationsfrom Chapter 5, the summer balance (bs) can be determined. For the 2013 study period,calving is responsible for just under 23% of the ice loss from Bridge Glacier. Calvinglosses are limited by a modest flow speed, small cross-sectional area at the terminus, andthe bathymetry of the lake. These factors serve to limit both the rate at which ice can besupplied to deep water, and the thickness of ice that will remain buoyant. Absent significantinter-annual changes in the flow speed or ice thickness, there is a relatively narrow rangeof potential calving losses. These constraints, coupled with a large ablation area, ensurethat surface melt is a consistently greater supplier of ice loss, and more substantial impacton the overall mass balance of Bridge Glacier.1096.3. Relative Importance and Long Term Implications6.3.2 Long-Term Mass Loss FluctuationsIn order to constrain the representativeness of the 2013 melt season, and better understandthe long-term variability of calving at Bridge Glacier, estimates of historical calving fluxesand surface melt are derived. Historical surface melt is estimated using ELA observationsand a fitted linear mass balance gradient (Shea et al., 2013). Below the snowline, the netbalance (where bn = bs) is estimated using the 2013 glacier hypsometry, wherebn(z) = b1(ELA− z) (6.1)and is calculated for the elevation of every point, z (m a.s.l.), below the ELA.The coefficient value (b1 = 6.62 m(w.e.)/m) taken from Shea et al. (2013) underes-timates the volume of ice loss during the 2013 melt season calculated using distributedenergy balance modelling (Chapter 3). The coeffiecient b1 is tuned using the mass balancegradient estimated from the DEBM (b1 = 9.07 m(w.e.)/m, Figure 6.2), and is used for allyears.The glacier area is determined from the end-of-season calving margin. All points thathave calved off are given an elevation of 1400 m (a.s.l.) and are considered in Equation6.1. Historical ELAs are measured from end-of-summer (between mid-September to mid-October) Landsat images from 1984 to 2013.Errors in ELA-derived mass balance calculations are estimated by assuming a 75 muncertainty in measuring the ELA, due to timing of available Landsat images, or 22%according to Shea et al. (2013), whichever is greater. The large error in ELA estimate is toaccount for errors that cannot be adequately quantified without additional historical data,such as the linearity of the mass balance gradient.Historical terminus positions, lake bathymetry, and estimates of ice thickness are madeusing field data from the 2013 field campaign, and methodology discussed in Section 5.2.Historical velocities were assumed to be equal to the 2013 summer flow speed (140 ma−1),and calving calculations are run with 70 ma−1 (50%) potential annual variability. A 60m uncertainty in measuring the terminus cross-section (x) (equal to 2 Landsat pixels) isapplied. The rate of change in terminus area (dAdt ) uncertainty is estimated as 7200 m2a−1(2 × 60 m × 60 m). The ice thickness uncertainty is estimated as 5.6% (from Section 5.3.4)plus an additional 10 m to account for changes in sedimentation, and ice thickness relativeto water depth. Before 1991, the terminus was not floating; therefore, an ice thicknessuncertainty of 60 m is estimated to account for a range of grounded terminus geometries.1106.3. Relative Importance and Long Term Implications1400 1600 1800 2000 2200 2400-8-6-4-202values(dem)values(bn_test)Elevation (m)Summer Balance (m w.e.)b1 = 6.62, Shea (2013)b1 = 9.07, Tuned ModelDistributed Energy Balance ModelFigure 6.2: Modelled mass balance gradients from Shea et al. (2013) and a tuned coefficient usingdistributed energy balance modelling from the 2013 melt season.Between 1991 and 2004, bathymetry has poor data coverage, and a water depth uncertaintyof 30 m is estimated (giving an additional ice thickness uncertainty of 33 m).During the study period, the volume of ice lost from surface melt is variable and showsonly a minor decrease over time. Between 1984 and 2013, the ELA varied from 1926 mto 2202 m; however, most years the ELA fell between 2050 m and 2150 m, leading toa volumetric ice loss standard deviation of 0.018 km3. The 2013 study season is aboveaverage, but within one standard deviation of the mean (x¯ = 0.107 km3 ) for the 30 yearperiod. Surface melt shows a minor decrease over time, which can be attributed to the lossof area in the lowest reaches of the glacier (z =1400 m) due to calving.Calving loses are minimal before 1991, due to the relative stability of the groundedterminus. From 1992 to 1994, the calving flux increases to 0.020 - 0.029 km3 (19 - 27% ofthe total annual ice loss) before a two year period of low flux (< 0.015 km3). From 1997to 2000, calving losses increase (0.023 - 0.052 km3), before settling into another period ofrelative stability in 2001 and 2002. The highest calving fluxes occur between 2003 to 20061116.3. Relative Importance and Long Term Implications1985 1990 1995 2000 2005 20100. Ice Loss (km3 )Surface Melt Calving FluxFigure 6.3: Historical ice loss from calving and surface melt, 1984-2013. Dark vertical line corre-sponds with terminus flotation.(0.030 - 0.084 km3) and again from 2008 to 2011 (0.036 - 0.100 km3) with a period ofstability in 2006 and 2007.Calving losses are characterized by several years of high flux, and periods of relativestability. The magnitude of the calving losses increased once the glacier achieved flotationin 1991. The calving flux increased again from 2005 to 2010, and the magnitude of surfacemelt losses decreased, making the calving flux a much larger contributor of ice loss. Thevolume of ice loss due to calving is roughly equal to the volume lost through surface meltin 2005, 2008 and 2010 (44 - 49% of total volumetric ice losses).The pattern of several high magnitude calving flux years followed by several low-fluxyears suggests that calving has the potential to contribute large volumetric ice losses oversingle seasons. However, high magnitude calving losses do not persist over more than a fewconsecutive seasons, suggesting that large calving fluxes are part of a transient stage in theglacier’s retreat. It is also possible that this pattern reflects the glacier’s response to largecalving events. Large calving fluxes could potentially restore stability to the terminus,as a result of a change in geometry, and a redistribution of buoyant forces acting on the1126.3. Relative Importance and Long Term Implicationsterminus. The built up stresses, released during a calving event, appears to take severalseasons to be fully re-established.6.3.3 Bridge Glacier Relative to Other Lacustrine Calving SystemsThis study above sheds light on the main controls and drivers of mass loss at BridgeGlacier. In order to draw broader conclusions about the nature and variability of mass lossin lacustrine glacier systems, Bridge Glacier is compared to other similar systems acrossthe globe (see Table 6.1, revisited from Table 1.2). Bridge Glacier falls in the middle of acontinuum of magnitude and frequency in lacustrine calving systems. The calving rate forBridge Glacier (281 ma−1 in 2013) is larger than that for smaller New Zealand glaciers, suchas Maug, Grey and Hooker (Warren and Kirkbride, 2003), as well as Mendenhall Glacierin Alaska at the end of the 1990s (Motyka et al., 2003; Boyce et al., 2007). Conversely,calving rates at Patagonian glaciers Leon and Ameghino are up to an order of magnitudegreater (Warren and Aniya, 1999).Bridge Glacier’s calving rate is also tempered by moderate water depth and flow speeds.Higher calving rates are associated with greater water depths and significantly larger ter-minus velocities. Large Patagonian and Icelandic glaciers have terminus velocities between258 and 1810 ma−1 (Haresign, 2004), up to an order of magnitude greater than what hasbeen measured at Bridge Glacier (140 ma−1). Conversely, smaller calving glaciers in NewZealand terminate in shallow lakes (< 50 m) and many have low flow speeds (< 70 ma−1).Lake temperatures also appear to play a role in the calving rate. Many Patagonianicefields terminate in large lakes where water temperatures can be up to 7.6◦C (Warren andAniya, 1999), significantly warmer than the well-mixed 1◦C water observed at Bridge Lake.This discrepancy is most likely related to the size of the lakes. Although Bridge Glacieris relatively large (6.3 km2 in 2013), it is dwarfed by the much longer lakes of SouthernPatagonia, where depths of over 300 m and large areas limit the cooling influence of glaciermelt.Bridge Glacier shares similar calving characteristics with both Tasman and MendenhallGlacers, which have undergone significant retreat as they transitioned from grounded tofloating termini (Boyce et al., 2007; Dykes et al., 2011). During this transition, terminusvelocities increased at Tasman from 69 ma−1 to 218 ma−1 (Dykes and Brook, 2010; Dykeset al., 2011), while the calving rates for both sites increased from 50 ma−1 to between 227and 431 ma−1 (Boyce et al., 2007; Dykes et al., 2011). These rates are consistent with1136.3. Relative Importance and Long Term ImplicationsTable 6.1: Characteristics of selected major lacustrine calving glacier studies with added BridgeGlacier variables. Dw is the mean water depth, Tw is the mean water (depth averaged or range)temperature, Uc is the calving rate, and UT is the terminus averaged flow speed . Citations: a:Boyce et al. (2007), b: Motyka et al. (2003), c: Warren and Kirkbride (2003), d: Dykes et al. (2011),e: Warren and Aniya (1999), f : Stuefer et al. (2007), g: Haresign (2004), h: this study.Location Year Dw (m) Tw (◦C) Uc (ma−1) UT (ma−1) SourceAlaskaMendenhall 1997 - 2004 45 - 52 1 - 3 12 - 431 45 - 55 a, bNew ZealandMaud, Grey,Godley, Ruth,Hooker1995 4 - 20 1.7 - 4.3 14 - 88 5 - 151 cTasman 1995 10 0.5 28 11 c2000 - 2006 50 1 - 10 78 69 d2006 - 2008 153 1 - 10 227 218 dPatagoniaUpsala West,Ameghino,Grey, Nef1995 153 - 325 2 - 7 355 - 2020 370 - 1620 ePerito Mereno 1995 - 2006 175 5.5 - 7.6 510 535 e, fLeon 2001 65 4.5 - 7.0 520 - 1770 520 - 1810 gIcelandFjallsjokull 2003 75 1.5 - 3.0 582 258 gCanadaBridge 2013 109 1.1 - 1.5 281 140 h1984 - 1990 61 30 70-210 h1991 - 2003 90 82 (0 - 351) 70-210 h2004 - 2012 102 237 70-210 hBridge Glacier.At both Tasman and Mendenhall Glaciers, water depth and buoyancy, control themagnitude of calving (Boyce et al., 2007; Dykes et al., 2011; Dykes, 2013). These findingsare corroborated by observations at Bridge Glacier, which suggests that the vast majorityof the ice discharged from the terminus is triggered by buoyant forces. While the rate andmagnitude of calving is greater at Bridge Glacier than in small, grounded glaciers, it isrelatively modest in relation to large floating lacustrine glaciers of Patagonia and Iceland.1146.3. Relative Importance and Long Term Implications6.3.4 Future ImplicationsSince 1991, Bridge Glacier has been experiencing large tabular calving events that haveresulted in an irregular retreat, defined by years of large retreat (and high calving flux), andyears of relative stability. This retreat is a product of buoyant terminus conditions, whichallow for large volumes of ice to be discharged in singular events as terminus crevassesconnect and ultimately fail due to buoyant torque. However, as Bridge Glacier retreatsfurther up-valley, it will terminate in shallower waters and will eventually cease to float(Figure 6.4).PFigure 6.4: Bridge Glacier from above, September 13, 2013, showing projected future groundingline.Once the terminus retreats to shallower water and becomes grounded, large tabularcalving events will no longer occur. Instead, the calving regime will transition from highmagnitude, irregular events, to low magnitude, higher frequency events. Observationsfrom this study indicate that small, high frequency ‘terminus force imbalance’ events area relatively small source of mass loss. This transition in the character of calving suggeststhat future calving should become an increasingly minor portion of mass loss at BridgeGlacier.Even when Bridge Glacier no longer experiences large calving losses characteristic of afloating terminus, it will still be responding to a legacy calving. High rates of calving, whichremove large volumes of ice from the terminus, also serve to reduce glacier backstress, and‘draw’ ice down from up-valley, leading to dynamic thinning (Naruse and Skvarca, 2000;Benn et al., 2007a). Although not yet quantified, this thinning is apparent in historicalsatellite images (see Figure 2.1). Dynamic thinning at Bridge Glacier would serve to1156.3. Relative Importance and Long Term Implicationsenhance and/or prolong future calving due to a smaller ice thickness, reducing the depthof water required to maintain flotation. This mechanism would serve to enhance futuredynamic thinning, and also enhance surface melt induced retreat once the glacier terminusretreats to dry land, due to a thinner terminus.116Chapter 7Conclusions7.1 Major FindingsThe primary findings from this study are summarized below:• During the 2013 melt season, surface ablation accounts for 77% of the ice loss. To-pographic and glaciological factors constrain the potential variability in calving flux,suggesting that surface melt is consistently the main driver of mass loss at BridgeGlacier.• Historical calving and surface melt rates show that calving does not significantlycontribute to mass loss before 1991. From 1991 to 2003 calving was more pronounced,while the calving flux was roughly equal to the volumetric ice loss from surface meltfor 2005, 2008 and 2010. Calving contributes significant mass loss in multi-annualperiods, alternating with years of low calving output.• Surface ablation at Bridge Glacier is driven primarily by solar radiation, particularlyover low albedo ice. Turbulent flux parameters are strongly affected, both spatiallyand in magnitude, by strong, persistent katabatic flow, and have a moderate effecton the volume, and spatial distribution of ice melt.• Calving events over a single season are found to be dissociated with major meteorolog-ical, limnological and glaciological variables. Calving events are triggered primarilyby buoyant torque, suggesting that variables that affect this force balance, such aswater level and ice thickness (melt rates) affect calving more through a cumulativeeffect, rather than short timescale perturbations.• As the retreat of Bridge Glacier continues, calving is expected to become an evenless important contributor of mass loss, particularly once the terminus retreats intoshallow water, and becomes grounded. However, the effect of calving on Bridge1177.2. Future WorkGlacier’s mass balance is not limited to the calving flux. Future retreat is expectedat Bridge Glacier due to a legacy of dynamic thinning brought about by its transientcalving phase.7.2 Future WorkThis study has explored the mechanics and relative importance of calving and ablationat Bridge Glacier and highlights further questions about mass loss in lacustrine calvingglaciers.One of the major limitations of this study is that high resolution mass balance data isonly taken from one melt season, and historical surface melt rates rely on coarse estimates.This broad investigation of historical melt rates limits the scope of the conclusions thatcan be reached. In order to better constrain the relative importance of calving in a floatingterminus lacustrine system, more work is needed to refine estimates of long-term surfacemelt rates and relate them to calving losses.Further work is also needed to better understand the transition from grounded tofloating terminus. A more robust understanding, either through modelling or observation,of the drivers of calving, and how they change over these transitions will improve our abilityto explain retreat rates of lacustrine calving glaciers.An examination and quantification of the dynamic thinning at Bridge Glacier wouldalso provide valuable insight into lacustrine glacier dynamics. The timing and rates ofthinning would allow for a better separation of the role of surface melt and calving, andthe interaction between the two main drivers of ice loss at Bridge Glacier.Finally, more work is needed to compare and contrast the history of Bridge Glacier withother lacustrine systems. While this study attempts to relate observations from the studyperiod with historical rates from other lacustrine glaciers (Section 6.3.3), these studiesonly represent snapshots in the life-cycle of a calving glacier. A more robust look at thecommonalities between lacustrine glaciers will elucidate the main drivers of calving, andallow for a clearer understanding of the relationship between climate and calving in glaciermass balance.118BibliographyAhn, Y. and Box, J. E. (2010). Glacier velocities from time-lapse photos: techniquedevelopment and first results from the Extreme Ice Survey (EIS) in Greenland. Journalof Glaciology, 56(198):723–734.Ahn, Y. and Howat, I. M. (2011). Efficient automated glacier surface velocitymeasurement from repeat images using multi-image/multichip and null exclusionfeature tracking. Geoscience and Remote Sensing, IEEE Transactions on,49(8):2838–2846.Allen, S. M. and Smith, D. J. (2007). Late Holocene glacial activity of Bridge Glacier,British Columbia Coast Mountains. Canadian Journal of Earth Sciences,44(12):1753–1773.Anslow, F. S., Hostetler, S., Bidlake, W. R., and Clark, P. U. (2008). Distributed energybalance modeling of South Cascade Glacier, Washington and assessment of modeluncertainty. Journal of Geophysical Research: Earth Surface, 113(F2).Barry, R. G. (2006). The status of research on glaciers and global glacier recession: areview. Progress in Physical Geography, 30(3):285–306.Bartholomew, I., Nienow, P., Mair, D., Hubbard, A., King, M. A., and Sole, A. (2010).Seasonal evolution of subglacial drainage and acceleration in a Greenland outletglacier. Nature Geoscience, 3(6):408–411.Beedle, M. J., Menounos, B., Luckman, B. H., and Wheate, R. (2009). Annual pushmoraines as climate proxy. Geophysical Research Letters, 36(20):L20501.Benn, D., Hulton, N., and Mottram, R. (2007a). ‘Calving laws’,‘sliding laws’ and thestability of tidewater glaciers. Annals of Glaciology, 46(1):123–130.Benn, D., Warren, C., and Mottram, R. (2007b). Calving processes and the dynamics ofcalving glaciers. Earth-Science Reviews, 82(3):143–179.Benn, D. I. and Evans, D. J. (2010). Glaciers and Glaciation. Hodder Education, NewYork, NY, USA, second edition.119BibliographyBoyce, E. S., Motyka, R. J., and Truffer, M. (2007). Flotation and retreat of alake-calving terminus, Mendenhall Glacier, southeast Alaska, USA. Journal ofGlaciology, 53(181):211–224.Braithwaite, R. (2002). Glacier mass balance: the first 50 years of internationalmonitoring. Progress in Physical Geography, 26(1):76–95.Braithwaite, R. J. and Olesen, O. B. (1989). Calculation of glacier ablation from airtemperature, West Greenland. Springer Netherlands.Brock, B., Willis, I., and Sharp, M. (2000). Measurement and parameterization of albedovariations at Haut Glacier d’Arolla, Switzerland. Journal of Glaciology,46(155):675–688.Brown, C. S., Meier, M. F., and Post, A. (1982). Calving speed of Alaska tidewaterglaciers, with application to Columbia Glacier. US Government Printing Office.Christy, J. R., Spencer, R. W., and Braswell, W. D. (2000). MSU tropospherictemperatures: Dataset construction and radiosonde comparisons. Journal ofAtmospheric and Oceanic Technology, 17(9):1153–1170.Collares-Pereira, M. and Rabl, A. (1979). The average distribution of solarradiation-correlations between diffuse and hemispherical and between daily and hourlyinsolation values. Solar Energy, 22(2):155–164.Corripio, J. G. (2004). Snow surface albedo estimation using terrestrial photography.International Journal of Remote Sensing, 25(24):5705–5729.Cowie, N. M., Moore, R. D., and Hassan, M. A. (2014). Effects of glacial retreat onproglacial streams and riparian zones in the Coast and North Cascade Mountains.Earth Surface Processes and Landforms, 39(3):351–365.Danielson, B. and Sharp, M. (2013). Development and application of a time-lapsephotograph analysis method to investigate the link between tidewater glacier flowvariations and supraglacial lake drainage events. Journal of Glaciology,59(214):287–302.Demuth, M., Munro, D., and Young, G. (2006). Peyto Glacier: One Century of Science.National Hydrology Research Institute Science Report # 8.Dykes, R., Brook, M., Robertson, C., and Fuller, I. (2011). Twenty-first century calvingretreat of Tasman Glacier, Southern Alps, New Zealand. Arctic, Antarctic, and AlpineResearch, 43(1):1–10.Dykes, R. C. (2013). A multi-parameter study of iceberg calving and the retreat ofHaupapa Tasman Glacier, South Island, New Zealand. PhD thesis, Massey University.120BibliographyDykes, R. C. and Brook, M. S. (2010). Terminus recession, proglacial lake expansion and21st century calving retreat of Tasman Glacier, New Zealand. New ZealandGeographer, 66(3):203–217.Dyurgerov, M. and Meier, M. (2005). Glaciers and the changing Earth system: a 2004snapshot. Institute of Arctic and Alpine Research, University of Colorado Boulder, CO.Dyurgerov, M., Meier, M., and Armstrong, R. L. (2002). Glacier mass balance andregime: data of measurements and analysis. Institute of Arctic and Alpine Research,University of Colorado Boulder, USA.Eiken, T. and Sund, M. (2012). Photogrammetric methods applied to Svalbard glaciers:accuracies and challenges. Polar Research, 31(1).Greuell, W. and Smeets, P. (2001). Variations with elevation in the surface energy balanceon Pasterze(Austria). Journal of Geophysical Research. D. Atmospheres, 106:31.Haresign, E. (2004). Glacio-limnological interactions at lake-calving glaciers. PhD thesis,University of St Andrews.Harrison, W., Echelmeyer, K., Cosgrove, D., and Raymond, C. (1992). Thedetermination of glacier speed by time-lapse photography under unfavorableconditions. Journal of Glaciology, 38(129):257–265.Hock, R. (1998). Modelling of glacier melt and discharge. PhD thesis, ETH Institute ofGeography.Hock, R. (2003). Temperature index melt modelling in mountain areas. Journal ofHydrology, 282(1):104–115.Hock, R. (2005). Glacier melt: a review of processes and their modelling. Progress inPhysical Geography, 29(3):362–391.Hock, R. and Holmgren, B. (2005). A distributed surface energy-balance model forcomplex topography and its application to storglacia¨ren, sweden. Journal ofGlaciology, 51(172):25–36.Idso, S., Jackson, R., Ehrler, W., and Mitchell, S. (1969). A method for determination ofinfrared emittance of leaves. Ecology, pages 899–902.Jiskoot, H. and Mueller, M. S. (2012). Glacier fragmentation effects on surface energybalance and runoff: field measurements and distributed modelling. HydrologicalProcesses, 26(12):1861–1875.Joughin, I., Abdalati, W., and Fahnestock, M. (2004). Large fluctuations in speed onGreenland’s Jakobshavn Isbrae Glacier. Nature, 432(7017):608–610.121BibliographyKaser, G., Cogley, J., Dyurgerov, M., Meier, M., and Ohmura, A. (2006). Mass balanceof glaciers and ice caps: Consensus estimates for 1961–2004. Geophysical ResearchLetters, 33(19).Kirkbride, M. P. and Warren, C. R. (1999). Tasman Glacier, New Zealand: 20th-centurythinning and predicted calving retreat. Global and Planetary Change, 22(1):11–28.Koppes, M., Conway, H., Rasmussen, L., and Chernos, M. (2011). Deriving mass balanceand calving variations from reanalysis data and sparse observations, Glaciar SanRafael, northern Patagonia. The Cryosphere, 5:791–808.Krimmel, R. M. and Vaughn, B. H. (1987). Columbia Glacier, Alaska: changes in velocity1977–1986. Journal of Geophysical Research: Solid Earth (1978–2012),92(B9):8961–8968.Letre´guilly, A. (1988). Relation between the mass balance of western Canadian mountainglaciers and meteorological data. Journal of Glaciology, 34(116):11–18.MacDougall, A. H. and Flowers, G. E. (2011). Spatial and temporal transferability of adistributed energy-balance glacier melt model. Journal of Climate, 24(5).Marshall, S., White, E., Demuth, M., Bolch, T., Wheate, R., Menounos, B., Beedle, M.,and Shea, J. (2011). Glacier water resources on the eastern slopes of the CanadianRocky Mountains. Canadian Water Resources Journal, 36(2):109–134.Meier, M. and Post, A. (1987). Fast tidewater glaciers. Journal of Geophysical Research,92(B9):9051–9058.Milner, A. M. and Bailey, R. (1989). Salmonid colonization of new streams in GlacierBay National Park, Alaska. Aquaculture Research, 20(2):179–192.Mokievsky-Zubok, O., Ommanney, C., and Power, J. (1985). NHRI glacier mass balance1964-1984 (Cordillera and Arctic). Glacier Section, Surface Water Division, NationalHydrology Research Institute Internal Report.Moore, R. D. (1983). On the use of bulk aerodynamic formulae over melting snow.Nordic Hydrology, 14(4):193–206.Moore, R. D. and Demuth, M. N. (2001). Mass balance and streamflow variability atPlace Glacier, Canada, in relation to recent climate fluctuations. HydrologicalProcesses, 15(18):3473–3486.Motyka, R. J., O’Neel, S., Connor, C. L., and Echelmeyer, K. A. (2003). Twentiethcentury thinning of Mendenhall Glacier, Alaska, and its relationship to climate, lakecalving, and glacier run-off. Global and Planetary Change, 35(1–2):93 – 112.122BibliographyMunro, D. S. (1989). Surface roughness and bulk heat transfer on a glacier: comparisonwith eddy correlation. Journal of Glaciology, 35(121).Munro, D. S. (2006). Linking the weather to glacier hydrology and mass balance at PeytoGlacier. In Demuth et al. (2006), page 278.Naruse, R. and Skvarca, P. (2000). Dynamic features of thinning and retreating GlaciarUpsala, a lacustrine calving glacier in southern Patagonia. Arctic, Antarctic, andAlpine Research, 32(4):485–491.Nick, F., van der Veen, C., Vieli, A., and Benn, D. (2010). A physically based calvingmodel applied to marine outlet glaciers and implications for the glacier dynamics.Journal of Glaciology, 56(199):781–794.Nussbaumer, S., Zumbu¨hl, H., and Steiner, D. (2007). Fluctuations of the ‘Mer deGlace’(Mont Blanc area, France) AD 1500–2050. Zeitschrift fu¨r Gletscherkunde undGlazialgeologie, 40(2005&2006):5–140.Nye, J. (1952). The mechanics of glacier flow. Journal of Glaciology, 2(12):82–93.Oerlemans, J. (1994). Quantifying global warming from the retreat of glaciers. Science,264:243–245.Oerlemans, J. (2005). Extracting a climate signal from 169 glacier records. Science,308(5722):675–677.Oerlemans, J. (2007). Estimating response times of Vadret da Morteratsch, Vadret daPalu, Briksdalsbreen and Nigardsbreen from their length records. Journal ofGlaciology, 53:357–362.Oerlemans, J. (2010). The Microclimate of Valley Glaciers. Igitur, Utrecht Publishing &Archiving Services.Oerlemans, J. (2011). Minimal Glacier Models. Igitur, Utrecht Publishing & ArchivingServices, second edition.Oerlemans, J. and Reichert, B. (2000). Relating glacier mass balance to meteorologicaldata by using a seasonal sensitivity characteristic. Journal of Glaciology, 46(152):1–6.Oke, T. (1988). Boundary Layer Climates. Routledge.Østrem, G. and Brugman, M. (1991). Glacier mass-balance measurements: a manual forfield and office work. Number 4. National Hydrology Research Institute Science Report.Paterson, W. (2000). The Physics of Glaciers. Butterworth-Heinemann, third edition.123BibliographyPebesma, E. J. (2004). Multivariable geostatistics in s: the gstat package. Computers &Geosciences, 30(7):683–691.Pelto, M. S. (1987). Mass balance of southeast Alaska and northwest British Columbiaglaciers from 1976 to 1984: methods and results. Annals of Glaciology, 9:189–194.Pelto, M. S. and Warren, C. R. (1991). Relationship between tidewater glacier calvingvelocity and water depth at the calving front. Annals of Glaciology, 15:115–118.Petersen, L. and Pellicciotti, F. (2011). Spatial and temporal variability of airtemperature on a melting glacier: Atmospheric controls, extrapolation methods andtheir effect on melt modeling, Juncal Norte Glacier, Chile. Journal of GeophysicalResearch: Atmospheres, 116(D23).Post, A., O’Neel, S., Motyka, R., and Streveler, G. (2011). A complex relationshipbetween calving glaciers and climate. Eos, Transactions American Geophysical Union,92(37):305.Quinn, P., Beven, K., Chevallier, P., and Planchon, O. (1991). The prediction of hillslopeflow paths for distributed hydrological modelling using digital terrain models.Hydrological processes, 5(1):59–79.Rasmussen, L. A. and Conway, H. (2004). Climate and glacier variability in westernNorth America. Journal of Climate, 17(9):1804–1815.Rignot, E. and Kanagaratnam, P. (2006). Changes in the velocity structure of theGreenland Ice Sheet. Science, 311(5763):986–990.Rignot, E., Koppes, M., and Velicogna, I. (2010). Rapid submarine melting of the calvingfaces of West Greenland glaciers. Nature Geoscience, 3(3):187–191.Rivera, A., Corripio, J., Bravo, C., and Cisternas, S. (2012). Glaciar Jorge Montt(Chilean Patagonia) dynamics derived from photos obtained by fixed cameras andsatellite image feature tracking. Annals of Glaciology, 53(60):147–155.Robertson, C., Douglas, I., Brook, M., Fuller, I., and Holt, K. (2012). Subaqueous calvingmargin morphology at Mueller, Hooker and Tasman glaciers in Aoraki/Mount CookNational Park, New Zealand. Journal of Glaciology, 58(212):1037.Rohl, K. (2006). Thermo-erosional notch development at fresh-water-calving TasmanGlacier, New Zealand. Journal of Glaciology, 52(177):203–213.Ryder, J. M. (1991). Geomorphological processes associated with an ice-marginal lake atBridge Glacier, British Colombia. Ge´ographie physique et Quaternaire, 45(1):35–44.124BibliographyShea, J. M. (2010). Regional-scale distributed modelling of glacier meteorology and melt,southern Coast Mountains, Canada. PhD thesis, University of British Columbia.Shea, J. M. and Marshall, S. J. (2007). Atmospheric flow indices, regional climate, andglacier mass balance in the Canadian Rocky Mountains. International Journal ofClimatology, 27(2):233–247.Shea, J. M., Menounos, B., Moore, R. D., and Tennant, C. (2013). An approach to deriveregional snow lines and glacier mass change from MODIS imagery, western NorthAmerica. The Cryosphere, 7(2):667–680.Shea, J. M. and Moore, R. D. (2010). Prediction of spatially distributed regional-scalefields of air temperature and vapor pressure over mountain glaciers. Journal ofGeophysical Research, 115.Shea, J. M., Moore, R. D., and Stahl, K. (2009). Derivation of melt factors from glaciermass-balance records in western Canada. Journal of Glaciology, 55(189):123–130.Skvarca, P., De Angelis, H., Naruse, R., Warren, C., and Aniya, M. (2002). Calving ratesin fresh water: new data from southern Patagonia. Annals of Glaciology, 34(1):379–384.Stahl, K., Moore, R. D., Shea, J., Hutchinson, D., and Cannon, A. (2008). Coupledmodelling of glacier and streamflow response to future climate scenarios. WaterResources Research, 44(2):W02422.Stuefer, M., Rott, H., and Skvarca, P. (2007). Glaciar Perito Moreno, Patagonia: climatesensitivities and glacier characteristics preceding the 2003/04 and 2005/06 dammingevents. Journal of Glaciology, 53(180):3–16.Sund, M., Eiken, T., and Rolstad Denby, C. (2011). Velocity structure, front positionchanges and calving of the tidewater glacier Kronebreen, Svalbard. The CryosphereDiscussions, 5(1):41–73.Tsutaki, S., Nishimura, D., Yoshizawa, T., and Sugiyama, S. (2011). Changes in glacierdynamics under the influence of proglacial lake formation in Rhonegletscher,Switzerland. Annals of Glaciology, 52(58):31–36.Van der Veen, C. (1996). Tidewater calving. Journal of Glaciology, 42(141):375–385.Van der Veen, C. (2002). Calving glaciers. Progress in Physical Geography, 26(1):96–122.Venteris, E. R. (1999). Rapid tidewater glacier retreat: a comparison between ColumbiaGlacier, Alaska and Patagonian calving glaciers. Global and Planetary Change,22(1):131–138.125Warren, C. and Aniya, M. (1999). The calving glaciers of southern South America.Global and Planetary Change, 22(1):59–77.Warren, C. and Sugden, D. (1993). The Patagonian Icefields: a glaciological review.Arctic and Alpine Research, 25(4):316–331.Warren, C. R. and Kirkbride, M. P. (2003). Calving speed and climatic sensitivity of NewZealand lake-calving glaciers. Annals of Glaciology, 36(1):173–178.Zwally, H. J., Abdalati, W., Herring, T., Larson, K., Saba, J., and Steffen, K. (2002).Surface melt-induced acceleration of Greenland ice-sheet flow. Science,297(5579):218–222.126Appendix ADistributed Energy BalanceModelling Data ProcessingA.1 Pressure Transducer Post-ProcessingIn order to obtain a high temporal resolution melt rate dataset, a pressure transducer wasinstalled in a water-filled borehole at Glacier AWS throughout the melt season. At thetime of installation, water was free-flowing in the vicinity, and the borehole was consistentlyrefilling to the surface. However, when the borehole was re-drilled July 19, water was onlyfilling to 25 cm from the surface.Raw water depths taken from the pressure transducer show a marked diurnal trendof emptying and filling. The magnitude of this diurnal trend is variable, with spreads ofseveral centimetres, up to 80 cm (see Figure A.1). This diurnal trend does not appearto be related to either air temperatures (at Glacier AWS), or modelled melt rates. Thespread does decrease, however, as the pressure transducer approaches the surface. Thisis due most likely to the fact that there is simply less daily filling and emptying that canoccur in a shallower hole.Although there is significant observed diurnal variation in water depth, taking themaximum water depth for each day produced no physically impossible results (i.e. anegative melt rate). There was, however, significantly higher variability in the magnitudeof daily melt rates relative to modelled results from both a distributed and point EnergyBalance model (see Figure A.2). While the DEBM consistently predicts melt between 4 -12 cm/day, measured results from the pressure transducer predict melt rates for some daysclose to 20 cm/day, substantially higher than would have otherwise been expected.The variability observed in the pressure transducer melt rate suggests that it is difficultto discern how the glacial hydrology refills and empties the borehole. While in aggregate itproduces meaningful results, higher frequency (sub-daily or daily) measurements become127A.2. Solar Geometry Calculations180 200 220 240 260020406080100Julian DayDaily Spread (cm)Figure A.1: Difference between maximum and minimum depth readings for each day for thepressure transducer. The two grey lines correspond to the days when the borehole was re-drilledand melted to the surface respectively.increasingly less reliable, and seem to produce daily melt rates that are up to 10 cm differentthan what would have been predicted from energy balance modelling. Conversely, thepressure transducer method seems to be a reliable method for measuring lower frequency(weekly or bi-weekly) melt rates, allowing for additional valuable low-cost information thatcould be used to inform and constrain modelled higher resolution melt rates.A.2 Solar Geometry CalculationsAll solar geometry calculations follow classical astronomical theory, and are outlined indetail by Oke (1988). The declination for every timestep is calculated asδ =0.006918− 0.399912 cos(γ) + 0.070257 sin(γ)− 0.006758 cos(2γ)+ 0.000907 sin(2γ)− 0.002697 cos(3γ) + 0.00148 sin(3γ) (A.1a)where γ is the fractional year. The fractional year is calculated as follows:γ = 2pi(DOY − 1)/365 (A.1b)128A.2. Solar Geometry Calculations0.00 0.05 0.10 0.15 Energy Balance Model Melt Rate (m/day)Measured Melt Rate (m/day)Figure A.2: Pressure transducer observed melt rate plotted against Distributed Energy BalanceModel melt rate, extracted for Glacier AWS location.and is a function of the Julian day of year (DOY ).The hour angle (h) is calculated ash = 15(12− t)) (A.2a)where t is the fractional hour. The fractional hour is converted to local apparent time(LAT ), which is calculated asLAT = LMST−60× 229.18(0.000075 + 0.001863 cos(γ)− 0.032077 sin(γ)− 0.014615 cos(2γ)− 0.040849 sin(2γ)) (A.2b)where LMST is the local mean standard time.The elevation angle is calculated ascosZ = cosϕ cos δ cosh+ sinϕ sin δ (A.3)and is a function of latitude (ϕ), declination (δ), and the hour angle (h). The azimuth129A.3. Modelled Data Gapsangle (Ω) is calculated as follows:cos Ω = (sin δ cosϕ− cosh cos δ sinϕ)(sinZ)−1. (A.4)The mean and instantaneous Sun-Earth distance ratio (RmR ) is calculated as(RmR)2 =1.00011 + 0.034221 cos(γ) + 0.001280 sin(γ)+ 0.000719 cos(2γ) + 0.000077 sin(2γ) (A.5)and is a function of the fractional year.A.3 Modelled Data GapsDuring the summer there were a couple instances when data was not collected from oneof the AWS sensors. To have a complete and continuous data record for energy balancemodelling, the data gaps were infilled using several statistical routines and assumptionsoutlined below.A.3.1 Incoming Longwave Radiation GapsBetween June 19 and 20, as well as August 17 for 4 hours, the pyrgeometer measuringincoming longwave radiation at Lake AWS malfunctioned, most likely due to a lack ofbattery power, and did not collect data for several timesteps. In order to infill the data,the Stefan-Botzmann equation was used, which was adapted following Oke (1988).Lin = εσ(T + 273.15)4(1− 0.22n2) (A.6)was assumed, and air temperatures were taken from the Ridge AWS. n is the fraction ofcloud cover during the period. Over the 12 hours before pyrgeometer failure, the incomingshortwave radiation was approximately 50% of the potential incoming radiation, leadingto an assumed 50% cloud cover during the time period.Emissivity for the period is calculated following Idso et al. (1969) whereεa(0) = 1− 0.261[−7.704(273.15− T )2] (A.7)130A.3. Modelled Data Gapsand T is the air temperature.A.3.2 Incoming Wind Speed and Direction GapsBetween July 14 and 19, Glacier AWS was tipped over due to a suspected combination ofuneven melt, high winds, and a crevasse opening up under one of the legs of the tripod.During this period, wind speed and direction data were unable to be measured. Otherdata measured from the Glacier AWS, such as temperature and humidity, were not foundto have varied in their correlations with both Ridge and Lake AWSs, suggesting those datawere not measurably affected.In order to fill wind speed data gaps, wind speed was correlated against ambient airtemperature, taken from Ridge AWS. The highest correlation was found with a 1 hourlag time (6 time steps), showing an adjusted r2 of 0.65, which is strongly statisticallysignificant.Wind direction was modelled based on a threshold difference between ambient andon-glacier air temperatures. If the difference was greater that 5◦C, the wind directionwas classified as downslope (255◦), while if the temperature difference was smaller, winddirection was classified as upslope (90◦).131Appendix BTime Lapse Camera Methodologyand WalkthroughB.1 Getting StartedTracker is a free Open Source modelling tool built in Java using Open Source Physics(OSP). It has a nice GUI, and allows an easy way to physically track features as theymove across successive pictures. It can be downloaded on Mac OS X, Windows, and Linuxfrom I have not used Linux, but bothMac and Windows options come with a nice installer, so all you need to do is follow theirdirections. Xuggle is the preferred video engine, and is automatically installed with yourTracker install (this can be disabled). Install FAQ and further documentation available at Importing PicturesTracker recognizes both videos (.mov format) and pictures. Images will have a visiblyhigher resolution than a video in all likelihood, due to video compression routines and arerecommended. In order to import time-lapse images, they need to be named successively,without gaps. An added wrinkle in this is that the numbering must all have the same char-acter length. For instance: nunatak 8.jpg, nunatak 9.jpg, nunatak 10.jpg will notwork because 9 and 10 are not the same character length; instead, they need to be namednunatak 08.jpg, nunatak 09.jpg, nunatak 10.jpg. This naming scheme should hap-pen automatically out of the camera, but if you are selecting only certain pictures andhave ‘jumps’ in numbering, that will need to be fixed before Tracker recognizes they aresuccessive frames.In all likelihood, the volume and resolution of your images will require added memory132B.2. Importing PicturesFigure B.1: Memory in useto properly display your images. When you get close to full, the banner in the top rightcorner will blink red and show you how close you are to full (see Figure B.1). When youare full, the program will have trouble displaying frames, and will display some as blank,leaving holes in your data. You can fix this by clicking on the Memory in use button(Figure B.1) and then navigating to the Runtime tab where you can manually increase thememory size (see Figure B.2).Figure B.2: Manually increase the memory sizeChoosing pictures may also end up being one of the most critical steps to getting gooddata:133B.2. Importing Pictures• Pick pictures with similar lighting. This will make tracking features easier, and moreaccurate.• Shoot roughly three times as many pictures as you may need. Inclement weather,such as fog, low clouds or rain/water on the lens can make pictures unusable.134


Citation Scheme:


Citations by CSL (citeproc-js)

Usage Statistics



Customize your widget with the following options, then copy and paste the code below into the HTML of your page to embed this item in your website.
                            <div id="ubcOpenCollectionsWidgetDisplay">
                            <script id="ubcOpenCollectionsWidget"
                            async >
IIIF logo Our image viewer uses the IIIF 2.0 standard. To load this item in other compatible viewers, use this url:


Related Items