High-resolution three-dimensional imaging and analysis of rock falls in Yosemite Valley, California

Yosemite Valley is a ∼1-km-deep, glacially carved canyon in the Sierra Nevada mountains of California that hosts some of the largest granitic rock faces in the world (Fig. 1). These steep walls are sculpted over time by rock falls that typically occur as exfoliation-type failures along surface-parallel sheeting joints (Matthes, 1930; Huber, 1987). Thick (≥100 m) talus accumulations flanking the cliffs record substantial rock-fall activity since the last glacier retreated from Yosemite Valley ca. 17 ka (Wieczorek and Jäger, 1996; Wieczorek et al., 1999). Over 700 rock falls and other slope movements have been documented in Yosemite National Park since A.D. 1857, and some of these events have resulted in several fatalities, numerous injuries, and damage to infrastructure (Wieczorek and Snyder, 2004). Recognition that rock falls pose a significant natural hazard and risk to the 3–4 million annual visitors has prompted detailed documentation and investigation of rock-fall triggering mechanisms, causative factors, and runout dynamics to better assess geological hazard and risk in Yosemite Valley (Wieczorek and Jäger, 1996; Wieczorek and Snyder, 1999, 2004; Wieczorek et al., 1998, 1999, 2000, 2008; Guzzetti et al., 2003; Stock and Uhrhammer, 2010).

Quantifying rock-fall events is a critical component of hazard analysis because (1) particle shapes, volumes, source area locations, and cliff surface geometry influence rock-fall trajectories (e.g., Okura et al., 2000; Guzzetti et al., 2003; Wieczorek et al., 2008), (2) computer programs that simulate rock-fall runout utilize these data (e.g., Jones et al., 2000; Agliardi and Crosta, 2003; Guzzetti et al., 2002, 2003; Lan et al., 2010), and (3) accurate and precise rock-fall volumes are needed to develop reliable probabilistic hazard assessments (Dussauge-Peisser et al., 2002; Dussauge et al., 2003; Guzzetti et al., 2003; Hantz et al., 2003; Malamud et al., 2004; Brunetti et al., 2009). However, in Yosemite Valley such quantification often is difficult due to the sheer scale of the rock faces, the relative inaccessibility of these faces, and the hazards associated with fieldwork in active rock-fall areas. Rock-fall volumes have traditionally been estimated using the product of detachment surface areas and an assumed failure depth or thickness; this value can be compared to volume estimates of fresh talus beneath the source area. Both techniques result in large volumetric uncertainties, typically on the order of ±20% and occasionally much larger (Wieczorek and Snyder, 2004). In Yosemite, these techniques are further limited by the inaccessibility of the vertical to overhanging ∼1-km-tall cliffs. Volume estimates from fresh talus may be imprecise because rock masses often fragment on impact, partially disintegrating into dust that can drift far from the impact area (Wieczorek et al., 2000), and because rock falls can mobilize talus from previous events (Wieczorek and Snyder, 2004). The resulting uncertainties propagate to additional uncertainties in the accuracy of hazard assessment.

New imaging tools offer opportunities for characterizing large rock faces and reducing uncertainties in their measurement and analysis. High-resolution (gigapixel) digital photography and airborne and ground-based terrestrial laser scanning (light detection and ranging [LiDAR]) are emerging remote sensing techniques that enable precise location, measurement, monitoring, and modeling of mass movement events (e.g., McKean and Roering, 2004; Derron et al., 2005; Lim et al., 2005; Lato et al., 2009; Sturznegger and Stead, 2009; Abellán et al., 2010; Lan et al., 2010). Terrestrial laser scanning is a particularly valuable tool for quantifying rock falls, identifying failure mechanisms, assessing slope stability, and monitoring vertical cliff faces, especially when repeat laser scans permit change detection (e.g., Rosser et al., 2005, Abellán et al., 2006, 2009; Jaboyedoff et al., 2007; Collins and Sitar, 2008; Rabatel et al., 2008; Oppikofer et al., 2008, 2009; Arnesto et al., 2009; Dunning et al., 2009). Importantly, these techniques provide noninvasive methods for rock-fall analysis in areas that are technically challenging or otherwise hazardous to access.

