首页 | 本学科首页   官方微博 | 高级检索  
相似文献
 共查询到20条相似文献,搜索用时 171 毫秒
1.
The detectability of an intermediate layer in a three-layer earth model in the time domain has been investigated. The calculations were made for the perpendicular loop (designated system II) and vertical-coplanar (designated system III) electromagnetic (EM) sounding systems. The primary excitation employed is a train of half-sinusoidal and square waveforms of alternating polarity. The time-domain response has been determined by Fourier transformation of the matched complex mutual coupling ratios into the time domain and by linear digital filtering. Top and bottom layers have equal resistivity. EM responses have been computed for conductive and resistive intermediate layer with a wide range of thickness and for two values (500 m and 1000 m) of loop-separation. For the detectability analyses, the root mean square (rms) difference between three-layer and homogeneous-earth responses is adapted. The threshold value for detectability is defined as an rms difference of 10% and the measurement error is arbitrarily assumed to be of the order of 3%. It is observed that the perpendicular-loop system is better than the vertical-coplanar system in detecting thin intermediate layers (either conductive or resistive). For a loop separation of 1000 m and half-sinusoidal pulse excitation, the detectable thickness ratio (h2/h1) is 0.10 by system II for the conducting middle layers; for square pulse excitation the corresponding thickness ratios are 0.06 for system II and 0.12 for system III. For a loop separation of 1000 m and half-sinusoidal pulse excitation in detecting the resistive intermediate layers, the corresponding thickness ratios are 0.9 for system II and 2.25 for system III; while for square pulse excitation the thickness ratios are 0.55 for system II and 1.55 for system III. Results in the frequency domain and time domain (for half-sinusoidal and square pulsed field) have also been presented for systems II and III for detecting conducting layers by considering an earth model where p1≠ p3 and p3 > p1 (p is the resistivity). The loop separa- tion used is 1000 m. Direct comparisons between the frequency domain and time-domain results clearly demonstrate the superiority of frequency-domain systems for detecting con- ducting intermediate layers.  相似文献   

2.
The digital linear filter method is used to compute the normalized vertical magnetic field for a circular loop in CFS system. Three-layer earth models with resistive and conductive basement are considered. The corresponding field expressions are suitably written, and the multifrequency response is computed and presented in convenient forms. Analysis of theoretical data indicates that for highly resistive basement, the variation in layer conductivity and intermediate layer thickness is well reflected on three-layer amplitude response curves at low frequencies and at high conductivity contrasts between first and second layers. This, however, is not true in the case of conductive basement, where the resolution of the intermediate layer is observed to be comparatively poor. The resolution of an intermediate conductive layer in a three-layer sequence is found to be satisfactory.  相似文献   

3.
Using approximate boundary conditions, expressions for electromagnetic fields have been derived for a thin, highly resistive layer lying between two homogeneous layers excited by an electric dipole grounded on the surface of the earth. The variations of the fields with the parameter T/T1 (ratio of the transverse resistance of the thin layer to the transverse resistance of the first layer) were studied in relation to frequency, time, the normalized separation source—receiver, and the angle between the source and the radius to the observation point. For a value of h2/h1 (ratio of thickness of second layer to the thickness of the first layer) approximately equal to 0.2, the general three-layer medium case gives the same results as this approach. It was found that the electric fields have a very strong dependence on the parameter T (transverse resistance) which characterizes the thin, highly resistive layer. However, the magnetic fields depend only very weakly on this parameter.  相似文献   

4.
Global optimization with very fast simulated annealing (VFSA) in association with joint inversion is performed for 1D earth structures. The inherent problems of equivalence and suppression in electromagnetic (EM) and direct current (DC) resistivity methods are studied. Synthetic phase data from multifrequency sounding using a horizontal coplanar coil system and synthetic apparent resistivity data from Schlumberger DC resistivity measurements are inverted individually and jointly over different types of layered earth structures. Noisy data are also inverted. The study reveals that global optimization of individual data sets cannot solve inherent equivalence or suppression problems. Joint inversion of EM and DC measurements can overcome the problem of equivalence very well. However, a suppression problem cannot be solved even after combination of data sets. This study reveals that the K-type earth structure is easiest to resolve while the A-type is the most difficult. We also conclude that the equivalence associated with a thin resistive layer can be resolved better than that for a thin conducting layer.  相似文献   

