ºÝºÝߣ

ºÝºÝߣShare a Scribd company logo
SCRIPT MATLAB
clc, clear all
%DEFINISI parameter PUSAT model hiposenter
xhipo=200; yhipo=400; zhipo=1000;
v_avg=7.8;
%tentukan koordinat stasiun pengukuran
x=[300 700 1350 1800];
y=[1200 200 450 650];
z=zeros(1,4);
%menghitung DATA t_obs
t_obs=zeros(length(x),1);
for i=1:length(x)
t_obs(i)=sqrt(((xhipo-x(i))^2 + (yhipo-y(i))^2 + (zhipo-z(i))^2))
/v_avg;
end
%lakukan PREDIKSI MODEL awal
iterasi=1; eps=1;
while eps >= 0.00001
if iterasi==1
xo_model=1600; yo_model=1600; zo_model=500;
else
xo_model=xo_pertu;
yo_model=yo_pertu;
zo_model=zo_pertu;
end
%menghitung DATA g_cal
t_cal=zeros(length(x),1);
for i=1:length(x)
t_cal(i)=sqrt(((xo_model-x(i))^2 + (yo_model-y(i))^2 + ...
(zo_model-z(i))^2)) / v_avg;
end
%menghitung MISFIT
dt_misfit=t_obs - t_cal;
eps=std(abs(dt_misfit));
e_plot(iterasi)=eps;
%membuat kondisi IF untuk inversi jacobi
if (std(abs(dt_misfit))) >= 0.00001
%membuat matriks JACOBI
for i=1:length(x)
derivative_x(i)=(xo_model-x(i)) / (v_avg*sqrt((xo_model-x(i))^2 +
(yo_model-y(i))^2 + (zo_model-z(i))^2));
derivative_y(i)=(yo_model-y(i)) / (v_avg*sqrt((xo_model-x(i))^2 +
(yo_model-y(i))^2 + (zo_model-z(i))^2));
derivative_z(i)=(zo_model-z(i)) / (v_avg*sqrt((xo_model-x(i))^2 +
(yo_model-y(i))^2 + (zo_model-z(i))^2));
end
J=ones(length(x),3);
J(:,1)=derivative_x';
J(:,2)=derivative_y';
J(:,3)=derivative_z';
%menghitung PERTURBASI MODEL
dm_perturbasi=inv(J'*J)*J'*dt_misfit;
xo_pertu=xo_model + dm_perturbasi(1);
yo_pertu=yo_model + dm_perturbasi(2);
zo_pertu=zo_model + dm_perturbasi(3);
iterasi=iterasi+1;
end
end
figure(1)
plot3(xo_model,yo_model,zo_model,'o',...

...

...
...
...
'MarkerFaceColor','r','MarkerEdgeColor','r')
grid on; set(gca,'zdir','reverse');
xlabel('Koord. X (m)'); ylabel('Koord. Y (m)');
zlabel('Kedalaman (m)'); title('Model Prediksi');
legend('Titik Hiposenter');
figure(2)
plot([1:1:(length(e_plot))],e_plot,'Color','m', ...
'LineStyle','-','LineWidth',2)
xlabel('Iterasi'); ylabel('Std Misfit'); title('Grafik Misfit')

maka dari script MATLAb diatas jika di-run akan memberikan output berupa grafik sebagai berikut,

More Related Content

Pendekatan Inversi Linier dengan Matriks Jacobi pada Kasus Perhitungan Hiposenter Gempa Sederhana

  • 1. SCRIPT MATLAB clc, clear all %DEFINISI parameter PUSAT model hiposenter xhipo=200; yhipo=400; zhipo=1000; v_avg=7.8; %tentukan koordinat stasiun pengukuran x=[300 700 1350 1800]; y=[1200 200 450 650]; z=zeros(1,4); %menghitung DATA t_obs t_obs=zeros(length(x),1); for i=1:length(x) t_obs(i)=sqrt(((xhipo-x(i))^2 + (yhipo-y(i))^2 + (zhipo-z(i))^2)) /v_avg; end %lakukan PREDIKSI MODEL awal iterasi=1; eps=1; while eps >= 0.00001 if iterasi==1 xo_model=1600; yo_model=1600; zo_model=500; else xo_model=xo_pertu; yo_model=yo_pertu; zo_model=zo_pertu; end %menghitung DATA g_cal t_cal=zeros(length(x),1); for i=1:length(x) t_cal(i)=sqrt(((xo_model-x(i))^2 + (yo_model-y(i))^2 + ... (zo_model-z(i))^2)) / v_avg; end %menghitung MISFIT dt_misfit=t_obs - t_cal; eps=std(abs(dt_misfit)); e_plot(iterasi)=eps; %membuat kondisi IF untuk inversi jacobi if (std(abs(dt_misfit))) >= 0.00001 %membuat matriks JACOBI for i=1:length(x) derivative_x(i)=(xo_model-x(i)) / (v_avg*sqrt((xo_model-x(i))^2 + (yo_model-y(i))^2 + (zo_model-z(i))^2)); derivative_y(i)=(yo_model-y(i)) / (v_avg*sqrt((xo_model-x(i))^2 + (yo_model-y(i))^2 + (zo_model-z(i))^2)); derivative_z(i)=(zo_model-z(i)) / (v_avg*sqrt((xo_model-x(i))^2 + (yo_model-y(i))^2 + (zo_model-z(i))^2)); end J=ones(length(x),3); J(:,1)=derivative_x'; J(:,2)=derivative_y'; J(:,3)=derivative_z'; %menghitung PERTURBASI MODEL dm_perturbasi=inv(J'*J)*J'*dt_misfit; xo_pertu=xo_model + dm_perturbasi(1); yo_pertu=yo_model + dm_perturbasi(2); zo_pertu=zo_model + dm_perturbasi(3); iterasi=iterasi+1; end end figure(1) plot3(xo_model,yo_model,zo_model,'o',... ... ... ... ...
  • 2. 'MarkerFaceColor','r','MarkerEdgeColor','r') grid on; set(gca,'zdir','reverse'); xlabel('Koord. X (m)'); ylabel('Koord. Y (m)'); zlabel('Kedalaman (m)'); title('Model Prediksi'); legend('Titik Hiposenter'); figure(2) plot([1:1:(length(e_plot))],e_plot,'Color','m', ... 'LineStyle','-','LineWidth',2) xlabel('Iterasi'); ylabel('Std Misfit'); title('Grafik Misfit') maka dari script MATLAb diatas jika di-run akan memberikan output berupa grafik sebagai berikut,