function [t] = bsp_20_fkt(x,y,z,r,t1)
% Determines the next t>s such that the sum of radii r(s)+r(t) < |P(s)-P(t)|
% INPUT: x   : x(t) symbolic function for x components
%        y   : y(t) symbolic function for y components
%        z   : z(t) symbolic function for z components
%        r   : r(t) symbolic function for radius of ball at P(t)
%        t1  : (numeric/symbolic) parameter t for given first point
% OUTPUT:  t : next (numeric) parameter t>s such that   r(t1)+r(t) < |P(t1)-P(t)| holds

syms s t

%  Quadrat des Eukldischen Abstands zwischen P(s) und P(t)
d_long(s,t) = (x(t)-x(s))^2 +(y(t)-y(s))^2 +(z(t)-z(s))^2;
%  Ausmultiplizieren und vereinfachen (vielleicht hilft das)
d(s,t) = simplify(expand(d_long));

% Bedingung d(s,t) >= (r(s)+r(t))^2 
% <==> g(s,t) := d(s,t) - (r(s)+r(t))^2  >=0 
g(s,t) = d - (r(s)+r(t))^2;

% Loese g(t1,t)==0
magic_param = 0.1;
interval_right = t1*[1    ,1+magic_param];
tsol = vpasolve(g(t1,t)==0, t, interval_right);

% Nimm das t aus tsol mit t-t1 = min(tsol-t1)
[~,idx] = min(tsol-t1);  % index of solution t with minimal distance to t1
t = tsol(idx);

t = double(t);           % symbolic --> numeric

end
