IRANIAN JOURNAL OF GEOPHYSICS
Year:2012 | Volume:6 | Issue:3|Papers Count:14

Horizontal and vertical gradients of the potential fields are used routinely to enhance the edge of the magnetic and gravity sources; furthermore, they are used as useful tools in interpreting and processing of magnetic and gravity data. In general, the derivatives of the potential fields are divided into horizontal and vertical derivatives, and they have always been significant tools in interpreting and processing of potential data. The derivatives can be determined in two procedures, direct measuring when the data are recorded, and calculation using mathematical and numerical methods. Many interpreting methods, that estimate the depth, location and the shape of a potential source, are based on using the gradients of potential fields. For example, both analytic signals and Euler Deconvolution methods, that have been widely applied, basically use the potential field derivatives. In these methods, different kinds of first order derivatives or derivatives of other positive integer orders are commonly used. In the basic equations of these methods, it is possible to use the derivatives of fractional orders in place of derivatives of other positive integer orders. Derivatives are high pass filters. They intrinsically amplify any noise and shallow anomalies present in the data. Therefore, using high order derivatives would be less common. Instead of using high order derivatives, one should use fractional order derivatives of the field. Besides, negative order derivatives are applicable in these kinds of methods and equations, and they can be considered as an interesting property of negative order derivation that acts as a low pass filter. In addition, horizontal fractional derivatives can be used instead of reduction to the pole at low latitudes to eliminate the instability of the reduced data. In this paper, the methods of the field gradient calculation, their alternation, and the application of the fractional order derivatives in analytic signals and reduction to the pole were inquired. To study the effects of the derivatives of different orders, the method was applied to synthetic data generated by various magnetic models such as a thin dike, and a horizontal cylinder. In the next stage, to simulate the real cases, the data was contaminated by random noise. To produce the synthetic data, the forward modeling was used. Finally, the method was applied to an aeromagnetic data set acquired over an area in Sweden. According to the geological studies in this region, there exists a granite intrusive body with certain fractures in which Diabase veins have penetrated. The results show that the fractional order derivatives as well as negative order ones are useful in data processing, and they can be considered as the principle of some of interpreting methods. All of the processing steps in this paper have been performed by using the code that we have written in Matlab.

Deconvolution is a longstanding problem in many areas of signal and image processing with applications in astronomy, remote-sensing imagery, medical imaging, and other fields working with imaging devices. It is also one of the major steps of seismic data processing and is studied in the framework of inverse problem theory. It is an ill-posed problem in the sense that the recovered solution (reflectivity series) is unstable and very sensitive to the presence of noise in the data. It is well known that the solution of an ill-posed problem is practically unusable unless taking into account some prior information about the original solution and the accuracy of such information highly affects the quality of the final regularized solution. Mathematically, usability of the prior information is ofgreat importance. The availability and usability of the prior information are two main concerns in solving inverse problems and hence deconvolution. Here, we introduce and develop some priors (in the category of heavy-tailed priors such as Cauchy and Laplace priors) that favor solutions having isolated spikes. One of the main advantages of such priors is that they have less penalization on large spikes corresponding to the true reflection coefficients while severely penalizing small spikes due to the noise and therefore results in a sparse reflectivity series. We then used the Bayes theory to incorporate the prior information into the formulation of deconvolution problems. Therefore, in this study, deconvolution was formulated in the framework of Bayes theory and the regularized solution of the problem was considered as a maximizer of the posterior probability distribution including the likelihood and the prior terms. In contrast to the conventional Wiener deconvolution which results in a minimum L2-norm solution, the methods presented in this paper recover the minimum structure or simple solutions. Sparse or simple solutions are more consistent with true earth reflectivity series, since the earth reflectivity series is simple in the sense that most of its coefficients are zero. The non-zero coefficients identify and quantify the impedance mismatches between different geological layers that are of great interest to the geophysicist.After formulating the deconvolution as a general cost function which can be convex or non-convex, we study an alternative method of determining its minimizer, as the limit of an Iteratively Re-weighted Least Squares (IRLS) algorithm. The IRLS algorithm benefits from simplicity and is easy to be coded. Furthermore, it is shown that its convergence to a local minimum from any initial guess is guaranteed and the convergence rate is superlinear. The main step of the proposed IRLS finds, for a given diagonal weight mtrix, the solution of a weighted zero order quadratic regularization where matrix is updated at each iteration. Moreover, using different priors about the reflectivity series results in only a simple change of the definition of matrix Q.Numerical experiments with synthetic and field data show that the proposed sparsity-based deconvolutions estimate the reflectivity with good resolution. Therefore, they can be used for accurate delineation of the thin layers in real poststack seismic data. The numerical results also show that the proposed methods perform much better than conventional Wiener deconvolution in the sense of the reconstruction error.

