First the code, some explanation later. The process is somewhat time-consuming, so the code was slighty optimized. Now it resembles a sharp selection, where very few candidates will pass all tests.
type
myType=Integer; {remember to change here to your own type,
an integers are good for an educational purpose}
point_2d=record xx,yy :myType end;
vector_2d=record vx,vy :myType end;
function bezier_has_a_loop(p0,p1,p2,p3:point_2d):boolean;
{-------------------------}
function f_vec_from_2_pts(pa,pb:point_2d):vector_2d;
begin
with Result do begin
vx:=pb.xx - pa.xx;
vy:=pb.yy - pa.yy end end;
{-------------------------}
function f_cross_prod(va,vb : vector_2d):myType;
begin
Result := va.vx*vb.vy - va.vy*vb.vx end;
{-------------------------}
function f_mult(m:myType; v:vector_2d):vector_2d;
begin
with Result do begin
vx := m*v.vx;
vy := m*v.vy end end;
{-------------------------}
function f_sum(va,vb : vector_2d):vector_2d;
begin
with Result do begin
vx := va.vx+vb.vx;
vy := va.vy+vb.vy end end;
{-------------------------}
function f_rotate(v, rotor : vector_2d):vector_2d;
var
m_sin,m_cos:myType; {only for a readability purpose}
begin
m_cos:=rotor.vx {
/sqrt(sqr(rotor.xx)+sqr(rotor.yy)) - unnecessary };
m_sin:=rotor.vy {
/sqrt(sqr(rotor.xx)+sqr(rotor.yy)) - unnecessary };
with Result do begin
vx := -v.vx * m_cos - v.vy *m_sin;
vy := v.vx * m_sin - v.vy *m_cos end end;
{-------------------------}
var
a,b,c, c1,c2,c3 :vector_2d;
ab,ac,bc:myType;
bb,s1,s2,delta,shift,t_1_2 : Double;
begin
a:=f_vec_from_2_pts(p0,p1);
b:=f_vec_from_2_pts(p1,p2);
c:=f_vec_from_2_pts(p2,p3);
ab:=f_cross_prod(a,b); {on the planar,
myType for a cross product is just fine}}
ac:=f_cross_prod(a,c);
bc:=f_cross_prod(b,c);
{==1== No inflection point(s) or cusp allowed}
if ac*ac >= 4*ab*bc then begin Result:=False; exit end;
c3:= f_sum( f_sum(a,c) , f_mult( -2,b ) );
if c3.vy<>0 then begin {Is the bag's handle horizontal?}
a := f_rotate(a,c3);
b := f_rotate(b,c3);
c := f_rotate(c,c3);
c3:= f_sum ( f_sum(a,c) , f_mult(-2,b) );
{Now the handle is forced to be horizontal.}
end;
c1:= f_mult ( 3,a );
c2:= f_sum ( f_mult(3,b) , f_mult(-3,a) );
{ Following 4 restrictions comes from a single caveats for a roots:
0<= t1<t2<=1}
bb:= -c1.vy / c2.vy; { always c2.vy<>0 }
{==2== A central point (t1+t2)/2 outside the limits}
if (bb<0) or (bb>2) then begin Result:=False; exit end;
s1:= c1.vx/c3.vx; { always c3.vx<>0 }
s2:= c2.vx/c3.vx;
delta := -bb*(4*s2+3*bb)-4*s1;
{==3== delta is to big}
if delta>1 then begin Result:=False; exit end;
shift:=sqrt(delta); { always delta>0 }
t_1_2:=bb-shift; {for readability purpose only,
one can omit this and write below: if shift>bb }
{==4== t1 is to low }
if t_1_2<0 then begin Result:=False; exit end;
t_1_2:=bb+shift; { no /2 here,beacause of *2 below}
{==5== t2 is to high}
if t_1_2>2 then Result:=False
else Result:=True { TA DA! Thank you for your patience. }
end;
Now some theory.
- Let take into account 4 points P0, P1, P2, P3. These points define (almost always) a cubic Bézier curve.
- Let establish the point H1 such as the vector P1_H1 is half of the vector P0_P1.
- Let establish the point H2 such as the vector P2_H2 is half of the vector P3_P2.
- At the end, let create the vetor h = H1_H2. (Personally I call this vector "the handle of the bag", guess why.)
Bezier curve and its handle
- No surprise here, when you start an isotropic scaling or a roation of P0, P1, P2, P3, then vector h will transform accordingly.
- Maybe for someone this will be a surprise, the vector h = [_h_x,_h_y] is proportional to the vector [c3x,c3y] created from the highest coefficients of the algebraic form of the cubic Bézier curve. The proportionality coefficient is (-1/2).
(In fact, when the points H1=H2 coincide, the h vector vanishes, c3x=c3y=0, thus the cubic Bézier curve reduces at least into the quadratic Bézier curve, created from P0,H1,P3 points.)
And now the clue: The right rotation can always turn the vector h horizontally (or verically), and this rotation will reduce c3y ( or c3x ) into null, however preserving the loop. In turn, one can reduce the stated problem into a trivial root-finding of a quadratic equation. (I suppose this is a decisive hint to find a solution.)
To measure the loop, I suggest take into account the variable delta from the code above.
0 < delta <= 1
When delta -> 0, the loop vanishes.
When delta -> 1, the loop becomes really pompous.
I'm not quite happy with this proposal, but it is still better than nothing.