# Using GMT for 2D and 3D Modeling of the Ryukyu Trench Topography, Pacific Ocean

Polina Lemenkova

► To cite this version:

Polina Lemenkova. Using GMT for 2D and 3D Modeling of the Ryukyu Trench Topography, Pacific Ocean. Miscellanea Geographica - Regional Studies on Development, 2021, 25 (4), pp.213-225. 10.2478/mgrsd-2020-0038 . hal-03020902v2

HAL Id: hal-03020902

https://hal.science/hal-03020902v2

Submitted on 24 Mar 2022

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Using GMT for 2D and 3D Modeling of the Ryukyu Trench Topography, Pacific Ocean

**Abstract**

This research focuses on the 2D and 3D geospatial analysis of the Ryukyu Trench, a deep-sea trench located in the western Pacific Ocean between Japan and Taiwan. The aim of the research is to visualize regional differences in the topography of the southern (S) and northern (N) parts of the trench. Technically, the methodology is based on using the Generic Mapping Tools (GMT) scripting toolset, for modelling the General Bathymetric Chart of the Oceans (GEBCO), and Earth Topography and Bathymetry dataset (ETOPO1) raster grids. The results demonstrated topographic differences in the two segments. The most frequent depths lie between -5,000 and -6,000 m. The N part has steeper gradient slopes and deeper bathymetry. Of the depth differences >-6,000 m, S has nine values with depths >-6,800 m while N shows 123 records (max -7,460 m). The submarine terraces of S have gentler slopes compared with the N segment. The technical approach presents GMT-based 2D and 3D cartographic modelling aimed at visualizing regional variations of the seafloor topography.

**Keywords**

Cartography • GMT • GEBCO • ETOPO1 • Ryukyu Trench • topography

© University of Warsaw – Faculty of Geography and Regional Studies

Polina Lemenkova <fnref n="1" />,<fnref n="2" />

Received: 11 December 2019

Accepted: 28 June 2020

**Introduction**

Seafloor topography is essential for a variety of applications, including theoretical and scientific purposes (modelling geological and oceanographic processes, ocean circulation), environmental studies (mapping marine habitats), practical issues (safety of navigation), as well as political and commercial purposes (mapping the limits of the continental shelf for national coastal bathymetry). Topographic seafloor mapping contributes towards ecological projects aimed at protecting the world's oceans. In addition to topographic data, geomorphic settings of the seafloor shape - for example, slope, aspect and curvature - are significant variables for modelling benthic species distribution (Wilson et al. 2007). Reliable seafloor maps provide the basis for habitat modelling to assist in environmental management. Since the world's oceans cover 71% of the Earth (Eakins and Sharman, 2010), cartographic modelling of the seafloor, based on high-resolution bathymetric data, also increases our understanding of the planet's shape.

Oceanic trenches are an important part of the ocean's seafloor since they serve as natural boundaries between the Earth's tectonic plates. Among other objects of the seafloor (such as continental shapes, slopes, mid-ocean ridges) oceanic trenches have unique and distinct geometric shapes (long, narrow depressions with irregular and often V-shaped cross-section profiles). At the same time, being the deepest parts of the ocean, the trenches are the least-studied topography due to the clear difficulty and high cost of direct observation. This necessarily places an emphasis on computer-based modelling of the seafloor topography, specifically for the deep-sea trenches.

The current paper presents the modelling of one of the oceanic trenches, the Ryukyu Trench, located in the western Pacific Ocean and stretching between the island of Taiwan and Japan (Kyushu Island). The importance of seafloor mapping, briefly described above, raises the question of which GIS software should be used for cartographic visualization and modelling. Many types of GIS software have been used for various applications in the geosciences - for example, QGIS as an example of open-source GIS for environmental and geological mapping (Nielsen et al. 2017; Molina-Navarro et al. 2018; Meyer & Riechert 2019). An example of commercial software includes ArcGIS in geological mapping (Suetova et al. 2005; Harris et al. 2014; Klaučo et al. 2013; Kuhn et al. 2006). GIS has functionalities for geological modelling through its standard menu (manipulating layers, manual digitizing, overlay, creating 2D and 3D models). However, the problems of automated digitizing, 2D and 3D modelling and shell scripting in cartography remain a challenge. Generic Mapping Tools (GMT) is an open-source command-line mapping toolset developed by Wessel & Smith (2006). It is a powerful cartographic toolset that proposes a solution for 2D and 3D visualization based on raster grids, vector layers and the processing of big data arrays (ASCII data, tables). In contrast with traditional GIS, GMT consists of modules and presents a suite of cartographic tools, based fully on a shell scripting console. Due to its module-based structure, GMT exceeds the mapping capabilities of classic traditional GIS through its applicability with UNIX utilities (Lemenkova 2019a), full support of bathymetric mapping (Gauger et al. 2007) and automated rapid modelling based on shell scripting (Lemenkova 2019b, Becker 2005; Lemenkova 2019c). Examples of the automated data mapping

<footnote n="1">¹College of Marine Geo-sciences, Ocean University of China, Qingdao, People's Republic of China</footnote>

<footnote n="2">²Analytical Center, Moscow, Russian Federation
e-mail: pauline.lemenkova@gmail.com</footnote>
exist (e.g., DeVasto et al. 2012; Schenke & Lemenkova 2008). GMT proposed further steps in machine learning methods of geodata processing with 2D and 3D outputs as print-quality maps, figures and graphics.

This research presents the application of GMT with the aim of analyzing differences in the topography of the northern (N) and southern (S) parts of the Ryukyu Trench based on high-quality datasets: the General Bathymetric Chart of the Oceans (GEBCO) and the 1 arc-minute global relief model of Earth's surface (ETOPO1). Data analysis in marine geology may be processed by means of a combination of GIS with statistical data modelling - for example, QGIS with Python and R (Lemenkova 2019d; 2019e). However, the advantages of GMT are the combination of both cartographic visualization and statistical data assessment. 2D and 3D spatial modelling and mapping were performed using the available documentation and examples (Wessel & Smith 2018; Lemenkova 2019f) and a case study of the N and S segments of the Ryukyu Trench.

## Research focus

