# NOAA Marine Geophysical Data and a GEBCO Grid for the Topographical Analysis of Japanese Archipelago by Means of GRASS GIS and GDAL Library

Polina Lemenkova

► To cite this version:

Polina Lemenkova. NOAA Marine Geophysical Data and a GEBCO Grid for the Topographical Analysis of Japanese Archipelago by Means of GRASS GIS and GDAL Library. Geomatics and Environmental Engineering, 2020, 14 (4), pp.25-45. 10.7494/geom.2020.14.4.25. hal-03011892

HAL Id: hal-03011892

https://hal.science/hal-03011892

Submitted on 18 Nov 2020

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
Polina Lemenkova<fnref n="1" />

# NOAA Marine Geophysical Data and a GEBCO Grid for the Topographical Analysis of Japanese Archipelago by Means of GRASS GIS and GDAL Library<fnref n="2" />

**Abstract:** This article analyzes topographical and geological settings in the Japan Archipelago for comparative raster data processing using GRASS GIS. Data include bathymetric and geological grids in NetCDF format: GEBCO, EMAG2, GlobSed, marine free-air gravity anomaly and EGM96. Data were imported to GRASS by gdalwarp utility of GDAL and projected via PROJ library. Methodology includes data processing (projecting and import), mapping and spatial analysis. Visualization was done by shell scripting using a sequence of GRASS modules: 'd.shade' for relief mapping, 'r.slope aspect'; for modelling based on DEM, 'r.contour' for plotting isolines, 'r.mapcalc' for classification, 'r.category' for associating labels, and auxiliary modules (d.vect, d.rast, d.grid, d.legend). The results of the geophysical visualization show that marine free-air gravitational anomalies vary in the Sea of Japan (–30 to above 40 mGal) reflecting density inhomogeneities of the tectonic structure, and correlating with the geological structure of the seafloor. Dominating values of geoid model are 30–45 m reflecting West Pacific rise, determined by deep density inhomogeneities associated with the mantle convention. Sediment thickness varies across the sea reflecting its geological development with density of 2 km in its deepest part and thinner in central part (Yamato Rise). The aspect map and reclassified map express gradient of the steepest descent.

**Keywords:** GRASS GIS, GEBCO, NOAA, Japan, topography, mapping

Received: 7 May 2020; accepted: 22 June 2020

© 2020 Author. This is an open access publication, which can be used, distributed and reproduced in any medium according to the Creative Commons CC-BY 4.0 License.

<footnote n="1">¹ Ocean University of China, College of Marine Geo-sciences, Qingdao, China, email: pauline.lemenkova@gmail.com
ORCID ID: https://orcid.org/0000-0002-5759-1089</footnote>

<footnote n="2">² This research was funded by China Scholarship Council (CSC), State Oceanic Administration (SOA), Marine Scholarship of China, Grant No. 2016SOA002, China</footnote>
# 1. Introduction

The current study presents GRASS GIS based research which entails the geoprocessing of several different raster grids and data analysis to visualize the regional geological and topographical settings in the Japanese archipelago: DEM based topographic modelling, geoid, marine free-air gravity, sediment thickness of the oceanic seafloor and magnetic anomalies. A special focus of this paper is placed on the methodological description of GRASS modules and snippets of scripting codes, together with visualized output maps demonstrating the cartographic functionality of GRASS GIS.

The main advantage of GRASS GIS is its varied functionalities: to process a wide variety of data which come in different formats (vector, raster netCDF, GRD, GeoTIFF). GIS plays a key role in the processing and visualizing of marine geospatial data, since it enables spatial analysis using advanced data processing and a variety of GIS approaches or statistical libraries with embedded mathematical equations, algorithms and models [1-3]. The proliferation of big data in the geosciences has brought a considerable amount of marine geological datasets to the fore and made it accessible for modelling and cartographic visualization.

Automatic data processing has actively been a question in cartography since the development of programming languages for the rapid processing of big data and methods for the automatization and optimization of data processing have constantly been progressing [4-6]. Therefore, another important asset of GRASS GIS consists in its applicability with Python and R, advanced high-level programming languages effectively used in geosciences [7, 8]. The use of geostatistical models, algorithms of data processing and approaches of data visualization has had a long history in spatial analysis [9-16]. The field of geospatial modeling had made certain advances by the early 1970s, along with rapid development of the IT industry. Advances in computational modelling methods applied in GIS have helped to bring the power of data modeling to classic cartography.

The main research problem consists in highlighting relationships established between the geospatial datasets indicating a correlation between geographic variables showing the complexity of the geologic formation that depends on a set of impact factors. For instance, the impact of the tectonic development occurs at the stage of the relief formation (faults, ridges, deep-se trenches, fracture zones) which finally results in the geomorphic signatures of the modern relief. The seafloor morphology, geomorphic signatures and topographic slope orientation of the Japanese archipelago show remarkably different features in geographically distant regions of the Sea of Japan, Japan Trench and archipelago. Hence, the research problem is to highlight the correlations between these variables, such as Earth magnetic anomaly and marine free-air gravity, combined with sediment thickness, GEBCO relief bathymetry and geoid. Spatial data integrity may reveal correlation of variable geomorphic signatures within and between the datasets overlaid with slope/aspect map and reclassified topographic models. In this way, it can indicate what are the regional differences
in geological, geophysical and tectonic settings and how does the topographic slope/aspect distribution of the modern relief mirrors variations in the geologic settings.

