Plotting the product of two functions.

2 views (last 30 days)
Hi,
I have plotted two offset Bessel functions of the first kind (blue and red dotted curves) - see figure below.
I would like to produce something similar to the black curve shown in the second attached figure below, but for the two curves I have plotted. Can someone please point me in the correct direction?
Thank you
lambda = 500e-6; % wavelength [m]
D = 0.0001 ; % diameter D of the primary aperture [m]
theta1 = -20:.01:20; % angular radius θ (as measured from the primary aperture) [Degrees]
u1 = (pi/lambda)*D*theta1;
theta2 = -20:.01:20;
u2 = (pi/lambda)*D*theta2;
x1 = u1./pi +0.61;
x2 = u2./pi -0.61;
I01 = 1.0;
I02 = 1.0;
i1 = real((besselj(1,u1)./(u1)).^2);
i2 = real((besselj(1,u2)./(u2)).^2);
I1 = I01*(i1./max(i1));
I2 = I02*(i2./max(i2));
colormap('default')
figure(1);
%clf;
hold on;
plot(x1,I1,'b:','LineWidth',1.5);
plot(x2,I2,'r:','LineWidth',1.5);
hold off;
grid on;
set(gca,'FontSize',14);
%axis([-4 4 0 1.05])
xlabel('\theta (\lambda/D)', 'fontsize', 16);
ylabel('Normalized Intensity', 'fontsize', 16);
legend('Source 1','Source 2','fontsize', 10);
legend;
figure(2);
clf;
[X,Y]=meshgrid(-10:.1:10);
%[X,Y]=meshgrid(x);
R=sqrt(X.^2+Y.^2);
z=(besselj(1,R)./R).^2;
Z=z./max(max(z));
surf(X,Y,Z,'EdgeColor','none');
colorbar
camlight left;lighting phong;
set(gca,'XTick',[], 'YTick', [], 'ZTick', [])
  1 Comment
John D'Errico
John D'Errico on 7 May 2020
Edited: John D'Errico on 7 May 2020
The black curve MIGHT be the sum of the two, or some linear combination thereof, but essentially NEVER the product.

Sign in to comment.

Accepted Answer

Ameer Hamza
Ameer Hamza on 7 May 2020
As pointed out by John, this is the sum of these two curves. Check the following code. I had to do a bit of interpolating because of the way you defined x1 and x2. I have added comments on the modified line. The graph looks quite similar to the graph you want
lambda = 500e-6; % wavelength [m]
D = 0.0001 ; % diameter D of the primary aperture [m]
theta1 = -20:.01:20; % angular radius θ (as measured from the primary aperture) [Degrees]
u1 = (pi/lambda)*D*theta1;
theta2 = -20:.01:20;
u2 = (pi/lambda)*D*theta2;
x1 = u1./pi +0.8;
x2 = u2./pi -0.8;
I01 = 1.0;
I02 = 0.5; % <<==== changed this parameter
i1 = real((besselj(1,u1)./(u1)).^2);
i2 = real((besselj(1,u2)./(u2)).^2);
I1 = I01*(i1./max(i1));
I2 = I02*(i2./max(i2));
% vvv Interpolated to make both vector equal
x = unique([x1 x2]);
I1_ = interp1(x1, I1, x);
I2_ = interp1(x2, I2, x);
colormap('default')
figure(1);
%clf;
hold on;
plot(x1,I1,'b:','LineWidth',1.5);
plot(x2,I2,'r:','LineWidth',1.5);
plot(x,I1_+I2_,'r:','LineWidth',1.5); % <<==== added this line
hold off;
grid on;
set(gca,'FontSize',14);
%axis([-4 4 0 1.05])
xlabel('\theta (\lambda/D)', 'fontsize', 16);
ylabel('Normalized Intensity', 'fontsize', 16);
legend('Source 1','Source 2','fontsize', 10);
legend;
figure(2);
clf;
[X,Y]=meshgrid(-10:.1:10);
%[X,Y]=meshgrid(x);
R=sqrt(X.^2+Y.^2);
z=(besselj(1,R)./R).^2;
Z=z./max(max(z));
surf(X,Y,Z,'EdgeColor','none');
colorbar
camlight left;lighting phong;
set(gca,'XTick',[], 'YTick', [], 'ZTick', [])
  2 Comments
Colin Ross
Colin Ross on 7 May 2020
Thank you Ameer! This is exactly what I was loking for.

Sign in to comment.

More Answers (0)

Categories

Find more on Graphics Object Programming in Help Center and File Exchange

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!