Lidar mapping of faults in Houston, Texas, USA

Accurate maps of Holocene surface faults help developers avoid building in zones of disturbed ground along the faults and provide insights on regional tectonics. These faults have deep-seated origins in the formation of the Gulf of Mexico (Sheets, 1971). They also appear to have been recently reactivated due to human withdrawal of groundwater and petroleum (Holzer, 1984).

The faults considered in this paper are part of a population of hundreds of faults that cut Pleistocene and Holocene sediments on the coastal plain between Beaumont and Victoria (Verbeek, 1979). Verbeek estimated that more than 10% of these faults have been active during the twentieth century. These faults belong to a regional system of down-to-the-basin faults along the Texas and Louisiana coasts.

Houston, Texas, lies within the broad passive margin of the Gulf of Mexico ocean basin (Fig. 1). Gulf of Mexico extension began with Triassic rifting (Salvador, 1991) followed by seafloor spreading during the middle Jurassic (Bird et al., 2005). Deposition in the northwestern Gulf Coast region has resulted in the progradation of the shelf margin into the Gulf of Mexico basin throughout the Cretaceous and Cenozoic (Winker, 1982). Paleogene deposition was predominantly in the South Texas area, but Neogene deposition has been focused on East Texas and South Louisiana (Williamson, 1959). The region of most active growth faulting occurs near the prograding shelf margin. The Oligocene prograding shelf margin was in the Houston area. Faulting occurs in Houston at present, but at a slower rate than at the shelf margin.