As briefly described above, the conceptualization of the seafloor mapping has been the subject of geographical research since the beginning of the GEBCO project in 1903 and has been actively developed thereafter (Smith & Sandwell 1997; Hall 2006; Heezen 1960). The rapid development of GIS tools in cartography has significantly assisted the development of precise oceanic mapping of such unreachable parts of the ocean as deep-sea trenches. Nevertheless, the combination of both GIS aspects and geographical concepts, specifically in the deep-sea trenches of the Pacific Ocean, has not yet been qualitatively explored in the existing literature, although some thematic overlaps can be found in works on the Ryukyu Trench (Ando 2012; Nishimura 2014). This paper aims to fill this gap and introduces results of both the technical application of GMT and a geographical regional overview of the Ryukyu Trench. The key questions addressed in this article are:

1. How can GMT technical tools be used for the mapping and 2D and 3D modelling of the ocean trench?

2. What are the topographic variations in spatially distinct parts of the trench?

3. How do these topographic variations correspond to the regional geological and geophysical settings?

Using raster and vector layers (GEBCO, ETOPO1, geology, gravity and geoid grids), this article attempts to perform a multi-source data analysis aimed at 2D and 3D regional topographic modelling. From the geographic point of view, the paper focuses on the regional analysis of the Ryukyu Trench: geological factors of its topographic variations in space and their correlations, as well as interaction between its bathymetry and movements of the tectonic plates. Technically, the paper uses GMT methodology, initially developed for terrestrial areas (Wessel & Smith 2006) and applied to submarine areas in this paper.

## Study area

### Geographic settings

The Ryukyu Trench stretches along the southeastern edge of the Ryukyu Islands Arc between the Philippine Sea, Taiwan, the Japanese island of Kyushu and the Pacific Ocean with a length of 1,400 km (Ando et al. 2018) (Figure 1). The system of the Ryukyu Trench and Ryukyu Island Arc serves as a natural border between two seas with highly contrasting bathymetric settings: the East China Sea (ECS) and the Philippine Sea. The ECS differs from other marginal seas by its highly developed shelf, which occupies most of its seafloor and extends to the

Yellow Sea. The shallow depths of the ECS clearly contrast with the depths of the Philippine Sea basin (Figure 1). The Ryukyu Island Arc is located on a submarine ridge with depth ranges of between 500 and 1000 m. The slopes of the ridge have a convex geomorphic profile and a stepped structure, formed by ancient folded areas bounded by Cenozoic and modern geosynclinal strata.

The Ryukyu Trench, located between the two seas with such contrasting bathymetry, extends along the foot of the outer slope of the Ryukyu Island Arc and reaches a maximal depth of 7,790 m (Litvin 1987).

### Bathymetry

The bathymetry of the Ryukyu Trench varies in its N and S segments, reflecting the correlation between its geomorphic shape, gradient steepness and local tectonic settings (Lemenkova 2019g). In its N part, the Philippine Sea Plate subduction is notable for shallow depths and low levels of steepness, reaching only 11° approximately. The Ryukyu Trench reaches depths of 7,000 m below sea level (BSL) (Figure 2). The N part of the Ryukyu Island Arc is divided into three parallel ridges. The forearc slope between the forearc island chain and the trench is between 130 and 160 km wide. The N group has Yaku-shima and Tanega-shima and the 200-km-long Tokara Channel gap between the island groups (Okamura et al., 2017). The depression of the Okinawa Trough is located north of the Ryukyu Island Arc (Figure 2). The continental crust decreases under the Okinawa Trough and increases again under the Ryukyu Island Arc. The basement surface of the Paleo-Mesozoic and Cenozoic geologic structures in the ECS is complicated by several sub-meridional depressions (with depths of >2 km) separated by uplifts (with depths of <1 km).

### Tectonics

The Ryukyu Trench was formed as a result of the subduction of the oceanic crust of the Philippine Sea Plate beneath the continental crust of the Eurasian Plate (Doo et al., 2018). The tectonic settings of the Ryukyu Trench reflect the evolutional processes of the building seamounts and structure of the Earth's crust. Therefore, they caused variations in the lithological patterns, the geological conditions of the rock formation, the stratigraphy and sedimentary settings of the ocean seafloor during the Jurassic, Mesozoic and Cretaceous (Lemenkova 2019h). The double subduction system of the Ryukyu Trench reflects complex trench motions as the expression of the interaction between the subducting plates and the underlying mantle of the Ryukyu Trench. The Pacific Plate subducts beneath the Philippine Sea Plate, which in turn subducts beneath the Eurasian plate along the Ryukyu Trench (Faccenna et al. 2018). The subduction of the Philippine Sea Plate causes the eruption of approximately 34 volcanoes (Figure 2) and earthquakes that produce tsunami.

Therefore, the Ryukyu arc-trench system is marked by active submarine margins of the continents located near northern Taiwan in the southern Ryukyu forearc (Figure 2). Southwestern Ryukyu is characterized by the edge regime influenced by the juxtaposing continental lithosphere (Kuo et al. 2012). Variations in the geomorphology of the trench follow the geometry of the slab dip which differs regionally: gentler in central (C) and southern (S) segments (40-50°) at a 7 km depth and steeper (70°) at a >7 km depth in the northern ones (Kodaira et al. 1996).

### Sedimentation

The sediment types of the Ryukyu Trench include mud diapirs, hemipelagic sediments, mud-supported breccia, fossil
<figure><img src="image_2.png" /><figcaption>Figure 1. Bathymetric map of the East China Sea, northern Philippine Sea, Ryukyu Arc and Ryukyu Trench. Bathymetry: GEBCO global terrain model, mapped using GMT</figcaption></figure>

Source: author

