Closed form not the same as the discrete form
1 回表示 (過去 30 日間)
古いコメントを表示
clear;
lambda = 3e-2;
x = 4 * pi/lambda * linspace(eps, 15, 100000);
T = 5e-3;
t = [0:0.001e-3:T] ; % 0.1:1e-3:0.1+T];
u = 3;
a = 4*pi/lambda * u;
for i = 1:length(x)
Z(i) = sum(-((cos(a.*t) - cos(x(i).*t)).^2 + (sin(a.*t) - sin(x(i).*t)).^2));
end
% Z1 = csc((a+x)/2) .* sin(1/2.*(2.*T + 1).*(a + x)) -1/2 .* csc(a) .* sin(2.*a.*T+a) - 1/2 .* csc(x) .* sin(2.*T.*x+x)
Z1 = csc((a-x)/2) .* sin((T+1/2).*(a-x))/(2*pi) - 2*T - 1;
figure; plot(x*lambda/(4 * pi), (Z));
hold on; plot(x*lambda/(4 * pi), (Z1), '*');
The above code does two things. one is a summation over `t` and stores it in `Z` for all `x` individually with a for loop. The second thing it is doing is trying an analytical function for the summation in `t` and and just use the vector `x` to find the same quantity as `Z`, but analytically. it is stored in `Z1`. Somehow they both are not the same. The expression of the analytical form can be found in https://mathworld.wolfram.com/DeltaFunction.html equation (41) on that page. There is also a plot representing the function there. However, the matlab plot is very different. However, it is interesting to notice that `Z` and `Z1` are the same at the limiting case `a == x`. The plot I get is shown below. It should have a maximum at `a = x`, that is the value it should take at that point is `Z = Z1 = 0` and for the rest of the `x`, it should be less than 0. However, the blue plot satisfies it, but not the red one.
0 件のコメント
採用された回答
Paul
2022 年 6 月 14 日
Hi Tworit,
I'm not quite following how the code is trying to implement equation 41 from the linked page. Here is one way to create a plot similar to that just above equation (41). Maybe you can modify this as needed? Also, where does the summation equation for Z come from?
n = 10:10:50; % different values of n
x = (-1:.001:1).';
x(x==0) = eps; % avoid divide by zero to keep things simple
Z1 = csc(x/2).*sin((n+1/2).*x)/2/pi;
figure
plot(x,Z1)
To shift the whole thing by 'a' :
a = 1; % for example
n = 10:10:50;
x = a + (-1:.001:1).';
x(x==a) = a+eps; % avoid divide by zero
Z1 = csc((x-a)/2).*sin((n+1/2).*(x-a))/2/pi;
figure
plot(x,Z1)
その他の回答 (0 件)
参考
製品
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!