find the intersection points?

Need help figuring out how to find the intersection points of the plot.

a0 = 1;

Ab = 5.21e6;

kenex = 1.05;

f = 0.1;

Lc = 1700;

T = 12.42*3600;

g = 9.81;

sigma = 2*pi/T;

Ck = 0.86;

% Hydraulic Curve

Ac = [1:1:2e3];

Rc = sqrt(Ac/20);

F = kenex + f*Lc./(4*Rc);

beta = F*Ab./(2*Lc*Ac)*a0;

alpha = (Lc*Ab./(g*Ac)).^0.5*sigma;

mu = 16*alpha.^2.*beta/(3*pi);

um_prime = ((((1-alpha.^2).^4 + mu.^2).^0.5 - (1-alpha.^2).^2)./(0.5*mu.^2)).^0.5;

um = um_prime*2*pi*a0/T.*Ab./Ac;

% Sedimentary Curve

P = (Ac/1.58e-4).^(1/0.95);

um2 = P*pi*Ck./(T*Ac);

figure(1);

semilogx(Ac,um)

xlabel('Channel Area (m^2)')

ylabel('u_m (m/s)')

grid on

hold on

semilogx(Ac,um2);

legend('Hydraulic Curve','Sedimentary Curve')

Answer 1

[Star Strider]

Add these ines:

curv_dif = um - um2;                                        % Subtract Curves
idx = find(diff(sign(curv_dif)));                           % Indicex Of Sign Changes
for k = 1:numel(idx)
    idxrng = [-1 0 1]+idx(k);
    xint(k) = interp1(curv_dif(idxrng), Ac(idxrng), 0);     % ‘x’ At Intersection
    yint(k) = interp1(Ac, um, xint(k));                     % ‘y’ At Intersection
end

so the complete revised code is now:

a0 = 1;
Ab = 5.21e6;
kenex = 1.05;
f = 0.1;
Lc = 1700;
T = 12.42*3600;
g = 9.81;
sigma = 2*pi/T;
Ck = 0.86;
% Hydraulic Curve
Ac = [1:1:2e3];
Rc = sqrt(Ac/20);
F = kenex + f*Lc./(4*Rc);
beta = F*Ab./(2*Lc*Ac)*a0;
alpha = (Lc*Ab./(g*Ac)).^0.5*sigma;
mu = 16*alpha.^2.*beta/(3*pi);
um_prime = ((((1-alpha.^2).^4 + mu.^2).^0.5 - (1-alpha.^2).^2)./(0.5*mu.^2)).^0.5;
um = um_prime*2*pi*a0/T.*Ab./Ac;
% Sedimentary Curve
P = (Ac/1.58e-4).^(1/0.95);
um2 = P*pi*Ck./(T*Ac);
figure(1);
semilogx(Ac,um)
xlabel('Channel Area (m^2)')
ylabel('u_m (m/s)')
grid on
hold on
semilogx(Ac,um2);

curv_dif = um - um2;                                        % Subtract Curves
idx = find(diff(sign(curv_dif)));                           % Indicex Of Sign Changes
for k = 1:numel(idx)
    idxrng = [-1 0 1]+idx(k);
    xint(k) = interp1(curv_dif(idxrng), Ac(idxrng), 0);     % ‘x’ At Intersection
    yint(k) = interp1(Ac, um, xint(k));                     % ‘y’ At Intersection
end
plot(xint, yint, 'xr')
hold off
legend('Hydraulic Curve','Sedimentary Curve','Intersections','Location','best')

EDIT — (12 Jun 2020 at 21:24)

Added plot figure —

.