planktonic foraminifera, and calcareous nannofossils (Ujjie
2000). The sediment composition of the Ryukyu Trench changes
spatially along the Ryukyu arc from smaller in the N segment
arc to heavier along the C and S parts. This is explained by the
dominance of the terrigenous and pelagic sedimentation along
the N and C Ryukyu arc (Shu et al. 2017). The deep-sea trench
presents a kind of structural trap in the continental margins of
the Pacific, as sediments are carried by the ocean currents
passing through the system of trenches of the Pacific Ocean
(Lemenkova 2018). At the same time, the sediment thickness
and composition of the trench seafloor strongly depend on its
location and curvature (slope steepness). According to the data
obtained from deep-sea drilling (Litvin 1987), the sedimentation
structure of the northwestern section of the Philippine Sea is
composed of Eocene-Pleistocene low carbonate silts and clays,
with volcanic tuffs, conglomerates, and basalts beneath them.
The carbonate content of precipitation increases markedly on
the geomorphic positive relief forms (small seamounts and
hills).

Environmental settings

Although the seafloor of the Ryukyu Trench is separated
vertically from the euphotic zone, it serves as a depositional
center for organic carbon accumulation. This is closely connected
with geomorphic factors: steep slopes enable lateral transport
from the continental margins (Itoh et al. 2011). As a consequence,
the spatial biodiversity of the Ryukyu Trench has lower values
at shallow depths, but these increase with water depth reaching
maximum levels between the trench floor and the abyssal plain
(Kitahashi et al. 2014). This proves that the distribution of the deep-
sea biomass of the Ryukyu Trench spatially differs in various
regions of the trench following geomorphic variations. Such
discrepancies are a result of differences in environmental factors:
primary productivity level. Moreover, migration patterns that have
been discovered in the Ryukyu Trench area show that reefs and
geomorphic landforms serve as land bridges for the migration
of terrestrial organisms from Okinawa-jima to Miyako-jima (Arai
et al. 2016).
<figure><img src="image_3.png" /><figcaption>**Figure 2.** Geologic map of the Ryukyu Trench. Bathymetric base map: ETOPO1 global terrain model, mapped using GMT
Source: author</figcaption></figure>

Methodology

Data acquisition

The dataset used in the present study includes the GEBCO global terrain model downloaded from the web (GEBCO Committee 2016). The GEBCO project was originally initialized by HSH Prince Albert I of Monaco in 1903 to enable studies of the physical shape, geomorphology and bathymetric characteristics of the seafloor (Weatherall et al. 2015). Since the initiation of the GEBCO project at the onset of the 20th century, it has been significantly improved due to the technical and cartographical

tools that have been developed, while the goal of mapping the
bathymetry of the world's oceans has been maintained until
recent times (Mayer 2006). The GEBCO project is focused on the
high-precision mapping of the ocean seafloor and is developed
and maintained by several oceanographic organizations
(Weatherall et al. 2015). The GEBCO_2019 grid used in the current
research is methodologically based on the development of the
previous versions of GEBCO (Schmitt & Weatherall 2014) and the
SRTM15_PLUS dataset (Olson et al. 2014). Besides, GEBCO is
collaborating with regional groups to continually improve its global
<figure><img src="image_4.png" /><figcaption>**Figure 3.** Model of the geoid in the Ryukyu Trench area. Raster data: Earth Gravitational Model EGM96, mapped by GMT
Source: author</figcaption></figure>

model: IBCAO and IBCSO gridded datasets (Jones 1994; Schmitt & Weatherall 2013; IOC, IHO & BODC 2003; GEBCO_2010). An ETOPO1 15 arc-sec grid was downloaded from the publicly available NOAA resources (Amante & Eakins 2009).

Geoid and marine free-air gravity field grids were based on Geosat and ERS-1 altimetry raster grids (Smith & Sandwell 1995; Sandwell et al. 2014). The geographic names were verified in the GEBCO Gazetteer of submarine features (IHO-IOC 2012). The vector lines on the maps (tectonic slabs, trench lines, coastlines, tectonic faults) were taken from the embedded layers of GMT (Wessel & Smith 1996). The data were converted to several projections: Mercator (Figure 1, Figure 6, C), Lambert Conic Conformal (Figure 2), Equidistant Conic (Figure 3), Cylindrical Equal-Area Gall-Peters (Figure 4), and a 3D model (Figure 5). A subset covering the area with geographic coordinates of 120°00′E 134°00′E 20°00′N 33°00′N was used for the study.

## Data processing

GMT is an open-source set of over 50 command-line cartographic commands suitable for manipulating big data both in 2D and 3D formats (arrays of XY and XYZ datasets in ASCII format). Shell scripts are executable on UNIX, Linux or Mac OS X. The GMT-integrated scripts provide a robust solution for generating high-accuracy maps based on GEBCO and ETOPO1 grids, extracting isolines, coastlines, relief, and plotting publication-quality maps. Eight GMT commands - integrated shell scripts - were written for plotting maps and graphics. Several GMT modules were used for cartographic visualizations.

For instance, 'grdcontour' was used to extract isolines from a GEBCO grid; 'psconvert', 'psscale' to visualize scale bars; 'pstext' to add annotations; 'psxy' to plot geologic lineaments and sample points; 'psbasemap' to map cartographic elements (isolines, geoid, bathymetry, net grids and titles), and the module
<figure><img src="image_5.png" /><figcaption>Figure 4. Model of the marine free-air gravity anomaly showing gravitational settings in the Ryukyu trench-arc area. Data: gravity grid, visualized by GMT</figcaption></figure>

Source: author

'grdtrack' to automatically extract data for plotting cross-section profiles from a GEBCO grid.

## GMT shell scripting

Each GMT module consists of a small line of code. Combined as a script, they were used to produce maps. The 'grdimage' module was used to map Figures 1, 3 and 4 (bathymetric grid, geoid, and gravity, respectively) using a sequence of the GMT codes. The 3D modelling was done using the 'grdview' module with the command: `grdview rt_relief.nc -J -R -JZ2i -CrainbowRT.cpt -p115/30 -Qsm -N-7500+glightgray -Wm0.07p -Wf0.1p,red -Wc0.1p,magenta -B4/4/2000:"Bathymetry and topography (m)":ESwZ -S5 -Y6.0c -O -K >> $ps`. The 3D map was overlaid with the 2D geoid model (Figure 5) presenting a combined 2D-3D modelling overlay approach. The trend curves were modelled using the 'trend1d' module. Statistical modelling and digitizing cross-section profiles were done using the 'pshistogram' and 'grdtrack' modules.

## GMT module approach

