Gridding of the NOS and NECOR soundings is accomplished using a combination of Generic Mapping Tool algorithms (Wessel and Smith, 1995), UNIX commands, and applications developed at NGDC specifically for merging gridded bathymetry with the USGS DEMs. The applications are linked together in a single Perl script command so that large volumes of digital bathymetric data can be gridded with minimal assistance from an operator. The script command is also versatile enough that it can be used to update individual grids as new data becomes available.
Each grid produced by the script encompasses an area of one degree of latitude by one degree of longitude. The horizontal resolution of the grid cells is identical to that of the USGS/NIMA DEMs; 3-arc-seconds, which is nominally 90 meters. The vertical resolution of the grid cells, on the other hand, is bimodal. Grid elevations of the onshore areas of the coastal zone are taken directly from the USGS/NIMA DEMs, which are resolved to one meter. But in the offshore areas, the grid elevations are resolved to one tenth of a meter.
This higher resolution was chosen at the urging of present users of the NOS hydrographic database who have found that the higher level of detail is not only supported by the sounding data, but reveals valuable morphologic information in near-shore areas and estuaries. To allow for these two vertical resolutions, the unit of elevation in the grids is one-tenth of a meter. Figure 1 is a flow chart of the procedure used to produce the grids. It diagrams the inputs to the Perl script and the various stages of the gridding process carried out by the script. These inputs and gridding stages are described in the sections that follow.
Gridding Inputs
Five input datasets must be assembled for each region of the US coastal zone to be gridded (Fig. 1). The first of these inputs is the compilation of all USGS/NIMA DEMs within the region. The USGS DEMs are acquired remotely via ftp from edcftp.cr.usgs.gov. The DEMs are then converted from their USGS ASCII format into a binary raster grid using a application called GS3 (written at NGDC). The final step is to convert this binary raster grid into a GMT grid using the GMT command grdraster. The second input dataset is the NOS medium-resolution (1:80,000) vector shoreline. The vector shoreline is referenced within a 3-arc-second resolution grid of the 1° x 1° region. The shoreline is then rasterized by flagging the grid cells that encompass the vector shoreline points as well as those cells that fall along straight lines connecting the points.
The third input dataset is an ASCII list containing the longitudes and latitudes of an arbitrary number of points that lie within the offshore area of the region to be gridded. These points are "seed" locations, which are used along with the rasterized shoreline to create a binary grid in which land areas have one value and ocean areas another. This land-sea mask and its role in the gridding process are described in detail below.
The fourth input dataset is a compilation of all the digital bathymetric soundings collected within the region to be gridded. These soundings come from three sources of data: (i) the NOS Hydrographic Soundings CD/DVD, (ii) the NGDC Multibeam database and (iii) recently digitized NOS soundings that have not yet been released on CD/DVD. The NOS Hydrographic Database CD/DVD does contain NOS multibeam data, but only five of the sixteen possible beams in a survey are included in this dataset. The Multibeam database contains not only the full resolution (i.e., all 16 beams) NOS multibeam surveys, but also the full resolution NECOR surveys archived at NGDC. As a result, all the multibeam bathymetry in the gridding process comes from this second database.
Soundings are extracted from all three datasets using the Search and Download software components of the NGDC GEODAS system (Sharman et al., 1998). Prior to gridding, as surveys are assimilated into the NOS Hydrographic Database, the soundings in each survey are manually checked against the soundings of adjacent and overlapping surveys to ensure depth consistency. In the future, the soundings will also all be corrected to the same vertical and horizontal datum (this issue was not addressed in constructing this first version of the gridded database). Currently, the soundings used to create the East Coast grids are referenced to two vertical datums and various horizontal datums. The vertical datums are mean low water (89% of the surveys) and mean lower low water (11%) of the surveys). The primary horizontal datums used were the NAD27 ellipsoid for soundings collected up until 1987-88, and the NAD83 ellipsoid for soundings collected since then. Despite their lack of uniformity, the different horizontal and vertical datums do not significantly alter the accuracy of the East Coast grids. This is because there is little difference in elevation between mean low water and mean lower low water in the gridding region, and the horizontal differences between the ellipsoids is less than the horizontal resolution of the cells within the grids (i.e.,<90m).
After being quality checked, all of the soundings are entered into a single, master sounding file. This file is an ASCII list of sounding longitude, latitude and depth. The fifth and final dataset is a list of all the 1° x 1° areas of longitude and latitude within the gridding region that contain bathymetric data.
Gridding Methodology and Outputs
The actual gridding is begun by entering into the Perl script command the paths of the directories to the five datasets above and the boundaries of the gridding region (Fig. 1). The first input dataset utilized by the Perl script is the list of 1° x 1° areas that contain bathymetric data. Using this list as a reference, every 1° x 1° area in the gridding region is classified as one of three grid types: land topography only, seafloor bathymetry only, or topography and bathymetry. Each area then undergoes a particular processing sequence based on its classification.
The simplest sequence is that for constructing grids of the areas of only land topography. In fact, the topography only grids are simply the USGS/NIMA DEMs converted to an internal format.
Slightly more complex is the sequence used to construct the seafloor bathymetry-only grids. In this sequence, two grids for each 1° x 1° area are produced. The first of these is the bathymetry-only elevation grid. A Perl script called trimXYZ is used to cull from the master sounding file all the soundings within an initial area of 1.2° x 1.2° (e.g., 67.9° - 69.1° W lon. by 37.9° - 39.1° N lat.), which is centered over the 1° x 1° area (68° - 69° W lon. by 38° - 39° N lat.). The culled soundings are input into the GMT application surface (Smith and Wessel, 1990) to create a grid of this larger. The final 1° x 1° grid is then extracted from this larger grid using the GMT application grdcut. The purpose of this process of gridding a larger area and then trimming the grid down to size is to ensure that the elevations along the boundaries of adjacent grids are identical. (See Fig. 2A).
The second grid produced for the bathymetry-only areas is a sounding-density grid (Fig. 2B). This grid has the same dimensions as the bathymetry grid, but instead contains the number of soundings in each cell of the bathymetry grid. Thus, the sounding-density grid provides the locations, spatial coverage and number of soundings upon which the depths of the bathymetry grids are derived. It also serves as an excellent reference for planning future bathymetric surveys, both to collect new soundings where data is lacking, and to document potential changes in areas where soundings have already been collected.
The sounding-density grid is created using the NGDC application grdDataLoc, which reads in all of the soundings within the data file created by trimXYZ, determines the grid cell within which each sounding falls, and then increments the value of that grid cell. Hence, cells containing for example three soundings have a value of three, while those with no soundings have a value of zero.
An additional step to constructing the bathymetry-only grids is performed for those grids that encompass the seaward limit of the hydrographic/bathymetric surveys. To prevent meaningless extrapolation of bathymetry into unsurveyed waters, an arbitrary limit is placed on the number of grid cells that the gridding process projects depths beyond a sounding location. Based on a trial-and-error evaluation of the maximum distance of extrapolation that preserves morphologic trends in the sounding data, a value of 110 grid cells was chosen.
The gridding sequence for areas encompassing both land and sea is the most complex. Again, both an elevation grid and a sounding-density grid are produced for each 1° x 1° area. The basic procedure for creating these grids is as follows:
1) The rasterized shoreline, the ocean "seed" locations, and the 1° x 1° area are input into the NGDC application grdSeaMask to generate a binary grid of land and ocean areas. The application creates the grid by first flagging the grid cells containing the seed locations as ocean area and those containing the rasterized shoreline as land area. The ocean area is then mapped out by flagging all the grid cells that are contiguous with the seed-location cells, but which are not part of the rasterized shoreline. When the mapping is complete, all the grid cells that have not been flagged as ocean area are flagged as land area, and the land-sea mask is complete.
2) The land-sea mask and the USGS DEM in the 1° x 1° land-topography-only grid area are then used in a boolean operation in the GMT application grdmath to create a third, land-topography/masked-bathymetry grid. In this third grid, the land areas in the land-sea mask are replaced with the DEM elevations, while the ocean areas remain identified by a single flag value. The purpose of this step is to restrict the seaward limit of the DEM to the rasterized NOS shoreline.
3) The land-topography-masked-bathymetry grid is converted into a list of longitude, latitude and topographic elevation using the GMT application grd2xyz. This list is then appended to the list of bathymetric soundings for the area extracted by trimXYZ. The result is a file containing all of the elevations for topography and bathymetry within a 1.2° x 1.2° region centered over the 1° x 1° area.
4) The GMT application surface is then run on the comprehensive elevation file to create an elevation grid of the 1.2° x 1.2° region. By using both land and seafloor elevations in the gridding algorithm morphological and/or structural trends are preserved across the shoreline.
5) The 1° x 1° grid is then extracted from the 1.2° x 1.2° grid using grdcut.
6) Finally, the land elevations in the 1° x 1° grid are replaced with the original USGS DEM values using a boolean operation in grdmath. This last step is taken to correct for slight alterations to the USGS DEM elevations by the gridding process and thus ensure that the land elevations remain the same.
7) The sounding-density grids for the combination topography/bathymetry grids are created in the same manner as those for the bathymetry-only grids with the exception that the land area in the grids are assigned a value of negative one.
The final sounding-density and elevation grids are stored in an internal format designed to minimize the size of the grids and make them compatible for use with viewing software described below. The NGDC grid format consists of a 128-byte header followed by the grid values, which are stored in rows that proceed from left to right and which are arranged from top to bottom (i.e., the grid origin is the upper left corner). The header contains 32 4-byte descriptors of the grid, which occur in the following order: version number, header length, data type (elevation grid, data-density grid, etc.), degrees of northernmost latitude, minutes of northernmost latitude, seconds of northernmost latitude, latitude dimension of grid cell (in arc-seconds), number of cells per grid row, degrees of westernmost longitude, minutes of westernmost longitude, seconds of westernmost longitude, longitude dimension of grid cell (in arc-seconds), number of cells per grid column, minimum grid elevation, maximum grid elevation, grid radius, elevation precision (meters or tenths of meters), empty grid-cell value (e.g., NaN), value type (how the data are stored, i.e. floating point, 2-byte integer, etc.), water datum (mean sea level or local datum), data limit ( maximum calculated value in grid cell), and 11 unused fields.