CN106908840A - Seismic data Hz noise automatic identification and drawing method based on principal component analysis - Google Patents

Abstract

The present invention relates to a kind of seismic data Hz noise automatic identification based on principal component analysis and drawing method, by calculating each road the average energy value of earthquake record, according to the average energy value otherness automatic identification Hz noise road, extract some neighboring track composition seismic interference road collection comprising interference way, according to the strong correlation of Hz noise noise, characteristic vector is reconstructed using principal component analysis technology reach the purpose of compacting Hz noise.The energy of useful signal is weakened while denoising due to introducing principal component analysis, is that this has carried out relative amplitude preserved processing again to reconstruct data, it is ensured that the trace consistency of earthquake record.Empirical tests, a kind of seismic data Hz noise automatic identification based on principal component analysis disclosed by the invention can realize the automatic identification to Hz noise in seismic data and compacting with drawing method, desert area seismic data quality is effectively increased, it is significant especially for the bad track for repairing the original seismic data caused by Hz noise.

Description

Seismic data Hz noise automatic identification and drawing method based on principal component analysis

Technical field：

The present invention relates to a kind of method of seismic data process denoising, the especially serious earthquake record denoising of Hz noise Method, using the average energy value identification interference way in each road of earthquake record and using principal component analysis technology compacting Hz noise Method.

Background technology：

The influence of Hz noise is suffered from during acquiring seismic exploration data, is often mixed with the earthquake record for collecting Industrial frequency noise, when Hz noise is very strong, and through one geological data all the time when, the seismic channel often shows as bad track, draws Seismic data quality is played to decline.By the analysis to a large amount of bad track data, many seismic channels do not stem from the strong of Near Ground Interference source, often from geophone problem in itself, the time depth with earthquake record is unrelated, currently to the treatment of bad track It is to reject bad track, this treatment reduces earthquake data quality, and the technology that this patent is provided, can effectively recognize due to power frequency The bad track for producing is disturbed, this provides possibility further to suppress Hz noise.Do not have very much also currently for Hz noise The solution of effect.Drawing method main at present has wave trap method, power frequency recurrence to subtract each other suppression method, the class of self-adaptive routing three. Article " automatic identification of Hz noise and compacting in seismic data " proposes that wave trap method filters industry disturbance, but the method is simultaneously The seismic signal of 50Hz frequencies can be filtered；Article " extraction of disturbance of industry frequency ripple and removing method " is proposed in time-domain with altogether Yoke gradient algorithm extracts mono-tone interference, but the party have ignored the changeability of work frequency in actual seismic data, and research shows Different zones, different time, Hz noise frequency is change；Article " force in compacting seismic data by the cosine of industrial noise The improvement and application of nearly method " proposes to filter industry disturbance based on cosine approach, but the method cannot realize that automatic identification is disturbed Road；Article " separating 50Hz industry disturbances in industry disturbance in geophone domain " proposes that the conversion of big gun collection filters industry disturbance, but the method efficiency is low simultaneously Cannot automatic identification interference way；Article " 50Hz industry disturbances denoising method and application based on Wiener filtering " proposes Wiener filtering Single-frequency noise is predicted in matching of the device to adjust the amplitude and phase of library track, but the method is not suitable for the power frequency of stationarity difference Interference；Article " the Hz noise technology for eliminating based on independent component analysis " proposition is isolated separate from mixed signal Each component of signal reaches the purpose for eliminating man-made noise, but the method lacks strict theoretical foundation by empirical；

" a kind of method of automatic identification based on dual factors with industry disturbance is removed " is proposed disclosed in CN104570118A Sine and cosine Weighted approximation method processes industry disturbance, but the method have ignored the frequency of power frequency, the not stationarity of phase and amplitude；

" method of alternating current disturbance signal in removal geological data " proposes that frequency domain is extracted disclosed in CN103630935A The method of work frequency, but the method does not account for the presence of industrial frequency harmonic；

" a kind of automatic identification and the method for eliminating seismic exploration industry electrical interference " proposes base disclosed in CN101907726A Determine Hz noise in self-correlation theory, but seismic data before first arrival time is depended on during the method automatic identification industrial frequency noise Quality.

The content of the invention：

The purpose of the present invention is that for above-mentioned the deficiencies in the prior art, there is provided a kind of earthquake based on principal component analysis Data Hz noise automatic identification and drawing method.

Main idea is that：Seismic data is often influenceed by Hz noise during earthquake-capturing, from And a large amount of bad tracks or strong jamming road in earthquake record are caused, and the quality of seismic exploration data is had a strong impact on, the present invention is to pass through Then automatic identification is isolated interference by the seismic channel of Hz noise by principal component analysis technology, then to the data after reconstruct Guarantor's width is carried out, automatic identification and the compacting of seismic data Hz noise is realized.