The design of GMT follows the UNIX approach in which multi-tasks are broken down into small and independent 'self-contained' commands and small sub-tasks. Consequently, taken together as a chain of such sub-tasks, individual GMT module commands were combined into a shell script for batch processing. A workflow of the GMT shell scripts contains a sequence of modules; each of them was used for a specific task. Visualizing raster grids as maps was performed using the 'grdimage' GMT module with local arguments set up for plotting grids in specific map projections and coordinate systems. For example, a command `-R120/134/20/33 -JM6i` means selecting a square region with 120°E 134°E 20°N 33°N coordinates and plotting a map in Mercator projection with a polygraphic (printing) dimension of the map layout of 6 inches. Applying color palettes was done using GMT module 'makecpt', cleaning data was performed using 'rm' from the UNIX utility. A map legend was plotted using the 'pslegend' module with annotating texts using
<figure><img src="image_6.png" /><figcaption>Figure 5. Composite overlay of the 3D topographic mesh model on the 2D geoid contour plot for topographic modelling of the Ryukyu trench-arc area</figcaption></figure>

Source: author

the 'pstext' module, generating a study area through 'cutting' selected squares ('grdcut'). Delineating isolines and adding bathymetric contour lines used to represent elevation from DEM was performed using the 'grdcontour' module.

### 3D and 2D modelling

3D mapping allowed the profiling of the topographic objects (cross-sections) in three dimensions. The map was technically based on a 'grdview' module which uses an advanced cartographic method for surface visualization and gathering information from a 3D scene. Extracting cross-section profiles was done by 'grdtrack' which enables linear graphs of the topographic curve from the start and end points of the cross-sectional segment to be plotted. The curve of the selected segment on a profile shows elevations perpendicular to the surface. Cross-section profiles

show the gradient, slope curvature and bathymetry represented on a 2D model. The 3D modelling of the plot with a mesh surface was done using a 'grdview' module. For 3D modelling, bathymetric data were exported as a NetCDF (Network Common Data Form) file with an .nc extension (@rt_relief.nc) and imported into GMT for raster processing, 3D visualization and evaluation. Bathymetric data were gridded at a 1 arc-minute resolution from the original ETOPO1 file (earth_relief_01m.grd).

### Bathymetric uncertainty

Error bars represent uncertainty in bathymetric data and the statistical range of data distribution. Bathymetric uncertainty can be also analyzed using QGIS plugins, a combination of R and Python, to implement methods and tools aimed at generating an estimate of the uncertainty (Schmitt et al., 2015). The estimation
of uncertainty depends on the methodology (in particular, on its accuracy) of the interpolation algorithms adapted for bathymetric use. For instance, interpolation algorithms of topographic isolines include inverse distance weighted interpolation (IDW), triangulated irregular network (TIN), B-Spline and Spline under tension (Charpy et al., 2015). Error bars (thin gray vertical lines on Figure 6, A and B) accompanying median values of the cross-section profile (red line on Figure 6, A and B) graphically present the variability of the bathymetric data used to indicate error or uncertainty in the input raster grid. A collection of GMT modules was used to perform mapping, projecting, rasterization, plotting isolines using techniques of interpolation and smoothing of raster grids (Wessel & Smith 2015).

The GMT 'psconvert' command was applied to convert PS files to JPEG (and can also be used for other desired extensions, such as EPS, PDF, and TIFF) through the following command: `gmt psconvert RTtrends2.ps -A20.0c -E720 -Tj -P -Z`, where 'RTtrends2.ps' is the name of the file, and -E720 is a 720 dpi resolution file in a JPG format (-Tj). Adding essential elements on the base map (title, legend, directional rose, grid, scale, time stamp) was technically possible using 'psbasemap'; calculating statistical mean and median, using 'pshistogram'. Error bars showing the standard deviation of the elevations were plotted using the `gmt grdtrack trenchRTn.txt -Grtrelief.nc -C400k/2k/20k -Sm+sstackRTn.txt > tableRTn.txt` command.

## Cartographic visualization

Each of the eight maps presented was created using eight GMT scripts, derived from Lemenkova (2019i) and modified to carry out this research. Within an individual GMT shell script, the attributes of each of the maps (coordinates of the map square, color palette, annotations, border ticks, legend and other elements) were strung together as a series of commands. Codes for text annotations were used as a template for iterations to repeat texts with modifications (coordinates, font size, color, rotation), as can be seen in Figure 2. The serious polygraphic advantage of GMT is its operation with output files. It creates print-quality outputs and graphics in a PostScript format, an industry standard for high-quality printing, which can be converted to other graphical formats (JPG, TIFF, BMP).

## Results

### 3D modelling

The geomorphology of the Ryukyu trench-arc seafloor results from the complex geodynamic settings and tectonic processes affecting its structure. These demonstrate traces of the continuous evolution and shape of the seabed surface as a system of mid-ocean ridges, trenches, troughs, transform faults, vast ocean plains complicated by ridges, plateaus, and such minor geomorphic forms as terraces, seamounts and hills (Figure 5). Perspective 3D scene visualization focused on the evaluation of such parameters as submarine topography and the comparison of the bathymetric ranges in different parts of the trench. Various factors affect variations in the submarine topography of the trench. Along with regional patterns, local variations are affected by the geological settings which include the character of the sedimentation processes. For example, the seafloor of the ECS is mainly covered by a thick continental crust with a thickness of 30 km approximately.

A raster layer was generated as a color-coded surface at the same resolution as for the visualized 3D model of the topographic slope, overlaying the geoid contour plot (Figure 5). Azimuthal rotation from a 115°/30° direction was added by rotating the initial 3D grid surface in an azimuth using the `grdview` module. Although GEBCO presents the highest quality for the global terrain model for ocean and land (15 arc-second intervals), the ETOPO1 raster grid

was selected as a base map for 3D modelling, since it has enough density to grid at an acceptable resolution (1 min). The ETOPO1 reduced the modelling and computation time for the 3D mesh and the limitations for data processing and visualization. The seafloor relief correlates with the core basement: the surface of the Ryukyu island arcs with steep ledges descends to the bottom of the Philippine Sea basins and Ryukyu Trench. Depths at the bottom of the seafloor generally exceed 5 km, in some places up to 6 km, and reach 7,790 m in the trench. Spatial differences of the Ryukyu Trench topography are mirrored in its sediment composition that changes systematically from lighter in the northern Ryukyu arc to heavier along the central part of the arc.