The technical research problem consists of the integration of multi-source thematic datasets to combine different raster data grids into one coherent GRASS GIS project using a scripting approach by importing and transforming GDAL utilities (r.in.gdal, gdalwarp, gdaloustate), and performing numerical spatial analysis by a set of GRASS GIS modules (r.slope-aspect, d.rast, d.shade, g.region etc). Aggregation of various high-resolution datasets undertaken in this study aimed to analyze geological, geophysical and topographic settings of the Japanese archipelago region for proper understanding of the deep-sea areas in the oceanic seafloor. There is a potential possibility that accurate orientation of the magnetic lineation zones, gravity anomalies, geologic fault lines, slabs, fracture zones striking in N-S, W-E slope orientation or various combinations of such orientations (e.g. NNW-SEE, NNE-SSW) correlate with geological lineations and their separation lines, and may help with geological prognosis and analysis. Accurate topographic modelling of the slope orientation (by r.slope-aspect and d.shade GRASS GIS modules) enables to compare distribution of the geophysical settings over study area. Different dynamic mechanisms of the seafloor spreading explain the phenomena of the variability of the sub-marine relief, continental rifting, slab pull and trench formation, during convergence of the Pacific Plate. According to gravimetric surveys and tectonic modelling, tectonic plate subduction plays an important role in trench formation and geophysical settings of the neighboring area. Therefore, a complex analysis of the multi-source datasets (by d.rast GRASS GIS module) clarifies the correlation between the geological and topographic variables.

The main goal of this study was to investigate if combination of diverse available spatial raster datasets including GEBCO terrain topography, Earth Magnetic Anomaly Grid EMAG2, GlobSed sediment thickness global grid, marine free-air gravity anomaly and geoid EGM96, processed by GRASS GIS allows to produce accurate thematic geological maps using scripting approach for complex geospatial analysis of the Japanese archipelago. Specifically, the study intended to produce a series of maps showing contemporary geologic settings in the Japanese Island area and the Sea of Japan at high spatial resolution based on high-resolution raster grids (one to five arc minute resolution) to allow consistent spatially explicit overlay of the topography, geopotential geoid grid, marine free-air gravity and magnetic anomalies, sediment thickness, and topographic modelling of slope and aspect based on raster data obtained from NOAA repository sources. To reach this goal, the study applied a scripting GRASS GIS algorithm to process several raster grids using GRASS code syntax, visualize and integrate various raster data layers in GRASS GIS environment through a sequence of commands. In the end, the study assessed topographic surface of the Japanese archipelago area through the reclassification of the initial topographic grid using algorithm of 'r.mapcalc' data processing and 'if-else' condition in the block of code followed by the 'r.category' operator for data analysis.
Finally, this work was meant to provide answers to the more specific research questions:

- How does the functionality of GRASS GIS enable to model and map geographic objects based on the netCDF, GRD and IMG raster numerical grids, aimed to visualize geographically continuous phenomena as xyz datasets (magnetic anomalies, sediment thickness, geopotential grids, gravimetric anomalies)?

- To what extent does it contribute to the development of the cartographic automatization with the ultimate aim of rapidly processing of the large datasets common in marine geology?

Answering these questions, the applied approach of GRASS GIS scripting algorithms enabled to create a series of the geological and topographical maps of the Japanese archipelago at the region scale, relating them to the phenomena of geospatial correlation between the topographical, geophysical and geological variables applying automatization techniques as advanced GIS methods.

## 2. Study Area

The study area is located in the Japanese archipelago and Sea of Japan within the square of coordinates 128°E-150°E, 30°N-46°N, Fig. 1. The Sea of Japan is characterized by a narrow shelf along the coasts of Korea, and Japan, with mean depths of 120–150 m [17]. The shelf break marks changes in continental slopes located at depths varying from 20 to 550 m (light blue, Fig. 1). The western coast of Japan is complicated by numerous underwater cliffs and small islands [18]. The continental slope in the northern part of the Sea of Japan is gentle, while along the coast of Korea is formed by a steep staircase-like structure with depths exceeding 3000 m. In the east of Korea there is a marginal plateau with depths of 1000–1200 m [19]. A submarine ridge extends in northward direction from the southern part of Honshu Island, reaching Yamato Rise in central part of Sea of Japan, with depths ranging from 300 to 1000 m. The seafloor of the Sea of Japan is flattened with occasional seamounts and slightly inclined to the north, which is contrasting with other seas of the Pacific Ocean: the Philippine Sea with complex and diverse seafloor relief, geological faults and fracture zones [20] and such prominent geomorphic features as the Philippine and Mariana trenches, containing the deepest point on the Earth [21]. The depths in the northern part of the basin exceed 3600 m, and in the southern part decrease up to 2600–2800 m.

The basement surface of the Mesozoic and Cenozoic structures of Japan descends from the coast to the seafloor of the basin as a steep slope with complex block division. The maximum depths of the basement seafloor in the Sea of Japan exceeds 5 km, and in its central part (Yamato Rise) 2–3 km. The geologic basement of the Japan Trench lies ca. 0.5–1.5 km deeper than its modern seafloor. The Japanese islands with the Japan Trench form a typical arc-trench system with developed continental structures. From the oceanic side, the Japanese archipelago is bordered by
a narrow shelf with depths of less than 100 m, which smoothly continues to a slope of the deep-sea Japan Trench [22] (Figs. 1, 2).

<figure><img src="image_2.png" /><figcaption>**Fig. 1.** Topographic map of the Japanese archipelago. Mapped using GRASS GIS (d.rast module), based on GEBCO raster grid imported via GDAL (r.in.gdal)</figcaption></figure>