Magnetic data is routinely presented as contours or color-shaded maps of the total magnetic intensity (TMI). An interpreter’s task is to identify features (anomalies) within the map and qualitatively and/or quantitatively interpret them into geologic structures at depth. An interpretation difficulty with TMI anomalies is that they are dipolar (anomalies having positive and negative components) such that the shape and the phase of the anomaly depends in part on the magnetic inclination and the presence of any remanent magnetization. This anomaly complexity makes the interpretation more difficult because the body and its edges do not necessarily coincide with the most obvious mapped feature (e.g., anomaly maxima). The reduction-to-the-pole (RTP) technique transforms TMI anomalies to anomalies that would be measured if the field were vertical. This RTP transformation makes the shape of magnetic anomalies more closely related to the spatial location of the source structure and makes the magnetic anomaly easier to interpret, as anomaly maxima will be located centrally over the body.Since the early 1970s a variety of automatic or semiautomatic methods, based on the use of the horizontal and/or the vertical gradients (derivatives) of potential-field anomalies, have been developed as efficient tools for the determination of the geometric parameters, such as the locations of the boundaries and the depths of the causative sources. The success of these methods results from the fact that quantitative or semi-quantitative solutions are found with no or few assumptions. To map the edges of the bodies, the horizontal derivatives of the RTP field are often used. The horizontal derivative will peak above a vertical contact. However, a dipping contact, an incorrect inclination used in the RTP transformation or presence of remanent magnetization, will tend to shift the anomaly maxima away from the true location of the contact.Vertical derivatives are used in the interpretation of potential field anomalies extensively. This filter enhances the details and sharpens anomalies. However, difficulty is that by this filter the noises are increased as signals are enhanced. This filter is normally used with the first and second orders. However, vertical derivatives are recently applied with a non–integer order, in order to produce a good equilibrium between signal and noise. The analytic signal for magnetic anomalies was initially defined as a “complex field deriving from a complex potential” (Nabighian, 1972). This function can be computed easily in the frequency domain; its real part is the horizontal derivative of the field and its imaginary part is the vertical derivative. Analytic signal processing and interpretation requires few initial assumptions regarding the source body geometry and magnetization and is particularly efficient at an early stage of the interpretation even if the constraints are not available. The amplitude of the 3-D analytic signal of the total magnetic field produces maxima over magnetic contact regardless of the magnetization direction. The “theta map” is a processing technique, derived from the analytic signal that highlights the magnetic contact in a TMI image. The method is equally valid for data that has been reduced to the pole or to the equator, but it was developed to process data gathered at low magnetic latitudes, where traditional reduction to the pole is not advisable. The theta map independently detects the edges of the strike and amplitude and is thus the most valuable at low latitudes where north-south-trending anomalies disappear in the TMI data. It can also be used to qualitatively estimate dip. The magnetic tilt angle is a normalized derivative based on the ratio of the vertical and horizontal derivatives of the RTP field. The tilt angle was first described by Miller and Singh (1994), before being further refined by Verduzco (2004) and Getech. The value of the tilt angle above the edges of the contact is zero. This suggests that contours of the magnetic tilt angle can identify the location (q=0o) of the contacts.In this study, we applied the described filters on synthetic and real data that gathered from the iron deposit in Robat Posht Badam in Yazd Province in Iran.

Tehran, Iran’s capital with more than 10 million population is located in the southern foothills of the Alborz collision zone. The Alborz active mountain range consists of several sedimentary and volcanic layers, EW trending mountain belt 100-km wide and 600-km long, is bounded by Talesh Mountains to the West and by the Kopet Dagh Mountains to the East.5±2 mm/yr shortening and 4±2 mm/yr left-lateral strike-slip motion in central Alborz implies a slip partitioning between strike-slip and reverse faults across Alborz. The city is bounded by active faults. Several of these faults have been mapped but their geometry at depth, their seismicity and kinematics are not precisely known. Historical earthquakes are associated with Mosha, Taleqan, Parchin and Garmsar faults, with the largest events on the Garmsar (Ms~7.6) and Taleqan (Ms~7.7) faults during the third and tenth centuries BC, respectively. Obtained information until now reveal that better understanding of the Alborz region needs more detailed studies in longer time intervals. Several questions about faults geometry, associated seismicity, their interactions and the mechanism of deformation in this region are remained unanswered. Considering the weak geological evidence of fault activity in some parts of Tehran, and rare calculated focal mechanisms for large earthquakes, moment tensor solution of small ones can help us with better understanding of fault behavior in this region.Combining the data recorded by 29 local seismic and accelerograph stations, October 17, 2009 Ray Earthquake, Mw 4.3, was located in the westernmost part of the Parchin fault in the south of Bibi Shahrbanoo Mountain, 35.57o latitude, 51.51o longitude and 15 km depth in south-east edge of Tehran mega city. Using first motion data, a reverse mechanism with a small component of the strike-slip motion was determined.Deviatoric moment tensor was inverted by using broadband data recorded by seven Iranian stations from National Seismic Network, INSN. We used ISOLA program (Sokos and Zahradnik, 2008) that is based on the multiple point-source representation and iterative deconvoloution method, similar to Kikuchi and Kanamori (1991) for teleseismic records, but here the full wavefield is considered, and Green functions are calculated by discrete wavenumber method of Bouchon (1981). Doing many tests, we selected the 0.06-0.095 Hz frequency range that resulted in the highest variance reduction. Besides, we examined the centroid-depth range between 5 and 23 km to find the best correlation. To calculate Green functions, we used the velocity model by Abbasi et. al. (2010) for the Southern flank of Alborz. Inversion with different data subsets verified the stability of the solution.The deviatoric moment tensor inversion for this earthquake by waveform modeling shows almost a pure reverse mechanism, 97% DC component, in northwest-southeast direction along Parchin fault and a centroid depth of 11 km. It is another evidence of dominant reverse mechanism in the southern edge of the Alborz region that implies the accommodation of deformation in Alborz by the slip partitioning. The estimated seismic moment for this earthquake was 3.096e15 Newton meter resulting in a 4.3 moment magnitude using Kanamori (1977) relation.

The availability of a large amount of the data recorded by the Iranian Seismic Telemetry Network (ISTN) has motivated this study to develop relations for the routine determination of ML scale for Central Alborz region of northern Iran. The ML is commonly used in engineering because it is determined within the frequency range (0.5-3 sec) of interest in most of such applications. For any comprehensive seismic hazard analysis, one needs a calibrated magnitude relationship as well as an earthquake catalog for the study region. It is a well-known fact that the regional geology has a great influence on magnitude relations. Therefore, for each seismic region a specific magnitude relation has to be developed. The ML scale is based on the arithmetic mean of horizontal components of the synthesized Wood–Anderson seismograms. We used both nonparametric and parametric methods for inversion. We used a large dataset of 3886 events including 62031 waveforms which recorded by Tehran, Semnan and Sari seismic networks during 02/03/1997 to 13/03/2011. These seismic networks comprise of 19 three-component stations. We calculated the associated synthesized Wood-Anderson seismogram for each SS-1 waveform which records the velocity. Based on Richter’s method, we used amplitudes which are arithmetic means of those of horizontal components.Richter’s ML formula first developed for southern California and Savage and Anderson introduced a nonparametric least-squares inversion method which has been used by others. In this method, the amplitudes recorded at arbitrary distances are linearly interpolated to yield values for the attenuation curve at some fixed distances. In this study, we used both methods.The resulting equations are -logA0=0.9819log (r/100)+0.0028 (r-100)+3.0 and-logA0=1.076log (r)+0.0029 (r)+0.5580 from parametric and non-parametric methods, respectively. Where r is hypocentral in kilometer and A0 is amplitude in millimeter. The two methods yielded very similar results. Unlike the parametric method, the nonparametric one does not impose any a priori assumption of the shape of the attenuation curve on the data and has the potential to detect hinges in the attenuation curve that are caused by structural boundaries such as Moho or geological variations affects on the attenuation curve. Thus the result obtained by nonparametric method was chosen as the final result.Bakun and Joyner (1984) give the following formula for the Q/f ratio: taking an average S-wave crustal velocity of VS=3.3 km/sec, the k value obtained by the non-parametric method, 0.0029, would imply a Q/f ratio of 150 in Central Alborz, Iran.