### 2D modelling

The comparison between two spatially different segments of the trench show that the N segment has higher depth values compared with S. This can also be noted in Figure 1 where the deepest values of the central part are clearly seen in the area extending seaward from the Miyako Fault (Figure 2). As can be seen when comparing Figures 6 A and B, the median stacked profile on the S part (Figure 6 A) demonstrates a slightly gentler slope shape with submarine terraces, which is not clear on the N part of the forearc slope (Figure 6 B). At the same time, the depth of the forearc slope increases in the N segment of the trench (Figure 6 B).

The geometric approximation of the trend curves (Figure 7) shows similarity between the segments in the applied coarse linear model and differences in geometric variations at applied polynomial functions.

### Statistical data processing

Statistical data processing (Figure 8) shows differences in the depth ranges of the Ryukyu Trench along the island arc. Depth and geometry of the slope vary considerably: maximal values of -5,500 to -5,700 m are notable for the N segment where 334 sample points were recorded (Figure 8 A).

The forearc basin at a depth of >2,000 m and <1,000 m BSL shows more samples in the N segment than in S (Figure 7). In contrast, a wide terrace with depths that are shallower than 25 m is generally to be found more frequently (161 samples) in the S part of the trench.

### Modelling marine free-air gravity and geoid

Variations of the slope geomorphology are highlighted by the contour of the increased marine free-air gravity anomaly, as a triangle shaped 'island' located at 126°E 24°N (Figure 4). Visualized gravity anomalies (Figure 4) show local variations of the gravity field around the Ryukyu arc-trench system and region of the East China and Philippine seas. Positive anomalies in the range of 40 to -240 mGal, (orange to red areas on Figure 4) demonstrate more gravity than predicted by the model which is explained by the presence of a sub-surface positive mass anomaly in the Ryukyu Island Arc and Taiwan Island. On the contrary, negative anomalies in the range of -20 to -240 mGal (colors blue to purple on Figure 4), clearly depicting the stripe of the Ryukyu Trench area showing lower values than predicted, which indicates a sub-surface mass deficit in the area.

