مجله ژئوفیزیک ایران
سال پانزدهم شماره 1 (پیاپی 50، بهار 1400)

One of the most important geophysical problems in exploration of mineral deposits is to estimate the depth of buried structure using observed gravity data. In this paper, we are trying to estimate mass anomaly depth by using one of the intelligence methods as Particle Swarm Optimization (PSO) with simple shapes as sphere, horizontal and vertical cylinder. In this modeling, two parameters of depth (z) and shape factor (q) were considered as particles and the maximum and minimum depth were used as the prior information. The method is tested for synthetic models with random noise. The method gives precise results for synthetic models contaminated with random noise which is quite acceptable and promising. This technique was also successfully applied to real data for mineral exploration. The applied real data belongs to an area with hilly topography located in the Fars province close to the Abadeh city where the barite deposit is under exploration. The method is used for a profile of real data that is provided from the residual anomalies and passed from the main detected positive anomaly in the area. The estimated depth from this method was 8.15 m which was in good agreement with the results obtained through Euler method and also outcrop‐scale observations.

Salt dome is a diapir shaped structure of salt that intrudes vertically through sediment layers and surrounding strata due to its low density. Salt area identification, determining its boundaries and its 3D modeling in seismic data is a crucial issue in the literature of the seismic data interpretation. Due to its high impermeability characteristic, it can form stratigraphic oil traps by sealing the hydrocarbon reservoirs and also could be used as underground storage for natural gas and disposal sites for hazardous waste such as isolation nuclear waste and creation of the compressed air reservoir. The steeply dipping complex-shaped structures related to the salt movement and significant difference in seismic wave propagation velocity inside the salt dome with the enclosed media, imposes significant challenges for seismic data processing and interpretation. Identification and delineation of salt body is a key step in seismic data processing and interpretation, which can help geophysicist to overcome aforementioned problems. In reflection seismic methods, salt boundaries are more often characterized by change of seismic character of the signal also called texture. There are several methods available for texture analysis in image processing that can be divided into seven classes which are statistical analysis, structural methods, transform based approaches, model-based methods, graph-based techniques, learning based strategies and entropy-based methods. Textural attributes characterize the spatial arrangement of neighboring amplitudes. Extraction of seismic texture attributes can be performed using spectral information of image such as Gabor filters and local 2D Fourier spectra. Dip, similarity and coherence are the common structural attributes which are used generally for textural analysis in seismic data. The most common rational approach to describe the texture in seismic image is to measure the statistical properties of the image. Gray Level Co-occurrence Matrix, chaos and variance are three conventional statistical seismic texture attributes used for this purpose. Due to the textural contrast of the salt dome with the surrounding layers and sediments, edge detection tools can also be used to determine the boundaries of textural changes and delineate the salt area. In this study, we used a new textural seismic attribute known as the gradient of texture to characterize the change of seismic character between the salt body and its surrounding geology. It calculates the texture gradient in two adjacent windows around a sample in different directions. It is supposed that different area in seismic image with different textural pattern will exhibit diverse gradient of texture. Thus, it will be appropriate for image segmentation for specific interpretation investigation. The gradient of texture attribute will differentiate desired area from the rest of the image through supervised classification and growth strategy in extending the selected classes. Efficiency of the introduced method for salt dome delineation and modeling in seismic data was investigated here by applying on a synthetic model and 3D seismic data from the Persian Gulf. Comparison between obtained results of the proposed method and conventional attributes revealed superiority of the 3D texture gradient in textural segmentation and salt dome modeling from seismic data.

