Topographic Model Research Articles

The evolving three-dimensional landscape of human adaptation.

Over the past 3 million years, humans have expanded their ecological niche and adapted to more diverse environments. The temporal evolution and underlying drivers behind this niche expansion remain largely unknown. By combining archeological findings with landscape topographic data and model simulations of the climate and biomes, we show that human sites clustered in areas with increased terrain roughness, corresponding to higher levels of biodiversity. We find a gradual increase in human habitat preferences toward rough terrains until about 1.1 million years ago (Ma), followed by a 300 thousand-year-long contraction of the ecological niche. This period coincided with the Mid-Pleistocene Transition and previously hypothesized ancestral population bottlenecks. Our statistical analysis further reveals that from 0.8 Ma onward, the human niche expanded again, with human species (e.g., H. heidelbergensis, H. neanderthalensis, and H. sapiens) adapting to rougher terrain, colder and drier conditions, and toward regions of higher ecological diversity.

2D Magnetotelluric Resistivity Structure Modeling Using Finite Element Method Based on Vector Triangular Grid and Its Application to Lembang Fault MT Data

Abstract Magnetotellurics is one of the geophysical exploration techniques that relies on the natural fluctuations of electromagnetic waves to delineate their influence on the Earth. The primary focus of this method is to reveal the structure of Earth’s subsurface with the value of resistivity. The application of numerical approaches in magnetotelluric modeling has proven to be an efficient method in various theoretical studies in the field of geophysics, particularly in the context of modeling two-dimensional structures. This research explains a 2D resistivity structure modeling using a vector finite element method. This approach utilizes the edges of elements as vector bases. The presented results include response values such as apparent resistivity and impedance phase at the surface. The study employs the standard model from COMMEMI as a reference to validate the modeling program. Furthermore, the results from this modeling program are compared with the outcomes of the modeling program developed by Weaver et al. The good results were obtained with error values for each model for layered and homogeneous Earth < 3%. Additionally, for the reference model COMMEMI, errors of 3.4393% and 1.4050% were obtained for TE and TM modes, respectively. Furthermore, apparent resistivity and impedance phase results closely approximated the reference values for the topography model. Subsequently, in the application to field data, specifically the Lembang Fault, errors were obtained for the TE and TM modes within the range of 1.16 – 9.16% for each MT data acquisition site.

Successive Tsunamigenic Events Near Sofu Seamount Inferred From High‐Frequency Teleseismic P and Regional T Waves

AbstractAn unexpected major tsunami from the region near Sofu Seamount was observed on 8 October 2023. Sofu Seamount is located approximately 600 km from the coast of Japan. Due to far epicentral distances and the successive occurrence of seismic events, the conventional seismic analysis to reveal the accompanying seismic sequence cannot work well. We investigated high‐frequency teleseismic P and regional T waves from the accompanying seismic sequence during the tsunamigenic events near Sofu Seamount. Envelope shapes of teleseismic P and regional T waves were similar, indicating that T‐wave envelopes also reflected source properties of seismic sequence. During seismic events near Sofu Seamount, observed regional envelopes were characterized by weak body waves and large amplitude T waves with durations of 39–68 s. According to numerical simulations of seismic wave propagation using a realistic topography model, characteristics of T waves exhibit weak slope‐angle and strong source‐depth dependencies. Strong T waves with durations less than 60 s only appeared in results with sources at depths ≤0.5 km below the seafloor. We concluded that high‐frequency radiation of the accompanying seismic sequence during the tsunamigenic events near Sofu Seamount possibly occurred at shallower depths just below the seafloor. If seismic and tsunami sources coincide, shallower source depths might cause tsunamigenic uplifts. The observed peak seafloor uplifts and T‐wave amplitudes during tsunamigenic events were scaled. This result suggests the possibility of tsunami forecasting based on T‐wave amplitudes from submarine volcanoes.

Mitigating Disparate Elevation Differences between Adjacent Topobathymetric Data Models Using Binary Code

Integrating coastal topographic and bathymetric data for creating regional seamless topobathymetric digital elevation models of the land/water interface presents a complex challenge due to the spatial and temporal gaps in data acquisitions. The Coastal National Elevation Database (CoNED) Applications Project develops topographic (land elevation) and bathymetric (water depth) regional scale digital elevation models by integrating multiple sourced disparate topographic and bathymetric data models. These integrated regional models are broadly used in coastal and climate science applications, such as sediment transport, storm impact, and sea-level rise modeling. However, CoNED’s current integration method does not address the occurrence of measurable vertical discrepancies between adjacent near-shore topographic and bathymetric data sources, which often create artificial barriers and sinks along their intersections. To tackle this issue, the CoNED project has developed an additional step in its integration process that collectively assesses the input data to define how to transition between these disparate datasets. This new step defines two zones: a micro blending zone for near-shore transitions and a macro blending zone for the transition between high-resolution (3 m or less) to moderate-resolution (between 3 m and 10 m) bathymetric datasets. These zones and input data sources are reduced to a multidimensional array of zeros and ones. This array is compiled into a 16-bit integer representing a vertical assessment for each pixel. This assessed value provides the means for dynamic pixel-level blending between disparate datasets by leveraging the 16-bit binary notation. Sample site RMSE assessments demonstrate improved accuracy, with values decreasing from 0.203–0.241 using the previous method to 0.126–0.147 using the new method. This paper introduces CoNED’s unique approach of using binary code to improve the integration of coastal topobathymetric data.

An overset-grid finite-difference algorithm for seismic wavefield propagations modelling in the polar coordinate system with a complex free-surface topography

Summary Nowadays, finite-difference method has been widely applied to simulate seismic wavefield propagation in a local seismic model with a complex topography by utilizing curvilinear grids and traction image method in Cartesian coordinates system. For a global seismic model in polar coordinate system, it is still a challenge for the conventional finite-difference method to deal with the grid singularity at the center and the topographic free surface at the top. In this study, we develop a finite-difference method in polar coordinates for seismic wavefield propagation in the 2D global model. In the proposed finite-difference method, the overset-grid algorithm is used to handle the grid singularity at the center, and the curvilinear grid technique as well as the traction image method are applied to implement free-surface boundary condition on the complex topography. The proposed finite-difference method is validated in flat and topographic free-surface models by comparing synthetic waveforms with reference solutions. The seismic wavefield propagation in a realistic Mars profile is computed by the proposed finite-difference method. The proposed finite-difference method is an efficient and accurate method for seismic wave modeling in the polar coordinate system with a complex free-surface topography.

Machine learning to predict morphology, topography and mechanical properties of sustainable gelatin-based electrospun scaffolds

Electrospinning is an outstanding manufacturing technique for producing nano-micro-scaled fibrous scaffolds comparable to biological tissues. However, the solvents used are normally hazardous for the health and the environment, which compromises the sustainability of the process and the industrial scaling. This novel study compares different machine learning models to predict how green solvents affect the morphology, topography and mechanical properties of gelatin-based scaffolds. Gelatin-based scaffolds were produced with different concentrations of distillate water (dH2O), acetic acid (HAc) and dimethyl sulfoxide (DMSO). 2214 observations, 12 machine learning approaches, including Generalised Linear Models, Generalised Additive Models, Generalised Additive Models for Location, Scale and Shape (GAMLSS), Decision Trees, Random Forest, Support Vector Machine and Artificial Neural Network, and a total of 72 models were developed to predict diameter of the fibres, inter-fibre separation, roughness, ultimate tensile strength, Young’s modulus and strain at break. The best GAMLSS models improved the performance of R2 with respect to the popular regression models by 6.868%, and the MAPE was improved by 21.16%. HAc highly influenced the morphology and topography; however, the importance of DMSO was higher in the mechanical properties. The addition of the morphological properties as covariates in the topographic and mechanical models enhanced their understanding.

UAV Modelling of Proglacial Landscapes, Axel Heiberg Island Nunavut

This project involved using a UAV to create digital models and products that enhance visualization and understanding of the moraines of White and Thompson Glaciers, Axel Heiberg Island, Nunavut. UAV surveying offers an effective alternative to high-resolution satellite imagery for documenting highly dynamic proglacial landscapes. Proglacial landforms and sediments provide insights into glacier characteristics and the behavior of glacier termini. UAV surveying allows for high-resolution examination of moraine features and the creation of digital elevation models (DEMs) that facilitate topographic modeling and the detection of changes. A Mavic 3M Drone was employed with UGCS flight planning software. UGCS flight planning software was selected for its ability to create orthomosaic flight plans that consider ground elevation, ensuring accurate Ground Sampling Distance (GSD). A 2m Arctic DEM was uploaded to UGCS prior to planning flight routes, allowing for more precise measurements. Post processing was done in Drone2Map and Agisoft Metashape to create both 2D orthomosaics and 3D models respectively. The results of this work will be used as a visual tool to study current aspects of the moraines, including extent, depositional processes, and present hydrodynamic energy environments. This information will be compared to historical air photos and previous studies to examine changes the moraines over a sixty-year period.

Seafloor Topography Recovery Using the Observation Data of Tiangong-2 Interferometric Imaging Radar Altimeter

Tens of years of altimetry sea surface height (SSH) data has continuously contributed to the development of global seafloor topography models because it can be used to precisely retrieve the gravity variation caused by seafloor topography. The development and application of altimeter has come to a new era from conventional nadir-looking altimeters to interferometric wide-swath altimeters as the successful launch of the Interferometric Imaging Radar Altimeter on board Chinese Tiangong 2 space laboratory (TG2 InIRA). Compared with conventional altimeters, the observation efficiency of TG2 InIRA has been greatly improved. This paper presents for the first time the results of seafloor topography derived from TG2 InIRA data in west Pacific Ocean and South China Sea regions. The root mean squared errors (RMSE) of the recovered topographies by gravity-geologic method (GGM) referring to the shipborne measurements in these two regions are 69.882 and 50.110 m, respectively. The recovered topographies are also evaluated by comparing with five well-known topography models. Furthermore, two subareas in west Pacific Ocean region are selected to conduct a supplemental experiment showing that more revisits can help achieve better recovery. The RMSE of the more revisited subarea is 66.808 m, while the RMSE of the less revisited subarea is 79.273 m.

Study of the Impact of Surface Topography on Wear Resistance

The surface texture parameter is a real reflection of the surface topography, which is closely related to the tribological properties of the surface, and the study of the correlation between the surface texture parameter and wear resistance is of great significance in revealing the tribological influence mechanism of the surface and realising the functional manufacturing of the surface. This paper takes the ball-end milling surface as the research object, establishes the three-dimensional simulation model of the surface topography, and analyses the surface topography and the surface texture parameter change rule. Based on the improved correlation analysis model, the correlation characteristics between the surface texture parameters, and between the surface texture parameters and the relative wear rate of per unit sliding distance KV, were investigated, and the prediction model of KV was established based on the surface texture parameters Sku, Sa, Sxp, Sp, and Ssk, and the correctness of the model was verified by experiments. The study in this paper provides a new idea to further reveal the relationship between surface topographical features and wear resistance and to guide the functional manufacturing of surfaces.

A New Precise Point Positioning with Ambiguity Resolution (PPP-AR) Approach for Ground Control Point Positioning for Photogrammetric Generation with Unmanned Aerial Vehicles

Unmanned aerial vehicles (UAVs) are now widely preferred systems that are capable of rapid mapping and generating topographic models with relatively high positional accuracy. Since the integrated GNSS receivers of UAVs do not allow for sufficiently accurate outcomes either horizontally or vertically, a conventional method is to use ground control points (GCPs) to perform bundle block adjustment (BBA) of the outcomes. Since the number of GCPs to be installed limits the process in UAV operations, there is an important research question whether the precise point positioning (PPP) method can be an alternative when the real-time kinematic (RTK), network RTK, and post-process kinematic (PPK) techniques cannot be used to measure GCPs. This study introduces a novel approach using precise point positioning with ambiguity resolution (PPP-AR) for ground control point (GCP) positioning in UAV photogrammetry. For this purpose, the results are evaluated by comparing the horizontal and vertical coordinates obtained from the 24 h GNSS sessions of six calibration pillars in the field and the horizontal length differences obtained by electronic distance measurement (EDM). Bartlett’s test is applied to statistically determine the accuracy of the results. The results indicate that the coordinates obtained from a two-hour PPP-AR session show no significant difference from those acquired in a 30 min session, demonstrating PPP-AR to be a viable alternative for GCP positioning. Therefore, the PPP technique can be used for the BBA of GCPs to be established for UAVs in large-scale map generation. However, the number of GCPs to be selected should be four or more, which should be homogeneously distributed over the study area.

Climate‐Driven Projections of Future Global Wetlands Extent

AbstractWetlands are crucial components of the Earth's system, interacting with various processes such as the hydrological cycle, energy exchanges with the atmosphere, and global nitrogen and carbon cycles. The future trajectory of wetlands is anticipated to be influenced not only by direct human activities, but also by climate change. Here we present our assessment of climate‐driven global changes in wetlands extent, focusing on the main wetland complexes. We used an approach based on the Topographic Hydrological model (TOPMODEL) and soil liquid water content projections from 14 models of the Coupled Model Intercomparison Project Phase 6 (CMIP6). Our analysis reveals a consistent decrease in wetlands extent in the Mediterranean, Central America, and northern South America, with a substantial loss of 28% in the western Amazon Basin for the end of the 21st century (2081–2100) under the SSP370 scenario. Conversely, Central Africa exhibits an increase in wetlands extent, except in the Congo Basin. Nevertheless, most of the areas studied (80%) present uncertain results, due to conflicting projections of changes between the models. Notably, we show that there is significant uncertainty among CMIP6 models regarding liquid soil water content in high latitudes. By narrowing our focus to 10 models, which seem to better represent the thawing of permafrost, we obtain a better inter‐model agreement. We then find a modest declines in the overall global area (<5%), but an average loss of 13% beyond 50°N. Specific areas like the Hudson Bay Lowlands experiencing a 21% decrease and the Western Siberian Lowlands a 15% decrease.

Formation Technique of a 3D Ground Topographic Numerical Model for Environmental Impact Prediction using Digital GIS Databases, Contour Lines and Constrained Delaunay Triangulation

Formation Technique of a 3D Ground Topographic Numerical Model for Environmental Impact Prediction using Digital GIS Databases, Contour Lines and Constrained Delaunay Triangulation

The Characteristics and Variation of the Golden Eagle Aquila chrysaetos Home Range

Satellite tracking allows for novel investigations into golden eagle home range characteristics. Understanding home range characteristics is important for conservation and for assessing the potential impact of landscape changes from forest planting, wind farms, etc. Small sample sizes, inconsistent definitions and methods restricted several previous studies. Our study involved 69 resident tagged eagles with over one year of data across five Scottish regions. Home range size was estimated from 95% isopleth contours extracted from Utilisation Distributions. Above a small threshold, estimated range size was not affected by the number of records but at least one year of data is required, largely because of the breeding and non-breeding seasonal differences. There were no significant range size differences between birds tagged as range holders and those previously tagged as nestlings. Across four regions, with considerable intra-regional variation, planar 95% isopleths did not differ (medians, km2): Argyll 58.9, Northwest Highlands 61.7, Northeast Highlands 89.3, South of Scotland 91.9. Ranges in the isolated Outer Hebrides region were exceptionally small, at 24.0 km2. Estimated range area was usually reduced to 70–80% of the planar area when restricted to usable habitat, as estimated by the Golden Eagle Topography (GET) model. Applying measures of known unsuitable habitat (closed-canopy commercial forest and wind turbines) further reduced usable open land. Loss of otherwise suitable habitat was substantially due to commercial forest. Larger ranges had larger extents of suitable habitat (according to GET), with no apparent optimum of preferred GET habitat. Range size was not different across a year between the sexes. Breeding ranges were smaller, and females’ breeding ranges were much smaller than those of males, but larger than males’ ranges in the non-breeding season. Breeding attempt duration was probably also influential. Our study provides novel insights into golden eagle home range characteristics and can guide further research and practical applications.

The Potential of Hybrid Data Integration Approaches Based UAV-SfM Photogrammetry in Identifying New Archaeological Features in The Archaeological UR Site

Structure-from-Motion (SfM) Photogrammetry-based unmanned aerial vehicles (UAV) are considered one of the most effective non-destructive techniques to detect features in archaeology. However, relying on only standalone-based photogrammetric approaches is not always valid in revealing new archaeological features and monuments. This study aims to verify the efficiency of the integration-based methodologies of SfM photogrammetry products to improve the detection of new archaeological features in the non-excavated regions of the archaeological site of Ur, south of Iraq. The methodology proposes an automatic workflow to extract 3D topographic models such as Digital Terrain Model (DTM) by Agisoft Metashape (MS) software and orthomosaic images as standard products for UAV-photogrammetry. These have delivered the following Structure-from-Motion and Multi-View Stereo (SfM-MVS) algorithms using an optimized flight planning, data capturing, and data processing pipeline. Stand-alone raster images derived from photogrammetry, such as (slope images, hill shade rasters, etc.), have limitations in archaeology compared to hybrid raster images in detecting covered and unrevealed features. The methodology focuses on fusing multiple rasters extracted from DTM using the QGIS software package analysis tools to detect features that may never be detected from single standalone image rasters. This is because combining multiple rasters of different parameters may highlight new potential features such as paths, temples, etc. The results findings revealed that the fusion methodology of raster images obtained from UAV-SfM Photogrammetry through its standard products and using high precision settings could successfully help to identify many new archaeological features in the undiscovered site in Ur, such as walls and houses, which may date back to the late periods of Ur. The results concluded that this enhances the role of photogrammetry and remote sensing techniques as non-destructive excavation tools to provide a rich database for archaeologists for further analysis and interpretation of undiscovered archaeological sites without destruction. Many new archaeological features were also revealed, and a general digital map of the site was prepared. It also prompts a review of the settlement theories in previous studies in the ancient city of UR and the documentation of excavated archaeological ruins.