In general, transition zones in the region from the ECS eastwards along the Philippine Sea towards the Mariana Islands are characterized by increased free-air gravity value (Faye's) anomalies of 0-25 mGal (green and light blue-greenish colors in Figure 4 which cover the majority of the area). A clear exception can be seen in the Ryukyu Trench: a purple-colored stripe with values below -42 mGal (Figure 4). Furthermore, Faye's anomalies on the shelf of the ECS are in the range of 20-30 mGal; in the Okinawa basin, they increase to 80-100; over the
<figure><img src="image_7.png" /><figcaption>Figure 6. Modelled stacked profiles for the southern (A) and northern (B) segments with the red solid line as a median, and error bars showing weighted standard deviation of the elevations; automatically digitized cross-section profiles in two segments (C), plotted in GMT</figcaption></figure>

Source: author
<figure><img src="image_8.png" /><figcaption>**Figure 7.** Modelling graphs of the trend curves of general geomorphological transects and slope gradients using various mathematic approximations: fitted regression models $y=f(x)+e$; Weighted Least Squares, polynomial and Fourier. Plotted by GMT
Source: author</figcaption></figure>

<figure><img src="image_9.png" /></figure>

island arc of Ryukyu, they decrease again to 20–40 mGal; and in the Philippine Sea basin, they reach 40–80 mGal (Figure 4). These results can be compared with previously measured Bouguer anomalies which are slightly positive on the shelf of the ECS, reach 160 mGal in the Okinawa Trough and 60–80 mGal over the Ryukyu Island Arc (Gaynanov 1980).

Being a joint structure, the Ryukyu Trench is located on the slopes of the depression of the surface of the upper mantle, at the junction of the subcontinental and oceanic crust types. The thickness of the Earth's crust in the Ryukyu Trench is 8.2 km (in comparison, it is 12.7 km in the Kuril-Kamchatka Trench). The thickness ratio of the sedimentary and basalt layers varies from 1:1 to 1:3 (Milashin 1971). The thickness of the basalt layer of the trench's seafloor exceeds the thickness of the same layer in the oceanic crust and is less than that of the subcontinental crust of anticline ridges. The mosaic distribution of crust types in such geosynclinal systems as oceanic trenches corresponds to the extreme values in gravity anomalies. Positive gravity anomalies (Bouguer reduction) are generally characterized by geosynclinal basins, which corresponds with the results presented. The belt of clearly expressed negative gravity anomalies along the Ryukyu Trench demonstrated corresponds with the general pattern of free-air gravity in deep-sea trenches (a rare exception is the Cayman Trench). The geoid model shows geoidal undulation in the SW part of the area east of Taiwan (Figure 3). Generally, a major part of the area has geoid values ranging from 25 to 35 mGal. Increased values (>35 to 85 mGal) can be seen on the SE part in the Philippine Sea, while lower values (15–25 mGal)

bordering the continental shelf stretch along the ECS, and the lowest (<15 m Gal) mainly occupy land areas.

## Discussion

In terms of cartographic output, using GEBCO and ETOPO1 raster grids takes full advantage of the application of such a professional toolset as GMT. The 3D modelling used for geological datasets enabled spatial differences in the submarine geomorphology to be highlighted using a three-dimensional scenic view. As mentioned before, the specifics of the submarine geomorphology are due to the remote location of the study object: the submerged forms of the Earth's relief can only be studied using modelling algorithms, remote sensing techniques, numerical approaches and computer-based software. Various factors affect submarine relief, which is reflected by complex structure and regional variations in the marine environment (Lemenkova 2019j): tectonic plate movements, magmatism, seismicity (Gutscher et al. 2016), volcanism and hydrothermal factors (Minami & Ohara 2018), thickness of the subduction channels, and earthquake rupture segments on the geomorphology of the oceanic landforms (Lemenkova 2019k). Understanding the complexity of these factors, their interrelationship and impact on the morphology of the deep-sea trench can be achieved by 2D and 3D modelling using the advanced cartographic solutions proposed by GMT.

## Conclusion

This paper demonstrated the cartographic application of the GMT modules used for 2D and 3D mapping for the topographic
<figure><img src="image_10.png" /><figcaption>Figure 8. Statistical histograms showing depth frequency in the cross-sections of the Ryukyu Trench, based on the bathymetric data: northern (A) and southern (B) parts</figcaption></figure>

Source: author

and geological dataset of the Ryukyu Trench. It aimed to analyze the topographic variations of the cross-section profiles of the oceanic trench in the context of its regional settings: geology, geomorphology, tectonics. Specifically, variations in the S and N segments of the trench were visualized. A multidisciplinary study combined cartographic visualization and geological analysis using high-resolution raster and vector datasets. The results obtained enabled topographic variations in the slope angle of the trench to be detected, showing steeper gradients and deeper bathymetric values in the N segment and shallower depths with gentler slopes in S.

Sequential use of shell scripting by GMT modules for thematic mapping and visualization of the raster grids and vector layers resulted in print-quality maps. The cartographic functionality of the GMT presented 2D maps (in four various cartographic projections), 3D models and cross-section graphs. This was implemented on the toolset level through the well-structured modular nature of GMT. The technical follow-up of this work presented snippets of the most important scripting codes.

The methodology demonstrated an approach for processing marine geological datasets using automated shell scripts comprising GMT and auxiliary UNIX utilities ('echo', 'rm', 'cat') for cartographic modelling combined with statistical analysis (modelled trends and histograms). Visualization of the thematic maps (geology, bathymetry, geoid, gravimetry, 3D modelling, 2D profile trends modelling, digitized cross-section transects) demonstrated the full technical functionality of GMT using 20 modules for mapping and plotting: gmtdefaults, grdcut, makecpt,

grdimage, grdcontour, psscale, psbasemap, psxy, grdtrack, convert, psconvert, pshistogram, psrose, set, pslegend, pscoast, grdimage, grd2cpt, pstext, and trend1d.

Effective cartographic visualization facilitates our better understanding of the seafloor geomorphology and contributes to the existing dataset on the Earth's topography. Advanced solutions of the GMT-integrated toolset make it possible to acquire accurate high-resolution 2D and 3D maps through scripting techniques. Using GEBCO and ETOPO1 topographic datasets for the generation of the bathymetric and geological maps is recommended as an effective data source for marine geological and bathymetric spatial-related studies.

GMT proved to be an effective cartographic instrument for 2D and 3D modelling, processing, visualization and producing print-quality polygraphic maps and graphics. As a special remark, it should be mentioned that GMT enabled big datasets and arrays (containing thousands of cells) to be processed in a rapid and precise way. This is especially important for geological studies which are mostly based on the processing of massive datasets obtained from field surveys. Needless to say, the production of high-quality press-ready maps is an important component of marine geological studies, due to difficult access to study areas such as deep-sea trenches. The increasing diversity of spatially distributed data on the observation of the Earth with advanced scripting cartographic methods enhances the potential of geospatial modelling using 2D and 3D approaches for regional studies of the submarine topography in selected areas of the World Ocean.
Acknowledgments

This research was funded by the China Scholarship Council, State Oceanic Administration, Marine Scholarship of China, Grant Nr. 2016SOA002. I cordially thank the two anonymous reviewers for their helpful comments on the manuscript.

ORCID

Polina Lemenkova https://orcid.org/0000-0002-5759-1089

References

Amante, C & Eakins, BW 2009, 'ETOPO1 1 arc-minute Global Relief Model: Procedures, Data Sources and Analysis', *NOAA Technical Memorandum NESDIS NGDC-24*.

Ando, M, Tu, Y, Kumagai, H, Yamanaka, Y & Lin, CH 2012, 'Very low frequency earthquakes along the Ryukyu subduction zone', *Geophysical Research Letters*, vol. 39.

Ando, M, Kitamura, A, Tu, Y, Ohashi, Y, Imai, T, Nakamura, M, Ikuta, R, Miyairi, Y, Yokoyama, Y & Shishikura, M 2018, 'Source of high tsunamis along the southernmost Ryukyu trench inferred from tsunami stratigraphy', *Tectonophysics*, vol. 722, pp. 265-276.

Arai, K, Matsuda, H, Sasaki, K, Machiyama, H, Yamaguchi, T, Inoue, T, Sato, T, Takayanagi, H & Iryu, Y 2016, 'A newly discovered submerged reef on the Miyako-Sone platform, Ryukyu Island Arc, Northwestern Pacific', *Marine Geology*, vol. 373, pp. 49-54.

Becker, NC 2005, 'Painting by numbers: A GMT primer for merging swath-mapping sonar data of different types and resolutions', *Computers & Geosciences*, vol. 31(8), pp. 1075-1077.

Charpy, C, Schmitt, T, Biscara, L, Maspataud, A, Avisse, L & Créach, R 2015, 'Précision et Performance des Méthodes d'Interpolation pour la Réalisation de Modèles Numériques de Terrain Bathymétriques' ['Precision and functionality of the interpolation methods for the realization of Digital Bathymetric Terrain Models'] in *Colloque merlGéo*, 24-26 novembre 2015, Brest, France.

DeVasto, MA, Czeck, DM & Bhattacharyy, P 2012, 'Using image analysis and ArcGIS® to improve automatic grain boundary detection and quantify geological images', *Computers & Geosciences*, vol. 49, pp. 38-45.

Doo, WB, Lo, CL, Wu, WN, Lin, JY, Hsu, SK, Huang, YS & Wang, HF 2018, 'Strength of plate coupling in the southern Ryukyu subduction zone', *Tectonophysics*, vol. 723, pp. 223-228.

Faccenna, C, Holt, AF, Becker, TW, Lallemand, S & Royden, LH 2018, 'Dynamics of the Ryukyu/Izu-Bonin-Marianas double subduction system', *Tectonophysics*, vol. 746, pp. 229-238.

Gauger, S, Kuhn, G, Gohl, K, Feigl, T, Lemenkova, P & Hillenbrand, C 2007 'Swath-bathymetric mapping'. *Reports on Polar and Marine Research*, vol. 557, pp. 38-45.

Gaynanov, AG 1980, *Gravimetric studies of the Earth's crust of the oceans*, MSU, Moscow.

GEBCO Committee 2016, *General Bathymetric Chart of the Oceans (GEBCO) - from the coast to the deepest trench*. Available from: <https://www.gebco.net/>. [21 April 2020].

GEBCO 2010, GEBCO_08 grid, version 20100927, *British Oceanographic Data Centre (BODC)*, Liverpool, UK.

Gutscher, MA, Klingelhoefer, F, Theunissen, T, Spakman, W, Berthet, T, Wang, TK & Lee, CS 2016, 'Thermal modeling of the SW Ryukyu forearc (Taiwan): Implications for the seismogenic zone and the age of the subducting Philippine Sea Plate (Huatung Basin)', *Tectonophysics*, vol. 692, pp. 131-142.

Hall, JK 2006, 'GEBCO Centennial Special Issue—Charting the secret world of the ocean floor: The GEBCO Project 1903–2003', *Marine Geophysical Research*, vol. 27, no. 1, pp. 1–5.

Harris, PT, Macmillan-Lawler, M, Rupp, J & Baker, EK 2014, 'Geomorphology of the oceans', Marine Geology, vol. 352, pp. 4–24.

Heezen, BC 1960, 'The rift in the ocean floor', Scientific American, vol. 203, no. 4, pp. 98–110.

IHO-IOC 2012, GEBCO Gazetteer of Undersea Feature Names. Available from: <www.gebco.net>. [21 April 2020].

IOC, IHO & BODC 2003, 'Centenary Edition of the GEBCO Digital Atlas, Published on CD-ROM on Behalf of the Intergovernmental Oceanographic Commission and the International Hydrographic Organization as Part of the General Bathymetric Chart of the Oceans', *British Oceanographic Data Centre*, Liverpool, UK.

Itoh, M, Kawamura, K, Kitahashi, T, Kojima, S, Katagiri, H & Shimanaga, M 2011, 'Bathymetric patterns of meiofaunal abundance and biomass associated with the Kuril and Ryukyu trenches, western North Pacific Ocean', *Deep-Sea Research Part I*, vol. 58, pp. 86-97.

Jones, MT 1994, *The GEBCO Digital Atlas*, pp. 17–20, NERC News, Swindon, UK.

Kitahashi, T, Kawamura, K, Kojima, S & Shimanaga, M 2014, 'Bathymetric patterns of a and b diversity of harpacticoid copepods at the genus level around the Ryukyu Trench, and turnover diversity between trenches around Japan', *Progress in Oceanography*, vol. 123, pp. 54–63.

Klaučo, M, Gregorová, B, Stankov, U, Marković, V & Lemenkova, P 2013, 'Determination of ecological significance based on geostatistical assessment: a case study from the Slovak Natura 2000 protected area', *Central European Journal of Geosciences*, vol. 5(1), pp. 28–42.

Kodaira, S, Iwasaki, T, Urabe, T, Kanazawa, T, Egloff, F, Makris, J & Shimamura, H 1996, 'Crustal structure across the middle Ryukyu trench obtained from ocean bottom seismographic data', *Tectonophysics*, vol. 263, pp. 39–60.

Kuhn, G, Hass, C, Kober, M, Petitat, M, Feigl, T, Hillenbrand, CD, Kruger, S, Forwick, M, Gauger, S & Lemenkova, P 2006, 'The response of quaternary climatic cycles in the South-East Pacific: development of the opal belt and dynamics behavior of the West Antarctic ice sheet', *Expeditions program Nr. 75 ANT XXIII/4, AWI Germany*.

Kuo, BY, Wang, CC, Lin, SC, Lin, CR, Chen, PC, Jang, JP & Chang, HK 2012, 'Shear-wave splitting at the edge of the Ryukyu subduction zone'. *Earth and Planetary Science Letters*, vol. 355–356, pp. 262–270.

Lemenkova, P 2018, 'R scripting libraries for comparative analysis of the correlation methods to identify factors affecting Mariana Trench formation'. *Journal of Marine Technology and Environment*, vol. 2, pp. 35–42.

Lemenkova, P 2019a, 'AWK and GNU octave programming languages integrated with generic mapping tools for geomorphological analysis', *GeoScience Engineering*, vol. 65 (4), pp. 1–22.

Lemenkova, P 2019b, 'Geomorphological modelling and mapping of the Peru-Chile Trench by GMT', *Polish Cartographical Review*, vol. 51(4), pp. 181–194.
Lemenkova, P 2019c, 'Automatic data processing for Visualising Yap and Palau Trenches by Generic Mapping Tools', *Cartographic Letters*, vol. 27(2), pp. 72–89.

Lemenkova, P 2019d, 'An empirical study of R applications for data analysis in Marine Geology'. *Marine Science and Technology Bulletin*, vol. 8, no. 1, pp. 1–9.

Lemenkova, P 2019e, 'Testing linear regressions by StatsModel library of Python for oceanological data interpretation', *Aquatic Sciences and Engineering*, vol. 34, pp. 51–60.

Lemenkova, P 2019f, 'Deep-Sea trenches of the Pacific Ocean: a comparative analysis of the submarine geomorphology by data modeling using GMT, QGIS, Python and R. Mid-Term PhD Thesis Presentation: Current Research Progress', *Presentation at OUC, College of Marine Geo-sciences*, Qingdao, China.

Lemenkova, P 2019g, 'Statistical analysis of the Mariana Trench geomorphology using R programming language', *Geodesy and Cartography*, vol. 45, no. 2, pp. 57–84.

Lemenkova, P 2019h, 'Numerical data modelling and classification in Marine geology by the SPSS statistics'. *International Journal of Engineering Technologies*, vol. 5, no. 2, pp. 90–99.

Lemenkova, P 2019i 'Topographic surface modelling using raster grid datasets by GMT: example of the Kuril-Kamchatka Trench, Pacific Ocean', *Reports on Geodesy and Geoinformatics*, vol. 108, pp. 9–22.

Lemenkova, P 2019j. 'GMT based comparative analysis and geomorphological mapping of the Kermadec and Tonga Trenches, Southwest Pacific Ocean'. *Geographia Technica*, vol. 14, no. 2, pp. 39–48.

Lemenkova, P 2019k, 'Processing oceanographic data by Python libraries NumPy, SciPy and Pandas', *Aquatic Research*, vol. 2, pp. 73–91.

Litvin, VM 1987, Morphostructure of the ocean floor, Nedra, Leningrad.

Mayer, LA 2006, 'Frontiers in seafloor mapping and visualization', *Marine Geophysical Research*, vol. 27, pp. 7–17.

Meyer, D, Riechert, M 2019, 'Open source QGIS toolkit for the Advanced Research WRF modelling system', *Environmental Modelling & Software*, vol. 112, pp. 166–178.

Milashin, AP 1971 'On the differences in the structure of the earth's crust of the seas and oceans' in *Marine Geology and Geophysics*, ed. AP Milashin, Nedra, Moscow, pp. 3–16.

Minami, H & Ohara, Y 2018, 'Detailed volcanic morphology of Daisan-Miyako Knoll in the southern T Ryukyu Arc', *Marine Geology*, vol. 404, pp. 97–110.

Molina-Navarro, E, Nielsen, A, Trolle, D 2018, 'A QGIS plugin to tailor SWAT watershed delineations to lake and reservoir waterbodies', *Environmental Modelling & Software*, vol. 108, pp. 67–71.

Nielsen, A, Bolding, K, Hu, F, Trolle, D 2017, 'An open source QGIS-based workflow for model application and experimentation with aquatic ecosystems', *Environmental Modelling & Software*, vol. 95, pp. 358–364.

Nishimura, T 2014, 'Short-term slow slip events along the Ryukyu Trench, southwestern Japan, observed by continuous GNSS', *Progress in Earth and Planetary Sciences*, vol. 1, no. 22.

Okamura, Y, Nishizawa, A, Oikawa, M & Horiuchi, D 2017, 'Differential subsidence of the forearc wedge of the Ryukyu (Nansei-Shoto) Arc caused by subduction of ridges on the Philippine Sea Plate', *Tectonophysics*, vol. 717, pp. 399–412.

Olson, CJ, Becker, JJ & Sandwell, DT 2014, 'A new global bathymetry map at 15 arcsecond resolution for resolving seafloor fabric: SRTM15_PLUS', in *Proceedings of the AGU Fall Meeting Abstracts 2014*, pp. 1-3, San Francisco, CA.

Pushkov, AN 1981, *Anomalies of the geomagnetic field and the deep structure of the Earth's crust*, Naukova Dumka, Kiev.

Sandwell, DT, Müller, RD, Smith, WHF, Garcia, E & Francis, R 2014, 'New global marine gravity model from CryoSat-2 and Jason-1 reveals buried tectonic structure', *Science*, vol. 346, issue 6205, pp. 65–67.

Schenke, HW & Lemenkova, P 2008, 'Zur Frage der Meeresboden-Kartographie: Die Nutzung von AutoTrace Digitizer für die Vektorisierung der Bathymetrischen Daten in der Petschora-See' ['To the question of seafloor mapping: the use of AutoTrace digitizer for the vectorization of bathymetric data in the Pechora Sea'], *Hydrographische Nachrichten*, vol. 81, pp. 16–21.

Schmitt, T & Weatherall, P 2014, 'GEBCO and EMODNet—Bathymetry hand in hand: Improving global and regional bathymetric models of European waters', *Abstract OS31B-0989* presented at 2014 Fall Meeting, AGU, San Francisco, California, US, 15–19 December.

Schmitt, T & Weatherall, P 2013, 'GEBCO and EMODnet Bathymetry hands in hands' in AGU Fall Meeting, San Francisco, US.

Schmitt, T, Penard, C & Waddle, J 2015, 'Uncertainty and bathymetric DEM—Developing an open source QGIS solution' in *Conference GEBCO Science day*, Kuala Lumpur, Malaysia.

Shu, Y, Nielsen, SG, Zeng, Z, Shinjo, R, Blusztajn, J, Wang, X & Chen, S 2017, 'Tracing subducted sediment inputs to the Ryukyu arc-Okinawa Trough system: Evidence from thallium isotopes', *Geochimica et Cosmochimica Acta*, vol. 217, pp. 462–491.

Smith, WHF & Sandwell, DT 1995, 'Marine gravity field from declassified Geosat and ERS-1 altimetry', *EOS Transactions American Geophysical Union*, vol. 76, Fall Mtng Suppl, F156.

Smith, WHF & Sandwell, DT 1997, 'Global sea floor topography from satellite altimetry and ship depth soundings', *Science*, vol. 277, issue 5334, pp. 1956–1962.

Suetova, IA, Ushakova, LA, & Lemenkova, P 2005, 'Geoinformation mapping of the Barents and Pechora Seas', *Geography and Natural Resources*, vol. 4, pp. 138–142.

Ujie, Y 2000, 'Mud diapirs observed in two piston cores from the landward slope of the northern Ryukyu Trench, northwestern Pacific Ocean', *Marine Geology*, vol. 163, pp. 149–167.

Weatherall, P, Marks, KM, Jakobsson, M, Schmitt, T, Tani, S, Arndt, JE, Rovere, M, Chayes, D, Ferrini, V & Wigley, R 2015, 'A new digital bathymetric model of the world's oceans', *Earth and Space Science*, vol. 2, issue 8, pp. 331–345.

Wessel, P & Smith, WHF 1996, 'A global self-consistent, hierarchical, high-resolution shoreline database', *Journal of Geophysical Research*, vol. 101, pp. 8741–8743.

Wessel, P & Smith, WHF 2006, 'New, improved version of the generic mapping tools released'. *EOS Transactions American Geophysical Union*, vol. 79, issue 47.

Wessel, P & Smith, WHF 2018, *The Generic Mapping Tools. Version 4.5.18 Technical Reference and Cookbook* [Computer software manual], US.