The present invention is achieved by the following technical solutions：

Seismic data Hz noise automatic identification and drawing method based on principal component analysis, comprise the following steps：

A, for pending original seismic data, choose the when window T that the record time is located at middle hypomere_W, calculate the earthquake Record is in T_WNei Ge roads the average energy value, such as formula

Wherein E_iFor when window in the i-th road of earthquake record the average energy value, x_i(t_j) for when window in the road of earthquake record i-th jth Individual sampling point value, i=1,2 ... n, n be the total road number of earthquake record, j=1,2 ... k, when k is in window earthquake record sampling number, When window T_WLength suggestion for summary journal lengthArriveIn time；

B, definition l roads characterization factorAccording to λ_lValue identification interference way, E_lFor earthquake record l roads exist T_WInterior the average energy value, q is T_WThe interior earthquake Taoist monastic name with minimum average B configuration energy value, wherein l=1,2 ... n recognize interference way Comprise the following steps that：

If b1, there is λ to all roads_l＜ 5, shows the earthquake record without notable Hz noise, without follow-up denoising Journey, otherwise performs step b2；

If b2, λ_l>=5, show that earthquake record has notable Hz noise, note l roads are interference way, are designated as U_l；

C, extraction and U_lIt is adjacent and including U_lContinuous 5 trace record constitute sub- seismic interference road collection, X is designated as, per together in X It is a row, obtains the average of X each columns, and the average of the row is subtracted with each element in X, obtains the sub- seismic channel after X treatment Collection, is designated as X^*；

D, calculating X^*Covariance matrix Γ, such as formula

Wherein b is X^*Columns, X^*TIt is X^*Transposed matrix, " " representing matrix multiplication；

E, eigenvalue matrix Λ and eigenvectors matrix R that covariance matrix Γ is calculated using singular value decomposition method, then There is formula

Γ=R Λ R^T (3)

Wherein Λ is characterized the diagonal matrix that value is arranged from big to small, each characteristic vectors for being classified as character pair value of R, R^T It is transposed matrix, meets R^TR=RR^T=E, E are unit matrix；

F, X are through R^TLinear Mapping, obtains principal component matrix, such as formula

Φ=R^T·X (4)

Wherein Φ is main component matrix；

G, Φ the first rows are set to 0, obtain Φ ', made

X₁ ^*=R Φ ' (5)

Wherein X₁ ^*It is restructuring matrix；

H, extraction X₁ ^*The data of middle correspondence interference way, are designated as U_l', to U_l' carry out relative amplitude preserved processing, such as formula

Wherein U_l ^*It is interference way signal, A after relative amplitude preserved processing₁It is the average value of interference way signal amplitude absolute value after denoising, A₂It is and U_lThe absolute value of amplitude average value sum of adjacent two track data it is average.

Beneficial effect：Through experiment, a kind of seismic data Hz noise based on principal component analysis disclosed by the invention from Dynamic identification can be realized suppressing the noise that industrial frequency interference source brings in seismic prospecting with drawing method, and being capable of automatic identification Interference way, the algorithm can repair the bad track that Hz noise causes, and effectively increase desert area seismic data quality, reduce and gather into This.

Brief description of the drawings：

Fig. 1 real seismic records

Earthquake record after Fig. 2 compacting Hz noises

Specific embodiment：

The present invention is described in further detail with reference to the accompanying drawings and examples：

Seismic data Hz noise automatic identification and drawing method based on principal component analysis, comprise the following steps：

A, using actual single shot record, this example T_W=700ms~800ms, calculates T_WInterior E_i, such as formula

Wherein E_iIt is the i-th road of earthquake record the average energy value, x_i(t_j) it is j-th sampled point unit in the road of earthquake record i-th Element, i=1 in this example, 2 ... n, n=48, j=1,2 ... k, k=100；

B, definition l roads characterization factorAccording to λ_lValue identification interference way, E_lFor earthquake record l roads exist T_WInterior the average energy value, q is T in this example_WInterior first, E₁=0.0042, wherein l=1,2 ... n, identification interference stage property Body step is as follows：

If b1, there is λ to all roads_l＜ 5, shows the earthquake record without notable Hz noise, without follow-up denoising Journey, meets λ in this example_lThe road of ＜ 5 is the 1st~7 road, the 9th~42 road, otherwise the 44th~48 road, execution step b2；

If b2, λ_l>=5, show that earthquake record has notable Hz noise, note l roads are interference way, are designated as U_l, it is full in this example Sufficient λ_l>=5 road is the 8th and the 43rd road, is designated as U₈And U₄₃, E₈=0.078, E₄₃=0.23, λ₈=18.57, λ₄₃=54.76；