Climate change is a major challenge for human society and the natural environment. Evidence suggests that human activities play a role in increasing temperature at various temporal-spatial scales. The effects of climate change can be assessed by analyzing air temperature trends. According to the latest IPCC report, global average temperatures will increase by 1.5 degrees Celsius by the end of this century. The main purpose of this study is to assess CMIP6 projected temperature over different climate zones of Iran and its trend in the future. The result of this study can be useful for a wide range of management areas, especially the study of water resources, snow reserves, agriculture, and tourism.    In this study, the mean annual temperature data of 43 synoptic stations of the Iran Meteorological Organization (IRIMO) were obtained from 1990 to 2009. To evaluate the anomaly and temperature trend in Iran until the end of the 21st century, the data of three models BCC-CSM2-MR, CAMS-CSM1-0, and MRI-ESM2-0 of CMIP6 models under two scenarios of SSP2.4-5 and SSP5.8-5 were used. We divided the period into four twenty-year periods, which are the first period (2020-2040), the second period (2041-2061), the third period (2061-2080), and the fourth period (2081-2100), respectively. To evaluate the air temperature which is the output of selected CMIP6 models, three statistical measures of RMSE, NSE and KGE were used. Using the Delta Change Factor (DCF) method, the bias of the models was corrected. Then, non-parametric Man-Kendall (MK) and Sen’s Slope Estimator tests were used to analyze trend analysis and trend slope in long-term data series.  The maximum temperature is seen on the northern coast of the Persian Gulf and the Sea of Oman and the minimum temperature is seen in the northwest following the heights of the Zagros and also in the north of Iran (Alborz Mountains). The minimum annual temperature of 10.80°C was calculated based on observed data for the period 1990-2009, and the maximum temperature was 27.90°C on the coasts of the Oman Sea in southeastern Iran and Khuzestan province on the shores of the Persian Gulf in southwestern Iran. The intensity of the increase in temperature in Iran in mountainous areas is mainly due to the increase in the minimum temperature rather than to the maximum one.Generally, the intensity of the warming in Iran is mostly projected in cold and temperate regions. There is also a tendency for the temperature to rise further at higher latitudes.  The projected temperature of CMIP6 models based on the SSP2.4-5 and SSP8.5-8 scenarios in Iran over four 20-year periods from 2020 to 2100 showed that the average slope of the temperature trend in Iran will reach 0.05°C per year, which shows an increase by a factor of 0.01 throughout Iran. As a general result, the annual trend of air temperature in Iran, based on observational data and projected output of CMIP6 models, shows warmer climate conditions for Iran. Temperature anomaly was not negative in any of the scenarios and periods; whereas positive anomaly is seen throughout the country. This increase in anomaly can be a major threat to the country's water resources.

Spherical harmonic transform (SHT) provides a standard tool for spectral analysis of data measured on a spherical manifold particularly in geodesy and geophysics, such as Earth’s gravity/magnetic field, topography, atmosphere, etc. Thanks to recent satellite missions such as SRTM, the spatial resolution of geoscience data decreases, leading to an increase in the degree/order of spherical harmonic. The accurate and fast computation of fully normalized associated Legendre functions (fnALF) is the main part of spherical harmonic computations.In polar region, for high degrees (n>2000), the arithmetic underflow may occur and Extended Range Arithmetic (ERA) can be applied to fix underflow problem. In Belikov method the four terms row-wise reccurrence formula is used for fnALF computation without any numerical issues. Previous studies have only addressed the efficiency of the ERA and Belikov method without mentioning to their performnace in spherical harmonic synthesis (SHS) and analysis (SHA). This study aims at investigating the effect of row-wise and colunm-wise fnALF on SHS/SHA for high harmonic degree/order. For this purpose, we developed parallel FORTRAN95 software that utilizes Belikov and ERA to perform fnALF in the row-wise and column-wise SHT, respectively. For SHA, exact Gauss-Legendre quadrature was used that allows one to determine the level error of SHT. In the first part of numerical results, the computation time and error of Belikov and ERA methods for calculating fnALF are compared in serial and parallel processing. Based on the results of the developed software, Belikov method is approximately 2 times faster than ERA for degrees 2160 to 64800 in both serial and parallel mode. Also for mid-latitudes, the errors for both methods are nearly the same. In polar region the error in ERA grows, while Belikov method has not latitudinal dependencies. To evaluate the performance of fnALF method on error in SHS/SHA, calculations were performed using random spherical harmonic coefficients (Cnm , Snm) between [-1, +1]. Using SHS and SHA in a closed cycle, the computed coefficients (C׳nm , S׳nm) were compared with original coefficients (Cnm , Snm). The RMS and maximum error for different degrees from N=2160 to 21600 were computed. The level error of two methods in terms of RMS is in the same level. Numerical results show that in serial mode, the row-wise SHS/SHA is approximately 1.8 times faster than column-wise method. However, in parallel computing with 10 threads, the speed up factor decreases to about 1.2.

