CN110673223B - SIP observation method without synchronous current acquisition and transmission - Google Patents
SIP observation method without synchronous current acquisition and transmission Download PDFInfo
- Publication number
- CN110673223B CN110673223B CN201911009751.9A CN201911009751A CN110673223B CN 110673223 B CN110673223 B CN 110673223B CN 201911009751 A CN201911009751 A CN 201911009751A CN 110673223 B CN110673223 B CN 110673223B
- Authority
- CN
- China
- Prior art keywords
- phase
- frequency
- current
- time
- length
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Expired - Fee Related
Links
- 238000000034 method Methods 0.000 title claims abstract description 57
- 230000005540 biological transmission Effects 0.000 title claims abstract description 44
- 230000001360 synchronised effect Effects 0.000 title claims abstract description 33
- 238000004364 calculation method Methods 0.000 claims abstract description 49
- 238000010183 spectrum analysis Methods 0.000 claims abstract description 22
- 238000010606 normalization Methods 0.000 claims abstract description 20
- 230000005284 excitation Effects 0.000 claims abstract description 18
- 230000009466 transformation Effects 0.000 claims description 15
- 238000005070 sampling Methods 0.000 claims description 14
- 239000011159 matrix material Substances 0.000 claims description 9
- 239000000126 substance Substances 0.000 claims description 3
- 238000004422 calculation algorithm Methods 0.000 description 16
- 238000005259 measurement Methods 0.000 description 15
- 238000012545 processing Methods 0.000 description 7
- 238000012937 correction Methods 0.000 description 6
- 230000000694 effects Effects 0.000 description 4
- 238000012360 testing method Methods 0.000 description 4
- 238000001514 detection method Methods 0.000 description 3
- 238000001914 filtration Methods 0.000 description 3
- 238000003672 processing method Methods 0.000 description 3
- 238000004088 simulation Methods 0.000 description 3
- 238000001228 spectrum Methods 0.000 description 3
- 238000004458 analytical method Methods 0.000 description 2
- 238000013480 data collection Methods 0.000 description 2
- 230000010363 phase shift Effects 0.000 description 2
- 230000010287 polarization Effects 0.000 description 2
- 230000008054 signal transmission Effects 0.000 description 2
- 230000003595 spectral effect Effects 0.000 description 2
- NAWXUBYGYWOOIX-SFHVURJKSA-N (2s)-2-[[4-[2-(2,4-diaminoquinazolin-6-yl)ethyl]benzoyl]amino]-4-methylidenepentanedioic acid Chemical compound C1=CC2=NC(N)=NC(N)=C2C=C1CCC1=CC=C(C(=O)N[C@@H](CC(=C)C(O)=O)C(O)=O)C=C1 NAWXUBYGYWOOIX-SFHVURJKSA-N 0.000 description 1
- 238000000540 analysis of variance Methods 0.000 description 1
- 238000006243 chemical reaction Methods 0.000 description 1
- 238000013461 design Methods 0.000 description 1
- 230000005611 electricity Effects 0.000 description 1
- 230000007613 environmental effect Effects 0.000 description 1
- 238000011835 investigation Methods 0.000 description 1
- 238000000691 measurement method Methods 0.000 description 1
- 230000002093 peripheral effect Effects 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V3/00—Electric or magnetic prospecting or detecting; Measuring magnetic field characteristics of the earth, e.g. declination, deviation
- G01V3/38—Processing data, e.g. for analysis, for interpretation, for correction
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V3/00—Electric or magnetic prospecting or detecting; Measuring magnetic field characteristics of the earth, e.g. declination, deviation
- G01V3/02—Electric or magnetic prospecting or detecting; Measuring magnetic field characteristics of the earth, e.g. declination, deviation operating with propagation of electric current
Landscapes
- Life Sciences & Earth Sciences (AREA)
- Engineering & Computer Science (AREA)
- Environmental & Geological Engineering (AREA)
- Geology (AREA)
- Remote Sensing (AREA)
- Physics & Mathematics (AREA)
- General Life Sciences & Earth Sciences (AREA)
- General Physics & Mathematics (AREA)
- Geophysics (AREA)
- Measurement Of Resistance Or Impedance (AREA)
Abstract
The invention relates to an SIP observation method without synchronous current acquisition and transmission, which is characterized in that an absolute phase calculation method is utilized to calculate the initial phase and amplitude normalization coefficient of all frequency currents which can be transmitted by a transmitter in advance, and an initial phase table file is generated on the basis of the initial phase and amplitude normalization coefficient; generating a radio frequency table and an initial phase table; the phase value of the excitation current with a certain frequency at any moment is obtained by utilizing the initial emission time and the acquisition time of the frequency voltage data according to a relative phase calculation method; and automatically constraining the voltage amplitude and the phase obtained by the full-phase frequency spectrum analysis according to the length, namely obtaining the complex resistivity at the frequency point. The invention can get rid of the dependence of the system on high-performance synchronous current acquisition and real-time transmission equipment or a high-precision frequency domain observation device, improves the flexibility of the SIP method, and realizes the smooth development of SIP method observation in complex environments such as weak GPS signals.
Description
Technical Field
The invention belongs to the field of geophysical electromagnetic detection, and particularly relates to an SIP observation method without synchronous current acquisition and transmission.
Background
The frequency domain induced polarization method (SIP method) is an important branch of the induced polarization method, and the observation target is the complex resistivity of a measurement area. The method can describe the electrical structure of the underground medium of the test area more thoroughly, and is an important detailed investigation means. The SIP method can realize high-density observation in both a time domain and a frequency domain, and has the advantages of multi-parameter measurement and comprehensive information output compared with other geophysical exploration means.
At present, the existing SIP data acquisition modes and corresponding observation systems are mainly divided into two types:
(1) directly collecting the amplitude and phase of the voltage signal and the current signal, and then calculating the complex resistivity of the measuring area by using the frequency domain parameters;
(2) and synchronously acquiring time domain data of the voltage signal and the current signal, then performing Discrete Fourier Transform (DFT), and solving the complex resistivity of the measurement area by using a frequency domain parameter result of the transformation.
The first method needs to realize high-precision frequency domain measurement, and the observation system mainly comprises a high-power transmitter, a high-precision frequency domain measurement receiver, electrodes and other devices. The receiver is composed of key modules such as an imaginary component channel, a real component channel, a relevant detection unit and the like. The existing phase measurement method has high requirements on a transmitter and a receiver. Firstly, the waveform distortion of the excitation current needs to be as small as possible, the zero point needs to be stable, and meanwhile, the receiver needs to have high zero point stability. However, electromagnetic signals observed in the field are generally weak, and the phase measurement accuracy is difficult to guarantee. In addition, although the amplitude measurement is easier than the phase measurement, when the environmental noise is large, the supply current needs to be increased to improve the signal-to-noise ratio when the high-precision amplitude measurement is completed. In general, the frequency domain data acquisition method requires a relatively high accuracy and stability of the observation system, which increases the observation cost and is not conducive to field work.
The second observation method comprises a high-power transmitter, a receiver, a current detection device, a synchronous current transmission device, electrodes, a matched cable and other related peripheral equipment. The system has high requirements on instrument hardware, particularly the synchronism and the real-time performance of the system, the noise resistance of an auxiliary device and the like.
Disclosure of Invention
The invention aims to simplify the traditional SIP observation system, improve the flexibility and the application range of the SIP method and provide a new SIP data acquisition and processing algorithm. The dependence of the system on high-performance synchronous current acquisition and real-time transmission equipment or a high-precision frequency domain observation device is eliminated, and the observation cost is reduced; the flexibility of the SIP method is improved, and the SIP method observation can be smoothly carried out in complex environments such as weak GPS signals.
In order to achieve the purpose, the technical scheme of the invention is as follows:
a SIP observation method without synchronous current acquisition and transmission is characterized in that: generating an initial phase table, calculating a relative phase and calculating an absolute phase; the SIP observation method comprises the following steps:
firstly, calculating initial phases and amplitude normalization coefficients of all frequency currents which can be transmitted by a transmitter in advance by using an absolute phase calculation method, and generating an initial phase table file, namely a radio frequency table and an initial phase table, on the basis of the initial phases and the amplitude normalization coefficients;
secondly, by using a relative phase calculation method, passing a certain frequency fnThe initial emission time and the voltage data acquisition time to obtain the certain frequency fnThe phase value of the excitation current at any time;
finally, according to the length, automatically constraining the full-phase frequency spectrum analysis to obtain the voltage amplitude and the phase, and obtaining the certain frequency fnThe complex resistivity of (d);
the length automatic constraint full-phase spectrum analysis comprises three parts of full-phase matrix transformation, length automatic constraint self-convolution windowing and Fourier transformation.
The absolute phase calculation method comprises the following steps:
firstly, constructing a full phase matrix through original data; again by the frequency f of the excitation signalinAnd corresponding sampling rate fsAdding a data length constraint length N to obtain a self-convolution window with the length of N; finally, the self-convolution window and the full phase matrix are subjected to windowed Fourier transform to obtain the absolute phase of the signal;
wherein, the data length constraint length N constrains the self-convolution window length of the full phase transformation to meet the integral period truncation, and N is the frequency f of the excitation signalinAnd corresponding sampling rate fsTo decide:
wherein d is1And d2Are all integers.
The calculation steps for obtaining the voltage amplitude and the phase according to the length automatic constraint full-phase spectrum analysis are as follows:
(1) looking up a table of the initial phase table according to the frequency information to find a current initial phase and an amplitude normalization coefficient corresponding to the current frequency;
(2) acquiring time offset deltat between the current data time and the initial transmitting time of the current frequency signal according to the time information provided by the transmitting frequency table and the time for starting to acquire the current voltage data; in addition, in the A/D acquisition process, because 128 data points are lost due to hardware starting intervals, the time deviation caused by the lost data points needs to be considered when calculating the delta t, and the value of the time deviation and the current sampling rate fsIn connection with this, the Δ t of any time of the data segment is calculated as follows:
wherein f isinIs the excitation signal frequency; (ii) a
nskipIs the starting position of the currently calculated data segment in the voltage data block;
tini(fin) Is the current starting emission time of the current frequency;
tnow(fin) Is the current voltage signal start acquisition time;
(3) calculating the phase and amplitude of the voltage by utilizing length automatic constraint full-phase spectrum analysis, wherein the DFT length N is the same as the length of a full-phase transformation self-convolution window of the current signal with the same frequency;
(4) the phase of the current signal, which is synchronized with the voltage signal being processed, is calculated as follows:
(5) and (3) calculating the phase and amplitude of the complex resistivity of the current frequency point by using the formula (8) and the formula (9):
wherein, Phaseρ(fin) And
respectively a phase value of the current complex resistivity and a phase value of the current voltage;
wherein the content of the first and second substances,
is the amplitude of the voltage that is applied,
is the amplitude normalization coefficient, I (f)in) Is the current amperage;
and finally, correcting the calculation result to obtain the actually measured complex resistivity information.
The transmission frequency table is a frequency-time table and comprises the following information: signal type, transmission frequency, current strength, and transmission duration.
Compared with the prior art, the invention has the substantive characteristics and remarkable effects:
the invention does not need a synchronous current acquisition and transmission and frequency domain measurement device, the phase of the current is calculated by combining an absolute phase calculation method and a relative phase calculation method with a frequency time table, and the amplitude is obtained by a normalization means. The phase calculation can obtain a high-precision phase calculation result without an additional correction means, and has the advantages of high precision and simple algorithm compared with the traditional DFT and windowed DFT. It is characterized in that:
1. the traditional SIP observation system is simplified, the dependence of the system on high-performance synchronous current acquisition and real-time transmission equipment or a high-precision frequency domain observation device is eliminated, and the observation cost is reduced;
2. the traditional complex resistivity calculation method is improved, and the phase calculation precision is improved by adopting the length automatic constraint full-phase spectrum analysis
3. The flexibility of the SIP method is improved, and SIP method observation can be smoothly carried out in complex environments such as weak GPS signals.
Drawings
Fig. 1 is an absolute phase calculation method.
Fig. 2 is a flow of data acquisition and processing by the SIP method of the present invention.
Fig. 3 illustrates an initial phase table generation method according to the present invention.
FIG. 4 is a complex resistivity calculation correction algorithm.
FIG. 5 is a comparison of calculated and theoretical complex resistivities obtained by different algorithms.
FIG. 6 is a comparison graph of calculated values and theoretical values of SIP complex resistivity under power frequency interference.
Fig. 7 shows the SIP field measurement result. The upper graph is the complex resistivity amplitude curve and the lower graph is the phase curve.
Fig. 8 is a SIP observation data processing software interface.
Detailed Description
The following detailed description of the invention refers to the accompanying drawings
Firstly, the SIP observation method without synchronous current acquisition and transmission proves that the real part and the imaginary part of the complex resistivity are equivalent to the amplitude and the phase through the existing theoretical analysis, and the real part and the imaginary part can be mutually converted. Accordingly, the complex resistivity of a certain frequency point can be obtained by using the voltage amplitude and phase of the frequency and the amplitude and phase of the current signal synchronized with the voltage data of the segment. In addition, the amplitude and phase of the voltage signal can be directly obtained by performing time-frequency transformation on the original sampling data. But because of the absence of synchronous current data, the spectral information of the current signal cannot be calculated in this way.
In this case, the present invention uses an absolute phase calculation method and a relative phase calculation method in combination with a frequency schedule to find the phase of the current, and the amplitude is obtained by a normalization means. The absolute phase calculation method is realized by utilizing length automatic constraint full-phase spectrum analysis.
The SIP observation method without synchronous current acquisition and transmission can obtain a high-precision phase calculation result without an additional correction means, and has the advantages of high precision and simple algorithm compared with the traditional DFT and windowed DFT.
An SIP observation method without synchronous current acquisition and transmission of the present invention utilizes an absolute phase calculation method, including a relative phase calculation, which refers to a phase shift associated with a current transmission time, and an absolute phase calculation, which refers to an initial phase of a current signal.
Firstly, calculating initial phase and amplitude normalization coefficients of all frequency currents which can be transmitted by a transmitter, and generating an initial phase table file on the basis of the initial phase normalization coefficients, namely generating a radio frequency table and an initial phase table;
secondly, by using a relative phase calculation method, passing a certain frequency fnThe initial emission time and the voltage data acquisition time to obtain the certain frequency fnThe phase value of the excitation current at any time;
and finally, automatically constraining the voltage amplitude and the phase obtained by the full-phase frequency spectrum analysis according to the length to obtain the complex resistivity at the frequency point.
For a SIP transmitter, a certain frequency current x (f) is outputin) Is constant, at least for a short period of time. Thus, the phase of the signal at any time can be calculated by equation (1).
phase_It(fin)＝mod(Δt×fin)×2π+phase_Iini(fin) (1)
Wherein: f. ofinAnd phase _ Iini(fin) Respectively, the frequency and the initial phase of the current, and deltat is the time difference from the moment when the current starts to transmit to the moment when the data is collected.
The above analysis shows that once the amplitude and frequency of the transmit current are fixed, its frequency domain parameters are known at any time. In addition, the frequency domain information of the voltage signal can be obtained by time-frequency conversion using the original data. Thus, the complex resistivity parameter can be calculated according to the formula (2) and the formula (3).
phase_ρ(fin)＝phase_U(fin)-Δphase_I(fin)-phase_Iini(fin)
＝phase_U(fin)-mod(Δt×fin)×2π-phase_Iini(fin) (2)
Amp_ρ(fin)＝Amp_U(fin)/Amp_Iini(fin). (3)
Formula (2) shows that the complex resistivity is at frequency finThe phase position is composed of a voltage phase _ U at the frequency point and a starting phase _ I of the frequency currentiniFrequency finAnd a transmission time Δ t (transmission start time t)iniAnd voltage data acquisition time tnThe difference therebetween) are determined in combination. While the amplitude can be obtained using a normalization method.
Thus, without a synchronous current data acquisition and transmission device, the following parameters are used: current frequency finCurrent intensity (amplitude) A, initial phase _ IiniAnd a transmission time deltat (in relation to the transmission start time), frequency domain information of the present excitation current signal can be determined. This is the theoretical basis for the present invention to design the real-time algorithm of the SIP data.
The key problem of the SIP observation method without synchronous current acquisition and transmission is phase calculation, including relative phase calculation and absolute phase calculation. The former refers to the phase shift associated with the current transit time, while the latter refers to the initial phase of the current signal.
In the case of a stable and known signal frequency, the complex resistivity phase calculation error may come from the voltage phase error phase _ Uerr(fin) Current initial phase error phase _ Iinierr(fin) And a time error Δ terr. The first two are determined by the accuracy of the algorithm, and the third term is controlled by the accuracy of the time recording. Both the current start phase and the voltage phase are obtained using an absolute phase calculation method.
The most basic and most common phase calculation method is the fourier transform, and for discrete signals the Discrete Fourier Transform (DFT). Both the direct DFT and the windowed DFT are affected to some extent by the effects of spectral leakage. When there is no frequency error or the frequency error is small, the windowing makes the spectrum leakage effect more obvious. In addition, the window function can inhibit side lobe leakage and simultaneously cause main lobe blurring, and the effect is particularly obvious in multi-frequency signal processing with dense frequency components.
In order to obtain more accurate spectrum parameters, particularly phases, the method adopts length automatic constraint full-phase spectrum analysis, and has the advantages of small calculation amount and high phase precision. Because the invention has very high requirement on the phase calculation precision, some constraint conditions are added to the full phase transformation.
The length automatic constraint full-phase spectrum analysis comprises three parts of full-phase matrix transformation, length automatic constraint self-convolution windowing and Fourier transformation.
As shown in fig. 1, the absolute phase calculation method of the SIP observation method without synchronous current collection and transmission of the present invention includes: firstly, constructing a full phase matrix through original data; again by the frequency f of the excitation signalinAnd corresponding sampling rate fsAdding a data length constraint length N to obtain a self-convolution window with the length of N points; and finally, obtaining the absolute phase of the signal by windowing Fourier transform from the convolution window and the full phase matrix.
The data length constraint length N is used for constraining the self-convolution window length of the full phase transformation to meet the requirement of the whole period truncation, and because the excitation signal of the SIP method is mostly a series of square waves with different frequencies, the data length constraint length N should be determined by the frequency f of the excitation signalinAnd corresponding sampling rate fsTo decide:
wherein d is1And d2Are all integers. d1The parameter d is required to satisfy the condition that the value in the symbol | | | | is an integer2To ensure that the window size is greater than a predetermined threshold σ, which is related to the sampling rate, e.g. table1 is shown.
TABLE 1 Window Length threshold without sampling Rate
Order to
Then d2Calculated according to the following formula:
where ceil () represents rounding up.
The specific implementation flow of the SIP observation method without synchronous current acquisition and transmission is as follows:
as shown in fig. 2, the SIP data acquisition and processing flow chart firstly generates a transmission frequency table and an initial phase table, and secondly establishes a "frequency-time table" (hereinafter referred to as "transmission frequency table"), which contains the following information: signal type, transmission frequency, current intensity, transmission time and the like. The transmit frequency table is then stored in the transmitter and acquisition station. After the observation system starts to work, the transmitter starts to circularly transmit current in sequence according to the transmission frequency table. The receiver can then start collecting voltage data at any time. This is less limited than the current common SIP data collection method. To ensure that voltage signals of all frequencies can be acquired, the total acquisition time cannot be less than the total time of one cycle. The calculation steps for obtaining the voltage amplitude and the phase according to the length automatic constraint full-phase spectrum analysis are as follows:
(1) looking up a table for an initial phase table according to the frequency information to find a current initial phase and an amplitude normalization coefficient corresponding to the current frequency;
(2) according to the time information provided by the transmitting frequency table and the time for starting to collect the current voltage data, the time offset between the current data time and the initial transmitting time of the current frequency signal can be obtainedThe quantity of the shift delta t is increased, in addition, in the A/D acquisition process, because 128 data points are lost due to the hardware starting interval, the time deviation caused by the lost data points needs to be considered when the delta t is calculated, and the value of the time deviation and the current sampling rate fsIn connection with this, the Δ t of any time of the data segment is calculated as follows:
wherein f isinIs the excitation signal frequency;
nskipis the starting position of the currently calculated data segment in the voltage data block;
tini(fin) Is the current starting emission time of the current frequency;
tnow(fin) Is the current voltage signal start acquisition time;
(3) calculating the phase and amplitude of the voltage by utilizing length automatic constraint full-phase spectrum analysis, wherein the DFT length N is the same as the length of a full-phase transformation self-convolution window of the current signal with the same frequency;
(4) the phase of the current signal, which is synchronized with the voltage signal being processed, is calculated as follows:
(5) and (3) calculating the phase and amplitude of the complex resistivity of the current frequency point by using the formula (8) and the formula (9):
wherein
And
current voltage and electricity respectivelyThe phase value of the stream.
Wherein the content of the first and second substances,
is the amplitude of the voltage that is applied,
is the amplitude normalization coefficient, I (f)in) Is the current intensity.
And finally, correcting the calculation result to obtain the actually measured complex resistivity information.
As shown in fig. 3, the initial phase table of the SIP observation method without synchronous current collection and transmission according to the present invention is generated by collecting signals of all frequencies generated by the SIP transmitter before observing by the SIP method, calculating their initial phase and normalized amplitude coefficient by using the length-constrained full-phase spectrum analysis, and generating an "initial phase-amplitude-frequency table" (hereinafter referred to as "initial phase table") based on the signals, wherein the table includes the following information: signal frequency, initial phase, amplitude normalization factor, and length constraint value N. The initial phase table file format is shown in table 2. Where the first column is the excitation current frequency and the second and third columns are the initial phase and amplitude normalization coefficients of the directly sampled data. The fourth and fifth columns are initial phase and amplitude normalization coefficients of the data processed with the harmonic filter. Because the field observation environment is complex, power frequency filtering is needed firstly when power frequency interference is strong, and a filter can cause certain influence on the amplitude and the phase of original sampling data, especially the phase, so that the frequency spectrum information after the emission current filtering is added into the initial frequency table to meet the requirement of field real-time data processing.
TABLE 2 "initial phase Table" file format
And (5) an initial phase table generation process. The data processing procedure in the two dashed boxes is identical. Table 3 is an example of a transmit frequency table. FD represents a square wave of the signal waveform. Each frequency signal is continuously transmitted circularly. The transmission end time refers to a difference between the current frequency signal transmission end time and the first frequency signal transmission start time. For example, assuming that the 128Hz signal in Table 3 begins to be transmitted at time 00:00:00, the frequency signal ends at 00:00:50, while the 64Hz signal begins to be transmitted to 00:01:40, and so on.
TABLE 3 Transmit frequency Table
Because the measurement result of the frequency domain electromagnetic method contains the system self response, the calculation result needs to be corrected through calibration data, and therefore the real complex resistivity information of the measurement area can be finally obtained.
As shown in fig. 4, the complex resistivity calculation results are corrected. Because the frequency values used in each observation are not necessarily identical, and the calibrated frequency points are limited, the correction coefficient needs to be calculated by an interpolation method. According to the signal frequency finLooking up the calibration data table to find finThen calculating f by interpolationinAnd finally, correcting the original calculation result by using the calculated correction coefficient. The interpolation algorithm adopted by the invention is parabolic interpolation, so that table lookup is firstly carried out to find three and finClosest frequency point: f. ofk-1、fkAnd fk+1. Then, the interpolation result is calculated according to the formula:
where r (f) represents the amplitude or phase calibration result at frequency f.
The complex resistivity correction formula is as follows:
φ(fin)＝φ(fin)-φc(fin)
ρ(fin)＝ρ(fin)/ρc(fin) (11)
wherein phi isc(fin) And ρc(fin) Respectively at frequency f for the observation systeminThe amplitude and phase of the self-response are calculated according to the formula (10).
The data acquisition and processing method flow is shown in fig. 1.
The comparison graph of the complex resistivity calculated value and the theoretical value can be obtained by different algorithms through simulation experiments.
As shown in fig. 5, the solid line represents the theoretical complex resistivity curve; the dotted lines represent the results obtained using the synchronous voltage and current data, and the discrete points represent the new method calculation results; "+," Δ ", and" □ "represent the results of the computation of the length auto-constrained full-phase spectral analysis, windowed DFT, and direct DFT, respectively. The upper graph is the complex resistivity amplitude and the lower graph is the phase.
In order to verify the reliability of the algorithm, simulated power frequency interference is added on the basis of the signal I to serve as excitation current. The power frequency interference consists of a set of sine waves with 50Hz and odd multiples thereof. The above simulation and data processing procedure is then repeated.
As shown in fig. 6, the complex resistivity calculation results are shown in fig. 6, and the algorithms represented by the different symbols and line types are consistent with fig. 5. The adopted data processing method is similar to that under the condition of no harmonic interference, and only power frequency filtering is added in the early stage.
The complex resistivity calculation accuracy of the above methods was discussed using an analysis of variance method. And the results are summarized in table 4. Simulation tests prove that under the condition of no power frequency interference, the calculation accuracy of the traditional SIP data processing method requiring synchronous current is higher than that of the novel algorithm provided by the invention, but under the condition of power frequency interference, the novel algorithm is superior to the traditional algorithm. For the new algorithm, the amplitude precision obtained by direct DFT, windowed DFT and length automatic constrained full-phase spectrum analysis is basically consistent, but the phase precision obtained by the length automatic constrained full-phase spectrum analysis is obviously higher than the former two.
Table 4 mean square error of different SIP data algorithms
As shown in fig. 7, it is a field test result of an SIP observation method without synchronous current collection and transmission according to the present invention, and three-channel voltage data is obtained by using table 3 as a transmission frequency table and data collection time about half an hour.
The complex resistivity amplitude and phase curve measured by the SIP method field test is smooth, three channels are basically consistent, and the phase of the high frequency band of the channel 1 is slightly lower than the measurement results of the other two channels. The results show that:
the SIP observation method without synchronous current acquisition and transmission can smoothly complete SIP observation under the condition of no synchronous current acquisition and transmission equipment and no frequency domain measuring device. In addition, the conventional SIP method requires that the voltage data and the current data must be acquired strictly and synchronously, and in comparison,
the SIP data acquisition method only acquires voltage data without synchronous current acquisition and transmission is more flexible and convenient.
As shown in fig. 8, the SIP observation data processing software interface shows the comparison between the observation results obtained by the SIP method and the full-phase spectrum analysis with the automatic length constraint.
In fig. 8, (a) and (b) are time domain waveforms and frequency domain waveforms of the original voltage sample data. (c) And (d) amplitude and phase calculations of complex resistivity. (e) Is a parameter setting and display section, where Array is a device type, "4" represents a symmetric quadrupole device; AB is the feed electrode distance; MN is the measured pole pitch. These parametersAre manually input according to the layout mode of the device before the measurement is started. K is the device coefficient, automatically calculated by the program from the input parameters. f. ofinIs the current signal frequency, and is judged according to the data acquisition time and the transmitting frequency table, and the sampling rate fsIs according to finAnd (4) determining. The solid line, the dotted line and the dashed line respectively represent the channel 1/2/3, "+" represents the calculation result of the length automatic constraint full-phase spectrum analysis of the invention, and "Δ" represents the calculation result of the windowed fourier transform method.
Claims (4)
1. A SIP observation method without synchronous current acquisition and transmission is characterized in that: generating an initial phase table, calculating a relative phase and calculating an absolute phase; the SIP observation method comprises the following steps:
firstly, calculating initial phases and amplitude normalization coefficients of all frequency currents which can be transmitted by a transmitter in advance by using an absolute phase calculation method, and generating an initial phase table file, namely a radio frequency table and an initial phase table, on the basis of the initial phases and the amplitude normalization coefficients;
secondly, by using a relative phase calculation method, passing a certain frequency fnThe initial emission time and the voltage data acquisition time to obtain the certain frequency fnThe phase value of the excitation current at any time;
finally, according to the length, automatically constraining the full-phase frequency spectrum analysis to obtain the voltage amplitude and the phase, and obtaining the certain frequency fnThe complex resistivity of (d);
the length automatic constraint full-phase spectrum analysis comprises three parts of full-phase matrix transformation, length automatic constraint self-convolution windowing and Fourier transformation.
2. The SIP observation method without synchronous current collection and transmission according to claim 1, wherein: the absolute phase calculation method comprises the following steps:
firstly, constructing a full phase matrix through original data; again by the frequency f of the excitation signalinAnd corresponding sampling rate fsAdding data length constraintsObtaining a self-convolution window with the length of N after the length of N is obtained; finally, the self-convolution window and the full phase matrix are subjected to windowed Fourier transform to obtain the absolute phase of the signal;
wherein, the data length constraint length N constrains the self-convolution window length of the full phase transformation to meet the integral period truncation, and N is the frequency f of the excitation signalinAnd corresponding sampling rate fsTo decide:
wherein d is1And d2Are all integers.
3. The SIP observation method without synchronous current collection and transmission according to claim 1, wherein: the calculation steps for obtaining the voltage amplitude and the phase according to the length automatic constraint full-phase spectrum analysis are as follows:
(1) looking up a table of the initial phase table according to the frequency information to find a current initial phase and an amplitude normalization coefficient corresponding to the current frequency;
(2) acquiring time offset deltat between the current data time and the initial transmitting time of the current frequency signal according to the time information provided by the transmitting frequency table and the time for starting to acquire the current voltage data; in addition, in the A/D acquisition process, because 128 data points are lost due to hardware starting intervals, the time deviation caused by the lost data points needs to be considered when calculating the delta t, and the value of the time deviation and the current sampling rate fsIn connection with this, the Δ t of any time of the data segment is calculated as follows:
wherein f isinIs the excitation signal frequency;
nskipis the starting position of the currently calculated data segment in the voltage data block;
tini(fin) Is the current starting emission time of the current frequency;
tnow(fin) Is the current voltage signal start acquisition time;
(3) calculating the phase and amplitude of the voltage by utilizing length automatic constraint full-phase spectrum analysis, wherein the DFT length N is the same as the length of a full-phase transformation self-convolution window of the current signal with the same frequency;
(4) the phase of the current signal, which is synchronized with the voltage signal being processed, is calculated as follows:
(5) and (3) calculating the phase and amplitude of the complex resistivity of the current frequency point by using the formula (8) and the formula (9):
wherein, Phaseρ(fin) And
respectively a phase value of the current complex resistivity and a phase value of the current voltage;
wherein the content of the first and second substances,
is the amplitude of the voltage that is applied,
is the amplitude normalization coefficient, I (f)in) Is the current amperage;
and finally, correcting the calculation result to obtain the actually measured complex resistivity information.
4. The SIP observation method without synchronous current collection and transmission according to claim 3, wherein: the transmission frequency table is a frequency-time table and comprises the following information: signal type, transmission frequency, current strength, and transmission duration.
Applications Claiming Priority (2)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201910059976.9A CN109581518A (en) | 2019-01-22 | 2019-01-22 | One kind is without synchronizing current acquisition transmission and the acquisition of frequency domain measurement SIP data and processing method |
| CN2019100599769 | 2019-01-22 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN110673223A CN110673223A (en) | 2020-01-10 |
| CN110673223B true CN110673223B (en) | 2021-06-15 |
Family
ID=65917724
Family Applications (2)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201910059976.9A Pending CN109581518A (en) | 2019-01-22 | 2019-01-22 | One kind is without synchronizing current acquisition transmission and the acquisition of frequency domain measurement SIP data and processing method |
| CN201911009751.9A Expired - Fee Related CN110673223B (en) | 2019-01-22 | 2019-10-23 | SIP observation method without synchronous current acquisition and transmission |
Family Applications Before (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201910059976.9A Pending CN109581518A (en) | 2019-01-22 | 2019-01-22 | One kind is without synchronizing current acquisition transmission and the acquisition of frequency domain measurement SIP data and processing method |
Country Status (1)
| Country | Link |
|---|---|
| CN (2) | CN109581518A (en) |
Families Citing this family (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN111707714B (en) * | 2020-06-24 | 2024-07-09 | 中国地质大学（武汉） | Soil organic chloride pollution rapid investigation device and method based on complex resistivity |
| CN115951836B (en) * | 2023-01-12 | 2023-09-05 | 上海奎芯集成电路设计有限公司 | DFI command and data channel universal controller and read-write method thereof |
Family Cites Families (9)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN1248015C (en) * | 2003-11-25 | 2006-03-29 | 中国地质大学（北京） | Broad spectrum current drive measuring system and measuring mode |
| CN2781399Y (en) * | 2005-02-01 | 2006-05-17 | 湖南继善高科技有限公司 | Chopper decoupling device for induced polarization measuring |
| CN100356194C (en) * | 2005-02-01 | 2007-12-19 | 湖南继善高科技有限公司 | Chopping decoupler with exciting polarizing measurement |
| CN101986098B (en) * | 2010-09-21 | 2012-02-22 | 东南大学 | Tricolor projection-based Fourier transform three-dimensional measuring method |
| US10209386B2 (en) * | 2012-08-30 | 2019-02-19 | Exxonmobil Upstream Research Company | Processing methods for time division CSEM data |
| US10241101B2 (en) * | 2014-01-24 | 2019-03-26 | Schlumberger Technology Corporation | Method and apparatus for determining permittivity of rock matrix |
| CN105301610B (en) * | 2015-09-17 | 2017-12-22 | 西安空间无线电技术研究所 | A kind of Novel GPS L5 signal quick catching methods of anti-symbol saltus step |
| CN106842332B (en) * | 2016-12-25 | 2018-06-22 | 中南大学 | A kind of comprehensive electrical prospecting method |
| CN109149968B (en) * | 2018-08-27 | 2020-07-28 | 重庆西南集成电路设计有限责任公司 | Synchronous rectifier diode and synchronous rectification control circuit |
-
2019
- 2019-01-22 CN CN201910059976.9A patent/CN109581518A/en active Pending
- 2019-10-23 CN CN201911009751.9A patent/CN110673223B/en not_active Expired - Fee Related
Also Published As
| Publication number | Publication date |
|---|---|
| CN109581518A (en) | 2019-04-05 |
| CN110673223A (en) | 2020-01-10 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| Ritzmann et al. | Comparison of measurement methods for 2–150-kHz conducted emissions in power networks | |
| Sleeman et al. | Three-channel correlation analysis: A new technique to measure instrumental noise of digitizers and seismic sensors | |
| CN101603985B (en) | Method for measuring sine signal with high accuracy | |
| CN101701984B (en) | Fundamental wave and harmonic wave detecting method based on three-coefficient Nuttall windowed interpolation FFT | |
| DK2406643T3 (en) | Method and system for monitoring of time domain interference | |
| CN105656485B (en) | A kind of multichannel time-interleaved AD C measures calibration method and device | |
| CN110673223B (en) | SIP observation method without synchronous current acquisition and transmission | |
| CN105137185A (en) | Frequency domain interpolation electric power harmonic wave analysis method based on discrete Fourier transform | |
| EP0107050A2 (en) | Signal generator for digital spectrum analyzer | |
| CN110244116B (en) | DC instantaneous power metering circuit and quasi-synchronous calculation method thereof | |
| CN110967658B (en) | Analog input merging unit calibrator tracing method based on digital differential method | |
| Liu et al. | Improved processing of harmonics and interharmonics by time-domain averaging | |
| WO2015154587A1 (en) | Measurement apparatus and method for ripple current | |
| CN110133564A (en) | A kind of synchronous detecting method and device of current sensor phase characteristic | |
| WO2007068972B1 (en) | Measuring electrical impedance at various frequencies | |
| CN105675126A (en) | Novel method for detecting sound pressure of multi-frequency multi-source complex stable sound field | |
| CN101718816B (en) | Fundamental wave and harmonic wave detection method based on four-item coefficient Nuttall window interpolation FFT | |
| CN104849569B (en) | Dielectric loss measuring method | |
| CN101308175A (en) | Phase spectrum analyzer | |
| CN101806835B (en) | Interharmonics measuring meter based on envelope decomposition | |
| CN109270404A (en) | A kind of voltage traveling wave accurate detecting method and device | |
| RU2695953C1 (en) | Method of estimating signal-to-noise ratio at input of receiving device for radio signal with digital amplitude modulation | |
| CN106603166B (en) | Vector measurement device and method for broadband modulation signal | |
| US9350470B1 (en) | Phase slope reference adapted for use in wideband phase spectrum measurements | |
| CN113341219B (en) | Method and device for measuring frequency modulation amplitude of constant-frequency alternating current power supply system |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant | ||
| CF01 | Termination of patent right due to non-payment of annual fee |
Granted publication date: 20210615 |