We collected airborne laser scanning (ALS) data, ground-based terrestrial laser scanning (TLS) data, and gigapixel photographs for the cliffs beneath Glacier Point in eastern Yosemite Valley (Fig. 1A) in September 2006, October 2007, and May 2008, respectively. Subsequently, two large rock falls occurred within the imaged area, the first at 13:30 Pacific Standard Time (PST) on 7 October 2008, and a second, larger rock fall from the same area at 05:55 PST on 8 October 2008 (Fig. 1B). Both rock falls occurred adjacent to an area of previous instability that had failed in August and September of 2001 (Wieczorek and Snyder, 2004). On both occasions in October 2008, the rock masses free fell ∼220 m, impacted a prominent east-dipping joint-controlled ledge, and fragmented into numerous boulders and smaller debris. Much of this debris traveled down a talus slope into Curry Village (Fig. 1B), causing minor injuries and damaging or destroying 25 buildings. In terms of structural damage, this was the most destructive rock fall in the history of Yosemite National Park. In the ten days following the second rock fall, we repeated the TLS and high-resolution photography of Glacier Point. The resulting pre– and post–rock-fall data, integrated to maximize data resolution, allow us to image the failed rock mass in three dimensions, quantify the rock-fall volume, evaluate the geologic structure that contributed to instability, and assess the likely failure mode.

Gigapixel photography is a digital mosaic approach used to achieve image resolution that far surpasses that of conventional digital photography, i.e., creation of individual images consisting of ≥1 billion pixels, 100 times the resolution of a standard 10-MP digital camera (e.g., Frenkel, 2010). Gigapixel photography uses a robotic control device that divides a field of view into several hundred positions that are shot with a telephoto lens and then stitched together to create a single high-resolution digital image.