Geomagnetism has always been at the forefront among the various branches of geophysics. In geomagnetics, we have different methods for estimating the depths and shapes of the magnetic bodies in data interpretation. One of the essential and significant methods to describe the geological complexity of earth’s crust is modeling of magnetic data by inversion.Modeling and inversion of the total magnetic field and its compliance with the basic principle to minimize the cross-sectional area of the source bodies are described. The software code, with an interactive graphical interface, operates in MATLAB environment. The code of the inversion procedure is based on a least-squares algorithm, according to a criterion of balancing the weight of the data inaccuracies and the compactness of the solution. The interpretation of magnetic data can sometimes involve two steps, calculation of the direct problem (forward modeling) and solution of the inverse problem (inversion). Forward modeling allows one to compute the theoretical response due to the magnetic source bodies, assuming some hypothesis on the shape and the volume of the magnetic body and the susceptibility contrast between the body and the hosting environment. The analytical solutions to these problems are usually given for bodies of simple shape and regular geometry. The inversion procedure considers the observed profile data or gridded data and using an optimization procedure, estimates the distribution of the susceptibility, shape and volume of the buried magnetic bodies. The principle of the compact inversion involves minimizing the area of the source body, which is the same as maximizing its compactness. Since most of the cases we deal with are slightly underdetermined problems, we need to solve the inversion problem using the least-squares method. The method consists of an iterative procedure in which the weighting matrices change in each iteration until a satisfying convergence of the solution is obtained. The input parameters for the inversion procedure are: (1) the maximum number of iterations, (2) the maximum allowed value of the susceptibility contrast, (3) the noise-over-signal ratio (The model resolution is greatly affected by the choice of the parameter of N/S). The choice of the iteration which offers the best fit is driven by the minimum norm between the calculated and experimental data. A quasi-automatic selection of the signal segments that could be considered as carrying information on the targets was proposed (mask signal). We derived the inversion operator on those parts of the signal that we called ‘‘useful signal’’, i.e. the main anomalies.In this study, we provided a method for magnetic data inversion to make 2D susceptibility models of an area with a suitable potential of Iron-ore. We made use of a 2D inversion method to study the magnetic data of Shahmirzad located in Semnan Province in Iran, to evaluate the hematite mineralization in the area. After data acquisition and processing, we applied an automatic 2D inversion to two profiles. This algorithm was based on the physical parameter distribution method. The subsurface beneath the profile was divided into a great number of infinitely long horizontal prisms with unknown susceptibilities. Solving an underdetermined system of equations in MATLAB resulted in a magnetic susceptibility distribution inside the earth which was related to the hematite content of the rocks. Inversion results on the selected profiles have shown some anomaly sources with trending east-west strike.Finally, after the modeling and inversion and applying the mask signal method to two profiles of Shahmirzad magnetic data, the models showed a steep anomaly in this region with an average thickness of 10 m, a depth of approximately 5 to 25 m and 100 m long. This iron-ore contains hematite mineral with a susceptibility of 0.05, located in the middle of the area of study between igneous intrusive masses and the sediments of limestone.

The full text of this article is the Persian language, please view the full text, please refer to the Persian language.Please click here to view the full text of this article.

