目錄
1、數(shù)據(jù)來源。
1）SEG? 系列。
2）OpenFWI 系列。
2、數(shù)據(jù)量，尺寸。
1） SEG 包含兩個數(shù)據(jù)集：SEGSaltData 和?SimulateData。
2）OpenFWI 包含12個數(shù)據(jù)集：
3、地震數(shù)據(jù)對應的觀測系統(tǒng)。
1） SEG系列
2）OpenFWI系列
4、顯示數(shù)據(jù)的源碼
5、正演的原理及源碼
6、我的疑問：
地震數(shù)據(jù)是非常寶貴的資源，很多真實數(shù)據(jù)并不是公開的，目前在網(wǎng)上流傳的都是合成數(shù)據(jù)。我將從以下角度來介紹合成數(shù)據(jù)：
1、數(shù)據(jù)來源。
1）SEG? 系列。
數(shù)據(jù)源自論文： Yang F, Ma J. Deep-learning inversion: A next-generation seismic velocity model building method DL for velocity model building[J]. Geophysics, 2019, 84(4): R583-R599.
作者基于unet架構，實現(xiàn)了端到端的全波形反演，從地震數(shù)據(jù)直接獲得速度模型，并在鹽體數(shù)據(jù)上測試性能。
論文代碼：GitHub - YangFangShu/FCNVMB-Deep-learning-based-seismic-velocity-model-building: Deep-learning inversion: A next-generation seismic velocity model building method
2）OpenFWI 系列。
數(shù)據(jù)源自論文：Deng C, Feng S, Wang H, et al. OpenFWI: Large-scale Multi-structural Benchmark Datasets for Full Waveform Inversion[J]. Advances in Neural Information Processing Systems, 2022, 35: 6007-6020.
數(shù)據(jù)和代碼網(wǎng)址：?https://openfwi-lanl.github.io/
它涵蓋了地球物理的多個領域(界面、斷層、CO2儲層等)，涵蓋了不同的地質地下構造(平面、曲線等)，包含了不同數(shù)量的數(shù)據(jù)樣本(2K - 67K)，還包括3D FWI數(shù)據(jù)集。此外，比較了三種二維FWI方法：InversionNet[20]提出了一種全卷積網(wǎng)絡來模擬地震反演過程；VelocityGAN[27]采用基于gan的模型求解FWI；UPFWI[31]將正演建模與CNN連接在一個回路中，實現(xiàn)無監(jiān)督學習，無需地面真值速度圖進行訓練。還有一種三維FWI方法：InversionNet3D[30]將InversionNet擴展到三維領域。為了減少內存占用和提高計算效率(即3D反演中最具挑戰(zhàn)性的兩個障礙)，該網(wǎng)絡在編碼器中使用了群卷積，并通過基于加性耦合的可逆層采用了部分可逆架構[46]。
2、數(shù)據(jù)量，尺寸。
1） SEG 包含兩個數(shù)據(jù)集：SEGSaltData 和?SimulateData。
| 數(shù)據(jù)集 | SEGSaltData | SimulateData |
| 數(shù)量（train+test） | 140（130+10） | 1700（1600+100） |
| 地震數(shù)據(jù)尺寸 | 400 x 301 | 201 x 301 |
| 速度模型尺寸 | 400 x 301 | 201 x 301 |
2）OpenFWI 包含12個數(shù)據(jù)集：
其中第二行的?B系列，比第一行的A系列有更高的難度。從左到右，依次為平面（FlatVel-A、FlatVel-B）、曲面（CurveVel-A?、CurveVel-B）、平面斷層（FlatFault-A、?FlatFault-B）、曲面斷層（CurveFault-A、CurveFault-B）、復雜的野外合成數(shù)據(jù)（Style-A、Style-B）、二氧化碳儲層（Kimberlina-CO2）、三維地震數(shù)據(jù)（3D Kimberlina-V1）。
在這里要注意，理解下表的尺寸時，如速度模型70x1x70, 中間為1表明是個二維速度模型。
3、地震數(shù)據(jù)對應的觀測系統(tǒng)。
1） SEG系列
?參考星移論文
2）OpenFWI系列
4、顯示數(shù)據(jù)的源碼
我的數(shù)據(jù)都存在.mat文件中：
1）顯示地震地震數(shù)據(jù)：
2）顯示速度模型：
5、正演的原理及源碼
原理及波動方程公式后續(xù)補上。
% 對文件夾內所有文件進行計算并存儲
% function [] = forward(vdir,vname,gdir,gname,start_num,end_num)
% for i = start_num:end_num
% vfile = [vdir,'/',vname,num2str(i),'.mat']
% gfile = [gdir,'/',gname,num2str(i),'.mat']
% calc(vfile,gfile);
% end
% end
% 對單個速度模型文件計算多炮數(shù)據(jù)
function [] = forward(vfile,outfile,sn) %sn 炮數(shù)
load(vfile,'vmodel'); %加載vfile文件中的vmodel變量
[nz,nx] = size(vmodel);% 201*301
Rec = ones(400,nx,sn);% 400*301 是地震數(shù)據(jù)的尺寸
%tic用來保存當前時間，而后使用toc來記錄程序完成時間，兩者往往結合使用
tic;
for i = 1:sn
sx = i*fix(nx/sn); % 這一步的意思是放炮的位置是均勻放置的么？
singleRec = calc(vmodel, sx);
Rec(:,:,i) = singleRec;
end
toc;
save(outfile,'Rec');
end
% 對單炮速度模型進行計算
function singleRec = calc(vmodel, sx)
[nz,nx] = size(vmodel); % x方向網(wǎng)格數(shù)（從左到右 和 z方向網(wǎng)格數(shù)（從上到下） 201*301
npmlz = 10; % 頂部和底部的PML層厚度
npmlx = 10; % 左邊和右邊的PML層厚度
% sx = 100; % 震源的x方向網(wǎng)格數(shù)
sz = 0; % 震源的z方向網(wǎng)格數(shù) 表明震源在地表
dx = 10; % x方向的網(wǎng)格大小 單位：米
dz = 10; % y方向的網(wǎng)格大小 單位：米
nt = 2000; % 計算的總時間步數(shù) 在方程的最后下采樣5，最后求的400 對應301*400的地震數(shù)據(jù)尺寸
dt = 1e-3; % 單位時間步長度 單位：秒
ampl = 1.0e0; % 震源子波的震幅
xrcvr = 1:1:nx; % 地面上x方向的接收器位置
nodr = 3; % half of the order number for spatial difference.
B = [1 zeros(1,nodr - 1)]';
A = NaN*ones(nodr,nodr);
for i = 1:1:nodr
A(i,:) = (1:2:2*nodr - 1).^(2*i - 1);
end
C = A\B;
Nz = nz + 2*npmlz;
Nx = nx + 2*npmlx;
rho = 1000*ones(Nz,Nx); % 密度; Unit: kg/m^3. 密度對結果有什么影響？ 此處密度1000
rho(fix(Nz/3):end,fix(Nx/2):end) = 500; % 通過修改密度，可以模擬介質中存在不同的物理特性，如不均勻的巖層、氣體或液體的存在等
vp = padarray(vmodel,[10,10],'replicate','both'); % 擴展邊界
% 邊界擴展是為了處理波場模擬中的邊界效應，確保在模擬過程中聲波傳播不會超出vmodel的范圍。
% 通過在vmodel的周圍添加額外的邊界值，可以避免波場反射和干擾。padarray函數(shù)的'replicate'選項表示將vmodel的邊界值復制并填充到擴展后的邊界上。
f0 = 25; % 震源頻率 單位 Hz
t0 = 1/f0; % the time shift of source Ricker wavelet; Unit: s; Suggest: 0.02 if fm = 50, or 0.05 if fm = 20.
t = dt*(1:1:nt);
src = (1 - 2*(pi*f0.*(t - t0)).^2).*exp( - (pi*f0*(t - t0)).^2); % the time series of source wavelet. 雷克子波
% The source wavelet formula refers to the equations (18) of Collino and Tsogka, 2001.
%% Perfectly matched layer absorbing factor PML層吸收因子設置
dpml0z = 3*max(vp(:))/dz*(8/15 - 3/100*npmlz + 1/1500*npmlz^2);
dpmlz = zeros(Nz,Nx);
dpmlz(1:npmlz,:) = (dpml0z*((npmlz: - 1:1)./npmlz).^2)'*ones(1,Nx);
dpmlz(npmlz + nz + 1:Nz,:) = dpmlz(npmlz: - 1:1,:);
dpml0x = 3*max(vp(:))/dx*(8/15 - 3/100*npmlx + 1/1500*npmlx^2);
dpmlx = zeros(Nz,Nx);
dpmlx(:,1:npmlx) = ones(Nz,1)*(dpml0x*((npmlx: - 1:1)./npmlx).^2);
dpmlx(:,npmlx + nx + 1:Nx) = dpmlx(:,npmlx: - 1:1);
% The PML formula refers to the equations (2) and (3) of Marcinkovich and Olsen, 2003.
%% Wavefield calculating 依據(jù)廣義胡克定律求的彈性系數(shù)
% rho1 和 rho2 是密度（rho）的副本。它們用于計算波場的系數(shù)，與波場更新方程中的密度項有關。
% Coeffi1 和 Coeffi2 是沿 x 軸和 z 軸方向的PML吸收因子的系數(shù)。它們根據(jù)PML吸收因子（dpmlx 和 dpmlz）和時間步長（dt）計算得出。這些系數(shù)在波場更新方程中用于考慮PML的吸收效果。
% Coeffi3 和 Coeffi4 是與密度（rho）和空間步長（dx 和 dz）相關的系數(shù)，用于考慮波場更新方程中的空間導數(shù)項。
% Coeffi5 和 Coeffi6 是與密度（rho）和速度（vp）的平方以及空間步長（dx 和 dz）相關的系數(shù)，用于考慮波場更新方程中的速度和應力項。
% 這些系數(shù)是根據(jù)彈性介質的性質和PML吸收層的影響，結合波場更新方程中的相關項，計算得出的。它們的計算方式是基于廣義胡克定律和對介質參數(shù)的離散化模擬。
% 通過使用這些系數(shù)，可以準確地更新波場的值，模擬彈性波在介質中的傳播和衰減行為。
rho1 = rho; % or = [(rho(:,1:end - 1) + rho(:,2:end))./2 (2*rho(:,end) - rho(:,end - 1))];
rho2 = rho; % or = [(rho(1:end - 1,:) + rho(2:end,:))./2; (2*rho(end,:) - rho(end - 1,:))];
Coeffi1 = (2 - dt.*dpmlx)./(2 + dt.*dpmlx);
Coeffi2 = (2 - dt.*dpmlz)./(2 + dt.*dpmlz);
Coeffi3 = 1./rho1./dx.*(2*dt./(2 + dt.*dpmlx));
Coeffi4 = 1./rho2./dz.*(2*dt./(2 + dt.*dpmlz));
Coeffi5 = rho.*(vp.^2)./dx.*(2*dt./(2 + dt.*dpmlx));
Coeffi6 = rho.*(vp.^2)./dz.*(2*dt./(2 + dt.*dpmlz));
% +++++++++++++++++++++++++++++++++++++ approximate coeffient ++++++++++++++++++++++++++++++++++++++
% Coeffi1 = 1 - dt.*dpmlx;
% Coeffi2 = 1 - dt.*dpmlz;
% Coeffi3 = 1./rho./dx.*dt;
% Coeffi4 = 1./rho./dz.*dt;
% Coeffi5 = rho.*(vp.^2)./dx.*dt;
% Coeffi6 = rho.*(vp.^2)./dz.*dt;
% --------------------------------------------------------------------------------------------------
NZ = Nz + 2*nodr; % All values of the outermost some columns are set to zero to be a boundary condition: all of wavefield values beyond the left and right boundary are null.
NX = Nx + 2*nodr; % All values of the outermost some rows are set to zero to be a boundary condition: all of wavefield values beyond the top and bottom boundary are null.
Znodes = nodr + 1:NZ - nodr;
Xnodes = nodr + 1:NX - nodr;
znodes = nodr + npmlz + 1:nodr + npmlz + nz;
xnodes = nodr + npmlx + 1:nodr + npmlx + nx;
nsrcz = nodr + npmlz + sz;
nsrcx = nodr + npmlx + sx;
Ut = NaN*ones(NZ,NX); % the wavefield value preallocation. 存儲波場的值，并在時間步進過程中進行更新
Uz = zeros(NZ,NX); % The initial condition: all of wavefield values are null before source excitation. 波場在z方向上的速度分量初始條件
Ux = zeros(NZ,NX); % The initial condition: all of wavefield values are null before source excitation. 波場在x方向上的速度分量初始條件
Vz = zeros(NZ,NX); % The initial condition: all of wavefield values are null before source excitation. 波場在z方向上的位移分量初始條件
Vx = zeros(NZ,NX); % The initial condition: all of wavefield values are null before source excitation. 波場在x方向上的位移分量初始條件
Psum = NaN*ones(Nz,Nx); % 存儲波場在不同時間步長上的壓力分量的累積和
U = NaN*ones(nz,nx,nt); % 用于存儲完整的波場數(shù)據(jù)
% tic;
for it = 1:1:nt
Ux(nsrcz,nsrcx) = Ux(nsrcz,nsrcx) + ampl*src(it)./2;
Uz(nsrcz,nsrcx) = Uz(nsrcz,nsrcx) + ampl*src(it)./2;
Ut(:,:) = Ux(:,:) + Uz(:,:);
U(:,:,it) = Ut(znodes,xnodes);
end
% toc;
syngram(:,:) = U(1,xrcvr,:);
singleRec = syngram';
singleRec = downsample(singleRec,5); % 下采樣
end
% 實驗
6、我的疑問：
1）評價指標，目前論文SSIM、MAE、RMSE,? ?地質勘探的角度，有什么評價指標？尤其在自然數(shù)據(jù)中，工程應用場景下，更希望通過FWI提供什么信息？
2）鹽體數(shù)據(jù)中，地層很薄，對這類的識別是否很重要？
3）地震數(shù)據(jù)的噪聲來源：采集設備、地理環(huán)境，哪些噪聲對地震數(shù)據(jù)的影響特別大？文章來源：http://www.zghlxwxcb.cn/news/detail-493573.html
4）地震數(shù)據(jù)采集后，到目前拿到手上的數(shù)據(jù)，中間是否經(jīng)過了別的處理？這種處理是否是加入噪聲，或者減少信息量的。文章來源地址http://www.zghlxwxcb.cn/news/detail-493573.html
到了這里，關于FWI 地震數(shù)據(jù)的認識的文章就介紹完了。如果您還想了解更多內容，請在右上角搜索TOY模板網(wǎng)以前的文章或繼續(xù)瀏覽下面的相關文章，希望大家以后多多支持TOY模板網(wǎng)！
