slice2D_toy_model.m 7.42 KB
Newer Older
1 2 3 4 5 6
%3D plot of model

clear all;
clc
close all

7 8 9 10 11 12
fname='toy.';
parameter='vp';
source_file='../par/sources/sources_toy.dat';
rec_file='../par/receiver/receiver_toy.dat';

iteration=60;
13

14
nx=160; ny=184; nz=160; %ny:vertical
15 16 17 18 19 20 21 22 23 24
outx=1; outy=1; outz=1; 
dh=0.8;
FW=10;
nx=nx/outx;ny=ny/outy;nz=nz/outz;

X=dh:dh*outx:nx*outx*dh;
Y=dh:dh*outy:ny*outy*dh;
Z=dh:dh*outz:nz*outz*dh;

% Slizes of the 3D-Model in m 
25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40
slx=56;
sly=24;
slz=52.8;

if strcmp(parameter,'vs')
caxis_value_1=3300;%vs
caxis_value_2=3900;%vs
end
if strcmp(parameter,'vp')
caxis_value_1=5700;%vp
caxis_value_2=6700;%vp
end


% No changes to be done after this line---------
%------------------------------------------------
41

42 43 44 45 46 47 48 49 50 51 52
%Input of receiver and source positions--

fp = fopen(source_file);
SRC=textscan(fp, '%f %f %f %f %f %f', 'delimiter','\t', 'commentstyle','%','MultipleDelimsAsOne',1);
fclose(fp);

fp = fopen(rec_file);
REC=textscan(fp, '%f %f %f ', 'delimiter','\t', 'commentstyle','#','MultipleDelimsAsOne',1);
fclose(fp);

% Lines to mark PML position
53 54 55 56 57 58
fontsize=14;
x_line=[FW FW nx-FW nx-FW FW]*dh;
y_line=[FW ny-FW ny-FW FW FW]*dh;
z_line=[FW nz-FW nz-FW FW FW]*dh;


59 60 61 62 63 64
%m --> model indices
slizex=int64(slx/dh/outx);
slizey=int64(sly/dh/outy);
slizez=int64(slz/dh/outz);


65 66 67
file_inp1=['../par/model/' fname parameter '.true'];
file_inp2=['../par/model/' fname parameter '_it' num2str(iteration)];

68 69 70 71



%Read 3D Models------------------------
72 73 74 75 76 77 78 79 80 81 82 83 84 85

fid=fopen(file_inp1,'r','ieee-le');
modelvec=zeros(ny/outy,nx/outx,nz/outz);
modelvec=fread(fid,(nx*ny*nz),'float');

model_true=reshape(modelvec,ny,nx,nz);

fid=fopen(file_inp2,'r','ieee-le');
modelvec=zeros(ny/outy,nx/outx,nz/outz);
modelvec=fread(fid,(nx*ny*nz),'float');

model=reshape(modelvec,ny,nx,nz);


86 87
%Creation of 2D Slices--------------

88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132
model_truex=zeros(ny,nz);
model_truey=zeros(nz,nx);
model_truez=zeros(ny,nx);

for y=1:1:ny
   for z=1:1:nz
    model_truex(y,z)=model_true(y,slizex,z);
   end
end

for x=1:1:nx
   for z=1:1:nz  
    model_truey(z,x)=model_true(slizey,x,z);
   end
end

for x=1:1:nx
   for y=1:1:ny  
    model_truez(y,x)=model_true(y,x,slizez);
   end
end


modelx=zeros(ny,nz);
modely=zeros(nz,nx);
modelz=zeros(ny,nx);

for y=1:1:ny
   for z=1:1:nz
    modelx(y,z)=model(y,slizex,z);
   end
end

for x=1:1:nx
   for z=1:1:nz  
    modely(z,x)=model(slizey,x,z);
   end
end

for x=1:1:nx
   for y=1:1:ny  
    modelz(y,x)=model(y,x,slizez);
   end
end

133 134
%Plotting-------------------------------------------

