本文通过地震数据二次谱求取地震子波谱，具体方法如下：

MATLAB代码实现如下：

function w = SndSpecExtWavelet(x, M)
% 功能：基于二次谱提取输入地震数据data的地震子波wavelet
%      Extracting Wavelet from Input Seismic Data Based on Second Spectrum
% 输入：
%       x     -  地震信号,Nt×Nz
%       M     -  提取地震子波的长度
% 输出：
%       w     -  提取的地震子波,1×M
% 时间：2023/10/14


[Nt,Nz] = size(x);
FxTwo = zeros(Nt,Nz);
% 多道平均法
for k = 1:Nz
    Fx = abs(fft(x(:,k)));
    FxTwo(:,k) = (fft(Fx));
end
FxTwoAvg = sum(FxTwo,2)/Nz;

cLow = zeros(Nt, 1); %低通滤波器
for i = 1:Nt
    cLow(i) = 0.5*(1+cos(2*pi*(i-1)/Nt));
end
FxTwoLow = FxTwoAvg .* cLow;
FxLow = (ifft(FxTwoLow));
we = real(ifft(FxLow));
wm = zeros(1,M);
wm(ceil(M/2):M) = we(1:ceil(M/2));
wm(1:ceil(M/2)) = we(Nt-floor(M/2):Nt);
p = max(abs(wm));
w = wm./p; % 归一化

end

自己写的代码，如有错误还请联系。创作不易，记得点赞支持哦，谢谢！