<figure><img src="image_3.png" /><figcaption>**Fig. 2.** Shaded relief map of the Japanese archipelago region. Mapping done using GRASS GIS (d.shade and r.slope(aspect modules), GEBCO raster grid</figcaption></figure>
The profile of the slope is convex, with steepness gradually increasing down-
wards. The Japan Trench has a flat seafloor, with a maximum depth of 8,046 m ac-
cording to the GMRT data [23] which is shallower when compared to the deepest
trenches of the Pacific Ocean, e.g. Tonga and Kermadec trenches with depths ex-
ceeding 10,000 m [24]. The convergent plate margin of the seaward edge of the con-
tinental plate is deformed by subduction of the oceanic Pacific Plate. The northern
Japan Trench is demarcated from the rigid continental plate as reported by a multi-
channel seismic survey [25].

The tectonic structure of the seafloor of the Sea of Japan, its continental shelf, Ja-
pan Trench and northwest segment of the Pacific Ocean are summarized in previous
works [26–31]. The seafloor structure of the Japanese archipelago and its geological
development have been detailed and reported using the datasets from regular cruise
expeditions published as series of volumes “Initial Reports of the Deep Sea Drilling
Project” [32–34].

**3. Materials and Methods**

The methodology of this work is based on the GRASS GIS (Geographic Re-
sources Analysis Support System), an open source GIS, continuously developed by
the GRASS Team since 1982 [35]. In contrast to traditional GIS, e.g. ArcGIS [36–39],
GRASS GIS is shell-scripting based GIS with advanced functionality in coding
and applicability to Python programming language and GDAL library [40]. GDAL
was used for re-projecting raster grids into cartographic projection (Equal Area Cy-
lindrical projection was selected for this research), warping, subsetting and convert-
ing data. The comparative analysis of several grids of datasets aims at highlighting
correlation between the geographical variables visualized on the maps.

Data selection was performed based on the data quality: reliability (NOAA)
and high-resolution. Spatial resolution of the input data is one of the key concepts
in GIS since the Earth's complexity restricts studies in full detail. Resolution is crit-
ical in determining a dataset's suitability for a given research problem, as it affects
the minuteness at which objects can be mapped and interpreted during spatial anal-
ysis. Technically, the resolution determines the volume of the processed geospatial
grids which affects computer memory and data storage volume.

Therefore, the following raster grids were selected as input data: General
Bathymetric Chart of the Oceans (GEBCO) with 15-second resolution for terrain
and 30 arc-second bathymetric grid for ocean seafloor [41], Earth Magnetic Anom-
aly Grid (EMAG2) with 2-arc-minute resolution [42], sediment thickness global
grid [43, 44], marine free-air gravity anomaly grid from Geosat and ERS 1 satellite
altimetry [45], geoid 15-minute grid EGM96 Geopotential Model grid [46].

The mapping is done using the raster grids visualized by GRASS GIS through
a sequence of modules using three main methodological blocks:
1) data preprocessing,

2) data visualization (mapping),

3) spatial analysis.

The coding methodology for each of these blocks is explained below stepwise.

## Data Preprocessing

Data preprocessing included the retrieval of information from the different raster grids aimed at highlighting coherency between geological, topographical and geophysical grids. Due to the differences between the original datasets in terms of content and consistency, it included the analysis of their metadata, spatial characteristics and structure: projection, datum, extent, range, spatial resolution.

The multi-source data origin resulted in their various resolutions: 30 arc-second GEBCO, 2 arc-minute EMAG2, 5 arc-minute GlobSed, 15 arc-minute grid EGM96, 1 arc-minute marine free-air gravity. Thus, data preprocessing in GRASS included the following four steps:

1. The first step involves cartographic projecting of the initial raster NetCDF in WGS 84 datum surface imported to GRASS GIS and warping it to an Equal Area Cylindrical projection by 'gralwarp' utility of GDAL using PROJ Cartographic Projections and Coordinate Transformations Library [47] via GMT [48]. The projecting was done using following code (here is the example for the NetCDF grid of topographic layer (Figs. 1, 2), repeated for other thematic grids (Figs. 3-6), respectively): `gdalwarp -t_srs '+proj=cea lat_ts=38 lon_0=138' jt_relief.nc jt_relief_EAC.nc`.

2. The second step involves the re-projecting of the raster grids read into the GRASS GIS project using 'r.in.gdal' module by following code: `r.in.gdal jt_relief_EAC.nc out=jt_relief_EAC title="Topography: GEBCO grid"`.

3. The third step involves applying date and time function of GRASS GIS: `r.timestamp map=jt_relief_EAC date='04 May 2020'`.

4. The list of the raster grids and information (projection, datum, min/max range) was received by modules 'g.list rast' which visualizes the available raster files in the current project and 'r.info'.

Vector data import from the GMT was done by 'v.in.ogr' command: `v.in.ogr trench.gmt out=trench`. The compatibility of the GRASS GIS with native GMT formats is an advantage of the GRASS GIS, as GMT is a functional cartographic toolset often used in geophysical mapping [49-51].

## Cartographic Mapping

Topographic dataset was processed based on the isolines generalization highlighting the particularities of the relief. Data at different levels of representation containing specific information (resolution and smoothness of magnetic and gravity anomalies, topography, geoid, sediment thickness) were integrated in a GRASS GIS project.
The resulting cartographic mapping and data visualization included following eight steps:

1. Displaying a raster map by initialized GRASS GIS monitor: d.mon wx0 and creating spatial subsets of the region of Japanese archipelago using the 'g.region module': g.region raster=jt_relief_EAC -p

2. Applying color palette for the maps via 'rcolors' GRASS module. Following color tables were selected for the maps: 'srtm_PLUS' (Fig. 1), 'haxby' [52] (Figs. 2, 6), 'roygbiv' (Fig. 3), 'celsius' (Fig. 4), 'rainbow' (Fig. 5), 'aspectcolr' (Fig. 7) and defined colors by 'rules' function for Figure 8. Color assignment was done by code: rcolors jt_relief_EAC col=haxby

3. Spatial analysis of slope, shaded relief and plotting isolines were performed using following modules: 'r.slope aspect', 'd.shade', 'r.contour' and 'd.vect', as follows in the GRASS code:
   r.slope aspect elevation=jt_relief_EAC aspect=jtAspect_EAC