135 136 137 138 139
figure(1)
imagesc(Z,Y,model_truex);
line(x_line,y_line,'LineStyle','--','Color','k')
colb=colorbar;
coll=get(colb,'xlabel');
140 141
set(coll,'String',[parameter ' in m/s'],'FontSize',fontsize);
colormap('jet');
142 143 144 145 146 147 148
caxis([caxis_value_1 caxis_value_2])
set(gca,'ydir','normal');
%set(gca,'xdir','reverse');
xlabel('z in m','FontSize',fontsize)
ylabel('y in m','FontSize',fontsize)
title(['true ' parameter '-model. Slize at x=' num2str(slx) 'm'],'FontSize',fontsize)
axis tight
149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164
hold on
%Plot sources if in model slice
for ii=1:length(SRC{1})
    if SRC{1}(ii)==slx;
        plot(SRC{3}(ii),SRC{2}(ii),'*blue','MarkerSize',15,'LineWidth',2.0);
    end
    
end
%Plot receiver if in model slice
for ii=1:length(REC{1})
    if REC{1}(ii)==slx;
        plot(REC{3}(ii),REC{2}(ii),'.blue','MarkerSize',15,'LineWidth',2.0);
    end
    
end
hold off
165

166 167
figure(2)
imagesc(Z,Y,modelx);
168 169 170 171
line(x_line,y_line,'LineStyle','--','Color','k')
colb=colorbar;
coll=get(colb,'xlabel');
set(coll,'String',[parameter ' in m/s'],'FontSize',fontsize); 
172
colormap('jet');
173 174 175
caxis([caxis_value_1 caxis_value_2])
set(gca,'ydir','normal');
%set(gca,'xdir','reverse');
176
xlabel('z in m','FontSize',fontsize)
177
ylabel('y in m','FontSize',fontsize)
178
title([ parameter '-model. Iteration ' num2str(iteration) ' Slize at x=' num2str(slx) 'm'],'FontSize',fontsize)
179
axis tight
180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195
hold on
%Plot sources if in model slice
for ii=1:length(SRC{1})
    if SRC{1}(ii)==slx;
        plot(SRC{3}(ii),SRC{2}(ii),'*blue','MarkerSize',15,'LineWidth',2.0);
    end
    
end
%Plot receiver if in model slice
for ii=1:length(REC{1})
    if REC{1}(ii)==slx;
        plot(REC{3}(ii),REC{2}(ii),'.blue','MarkerSize',15,'LineWidth',2.0);
    end
    
end
hold off
196

197 198 199
figure(3)
imagesc(X,Z,model_truey);
line(x_line,z_line,'LineStyle','--','Color','k')
200 201 202
colb=colorbar;
coll=get(colb,'xlabel');
set(coll,'String',[parameter ' in m/s'],'FontSize',fontsize); 
203
colormap('jet');
204 205
caxis([caxis_value_1 caxis_value_2])
set(gca,'ydir','normal');
206 207 208
xlabel('x in m','FontSize',fontsize)
ylabel('z in m','FontSize',fontsize)
title(['true ' parameter '-model. Slize at x=' num2str(sly) 'm'],'FontSize',fontsize)
209
axis tight
210 211 212 213 214 215 216 217 218 219 220 221 222 223 224 225
hold on
%Plot sources if in model slice
for ii=1:length(SRC{1})
    if SRC{2}(ii)==sly;
        plot(SRC{1}(ii),SRC{3}(ii),'*blue','MarkerSize',15,'LineWidth',2.0);
    end
    
end
%Plot receiver if in model slice
for ii=1:length(REC{1})
    if REC{2}(ii)==sly;
        plot(REC{1}(ii),REC{3}(ii),'.blue','MarkerSize',15,'LineWidth',2.0);
    end
    
end
hold off
226 227 228 229 230 231 232 233


