polyxpoly to find intersection

Question

0 个投票

i have two curves as given in the code 'p' and 'q' for diffeent values of 'x' the curves intersect at different points. Can anyone tell me why is my code not running?

clc
clf
b=linspace(0,1);
x=8;
p= x.*(sqrt(1-b)).*(((besselj(1,x.*sqrt(1-b)))./(besselj(0,x.*sqrt(1-b)))));
q = x.*(sqrt(b)).*((besselk(1,x.*sqrt(b)))./(besselk(0,x.*sqrt(b))));
plot(b,p)
hold on
plot(b,q)
[xx,yy] = polyxpoly(b,p,b,q)

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

请先登录，再进行评论。

请先登录，再回答此问题。

Follow Question

Answer 1

KSSV 2020-4-24

编辑：KSSV 2020-4-25

在 MATLAB Online 中打开

0 个投票

You may consider this https://in.mathworks.com/matlabcentral/fileexchange/22441-curve-intersections?focused=5165138&tab=function to calculate the intersection points.

Save the above function in your working folder and use the below code.

clc
clf
b=linspace(0,1);
x=8;
p= x.*(sqrt(1-b)).*(((besselj(1,x.*sqrt(1-b)))./(besselj(0,x.*sqrt(1-b)))));
q = x.*(sqrt(b)).*((besselk(1,x.*sqrt(b)))./(besselk(0,x.*sqrt(b))));
L1 = [b ;p] ;
L2 = [b ;q] ; 
P = InterX(L1,L2) ; % get intersection points using the function 
plot(b,p)
hold on
plot(b,q)
plot(P(1,:),P(2,:),'*r')

13 个评论
显示 11更早的评论隐藏 11更早的评论

vipul kumar 2020-4-25

clc

clf

b=linspace(0,1);

x=2;

p= x.*(sqrt(1-b)).*(((besselj(1,x.*sqrt(1-b)))./(besselj(0,x.*sqrt(1-b)))));

q = x.*(sqrt(b)).*((besselk(1,x.*sqrt(b)))./(besselk(0,x.*sqrt(b))));

plot(b,p,b,q,P(1,:),P(2,:),'ro')

function P = InterX(L1,varargin)

%INTERX Intersection of curves

% P = INTERX(L1,L2) returns the intersection points of two curves L1

% and L2. The curves L1,L2 can be either closed or open and are described

% by two-row-matrices, where each row contains its x- and y- coordinates.

% The intersection of groups of curves (e.g. contour lines, multiply

% connected regions etc) can also be computed by separating them with a

% column of NaNs as for example

%

% L = [x11 x12 x13 ... NaN x21 x22 x23 ...;

% y11 y12 y13 ... NaN y21 y22 y23 ...]

%

% P has the same structure as L1 and L2, and its rows correspond to the

% x- and y- coordinates of the intersection points of L1 and L2. If no

% intersections are found, the returned P is empty.

%

% P = INTERX(L1) returns the self-intersection points of L1. To keep

% the code simple, the points at which the curve is tangent to itself are

% not included. P = INTERX(L1,L1) returns all the points of the curve

% together with any self-intersection points.

%

% Example:

% t = linspace(0,2*pi);

% r1 = sin(4*t)+2; x1 = r1.*cos(t); y1 = r1.*sin(t);

% r2 = sin(8*t)+2; x2 = r2.*cos(t); y2 = r2.*sin(t);

% P = InterX([x1;y1],[x2;y2]);

% plot(x1,y1,x2,y2,P(1,:),P(2,:),'ro')

% Author : NS

% Version: 3.0, 21 Sept. 2010

% Two words about the algorithm: Most of the code is self-explanatory.

% The only trick lies in the calculation of C1 and C2. To be brief, this

% is essentially the two-dimensional analog of the condition that needs

% to be satisfied by a function F(x) that has a zero in the interval

% [a,b], namely

% F(a)*F(b) <= 0

% C1 and C2 exactly do this for each segment of curves 1 and 2

% respectively. If this condition is satisfied simultaneously for two

% segments then we know that they will cross at some point.

% Each factor of the 'C' arrays is essentially a matrix containing

% the numerators of the signed distances between points of one curve

% and line segments of the other.

%...Argument checks and assignment of L2

error(nargchk(1,2,nargin));

if nargin == 1,

L2 = L1; hF = @lt; %...Avoid the inclusion of common points

else

L2 = varargin{1}; hF = @le;

end

%...Preliminary stuff

x1 = L1(1,:)'; x2 = L2(1,:);

y1 = L1(2,:)'; y2 = L2(2,:);

dx1 = diff(x1); dy1 = diff(y1);

dx2 = diff(x2); dy2 = diff(y2);

%...Determine 'signed distances'

S1 = dx1.*y1(1:end-1) - dy1.*x1(1:end-1);

S2 = dx2.*y2(1:end-1) - dy2.*x2(1:end-1);

C1 = feval(hF,D(bsxfun(@times,dx1,y2)-bsxfun(@times,dy1,x2),S1),0);

C2 = feval(hF,D((bsxfun(@times,y1,dx2)-bsxfun(@times,x1,dy2))',S2'),0)';

%...Obtain the segments where an intersection is expected

[i,j] = find(C1 & C2);

if isempty(i),P = zeros(2,0);return; end;

%...Transpose and prepare for output

i=i'; dx2=dx2'; dy2=dy2'; S2 = S2';

L = dy2(j).*dx1(i) - dy1(i).*dx2(j);

i = i(L~=0); j=j(L~=0); L=L(L~=0); %...Avoid divisions by 0

%...Solve system of eqs to get the common points

P = unique([dx2(j).*S1(i) - dx1(i).*S2(j), ...

dy2(j).*S1(i) - dy1(i).*S2(j)]./[L L],'rows')';

function u = D(x,y)

u = bsxfun(@minus,x(:,1:end-1),y).*bsxfun(@minus,x(:,2:end),y);

end

vipul kumar 2020-4-25

ohh, you edited the parent question. My bad. Got it bruh, thanks a lot.

vipul kumar 2020-4-25

There is one more thing, if you lookt at it. Since the function goes to infinity so the intersection point generated in the graph are more i.e. if you look from the right of the graph, only alternate points are taken into account. First relevent point is the rightmost and then alternate are to be taken.

I verified this with the data given but can not find the exact explanation why does it give one extra point after one relevant one?

请先登录，再进行评论。

polyxpoly to find intersection

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

采纳的回答

13 个评论
显示 11更早的评论隐藏 11更早的评论

更多回答（0 个）

类别

标签

Community Treasure Hunt

polyxpoly to find intersection

0 个评论 显示 -2更早的评论 隐藏 -2更早的评论

采纳的回答

13 个评论 显示 11更早的评论 隐藏 11更早的评论

更多回答（0 个）

类别

标签

另请参阅

Community Treasure Hunt

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

13 个评论
显示 11更早的评论隐藏 11更早的评论