Remote sensing is a maturing discipline that calls for a wide range of specialties and crosses boundaries between traditional scientific and technological disciplines. Its multidisciplinary nature requires its practitioner to have a good basic knowledge in many areas of science and requires interactions with researchers in a wide range of areas such as electromagnetic theory, spectroscopy, applied physics, geology, atmospheric sciences, agronomy, oceanography, plasma physics, electrical engineering, and optical engineering. The scattering of electromagnetic waves by precipitation particles and their propagation through precipitation media are of fundamental importance in understanding the signal returns from dual-polarized, Doppler weather radars.The main advantage of using radars for precipitation estimation is that they can provide measurements over large areas (about 10 000 km2) with fairly high temporal and spatial resolutions. Installing just one guage for each radar spatial sample (150 m resolution in range and one-degree resolution in azimuth) would require more than one-quarter of a million guages over a 150-km radius. These measurements are sent to a central location at the speed of light by “natural” networks. In addition, radars can provide fairly rapid updates of the three-dimensional structure of precipitation.The use of the radar data to detect atmospheric phenomena with suitable spatial and temporal resolutions has become one of the main methods to improve the performance of numerical weather prediction models. The effects of assimilating radar data to the ARPS numerical model on short time rain forecasts were investigated for a region covering parts of Tehran and Qom Provinces. The investigation was carried out for a synoptic system that affected the central and southern regions of Iran on March 31, 2009. The result of the juxtaposition of a Sudanese low and a strong Siberian high, the synoptic system led to remarkable rainfall in the main parts of the region of interest while leaving the southern flanks of Eastern Alborz with little rain.The ARPS numerical model was ran in two different ways: first, with the GFS (Global Forecast System) data in 3-hour time intervals; second, using the same GFS data together with the assimilation of the data of Tehran's meteorological radar. The results of the latter two applications were compared with the actual observed rainfall accumulated over 6-hour and 24-hour intervals on March 31, 2009.The results demonstrated the usefulness of assimilating radar data to improve the rainfall forecast, both quantitatively and qualitatively. The effects of the radar data are felt more strongly at the final hours of the model run. This is due to the fact that the last part of the radar data was assimilated to the model at 10: 30 UTC. The usefulness of the radar data assimilation is less felt in the high-altitude parts where the rain forecast critically depends on the particular cloud and the scheme used in convection parameterization. For the same reason, the rainfall forecast error is usually larger in the high-altitude parts of the region.

Temporal variations of radon concentration in soil and groundwater might be one of the few promising precursors for earthquake prediction. In this study, the relation between radon concentration and aftershocks of Bam Earthquake (26/12/2003, Ms=6.8) has been investigated. The radon monitoring station was located at 29oN and 58.4oE, precisely on Bam Fault where there have been high occurrences of seismic activities. The study was carried out using an active method involving an Alpha Gurad PQR2000, Alpha Pump and relative accessories which is a device capable of accurately measuring radon concentrations every 10 minutes. Air was being pumped from ground to the measuring system with a flux of 1 L/min. Forced air suction was chosen in order to avoid stratification effects, very common for radon, due to its elevated weight. Radon-monitoring sites are usually chosen in the areas where higher concentrations of radon in the surface soil layer can be expected. For this propose, the radon monitoring site was placed exactly on Bam Fault, which was placed between Bam and Baravat Cities. Radon concentration monitoring data was collected in soil at 90 cm depth exposed for a period of 90 days, every 10 minutes. Radon concentration changes are not only controlled by an earthquake, but they are also controlled by meteorological parameters at the radon monitoring site such as rainfall, soil moisture, temperature and atmospheric pressure. Therefore, in order to use radon variations as a reliable earthquake precursor, we must be able to differentiate changes that are due to earthquake from those which are not.In recent years, artificial neural networks have become very powerful, intelligent tools, used widely in signal processing, pattern recognition and other applications. The main advantages of the method are the learning capability for developing new solutions to problems that are not well defined, an ability to deal with computational complexities, a facility of carrying out quick interpolative reasoning, and finding functional relationships between sets of data.We have used a modified Adaline structure to estimate the temporal variation of radon concentration related to environmental parameters. This enables us to differentiate the changes due to phenomena in the earth such as earthquakes from those of environmental parameters. Radon concentration data obtained from our site and meteorological parameters measured in meteorological station of Bam were processed by the adaptive linear neural network, Adaline. It was indicated that the linear neural network was able to differentiate linear variations of radon concentration caused by the meteorological parameters from those arose from anomaly phenomena due to the aftershocks.

