CN104656140B - Median filtering method and system based on vertical earthquake attenuation laws - Google Patents
Median filtering method and system based on vertical earthquake attenuation laws Download PDFInfo
- Publication number
- CN104656140B CN104656140B CN201410718012.8A CN201410718012A CN104656140B CN 104656140 B CN104656140 B CN 104656140B CN 201410718012 A CN201410718012 A CN 201410718012A CN 104656140 B CN104656140 B CN 104656140B
- Authority
- CN
- China
- Prior art keywords
- data
- sequence
- section
- vsp
- amplitude
- 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.)
- Active
Links
Landscapes
- Geophysics And Detection Of Objects (AREA)
Abstract
The invention provides a median filtering method and system based on vertical earthquake attenuation laws. The method comprises the following steps: collecting VSP (Vertical Seismic Profiling) data; selecting single component data according to the VSP data; carrying out zero-phasing treatment on the single component data, so as to obtain a profile treated through zero-phasing treatment; obtaining a time curve according to the profile treated through zero-phasing treatment; determining spherical spreading compensation TAR factors according to the profile treated through zero-phasing treatment and the time curve; pulling the profile treated through zero-phasing treatment to be homogeneous along the time curve, and obtaining a homogeneously pulled profile; carrying out directional weight median filtering on the homogeneously pulled profile according to the TAR factors of the VSP data. The median filtering method and system, provided by the invention, solve the problem of no reliable wave field separation method in VSP data processing used for well control.
Description
Technical field
The present invention with regard to technical field of geophysical exploration, at vertical seismic data in geophysical exploration
Reason neck technology, is concretely a kind of median filter method based on vertical seismic attenuation law and system.
Background technology
Wave field abundant information is an important feature of vertical seismic profiling (VSP) VSP data, can be seen in section various
The phenomenon that wave field is interweaved.Generally, VSP data do not possess the higher degree of covering of surface seismic again, so it is difficult to pass through to fold
Plus method reinforcing or a certain wave field of decaying.Therefore, application of the wave field separation method in the process of VSP data just shows especially
For important.
In prior art, not targetedly wave field separation method in the processing scheme of VSP data, its wave field separation hand
Section is all the method that surface seismic data is processed that borrows, such as frequency-wavenumber filtering (FK filtering), mean filter and medium filtering
Deng.Wherein, FK filtering capacities are also easy to produce alias, and the application in the VSP data wave fields that degree of covering is low, quality requirement is high extremely has
Limit；Mean filter and median filter method belong to the multiple tracks pie slice of time-domain, and the anti-noise ability of medium filtering is better than equal
Value filtering, applies more, and derives self adaptation and go the median filter methods such as exception, border lane change number in actual production.
In VSP Data processings, although the characteristics of these wave field separation methods do not consider VSP data itself, but traditional
In handling process, wave field separation step is carried out after amplitude compensation, and now section amplitude energy is relatively uniform, can be
Wave field separation is carried out using median method.
In recent years, VSP data are developed rapidly in terms of well control technique, and its amplitude attenuation parameter, spherical diffusion are mended
Repay parameter (TAR), formation absorption compensating parameter (Q) etc. to widely apply in surface seismic process, serve good well control about
Shu Zuoyong.However, the extraction of these well control parameters is required for being completed under the premise of the relative fidelity of VSP data amplitudes, it is desirable to shaking
Wave field separation is carried out in the section that width exponentially changes with VSP data record depth.Generally, the diverse location of a lineups is deposited
In very big amplitude energy difference, this is by the reliability of strong influence median filter method.
Therefore, how to propose a kind of new filtering method, and then solve not have in processing for the VSP Data Datas of well control
Reliable wave field separation is this area technical barrier urgently to be resolved hurrily.
The content of the invention
It is used for asking without reliability wave field separation method during the VSP Data Datas of well control are processed in order to solve prior art
Topic, the present invention is proposed a kind of based on vertical seismic decay on the basis of self adaptation removes abnormal, border lane change number medium filtering
The median filter method and system of rule, is a kind of scheme of the accurate medium filtering based on vertical seismic attenuation law, by
VSP data obtain simple component data, and to simple component data dephasing process is carried out, and determine that spherical diffusion is mended with reference to time curve
(TAR) factor is repaid, and section is pulled together along time curve, be oriented Weighted median filtering, greatly improve medium filtering can
By property.
It is an object of the invention to provide a kind of median filter method based on vertical seismic attenuation law, including：Adopt
Collection vertical seismic profiling (VSP) VSP data；Simple component data are sub-elected according to described VSP data；Described simple component data are entered
The process of row dephasing, obtains the section after dephasing is processed；According to the dephasing process after section obtain when away from
Curve；The spherical diffusion of the VSP data is determined according to the section and described time curve after dephasing process
The compensation TAR factors；Section after dephasing is processed is pulled together along the time curve, the section after being pulled together；According to institute
Section after pulling together described in the TAR factor pairs for stating VSP data is oriented Weighted median filtering.
It is an object of the invention to provide a kind of system of the medium filtering based on vertical seismic attenuation law, bag
Include：Profile collection device, for gathering vertical seismic profiling (VSP) VSP data；Simple component data sorting unit, for according to described
VSP data sub-elect simple component data；Dephasing processing meanss, for carrying out at dephasing to described simple component data
Reason, obtains the section after dephasing is processed；Time curve determining device, for according to the section after dephasing process
Obtain time curve；TAR factor determining devices, for according to the section after dephasing process and it is described when away from song
Line determines the spherical diffusion compensation TAR factor of the VSP data；Device is pulled together, for the section edge after dephasing is processed
The time curve is pulled together, the section after being pulled together；Medium filtering device, for according to the TAR factor pairs of the VSP data
It is described pull together after section be oriented Weighted median filtering.
The beneficial effects of the present invention is, on the basis of self adaptation removes abnormal, border lane change number medium filtering, there is provided
A kind of median filter method and system based on vertical seismic attenuation law is a kind of accurately based on vertical seismic decay rule
The scheme of the medium filtering of rule, by VSP data simple component data are sub-elected, and simple component data are carried out with dephasing process, knot
Close time curve and determine the TAR factors, and section is pulled together along time curve, Weighted median filtering is oriented, in greatly improving
The reliability of value filtering, to provide reliable wave field separation scheme for the process of the VSP Data Datas of well control.
It is that the above and other objects, features and advantages of the present invention can be become apparent, preferred embodiment cited below particularly,
And coordinate institute's accompanying drawings, it is described in detail below.
Description of the drawings
In order to be illustrated more clearly that the embodiment of the present invention or technical scheme of the prior art, below will be to embodiment or existing
The accompanying drawing to be used needed for having technology description is briefly described, it should be apparent that, drawings in the following description are only this
Some embodiments of invention, for those of ordinary skill in the art, on the premise of not paying creative work, can be with
Other accompanying drawings are obtained according to these accompanying drawings.
Fig. 1 is a kind of flow process of median filter method based on vertical seismic attenuation law provided in an embodiment of the present invention
Figure；
Fig. 2 is the particular flow sheet of step S104 in Fig. 1；
Fig. 3 is the particular flow sheet of step S105 in Fig. 1；
Fig. 4 is the particular flow sheet of step S107 in Fig. 1；
Fig. 5 is the particular flow sheet of step S404 in Fig. 4；
Fig. 6 is the particular flow sheet of step S406 in Fig. 4；
Fig. 7 is a kind of structural frames of medium filtering system based on vertical seismic attenuation law provided in an embodiment of the present invention
Figure；
Fig. 8 be in a kind of medium filtering system based on vertical seismic attenuation law provided in an embodiment of the present invention when away from
The concrete structure block diagram of curve determining device；
Fig. 9 is the TAR in a kind of medium filtering system based on vertical seismic attenuation law provided in an embodiment of the present invention
The concrete structure block diagram of factor determining device；
Figure 10 be in a kind of medium filtering system based on vertical seismic attenuation law provided in an embodiment of the present invention in
The concrete structure block diagram of value filtering device；
Figure 11 is adding in a kind of medium filtering system based on vertical seismic attenuation law provided in an embodiment of the present invention
The concrete structure block diagram of power data sequence determining module；
Figure 12 is going in a kind of medium filtering system based on vertical seismic attenuation law provided in an embodiment of the present invention
The concrete structure block diagram of outlier processing module；
Original VSP data Z component schematic diagram in the specific embodiment that Figure 13 is provided for the present invention；
Seismic wavelet dephasing schematic diagram in the specific embodiment that Figure 14 is provided for the present invention；
VSP data first break pickup schematic diagrames in the specific embodiment that Figure 15 is provided for the present invention；
The result of calculation figure of the TAR factors calculated using VSP data in the specific embodiment that Figure 16 is provided for the present invention；
The calculating difference curve figure of TAR in the specific embodiment that Figure 17 is provided for the present invention；
Figure 18 in the specific embodiment that provides of the present invention along down going wave time curve even up after VSP schematic diagram datas；
Sampling point sequence diagram in the specific embodiment that Figure 19 is provided for the present invention in a window；
According to the weight coefficient sequence diagram of VSP data rules in the specific embodiment that Figure 20 is provided for the present invention；
Sampling point sequence diagram in the specific embodiment that Figure 21 is provided for the present invention in latter window of weighting；
The sampling point sequence diagram that self adaptation is gone after exception in the specific embodiment that Figure 22 is provided for the present invention；
Figure 23 in the specific embodiment that provides of the present invention directional weight, self adaptation go it is abnormal, border lane change number
VSP schematic diagram datas 1 after medium filtering；
Figure 24 in the specific embodiment that provides of the present invention directional weight, self adaptation go it is abnormal, border lane change number
VSP schematic diagram datas 2 after medium filtering.
Specific embodiment
Below in conjunction with the accompanying drawing in the embodiment of the present invention, the technical scheme in the embodiment of the present invention is carried out clear, complete
Site preparation is described, it is clear that described embodiment is only a part of embodiment of the invention, rather than the embodiment of whole.It is based on
Embodiment in the present invention, it is every other that those of ordinary skill in the art are obtained under the premise of creative work is not made
Embodiment, belongs to the scope of protection of the invention.
The present invention relates to vertical seismic data processing field in geophysical exploration, belongs to seism processing means
The category of wave field separation is carried out, is that a kind of on the basis of self adaptation goes exception and border lane change number medium filtering, application is vertical
The median filter method that geological data attenuation law is weighted.
Without reliability wave field separation method, the present invention in order to solve the problems, such as the VSP Data Datas process for well control
On the basis of self adaptation removes abnormal, border lane change number medium filtering, propose that one kind carries out adding using VSP data attenuation law
The brand-new median filter method of power.Fig. 1 is that a kind of intermediate value based on vertical seismic attenuation law provided in an embodiment of the present invention is filtered
The flow chart of wave method, as shown in Figure 1, described method includes：
S101：Collection vertical seismic profiling (VSP) VSP data.
S102：Simple component data are sub-elected according to described VSP data.
In the particular embodiment, VSP acquired original data are arranged using existing method and seismic processing software,
Acquisition is ready for the simple component data of wave field separation, as shown in figure 13, original VSP in the specific embodiment provided for the present invention
Data Z component schematic diagram.
Described data preparation does not include that amplitude compensation is processed, trace equalization is processed, automatic gain (AGC) is processed or other can
The process step of amplitude or energy relativeness can be changed.
S103：Dephasing process is carried out to described simple component data, the section after dephasing is processed is obtained.
In the particular embodiment, dephasing process is carried out using existing software.Dephasing process purpose be by
The complicated wavelet of mixed-phase is converted into the wavelet of zero phase, makes the crest or trough of wavelet corresponding with bed boundary, such as schemes
Shown in 14, seismic wavelet dephasing schematic diagram in the specific embodiment provided for the present invention.
S104：Time curve is obtained according to the section after dephasing process.
Fig. 2 is the particular flow sheet of step S104, as shown in Figure 2, is obtained according to the section after dephasing process
Time curve includes：
S201：Using the crest location of descending direct wave in the section after dephasing process as take-off location；
S202：The descending direct wave is picked up with the take-off location, time curve is obtained.
In a particular embodiment, descending direct wave first arrival, m- depth when obtaining are picked up using existing seismic processing software
Degree relational sequence, the sequence constitutes in the plane the relation curve of a time and distance, hereinafter referred to as time curve, such as schemes
Shown in 15, VSP data first break pickup schematic diagrames in the specific embodiment provided for the present invention.
The first break pickup position is the crest location of descending direct wave in section.
The down going wave time curve is simultaneously as the time curve of follow-up spherical diffusion compensation (TAR) factor.
As shown in Figure 1, the method also includes：
S105：The VSP data are determined according to the section and described time curve after dephasing process
The spherical diffusion compensation TAR factor.Fig. 3 is the particular flow sheet of step S105, from the figure 3, it may be seen that the step is specifically included：
S301：In section after dephasing process, the amplitude of the time curve is picked up, it is m- when obtaining
Amplitude relation sequence；
S302：To it is described when m- amplitude relation sequence in time and amplitude ask for natural logrithm respectively, obtain new sequence
Row；
S303：Linear fit is carried out to described new sequence, fitting a straight line is obtained；
S304：Determine the slope of the fitting a straight line；
S305：Opposite number is taken to the slope, the spherical diffusion compensation TAR factor of the VSP data is obtained.
In a particular embodiment, in the section obtained by step S103, step S104 time curve position is picked up
Amplitude, m- amplitude relation sequence during acquisition asks for respectively natural logrithm, then to institute to time in the sequence and amplitude
The new sequence for obtaining carries out linear fit, is calculated the slope of fitting a straight line, and takes opposite number to the slope value, that is, hung down
The TAR factors of straight earthquake, as shown in figure 17, the TAR factors calculated using VSP data in the specific embodiment provided for the present invention
Result of calculation figure.
As shown in Figure 1, the method also includes：
S106：Section after dephasing is processed is pulled together along the time curve, the section after being pulled together.
In a particular embodiment, using existing software by the data obtained by step S103 obtained by step S104
Time curve pull together, obtain and prepare to separate the section after the pulling together of down-going wave fields, as shown in figure 18, for the tool that the present invention is provided
In body embodiment along down going wave time curve even up after VSP schematic diagram datas.
S107：Section after pulling together according to the spherical diffusion compensation TAR factor pair of the VSP data is oriented and adds
Power medium filtering.Fig. 4 is the particular flow sheet of step S107, and as shown in Figure 4, the step is specifically included：
S401：Seismic channel in section after described pulling together is recorded to most deep one-level wave detector note from most shallow one-level wave detector
Record sequence；
S402：Select the seismic channel number for participating in computing in the section after described pulling together every time, referred to as medium filtering window.It is described
Medium filtering window progressively reduces at the border of the section after pulling together, to ensure to occur without sky sampling point in medium filtering window
Value so that subsequent treatment result is with respect to fidelity.
S403：Amplitude in the medium filtering window is extracted, discrete data sequences are constituted；
In a particular embodiment, determine medium filtering window (span), i.e., the maximum number of samples of computing per treatment,
The seismic channel number for participating in statistical computation each time is namely selected on VSP data sections.
Sample value is carved while on these seismic channels and constitutes discrete data sequence.The discrete data sequences, it is desirable to from
Most shallow one-level wave detector correspondence seismic channel is to the corresponding seismic channel direction sequencing of most deep one-level wave detector.
In a particular embodiment, it is determined that the maximum number of samples for median filter process is set to n, that is, select and participate in every
The seismic channel number of statistical computation, carves sample value and constitutes discrete data sequence X, such as (1) formula while on these seismic channels
It is shown.
X={ a1,a2,a3,…an} (1)
In formula：X is discrete amplitude sequence；A represents the amplitude of single sampling point；N is maximum number of samples.
From when deep relation, this n sampling point corresponding time series T such as following formula：
T={ t1,t2,t3,…tn} (2)
In formula：T is discrete time series；T represents the time of single sampling point, and unit is for the second (s)；N maximum number of samples.
In a particular embodiment, the seismic channel number span for statistical computation is between 7 and 15
Odd number, the determination of the parameter wants two aspects of comprehensive analysis wave field separation effect and data fidelity, concrete numerical value to regard letter section
State of signal-to-noise and determine, it is higher in signal to noise ratio from 7 as number of samples in signal to noise ratio compares relatively low geological data
In the case of from 15 as number of samples, odd number 9,11 or 13 is selected when signal to noise ratio falls between.The value is pulled together close
Progressively reduce during the border of section afterwards, to ensure to occur without sky sample value in window so that subsequent treatment result is with respect to fidelity.
S404：Weighted number is obtained according to the spherical diffusion compensation TAR factor and the medium filtering window of the VSP data
According to sequence.In a particular embodiment, using the calculated spherical diffusion compensation TAR factor of institute in step S105, and
Sliding window length determined by step S402, designs weighting function, obtains incremental weighted data sequence.
Fig. 5 is the particular flow sheet of step S404.As shown in Figure 5, the step is specifically included：
S501：The TAR factors in any one medium filtering window are determined according to the TAR factors of the VSP data.Using step
Calculated TAR factors k of institute, can obtain the sphere of any one moment t in the section of the geological data in rapid S105
Diffusion compensation, is shown below：
At=A0tk (3)
In formula：AtFor the amplitude at any one moment；A0For the amplitude of zero moment；T represents the moment that sampling point is located, unit
For the second (s)；K is the TAR factors.
S502：The discrete data sequences pair according to the TAR factor pairs in any one medium filtering window in geological data
The sample value in medium filtering window answered carries out relative amplitude relationships compensation.
In a particular embodiment, due to this invention address that the original amplitude information of protection, therefore can not be based on public affairs
Formula (3) carries out the compensation of the entire profile, but can carry out relative amplitude relationships to the sampling point in step S502 window using the formula
Compensation.Concrete grammar is select window thePoint carries out compensation deals on the basis of individual sampling point, as shown in (4) formula.Institute's sampling
Point is in initial data the sampling point for preparing to replace.The step is the committed step of the present invention, and original amplitude of vibration can be made to obtain most
Good protection.
Arrangement is carried out to (4) formula to obtain：
S503：Datum mark normalized is carried out to the sample value after relative amplitude relationships compensation.
In a particular embodiment, withIndividual sampling point amplitude is obtained to carry out datum mark normalized：
In (4) to (6) formula：AiIt is with windowThe theoretical com-pensation amplitude put on the basis of individual sampling point；A′iFor
AiRelative to the value after datum mark normalization；T represents the moment that sampling point is located, and unit is the second (s)；K is the TAR factors；N is in window
Maximum number of samples.
S504：The sample value obtained according to datum mark normalized determines weighted data sequence.
Incremental weighted data sequence Y, can be represented with following formula.
In formula, Y is amplitude weight sequence；T represents the moment that sampling point is located；K is the TAR factors；N is maximum number of samples.
As shown in Figure 4, the step also includes：
S405：Weighting is oriented to described discrete data sequences and weighted data sequence, after obtaining directional weight
Data sequence.
In a particular embodiment, sampling point correspondence is carried out to sequence X and sequence Y to be multiplied, that is, complete directional weight calculating,
New data sequence Z is obtained, as shown in formula (8).
Represent every in (8) formula with A ', be denoted as：
Z={ A '1,A′2,A′3,…A′n} (9)
In (8) formula and (9) formula：Z is amplitude sequence new after weighting；X is original amplitude sequence；Y is amplitude weight
Sequence；A represents the amplitude of single sampling point in original series；T represents the moment that sampling point is located；A ' represents single sampling point after weighting
Amplitude；N is maximum number of samples.
S406：Outlier processing is carried out to the data sequence after directional weight, the data after outlier processing are obtained
Sequence.Fig. 6 is the particular flow sheet of step S406, it will be appreciated from fig. 6 that the step is specifically included：
S601：Determine the mean value of the absolute value of each sampling point in the data sequence after directional weight；
S602：Determined according to described mean value and remove abnormal threshold value.In a particular embodiment, peek row in each
The absolute value of sampling point mean value, and the value plus-minus 20% as a pair is removed into abnormal threshold value.
S603：Absolute value is filtered out in each sampling point of data sequence from after directional weight and removes abnormal threshold value described
Sampling point in addition, as abnormal spots.In a particular embodiment, by the absolute value of each sampling point in former ordered series of numbers with go exception
Threshold value is compared one by one, filter out absolute value it is described remove abnormal threshold value beyond sampling point, as abnormal spots.
S604：The abnormal sample is replaced with the mean value of the absolute value of each sampling point in the data sequence after directional weight
Point, obtains the data sequence after outlier processing.
Namely the present invention removes outlier processing, concrete grammar is the absolute value of each sampling point mean value in peek row, and
The value plus-minus 20% as a pair is removed into abnormal threshold value, then the absolute value of each sampling point in former ordered series of numbers and threshold value are carried out
Compare one by one, and the sample value beyond threshold value scope is replaced with mean value.As shown in formula (10).
Now, the amplitude sequence Z ' that self adaptation is gone after exception is：
Z '=A "1,A″2,A″3,…A″n} (11)
In (10) formula and (11) formula：The amplitude sequence that Z ' is gone after exception for self adaptation；A ' represents single sampling point after weighting
Amplitude；A " represents that self adaptation removes the amplitude of single sampling point after exception；N is maximum number of samples.
The weighting of a sampling point window, self adaptation are illustrated in Figure 19 to Figure 22 and goes abnormal process.Wherein, Figure 19 is this
Sampling point sequence diagram in the specific embodiment of bright offer in a window, Figure 20 in the specific embodiment that provides of the present invention according to
The weight coefficient sequence diagram of VSP data rules, latter window of weighting in the specific embodiment that Figure 21 is provided for the present invention
Interior sampling point sequence diagram, the sampling point sequence that self adaptation is gone after exception in the specific embodiment that Figure 22 is provided for the present invention is illustrated
Figure.
The threshold value is set as the 20% of average, and the value is gone in abnormal method far below other geological data self adaptations
Threshold value.This is that reception environment is relatively quiet, thus abnormal high level or low value occurs because VSP data wave detectors go deep into down-hole
Earthquake outlier situation it is less；And after sampling point weighting is processed, amplitude becomes more balanced in range of normal value.
In other embodiments of the invention, when window slides into boundary position, automatically reduction window is long, to meet calculating need
Will, avoid the abnormality processing effect on border, the technology is known technology, during many business seismic processing softwares have been applied to.
As above preparation is completed, every bit just corresponds to a window in section.
As shown in Figure 4, the step also includes：
S407：Data sequence after to removing outlier processing carries out median filter process, obtains directional weight medium filtering
Processing data.
In a particular embodiment of the present invention, the specific algorithm of the step is according to from big to small by all sampling points in window
Rule rearrangement, and obtain the numerical value of one sampling point in centre of new sequence, to replace section in the point, so as to complete from
Adaptation goes directional weight median filter process abnormal, border lane change number, according to VSP data rules, such as Figure 23 to be this
In the specific embodiment of bright offer directional weight, self adaptation go VSP data after abnormal, border lane change number medium filtering to show
Be intended to 1, Figure 24 in the specific embodiment that provides of the present invention directional weight, self adaptation gone in abnormal, border lane change number
VSP schematic diagram datas 2 after value filtering.
The handling process of circulation step S401 to step S407, can complete multigroup down going wave median filter process.To up
When ripple carries out median filter process, pick up upgoing wave time curve, while the weighting sequence than successively decreasing such as use, the change be by
Seimic wave propagation feature is determined in VSP data methods, and its upgoing wave amplitude reduces with shoaling for wave detector depth is received,
But attenuation law is constant, it is not necessary to recalculate the TAR factors.
As above it is of the invention on the basis of self adaptation removes abnormal, border lane change number medium filtering, proposes a kind of profit
The brand-new median filter method being weighted with VSP data attenuation laws, solve for well control VSP Data Datas process in do not have
There is the problem of reliable wave field separation method.
Fig. 7 is that a kind of medium filtering decayed based on vertical seismic provided in an embodiment of the present invention is decayed based on vertical seismic
The structured flowchart of the system of the medium filtering of rule, as shown in Figure 7, described system includes：
Profile collection device 101, for gathering vertical seismic profiling (VSP) VSP data.
Simple component data sorting unit 102, for sub-electing simple component data according to described VSP data.
In the particular embodiment, VSP acquired original data are arranged using existing method and seismic processing software,
Acquisition is ready for the simple component data of wave field separation, as shown in figure 13, original VSP in the specific embodiment provided for the present invention
Data Z component schematic diagram.
Described data preparation does not include that amplitude compensation is processed, trace equalization is processed, automatic gain (AGC) is processed or other can
The process step of amplitude or energy relativeness can be changed.
Dephasing processing meanss 103, for carrying out dephasing process to described simple component data, obtain zero phase
Section after change process.
In the particular embodiment, dephasing process is carried out using existing software.Dephasing process purpose be by
The complicated wavelet of mixed-phase is converted into the wavelet of zero phase, makes the crest or trough of wavelet corresponding with bed boundary, such as schemes
Shown in 14, seismic wavelet dephasing schematic diagram in the specific embodiment provided for the present invention.
Time curve determining device 104, for obtaining time curve according to the section after dephasing process.
Fig. 8 is the concrete structure block diagram of time curve determining device, and as shown in Figure 8, the device includes：
Take-off location determining module 201, for the crest of descending direct wave in the section after the dephasing is processed
Position is used as take-off location；
Time curve determining module 202, for picking up the descending direct wave with the take-off location, away from song when obtaining
Line.
In a particular embodiment, descending direct wave first arrival, m- depth when obtaining are picked up using existing seismic processing software
Degree relational sequence, the sequence constitutes in the plane the relation curve of a time and distance, hereinafter referred to as time curve, such as schemes
Shown in 15, VSP data first break pickup schematic diagrames in the specific embodiment provided for the present invention.
The first break pickup position is the crest location of descending direct wave in section.
The down going wave time curve is simultaneously as the time curve of the follow-up TAR factors.
As shown in Figure 7, the system also includes：
TAR factors determining device 105, for according to the section and described time curve after dephasing process
Determine the spherical diffusion compensation TAR factor of the VSP data.Fig. 9 is the concrete structure block diagram of TAR factor determining devices, by Fig. 9
Understand, the device is specifically included：
When m- amplitude relation sequence determining module 301, for the dephasing process after section in, pick up institute
The amplitude of time curve is stated, m- amplitude relation sequence when obtaining；
New sequence determining module 302, for it is described when m- amplitude relation sequence in time and amplitude ask for respectively from
So logarithm, obtains new sequence；
Fitting a straight line determining module 303, for carrying out linear fit to described new sequence, obtains fitting a straight line；
Fitting a straight line slope determination module 304, for determining the slope of the fitting a straight line；
TAR factors determining module 305, for taking opposite number to the slope, obtains the TAR factors of the VSP data.
In a particular embodiment, in the section obtained by step S103, step S104 time curve position is picked up
Amplitude, m- amplitude relation sequence during acquisition asks for respectively natural logrithm, then to institute to time in the sequence and amplitude
The new sequence for obtaining carries out linear fit, is calculated the slope of fitting a straight line, and takes opposite number to the slope value, that is, hung down
The TAR factors of straight earthquake, as shown in figure 17, the TAR factors calculated using VSP data in the specific embodiment provided for the present invention
Result of calculation figure.
As shown in Figure 7, the system also includes：
Device 106 is pulled together, is pulled together along the time curve for the section after dephasing is processed, after being pulled together
Section.
In a particular embodiment, using existing software by the data obtained by step S103 obtained by step S104
Time curve pull together, obtain and prepare to separate the section after the pulling together of down-going wave fields, as shown in figure 20, for the tool that the present invention is provided
In body embodiment along down going wave time curve even up after VSP schematic diagram datas.
Medium filtering device 107, it is fixed to carry out for the section after pulling together according to the TAR factor pairs of the VSP data
To Weighted median filtering.Figure 10 is the concrete structure block diagram of medium filtering device, and as shown in Figure 10, the device is specifically included：
Seismic channel order module 401, for seismic channel in the section after described pulling together to be recorded from most shallow one-level wave detector
To most deep one-level wave detector record ordering；
Medium filtering window determining module 402, for selecting described pulling together after section in participate in the seismic channel of computing every time
Number, referred to as medium filtering window.The medium filtering window progressively reduces at the border of the section after pulling together, to ensure intermediate value
Sky sample value is occurred without in spectral window so that subsequent treatment result is with respect to fidelity.
Discrete data sequences determining module 403, for the amplitude in the medium filtering window to be extracted, constitute from
Scattered data sequence；
In a particular embodiment, determine medium filtering window (span), i.e., the maximum number of samples of computing per treatment,
The seismic channel number for participating in statistical computation each time is namely selected on VSP data sections.
Sample value is carved while on these seismic channels and constitutes discrete data sequence.The discrete data sequences, it is desirable to from
Most shallow one-level wave detector correspondence seismic channel is to the corresponding seismic channel direction sequencing of most deep one-level wave detector.
In a particular embodiment, it is determined that the maximum number of samples for median filter process is set to n, that is, select and participate in every
The seismic channel number of statistical computation, carves sample value and constitutes discrete data sequence X, such as (1) formula while on these seismic channels
It is shown.
From when deep relation, this corresponding time series T of n sampling point such as (2) formula.
In a particular embodiment, the seismic channel number span for statistical computation is between 7 and 15
Odd number, the determination of the parameter wants two aspects of comprehensive analysis wave field separation effect and data fidelity, concrete numerical value to regard letter section
State of signal-to-noise and determine, it is higher in signal to noise ratio from 7 as number of samples in signal to noise ratio compares relatively low geological data
In the case of from 15 as number of samples, odd number 9,11 or 13 is selected when signal to noise ratio falls between.The value is pulled together close
Progressively reduce during the border of section afterwards, to ensure to occur without sky sample value in window so that subsequent treatment result is with respect to fidelity.
Weighted data sequence determining module 404, for according to the TAR factors and the medium filtering of the VSP data
Window obtains weighted data sequence.In a particular embodiment, using the calculated TAR factors of institute in step S105, and
Sliding window length determined by step S402, designs weighting function, obtains incremental weighted data sequence.Figure 11 is weighted number
According to the concrete structure block diagram of sequence determining module.As shown in Figure 11, the module is specifically included：
TAR factor specifying units 501, for determining any one medium filtering window according to the TAR factors of the VSP data
The interior TAR factors.Using calculated TAR factors k of institute in step S105, can obtain appointing in the section of the geological data
The spherical diffusion compensation of one moment t of meaning, as shown in (3) formula：
In formula：AtFor the amplitude at any one moment；A0For the amplitude of zero moment；T represents the moment that sampling point is located, unit
For the second (s)；K is the TAR factors.
Relative amplitude relationships compensating unit 502, for according to the TAR in any one medium filtering window in geological data because
Son carries out relative amplitude relationships compensation to the sample value in the corresponding medium filtering window of the discrete data sequences.
In a particular embodiment, due to this invention address that the original amplitude information of protection, therefore can not be based on public affairs
Formula (3) carries out the compensation of the entire profile, but can carry out relative amplitude relationships to the sampling point in step S502 window using the formula
Compensation.Concrete grammar is select window thePoint carries out compensation deals on the basis of individual sampling point, as shown in (4) formula.Institute's sampling
Point is in initial data the sampling point for preparing to replace.The step is the committed step of the present invention, and original amplitude of vibration can be made to obtain most
Good protection.
Arrangement is carried out to (4) formula and obtains (5) formula.
Normalized unit 503, is carried out at datum mark normalization for the sample value after compensating relative amplitude relationships
Reason.
In a particular embodiment, withIndividual sampling point amplitude is obtained to carry out datum mark normalized
(5) formula.
In (4) to (6) formula：AiIt is with windowThe theoretical com-pensation amplitude put on the basis of individual sampling point；A′iFor
AiRelative to the value after datum mark normalization；T represents the moment that sampling point is located, and unit is the second (s)；K is the TAR factors；N is in window
Maximum number of samples.
Weighted data sequence determination unit 504, the sample value for being obtained according to datum mark normalized determines weighting
Data sequence.
Incremental weighted data sequence Y, can be represented with (7) formula.
As shown in Figure 10, the module also includes：
Directional weight module 405, for being oriented weighting to described discrete data sequences and weighted data sequence,
Obtain the data sequence after directional weight.
In a particular embodiment, sampling point correspondence is carried out to sequence X and sequence Y to be multiplied, that is, complete directional weight calculating,
New data sequence Z is obtained, as shown in formula (8).
Represent every in (8) formula with A ', be denoted as (9).
In (8) formula and (9) formula：Z is amplitude sequence new after weighting；X is original amplitude sequence；Y is amplitude weight
Sequence；A represents the amplitude of single sampling point in original series；T represents the moment that sampling point is located；A ' represents single sampling point after weighting
Amplitude；N is maximum number of samples.
Outlier processing module 406 is gone, for carrying out outlier processing to the data sequence after directional weight, is gone
Data sequence after outlier processing.Figure 12 is the concrete structure block diagram of outlier processing module, as shown in Figure 12, the module
Specifically include：
Mean value determining unit 601, for determining directional weight after data sequence in each sampling point absolute value it is flat
Average；
Abnormal threshold value determining unit 602 is gone, for determining according to described mean value abnormal threshold value is removed.Specific
In embodiment, the absolute value of each sampling point mean value in peek row, and the value plus-minus 20% as a pair is removed into abnormal threshold
Value.
Abnormal spots determining unit 603, for filtering out in each sampling point of the data sequence from after directional weight definitely
Value it is described remove abnormal threshold value beyond sampling point, as abnormal spots.In a particular embodiment, by former ordered series of numbers each
The absolute value of sampling point is compared one by one with abnormal threshold value is removed, filter out absolute value it is described remove abnormal threshold value beyond sample
Point, as abnormal spots.
Abnormal spots replacement unit 604, for in the data sequence after directional weight the absolute value of each sampling point it is flat
Average replaces the abnormal spots, obtains the data sequence after outlier processing.
Namely the present invention removes outlier processing, concrete grammar is the absolute value of each sampling point mean value in peek row, and
The value plus-minus 20% as a pair is removed into abnormal threshold value, then the absolute value of each sampling point in former ordered series of numbers and threshold value are carried out
Compare one by one, and the sample value beyond threshold value scope is replaced with mean value.As shown in formula (10).
Now, the amplitude sequence Z ' that self adaptation is gone after exception is shown in (11).
In (10) formula and (11) formula：The amplitude sequence that Z ' is gone after exception for self adaptation；A ' represents single sampling point after weighting
Amplitude；A " represents that self adaptation removes the amplitude of single sampling point after exception；N is maximum number of samples.
The weighting of a sampling point window, self adaptation are illustrated in Figure 19 to Figure 22 and goes abnormal process.Wherein, Figure 19 is this
Sampling point sequence diagram in the specific embodiment of bright offer in a window, Figure 20 in the specific embodiment that provides of the present invention according to
The weight coefficient sequence diagram of VSP data rules, latter window of weighting in the specific embodiment that Figure 21 is provided for the present invention
Interior sampling point sequence diagram, the sampling point sequence that self adaptation is gone after exception in the specific embodiment that Figure 22 is provided for the present invention is illustrated
Figure.
The threshold value is set as the 20% of average, and the value is gone in abnormal method far below other geological data self adaptations
Threshold value.This is that reception environment is relatively quiet, thus abnormal high level or low value occurs because VSP data wave detectors go deep into down-hole
Earthquake outlier situation it is less；And after sampling point weighting is processed, amplitude becomes more balanced in range of normal value.
In other embodiments of the invention, when window slides into boundary position, automatically reduction window is long, to meet calculating need
Will, avoid the abnormality processing effect on border, the technology is known technology, during many business seismic processing softwares have been applied to.
As above preparation is completed, every bit just corresponds to a window in section.
As shown in Figure 10, the module also includes：
Medium filtering module 407, for the data sequence after to removing outlier processing median filter process is carried out, and is determined
To Weighted median filtering processing data.
In a particular embodiment of the present invention, the specific algorithm of the step is according to from big to small by all sampling points in window
Rule rearrangement, and obtain the numerical value of one sampling point in centre of new sequence, to replace section in the point, so as to complete from
Adaptation goes directional weight median filter process abnormal, border lane change number, according to VSP data rules, such as Figure 23 to be this
In the specific embodiment of bright offer directional weight, self adaptation go VSP data after abnormal, border lane change number medium filtering to show
Be intended to 1, Figure 24 in the specific embodiment that provides of the present invention directional weight, self adaptation gone in abnormal, border lane change number
VSP schematic diagram datas 2 after value filtering.
The handling process of circulation step S401 to step S407, can complete multigroup down going wave median filter process.To up
When ripple carries out median filter process, upgoing wave time curve is picked up；The weighting sequence than successively decreasing such as use simultaneously, the change be by
Seimic wave propagation feature is determined in VSP data methods, and its upgoing wave amplitude reduces with shoaling for wave detector depth is received,
But attenuation law is constant, it is not necessary to recalculate the TAR factors.
As above it is of the invention on the basis of self adaptation removes abnormal, border lane change number medium filtering, proposes a kind of profit
The brand-new medium filtering system being weighted with VSP data attenuation laws, solve for well control VSP Data Datas process in do not have
There is the problem of reliable wave field separation method.
With reference to specific embodiment, technical scheme is discussed in detail.The purpose of the present invention is to propose to
A kind of median filter method being weighted for VSP data attenuation laws, solves not have for the VSP Data processings of well control
The problem of reliable wave field separation method.It is on the basis of self adaptation removes abnormal, border lane change number medium filtering, using VSP
The attenuation law of data is weighted the brand-new median filter method of process.
In a particular embodiment, the program is specifically included：
1) data prepare：
VSP data acquired original data are arranged using existing method and seismic processing software, is ready for
The simple component data of wave field separation, as shown in figure 13, original VSP data Z component shows in the specific embodiment provided for the present invention
It is intended to.In order to show VSP data amplitudes relations, the data do not add gain to show (similarly hereinafter).
Described data preparation does not include that amplitude compensation is processed, trace equalization is processed, automatic gain (AGC) is processed or other can
The process step of amplitude or energy relativeness can be changed.
2) dephasing is processed：
Dephasing process is carried out using existing software.The purpose of dephasing process is by the complicated wavelet of mixed-phase
The wavelet of zero phase is converted into, makes the crest or trough of wavelet corresponding with bed boundary, be offer of the present invention as shown in figure 14
Specific embodiment in seismic wavelet dephasing schematic diagram.
3) m- depth relationship pickup when：
Descending direct wave first arrival is picked up using existing seismic processing software, m- depth relationship sequence when obtaining, the sequence exists
The relation curve of one time of composition and distance in plane, hereinafter referred to as time curve, are offer of the present invention as shown in figure 15
Specific embodiment in VSP data first break pickup schematic diagrames.
The first break pickup position is the crest location of descending direct wave in section.
4) the spherical diffusion compensation factor (the TAR factors) is calculated：
In step 2) obtained by section in, pick up step 3) amplitude of time curve position, m- amplitude is closed during acquisition
It is sequence, natural logrithm is asked for respectively to time in the sequence and amplitude, then linear fit is carried out to resulting new sequence,
The slope of fitting a straight line is calculated, and opposite number is taken to the slope value, that is, obtain the TAR factors of vertical seismic, such as Figure 18 is
The result of calculation figure of the TAR factors calculated using VSP data in the specific embodiment that the present invention is provided, Figure 17 is provided for the present invention
Specific embodiment in TAR calculating difference curve figure.For monitoring the precision of TAR calculating.
5) time curve is pulled together：
Pulled together using time curve of the existing software by the data obtained by 2) step obtained by 3) step, prepared
The geological data of down-going wave fields is separated, as shown in figure 18, along down going wave time curve in the specific embodiment provided for the present invention
VSP schematic diagram datas after evening up.
6) median filter process is carried out in the following manner：
(1) determine maximum number of samples n for median filter process, that is, select the seismic channel for participating in statistical computation each time
Number, carves sample value and constitutes discrete data sequence X, as shown in (1) formula while on these seismic channels.
From step 3) when deep relation, this corresponding time series T of n sampling point is as shown in (2) formula.
The seismic channel number span for statistical computation is the odd number between 7 and 15, the determination of the parameter
Two aspects of comprehensive analysis wave field separation effect and data fidelity, concrete numerical value are wanted to determine depending on believing section state of signal-to-noise,
Signal to noise ratio compares in relatively low geological data from 7 as number of samples, in good signal to noise situations from 15 as sample
Points, select odd number 9,11 or 13 when signal to noise ratio falls between.
The discrete data sequences, it is desirable to corresponding to most deep one-level wave detector from most shallow one-level wave detector correspondence seismic channel
Seismic channel direction sequencing.
(2) using step 4) middle calculated TAR factor k, the sphere of any one moment t in section can be obtained
Diffusion compensation, as shown in (3) formula.
Due to this invention address that the original amplitude information of protection, therefore the entire profile can not be carried out based on formula (3)
Compensation, but can utilize the formula that the compensation of relative amplitude relationships is carried out to the sampling point in (1) step window.Concrete grammar is to select window
Point carries out compensation deals on the basis of individual sampling point, as shown in (4) formula.Institute's sampling point is in initial data and prepares to replace
The sampling point for changing.The step is the committed step of the present invention, and original amplitude of vibration can be made to obtain best protection.
Arrangement is carried out to (4) formula and obtains (5) formula, then withIndividual sampling point amplitude is to carry out at datum mark normalization
Reason, obtains (6) formula.In (4) to (6) formula：AiIt is with windowThe theoretical com-pensation amplitude put on the basis of individual sampling point；
A′iFor AiRelative to the value after datum mark normalization；T represents the moment that sampling point is located, and unit is the second (s)；K is the TAR factors；N is
Maximum number of samples in window.
In weighted data sequence Y being incremented by, as shown in (7) formula.
(3) sampling point correspondence is carried out to sequence X in (1) step and sequence Y in (2) step to be multiplied, that is, complete directional weight meter
Calculate, new data sequence Z is obtained, as shown in formula (8).Represent every in (8) formula with A ', be denoted as (9) formula.
In (8) formula and (9) formula：Z is amplitude sequence new after weighting；X is original amplitude sequence；Y is amplitude weight
Sequence；A represents the amplitude of single sampling point in original series；T represents the moment that sampling point is located；A ' represents single sampling point after weighting
Amplitude；N is maximum number of samples.
(4) outlier processing is carried out to the sequence obtained by previous step, concrete grammar is that each sampling point is put down in peek row
The absolute value of average, and using value plus-minus 20% as removing abnormal threshold value for a pair, then by former ordered series of numbers each sampling point it is exhausted
Value is compared one by one with threshold value, and the sample value beyond threshold value scope is replaced with mean value.As shown in formula (10).
Now, the amplitude sequence Z ' that self adaptation is gone after exception is as shown in (11) formula.
In (10) formula and (11) formula：The amplitude sequence that Z ' is gone after exception for self adaptation；A ' represents single sampling point after weighting
Amplitude；A " represents that self adaptation removes the amplitude of single sampling point after exception；N is maximum number of samples.
The weighting of a sampling point window, self adaptation are illustrated in Figure 19 to Figure 22 and goes abnormal process.Figure 19 is the present invention
Sampling point sequence diagram in the specific embodiment of offer in a window, Figure 20 in the specific embodiment that provides of the present invention according to
The weight coefficient sequence diagram of VSP data rules, latter window of weighting in the specific embodiment that Figure 21 is provided for the present invention
Interior sampling point sequence diagram, the sampling point sequence that self adaptation is gone after exception in the specific embodiment that Figure 22 is provided for the present invention is illustrated
Figure.
The threshold value is set as the 20% of average, and the value is gone in abnormal method far below other geological data self adaptations
Threshold value.This is that reception environment is relatively quiet, thus abnormal high level or low value occurs because VSP data wave detectors go deep into down-hole
Earthquake outlier situation it is less；And after sampling point weighting is processed, amplitude becomes more balanced in range of normal value.
(5) when window slides into boundary position, automatically reduction window is long, to meet calculating needs, avoids the exception on border
Reason effect, the technology is known technology, during many business seismic processing softwares have been applied to.
(6) preparation of step (1) to step (5) is completed, every bit just corresponds to a window in section, final step is
Medium filtering calculating is carried out, specific algorithm is that all sampling points in window are resequenced according to rule from big to small, and obtains new
The numerical value of one sampling point in centre of sequence, to replace section in the point, go abnormal, border lane change so as to complete self adaptation
Directional weight median filter process several, according to VSP data rules, orients in the specific embodiment that Figure 23 is provided for the present invention
Weight, self adaptation goes VSP schematic diagram datas 1 after abnormal, border lane change number medium filtering, Figure 24 to provide for the present invention
In specific embodiment directional weight, self adaptation remove VSP schematic diagram datas 2 after abnormal, border lane change number medium filtering.
7) separation of multigroup wave field is completed：
Circulation step 3) to step 6) handling process, multigroup down going wave median filter process can be completed.If to up
Ripple carries out median filter process, then need in step 3) in, pick up upgoing wave time curve；Simultaneously in step 6) used in wait ratio
The weighting sequence for successively decreasing, the change determined by seimic wave propagation feature in VSP data methods, and its upgoing wave amplitude is with connecing
Shoaling and reducing for wave detector depth is received, but attenuation law is constant, it is not necessary to recalculate the TAR factors.
In sum, on the basis of self adaptation removes abnormal, border lane change number medium filtering, the invention provides a kind of
Median filter method and system based on vertical seismic attenuation law, it is a kind of accurately based in vertical seismic attenuation law
The scheme of value filtering, by VSP data obtain simple component data, dephasing process is carried out to simple component data, with reference to when away from song
Line determines the TAR factors, and pulls together section along time curve, is oriented Weighted median filtering, greatly improves medium filtering
Reliability, to provide reliable wave field separation scheme for the process of the VSP Data Datas of well control.
One of ordinary skill in the art will appreciate that realizing all or part of flow process in above-described embodiment method, Ke Yitong
Computer program is crossed to instruct the hardware of correlation to complete, described program can be stored in general computer read/write memory medium
In, the program is upon execution, it may include such as the flow process of the embodiment of above-mentioned each method.Wherein, described storage medium can be magnetic
Dish, CD, read-only memory (Read-Only Memory, ROM) or random access memory (Random Access
Memory, RAM) etc..
Those skilled in the art will also be appreciated that the various functions that the embodiment of the present invention is listed are by hardware or soft
Part is realizing depending on the design requirement of specific application and whole system.Those skilled in the art can be specific for every kind of
Using, it is possible to use various methods realize described function, but this realization is understood not to be protected beyond the embodiment of the present invention
The scope of shield.
Apply specific embodiment in the present invention to be set forth the principle and embodiment of the present invention, above example
Explanation be only intended to help and understand the method for the present invention and its core concept；Simultaneously for one of ordinary skill in the art,
According to the thought of the present invention, will change in specific embodiments and applications, in sum, in this specification
Appearance should not be construed as limiting the invention.
Claims (12)
1. a kind of median filter method based on vertical seismic attenuation law, is characterized in that, described method includes：
Collection vertical seismic profiling (VSP) VSP data；
Simple component data are sub-elected according to described VSP data；
Dephasing process is carried out to described simple component data, the section after dephasing is processed is obtained；
Time curve is obtained according to the section after dephasing process；
Determine that the spherical diffusion of the VSP data is mended according to the section and described time curve after dephasing process
Repay the TAR factors；
Section after dephasing is processed is pulled together along the time curve, the section after being pulled together；
Section after pulling together according to the spherical diffusion compensation TAR factor pair of the VSP data is oriented weighted median filter
Ripple.
2. method according to claim 1, is characterized in that, away from song when being obtained according to the section after dephasing process
Line includes：
Using the crest location of descending direct wave in the section after dephasing process as take-off location；
The descending direct wave is picked up with the take-off location, time curve is obtained.
3. method according to claim 2, is characterized in that, according to the section and described after dephasing process
Time curve determines that the spherical diffusion compensation TAR factor of the VSP data includes：
In section after dephasing process, the amplitude of the time curve, m- amplitude relation sequence when obtaining are picked up
Row；
To it is described when m- amplitude relation sequence in time and amplitude ask for natural logrithm respectively, obtain new sequence；
Linear fit is carried out to described new sequence, fitting a straight line is obtained；
Determine the slope of the fitting a straight line；
Opposite number is taken to the slope, the spherical diffusion compensation TAR factor of the VSP data is obtained.
4. method according to claim 3, is characterized in that, according to the spherical diffusion compensation TAR factor pair of the VSP data
It is described pull together after section be oriented Weighted median filtering and process and include：
Seismic channel in section after described pulling together is recorded to most deep one-level wave detector record ordering from most shallow one-level wave detector；
Select the seismic channel number for participating in computing in the section after described pulling together every time, referred to as medium filtering window；
Amplitude in the medium filtering window is extracted, discrete data sequences are constituted；
Weighted data sequence is obtained according to the spherical diffusion compensation TAR factor and the medium filtering window of the VSP data；
Weighting is oriented to described discrete data sequences and weighted data sequence, the data sequence after directional weight is obtained
Row；
Outlier processing is carried out to the data sequence after directional weight, the data sequence after outlier processing is obtained；
Data sequence after to removing outlier processing carries out median filter process, obtains directional weight median filter process data.
5. method according to claim 4, is characterized in that, according to the spherical diffusion compensation TAR factor of the VSP data with
And the medium filtering window obtains weighted data sequence and includes：
Determine that the spherical diffusion in any one medium filtering window is mended according to the spherical diffusion compensation TAR factor of the VSP data
Repay the TAR factors；
The corresponding intermediate value of discrete data sequences according to the spherical diffusion compensation TAR factor pair in any one medium filtering window
Sample value in spectral window carries out relative amplitude relationships compensation；
Datum mark normalized is carried out to the sample value after relative amplitude relationships compensation；
The sample value obtained according to datum mark normalized determines weighted data sequence.
6. method according to claim 4, is characterized in that, the data sequence after directional weight is carried out at exceptional value
Reason, obtaining the data sequence after outlier processing includes：
Determine the mean value of the absolute value of each sampling point in the data sequence after directional weight；
Determined according to described mean value and remove abnormal threshold value；
Filter out in each sampling point of data sequence from after directional weight absolute value it is described remove abnormal threshold value beyond sample
Point, as abnormal spots；
The abnormal spots are replaced with the mean value of the absolute value of each sampling point in the data sequence after directional weight, exception is obtained
Data sequence after value process.
7. a kind of medium filtering system based on vertical seismic attenuation law, is characterized in that, described system includes：
Profile collection device, for gathering vertical seismic profiling (VSP) VSP data；
Simple component data sorting unit, for sub-electing simple component data according to described VSP data；
Dephasing processing meanss, for carrying out dephasing process to described simple component data, obtain dephasing process
Section afterwards；
Time curve determining device, for obtaining time curve according to the section after dephasing process；
Spherical diffusion compensation TAR factor determining device, for according to the section after dephasing process and it is described when
Determine the spherical diffusion compensation TAR factor of the VSP data away from curve；
Device is pulled together, is pulled together along the time curve for the section after dephasing is processed, the section after being pulled together；
Medium filtering device, for the section after pulling together according to the spherical diffusion compensation TAR factor pair of the VSP data
It is oriented Weighted median filtering.
8. system according to claim 7, is characterized in that, described time curve determining device includes：
Take-off location determining module, for using the dephasing process after section in descending direct wave crest location as
Take-off location；
Time curve determining module, for picking up the descending direct wave with the take-off location, obtains time curve.
9. system according to claim 8, is characterized in that, described spherical diffusion compensation TAR factor determining device includes：
When m- amplitude relation sequence determining module, for the dephasing process after section in, pick up it is described when away from song
The amplitude of line, m- amplitude relation sequence when obtaining；
New sequence determining module, for it is described when m- amplitude relation sequence in time ask for natural logrithm respectively with amplitude,
Obtain new sequence；
Fitting a straight line determining module, for carrying out linear fit to described new sequence, obtains fitting a straight line；
Fitting a straight line slope determination module, for calculating the slope of the fitting a straight line；
Spherical diffusion compensation TAR factor determining module, for taking opposite number to the slope, obtains the sphere of the VSP data
The diffusion compensation TAR factor.
10. system according to claim 9, is characterized in that, described medium filtering device includes：
Seismic channel order module, for seismic channel in the section after described pulling together to be recorded to most deep one from most shallow one-level wave detector
Level wave detector record ordering；
Medium filtering window determining module, for selecting described pulling together after section in participate in the seismic channel number of computing every time, referred to as
Medium filtering window；
Discrete data sequences determining module, for the amplitude in the medium filtering window to be extracted, constitutes discrete data
Sequence；
Weighted data sequence determining module, for according to the spherical diffusion compensation TAR factor and the intermediate value of the VSP data
Spectral window obtains weighted data sequence；
Directional weight module, for being oriented weighting to described discrete data sequences and weighted data sequence, is determined
To the data sequence after weighting；
Outlier processing module is gone, for carrying out outlier processing to the data sequence after directional weight, exceptional value is obtained
Data sequence after process；
Medium filtering module, for the data sequence after to removing outlier processing median filter process is carried out, and obtains directional weight
Median filter process data.
11. systems according to claim 10, is characterized in that, described weighted data sequence determining module includes：
Spherical diffusion compensation TAR factor specifying unit, for being determined according to the spherical diffusion compensation TAR factor of the VSP data
The spherical diffusion compensation TAR factor in any one medium filtering window；
Relative amplitude relationships compensating unit, for according to the spherical diffusion compensation TAR factor pair in any one medium filtering window
Sample value in the corresponding medium filtering window of the discrete data sequences carries out relative amplitude relationships compensation；
Normalized unit, for the sample value after compensating relative amplitude relationships datum mark normalized is carried out；
Weighted data sequence determination unit, the sample value for being obtained according to datum mark normalized determines weighted data sequence
Row.
12. systems according to claim 10, is characterized in that, described outlier processing module of going includes：
Mean value determining unit, for determining directional weight after data sequence in each sampling point absolute value mean value；
Abnormal threshold value determining unit is gone, for determining according to described mean value abnormal threshold value is removed；
Abnormal spots determining unit, for filtering out absolute value described in each sampling point of the data sequence from after directional weight
The sampling point gone beyond abnormal threshold value, as abnormal spots；
Abnormal spots replacement unit, for being replaced with the mean value of the absolute value of each sampling point in the data sequence after directional weight
The abnormal spots, obtain the data sequence after outlier processing.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201410718012.8A CN104656140B (en) | 2014-12-01 | 2014-12-01 | Median filtering method and system based on vertical earthquake attenuation laws |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201410718012.8A CN104656140B (en) | 2014-12-01 | 2014-12-01 | Median filtering method and system based on vertical earthquake attenuation laws |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN104656140A CN104656140A (en) | 2015-05-27 |
| CN104656140B true CN104656140B (en) | 2017-05-10 |
Family
ID=53247487
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201410718012.8A Active CN104656140B (en) | 2014-12-01 | 2014-12-01 | Median filtering method and system based on vertical earthquake attenuation laws |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN104656140B (en) |
Families Citing this family (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN105158801B (en) * | 2015-07-29 | 2017-10-17 | 中国石油天然气集团公司 | The compression process method and apparatus of optical cable coupled noise |
| CN107315196B (en) * | 2016-04-26 | 2020-08-25 | 中国石油化工股份有限公司 | Improved VSP wave field separation method and device based on median filtering |
| CN107255836B (en) * | 2017-07-13 | 2019-02-05 | 中国石油天然气集团公司 | Multi-component earthquake data filtering method based on vector mixing distance-taxis |
| CN108035708B (en) * | 2017-11-20 | 2021-04-30 | 中国石油天然气股份有限公司 | Method and device for removing stratum interface reflected waves |
| CN112346129A (en) * | 2020-10-29 | 2021-02-09 | 中国石油天然气集团有限公司 | Logging curve extraction and synthetic seismic record manufacturing method and device |
Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US6292754B1 (en) * | 1999-11-11 | 2001-09-18 | Bp Corporation North America Inc. | Vector recomposition of seismic 3-D converted-wave data |
| CN102323617A (en) * | 2011-06-13 | 2012-01-18 | 中国石油化工股份有限公司 | Merging processing method of 2D seismic data of complex surfaces |
| CN104090299A (en) * | 2014-07-16 | 2014-10-08 | 中国石油集团川庆钻探工程有限公司地球物理勘探公司 | Surface seismic data amplitude compensation method based on VSP primary waves |
Family Cites Families (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| FR2831961B1 (en) * | 2001-11-07 | 2004-07-23 | Inst Francais Du Petrole | METHOD FOR PROCESSING SEISMIC DATA OF WELLS IN ABSOLUTE PRESERVED AMPLITUDE |
-
2014
- 2014-12-01 CN CN201410718012.8A patent/CN104656140B/en active Active
Patent Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US6292754B1 (en) * | 1999-11-11 | 2001-09-18 | Bp Corporation North America Inc. | Vector recomposition of seismic 3-D converted-wave data |
| CN102323617A (en) * | 2011-06-13 | 2012-01-18 | 中国石油化工股份有限公司 | Merging processing method of 2D seismic data of complex surfaces |
| CN104090299A (en) * | 2014-07-16 | 2014-10-08 | 中国石油集团川庆钻探工程有限公司地球物理勘探公司 | Surface seismic data amplitude compensation method based on VSP primary waves |
Non-Patent Citations (4)
| Title |
|---|
| Sub-salt imaging on VSP multi-component data:A case study in DB1 well-site,Tarim Basin,China;Yan Yousheng et al.;《Expanded Abstracts CPS/SEG 2004 International Geophysical Conference （Volume Ⅱ）》;20041231;第626-630页 * |
| 一种改进的中值滤波方法;袁勇等;《成都理工大学学报（自然科学版）》;20130430;第40卷(第2期);第125-129页 * |
| 井控处理中的真振幅恢复与Q补偿方法及应用;王正和等;《物探与化探》;20080831;第32卷(第4期);第434-437页 * |
| 面向油田开发的井控地震资料处理技术及应用;康有元等;《石油物探》;20130731;第52卷(第4期);第402-408页 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN104656140A (en) | 2015-05-27 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN104656140B (en) | Median filtering method and system based on vertical earthquake attenuation laws | |
| WO2017024523A1 (en) | Inversion method for ray elastic parameter | |
| CN107783185B (en) | A kind of processing method and processing device of tomographic statics | |
| CN105607124B (en) | Seismic wave near surface interval quality factors compensation method and device | |
| CN104570076B (en) | Automatic seismic wave first-arrival picking method based on dichotomy | |
| CN104280777B (en) | Method for suppressing interference of seismic data multiples on land | |
| CN105863628B (en) | A kind of phase of development subtle hydrocarbon reservoir method | |
| CN109085663A (en) | A kind of tight sandstone reservoir stratification seam recognition methods | |
| CN109254324A (en) | Full range protects width seismic data processing technique and device | |
| CN109884707A (en) | Near surface is layered time-depth curve static correcting method | |
| CN111381275A (en) | First arrival picking method and device for seismic data | |
| CN105629300B (en) | The method for improving complicated structure offset data signal-to-noise ratio | |
| CN104570116A (en) | Geological marker bed-based time difference analyzing and correcting method | |
| CN103744114B (en) | Method and device for estimating quality factor on basis of zero offset VSP (vertical seismic profile) data | |
| CN112253087A (en) | Biological disturbance reservoir physical property calculation method based on multi-source logging data | |
| CN104330821B (en) | Poststack signal-to-noise ration estimation method and device | |
| CN112183407B (en) | Tunnel seismic wave data denoising method and system based on time-frequency domain spectral subtraction | |
| CN107065007B (en) | A kind of seismic data amplitude method of adjustment and device | |
| Button et al. | Distribution of mafic sills in the Transvaal Supergroup, northeastern South Africa | |
| CN107942389A (en) | For suppressing method, system and the computer-readable medium of adjacent big gun interference | |
| CN109298451B (en) | Method for automatically picking S wave seismic phase by improving skewness | |
| CN107643541B (en) | Normal-moveout spectrum means of interpretation based on rate pattern | |
| CN104199107B (en) | Depth prediction approach and system before brill based on the many wave datum of vertical seismic | |
| CN106772596A (en) | A kind of method and device for determining pre-stack time migration velocity field | |
| CN108535778A (en) | A kind of seismic first breaks detection method and system |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| C06 | Publication | ||
| PB01 | Publication | ||
| C10 | Entry into substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |