Feed aggregator

Some Features of Interacting Solar Wind Disturbances

Geomagnetism and Aeronomy - Thu, 08/01/2024 - 00:00
Abstract

The updated database of Forbush effects and interplanetary disturbances (https://tools.izmiran.ru/feid) is used for an extensive analysis of various characteristics of events caused by the influence of interacting solar wind disturbances on the near-Earth space. In particular, the cases of different combinations of the pair interaction of high-speed streams from coronal holes and coronal mass ejections over the long period from 1995 to 2022 are considered. Variations in the flux of galactic cosmic rays (with a rigidity of 10 GV) and changes in the parameters of the interplanetary medium and geomagnetic activity are described. It is shown that the degree of mutual influence depends on the time between the detection of neighboring events; with the most pronounced changes in various parameters for events whose interaction occurred before reaching Earth’s orbit. It is also established that in interacting solar wind disturbances not only the extremes of the parameters of cosmic rays, interplanetary medium, and geomagnetic activity but also their time profiles are subject to changes.

Long-Term Trends in the Height of the Ionospheric F2 Layer Peak

Geomagnetism and Aeronomy - Thu, 08/01/2024 - 00:00
Abstract

Long-term variations (trends) in the height of the ionospheric F2 layer peak hmF2 is analyzed based on the data of Moscow and Juliusruh stations. The near-noon LT hours and two winter months (January and February) and two summer months (June and July) are considered for a period of 1996–2023. Well-pronounced and statistically significant negative hmF2 trends are found both in summer and winter. Overall, the F2 layer height decreased during the analyzed period by 0.5–1 km per year. The “Delta” method developed and published by the authors earlier is applied to the same data. The results confirm a systematic decrease in the hmF2 value in the past two decades. It is found that the F2 layer height has decreased in recent years more rapidly than in the earlier years.

HASPPP: an open-source Galileo HAS embeddable RTKLIB decoding package

GPS Solutions - Thu, 08/01/2024 - 00:00
Abstract

Galileo supports global precise point positioning (PPP) service by delivering precise products via satellites, which is known as high accuracy service (HAS). To improve reception efficiency, the Galileo HAS employs a high-parity vertical reed-solomon encoding scheme, increasing the complexity of HAS corrections recovery. To promote research and application of the HAS, an open-source C/C + + decoding package, HASPPP, has been developed for seamless embedding into the prevalent C/C++-based software, such as RTKLIB. HASPPP provides interfaces for effortless decoding support of raw HAS binary data from various manufacturers. The HASPPP manual offers code analysis along with simple to complex examples to aid users in swiftly mastering decoding. HASPPP has undergone rigorous validation of its decoding accuracy and reliability.

Compact RTK for expanded area (COREA): a new method for carrier-phase-based satellite augmentation system

GPS Solutions - Thu, 08/01/2024 - 00:00
Abstract

This study proposes a new concept of carrier-phase-based satellite augmentation system named “Compact Real-time Kinematic for Expanded Area (COREA),” which provides centimeter-level positioning services across nationwide coverage. The proposed system’s architecture is very similar to that of the satellite-based augmentation system (SBAS), a meter-level aviation safety service. While network real-time kinematic and precise-point-positioning-RTK (PPP-RTK) rely on several densely positioned reference stations, COREA provides carrier-phase-based corrections using a few reference stations with a distance of 400–1000 km. Furthermore, the COREA corrections can be transmitted by satellite signals with extremely low-speed data links of 250 bps, similar to SBAS. This study focused on the generation method for satellite code/phase clock (CPC) corrections, which is the most significant part of the system. We analyzed the user performance of the COREA system constructed in the Midwest and South of the United States with six reference stations. Consequently, satellite CPC corrections are resilient to communication failures and highly accurate in identifying user integer ambiguity. The 95% position accuracy is approximately 1.8 cm horizontally and 7.1 cm vertically, with an average convergence time of 1–3 min using only GPS triple-frequency signals. In summary, the COREA system preserves the hardware architecture of the legacy SBAS while providing centimeter-level services with fast convergence times by utilizing extremely low-speed satellite data links across the country.

Dynamic evolution of West Gondwana inferred from crustal anisotropy of the South American platform

Geophysical Journal International - Thu, 08/01/2024 - 00:00
SummaryThe amalgamation and breakup of the West Gondwana shaped the South American platform. The dynamics during the processes can be reflected by crust anisotropy of the platform, but there are no specialized crustal anisotropic measurements yet. Splitting analysis of Moho-converted shear waves in P-wave receiver functions (Pms) can reveal crustal-scale anisotropy, which is important for understanding the dynamic evolution of the crust and for the interpretation of mantle anisotropy from splitting analysis of core–mantle refracted shear waves (XKS phases). This study measured crustal anisotropy for the old and stable South American platform by Pms splitting analysis. The splitting times vary mainly in the range of 0–0.5 s, with a regional mean of 0.2 s, slightly lower than that observed in tectonically active regions. The detected crustal anisotropy shows distinct characteristics and spatial zoning, providing insights into tectonic processes. (1) Fast polarization directions at stations close to the Transbrasiliano Lineament (TBL) are oriented NNE–SSW, generally consistent with the strike of the TBL but inconsistent with the maximum horizontal compressive stress, implying that they might be formed by dynamic metamorphism during the formation of the TBL. (2) Crustal anisotropy along the passive continental margin in the east and northeast is weak. Still, the fast polarization directions tend to be oriented along the margin, implying the existence of fossil extensional crustal fabrics formed during the continental rifting of West Gondwana. (3) The Paraná Basin, one of the world's largest Large Igneous Provinces (LIP) covered by continental flood basalts, shows distinctively strong anisotropy, with fast polarization directions highly aligned with mantle anisotropy, implying that synchronous crust–mantle deformation occurred in these regions as a result of magmatism during the breakup of West Gondwana.

A comparative investigation of geomagnetic jerks across the SAA during the period 2000-2020

Geophysical Journal International - Thu, 08/01/2024 - 00:00
SummaryGeomagnetic field data from six magnetic observatories in and adjacent to the South Atlantic Anomaly were individually analysed to detect abrupt secular variation changes occurring on time scales of less than a year and to explore any correlation with the evolution of the South Atlantic Anomaly. After applying external field corrections by means of the CHAOS-7 model, twelve-month differences of the respective observatory monthly mean of the eastward component Y revealed evidence of several geomagnetic jerks with varying amplitudes during the period 2000-2020. These observations were compared to the CHAOS-7 spherical harmonic model and previous studies of the South Atlantic Anomaly’s evolution. It emerged from this study that global field models like CHAOS are not very effective in identifying rapid localised changes in the geomagnetic field, highlighting the importance of using observatory data in conjunction with satellite data when studying geomagnetic jerks.

Interpolated Fast Damped Multichannel Singular Spectrum Analysis for Deblending of Off-the-Grid Blended Data

Surveys in Geophysics - Thu, 08/01/2024 - 00:00
Abstract

Blended acquisition offers significant cost and period reduction in seismic data acquisition. However, fired blended sources are usually deployed at off-the-grid (OffG) samples due to obstacle limitation and economic cost considerations. The irregular distribution of coordinates, along with the blending noise, has a detrimental effect on the performance of subsequent seismic processing and imaging. The interpolated multichannel singular spectrum analysis (I-MSSA) algorithm effectively provides on-the-grid deblended results by employing an interpolator, in conjunction with a projected gradient descent strategy. However, the deblending accuracy and computational efficiency of the I-MSSA are still a concern due to the limitations of the traditional singular value decomposition (SVD). To address these limitations, we propose an interpolated fast damped multichannel singular spectrum analysis (I-FDMSSA) rank-reduction algorithm. The proposed algorithm incorporates the damping operator, the randomized SVD (RSVD) and the fast Fourier transform (FFT) strategy. The damping operator can further attenuate the remaining noise in the estimated signal obtained from the truncated SVD, resulting in an improved deblending performance. The RSVD accelerates the rank-reduction process by shrinking the size of the Hankel matrix. To expedite the rank-reduction and anti-diagonal averaging stages without explicitly constructing large-scale block Hankel matrices, the FFT strategy is employed. By incorporating a 2D separable sinc interpolator, the I-FDMSSA enables an efficient and accurate deblending of 3D OffG blended data. The deblending performance and operational efficiency improvements of the proposed I-FDMSSA algorithm over the traditional I-MSSA algorithm are demonstrated through OffG synthetic and field blended data examples.

Categories:

Constructing Priors for Geophysical Inversions Constrained by Surface and Borehole Geochemistry

Surveys in Geophysics - Thu, 08/01/2024 - 00:00
Abstract

Prior model construction is a fundamental component in geophysical inversion, especially Bayesian inversion. The prior model, usually derived from available geological information, can reduce the uncertainty of model characteristics during the inversion. However, the prior geological data for inferring a prior distribution model are often limited in real cases. Our work presents a novel framework to create 3D geophysical prior models using soil geochemistry and borehole rock sample measurements. We focus on the Bayesian inversion, which enables encoding of knowledge and multiple non-geophysical data into the prior. The new framework developed in our research comprises three main parts, namely correlation analysis, prior model reconstruction, and Bayesian inversion. We investigate the correlations between surface and subsurface geochemical features, as well as the correlation between geochemistry and geophysics, using canonical correlation analysis for the surface and borehole geochemistry. Based on the resulting correlations, we construct the prior susceptibility model. The informed prior model is then tested using geophysical forward modeling and outlier detection methods. In this test, we aim to falsify the prior model, which happens when the model cannot predict the field geophysical observation. To obtain the posterior models, the reliable prior models are incorporated into a Bayesian inversion framework. Using a real case of exploration in the Central African Copperbelt, we illustrate the workflow of constructing the high-resolution 3D stratigraphic model conditioned on soil geochemistry, borehole data, and airborne geophysics.

Categories:

Stress-Dependent PP-Wave Reflection Coefficient for Fourier-Coefficients-Based Seismic Inversion in Horizontally Stressed Vertical Transversely Isotropic Media

Surveys in Geophysics - Thu, 08/01/2024 - 00:00
Abstract

The subsurface in situ stress fields significantly influence the elastic and anisotropic properties of rocks, yet traditional linear elastic theories often overlook the impact of stress on seismic response characteristics. Nonlinear acoustoelastic theory integrates third-order elastic constants (TOECs) to elucidate the influence of stress on changes in elastic and anisotropic properties of stressed rocks. A comprehensive examination of recent scholarly investigations on nonlinear acoustoelastic phenomena precedes the introduction of an innovative stress-dependent equation for the PP-wave reflection coefficient. This equation delineates the dependence of azimuthal seismic response on horizontal uniaxial stress in inherently vertical transversely isotropic (VTI) media, or those VTI formations induced by a single set of horizontal aligned fractures. Emphasis is placed on delineating stress-induced anisotropy and elucidating azimuthal PP-wave reflection characteristics in horizontally uniaxially stressed VTI media. Additionally, this discourse extends to more intricate scenarios involving horizontally biaxially and triaxially stressed VTI media, as delineated by nonlinear acoustoelastic theory. Subsequently, the reflection coefficient of horizontally uniaxially stressed VTI media is expressed in terms of azimuthal Fourier coefficients (FCs), revealing that the unstressed VTI background exhibits heightened sensitivity to zeroth-order FC, while the stress-induced anisotropy manifests greater sensitivity to second-order FC. Through the application of azimuthal FCs-based amplitude versus offset and azimuth (AVOAz) inversion method to both synthetic and field datasets, the proposed model and approach offer promising avenues for reservoir characterization in VTI media subject to horizontal uniaxial stress conditions.

Categories:

Multiscalar Integration of Dense and Sparse Spatial Data: an Archaeological Case Study with Magnetometry and Geochemistry

Surveys in Geophysics - Thu, 08/01/2024 - 00:00
Abstract

Integration of different kinds of data is an important issue in archaeological prospection. However, the current methodological approaches are underdeveloped and rarely use the data to their maximum potential. Common approaches to integration in the geophysical sciences are mostly just various forms of comparison. We argue that true integration should involve the mathematical manipulation of input data such that the original values of the input data are changed, or that new variables are produced. To address this important research gap, we present an innovative approach to the analysis of geochemical and geophysical datasets in prospection-focused disciplines. Our approach, which we refer to as “multiscalar integration” to differentiate it from simpler methods, involves the application of mathematical methods and tools to process the data in a unified way. To demonstrate our approach, we focus on integrating geophysical data (magnetometry) with geochemical data (elemental content). Our approach comprises three main stages: Quantification of the data deviation from random distributions, linear modelling of geophysical and geochemical data and integration based on weighting of the different elements derived in previous steps. All the steps of the workflow can be also applied separately and independently as needed or preferred. Our approach is implemented in the R environment for statistical computing. All data, functions and scripts used in the work are available from open access repositories (Zenodo.org and Github.com) so that others can test, modify and apply our proposed methods to new cases and problems. Our approach has the following advantages: (1) It allows the rapid exploration of multiple data sources in an unified way; (2) it can increase the utility of geochemical data across diverse prospection disciplines; (3) it facilitates the identification of links between geochemical and geophysical data (or generally, between point-based and raster data); (4) it innovatively integrates various datasets by weighting the information provided by each; (5) it is simple to apply following a step-by-step framework; (6) the code and workflow is fully open to allow for customization, improvements and additions.

Categories:

Meta-Processing: A robust framework for multi-tasks seismic processing

Surveys in Geophysics - Thu, 08/01/2024 - 00:00
Abstract

Machine learning-based seismic processing models are typically trained separately to perform seismic processing tasks (SPTs) and, as a result, require plenty of high-quality training data. However, preparing training data sets is not trivial, especially for supervised learning (SL). Despite the variability in seismic data across different types and regions, some general characteristics are shared, such as their sinusoidal nature and geometric texture. To learn the shared features and thus, quickly adapt to various SPTs, we develop a unified paradigm for neural network-based seismic processing, called Meta-Processing, that uses limited training data for meta learning a common network initialization, which offers universal adaptability features. The proposed Meta-Processing framework consists of two stages: meta-training and meta-testing. In the former, each SPT is treated as a separate task and the training dataset is divided into support and query sets. Unlike conventional SL methods, here, the neural network (NN) parameters are updated by a bilevel gradient descent from the support set to the query set, iterating through all tasks. In the meta-testing stage, we also utilize limited data to fine-tune the optimized NN parameters in an SL fashion to conduct various SPTs, such as denoising, interpolation, ground-roll attenuation, image enhancement, and velocity estimation, aiming to converge quickly to ideal performance. Extensive numerical experiments are conducted to assess the effectiveness of Meta-Processing on both synthetic and real-world data. The findings reveal that our approach leads to a substantial improvement in the convergence speed and predictive performance of the NN.

Categories:

Detectability of Seamount Eruptions Through a Quantum Technology Gravity Mission MOCAST+: Hunga Tonga, Fani Maoré and Other Smaller Eruptions

Surveys in Geophysics - Thu, 08/01/2024 - 00:00
Abstract

Seamount eruptions alter the bathymetry and can occur undetected due to lack of explosive character. We review documented eruptions to define whether they could be detected by a future satellite gravity mission. We adopt the noise level in acquisitions of multi-satellite constellations as in the MOCAST+ study, with a proposed payload of a quantum technology gradiometer and clock. The review of underwater volcanoes includes the Hunga Tonga Hunga Ha’apai (HTHH) islands for which the exposed surface changed during volcanic unrests of 2014/2015 and 2021/2022. The Fani Maoré submarine volcanic eruption of 2018–2021 produced a new seamount 800 m high, emerging from a depth of 3500 m, and therefore not seen above sea surface. We review further documented submarine eruptions and estimate the upper limit of the expected gravity changes. We find that a MOCAST+ type mission should allow us to detect the subsurface mass changes generated by deep ocean submarine volcanic activity for volume changes of 6.5 km3 upwards, with latency of 1 year. This change is met by the HTHH and Fani Maoré volcanoes.

Categories:

Pressure Effects on Plane Wave Reflection and Transmission in Fluid-Saturated Porous Media

Surveys in Geophysics - Thu, 08/01/2024 - 00:00
Abstract

The wave reflection and transmission (R/T) coefficients in fluid-saturated porous media with the effect of effective pressure are rarely studied, despite the ubiquitous presence of in situ pressure in the subsurface Earth. To fill this knowledge gap, we derive exact R/T coefficient equations for a plane wave incident obliquely at the interface between the dissimilar pressured fluid-saturated porous half-spaces described by the theory of poro-acoustoelasticity (PAE). The central result of the classic PAE theory is first reviewed, and then a dual-porosity model is employed to generalize this theory by incorporating the impact of nonlinear crack deformation. The new velocity equations of generalized PAE theory can describe the nonlinear pressure dependence of fast P-, S- and slow P-wave velocities and have a reasonable agreement with the laboratory measurements. The general boundary conditions associated with membrane stiffness are used to yield the exact pressure-dependent wave R/T coefficient equations. We then model the impacts of effective pressure on the angle and frequency dependence of wave R/T coefficients and synthetic seismic responses in detail and compare our equations to the previously reported equations in zero-pressure case. It is inferred that the existing R/T coefficient equations for porous media may be misleading, since they lack consideration for inevitable in situ pressure effects. Modeling results also indicate that effective pressure and membrane stiffness significantly affect the amplitude variation with offset characteristics of reflected seismic signatures, which emphasizes the significance of considering the effects of both in practical applications related to the observed seismic data. By comparing the modeled R/T coefficients to the results computed with laboratory measured velocities, we preliminarily confirm the validity of our equations. Our equations and results are relevant to hydrocarbon exploration, in situ pressure detection and geofluid discrimination in high-pressure fields.

Categories:

Joint Inversion Method of Gravity and Magnetic Data with Adaptive Zoning Using Gramian in Both Petrophysical and Structural Domains

Surveys in Geophysics - Thu, 08/01/2024 - 00:00
Abstract

Different observation data are utilized to obtain a unified geophysical model based on the correlations of underground geological bodies in joint inversions. By specifying a type of Gramian constraints, Gramian as a coupling term can link geophysical models through relationships of physical properties or structural similarities. Considering the complex relationships of physical properties of underground geological bodies, we proposed an adaptive zoning method to automatically divide the whole inversion area into subregions with different relationships of physical properties and to determine the number and range of subregions that utilized correlation between geophysical data before joint inversions. On this basis, we considered the use of a combination of Gramian coupling terms rather than one term to link petrophysical and structural domains during joint inversions. Synthetic tests showed that the algorithm is capable of having a robust estimate of the spatial distribution and relationships between density and magnetization intensity of geological bodies. The idea was also applied to the ore concentration area in the middle and lower reaches of the Yangtze River to obtain the three-dimensional (3-D) distribution model of magnetite-bearing rocks within 5 km underground, which corresponds well with the existing shallow ore sites and demonstrates the existence of available deep resources in the study area.

Categories:

Photonic Seismology: A New Decade of Distributed Acoustic Sensing in Geophysics from 2012 to 2023

Surveys in Geophysics - Thu, 08/01/2024 - 00:00
Abstract

This paper delivers an in-depth bibliometric analysis of distributed acoustic sensing (DAS) research within the realm of geophysics, covering the period from 2012 to 2023 and drawing on data from the Web of Science. By employing bibliographic and structured network analysis methods, including the use of Bibliometrix and VOSviewer®, the study highlights the most influential scholars, leading institutions, and pivotal research contributions that have significantly shaped the field of DAS in geophysics. The research delves into key collaborative dynamics, unraveling them through co-authorship network analysis, and delves into thematic developments and trajectories via comprehensive co-citation and keyword co-occurrence network analyses. These analyses elucidate the most robust and prominent areas within DAS research. A critical insight gained from this study is the rise of ‘photonic seismology’ as an emerging interdisciplinary domain, exemplifying the fusion of photonic sensing techniques with seismic science. This paper also discusses certain limitations inherent in the study and concludes with implications for future research.

Categories:

Feasibility Study of Anisotropic Full-Waveform Inversion with DAS Data in a Vertical Seismic Profile Configuration at the Newell County Facility, Alberta, Canada

Surveys in Geophysics - Thu, 08/01/2024 - 00:00
Abstract

As an emerging seismic acquisition technology, distributed acoustic sensing (DAS) has drawn significant attention in earth science for long-term and cost-effective monitoring of underground activities. Field seismic experiments with optical fibers in a vertical seismic profile (VSP) configuration were conducted at the Newell County Facility of Carbon Management Canada in Alberta, Canada, for \({\text{CO}}_2\) injection and storage monitoring. Seismic full-waveform inversion (FWI) represents one promising approach for high-resolution imaging of subsurface model properties. In this study, anisotropic FWI with variable density is applied to the DAS-recorded walk-away VSP data for characterizing the subsurface velocity, anisotropy, and density structures, serving as baseline models for future time-lapse studies at the pilot site. Synthetic inversion experiments suggest that, without accounting for anisotropy, the inverted density structures by isotropic FWI are damaged by strong trade-off artifacts. Anisotropic FWI can provide more accurate P-wave velocity, density, and valuable anisotropy models. Field data applications are then performed to validate the effectiveness and superiority of the proposed methods. Compared to the inversion outputs of isotropic FWI, the inverted P-wave velocity by anisotropic FWI matches trend variation of the well log more closely. In the inverted density model, the \({\text{CO}}_2\) injection formation can be clearly resolved. The inverted anisotropy parameters provide informative references to interpret the structures and lithology around the target \({\text{CO}}_2\) injection zone.

Categories:

Multiple Satellite Observations of the High‐Latitude Cusp Aurora During Northward IMF Conditions

JGR:Space physics - Wed, 07/31/2024 - 20:29
Abstract

Cusp auroras poleward of the typical auroral oval are ascribed to high-latitude lobe reconnection when the Interplanetary Magnetic Field (IMF) B z is predominantly northward. In this study, we further investigate the ionospheric characteristics of a unique high-latitude cusp region employing multiple satellite observations. A cusp aurora event with wide spatial spread was observed in the postnoon polar cap region. It was found to be associated with northward IMF B z and positive B y components. The cusp aurora was located from 68° to 86° in magnetic latitude and within 15–17 hr in magnetic local time. This broad coverage in the polar cap indicates direct precipitating particles from the magnetosheath. Particle energy is different between the equatorward and poleward edges of the cusp aurora. The precipitating ions at the equatorward side maintain magnetosheath particle characteristics as expected, while ions with higher energies occurred in the poleward side. Further, the poleward edge of the cusp aurora was nearly situated in the center of a convection shear and was associated with an upward field-aligned current. These observations suggest a lobe cell circulation, hence we attribute the formation of the cusp aurora to the high-latitude lobe reconnection. Simultaneous observations in the southern hemisphere indicate the absence of cusp aurora. The auroral presence only in the northern hemisphere is probably due to the combination of large dipole tilt angle and positive IMF B z , which facilitates the lobe reconnection.

Formation Mechanism of Fingers That Protrude Eastward From the Io Plasma Disk During the Interchange Instability

JGR:Space physics - Wed, 07/31/2024 - 20:25
Abstract

The solar wind-magnetosphere-ionosphere interaction at Jupiter is reproduced numerically adopting the nine-component magnetohydrodynamic simulation. Calculations take into account the magnetosphere-ionosphere coupling, Jovian rotation, and Io plasma source. High-speed rotating plasma inside restricted magnetospheric space causes expansion and contraction of magnetic field, forming super-rotation at radial distance 20∼30 Rj and co-rotation breakdown further outside. Field-perpendicular current that restores co-rotational delay beyond 30 Rj is connected via field-aligned current to the main oval in the ionosphere. Inside 20 Rj, there is almost co-rotation region (deviation from co-rotation less than 20 km/s). Particularly within 10 Rj, the deviation from co-rotation is less than 2 km/s. In the nearly co-rotating region, the Io plasma forms a disk structure through field-aligned redistribution. The interchange instability occurs near the outer wall of the Io plasma disk, and instability flow develops to vortex. Through this instability, a part of the centrifugal drift current supporting the Io plasma disk is connected to low-latitude field-aligned current that generates beads-like spots on the lower latitude side of the main oval. Resulting interchange instability comes to satisfy the structure of convection and enables further development of vortex. The Coriolis force acting on eastward flow inside the developing vortex makes this flow protrude further outward, forming eastward bending fingers. Inside 10 Rj, Io plasma transport by the interchange instability becomes slower, despite the center of the disk. Io plasma escapes from the inner magnetosphere with a time constant of 20 days if this slow transport is taken into account.

Precise point positioning (PPP) based on the machine learning-based ionospheric tomography

Publication date: Available online 22 July 2024

Source: Advances in Space Research

Author(s): Pengxiang Chen, Dunyong Zheng, Wenfeng Nie, Fei Ye, Sichun Long, Changyong He, Mengguang Liao, Jian Xie

Evaluation and refinement of ERA5-land 2 m atmospheric temperature in GNSS precipitable water vapor

Publication date: Available online 20 July 2024

Source: Advances in Space Research

Author(s): Caiya Yue, Hu Wang, Liya Hu, Yamin Dang, Yafeng Wang

Theme by Danetsoft and Danang Probo Sayekti inspired by Maksimer