C, with U₈As a example by, extract the 6th, 7,8,9,10 road and constitute sub- seismic interference road collection, X is designated as, every in X is together one Row, obtain the average of X each columns, and the average of the row is subtracted with each element in X, obtain the sub- seismic channel set after X treatment, note It is X^*；

D, calculating X^*Covariance matrix Γ, such as formula

Wherein b is X^*Columns, X^*TIt is X^*Transposed matrix, " " representing matrix multiplication, b=5 in formula；

E, eigenvalue matrix Λ and eigenvectors matrix R that covariance matrix Γ is calculated using singular value decomposition method, then There is formula

Γ=R Λ R^T (3)

Wherein Λ is characterized the diagonal matrix that value is arranged from big to small, each characteristic vectors for being classified as character pair value of R, R^T It is transposed matrix, meets R^TR=RR^T=E, E are unit matrix；

F, X are through R^TLinear Mapping, obtains principal component matrix, such as formula

Φ=R^T·X (4)

Wherein Φ is main component matrix；

G, Φ the first rows are set to 0, obtain Φ ', made

X₁ ^*=R Φ ' (5)

Wherein X₁ ^*It is restructuring matrix；

H, extraction X₁ ^*The data of middle correspondence interference way, are designated as U_l', to U_l' carry out relative amplitude preserved processing, such as formula

Wherein U_l ^*It is interference way signal, A after relative amplitude preserved processing₁It is the average value of interference way signal amplitude absolute value after denoising, A₂It is and U_lAverage, the A in this example of the absolute value of amplitude average value sum of adjacent two track data₂=0.0689, A₁=0.0057.

Claims (1)

1. a kind of seismic data Hz noise automatic identification and drawing method based on principal component analysis, comprise the following steps：

A, for pending original seismic data, choose the when window T that the record time is located at middle hypomere_W, calculate the earthquake record and exist T_WNei Ge roads the average energy value, such as formula

$E_i=\frac{1}{k}\underset{j=1}{\overset{k}{\Σ}}{x}_i{\left(t_j\right)}^2---\left(1\right)$

Wherein E_iFor when window in the i-th road of earthquake record the average energy value, x_i(t_j) for when window in j-th of the road of earthquake record i-th adopt Sample value, i=1,2 ... n, n be the total road number of earthquake record, j=1,2 ... k, when k is in window earthquake record sampling number, when window T_WLength suggestion for summary journal lengthArriveIn time；

B, definition l roads characterization factorAccording to λ_lValue identification interference way, E_lIt is earthquake record l roads in T_WIt is interior The average energy value, q is T_WThe interior earthquake Taoist monastic name with minimum average B configuration energy value, wherein l=1,2 ... n, identification interference stage property Body step is as follows：

If b1, there is λ to all roads_l＜ 5, shows the earthquake record without notable Hz noise, no without follow-up denoising process Then perform step b2；

If b2, λ_l>=5, show that earthquake record has notable Hz noise, note l roads are interference way, are designated as U_l；

C, extraction and U_lIt is adjacent and including U_lContinuous 5 trace record constitute sub- seismic interference road collection, X is designated as, per being together one in X Row, obtain the average of X each columns, and the average of the row is subtracted with each element in X, obtain the sub- seismic channel set after X treatment, note It is X^*；

D, calculating X^*Covariance matrix Γ, such as formula

$\mathrm{\Γ}=\frac{X^{*}\·X^{{*}^T}}{b-1}---\left(2\right)$

Wherein b is X^*Columns, X^*TIt is X^*Transposed matrix, " " representing matrix multiplication；

E, eigenvalue matrix Λ and eigenvectors matrix R that covariance matrix Γ is calculated using singular value decomposition method, then be present Formula

Γ=R Λ R^T (3)

Wherein Λ is characterized the diagonal matrix that value is arranged from big to small, each characteristic vectors for being classified as character pair value of R, R^TTo turn Matrix is put, R is met^TR=RR^T=E, E are unit matrix；

F, X are through R^TLinear Mapping, obtains principal component matrix, such as formula

Φ=R^T·X (4)

Wherein Φ is main component matrix；

G, Φ the first rows are set to 0, obtain Φ ', made

X₁ ^*=R Φ ' (5)

Wherein X₁ ^*It is restructuring matrix；

H, extraction X₁ ^*The data of middle correspondence interference way, are designated as U_l', to U_l' carry out relative amplitude preserved processing, such as formula

${U_l}^{*}={U_l}^{\′}\·\frac{A_2}{A_1}---\left(6\right)$

Wherein U_l ^*It is interference way signal, A after relative amplitude preserved processing₁It is the average value of interference way signal amplitude absolute value after denoising, A₂Be with U_lThe absolute value of amplitude average value sum of adjacent two track data it is average.