5.
The objective of this study is to generate the separation-distance-domain (r-domain) transformation of the theoretically calculated wave number domain (m-domain) electromagnetic induction field component Bz(m, ω) of a stratified medium and to search for interpretive information which has been absent in the previously achieved numerical solutions of the problem. The r-domain kernel R?(r, ω) function defining the induction field appears to adequately reflect the layering and electrical properties of the medium if it is expressed as a function of the frequency if the source-receiver separation r is small with respect to the thickness of the first layer. However, exact values of the conductivity cannot be distinguished from those of the neighboring values unless a resistive basement layer is present. This feature is the result of the truncation in series representation of the kernel function R?(m, ω). However, this truncation is regarded as significant in the case of a conductive first layer. In m-domain static-zone studies, a conductive first layer slightly influences its r-domain correspondent. Although the computational cost of obtaining the kernel B(r, ω) by evaluation of the convolution in a cylindrical coordinate system is high, this semi-analytic solution is still superior to those based on the asymptotic assumptions.  相似文献   

6.
The transient response of a conductive shell-shell model in the one-loop version was obtained analytically. The results indicate that four zones, namely early, late early, intermediate, and late zone can be identified in the total transient characteristic of the model. In case the measurements are carried out in the late early zone, a conductive target appears as a resistive one. It is suggested that the optimum time of measurement should be so selected as to fall in the intermediate zone.  相似文献   

7.
The in-loop pulse electromagnetic response of a stratified earth has been expressed in terms of an apparent resistivity- time plot using the PEM response over a homogeneous half-space which is typically unipolar with monotonic decay. This half-space response characteristic provides a unique relationship between Crone PEM channel amplitude and the apparent half-space resistivity. The possibility to resolve a thin intermediate conductive and resistive layer with the in-loop PEM system has been investigated. The system is well in shallow geoelectric mapping.  相似文献   

8.
A complex transfer function c (or generalized skin depth) can be derived from data for the secondary magnetic field measured by a dipole system with small coil spacing at height h above the ground. This function has a useful property: For a uniform or layered ground, the real part of c yields the‘ centroid depth’z* of the in-phase current system as a function of frequency. This parameter can be combined with the apparent resistivity ρa derived by conventional methods. The function ρa(z*), if known over a broad frequency range, yields a smoothed approximation of the true distribution ρ(z) without an initial model. The relations between ρa(z*) and ρ(z) are studied for a number of multilayer models. An example of the application of the ρa*) algorithm to data from a groundwater survey is given.  相似文献   

9.
Grounded-source TEM modelling of some deep-seated 3D resistivity structures   总被引:2,自引:0,他引:2  
Long-offset transient electromagnetics (LOTEM) is now regarded as a suitable electrical method for deep exploration along with magnetotellurics (MT). In this method, the vertical magnetic-field impulse response and, occasionally, the horizontal electric-field step response of a grounded-wire source on the surface of the earth are measured. Here, these two responses are computed for 3D models of three deep resistivity structures of interest in hydrocarbon exploration: (i) a faulted graben in a resistive basement rock at a depth of 4 km beneath a conductive overburden; (ii) a facies change in a resistive layer buried at a depth of 2 km in the conductive overburden above a resistive basement; and (iii) an anticlinal uplift of a resistive layer at a depth of 1 km in the conductive overburden above a resistive basement. The results show that the sensitivity of the electric-field response to model perturbation is generally greater than that of the magnetic-voltage response. Further, the electric-field sensitivity is confined to early and intermediate times while that of the magnetic-voltage response is confined to intermediate and late times. Hence it is recommended that both electric and magnetic recordings are made in a LOTEM survey so that the final results can be presented as apparent-resistivity curves derived from the two responses jointly as well as separately.  相似文献   

10.
For the computation of the vertical component Hz of the magnetic field of a horizontal A.C. dipole lying on the earth's surface, a recurrence formula is presented for a horizontally stratified half space, to obtain the (n+ 1)-layer case from the w-layer case. By means of several computed diagrams for the two-layer case, Hz can be determined for different ratios of conductivity of the subsoil and that of the overburden. Thereby the distance from the dipole as well as the layer thickness h are expressed in terms of the wave length A of a plain wave in the overburden. Assuming a sufficiently large conductivity difference, the results show that evidence about the subsurface conditions can be obtained if the distance between the measuring coil and the dipole is of the order of A/3, and if the thickness h of the layer varies within the range A/100 < h < A/6. As an example for the 3-layer case, a nonconducting intermediate layer is assumed.  相似文献   

