EFFECTS OF COASTAL LANDFORM CHANGES ON STORM SURGE ALONG THE HATTERAS ISLAND BREACH AREA

Storm surge is known to be the source of many coastal problems including dune erosion, overwash deposits and inlet breaching. The work discussed in this paper arises out of an effort to develop geospatial tools and workflows that integrate with process based numerical models in order to improve the understanding of the importance of landform and complex landform changes in these processes building on the use of existing models. The effects of landform changes on calculation of storm surge characteristics (elevation and velocity) are investigated for three different cases created using the geospatial workflows coupled with ADCIRC modeling efforts.


INTRODUCTION
Impacts to coastal infrastructure usually occur in developed coastal areas with complex coastal hydraulics such as inlets and in areas with vulnerable landform features such as low-lying barrier island systems.Barrier island systems are susceptible to remarkable landform changes during extreme events that can lead to social and economic losses in the area.The most damaging extreme events are hurricanes due to the large wind waves and storm surge that can be produce.Water level increase caused by the storm surge creates the base for the storm waves to run high on the dunes which may lead to loss of the dune and wash over fans, starting the critical landform change in the form of dune erosion process.While the effect of storm surge on coastal landform change is perceptible, the opposite effect is usually not investigated.In this study we aim to report our findings on the effects of coastal landform changes on storm surge including water elevations and current velocities since both parameters are important to in determining the conditions that will cause the system to fail.

Study Location
The area selected for this study is the Hatteras hotspot area which is between Buxton and Hatteras Village North Carolina, USA on the southeasterly facing portion of Hatteras Island and is part of the Cape Hatteras National Seashore (Figure 1).The area is considered a "hotspot" due to the narrow width (~130 m) of the barrier and low dune elevation (~2.5 m NAVD88) at this location.Hurricane Isabel breached the Hatteras hotspot in September, 2003.The inlet developed multiple channels and eroded to a depth of approximately 6 m (NAVD88) in the main channel (Wamsley and Hathaway 2004).In the months following the storm, the inlet was artificially filled restoring it to approximately the pre-storm topography (Wutkowski 2004).This recent event and the availability of the LIDAR data sets of the area just before and after the Hurricane Isabel make it a robust test bed for studying the effects of coastal landform changes on storm surges.

GIS-BASED WORKFLOW
Over the last decade, LIDAR (Light Detection and Ranging) has become very popular technique for collecting topographic data.LIDAR data has advantages over traditional photogrammetric techniques due to the high point density and ease of data processing.The advance in LIDAR technology and the nearly annual data availability enables us to perform time series analysis of coastal terrain evaluation based on GIS-based workflows (Overton et al. 2008, Mitasova et al., 2009, 2010).Time series analysis, basically, involves data acquisition, systematic error correction, interpolation of digital surface models (DSMs) and per cell analysis of the DSMs.

Data Acquisition and systematic error correction
The lidar data used in the study spans approximately ten years capturing the evolution of the area before Isabel and a five year post-storm recovery period.The lidar data used in the study were downloaded from the National Oceanic and Atmospheric Administration (NOAA) distribution website Digital Coast.The data is the NC State Plane North American Datum 1983 coordinate system with North American Vertical Datum 1988 (NAVD88) and units in meters.The published vertical accuracy of most of the data is 0.15-0.2m and the horizontal accuracy is 2 m and better.The data is provided to the distribution website from various missions carried out by different agencies which include the North Carolina Floodplain Mapping Program (NCFMP), U.S. Army Corps of Engineers (USACE) National Coastal Mapping Program, National Aeronautics and Space Administration (NASA), the NASA / USGS Experimental Advanced Airborne Research Lidar (EAARL), the Airborne Lidar Assessment of Coastal Erosion (ALACE) project by the National Oceanic and Atmospheric Administration (NOAA), and U.S. Geological Survey (USGS) (Table 1).To correct the systematic errors, the NCDOT benchmarks located along NC12 highway centerline are taken as true elevations and the elevations from the DSMs were sampled at these points at the study site.The differences between the DSM and NCDOT elevations and the benchmarks were then calculated.Although the values of the mean and median differences were close to each other, median was considered more appropriate for systematic error correction because of its lower sensitivity to outliers present in some of the data sets.The median difference estimated for each year was then applied to the entire DSM, assuming uniform distribution of the systematic error.The improvement gained by the correction was verified by computing the differences between the corrected DSMs and the NCDOT benchmarks (Mitasova et al. 2010).

Per Cell Analysis
Raster based measures characterizing spatial patterns in terrain dynamics were derived by applying per cell statistics based on time series analysis of the multi-temporal DSMs.The derived measures include the minimum, maximum, coefficient of determination and linear regression rate of elevation change.These derived measures are explained in detail by Mitasova et al. 2009 and2010.In this study, the minimum surface is used as one of the different topographies used to investigate the effects of landform changes on storm surge.The minimum elevation surface is called the core surface and is created by recording the minimum elevation at each 0.5m by 0.5m cell over the study area for each year and compiling the minimum elevation cells together to create the core surface.The core surface indicates the volume of the study area that has never changed during the study period.The resulting DSM is illustrated in Figure 2 below.The maximum surface is developed in a similar manner using the temporal maximum elevation per grid cell (green surface).The core surface in the figure below is representative and is created by using the whole data set mentioned in Table 1.The core surface created for this study is modified as explained below.