figure(4)
imagesc(X,Z,modely);
line(x_line,z_line,'LineStyle','--','Color','k')
colb=colorbar;
coll=get(colb,'xlabel');
set(coll,'String',[parameter ' in m/s'],'FontSize',fontsize); 
234
colormap('jet');
235 236 237 238 239 240
caxis([caxis_value_1 caxis_value_2])
set(gca,'ydir','normal');
xlabel('x in m','FontSize',fontsize)
ylabel('z in m','FontSize',fontsize)
title([parameter '-model. Iteration ' num2str(iteration) ' Slize at x=' num2str(sly) 'm'],'FontSize',fontsize)
axis tight
241 242 243 244 245 246 247 248 249 250 251 252 253 254 255 256
hold on
%Plot sources if in model slice
for ii=1:length(SRC{1})
    if SRC{2}(ii)==sly;
        plot(SRC{1}(ii),SRC{3}(ii),'*blue','MarkerSize',15,'LineWidth',2.0);
    end
    
end
%Plot receiver if in model slice
for ii=1:length(REC{1})
    if REC{2}(ii)==sly;
        plot(REC{1}(ii),REC{3}(ii),'.blue','MarkerSize',15,'LineWidth',2.0);
    end
    
end
hold off
257

258 259 260 261 262 263
figure(5)
imagesc(X,Y,model_truez);
line(x_line,y_line,'LineStyle','--','Color','k')
colb=colorbar;
coll=get(colb,'xlabel');
set(coll,'String',[parameter ' in m/s'],'FontSize',fontsize); 
264
colormap('jet');
265 266 267 268 269 270 271 272 273 274 275 276 277 278 279 280 281 282 283 284 285 286 287
caxis([caxis_value_1 caxis_value_2])
set(gca,'ydir','normal');
%set(gca,'xdir','reverse');
xlabel('x in m','FontSize',fontsize)
ylabel('y in m','FontSize',fontsize)
title(['true ' parameter '-model. Slize at x=' num2str(slz) 'm'],'FontSize',fontsize)
axis tight
hold on
%Plot sources if in model slice
for ii=1:length(SRC{1})
    if SRC{3}(ii)==slz;
        plot(SRC{1}(ii),SRC{2}(ii),'*blue','MarkerSize',15,'LineWidth',2.0);
    end
    
end
%Plot receiver if in model slice
for ii=1:length(REC{1})
    if REC{3}(ii)==sly;
        plot(REC{1}(ii),REC{2}(ii),'.blue','MarkerSize',15,'LineWidth',2.0);
    end
    
end
hold off
288 289 290 291 292 293

figure(6)
imagesc(X,Y,modelz);
line(x_line,y_line,'LineStyle','--','Color','k')
colb=colorbar;
coll=get(colb,'xlabel');
294
colormap('jet');
295 296 297 298 299 300 301 302
set(coll,'String',[parameter ' in m/s'],'FontSize',fontsize); 
caxis([caxis_value_1 caxis_value_2])
set(gca,'ydir','normal');
%set(gca,'xdir','reverse');
xlabel('x in m','FontSize',fontsize)
ylabel('y in m','FontSize',fontsize)
title([parameter '-model. Iteration ' num2str(iteration) ' Slize at x=' num2str(slz) 'm'],'FontSize',fontsize)
axis tight
303 304 305 306 307 308 309 310 311 312 313 314 315 316 317 318
hold on
%Plot sources if in model slice
for ii=1:length(SRC{1})
    if SRC{3}(ii)==slz;
        plot(SRC{1}(ii),SRC{2}(ii),'*blue','MarkerSize',15,'LineWidth',2.0);
    end
    
end
%Plot receiver if in model slice
for ii=1:length(REC{1})
    if REC{3}(ii)==sly;
        plot(REC{1}(ii),REC{2}(ii),'.blue','MarkerSize',15,'LineWidth',2.0);
    end
    
end
hold off
319 320


321 322
 %exportfig(2, [fname parameter '.eps'],'bounds','tight', 'color','rgb', ...
 % 'preview','none', 'resolution',200, 'lockaxes',1);