How to plot for this code?
3 Ansichten (letzte 30 Tage)
Ältere Kommentare anzeigen
clc
clear all
close all
syms k1 k2
eta0 = 1; % Maximum displacement(m)
L = 100; % Propagate length (km)
W = 100; % Source width (km)
v = 0.14; % Ruptuer velocity to obtain maximum surface amplitude(km/s)
h = 2; % Water depth (km)
g = 0.0098; % Acceleration due to gravity (km/s^2)
t1 = 50/v; % (s)
t = (t1); % (s)
k = sqrt(k1.^2+k2.^2);
w = sqrt(g.*k.*tanh(k.*h));
A = ((sin(w*t))./(w.*cosh(k*h))).*((eta0*v)./(2*L));
B = (1-exp(-1i*100*k1))/(1i*k1);
C = ((exp(-1i*100*k1))./(1-((50*k1)./pi).^2)).*(1i*k1).*((50./pi).^2).*((exp(1i*100*k1))-1);
D = (exp(1i*150*k2)-exp(1i*50*k2))./(1i*k2);
E = (1./(1-(((100*k2)./pi).^2))).*(1i*k2).*((100./pi)^2).*(exp(1i*50*k2)+exp(1i*150*k2));
F = (4*sin(50*k2))./k2;
G = (exp(-1i*50*k2)-exp(-1i*150*k2))./(1i*k2);
H = (1./(1-(((100*k2)./pi).^2))).*(1i*k2).*((100./pi).^2).*(exp(-1i*150*k2)+exp(-1i*50*k2));
K1 = 0:100;
K2 = -150:150;
for x1 = 1:101;
for x2 = 1:301;
k1 = K1(x1);
k2 = K2(x2);
if (k==0)
eta_(x1,x2)= 200*eta0*v*t;
else
if (k1 == 0)
eta_(x1,x2)= (A.*(100-C)).*((D-E)+F+(G-H));
else
eta_(x1,x2)= (A.*(B-C)).*((D-E)+F+(G-H));
end
if (k2 == 0)
eta_(x1,x2)= (A.*(B-C)).*((100-E)+200+(100-H));
else
eta_(x1,x2)= (A.*(B-C)).*((D-E)+F+(G-H));
end
end
end
end
ETA = ifft2(eta_)
X = abs(ifft2(ETA,K1,K2));
X = ifftshift(X);
F = [-K1/2:K1/2-1]/K1;
plot(F,X)
0 Kommentare
Antworten (1)
Image Analyst
am 21 Feb. 2021
Since youi're calline fft2(), you're dealing with images not 1-D signals. So call imshow() instead of plot(), unless you want to display some cross sections (profiles) through the image.
0 Kommentare
Siehe auch
Kategorien
Mehr zu Calculus finden Sie in Help Center und File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!