STORM SURGE MODELLING
For storm-surge modeling, advanced circulation model for oceanic, coastal and estuarine waters (ADCIRC) is used.ADCIRC requires, as input, a finite element grid built by nodes with corresponding spherical coordinates and elevations.The ADCIRC runs are built on the North Carolina Floodplain ADCIRC grid as described in Blanton and Luettich 2008.This grid is used in the "Original grid" run which served as the reference set for comparisons.In order to study the effects of changing landform on storm surge, two different DSMs obtained by following the geospatial workflow discussed above are used to modify the grid node elevations and create two different ADCIRC grids.These are: the "Core grid" and the "Breach grid".The core surface DSM is created by the methodology described above with the datasets listed in Table 1 except the post Isabel dataset.This is important since the post Isabel dataset has a gap in the data at the breach location.The breach DSM is created by interpolating the USACE bathymetry survey data taken after the storm (Wamsley and Hathaway 2004) into the existing original DSM.

Grid modification
The original grid is directly modified to create different case studies.Grid modifications are carried out by updating node elevations and attributes from the derived DSMs.Using a python code, this approach is efficient since it uses the raster surfaces as the source of information, making the density and the number of points a minor consideration to the modification that will be carried out.The "core grid" is simply created by running the code to update the grid with core surface DSM elevations along the outer banks barrier island system.
For the breach case, the USACE bathymetry data (Figure 3a) is used to create the breach DSM (Figure 3b).Using the survey extent, an area is carved out of the original DSM with a buffer left around the survey extent in order to ensure seamless integration of the DSM and the survey data.The survey data is interpolated using spline tension method and integrated in the DSM (Figure 3b).The ADCIRC nodes that need to be updated for the grid modification are essentially overlaid (Figure 3c) on the new DSM and their elevation information is updated using the routine in the python code.During this process, elevation change is monitored and node attributes are modified for nodes that have new dry/wet status.The final product, the breach grid, is illustrated in Figure 3d.

Model Runs
Three ADCIRC model runs were carried out using the three grids presented above.The grids have approximately five hundred thousand nodes and over a million elements connecting them.45 days of tide ramping runs were carried out for every case.After ramping, the model runs have been hot started from day 45 which corresponds to September 14 th , 2003 midnight EST and the simulation continued for another 5.5 days with the hurricane Isabel winds that were modeled as reported in Vickery and Blanton, 2008.Water elevations and current velocities are recorded at 30 minutes intervals.The results are compared at September 18 th , 2003 5pm EST when the first overwash is spotted on the core grid.This time also coincides with the maximum elevation reached at the location.

NUMERICAL RESULTS OF CHANGING LANDFORMS Storm Surge
To study the effect of landform changes on storm surge, the water elevation values are extracted for all three cases and for the first set of nodes that below sea level in front of the breach location (Figure 4).The first row of along shore nodes has an approximate depth of 1.5 meters.The second row (parallel to the first row) of along shore nodes has an approximate depth of 5 meters.For comparison purposes, the extracted information is then plotted at a time step corresponding to 2pm EST which is before the storm surge reaches the location (Figure 5) and at the time step corresponding to 5pm EST (Figure 6). Figure 5 and Figure 6 include 4 subplots: the 3D view, top view, along shore and cross shore cross section view referred to as a, b, c and d respectively.Cross sectional plots are sampled at locations along transects visualized in the 3D view shown in Figure 5a and Figure 6a.The green color represents the original run, blue color represents the core run and the red color represents the breach run.The data presented in Figure 5 corresponds to a single time step (September 18 th , 2003 2pm EST).At this time, the storm surge starts to pile up on the shore and the water level keeps increasing until 6pm then decreases back to its normal level.Figure 5a illustrates the 3D view of water level (surge + tide) in front of the breach location for the three model runs carried out.Looking at Figure 5b, since the elevation difference magnitudes (~1 cm) are not significant (Figure 5c), we can say that the water level distribution does not have any significant spatial distribution.Looking at Figure 5c, we see that water level at the breach channel entrance is approximately 7 cm lower than the other two cases.One possible explanation is breach bathymetry integrated into the grid is preventing the water level to reach higher levels by letting the surge flow through the inlet.Figure 5d shows that the breach run water levels are getting higher when approaching from off shore but the magnitude of the difference between all three runs is negligibly small (less than 1 cm).We see overall similar results at the time step when the first overwash occurs on the core surface at 5pm EST (Figure 6).A notable and expected result at this time step is that the original grid run water levels are higher than the other two cases in the domain we are plotting (Figure 6a, b) the water levels.We suspect that the topography prevents overwashing thus the surge piles up.Although, this is visible in the plot, looking at Figure 6c shows that the water level differences are negligibly small (average 5 cm) for the original and core grid runs.On the other hand, the difference between these two runs (original grid and core grid) and the breach run is in the order of 35 to 40 cm.