We obtained baseline gigapixel photography for Yosemite Valley by simultaneously photographing the rock faces from 20 locations along the valley rim (http://www.xrez.com/case-studies/national-parks/yosemite-extreme-panoramic-imaging-project/). At each location we collected ∼400–700 high-resolution overlapping digital photographs using a motion-controlled camera tripod (GigaPan™ unit) and a Canon G9 camera with a 2× extender (effective 300 mm focal length). We then aligned the overlapping photographs, stitched them together using PTGui™ software, and rendered the stitched images to create 20 individual gigapixel panoramic images of Yosemite Valley. Images can be viewed online at http://gigapan.org/gigapans/most_recent/?q=xrez and at http://www.xrez.com/yose_proj/Yose_result.html. As part of this process, we captured a 3.7-gigabyte panoramic photograph (63,232 × 20,224 pixel resolution) of the cliffs beneath Glacier Point on 28 May 2008 from a position on the opposite valley wall (Figs. 1 and 2). Two days after the 8 October 2008 rock fall, we repeated the high-resolution photography of Glacier Point, creating a detailed panoramic image of the rock-fall source area and surrounding cliff before and after the rock falls occurred (Fig. 2).

In September of 2006, the National Center for Airborne Laser Mapping (NCALM), in collaboration with the National Park Service, collected ALS data for Yosemite Valley and vicinity (Fig. 1A), an area of ∼43 km². ALS data were collected with an Optech 1233 ALTM scanner mounted in a turbocharged twin engine Cessna 337. Above–ground-level flying heights varied from less than 100 m to over 2 km, with an average range of 1050 m. The Glacier Point area (blue box in Fig. 1) consists of ∼9.1 million data points, corresponding to a horizontal plane point spacing of ∼60 cm.

The ALS point-cloud data and digital elevation models (DEMs) derived from them resolve the lower angle topography and upward-facing surfaces in the Glacier Point area in high resolution. However, because ALS is a downward scanning system it cannot fully resolve topographic overhangs where there may be three or more surface measurements for a given x-y coordinate: the uppermost surface, the overhang surface, and the cliff face underneath the overhang. Typically, DEMs created from ALS point-cloud data only use the uppermost position, and are therefore unable to image and characterize vertical to overhanging surfaces, such as the 2008 rock-fall detachment surface. In order to fully resolve the complex topography of the Glacier Point area, we expanded the spatial resolution by also collecting ground-based TLS data.

We collected TLS data for the cliffs below Glacier Point from a position on the northern edge of Stoneman Meadow, 1.23 km line-of-sight distance from the rock-fall detachment surface (Fig. 1). We used an Optech ILRIS-3D^ER extended range scanner to collect pre–rock-fall TLS data on 11 October 2007, and post–rock-fall TLS data on 18 October 2008. The scan area of ∼612,000 m² covered nearly the entire vertical extent of the cliff (728 m; Fig. 3). Glacier Point is a challenging TLS target due to high incident scanning angles from the valley floor, dark rock surface staining, and smooth granitic surface properties that can reflect the laser signal away from the scanner, all of which tend to yield low signal return. Furthermore, laser attenuation reduced the signal from surfaces >1 km from the scanner. The TLS scans consisted of ∼8 million data points, corresponding to a vertical plane point spacing at the rock-fall detachment area of ∼50 cm. Because of the short time interval between the 7 and 8 October rock falls (16.5 h), we were unable to repeat the TLS surveys until after the second event; thus, the rock-fall volume and other metrics reported here are cumulative for both events.

To aid our three-dimensional characterization of the rock-fall detachment surface and adjacent cliff area, we integrated the gigapixel photographs with the LiDAR data using 3D point-cloud meshing and 3D animation software. As described above, ALS point-cloud data do not accurately resolve vertical to overhanging cliffs because they cannot capture surfaces with multiple z values for one x-y position. In contrast, TLS data accurately resolve vertical and overhanging surfaces, but our TLS point-cloud data have large gaps due to shielding of lower angle surfaces from the scanning position on the valley floor, and also due to laser attenuation beginning at distances >1 km (Fig. 3). To compensate for these data gaps, we merged the ALS and TLS point clouds into a single point cloud using a best-fit alignment algorithm described below (Fig. 4). By merging the ALS and TLS point clouds, we are able to resolve the complex topographic surfaces of the Glacier Point area in three dimensions.

Once the ALS and TLS data were merged, we projected the gigapixel photograph of Glacier Point onto an interpolated surface created from the point clouds. To do this, we imported the merged point cloud into VRMesh™ 5.0 and interpolated surfaces from the point cloud as described below (Fig. 5A). We then imported the surface model in Maya™ 3D animation software as an exported object file from VRMesh and properly scaled and positioned it. We chose a smoothing angle to create shading along the surface, then applied texture coordinate mapping onto the surface. Finally, we projected a texture map derived from the gigapixel panoramic photograph onto the surface, yielding a three-dimensional form of the photographic imagery that reveals the morphology, structure, and texture of the cliff below Glacier Point in high resolution (Fig. 5B).

Repeat TLS scans before and after the rock falls occurred allow us to precisely calculate the volume of the failed rock mass. To do so, we aligned the 2007 and 2008 TLS point-cloud data using InnovMetric PolyWorks™ software (InnovMetric, 2010), first using manual point-pair matching and then using the surface-to-surface iterative closest point algorithm IMAlign. This routine creates a best-fit surface through the point-cloud data for each time period and then uses a least-square inversion approach to minimize the misfit between the two epochs. We aligned on common cliff-face data points outside of the rock-fall detachment area.

As described above, variations in reflective properties of the cliff face and laser attenuation produced data gaps in the TLS point clouds (Fig. 3B). To fill these gaps, we generated surface models using three different interpolation methods: kriging, triangular irregular network (TIN), and inverse distance to power (IDP). Kriging and IDP use a weighted distribution of neighboring points to project the cliff surface across data gaps. For both methods, we used a linear search of 10 m and 0.25 m spot spacing; for kriging we used zero power and for IDP we used a power of 2, both calculated with Surfer™ Version 8 software (Golden Software, 2010). We also created a TIN model, which connects triangular surfaces between adjacent points. All three approaches yield reasonable rock surface models because there are very few data points associated with vegetation on the cliff face; however, kriging and IDP more effectively model data gaps by incorporating neighboring points. To take advantage of the interpolation routines and still account for the overhanging character of the rock-fall detachment area, we rotated the xyz coordinate system from x into the cliff, y along the cliff face, and z is up, to a coordinate system where x-y are in the cliff face and z is perpendicular to the cliff face.

Once these surface models were created, we calculated the volume change at the rock-fall source area between the 2007 and 2008 models using Applied Imagery Quick Terrain Modeler™ (Applied Imagery, 2010) and visualized and assessed the data misfit with LiDARViewer (Kreylos et al., 2008). Although the alignment precision for the two TLS data sets is on the order of centimeters, there is additional uncertainty in the volumetric change associated with the interpolation methods. We assessed this uncertainty by calculating volume changes for five circular areas (55.4 m diameter, the diameter of a circle with the same approximate area as the rock-fall detachment surface) on the cliff adjacent to, but outside of, the rock-fall detachment area (Fig. 6A). The high-resolution photographs confirm that these areas did not experience rock falls, and thus did not change volume, between the 2007 and 2008 scans. Determining volume changes for these areas provides a measure of volumetric uncertainty associated with each interpolation method (Table 1).

Comparison of the 2007 and 2008 TLS-based surface models (Fig. 6) reveals that the rock-fall detachment surface is 69.0 m along its longest axis (A–A′; Figs. 6B and 7A), and has a total surface area of 2750 m². The failed rock mass, which was approximately lens-shaped, had a mean thickness of 2.1 m and a maximum thickness of 7.1 m near the upper, eastern corner (Figs. 6B and 7A); these thickness calculations are confirmed by the actual measured thicknesses of fresh boulders on the talus slope (Fig. 7B). Cross sections through the slab reveal that it was of relatively uniform thickness across most of its width, with the detachment surface remarkably parallel to the pre-failure cliff surface (Fig. 7A, Supplemental File 1¹, Animation 1). Kriging and TIN interpolation methods both yield rock-fall volumes of 5667 m³, although the kriging method has greater uncertainty (±36 m³ versus ± 27 m³, respectively; Table 1). The IDP method yields a slightly smaller volume of 5658 ± 44 m³. Our best estimate of the cumulative rock-fall volume is 5663 ± 36 m³ (the error-weighted mean and uncertainty of the three interpolation methods). Photographs of the source area taken immediately after the 7 October rock fall suggest that ∼20% of the total volume (∼1133 m³) is attributable to this first event and the remainder (∼4530 m³) attributable to the subsequent 8 October rock fall (Fig. 2B). Notably, observation-based estimates of the cumulative volume made immediately after the events underestimated the actual volume by roughly a factor of two, primarily because the relatively thin (≤1.0-m) overhangs at the top of the detachment area proved to be a poor indicator of the mean slab thickness (4.1 m). This highlights the difficulty of attaining accurate rock-fall volumes without quantitative topographic data. Our results also highlight the uncertainty associated with estimating rock-fall volumes from fresh talus deposits. Talus resulting from the 7 and 8 October 2008 rock falls is spread over an area of ∼118,000 m² on the talus slope beneath Glacier Point. The largest boulder on the talus slope accurately records the original slab thickness (Fig. 7B), but represents just 2.5% of the total rock-fall volume. Furthermore, some fresh-appearing talus was actually older debris on the talus slope that was remobilized by the event. Estimating the cumulative rock-fall volume from fresh talus deposits alone would likely have also led to a substantial underestimation of the actual volume. These discrepancies illustrate the importance of repeat high-resolution topographic data for accurately determining rock-fall volumes for individual events.

High-resolution, three-dimensional imaging helps to evaluate the geologic structure that contributed to failure. The October 2008 rock falls were clear examples of granitic exfoliation along a sheeting joint (Matthes, 1930; Huber, 1987). Images of the detachment surface and adjacent cliff reveal that the dominant structural feature controlling detachment was a vertically oriented, near-planar sheeting joint, which is part of the surface-parallel J₁ joint set (Wieczorek and Snyder, 1999; Wieczorek et al., 2008). Sheeting joints are common in Yosemite Valley, and they often form detachment surfaces for rock falls (Matthes, 1930; Huber, 1987; Wieczorek and Snyder, 2004). As determined by plane-fitting to the TLS data, the primary J₁ joint-controlled detachment surface of the October 2008 rock falls is oriented 027°/89° (dip direction/dip angle), identical to the orientation of the cliff face prior to failure (Fig. 7A, Supplemental File 1 [see footnote 1], Animation 1). In fact, the detachment surface very closely mirrors the pre-failure cliff surface, not only in its overall orientation but also in the location and magnitude of surface convexities (Fig. 7A, Animation 1). This tends to support the suggestion that cliff surface morphology, in particular the degree of curvature, may strongly control the development of sheeting joints (Martel, 2006) and thus influence rock-fall susceptibility.

The failed slab was further bounded on its upper and lower edges by several predominantly east-dipping joints, part of a pervasive joint set, termed J₂, which is prominent throughout the Glacier Point area (Wieczorek and Snyder, 1999; Wieczorek et al., 2008; Fig. 5B). The dominant J₂ joint exposed directly above the detachment surface (Fig. 5B) has a dip direction/dip angle of 094°/30°. The detachment surface was further bounded on its upper western edge by a series of subvertical fractures (Fig. 2B). The lower western edge of the failed slab consisted of an overhang that resulted from earlier rock falls occurring on 14 and 25 September 2001 (Wieczorek and Snyder, 2004) (Fig. 2B). These bounding features provided structural weaknesses that likely contributed to instability. A prominent light-colored dike extending across the upper portion of the detachment surface (Fig. 8C) may have contributed to the greater thickness of the failed slab in the upper northeast portion of the detachment area (Figs. 6B and 7A).

The high-resolution data provided by gigapixel photography, ALS, and TLS help to clarify the likely failure mode of the October 2008 rock falls. Typically, rock falls from steep cliffs fail by one of two modes, shear (sliding along the cliff surface) or tension (rotation away from the cliff surface, also known as toppling; Goodman, 1989). These can be analyzed using various methods (limit equilibrium, constitutive modeling, fracture mechanics, etc.) Here we explore the use of limit equilibrium, where the likelihood of one or another failure mode is determined by the relationship between driving forces (the geometry and volume of the rock mass and external forces such as cleft pressures) and resisting forces (the shear and tensile strengths of the rock material) (e.g., Norrish and Wyllie, 1996; Wyllie and Mah, 2004).

Our use of the limit equilibrium method provides a first-order failure analysis of rock slab detachment, as is commonly performed in conventional soil and rock slope stability analysis. Using these methods, we compared the driving forces and moments of the failed slab, based on the measured volume and detachment surface geometry, to the resisting forces and moments of the slab (Figs. 8A and 8B). We assumed typical values for granite rock unit weight (26.5 kPa, Goodman, 1989) and low values for strength (friction angle = 31°, Jaeger et al., 2007; cohesion = 25,000 kPa, assumed value; overall shear strength = 13,445 kPa, West, 1995; tensile strength = 6688 kPa, West, 1995); these values were verified by comparison with site-specific values for Sierra Nevada granodiorite (Krank and Watters, 1983). Further, we assume that there is only minimal variation between strength parameters for granite and the light-colored dike that runs through the detachment surface (Fig. 8C), and thus used uniform parameters for both rock types. Analysis of high-resolution imagery indicates that freshly broken surfaces were distributed over seven areas on the detachment surface, encompassing ∼26% of the total surface area; the remaining areas display slight staining or weathering, suggesting earlier detachment of these areas (Fig. 8C). Based on these observations, we assume that prior to failure the majority of the slab was detached from the main cliff surface, and we calculated strength contributions only for those portions of the slab that were previously attached at these freshly broken surfaces. We analyzed shear failure based on both an inclined plane formulation dipping 89° using Mohr-Coulomb shear strength parameters, and a static vertical analysis using an overall shear strength parameter (Fig. 8A; West, 1995) to calculate the factor of safety for sliding (i.e., the ratio of shear strength to rock slab gravitational stress). We analyzed tensile failure based on a moment analysis of the slab acting at a point 2.0 m from the detachment surface and rotating about a pivot point located at the bottom of the lowest attached surface (Figs. 8B and 8C), held in place by internal tensile strength acting along the fresh areas of broken rock. We calculated the acting location of tensile forces through general mechanical analysis of a composite tensile strength centroid, resulting in an acting tensile strength vector positioned 31.1 m above the bottom of the lowest attachment surface (Fig. 7C). This provided a factor of safety for tensile failure (i.e., the ratio of tensile strength to rock slab outward rotational moment). We also considered a third potential failure mode, that of lateral shearing (i.e., tearing) along the detachment surface, but did not have sufficient information on the likely point of rotation to develop a meaningful analysis.

The results of these calculations imply that shear (sliding) failure was the more likely failure mode because the calculated factor of safety for shear failure is nearly five times lower than that for moment-driven (tensile rotation or toppling) failure. However, these analyses highlight the limitations of using limit equilibrium methods to back-calculate the stability of exfoliating rock slabs, in that both the shear and moment driven analyses yielded factors of safety that were considerably higher than those defining instability (i.e., safety factors >>1). This suggests additional driving forces acted to initiate shear failure, but the static limit equilibrium methods used in this analysis are not capable of specifically identifying these forces. Observations indicating that the detachment surface was dry at the time of failure suggest that water pressures did not act as an additional driving force; however, we cannot rule out the possibility of increased cleft pressures immediately prior to failure. Stress concentration at crack tips and resulting fracture propagation has been previously proposed as a driving force for rock failures (e.g., Martel, 2004; Ishikawa et al., 2004), including exfoliation events in Yosemite Valley (Bahat et al., 1999; Wieczorek and Snyder, 1999, 2004), and fracture propagation likely played a role in the 7 and 8 October 2008 rock falls. Fracture propagation is consistent with reports of cracking sounds preceding the 7 October 2008 rock fall by a few hours. We tentatively suggest that a change in the equilibrium of the rock slab, for example, stress redistribution on fracture surface asperities that exceeded the shear or tensile strength in the areas of attachment, likely led to fracture propagation. It is possible that the earlier rock falls in August and September of 2001, also dry failures without recognized triggering mechanisms (Wieczorek and Snyder, 2004), served to destabilized the rock slab, with fracture propagation occurring sporadically until the October 2008 failures.

Shear or tensile failure of a rock slab along surface parallel sheeting joints highlights a remaining challenge for stability assessments based on laser scanning data. Because they form within the rock mass roughly parallel to topographic surfaces, sheeting joints are often difficult or impossible to detect even with high-resolution imaging. Laser scanning has great potential for identifying rock mass configurations susceptible to failure (e.g., Jaboyedoff et al., 2007; Lato et al., 2009; Sturznegger and Stead, 2009), but more work is needed to remotely detect and monitor surface-parallel sheeting joints, which form the detachment surfaces for many rock falls in Yosemite Valley.

Integrated repeat gigapixel photography and airborne and ground-based terrestrial laser scanning provide accurate and precise volume calculations and three-dimensional geologic characterization of the October 2008 rock-fall source area that could not be otherwise attained by traditional assessment methods. They also provide a means of quantifying uncertainties associated with these calculations. Our results obtained with a long-range laser scanner demonstrate that high-precision topographic data are attainable for vertical rock faces >1.2 km distant, and that rock-fall volumes can be computed at these distances with uncertainties of <1%. Merging ALS and TLS point clouds, and integrating with gigapixel photography, provides unprecedented imaging capabilities for rock-fall analysis in areas that are technically challenging or otherwise hazardous to access. Volumetric, structural, and other geologic and topographic data pertaining to rock falls are critical for deriving accurate and precise hazard assessment based on probabilistic or deterministic methods, and allow for evaluating potential failure modes. The vast improvements in quantitative analyses for tall cliffs resulting from integrated high-resolution imaging techniques should lead to a reduction of related uncertainties in rock-fall hazard assessments that rely on these analyses.

We thank David Haddad and Michel Jaboyedoff for constructive manuscript reviews, and Steve Martel and Jonathan Stock for commenting on an earlier draft. We appreciate the assistance of the many volunteer photographers who helped to acquire the high-resolution panoramic images. ALS data were acquired by the National Center for Airborne Laser Mapping (NCALM) at the University of Houston and processed by NCALM at the University of California, Berkeley. NCALM is funded by the National Science Foundation. This work was supported in part by funding from the Yosemite Conservancy.