d.shade shade=jtAspect_EAC color=jt_relief_EAC (the 'd.shade' module resulted in overlay of the topographic DEM and aspect models, in which raster cell values are equally spaced, gridded according to a matrix structure of topography and aspect which gives the effect of shaded relief visualization (Fig. 2).
r.contour jt_relief_EAC out=TopoJT_EAC step=2000 -overwrite
d.vect TopoJT_EAC color='75:23:3' width=0.

4. The next step included plotting a vector border box around the study area (red line on Figure 1) by module 'v.in.region' using following code:
v.in.region output=jt_EAC_bbox
v.info map=jt_EAC_bbox
d.vect jt_EAC_bbox color=red width=3 fill_color="none".

5. Cartographic projected coordinate grids along the map borders were plotted using following GRASS code:
d.grid size=5 color='172:219:250' border_color=yellow width=0.1 fontsize=8 text_color=red bgcolor=white.

6. Legend to the left of each map was plotted using following code where self-explaining cartographic attributes, data range and visual settings of the legend are plotted by flags, such as font size, color, border, background color (here is example for the topographic map (Fig. 1):
d.legend raster=jt_relief_EAC range=-9759.701,3700.742 title=Topography,m title_fontsize=8 font= Arial fontsize=7 -t -b bgcolor=white label_step=1000 border_color=gray thin=8.

7. Annotating textual information was done by following code example:
d.text text="Honshu Island" color=green size=2.5 font="Trebuchet MS".

8. Adding the title of the maps was piped by Unix redirection combining 'd.title' and 'd.text' GRASS modules using following code example: d.title map=jt_relief_EAC | d.text text="Shaded topography: Japan Islands" color="green" size=3.
## Spatial Analysis

After the raster grids were visualized (Figs. 1–6), the next step included topographic spatial analysis (Figs. 7, 8). The spatial analysis performed on the topographic dataset can be considered as a derivation process of the initial reference grid (GEBCO) aimed to retrieve information on the slope and aspect of the topography. The spatial analysis was made by following steps:

1. The compass orientation was generated and classified by four major directions (N, E, S, W) using adopted code [53] developed based on the algorithms of data analysis [54–56]:

```python
r.mapcalc "aspect_4 directions_Japan = eval(\ 
    compass=(450 - jtAspect_EAC)% 360, \ 
    if compas >=0. && compass < 45., 1) \ 
    + if compas >=45. && compass < 135., 2) \ 
    + if compas >=135. && compass < 225., 3) \ 
    + if compas >=225. && compass < 315., 4) \ 
    + if compas >=315., 1) \ 
    )"
```

2. The next step included assigning of text labels using 'r.category' module and condition encountered by a EOF (end-of-file) operator:

```python
r.category aspect_4 directions_Japan separator=comma rules=- << EOF
1,North 2,East 3,South 4,West EOF
```

3. Assignment of color table was done using previously explained 'rcolors' GRASS GIS module using RGB triplets: `rcolors aspect_4 directions rules=- << EOF 1 255,234,0 2 136,72,152 3 230,0,51 4 195,216,37 EOF`

The resulting map re-classified from Figure 7 is presented on Figure 8.

## 4. Results

The topographic map (Fig. 1) and derived shaded relief map (Fig. 2) are based on the GEBCO global bathymetry and topography grid with a spatial resolution interval of 15 arc seconds where bathymetry is produced using a combination of shipboard and satellite altimetry predicted data.

The map of the geoid in the Sea of Japan (Fig. 3) is based on the visualized EGM96 raster grid. It shows the transition zone of geoid undulations from the Pacific Ocean to the continental area. Composition relationships in geophysical datasets are based on the analysis of several sub-entities (marine free-air gravity anomaly values distribution correlated with geoid and topographic isolines over the study area). Thus, the region of the Japan Archipelago is located in the vast West Pacific rise of the geoid (the long-wave part of the spectrum of the geoid).
<figure><img src="image_4.png" /><figcaption>**Fig. 3.** Geoid model of the Japanese archipelago.
Mapping done using GRASS GIS, based on EGM96 raster grid processed by GDAL (gdalwarp, gdaloust translate)</figcaption></figure>

Therefore, as can be seen in Figure 2, the dominating values are 30–45 m (orange to red colors on the map), followed by 26–30 m (yellow colors). Lower geoid values (17–25 m, green and 10–17, cyan) are mostly repeating topographic depressions of the Pacific Ocean and the lowest –5 to 10 (blue to purple) corresponding to the Japan Trench topography (Fig. 3). The nature of the geoid is determined by deep density inhomogeneities that may be associated with the mantle convention.

Map of the marine free-air gravity (Fig. 4) shows gravimetric patterns in the basin of the Sea of Japan. There are weak positive anomalies (0–5 mGal, cyan colors, Fig. 4) or negative anomalies (0 to –17 mGal, dark blue colors, Fig. 4) and in some depressions, larger negative values (–17 up to –30 mGal, purple colors, Fig. 4) are observed. The Yamato Rise (in the central part in the Sea of Japan) is distinguished by positive Faye anomalies (>20 mGal, orange to red colors, Fig. 4). Above the Japanese Islands, the marine free-air gravity anomalies reach values above 40 mGal and higher (grey-pink colors, Fig. 4). They sharply decrease to a deep-sea trench, up to strongly negative values (below –30 mGal, white color, Fig. 4).

The marine free-air gravitational anomalies field varies significantly in the Sea of Japan (from –30 mGal to above 40 mGal) reflecting density inhomogeneities of the tectonic structure, but in general they correlate with geological structures
of the seafloor. On the shelf they are mostly positive (20–40 mGal, orange to red colors). The shelf area is characterized by linearly elongated anomalies parallel to the shelf isolines. In the Japan Trench, they reach values lower than –30 mGal.

<figure><img src="image_5.png" /><figcaption>Fig. 4. Marine free-air gravity map of the Japanese archipelago region. Mapping done using GRASS GIS, based on 'grav_27.1.img' gravity raster grid (SIO)</figcaption></figure>

An accurate and consistent visualization of the marine free-air gravity field in the Japan region is crucial for geophysical studies of western Pacific Ocean and thus also for a better understanding of the geological structures of the Japanese archipelago. The map of the magnetic anomalies of Japan (Fig. 5) shows the magnetic anomalies zone, which is one of the three magnetic anomalies zones of the Pacific Ocean forming a system of ancient magnetic anomalies which reflects the basin development in Mesozoic period (the other two are Phoenix and Hawaiian). Theoretical basis of the magnetic anomalies of the oceans is described in the existing studies [57, 58] and explained by regional particularities of the magnetic lineament anomalies [59] or anisotropy of magnetic susceptibility [60]. The geometric relationships analyzed from the comparison of datasets on magnetic anomalies, topography, geoid, and marine free-air gravity are meant to translate spatial relationships between the geographical and geophysical entities. These are especially concerned with the notions of the geometric similarities notable on data value distribution (depicted as colored areas) and isolines (depicted as lines), the intersection as overlay in a GRASS GIS project and spatial connection between the data.
<figure><img src="image_6.png" /><figcaption>**Fig. 5.** Magnetic anomaly grid model in the Japanese Archipelago region. Data: EMAG2v3 2-arc-min resolution raster grid, NOAA. Model mapping using GRASS GIS (modules: r.in.gdal, d.rast). Data subset projection: GDAL</figcaption></figure>

Thus, the Japanese system of the magnetic anomalies occupies the western and northern parts of the Northwest Pacific Basin which correlates with the topographic subdivision of the basin based on the elevation distribution. The magnetic anomalies here are oriented as strips in a northeast direction (Fig. 5). The magnetic anomalies of Japan region are separated from the Hawaiian system by the Shatsky Rise which is the 3rd Earth's largest oceanic basaltic plateau formed during a period of geomagnetic reversals [61, 62]. The magnetic anomalies are then 'absorbed' by the Japanese Trench. The amplitudes of the magnetic anomalies exceed 400 gammas in their highest values (Fig. 5, red colors).

The correspondence between the dataset of sediment thickness and topographic relief represent two closely correlated phenomena of marine geological settings: sedimentation and bathymetry. Map of sediment thickness (Fig. 6) shows distribution of the marine sediments as important variables in geological datasets. Dependencies in sedimentation distribution can be clearly observed comparing two maps (Figs. 1, 6) showing higher sedimentation values in shallow shelf areas and river mouth estuaries and almost absent or low values in the open ocean basin.

Visualization of sediments reflects seafloor sedimentation records the Earth's tectonic movements during geological history, and may assist in reconstructing palaeoclimate data. Visualized sediment maps can assist in highlighting correlations
of oceanological parameters, such as ocean current directions and strength, changes
in sea level and sea temperature. Besides, seafloor sedimentation is important for
studies of the marine biological productivity and variations in environmental con-
ditions, biological evolution and ecology since it indirectly reflects the complexity
of the climatic and geological processes [63]. Therefore, it has always been important
to visualize seafloor sediment coverage for geological studies.

<figure><img src="image_7.png" /><figcaption>Fig. 6. Sediment thickness map of the Japanese Archipelago region.
Mapping done using GRASS GIS, based on GlobSed 5-arc-min raster grid V.3</figcaption></figure>

Sediment coverage in the Japan Trench region is contrasting with other trenches
of the Pacific Ocean, e.g. Middle America Trench having Quaternary trench-filling
turbidites, fault scarp in lower Miocene chalk, pelagic clays and diatom ooze [64], or
Mariana Trench with sedimentary cover of basalt of the leitic composition, a com-
plex of parallel dykes of diabases, an isotropic gabbro, a banded gabbro-ultrabasic
complex [65]. In the deepest part of the basin of the Sea of Japan, the sediment thick-
ness represented by Meso-Cenozoic sediments is about 2 km. In its central segment
on the Yamato Rise and other the sediment coverage is thinner and has a sharply
variable thickness: it decreases to a minimum on the individual peaks and steep
slopes and increases in local depressions and saddles [66]. Some data received from
drilling in the Yamato Rise discovered 529 m thick layer of siliceous-siliceous silts,
gradually changing to diatomaceous silts of the Upper Miocene, followed by volca-
nic sand and tuffs. Further reading regarding regional patterns of the sedimentation
at the seafloor of the Pacific Ocean and theoretical processes of sediment deposition through modelling and observations are presented recently [67, 68] and in earlier works [69].

The dependency relationships were established based on the analysis of the datasets visualized as several raster grids. They were used to compute the attribute values of the marine free-air gravity, topographic grids, sedimentation and geoid data with slope/aspect reclassified map. For geographic databases, the geometry was considered as an attribute for data comparison. The topographic spatial analysis is presented on Figures 7 and 8. Figure 8 shows the results of the reclassification of a continuous field map showing aspect of the topographic slope in the Japanese archipelago (Fig. 7).

<figure><img src="image_8.png" /><figcaption>Fig. 7. Aspect map of the slope topography: Japanese Archipelago region. Mapping done using GRASS GIS, based on GEBCO elevation raster grid</figcaption></figure>

The initial aspect angle ranged from 0 to 360 degrees with the origin in the east and a counter-clockwise rotation. The reclassified map (Fig. 8) shows the aspect map with values grouped into four aspect clusters reflecting the terrain compass orientation (North-South-West-East). The topographic parameters describe the geometry of the terrain and seafloor surface in a given point, such as slope and aspect. Estimating and visualizing a slope aspect and steepness is useful for the analysis of such hazardous geological parameters as soils creep, debris avalanching and landslides [70]. Transport processes on the slopes include surface processes (e.g. flow)
and subsurface process (leaching) for middle-scale spatial analysis and large mass movements that are the subject of gravity tectonics, for a small-scale analysis.

<figure><img src="image_9.png" /><figcaption>Fig. 8. Reclassified aspect map of the topography. Mapping done using GRASS GIS (modules r.slope aspect, r.mapcalc, r.category), based on GEBCO raster grid</figcaption></figure>

A geometrical arrangement of patterned ground on the slopes is represented by geomorphological landforms such as polygons, local depressions, circles, steps etc. Slope analysis is essential for the geomorphological viewpoint, since both terrain and submarine ground surfaces are composed of the landform elements related to slope: curved geometric shapes, related to upslope, downslope, and orientation. Hence, altitude, extent, slope, curvature, ruggedness, and aspect are essential parts of the spatial analysis based on DEM [71]. The slope of the terrain surface allows gravity to induce the flow of materials (debris, rocks) and water. Therefore, it lies at the core of many geological models.

The aspect of the slope is a compass direction of the maximum rate of variation. Continuous classes in the topographic aspect modelling were first visualized (Fig. 7) and then modeled based on the 'r.mapcalc' re-classification algorithm applied to the initial dataset of GEBCO, merged in four classes in the derived grid by means of the attribute re-calculation (Fig. 8). The reclassified maps express the orientation of the gradient of the steepest descent, converted to a compass attributes N-S-W-E. Comparison of the map on Fig. 8 with topographic DEM (Figs. 1, 2) shows that slope aspect in the Japanese archipelago region corresponds to the general topography.
# 5. Conclusion

A unique aspect of the submarine geomorphology comparing to the terrestrial studies lies in the fact that seabed is not available for direct observations. Marine geological datasets could only be received by the remote sensing methods of the seafloor survey and data post-processing using GIS. In view of this, the selection of the advanced GIS tool for data processing is clear. The GRASS GIS, as demonstrated in this paper, presents an excellent tool for the processing of the raster datasets enabling to visualize, perform spatial topographic analysis and map thematic grids.

The scripting approach provided by GRASS GIS affords the opportunity to combine a number of thematic raster grids and visualize layers using selected code lines of GRASS GIS syntax. Such a collection of codes can be then saved as shell scripts. Using shell scripts in GIS, rather than plotting maps through a Graphical User Interface (GUI), provides significant editing and analytical possibilities due to the principle of repeatability. Shell scripting presents a new way in cartographic automatization through visualizing spatial big datasets in a rapid and effective way. In this sense, this paper presented a contribution to the development of the cartographic methods with a case study of complex analysis of the marine geological settings of the Japanese archipelago by raster datasets.

The GRASS GIS scripting-based mapping relied on using a specific code syntax used for processing a series of raster datasets, as demonstrated in this research. Such a technical approach may use the advantage of the repeatability of scripts and perform a series of maps with identical spatial extend and coverage reflecting correlations between topographic setting, geophysical anomalies and geologic conditions of the study area. As one of the goals driving this research was the visualization of the series of the thematic maps by scripts written on GRASS GIS syntax, such a paradigm demanded a careful evaluation and standardization of the input multi-source data (projections, resolution, spatial extent, coverage) to reduce the related uncertainty of the output spatial series of the geological, geophysical and topographical maps.

# References

[1] Lemenkova P.: GMT Based Comparative Geomorphological Analysis of the Vityaz and Vanuatu Trenches, Fiji Basin. Geodetski List, vol. 74 (97), no. 1, pp. 19–39.

[2] Klaučo M., Gregorová B., Stankov U., Marković V., Lemenkova P.: Landscape metrics as indicator for ecological significance: assessment of Sitno Natura 2000 sites, Slovakia. Ecology and Environmental Protection, 2014, pp. 85–90.

[3] Lemenkova P.: Testing Linear Regressions by StatsModel Library of Python for Oceanological Data Interpretation. Aquatic Sciences and Engineering, vol. 34, 2019, pp. 51–60.
[4] Schenke H.W., Lemenkova P.: Zur Frage der Meeresboden-Kartographie: Die Nutzung von AutoTrace Digitizer für die Vektorisierung der Bathymetrischen Daten in der Petschora-See. Hydrographische Nachrichten, vol. 81, 2008, pp. 16–21.

[5] Farquharson C.G., Oldenburg D.W.: A comparison of automatic techniques for estimating the regularization parameter in non-linear inverse problems. Geophysical Journal International, vol. 156, no. 3, 2004, pp. 411–425.

[6] Lemenkova P.: Automatic Data Processing for Visualising Yap and Palau Trenches by Generic Mapping Tools. Cartographic Letters, vol. 27, no. 2, 2019, pp. 72–89.

[7] Cockett R., Kang S., Heagy L.J., Pidlisecky A., Oldenburg D.W.: SimPEG: An open source framework for simulation and gradient based parameter estimation in geophysical applications. Computers & Geosciences, vol. 85, part A, 2015, pp. 142–154.

[8] Lemenkova P.: Geospatial Analysis by Python and R: Geomorphology of the Philippine Trench, Pacific Ocean. Electronic Letters on Science and Engineering, vol. 15, no. 3, 2019, pp. 81–94.

[9] Deutsch C.V., Journel A.G.: GSLIB: Geostatistical software library and user's guide. 2nd ed. Oxford University Press, 1998.

[10] Lemenkova P.: Processing oceanographic data by Python libraries NumPy, SciPy and Pandas. Aquatic Research, vol. 2, 2019, pp. 73–91.

[11] Lemenkova P.: 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, 2018, pp. 35–42.

[12] Lemenkova P.: Numerical Data Modelling and Classification in Marine Geology by the SPSS Statistics. International Journal of Engineering Technologies, vol. 5, no. 2, 2019, pp. 90–99.

[13] Lindh P.: Compaction- and strength properties of stabilised and unstabilised fine-grained tills. Lund 2004 [Ph.D. thesis].

[14] Klaučo M., Gregorová B., Stankov U., Marković V., Lemenkova P.: Interpretation of Landscape Values, Typology and Quality Using Methods of Spatial Metrics for Ecological Planning. [in:] Poster at the 54th International Conference Environmental & Climate Technologies (Riga Technical University, October 14, 2013. Riga, Latvia). https://doi.org/10.13140/RG.2.2.23026.96963.

[15] Lemenkova P.: GMT Based Comparative Geomorphological Analysis of the Vityaz and Vanuatu Trenches, Fiji Basin. Geodetski List, vol. 74 (1), 2020, pp. 19–39.

[16] Salles C., Poesen J., Govers G.: Statistical and physical analysis of soil detachment by raindrop impact: rain erosivity indices and threshold energy. Water Resources Research, vol. 36, 2000, pp. 2721–2729.

[17] Udintsev G.B.: Geomorfologiya i tektonika dna Tikhogo okeana. Nauka, Moskva 1972 [Удинцев Г.Б.: Геоморфология и тектоника дна Тихого океана. Наука, Москва 1972].

[18] Lemenkova P.: Variations in the bathymetry and bottom morphology of the Izu-Bo nin Trench modelled by GMT. Bulletin of Geography. Physical Geography Series, vol. 18 (1), 2020, pp. 41–60.
[19] Litvin V.M.: Morfostruktura dna okeanov. Nedra, Leningrad 1987 [Литвин В.М.: Морфоструктура дна океанов. Недра, Ленинград 1987].

[20] Lemenkova P.: Visualization of the geophysical settings in the Philippine Sea margins by means of GMT and ISC data. Central European Journal of Geography and Sustainable Development, vol. 2, no. 1, 2020, pp. 5–15.

[21] Lemenkova P.: Statistical Analysis of the Mariana Trench Geomorphology Using R Programming Language. Geodesy and Cartography, vol. 45, no. 2, 2019, pp. 57–84.

[22] Nakanishi M.: Bending-related topographic structures of the subducting plate in the northwestern Pacific Ocean. [in:] Ogawa Y., Anma R., Dilek Y. (eds.), Accretionary Prisms and Convergent Margin Tectonics in the Northwest Pacific Basin, Springer Science Business Media B.V, Berlin 2011, pp. 1–38.

[23] Ryan W.B.F., Carbotte S.M., Coplan J.O., O'Hara S., Melkonian A., Arko R., Weissel R.A., Ferrini V., Goodwillie A., Nitsche F., Bonczkowski J., Zemsky R.: Global Multi-Resolution Topography synthesis. Geochemistry, Geophysics, Geosystems, vol. 10, no. Q03014, 2009. https://doi.org/10.1029/2008GC002332.

[24] Lemenkova P.: GMT Based Comparative Analysis and Geomorphological Mapping of the Kermadec and Tonga Trenches, Southwest Pacific Ocean. Geographia Technica, vol. 14, no. 2, 2019, pp. 39–48.

[25] Tsuru T., Park J.-O., Takahashi N., Kodaira S., Kido Y., Kaneda Y., Kono Y.: Tectonic features of the Japan Trench convergent margin off Sanriku, northeastern Japan, revealed by multichannel seismic reflection data. Journal of Geophysical Research, vol. 105, 2000, pp. 16403–16413.

[26] Belousov V.V.: Zemnaya kora i verkhyaya mantiya okeanov. Nauka, Moskva 1968 [Белоусов В.В.: Земная кора и верхняя мантия океанов. Наука, Москва 1968].

[27] Lemenkova P.: GMT-based geological mapping and assessment of the bathymetric variations of the Kuril-Kamchatka Trench, Pacific Ocean. Natural and Engineering Sciences, vol. 5, no. 1, 2020, pp. 1–17.

[28] Miura S., Takahashi N., Nakanishi A., Tsuru T., Kodaira S., Kaneda Y.: Structural characteristics off Miyagi forearc region, the Japan Trench seismogenic zone, deduced from a wide-angle reflection and refraction study. Tectonophysics, vol. 407, 2005, pp. 165–188.

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

[30] Kita S., Okada T., Hasegawa A., Nakajima J., Matsuzawa T.: Anomalous deepening of a seismic belt in the upper-plane of the double seismic zone in the Pacific slab beneath the Hokkaido corner: possible evidence for thermal shielding caused by subducted forearc crust materials. Earth and Planetary Science Letters, vol. 290, 2010, pp. 415–426.
[31] Pushcharovski Yu.M.: *Vvedeniye v tektoniku tikhookeanskogo segmenta Zemli*. Nauka, Moskva 1972 [Пущаровский Ю.М.: *Введение в тектонику тихоокеанского сегмента Земли*. Наука, Москва 1972].

[32] Shipboard Scientific Party: *Initial Reports of the Deep Sea Drilling Project*. 56. U.S. Government Printing Office, Washington 1980.

[33] Shipboard Scientific Party: *Site 584*. [in:] Kagami H., Karig D.E., Coulbourn W.T. (eds.), *Initial Reports of Deep Sea Drilling Project*, 87. U.S. Government Printing Office, Washington 1986, pp. 16–17.

[34] Shipboard Scientific Party: *Site 1150*. [in:] Sacks I.S., Suyehiro K., Acton G.D. (eds.), *Proceedings of the Ocean Drilling Program Initial Reports*, vol. 186, TX 77845-9547, USA, College Station 2000.

[35] GRASS Development Team: *Geographic Resources Analysis Support System (GRASS GIS) Software*. Open Source Geospatial Foundation, USA, 2018. https://grass.osgeo.org [access: 1.05.2020].

[36] Suetova I.A., Ushakova L.A., Lemenkova P.: *Geoinformation mapping of the Barents and Pechora Seas*. *Geography and Natural Resources*, vol. 4, 2005, pp. 138–142.

[37] Basista I.: *Application of GIS Tools to Describe the Location of New Registered Parcels*. Geomatics and Environmental Engineering, vol. 14, no. 1, 2020, pp. 5–13.

[38] Klaučo M., Gregorová B., Stankov U., Marković V., Lemenkova P.: *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, no. 1, 2013, pp. 28–42.

[39] Lemenkova P.: *Seagrass Mapping and Monitoring Along the Coasts of Crete, Greece*. University of Twente, Enschede 2011 [M.Sc. thesis].

[40] GDAL/OGR contributors: GDAL/OGR Geospatial Data Abstraction Software Library. Open Source Geospatial Foundation, 2020. https://gdal.org [access: 2.05.2020].

[41] GEBCO Compilation Group: GEBCO 2020 Grid. https://doi.org/10.5285/ a29c5465-b138-234d-e053-6c86abc040b9.

[42] Meyer B., Saltus R., Chulliat A.: *EMAG2v3: Earth Magnetic Anomaly Grid (2-arc-minute resolution)*. National Centers for Environmental Information, NOAA, 2017. https://doi.org/10.7289/V5H70CVX.

[43] Wobbe F., Lindeque A., Gohl K.: *Total sediment thickness grid of the Southern Pacific Ocean off West Antarctica, links to NetCDF files*. Pangea, 2014. https:// doi.org/10.1594/PANGAEA.835589.

[44] Divins D.: *Total sediment thickness of the world's oceans and marginal seas*. NOAA National Geophysical Data Center, Boulder 2003.

[45] Smith W.H.F., Sandwell D.T.: *Marine gravity field from declassified Geosat and ERS-1 altimetry*. EOS, vol. 76 (46), 1995, pp. F156.
[46] Lemoine F.G., Smith D.E., Kunz L., Smith R., Pavlis E.C., Pavlis N.K., Klosko S.M., Chinn D.S., Torrence M.H., Williamson R.G., Cox C.M., Rachlin K.E., Wang Y.M., Kenyon S.C., Salman R., Trimmer R., Rapp R.H., Nerem R.S.: *The development of the Joint NASA GSFC and the National Imagery and Mapping Agency (NIMA) Geopotential Model EGM96*. [in:] *Gravity, Geoid and Marine Geodesy*, International Association of Geodesy Symposia (IAG SYMPOSIA), vol. 117, Springer-Verlag Berlin Heidelberg, 1997, pp. 461–469.

[47] PROJ contributors: PROJ coordinate transformation software library. Open Source Geospatial Foundation, 2020. https://proj.org/ [access: 7.05.2020].

[48] Wessel P., Smith W.H.F., Scharroo R., Luis J., Wobbe F.: *Generic mapping tools: Improved version released*. Eos, vol. 94, no. 45, 2013, pp. 409–410.

[49] Lemenkova P.: *Geomorphological modelling and mapping of the Peru-Chile Trench by GMT*. Polish Cartographical Review, vol. 51, no. 4, 2019, pp. 181–194.

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

[51] Lemenkova P.: *AWK and GNU Octave Programming Languages Integrated with Generic Mapping Tools for Geomorphological Analysis*. GeoScience Engineering, vol. 65, no. 4, 2019, pp. 1–22.

[52] Haxby W.F., Karner G.D., LaBrecque J.L., Weissel J.K.: *Digital images of combined oceanic and continental data sets and their use in tectonic studies*. Eos, vol. 64, no. 52, 1983, pp. 995–1004.

[53] GRASS Development Team: GRASS GIS 7.8.4dev Reference Manual. https://grass.osgeo.org/grass78/manuals/ [access: 7.05.2020].

[54] Horn B.K.P.: *Hill Shading and the Reflectance Map*. Proceedings of the IEEE, vol. 69, no. 1, 1981, pp. 14–47.

[55] Mitasova H.: *Cartographic aspects of computer surface modeling*. Slovak Technical University, Bratislava 1985 [Ph.D. thesis].

[56] Hofierka J., Mitasova H., Neteler M.: *Geomorphometry in GRASS GIS*. [in:] Hengl T., Reuter H.I. (eds.), *Geomorphometry: Concepts, Software, Applications*, Developments in Soil Science, vol. 33, Elsevier, 2009, pp. 387–410. http://geomorphometry.org/book [access: 7.05.2020].

[57] Mirlin Ye.G.: *Magnitnyye anomalii okeanov i novaya global'naya tektonika*. Nauka, Moskva 1981 [Мирлин Е.Г.: Магнитные аномалии океанов и новая глобальная тектоника. Наука, Москва 1981].

[58] Kitamura Y., Kanamatsu T., Zhao X.: *Structural evolution in accretionary prism toe revealed by magnetic fabric analysis from IODP NanTroSEIZE Expedition 316*. Earth and Planetary Science Letters, vol. 292, 2010, pp. 221–230.

[59] Kobayashi K.: *Preliminary report of the Hakuho Maru cruise KH 90-t*. Ocean Research Institute, University of Tokyo, Tokyo 1991.
[60] Yang T., Mishima T., Ujiie K., Chester F.M., Mori J.J., Eguchi N., Toczko S.: *Strain decoupling across the décollement in the region of large slip during the 2011 Tohoku-Oki earthquake from anisotropy of magnetic susceptibility.* Earth and Planetary Science Letters, vol. 381, 2013, pp. 31–38.

[61] Geldmacher J., van den Bogaard P., Heydolph K., Hoernle K.: *The age of Earth's largest volcano: Tamu Massif on Shatsky Rise (northwest Pacific Ocean).* International Journal of Earth Sciences, vol. 103, no. 8, 2014, pp. 2351–2357.

[62] Nakanishi M., Sager W.W., Klaus A.: *Magnetic lineations within Shatsky Rise, northwest Pacific Ocean: Implications for hot spot-triple junction interaction and oceanic plateau formation.* Journal of Geophysical Research: Solid Earth, vol. 104, no. B4, 1999, pp. 7539–7556.

[63] Kuhn G., Hass C., Kober M., Petitat M., Feigl T., Hillenbrand C.D., Kruger S., Forwick M., Gauger S., Lemenkova P.: *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.* [in:] Gohl K. (ed.), *Expeditions program, no. 75 ANT XXIII/4*, 2006.

[64] Lemenkova P.: *Geophysical modelling of the Middle America Trench using GMT.* Annals of Valahia University of Targoviste. Geographical Series, vol. 19, no. 2, 2019, pp. 73–94.

[65] Lemenkova P.: *An Empirical Study of R Applications for Data Analysis in Marine Geology.* Marine Science and Technology Bulletin, vol. 8, no. 1, 2019, pp. 1–9.

[66] Levin L.E.: *Geologiya osadochnogo chekhla dna morey i okeanov.* Nedra, Moskva 1984 [Левин Л.Э.: *Геология осадочного чехла дна морей и океанов.* Недра, Москва 1984].

[67] Storti F., McClay K.R.: *Influence of syntectonic sedimentation thrust wedges in analogue models.* Geology, vol. 23, 1995, pp. 999–1002.

[68] Tanioka Y., Seno T.: *Sediment effect on tsunami generation of the 1986 Sanriku tsunami earthquake.* Geophysical Research Letters, vol. 28, 2001, pp. 3389–3392.

[69] Von Herzen R., Maxwell A.E.: *The measurement of thermal conductivity of deep-sea sediments by a needle-probe method.* Journal of Geophysical Research, vol. 64, 1959, pp. 1557–1563.

[70] Lemenkova P., Promper C., Glade T.: *Economic Assessment of Landslide Risk for the Waidhofen a.d. Ybbs Region, Alpine Foreland, Lower Austria.* [in:] Eberhardt E., Froese C., Turner A.K., Leroueil S. (eds.), *11th International Symposium on Landslides and the 2nd North American Symposium on Landslides & Engineered Slopes (NASL). Protecting Society through Improved Understanding, June 2–8, 2012, Banff, AB, Canada*, 2012, pp. 279–285.

[71] Giles P.T.: *Geomorphological signatures: classification of aggregated slope unit objects from digital elevation and remote sensing data.* Earth Surface Processes and Landforms, vol. 23, 1998, pp. 581–594.