11.
The rotating current EM method has been applied to the delineation of two conductive orebodies, Elura near Cobar, NSW, and Thalanga near Charter's Towers, Queensland. The field data were collected in the form of observations of the vertical magnetic field strength ratio and phase difference using a Turam-style receiver with twin vertical coils. By reconstituting this data back to the ring source field and phase, i.e. the observed Hz, phasor, it is possible to present contoured maps of the EM field. Anomaly phasors are obtained by subtracting theoretical phasors from the observed phasors in the complex plane of the Hz phasor. The theoretical phasors for the finite source are based on horizontally layered, half-space earth models, computed at each point of the survey grids, then normalized to selected points of the observed fields. Use is made of the intrinsic circular symmetry of the method in X–Y plots of field versus source-receiver distance to ascertain geoelectric parameters for the earth models. A steel picket fence at Thalanga is modelled by a line source grounded at each end and its Hz, phasor is removed by the same process. A considerable improvement in anomaly delineation is gained over previous Turam-style anomalies and the two survey examples illustrate the limitations of the method in the presence of a conductive overburden (Elura) and its abilities in the absence of a conductive overburden (Thalanga).  相似文献   

12.
For a thin highly-conducting layer with given longitudinal conductance the recurrence formulae for an n-fold horizontally stratified subsoil are established for d.c. resistivity and magnetotelluric soundings. Similarly, a thin low conductivity layer with given transverse resistance is treated in the d.c. case and a non-conducting intermediate bed in magnetotellurics. Model curves for a thin high- or low-conductivity intermediate layer in the three-layer case have been carried out, which may serve as an extension of the well-known three-layer diagrams for a Schlumberger configuration. The corresponding model curves in magnetotellurics are given. By numerical comparison of these curves with real three-layer curves some diagrams have been developed to show the allowed thicknesses of the intermediate layer in the Schlumberger case and in the case of magnetotelluric sounding.  相似文献   

13.
14.
An infinitely resistive/conductive horizontal bed is assumed in an otherwise homogeneous and isotropic half space. Schlumberger, three electrode, and unipole profiles are computed at right angles to the strike of the bed. The Schwarz-Christoffel method of conformal transformation and numerical methods of solving non-linear differential equations are used to solve the boundary value problem. It is observed that (i) the three electrode system is the most sensitive gradient electrode configurations for electrical profiling, (ii) the apparent resistivities for Schlumberger, three electrode, and unipole methods become maximum when the depth of the bed is 0.06 L, 0.1 L, and 0.055 L for a resistive bed and minimum when depths are 0.085 L, 0.04 L-0.02 L and indeterminate for conductive beds, respectively, (iii) the limiting depths of detection (defined in the text) by Schlumberger, three electrode, and unipole configurations are respectively 0.9 L, 6.6 L and 2.0 L for resistive beds and 0.58 L, 1.17 L and 1.5 L for conductive beds. The electrode separation L is the distance between the two farthest active electrodes.  相似文献   

15.
Abstract

The generation of eddies by a large-scale flow over mesoscale topography is studied with the help of two- and three-layer nonlinear quasi-geostrophic models of the open ocean. The equations are integrated forward in time with no eddies present initially. For a given time, the displacement of the interface between layers two and three (ζ) tends to a well-defined limit (function of the horizontal spatial coordinates) as ρ 3 - ρ 20 (ρr is the density of layer r). Even for values of α[= (ρ 3 - ρ 2)/(ρ 2 - ρ 1)] as small as 0.01 the potential energy due to ζ is not negligible and it can reach, in some cases, a considerable fraction of the total eddy energy.  相似文献   

16.
The problem of equivalence in direct current (DC) resistivity and electromagnetic methods for a thin resistive and conducting layer is well‐known. Attempts have been made in the past to resolve this problem through joint inversion. However, equivalence still remains an unresolved problem. In the present study, an effort is made to reduce non‐uniqueness due to equivalence using global optimization and joint inversion by successive refinement of the model space. A number of solutions derived for DC resistivity data using very fast simulated annealing global inversion that fits the observations equally well, follow the equivalence principle and show a definite trend. For a thin conductive layer, the quotient between resistivity and thickness is constant, while for a resistive one, the product between these magnitudes is constant. Three approaches to obtain very fast simulated annealing solutions are tested. In the first one, layer resistivities and thicknesses are optimized in a linear domain. In the second, layer resistivities are optimized in the logarithmic domain and thicknesses in the linear domain. Lastly, both layer resistivities and thicknesses are optimized in the logarithmic domain. Only model data from the mean models, corresponding to very fast simulated annealing solutions obtained for approach three, always fit the observations. The mean model defined by multiple very fast simulated annealing solutions shows extremely large uncertainty (almost 100%) in the final solution after inversion of individual DC resistivity or electromagnetic (EM) data sets. Uncertainty associated with the intermediate resistive and conducting layers after global optimization and joint inversion is still large. In order to reduce the large uncertainty associated with the intermediate layer, global optimization is performed over several iterations by reducing and redefining the search limits of model parameters according to the uncertainty in the solution. The new minimum and maximum limits are obtained from the uncertainty in the previous iteration. Though the misfit error reduces in the solution after successive refinement of the model space in individual inversion, it is observed that the mean model drifts away from the actual model. However, successive refinement of the model space using global optimization and joint inversion reduces uncertainty to a very low level in 4–5 iterations. This approach works very well in resolving the problem of equivalence for resistive as well as for conducting layers. The efficacy of the approach has been demonstrated using DC resistivity and EM data, however, it can be applied to any geophysical data to solve the inherent ambiguities in the interpretations.  相似文献   

17.
The time-domain controlled source electromagnetic method is a geophysical prospecting tool applied to image the subsurface resistivity distribution on land and in the marine environment. In its most general set-up, a square-wave current is fed into a grounded horizontal electric dipole, and several electric and magnetic field receivers at defined offsets to the imposed current measure the electromagnetic response of the Earth. In the marine environment, the application often uses only inline electric field receivers that, for a 50% duty-cycle current waveform, include both step-on and step-off signals. Here, forward and inverse 1D modelling is used to demonstrate limited sensitivity towards shallow resistive layers in the step-off electric field when transmitter and receivers are surrounded by conductive seawater. This observation is explained by a masking effect of the direct current signal that flows through the seawater and primarily affects step-off data. During a step-off measurement, this direct current is orders of magnitude larger than the inductive response at early and intermediate times, limiting the step-off sensitivity towards shallow resistive layers in the seafloor. Step-on data measure the resistive layer at times preceding the arrival of the direct current signal leading to higher sensitivity compared to step-off data. Such dichotomous behaviour between step-on and step-off data is less obvious in onshore experiments due to the lack of a strong overlying conductive zone and corresponding masking effect from direct current flow. Supported by synthetic 1D inversion studies, we conclude that time-domain controlled source electromagnetic measurements on land should apply both step-on and step-off data in a combined inversion approach to maximize signal-to-noise ratios and utilize the sensitivity characteristics of each signal. In an isotropic marine environment, step-off electric fields have inferior sensitivity towards shallow resistive layers compared to step-on data, resulting in an increase of non-uniqueness when interpreting step-off data in a single or combined inversion.  相似文献   

18.
航空电磁系统校准是开展实际测量工作的基础,校准情况直接影响数据处理和解释.传统校准方法通常假设在自由空间中进行,忽略导电大地耦合影响.然而,实际工作中很难找到绝对高阻的校准场地,导电大地对系统校准和观测数据的影响无法忽视.本文以频率域航空电磁系统为例,对导电大地上航电系统校准技术和校准误差改正方法进行研究.我们首先推导了层状导电大地上水平共面和直立共轴线圈系统的校准公式,结果表明导电大地对航电系统校准尤其是水平共面装置的高频信号影响很大.针对校准过程中大地电导率已知的情况,本文采用非线性方程求解技术一次性确定校准线圈位置和Q值;在没有任何辅助信息情况下,也可直接利用实测数据计算校正因子进行迭代求解.测试结果表明该方法快速、准确、有效.考虑到系统相位和增益调整直接影响观测数据,本文提出了航空电磁数据校准误差的改正算法.实测数据误差改正结果表明,导电大地对高频信号影响严重,校准误差改正后的航空电磁数据与实际地质资料更好吻合.  相似文献   

