How to find the intersection values of line and curve?
    6 views (last 30 days)
  
       Show older comments
    
How to find the intersection values  of line(black) and the curve(blue). And I need the intersection points in ascending order in array. Can you please help.....
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% This is a program for 1-D Photonic crystal
% For Photonic band structure calculation and Transmission Spectrum
% This program uses Transfer Matrix Method
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
clc
clear all
close all
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Initialization of parameters
dB=7778e-9;                         % Thickness of Second Layer
dA=0.1*5303e-9;                         % Thickness of First Layer
nA=3.3;                             % Refractive index of First Layer
nB=2.25;                            % Refraive index of Second Layer
n0=1;                               % Refractive index of Incident and transmission medium
c=3e8;                              % Velocity of Light (m/s)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Frequency Range
frequncy=0.1:0.1:40;                    % Frequncy in THz
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
for loop=1:length(frequncy)
    % angular frequency
    w=2*pi*frequncy*1e12;
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    d=dA+dB;    DA=((w(loop))/c)*dA*nA;    DB=((w(loop))/c)*dB*nB;
    m11=cos(DA);   m12=-1i*sin(DA)/nA;   m21=-1i*nA*sin(DA);   m22=cos(DA);   mA=[m11 m12; m21 m22];
    l11=cos(DB);   l12=-1i*sin(DB)/nB;   l21=-1i*nB*sin(DB);   l22=cos(DB);   mB=[l11 l12; l21 l22];
      m=(mA*mB);
    RHS(loop)=(m(1,1)+m(2,2))/2;
end
plot(frequncy,RHS)
hold on
yline(-1)
yline(1)
hold off
4 Comments
  Dyuman Joshi
      
      
 on 29 Aug 2023
				"Here i'm using for loop."
That does not make any difference. 
You can still apply the method suggested in that answer.
"I'm confusing to use that previous method. "
How are you confused?
  Bruno Luong
      
      
 on 29 Aug 2023
				If you are confused the iI recommend asking specific question focusing on what confuses you, instead of ignoring it and ask again the same question. You'll end of with a tone of things that confuse you even after 10 year of MATLAB practicing.
Accepted Answer
  Bruno Luong
      
      
 on 29 Aug 2023
        OK I'll be nice with you because of your laziness I just copy the solution of another thread here; wirh a little adaptation to sort them 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% This is a program for 1-D Photonic crystal
% For Photonic band structure calculation and Transmission Spectrum
% This program uses Transfer Matrix Method
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%clc
%clear all <====== BAD programming practice
%close all
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Initialization of parameters
dB=7778e-9;                         % Thickness of Second Layer
dA=0.1*5303e-9;                         % Thickness of First Layer
nA=3.3;                             % Refractive index of First Layer
nB=2.25;                            % Refraive index of Second Layer
n0=1;                               % Refractive index of Incident and transmission medium
c=3e8;                              % Velocity of Light (m/s)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Frequency Range
frequncy=0.1:0.1:40;                    % Frequncy in THz
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
for loop=1:length(frequncy)
    % angular frequency
    w=2*pi*frequncy*1e12;
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    d=dA+dB;    DA=((w(loop))/c)*dA*nA;    DB=((w(loop))/c)*dB*nB;
    m11=cos(DA);   m12=-1i*sin(DA)/nA;   m21=-1i*nA*sin(DA);   m22=cos(DA);   mA=[m11 m12; m21 m22];
    l11=cos(DB);   l12=-1i*sin(DB)/nB;   l21=-1i*nB*sin(DB);   l22=cos(DB);   mB=[l11 l12; l21 l22];
      m=(mA*mB);
    RHS(loop)=(m(1,1)+m(2,2))/2;
end
allyx = [-1 1];
xintersect = [];
yintersect = [];
for yx = allyx
    ys = RHS-yx;
    i = find(ys(1:end-1).*ys(2:end) <= 0);
    % linear interpolation
    i1 = i; ss1 = ys(i1);
    i2 = i1 + 1; ss2 = ys(i2);
    w = ss2./(ss2-ss1);
    xx = w.*frequncy(i1) + (1-w).*frequncy(i2);
    yy = yx+zeros(size(xx));
    xintersect = [xintersect, xx];
    yintersect = [yintersect, yy];
end
[xintersect, is] = sort(xintersect);
yintersect= yintersect(is);
figure
plot(frequncy,RHS)
hold on
yline(allyx)
plot(xintersect, yintersect, 'or')
hold off
0 Comments
More Answers (1)
  Rawad SADAH
 on 29 Aug 2023
        Try this code
% Finding intersection points
threshold = 0.01;  % Adjust this threshold based on your data and precision needs
intersection_points = [];
for loop = 2:length(frequncy)
    if (RHS(loop) > -1 && RHS(loop - 1) < -1) || (RHS(loop) < 1 && RHS(loop - 1) > 1)
        % Perform linear interpolation to find the approximate intersection point
        t = (1 - RHS(loop - 1)) / (RHS(loop) - RHS(loop - 1));
        intersection_freq = frequncy(loop - 1) + t * (frequncy(loop) - frequncy(loop - 1));
        % Add the intersection point to the array
        intersection_points = [intersection_points, intersection_freq];
    end
end
% Sort the intersection points in ascending order
sorted_intersection_points = sort(intersection_points);
disp("Intersection Points:");
disp(sorted_intersection_points);
% Plotting the intersection points on the blue curve
hold on;
plot(sorted_intersection_points, ones(size(sorted_intersection_points)), 'ro');
hold off;
See Also
Categories
				Find more on Photonics 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!