Many active surface faults in the Houston area have been mapped prior to our study (e.g., Clanton and Amsbury, 1975; Heuer, 1979; Verbeek et al., 1979; O'Neill and Van Siclen, 1984; Shaw and Lanning-Rush, 2005). Figure 2 shows active surface faults that we investigated using lidar data. Rates of local movement on individual faults can be as much as three centimeters per year (Norman, 2005). Because there are no recorded earthquake epicenters in the Houston area, fault motion by aseismic creep is inferred. Houston fault motion exhibits both spatial and temporal variability (Mastroianni, 1991). Movement on active faults has caused damage to buildings, pipelines, petroleum and water wells, sewer lines, power lines, and roads. In some cases, fault locations have remained unknown until accumulated slip resulted in significant damage. Ongoing building maintenance has, in some cases, minimized damage from active faults.

The locations for many active surface faults in the Houston area are not published, mainly because some data are proprietary (i.e., site studies produced by geotechnical consultants for clients). Also, the threat of lawsuits limits publication. The principal legal liability is that the absence of a fault from a map may be incorrectly interpreted as indicating that a fault is not present. Unlike California's Alquist-Priolo Act (Hart and Bryant, 1997), Texas and Houston lack regulations concerning construction in fault zones. In some areas where fault locations have been recognized, the land along the fault has been left open, often as a park or as a storm water detention pond.

Many of Houston's surface faults have been correlated with subsurface faults (Heuer, 1979; Sheets, 1979; Van Siclen, 1967). Verbeek et al. (1979) recognized that these faults are growth faults. Down-to-the-coast faults, representing extension into the basin, dominate and define a regional southwest to northeast trend (Sheets, 1971). Antithetic surface faults in Houston are found opposite the most active sections of the down-to-the-coast faults (Norman, 2005). Many of these antithetic faults are currently active.

Most (80%) of the faults in the Houston area occur over salt domes (Norman, 2005). Many are radial faults (Verbeek and Clanton, 1978) and are typically short and commonly form grabens. Salt domes and directly associated faults predominate in the southeastern part of the region. The association of radial faults and salt domes is well known and readily understood (Schultz-Ela et al., 1994).

The faults in the northwestern Houston area express the regional trend and are the focus of our study. There are three main fault systems in this area—the Hockley-Conroe Fault System, Addicks Fault System, and Long Point-Eureka Heights Fault System (Figs. 2 and 3). The Hockley-Conroe Fault System extends well outside of Harris County, where we have data. Interference of drainage patterns with fault scarps precludes clear identification of a northeast continuation of the fault within Harris County on lidar images. The Addicks Fault System extends from the Barker Reservoir north-eastward toward Bush Intercontinental Airport. The active antithetic Lee Fault is southeast of the airport, suggesting that an active continuation of the Addicks Fault System runs through the northwest corner of the airport. The Long Point Fault is one of the more active faults in the Houston area and is probably the most studied fault in the region.

Salt tectonics may also play a significant role in the fault activity in the northwestern Houston area. In this model (Jackson et al., 2003), the withdrawal of salt causes additional accommodation space in the hanging wall, with salt welds permitting gravity gliding.

Surface subsidence in Houston can be related to several causes including water withdrawal, sediment compaction, and surface faulting. The abandonment of the Brownwood subdivision illustrates the significance of subsidence in the Houston area (Coplin and Galloway, 1999). The Brownwood subdivision was built in the late 1930s. Initially, elevations were approximately 3 m above sea level, but by the late 1970s subsidence of more than 2.5 m had occurred. The subdivision was subject to frequent flooding. Subsidence in this area was related to groundwater withdrawal for petrochemical plants along the Houston Ship Channel and for the city of Baytown. Recent InSAR (interferometric synthetic aperture radar) measurements (Stork and Sneed, 2002) confirm the abatement of subsidence in the Houston Ship Channel area and ongoing subsidence in the Addicks reservoir area of west Houston.

The interaction between Houston surface faults and subsidence is complex and not well understood. Holzer and Gabyrsch (1987) find a temporal correlation between groundwater withdrawal and the amount of fault slip. Kreitler (1976) finds that the faults compartmentalize subsidence. The association of kilometer-scale subsidence depressions in the Houston area with active faults has not been well established (O'Neill and Van Siclen, 1984). Paine (1993) attributes increased subsidence since the Pleistocene throughout the Texas coastal region to oil and gas withdrawal since groundwater use is absent in many areas. Glacial-isostasy also plays a role (Gonzalez and Tornqvist, 2006). Dokka (2006) attributes subsidence to tectonic factors in coastal Louisiana.

Lidar is an acronym for light detection and ranging (analogous to radar, but with laser light as a source). A laser beam is directed toward the ground as the survey aircraft flies over the area (Fowler, 2001). The laser beam is pulsed so that the time between the pulse and the returned reflection can be measured to determine the distance between the instrument and the ground. Inertial guidance and global positioning system (GPS) units track and record the position and orientation of the airborne platform. The two-way laser travel time, inertial measurement unit (IMU), and GPS are integrated to determine the x, y, z position of each reflection (Wehr and Lohr, 1999). Reflections are returned from the ground, vegetation, and buildings. Many lidar systems acquire multiple returns from a single outgoing laser pulse, which permits penetration of the vegetation canopy. One of the key processing tasks is to remove returns from buildings and vegetation, yielding what is known as a bare-earth model. Different filters have been used to separate ground returns from vegetation and human constructions (Sithole and Vosselman, 2004). Frequently, other data such as GIS locations of buildings and aerial photos are used to aid in classifying returns.

Lidar allows for rapid and accurate generation of digital elevation models (DEMs). A DEM is a grid in which the elevations are equally sampled in both x and y directions (Maune, 2001). Typically the grid is formed by generating a triangulated irregular network (TIN) from the filtered points. The TIN is then interpolated to create the DEM grid. Lidar-derived DEMs and hillshading techniques have been used for detecting surface faults (e.g., Haugerud, 2003; Wieczorek et al., 2004; Egnew, 2005), primarily in areas with significant geologic structure and rarely in urban environments.

The data used in this study were acquired as part of the Tropical Storm Allison Recovery Project (TSARP). Tropical Storm Allison produced up to one meter of rainfall over a five-day period in June of 2001 (TSARP, 2004). The study was sponsored and funded by the Federal Emergency Management Administration (FEMA) and the Harris County Flood Control District (HCFCD) (Meyer, 2002) to update floodplain maps. Terrapoint gathered and processed the lidar data in late 2002. The survey was done with fixed-wing aircraft flown at an elevation of 1000 m with a swath width of 700 m (Quarles et al., 2002). The Terrapoint sensor employs a rotating mirror and gathers 20,000 points per second. Up to four returns may be obtained. Orthophotos were obtained during the same time period (TSARP, 2004) along with detailed topographic cross sections along stream channels (Meyer and Corbley, 2004). More than 10,000 cross sections were taken and used to correct the raw lidar data in producing the DEM. The cross sections began and ended outside of the channel bank and included measuring water depth. GPS was used to tie each cross section to the map with conventional survey methods employed to measure the elevation changes, referenced to North American Vertical Datum (NAVD) 1988 (2000 adjustment). A survey of 1500 monuments used by HCFCD was also done, with work completed in April, 2003 (TSARP, 2003). Determination of non-ground returns was aided by examination of the aerial photos and using GIS data showing the location of buildings (TSARP, 2001). The lidar elevations were linearly interpolated from TINs after removal of non-ground returns. Examination of the resulting DEM shows a few artifacts that should have been removed (Fig. 3). Vertical accuracy was checked for each ground cover class. The overall root mean square (RMS) elevation error was found to be 11.6 cm (Meyer, 2002).

The bare-Earth DEM grid obtained from lidar is 23,000 by 17,000 nodes (1.5 gigabytes). This lidar grid has a 5-m cell size; but, the density of the raw data (30,000,000 points per quarter quadrangle) permits a finer grid spacing, such as 1.5 m. Increased resolution allows us to more accurately map scarps. The data consist of the bare-Earth DEM and the raw point cloud. We developed a grid-refinement algorithm (described below) to achieve higher resolution and avoid refiltering of the raw data. This algorithm is more computationally efficient than starting with the raw data, but it may also produce smoother results owing to the input grid resolution. Grid refinement also accommodates the survey adjustments to the DEM (noted above). Refinement was carried out for selected quarter quadrangles. In this process, the raw data were gridded along with the existing coarse (5-m cell size) DEM to generate a finer (1.5-m cell size) DEM.

Our grid-refinement algorithm compares the elevation of each raw point with the interpolated value of the original DEM. A trapezoidal filter weights the points according to this elevation difference. If the difference was less than a specified minimum threshold, then the raw point was kept with a weight of one. If the difference exceeded a maximum threshold, the point was discarded. Points whose difference fell between the two thresholds received a linearly interpolated weight. The weighting was one at the minimum threshold and reduced to zero at the maximum threshold. Experimentation and consideration of non-ground returns resulted in choosing threshold values of one and two meters for refinement. Our basic algorithm is outlined below, with z for the raw data value, dz the absolute difference between the raw value and input DEM value, dmin the minimum threshold, dmax the maximum threshold, and w the resultant weight:

• if dz < dmin

• w = 1

• else if dz > dmax

• w = 0

• else

• w = 1 – (dz − dmin)/(dmax − dmin).

Weighted bilinear interpolation is performed for the raw points. The weighted elevations and weights were summed in the four adjacent output grid nodes with further weighting based upon the position of the point within the grid cell. Values for each node were computed from the summed weighted values normalized by the summed weights. Bilinear interpolation of the original grid was used where necessary to fill out the refined grid.

Hillshading involves illuminating the surface from a particular sun position. Surfaces inclined toward the light are lighter than average, whereas surfaces inclined away are darker. We tried other techniques such as contouring, calculation of derivative products such as slope, aspect, various overlays of different derived grids, and direct examination of the DEM. These did not prove nearly as effective in scarp identification, with the results being harder to interpret (Fig. 4),⁴⁰². Other workers (e.g., Ganas et al., 2005) also found hillshading to be the most effective method for visualizing surface faults. Hillshading is controlled by both light-source elevation and azimuth. Low-elevation angles were noisy, highlighting streets as well as scarps; whereas high angles resulted in less contrast. A medium angle, such as 45°, proved best for this study. Choice of azimuth influences slope recognition. With the regional SW-NE fault trend (Fig. 4B),⁴⁰² an azimuth of 315° highlights the scarps in northwestern Houston. An azimuth of 225° did not highlight any additional scarps in northwestern Houston. However, in southeastern Houston, some scarps were only highlighted by an azimuth of 225°. Other azimuth angles did not yield results that could not be seen with 315 and 225° azimuths.

Visual interpretation of the hillshaded images was used to identify potential fault scarps. Map scales of 1:10,000 or larger were used (Fig. 5). In this location, the scarp elevation contrast in the DEM is prominent, but in many other locations, the elevation contrast is hard to observe. The hillshaded image provides a clear indication of the scarp. Scarps were visible in other locations only with hillshading.

The locations of published faults were posted on the hillshaded image and used as an initial guide for fault interpretation. Two fault maps were used—a recent digital version extracted from the Geologic Map of Texas (Fisher, 1982) and georeferenced versions from scanned older maps of faults (Verbeek and Clanton, 1978; Verbeek and Clanton, 1979; Verbeek et al., 1979). The digital map is more complete but less accurate than the older scanned maps. Our new maps are more accurate than both, with scarp locations mapped within a couple of meters when using refined grids. In some cases, the old maps were in error by more than 100 m (Fig. 6).

Lidar imagery emphasizes elevation differences, such as scarps. Aerial photos, in contrast, often show differences in vegetation on the upthrown and downthrown sides of the scarp because the downthrown side is typically wetter. The Long Point Fault can be clearly seen in the lidar image (Fig. 7A), but it cannot be seen in the orthophoto (Fig. 7B).

In some places, the elevation transition across the scarp is broad and gently dipping. Figure 8 shows one model for this geometry. Erosion removes material from the footwall and deposits it on the hanging wall (Verbeek and Clanton, 1979). Deformation near the surface results in a zone of shear in the unconsolidated sediments, rather than a distinct fault (Thomas and Nunn, 2006). Erosional leveling and distributed deformation can both contribute to scarp morphology.

Some maps have categorized the Renn Scarp in southwest Houston as a fault. Drilling has confirmed that the scarp is the cutbank of an ancient stream channel. The scarp was subsequently obliterated by recent urban development. Ground observations are required to confirm whether a particular feature corresponds to a fault or has some other cause.

Many visits were made to field locations of scarps identified by the lidar. In some cases, there is a clearly visible scarp, whereas in others, the elevation changes detected by the lidar are subtle and unlikely to be readily noticed on the ground. Photographs of the sites were taken to record their character. Scarps were visited where access was permitted. Access along roads proved easiest, although roads and other constructions alter terrain. Recent modifications tend to smooth out and flatten the scarps. Scarps in fields are often difficult to find because of vegetation cover. In areas where the fault has been recently active, there is sometimes cracked pavement as well as a noticeable scarp. Deformation of manmade features along the trace of an extensive scarp is a strong indication of an active fault (e.g., Fig. 9).

Our analysis of the lidar data confirmed the location of many of the previously mapped active faults (Fig. 2) and found several additional scarps that may be faults (Fig. 3). Field visits did not provide a clear indication of faulting at these locations. Subsurface studies are needed to help decide whether such features are likely to be faults. Possible tools include shallow borehole well logging, ground penetrating radar (GPR), and high-resolution seismic reflection profiling.

Fault scarps sometimes have kinks and bends. In many cases, local kinks are due to construction of a parking lot or leveling of property for a house (Fig. 10). A field visit shows the scarp on the east side of the street is along the property line between two houses. This strongly indicates that the lot on the south (downthrown) side was leveled by grading. Continued fault activity may damage this house.

Normal faults commonly bifurcate (Shelton, 1984). One prominent Houston-area example occurs toward the eastern end of the Long Point Fault, north of where it passes over Long Point Road (Fig. 11). We found two to three overlapping arcuate scarps here with the lidar imagery. The orientation and shape of the scarps suggest they are not manmade features. Field observations show pavement deformation, indicating these may be fault scarps. Our interpretation is supported by an observation of Ambs (1980), who found bifurcation on the Woodgate fault.

Regional subsidence depressions do not show associated systematic fracture patterns using the lidar data. Apparently the deformation is fairly evenly distributed from the depression edges toward the depression center. Figure 12 shows lidar imagery covering the area of the Addicks subsidence depression (Stork and Sneed, 2002). The subsidence depression does not form a prominent topographic feature, even after removing a linear trend. The only prominent geomorphic features are stream courses, ditches, and the major faults. We removed the linear trend with a custom program that computed a least square surface (Davis, 1973), which was subtracted from the DEM grid to produce the residual grid. The locations of surface faults in northwest Houston do not correspond well with regional subsidence contours for other periods either. However, Pratt and Johnson (1926) reported fissures resulting from the subsidence due to oil pumping from the Goose Creek oil field. Fluid withdrawal at Goose Creek may have been more rapid and localized than in the Addicks subsidence depression. Goose Creek overlies a salt dome; therefore, oil extraction may have activated existing faults.

Scarps must be evaluated to determine whether they are caused by faulting or some other process. Correlation of subsur-face horizons between boreholes on opposite sides of a fault provides a clear indication of displacement and also confirms downthrown side thickening of strata that is characteristic of growth faults. However, boreholes are expensive to drill and log, and also require permission from landowners. Initial geophysical studies at two locations on the Long Point Fault show a strong indication of faulting on GPR data (Engelkemeir and Khan, 2007). The GPR survey of the possible Long Point splinter fault does not show evidence of faulting. Saribudak and Van Nieuwenhuise (2007) used several geophysical techniques to image anomalies across a Houston area fault.

Additional investigation is required to determine whether a particular fault is active or has been recently active in a particular area. The association of pavement and building cracks with scarps indicates recent displacement. Fault displacements have been measured for more than 50 locations on 35 Houston faults (Mastroianni, 1991; Norman, 2000) with a leveling system set up at measurement locations and repeated over at least a few years. One possible technique to investigate current fault activity might be Interferometric Synthetic Aperture Radar (InSAR). The vertical component of deformation resolved by InSAR (Hensley et al., 2001) appears to be adequate for the slip rates that have been observed in Houston. Recent InSAR work provides indication of differential deformation across the Long Point Fault (Buckley et al., 2003; Buckley, 2005). Uplift on the salt domes in southeastern Houston should provide corroborating evidence for the model mentioned in section 1.2.

Lidar has been shown to be a highly effective tool for accurately mapping faults in the Houston area. Grid refinement is an efficient technique for improved resolution without the need to refilter the raw data. More accurate and current fault maps will aid developers in building safely. Published paper maps are either at too small a scale for siting requirements or are specific to a particular location. Digital maps do not suffer from these limitations and can easily be used in GIS or computer-aided design (CAD) software in association with GIS data such as street locations. Similar techniques should be applicable to study faults on other parts of the Texas coast and in other active extensional salt basin margins.

The fault locations do not correspond closely with the location of subsidence depressions, although fault activation by shallow ground-water withdrawal remains an open issue. They instead appear to be manifestations of regional and salt dome tectonics.

This and other studies of these faults contribute to broadening our understanding of Houston and Gulf Coast neotectonics. Better mapping of fault locations helps in defining their overall patterns.

We are grateful to Jean Parcher of the U.S. Geological Survey for providing lidar and associated GIS data. Khan thanks Environmental Institute of Houston for providing partial funding. Engelkemeir is supported by Schlumberger for his education. We appreciate the careful review of the manuscript by Carl Norman, Kevin Burke, William Dupre, Christopher Crosby, and Alan Nelson.