Current Velocity
The velocity and water elevation values for the single time step (September 18 th 2003, 5pm EST) are extracted and presented below.Figure 7 illustrates the velocity comparison at the Hatteras inlet and the breach location.Figure 7a illustrates the result on the original grid run.No overwash occurs during the original grid run.Also, the current velocities at the location followed regular flow patterns consistent with no overwashing of the dunes.At the breach location, at a depth approximately 1.5 meters, the current velocities had an average of 0.322 meters/second in the southwestwest direction (along shore).The average velocity in the same direction at 5 meters depth was 0.8 meters/second.Figure 7b illustrates the results for the core grid run.The average elevation of the dune crest in the original topography is 2.6 meters above MSL at the breach location.At the same location, the average elevation of the dune crest from the core surface is 1.6-2 meters above MSL.Given a maximum water elevation of 1.8 meters at this location, overwash of the core surface occurs.The current flow pattern is also different than the original grid run.The current velocities start to point towards the breach location at and around the shoreline well before the first overwash is observed.The overwash only occurs when the water level reaches a level sufficient to allow the water to flow over the dunes.At 1.5 meters depth, the current velocities had an average of 0.35 meters/second along the shore and 0.9 meters/second at 5 meters depth.Although the current velocity differences with the original grid run are small, the real difference in this case is that the average current velocity overwashing the location is approximately 1.5 meters/second.The third and the last case illustrated in Figure 7c is the breach grid case.At 1.5 meters depth, the current velocities reached an average of 1.9 meters/second along the shore and 1.4 meters/second at 5 meters depth.Since this grid had the breach bathymetry implemented in it, as expected, the current velocities inside the breach channels had higher velocities.The average velocity in the breach channel is approximately 3.2 meters/second.The current velocities are summarized for the three cases in Table 2.The breach bathymetry allows the surge current to flow in to the sound side of the barrier islands.While doing so, the current velocities increase not only at the breach location but also in front of the breach towards offshore.In addition, because of this the storm surge does not pile up thus creating a lower water elevation of approximately 40 cm as compared to the original and core grid run as mentioned in the section above.

SUMMARY OF RESULTS AND FUTURE WORK
The objective of the study was to investigate the effects of landform changes on storm surge while developing workflows that integrate geospatial tools with process based numerical models in order to improve the simulation of complex landform changes.Geospatial processing tools provide methods to efficiently incorporate updated topography and bathymetry into existing grids to create real and hypothetical case study scenarios.Predictions based on model runs depend on the grid accuracy and the grid attributes.It is very important that elevation data as well as associated attribute data (e.g.vegetation cover) is frequently collected at the coastal regions for this kind of geospatial tools to prove useful to help solving the coastal zone problems.Also, the grid density plays an important role on the model output accuracy.Future work, include modification of the ADCIRC grid density in the vicinity of the channels and the dune ridges in order to capture the spatial resolution and detail of these spatial features.
Clearly the process of dune erosion, which is not modeled in this ADCIRC simulation, would be expected as the storm surge rises allowing the waves to reach the dune face.The wave action on the dunes would serve to erode the dune and flatten the profile providing the antecedent topography for breaching.Multiple ADCIRC runs, with changing topography (simulating change in landforms) based on our geospatial analysis illustrated above, enabled an analysis of storm surge as a function of changing topography as well as a study of the post storm tidal hydrodynamics that develop as a function of the breach.The results show that the changing topography does not notably affect the storm surge levels outside the area of dune collapse.On the other hand, the change in current velocities is remarkable, particularly for the core grid run.Even though a full breach is not implemented for this case, current flow patterns start to focus at the breach location before the maximum surge hits the location.Further studies can be carried out to look for these current patterns around the barrier island systems in order to characterize vulnerable locations and develop empirically based dune erosion models.Future work will include comparison of this suggested empirical approach to more robust dune erosion models such as XBeach.

Figure 2 .
Figure 2. The core surface is illustrated as the red layer beneath the higher elevation green surfaces at three transects.The cross section in the area that was breached during the Hurricane Isabel 2003 has no core surface as developed from the subaerial LIDAR (adapted from Mitasova et al. 2010).

Figure 3
Figure 3. a) Bathymetry survey lines by USACE b) Breach DEM interpolated and integrated in the original DSM c) ADCIRC grid nodes overlaid on the DSM d) Modified ADCIRC grid with breach topography and bathymetry.

Figure 5 .
Figure 5. Storm surge distribution at the selected nodes (2pm EST) a) 3D View b) Top View c) Along Shore Section d) Cross shore section.

Figure 6 .
Figure 6.Storm surge distribution at the selected nodes (5pm EST) a) 3D View b) Top View c) Along Shore Section d) Cross shore section.

Figure 7 .
Figure 7. Velocity comparison at the breach location and Hatteras inlet at 5pm, September 18 th ,2003.a) Original grid b) Core grid c) Breach grid.