19.
I. Introduction In this section the problem is stated, its physical and mathematical difficulties are indicated, and the way the authors try to overcome them are briefly outlined. Made up of a few measurements of limited accuracy, an electrical sounding does not define a unique solution for the variation of the earth resistivities, even in the case of an isotropic horizontal layering. Interpretation (i.e. the determination of the true resistivities and thicknesses of the ground-layers) requires, therefore, additional information drawn from various more or less reliable geological or other geophysical sources. The introduction of such information into an automatic processing is rather difficult; hence the authors developped a two-stage procedure:
  • a) the field measurements are automatically processed, without loss of information, into more easily usable data;
  • b) some additional information is then introduced, permitting the determination of several geologically conceivable solutions.
The final interpretation remains with the geophysicist who has to adjust the results of the processing to all the specific conditions of his actual problem. II. Principles of the procedure In this section the fundamental idea of the procedure is given as well as an outline of its successive stages. Since the early thirties, geophysicists have been working on direct methods of interpreting E.S. related to a tabular ground (sequence of parallel, homogeneous, isotropic layers of thicknesses hi and resistivities ρi). They generally started by calculating the Stefanesco (or a similar) kernel function, from the integral equation of the apparent resistivity: where r is the distance between the current source and the observation point, S0 the Stefanesco function, ρ(z) the resistivity as a function of the depth z, J1 the Bessel function of order 1 and λ the integration variable. Thicknesses and resistivities had then to be deduced from S0 step by step. Unfortunately, it is difficult to perform automatically this type of procedure due to the rapid accumulation of the errors which originate in the experimental data that may lead to physically impossible results (e.g. negative thicknesses or resistivities) (II. 1). The authors start from a different integral representation of the apparent resistivity: where K1 is the modified Bessel function of order I. Using dimensionless variables t = r/2h0 and y(t)=ζ (r)/ρ1 and subdividing the earth into layers of equal thicknesses h0 (highest common factor of the thicknesses hi), ø becomes an even periodic function (period 2π) and the integral takes the form: The advantage of this representation is due to the fact that its kernel ø (function of the resistivities of the layers), if positive or null, always yields a sequence of positive resistivities for all values of θ and thus a solution which is surely convenient physically, if not geologically (II.3). Besides, it can be proved that ø(θ) is the Fourier transform of the sequence of the electric images of the current source in the successive interfaces (II.4). Thus, the main steps of the procedure are: a) determination of a non-negative periodic, even function ø(θ) which satisfies in the best way the integral equation of apparent resistivity for the points where measurements were made; b) a Fourier transform gives the electric images from which, c) the resistivities are obtained. This sequence of resistivities is called the “comprehensive solution”; it includes all the information contained in the original E.S. diagram, even if its too great detail has no practical significance. Simplification of the comprehensive solution leads to geologically conceivable distributions (h, ρ) called “particular solutions”. The smoothing is carried out through the Dar-Zarrouk curve (Maillet 1947) which shows the variations of parameters (transverse resistance Ri= hii–as function of the longitudinal conductance Ci=hii) well suited to reflect the laws of electrical prospecting (principles of equivalence and suppression). Comprehensive and particular solutions help the geophysicist in making the final interpretation (II.5). III. Computing methods In this section the mathematical operations involved in processing the data are outlined. The function ø(θ) is given by an integral equation; but taking into account the small number and the limited accuracy of the measurements, the determination of ø(θ) is performed by minimising the mean square of the weighted relative differences between the measured and the calculated apparent resistivities: minimum with inequalities as constraints: where tl are the values of t for the sequence of measured resistivities and pl are the weights chosen according to their estimated accuracy. When the integral in the above expression is conveniently replaced by a finite sum, the problem of minimization becomes one known as quadratic programming. Moreover, the geophysicist may, if it is considered to be necessary, impose that the automatic solution keep close to a given distribution (h, ρ) (resulting for instance from a preliminary interpretation). If φ(θ) is the ø-function corresponding to the fixed distribution, the quantity to minimize takes the form: where: The images are then calculated by Fourier transformation (III.2) and the resistivities are derived from the images through an algorithm almost identical to a procedure used in seismic prospecting (determination of the transmission coefficients) (III.3). As for the presentation of the results, resorting to the Dar-Zarrouk curve permits: a) to get a diagram somewhat similar to the E.S. curve (bilogarithmic scales coordinates: cumulative R and C) that is an already “smoothed” diagram where deeper layers show up less than superficial ones and b) to simplify the comprehensive solution. In fact, in arithmetic scales (R versus C) the Dar-Zarrouk curve consists of a many-sided polygonal contour which múst be replaced by an “equivalent” contour having a smaller number of sides. Though manually possible, this operation is automatically performed and additional constraints (e.g. geological information concerning thicknesses and resistivities) can be introduced at this stage. At present, the constraint used is the number of layers (III.4). Each solution (comprehensive and particular) is checked against the original data by calculating the E.S. diagrams corresponding to the distributions (thickness, resistivity) proposed. If the discrepancies are too large, the process is resumed (III.5). IV. Examples Several examples illustrate the procedure (IV). The first ones concern calculated E.S. diagrams, i.e. curves devoid of experimental errors and corresponding to a known distribution of resistivities and thicknesses (IV. 1). Example I shows how an E.S. curve is sampled. Several distributions (thickness, resistivity) were found: one is similar to, others differ from, the original one, although all E.S. diagrams are alike and characteristic parameters (transverse resistance of resistive layers and longitudinal conductance of conductive layers) are well determined. Additional informations must be introduced by the interpreter to remove the indeterminacy (IV.1.1). Examples 2 and 3 illustrate the principles of equivalence and suppression and give an idea of the sensitivity of the process, which seems accurate enough to make a correct distinction between calculated E.S. whose difference is less than what might be considered as significant in field curves (IV. 1.2 and IV. 1.3). The following example (number 4) concerns a multy-layer case which cannot be correctly approximated by a much smaller number of layers. It indicates that the result of the processing reflects correctly the trend of the changes in resistivity with depth but that, without additional information, several equally satisfactory solutions can be obtained (IV. 1.4). A second series of examples illustrates how the process behaves in presence of different kinds of errors on the original data (IV.2). A few anomalous points inserted into a series of accurate values of resistivities cause no problem, since the automatic processing practically replaces the wrong values (example 5) by what they should be had the E.S. diagram not been wilfully disturbed (IV.2.1). However, the procedure becomes less able to make a correct distinction, as the number of erroneous points increases. Weights must then be introduced, in order to determine the tolerance acceptable at each point as a function of its supposed accuracy. Example 6 shows how the weighting system used works (IV.2.2). The foregoing examples concern E.S. which include anomalous points that might have been caused by erroneous measurements. Geological effects (dipping layers for instance) while continuing to give smooth curves might introduce anomalous curvatures in an E.S. Example 7 indicates that in such a case the automatic processing gives distributions (thicknesses, resistivities) whose E.S. diagrams differ from the original curve only where curvatures exceed the limit corresponding to a horizontal stratification (IV.2.3). Numerous field diagrams have been processed (IV. 3). A first case (example 8) illustrates the various stages of the operation, chiefly the sampling of the E.S. (choice of the left cross, the weights and the resistivity of the substratum) and the selection of a solution, adapted from the automatic results (IV.3.1). The following examples (Nrs 9 and 10) show that electrical prospecting for deep seated layers can be usefully guided by the automatic processing of the E.S., even when difficult field conditions give original curves of low accuracy. A bore-hole proved the automatic solution proposed for E.S. no 10, slightly modified by the interpreter, to be correct.  相似文献   

20.
Following up our recent study of an indirect procedure for the practical determination of the maximum frequency-effect, defined as fe = 1 ? pρdc with ρ the resistivity at infinite frequency, we show at first how, through the Laplace transform theory, ρ can be related to stationary field vectors in the simple form of Ohm's law. Then applying the equation of continuity for stationary currents with a suitable set of boundary conditions, we derive the integral expression of the apparent resistivity at infinite frequency ρ,a in the case of a horizontally layered earth. Finally, from the definition of the maximum apparent frequency-effect, analytical expressions of feα are obtained for both Schlumberger and dipole arrays placed on the surface of the multi-layered earth section in the most general situation of vertical changes in induced polarization together with dc resistivity variations not at the same interfaces. Direct interpretation procedures are suggested for obtaining the layering parameters directly from the analysis of the sounding curves.  相似文献   

设为首页 | 免责声明 | 关于勤云 | 加入收藏

Copyright©北京勤云科技发展有限公司  京ICP备09084417号