clc;
clear all;
%path_root = '.\marmousi2_cut';
path_vp = ('E:\yfmmm\FWI源码和可执行程序\FWI可执行文件\Marmousi2\model\pup2_vp.dat');
path_vs = ('E:\yfmmm\FWI源码和可执行程序\FWI可执行文件\Marmousi2\model\pup2_vs.dat');
path_rho = ('E:\yfmmm\FWI源码和可执行程序\FWI可执行文件\Marmousi2\model\pup2_rho.dat');
nx=175;
ny=250;
fid = fopen(path_vp,'rb');
vp= fread(fid,[ny,nx],'float32');
fclose(fid);
fid = fopen(path_vs,'rb');
vs= fread(fid,[ny,nx],'float32');
fclose(fid);
fid = fopen(path_rho,'rb');
rho= fread(fid,[ny,nx],'float32');
fclose(fid);
figure
imagesc([0:175]*4,[0:250]*4,vp);
colorbar;
xlabel('x / m');
ylabel('z / m');
axis image
figure
imagesc([0:175]*4,[0:250]*4,vs);
colorbar;
xlabel('x / m');
ylabel('z / m');
axis image
figure
imagesc([0:175]*4,[0:250]*4,rho);
colorbar;
xlabel('x / m');
ylabel('z / m');
axis image
vpp=vp;
vss=vs;
rhoo=rho;
f=fspecial('gaussian',[60 60],50);%%%%%%%%%%%%生成一个滤波器f
vp_smooth=imfilter(vp,f,'replicate');
figure
imagesc([0:175]*4,[0:250]*4,vp_smooth);
colorbar;
xlabel('x / m');
ylabel('z / m');
axis image
f=fspecial('gaussian',[60 60],50);%%%%%%%%%%%%生成一个滤波器f
vs_smooth=imfilter(vs,f,'replicate');
figure
imagesc([0:175]*4,[0:250]*4,vs_smooth);
colorbar;
xlabel('x / m');
ylabel('z / m');
axis image
%
f=fspecial('gaussian',[60 60],50);%%%%%%%%%%%%生成一个滤波器f
rho_smooth=imfilter(rho,f,'replicate');
figure
imagesc([0:175]*4,[0:250]*4,rho_smooth);
colorbar;
xlabel('x / m');
ylabel('z / m');
axis image
%%%%%%%%%%%%%%%%当只有vp模型时，利用公式对VS和rho进行计算
% vs_smooth=vp_smooth/(sqrt(3.0));
% vs_smooth(1:23,:)=0;%%%%%%%%%%%%由于之前的模型是marmousi，因此会有一段深度为海水，横波为0，速度固定为1000
% rho_smooth=0.31*1000*vp_smooth.^(1/4);
% rho_smooth(1:23,:)=1000;
% vs_smooth=imfilter(vs,f,'replicate');
% figure
% imagesc(vs_smooth);
% colorbar;
% axis image
%
%
% rho_smooth=imfilter(rho,f,'replicate');
% figure
% imagesc(rho_smooth);
% colorbar;
% axis image
fid = fopen('E:\yfmmm\FWI源码和可执行程序\FWI可执行文件\Marmousi2\model\pup2_init_vp.dat','wb');
fwrite(fid,vpp,'float32');
fclose(fid);
fid = fopen('E:\yfmmm\FWI源码和可执行程序\FWI可执行文件\Marmousi2\model\pup2_init_vs.dat','wb');
fwrite(fid,vss,'float32');
fclose(fid);
fid = fopen('E:\yfmmm\FWI源码和可执行程序\FWI可执行文件\Marmousi2\model\pup2_init_rho.dat','wb');
fwrite(fid,rhoo,'float32');
fclose(fid);