A precise earthquake location and location error estimation is a crucial element in many seismological applications such as local earthquake tomography, seismicity and seismic hazard assessment. Location error estimates may also be crucial to establish whether the hypocenter trend of an earthquake sequence really marks the seismogenic structure or simply reflects ill-conditioning of the location process.So far many methods have been introduced to locate earthquakes. Earthquake location methods have undergone many changes by Geiger’s (1912) principles. One of the first programs based on Geiger’s principles is Hypo71 (Lee and Lahr, 1972), which has already been used in many studies. The basic theory of Geiger (1912) is using Taylor series expansion of the travel time function of source to station. In order to simplify the earthquake location problem, Geiger used only the first term of Taylor series expansion that led to a straight-line equation. Therefore, they are known as linearized relationships. Using the linearized relationships results in a decrease in the accuracy of earthquake location due to losing the higher terms of Taylor series; it may lead to failure in determining the location of earthquakes using a suboptimal network, e.g. where earthquakes are located outside the seismic network. Because of the nonlinearity of the earthquake location problem, all of the algorithms and methods based on theses linearized relationships solve the earthquake location equation iteratively. Thurber (1985) showed that when the depth of earthquake is smaller than the closest distance to station, determining the focal depth is not possible in the linearized methods. Furthermore, for using the higher terms of Taylor series, it is necessary to calculate higher degree derivatives, which are very complex, and sometimes impossible, using a three-dimensional velocity model. However, the non-linear earthquake location problem can also be solved directly by a range of probabilistic algorithms (Tarantola and Valet, 1982). Tarantola and Valet (1982) presented a method that determined the location of earthquakes with fully non-linear relationships without any need to calculate the partial derivatives. The basic theory of nonlinear probabilistic method to determine the location of the earthquakes was introduced by Tarantola and Valet (1982) and Tarantola (1987).Reporting a reliable uncertainty for the location of an earthquake is one of the most important parts of earthquake location, so that presenting the epicenter and depth for an earthquake without the uncertainty is completely meaningless. Moreover, knowing the uncertainty of a location is very important in many other studies such as seismicity and tomography. Thus all the methods and algorithms designed to earthquake location; present an uncertainty for the depth and epicenter of the location. Calculation of uncertainty in an entire earthquake location problem, such as Hypo71 (Lee and Lahr, 1972) based on Geiger’s principles and NonLinLoc (Lomax et al., 2000) is by calculation of a covariance matrix. The basic premise in these methods is that the uncertainties of the observed arrival times and their relationship with the predicted travel times are assumed to be Gaussian (bell-shaped). A bell-shaped error in the time of receipt will be achieved only if the error is observed at the time and is calculated from a random and independent model. However, apart from errors that result from picking the seismic phases (in arrival times); the biggest error in an earthquake location is given by the seismic network. Bondar et al. (2004) identified four main network criteria for epicenter accuracy: (1) the number of phases used in per location; (2) the distance to the closest station; (3) the azimuthal gap; and (4) the secondary azimuthal gap.Thus, many studies are done to find optimal conditions for the use of a network station, e.g. Chatelain et al. (1980), Kissling et al. (1988), Gomberg et al. (1990). Based on the relocation of explosions, Bondar et al. (2004) introduced four characteristics for an optimal seismic network to achieve a location within a 95% confidence level and under 5 km error in depth and epicenter: (1) there are 10 or more stations, all within 250 km, (2) an azimuthal gap of less than 110o, (3) a secondary azimuthal gap of less than 160o, and (4) at least one station within 30 km.Another source of related errors is to use an inappropriate velocity model of the Another source of related errors is to use an inappropriate velocity model of the seismic waves to predict the travel times from source to stations.In this work, to investigate the calculation of uncertainty in different location methods, we compared the performances of nonlinear and linear earthquake location methods with synthetic data by simulation of three clusters of earthquakes in Central Alborz region where the location problem was ill-conditioned. Comparisons were made between the non-linear probabilistic algorithm named NonLinLoc and linear location method known as Hypo71. We studied the performance of these algorithms under different suboptimal network conditions including primary and secondary (largest azimuthal gap by removing single station) azimuthal gaps, an inappropriate velocity model, phase-reading error and the distance to the nearest station using various synthetic tests in the same network-geometry conditions of the real earthquake sequences in Mosha, Firuzkuh and Qom regions.We found out that in the suboptimal network conditions, the location error estimates from Hypo71 were, in general, less accurate than NonLinLoc's and NonLinLoc solutions were more reliable. For earthquakes occurred inside a dense seismic network, we concluded that linearized methods produced lower quality location error estimates with no overall bias in the hypocentral coordinates compared to non-linear methods.