Comparative Study of Seafloor Topography Prediction from Gravity–Geologic Method and Analytical Algorithm

Seafloor topography prediction can fill in sea areas without ship sounding data. However, the dependence of various topographic prediction algorithms on ship soundings varies significantly. Hence, this study explores the impact of the number and distributions of ship soundings on topographic prediction using the gravity–geologic method (GGM) and an analytical algorithm. Firstly, this study investigates the influence of ship sounding coverage on the two algorithms. The simulation results demonstrate that increasing coverage from 5.40% to 31.80%, coupled with more uniform distributions across the study area, substantially reduces the RMS error of the GGM. Specifically, the RMS error decreases from 238.68 m to 42.90 m, an improvement of 82.03%. The analytical algorithm maintains a consistent RMS error of 40.39 m because it does not depend on ship soundings. Furthermore, we select a 1° × 1° sea area (134.8°–135.8°E, 30.0°–31.0°N), and the ship soundings are divided into two control groups, Part I and Part II, with coverages of 8.19% and 33.19%, respectively. When Part II is used for calculation, the RMS error of the GGM decreases from 204.17 m to 126.95 m compared to when Part I is used, while the analytical algorithm exhibits an RMS error of 167.94 m. The findings indicate that the prediction accuracy of the GGM is significantly affected by ship soundings, whereas the analytical algorithm is more stable and independent of ship soundings. Based on simulation experiments and realistic examples, when the effective ship soundings coverage exceeds 30%, the GGM may have more advantages. Conversely, the analytical algorithm may be better. This suggests that effectively combining and utilizing different algorithms based on the ship sounding coverage can improve the accuracy of topographic prediction. This will provide a basis for integrating multiple algorithms to construct a global seafloor topography model.

Refined Bathymetric Prediction Based on Feature Extraction of Gravity Field Signals: BATHY‐FE

AbstractThe resolution of altimetry‐derived gravity signals renders traditional bathymetry inversion methods inadequate for detecting short‐wavelength bathymetric features. This study presents a new global seafloor topography model by merging regional models constructed via a deep learning approach recently developed by the authors. Feature extraction techniques were used to refine each regional bathymetric model, resulting in sharpened features that are fuzzily revealed (or not seen) in global bathymetric models such as GEBCO_2023 and SRTM15+V2, especially in areas close to the poles. This is proof that, with regard to current gravity field products, conventional bathymetry inversion methods are weaker at generalizing to uncharted locations of the global ocean floor. At test points, the global seafloor model has mean error and error standard deviation of 1.75 ± 81.15 m. The refined seafloor topography can be used to supplement existing seafloor maps, especially in the polar regions. It is useful for studying seafloor geography, and also for studies in which seafloor ruggedness is required.

Application of active MASW study on irregular and undulated sloping terrain: a case study

ABSTRACT Multichannel Analysis of Surface Waves (MASW) is a geophysical technique used for deciphering subsurface profiling by the determination of shear wave velocity. It is increasingly utilized in sloping terrain; however, its accuracy is limited by the presence of undulating or irregular topography. The study analysed the impact of undulating topography on sloping terrain using models of depressed and elevated topography. The surface waves lose more energy as depth, width, and steepness of the undulated topography increased. The study revealed that the elevated topography had a greater influence than depressed topography. The undulating topography in close proximity to the source resulted in a decrease in resolution dispersion, however the terrain located further away did not have a substantial impact. The study aims to identify suitable configurations that must be taken into account while identifying the subsurface characterization of a sloping terrain using the MASW.

Topography Modeling of Surface Grinding Based on Random Abrasives and Performance Evaluation