In water resources management, correct identification of bedrock characteristics, such as type, depth and structure is useful and various geophysical methods are used for this purpose. In geophysical and subsurface investigations, solving inverse problems using physics of the problem in partial differential equation (PDE) system is very important. Inverse modeling is one of the useful solutions which creates a logical model between observed and measured values. The limitations and ambiguity of individual techniques can be significantly reduced by adopting joint inversion schemes. Using inversion approaches which combine two or more geophysical methods, allows substantial improvements in modeling processes. The joint inversion of different electromagnetic and resistivity data sets in hydrogeology, structural geology and mining has been used in some studies but the gravity method is less used in hydrogeological investigations. One example of using gravity method in previous studies is an algorithm for joint inversion of gravity and resistivity data, where the interfaces corresponding to changes in the bulk density were interpreted as interfaces of porosity and water content changes, i.e. interfaces of electrical resistivity changes.    In this paper, we present a framework for joint inversion of gravity and geoelectric data solving under-determined inverse problems. It can be used in a wide range of physical systems governed by PDEs. The integrated approach is based on the connection between density and resistivity of a sedimentary sequence through the porosity. Also, we present a general adjoint state formulation which may be used in this framework. It increases the calculation speed of sensitivity matrices in a variety of commonly encountered under-determined problems. There are two steps in this research. First, 2D joint inversion of gravity and geoelectric data is run and validated in COMSOL multi-physics software using one synthetic model and synthetic data in a forward modeling process. Afterwards, using real gravity and geoelectric data surveyed along a cross-section in part of the Qotrum plain in the SE Yazd, the lateral structure of bedrock is estimated. Because Qotrum basin is already well surveyed, it offers opportunities for evaluation of accuracy and reliability in joint inversion approaches. Fifty-one gravity stations and eleven vertical electrical soundings were carried out across the central part of this basin to apply the proposed joint inversion. The results indicated that gravity-geoelectric joint inversion is more accurate than the individual interpretation of each of these methods and significantly improves the solution by decreasing the ambiguity of models. Furthermore, this method, while high in computational speed, can be used in modeling of a wide range of geophysical, hydrogeological and physical systems governed by the partial differential equation laws.

Precipitation is a favorable phenomenon for human due to the water resources for agriculture, drinking water and health purposes. Because of the importance of precipitation, understanding the mechanism of this phenomenon in synoptic, dynamics and climatology has always been of interest for researchers. Although heavy rainfall is useful for water supply, it has destructive effects, such as floods. Heavy rains show the mechanism of precipitation more clearly because the pressure generating systems are more clearly detectable.    In this study, mean sea level pressure pattern of severe precipitation in southwestern Iran over a period of 30 years was identified and classified. Moreover, the area and time of research were defined, and using the Climate Prediction Center (CPC) precipitation data, the 30-years extreme precipitation with 99th percentiles was determined. There were 64 cases of extreme precipitation from 1989 to 2018 with 99th percentiles.  Classification of extreme precipitation systems based on mean sea level pressure field was performed with a spatial resolution of 2.5° data from the NCEP-NCAR center. Extreme precipitation in southwestern Iran at warm time of year was classified using factor analysis technique. In this study, it was shown that the pressure systems that cause extreme precipitation in southwestern Iran during the warm season are classified into five main patterns: (A) The first dominant pattern is created by the Red Sea trough and the low-pressure center of Saudi Arabia. It should be noted that in this pattern, the role of the Red Sea trough is more prominent than the low-pressure Saudi Arabia. In some cases, Saudi Arabia has formed along the Red Sea trough, reflecting the influence of the Red Sea on creation and development of the Saudi Arabian low-pressure center. (B) The second dominant pattern is formed by the-low pressure of Saudi Arabia and the low-pressure center in eastern Iran. The low-pressure center in the eastern part of Iran has an important role in the precipitation of Iran. It is present in all patterns, but first one. The presence of this low-pressure center has created atmospheric fronts in Iran. The rainfall region of this pattern shows that the cold front activity of this low-pressure center causes precipitation. (C) The third pattern is similar to the second one, but a high-pressure center located in the north and northwest of Iran. The presence of this high-pressure center indicates the cold air advection behind the low-pressure cold front of eastern Iran. (D) The low-pressure center of the eastern part of Iran and the high-pressure center in the north and northwest of Iran exist in the fourth pattern. The difference between this pattern and the third one is elimination of Saudi low-pressure center in this pattern. Therefore, the dynamic effects are similar to the third pattern and the extreme precipitation caused by this pattern is due to the cold front of the low-pressure system of the eastern part of Iran. (E) The fifth pattern, like the third one, consists of low-pressure systems in eastern Iran and Saudi Arabia. The difference between these patterns is presence of low-pressure of Cyprus in the fifth pattern. The low-pressure of Cyprus is cause of classical precipitation pattern (formation of the occluded front) in this pattern. The occluded front in this pattern has expanded from low-pressure eastern of Iran to low-pressure in Cyprus and converted into a secondary cold front over several stages of the frontogenesis process.  In this research, the southwestern precipitation systems of Iran during the warm season are classified into five patterns. The most important factor of precipitation in the southwestern part of Iran is the Red Sea trough with 48% frequency of occurrence. Therefore, it can be concluded that almost half of the extreme precipitation of southwestern Iran during the warm season is related to this area that extends from the Red Sea trough to the west of Iran.