The surface morphology and roughness of a workpiece are crucial parameters in grinding processes. Accurate prediction of these parameters is essential for maintaining the workpiece’s surface integrity. However, the randomness of abrasive grain shapes and workpiece surface formation behaviors poses significant challenges, and accuracy in current physical mechanism-based predictive models is needed. To address this problem, by using the random plane method and accounting for the random morphology and distribution of abrasive grains, this paper proposes a novel method to model CBN grinding wheels and predict workpiece surface roughness. First, a kinematic model of a single abrasive grain is developed to accurately capture the three-dimensional morphology of the grinding wheel. Next, by formulating an elastic deformation and formation model of the workpiece surface based on Hertz theory, the variation in grinding arc length at different grinding depths is revealed. Subsequently, a predictive model for the surface morphology of the workpiece ground by a single abrasive grain is devised. This model integrates the normal distribution model of abrasive grain size and the spatial distribution model of abrasive grain positions, to elucidate how the circumferential and axial distribution of abrasive grains influences workpiece surface formation. Lastly, by integrating the dynamic effective abrasive grain model, a predictive model for the surface morphology and roughness of the grinding wheel is established. To examine the impact of changing the grit size of the grinding wheel and grinding depth on workpiece surface roughness, and to validate the accuracy of the model, experiments are conducted. Results indicate that the predicted three-dimensional morphology of the grinding wheel and workpiece surfaces closely matches the actual grinding wheel and ground workpiece surfaces, with surface roughness prediction deviations as small as 2.3%.

A hp-adaptive finite element approach for 3D controlled source electromagnetic forward modeling

Accurate and efficient 3D controlled-source electromagnetic (CSEM) forward modeling is crucial for deciphering CSEM responses in complex geological contexts and for assessing the efficacy of survey designs. It also serves as the foundation for 3D inversion, a critical technique for interpreting CSEM data. We present a novel hp-adaptive finite element method for 3D CSEM forward modeling. This method enhances the traditional h-adaptive approach, which only refines the mesh while keeping the polynomial degree of the basis functions constant, by enabling simultaneous refining mesh and increasing polynomial degree of basis functions. High-order elements yield superior approximations when the solution exhibits sufficient smoothness, while mesh refinement enhances accuracy in regions with less smooth solutions. The adaptive refinement process is guided by a goal-oriented error estimator that relies on the residual of the governing equation and the continuity of electromagnetic fields. The hp-adaptive strategy, which determines whether to refine elements or increase their polynomial degree, is based on the refinement levels of the elements. Constraints are applied to hanging nodes, introduced during the refinement process, to ensure the continuity of the finite element space. We validate the accuracy, efficiency, and convergence of the hp-adaptive method through tests on two layered models. Further simulations on a topographic model and a realistic model corroborate its effectiveness. Our findings indicate that the hp-adaptive method outperforms the traditional h-adaptive method in terms of accuracy, convergence rate and computational cost, making it a promising technique for 3D CSEM forward modeling.

Formation of orientation-dependent surface morphologies on tungsten after low energy helium plasma exposure: multiscale characterization and new insights

In ITER, the helium (He) impurity produced by the deuterium-tritium reaction will bombard the tungsten (W) divertor armor at the strike points. Consequently, strong interaction occurs therein that both impact the performance of the plasmas and the lifetime of the divertor. Despite an ever-increasing understanding of this interaction, some experimental phenomena remain mysterious, especially the formation of orientation-dependent surface morphologies. Here, we combine multiscale experimental characterization and theoretical models to shed new light on this problem. After low-energy He plasma exposure in a linear plasma generator, the polycrystalline W surface developed various morphologies. Through electron backscatter diffraction analysis, we found that the {111} grains developed cube-corner structures, the {110} grains developed ripple structures, whereas the {100} grains remained smooth. Then, electron-transparent lamellae were extracted from such grains to observe the subsurface He bubbles by transmission electron microscopy. The volume density, size distribution, and depth range of the He bubbles weakly depend on the crystallographic orientation, suggesting that the migration of W atoms causes the morphology variety. Accordingly, we proposed a two-stage formation mechanism. First, W atoms generated by over-pressurized He bubbles glide on the slip plane and in the slip direction to reach the surface, forming characteristic patterns that are enclosed by the slip traces. Second, morphological instability drives the evolution of the surface patterns, in which the initial surface structure and surface self-diffusion kinetics mediate. The proposed mechanism has been incorporated into a topographical instability model to enable asemi-quantitative analysis. The obtained new insights are valuable to the impurity control of the core plasmas and the lifetime analysis of the divertor for ITER.