The goal of this study is to invert the seismic data directly to porosity as well as to quantify the associated uncertainty in one of the carbonate reservoirs located in southwestern Iran. Probabilistic inverse methods are able to present the model parameters as a posterior probability distribution function by combining the probability density function of the model prior information and the likelihood model. The likelihood density function is defined based on the noise model. The answer to the probabilistic inverse problem is a posterior distribution function that is not only consistent with the prior model but also is constrained to the seismic data. In this study, one of the sampling methods based on Markov chain Monte Carlo is used, which is able to generate realizations of the desired model parameter by sampling from the posterior distribution function. Unlike deterministic inverse methods, which provide only one answer for the model parameters, in probabilistic inverse methods the realizations generated from the posterior distribution function allow statistical analysis and model uncertainty quantification. The results of the implementation of the proposed method on synthetic seismic data showed that the estimation of noise variance has a significant effect on the results of probabilistic inversion and the uncertainty of the model realizations. The underestimation of the noise variance leads to fitting the noise on the data and subsequently generates artifacts on the output realizations. The overestimation of the noise variance provides smooth realizations with higher uncertainty. In the latter case, the reference porosity model is in the 95% confidence interval in contrast to the former case. Therefore, care should be taken in estimating the noise variance in probabilistic inverse methods. Considering a calibrated rock physics model for the carbonate reservoir under study, which is the main core in a direct inversion approach, the proposed algorithm was applied to the seismic traces adjacent to the four well logs. The uncertainty of the porosity was quantified in each well location. The correlation coefficient of the mean of the porosity realizations and the true porosity in four well locations were approximated about 79%, 63%, 51% and 67%. The consistency of the results obtained from the inversion with the observed porosity at the well locations indicates the good performance of the algorithm in estimating the porosity and its associated uncertainty. Due to the ability of the probabilistic inverse methods in a direct inverse of seismic data to the petrophysical properties, and their applicability in being performed in a parallel structure in processing clusters, these algorithms can be used in reservoir characterization of 2D and 3D data.

Petrophysical reservoir properties are usually estimated from elastic properties such as acoustic impedance (AI) using petro-elastic models. AI data are only available at sparse well locations, especially in early exploration stages of oil fields and therefore to quantify the spatial reservoir properties, it is necessary to estimate the AI for the rest of the reservoir area. Alongside the deterministic seismic inversion methods, multi-realization geostatistical simulation approaches can be used for the parameter estimations and uncertainty quantification. Geostatistical simulation methods allow the production of two- or three dimensional models of reservoir properties using existing well log data. Conventionally, sequential Gaussian simulation (SGS) method is used because of its simplicity. However, its accuracy is not always guaranteed. Nowadays, turning bands simulation method (TBSim) has received much attention because of its capability to reproduce the statistical properties of the original data. The main principle behind the TBSim is simplifying a multi-dimensional geostatistical simulation problem into a set of fast 1D simulation problems. The simulations are performed along the uniformly distributed lines spanning a unit sphere using sinusoid functions.  In this study, we are going to generate 3D AI models in an oil field in SW of Iran via the TBSim. We will also compare the results with the AI data computed from seismic inversion. Geostatistical modeling has less computational complexity than seismic inversion. Although AI modeling cannot be used as a definite alternative to seismic inversion, it can be used as a primary method for estimating AI, especially in areas without seismic data. It also has the capability to be involved with some seismic inversion methods known as stochastic seismic inversion methods. The dataset used in this study includes AI logs of seven wells, one of which (known as test well) has been excluded from calculations to check the accuracy of the results. The study reservoir zone includes the Ahwaz sandstone member in the upper part of Asmari formation and also the lower carbonates. The lower carbonates of Asmari formation are separated by an unconformity with Jahrum formation. After constructing the structural model and upscaling the AI logs in it, 50 realizations of AI are generated in the gridded model. The results in the test well indicate a high correlation between the modeled values and real AI data as well as the model obtained by seismic inversion. The maximum correlation of AI values of different realizations with real values in the test well equals to 78.8%. The correlation coefficient achieves to 82.9 % for mean of realizations and is 72.8% for seismic inversion data. The cube of mean of realizations is also in good agreement with the seismic inversion cube and reproduces the dominant trend of vertical and lateral AI variations. Comparison of results histogram with real data as well as the seismic inversion data reveals the capability of geostatistical AI modeling in reproducing the statistical properties of original data. In addition, the results of the uncertainty analysis of produced models also confirm the reliability of these models, especially in the test well. Therefore, we would recommend the TBSim as a powerful method for AI modeling during reservoir characterization.