Unsatisfying Point Type Change From Corner To Curve

[linusromer]

The function SPChangePointType() in splinechar.c has a unsatisfying behaviour that I will try to describe in two GIFs:

1. If you change a single-linked corner point to a curve point, it will head to the next on-curve-node in the distance of the nice proportion = .39:
2. If you change a double-linked corner point, the behaviour is basically the same but the direction is defined by the vector from the next and the previous on-curve-node:

This problem occurs both with the 2020-11-07 release as well as the current main branch of FontForge.

I have tried to fix the problem by changing line 699 of splinechar.c to something like

if( pointtype==pt_curve ) { 
    if ( sp->pointtype==pt_corner ) { 
        if ( sp->prev!=NULL && sp->next!=NULL && prevlen!=0 && nextlen!=0 ) {
            unitnext.x = .5*(unitnext.x-unitprev.x); /* average direction */
            unitnext.y = .5*(unitnext.y-unitprev.y); 
            sp->nextcp.x = sp->me.x + unitnext.x*nextlen;
            sp->nextcp.y = sp->me.y + unitnext.y*nextlen;
            sp->prevcp.x = sp->me.x - unitnext.x*prevlen; /* prev is antidirectional to next */
            sp->prevcp.y = sp->me.y - unitnext.y*prevlen;
        }
        if ( ( sp->prev==NULL && sp->next!=NULL && nextlen!=0 ) || (sp->next==NULL && sp->prev!=NULL && prevlen!=0) ) 
            makedflt = false; /* just do nothing */
    } 
    else makedflt = true; /* original behaviour */
}   

to change this behaviour but it did not work as desired. Maybe we should consider changing the function SplineRefigureFixup() as well? If someone has a good hint, I will try to figure out a patch and make a pull request.

[skef]

Huh, I thought I had filed an issue on this a few years ago but apparently not.

It depends on what you mean by "did not work as desired". Is it not taking effect or is the change just not desirable?

I would say:

1. When converting from corner to curve, if the point only has one defined side the control point should stay in the same position and the "hanging" control point, if any, should keep the same distance while moving into alignment with it.
2. When converting from corner to curve and both sides are defined, you should consider averaging the control point angles (that is, changing the angle of each an equal distance towards colinear) and then fitting each side at that slope to the previous spline. That way you get as close to the previous curve as possible with the new slope. (Your more accurate brute-force fitter would be good for this.) Or at least that's a pretty good starting heuristic.
3. Another potential option would be that if one control point is selected, keep that one as it is and move the other to alignment with it. That way you have more control/options about what happens.

(Hmm -- should one select the one that won't change or the one that will change? Hmm.)

Tangent point conversions are messed up as well. I would favor:

1. If there is only one (non-colocated) control point, keep that side.
2. Otherwise if one control point is selected, and keep that one while changing the other side to a line (or vice versa).
3. Otherwise raise a dialog saying one control point has to be selected.

That way people who change point types this way gradually learn how to select the side (although other reviewers may have different opinions). The dialog might be overkill when making a "mass change", although the current alternative is a terrible one. (I understand the rationale for why the tangent points behave as they do, but converting to one with both control points extended is bad.)

[linusromer]

For the sake of simplicity and consistency (with "mass change") I would not make it dependable from selected control points, therefore I would keep your suggestions 1 and 2 (but not 3) for the conversion to smooth points, the average direction should work as I describe before with

unitnext.x = .5*(unitnext.x-unitprev.x); 
unitnext.y = .5*(unitnext.y-unitprev.y);

Keeping this logic for the conversion to tangent points, I would suggest keeping your suggestion1 but for the generic case calculate a measure for the "linearity" of the two adjacent splines and change the more linear spline to a line.

[skef]

All this sounds reasonable.

The very simple measures of linearity I can think of have foibles. A reasonable not-quite-simple measure is the max distance from the line formed between the start and end divided by the distance between the start and end. If there's not a handy way of calculating that directly you could rotate to make the line vertical and then use the "accurate bounding box" routine. Or maybe there's something better. (One rarely does changes hundreds or even dozens of points so I think it's OK if the heuristics are a little "expensive".)

[linusromer]

A fast computable measure for linearity could be the following:

1. Norm the spline: Scale and rotate the spline such that one end lies on (0,0) and the other on (1,0). We denote the signed angle of the start handle to the x-axis with α and the signed angle of the end handle to the x-axis with β. The lengths of the handles are denoted by a resp. b.
2. Calculate the second moment about the x-axis (along the y-axis) of this normed spline. This is the same as the volume of revolution about the x-axis except for the factor 1/2 rather than π:

   1/2 \int y^2 dA = 1/2 \int y^2 dA/dt dt
   = 1/2*integrate( (3*(1-t)^2*t*a*sin(α)+3*(1-t)*t^2*sin(β)*b)^2*(-9*(t-1)*t*(b*sin(β)*t-a*sin(α)*t+a*sin(α))*(3*b*cos(β)*t^2+3*a*cos(α)*t^2-2*t^2-2*b*cos(β)*t-4*a*cos(α)*t+2*t+a*cos(α))),t,0,1);
   = -(9*((9*a*cos(α)-14)*b^3*sin(β)^3+(9*a*sin(α)*b^3*cos(β)+(15*a^2*cos(α)-30*a)*sin(α)*b^2)*sin(β)^2+(15*a^2*sin(α)^2*b^2*cos(β)+(9*a^3*cos(α)-30*a^2)*sin(α)^2*b)*sin(β)+9*a^3*sin(α)^3*b*cos(β)-14*a^3*sin(α)^3))/6160

   Of course, we can drop the division by 6160 and the formula can somewhat be optimized by partial factorization. I think this is a super fast and super exact sum of square error for linearity.

[skef]

Is that really faster than a bounding box in 1d? (That is, rotating to vertical wrt on-curve points, finding the roots of the x-dimension derivative, and testing them for left or right positions (adding if needed), and then dividing to get the ratio?) Or is the metric better in some way?

[linusromer]

Speaking about speed, the second moment is probably comparable to the bounding box in 1d. However, the properties of the second moment might indeed be more desirable. Consider the two following splines: The upper spline would have the same bounding box in 1d as the lower spline but the bigger second moment because more points lie far from the chord. So the lower spline is more linear.

If the bounding box is calculated using the derivative, cusps need special treatment since the derivative is not defined there. But cusps are no problem for the second moment.

[linusromer]

So, I have invested some time on this again. With replacing the line 699 of splinechar.c by

if( pointtype==pt_curve ) { /* added by Linus Romer */
        if ( sp->pointtype==pt_corner ) { 
            makedflt = false; /* just do nothing (later )*/ 
            if ( sp->prev!=NULL && sp->next!=NULL) {
                if ( prevlen!=0 && nextlen!=0 ) { /* take the average direction */
                    sp->nextcp.x = sp->me.x + .5*(unitnext.x-unitprev.x)*nextlen; 
                    sp->nextcp.y = sp->me.y + .5*(unitnext.y-unitprev.y)*nextlen;
                    sp->prevcp.x = sp->me.x - .5*(unitnext.x-unitprev.x)*prevlen; /* prev is antidirectional to next */
                    sp->prevcp.y = sp->me.y - .5*(unitnext.y-unitprev.y)*prevlen;
                }
                else if ( prevlen!=0 ) { /* and therefore nextlen==0 */
                    unitprev.x /= prevlen; unitprev.y /= prevlen; /* make it really a unit vector */
                    nextlen = BPNorm(BPSub(sp->next->to->me,sp->me));
                    sp->nextcp.x = sp->me.x - unitprev.x*.39*nextlen;
                    sp->nextcp.y = sp->me.y - unitprev.y*.39*nextlen;
                    nextlen = 0; /* just to keep the former information */
                }
                else if ( nextlen!=0 ) { /* and therefore prevlen==0 */
                    unitnext.x /= nextlen; unitnext.y /= nextlen; /* make it really a unit vector */
                    prevlen = BPNorm(BPSub(sp->prev->from->me,sp->me));
                    sp->prevcp.x = sp->me.x - unitnext.x*.39*prevlen;
                    sp->prevcp.y = sp->me.y - unitnext.y*.39*prevlen;
                    prevlen = 0; /* just to keep the former information */
                }
                else makedflt = true;               
            }
        } 
        else makedflt = true; /* original behaviour */
    }   

and including utanvec.h, I have the desired property of the conversion from corner to curve, see splinechar.c.txt.
But then the conversion back does not work anymore. Probably because I did not call SplineRefigureFixup(), where some sanity checks are done.

Can anybody give me a hint? Or should I make a pull request in this very early stage?

[linusromer]

Stupid me! I have overseen sp->pointtype = pointtype; (strange that it still worked partially). I have now restored int oldpointtype = sp->pointtype; and I am again on my way. I will keep you informed.

[linusromer]

After working longer on this, I have changed my mind on two things:

• When I convert a corner point to a smooth point, the handles of the corner are already visually parallel. Then the simple "averaging angle algorithm" with projection of the handles to the new direction is very close to what I expect - there is no need for a fitting algorithm. If the handles of the corner are not parallel at all I don't expect FontForge to search for a super close solution because I would fiddle the handles afterwards anyway.
• The second moment is a cool thing but my calculation of 1/2 \int y^2 dA suffers from the fact, that my dA is signed. Hence, the formula given above is useless for a measure of linearity. I make use of skef's suggestion: "A reasonable not-quite-simple measure is the max distance from the line formed between the start and end divided by the distance between the start and end".

The code for the measure of linearity is not super short but fast (at least faster than a linear transformation and a bounding box I guess):

/* This is a fast computable measure for linearity. */
/* First, we norm the spline: Scale and rotate the spline such that */
/* one end lies on (0,0) and the other on (1,0). We denote the signed */
/* angle of the start handle to the x-axis with alpha and the signed */
/* angle of the end handle to the x-axis with beta. The lengths of */
/* the handles are denoted by a resp. b. */
static bigreal Linearity(Spline *s) { /* added by Linus Romer */
    if ( s->islinear ) return 0;  
    BasePoint ftunit = BPSub(s->to->me, s->from->me);
    bigreal ftlen = BPNorm(ftunit);
    if ( ftlen==0 ) return -1; /* flag for error, no norming possible */
    BasePoint fromunit = BPSub(s->from->nextcp, s->from->me);
    BasePoint tounit = BPSub(s->to->prevcp, s->to->me);
    bigreal a = BPNorm(fromunit)/ftlen;
    bigreal b = BPNorm(tounit)/ftlen;
    ftunit = BPScale(ftunit, 1/ftlen);
    if ( a==0 && b==0 ) return 0;
    if ( a==0 ) fromunit = BPSub(s->to->prevcp, s->from->me);
    if ( b==0 ) tounit = BPSub(s->from->nextcp, s->to->me);
    fromunit = NormVec(fromunit);
    tounit = NormVec(tounit);
    bigreal sinalpha = BPCross(ftunit, fromunit); 
    bigreal sinbeta = BPCross(ftunit, tounit); 
    /* just run t from 0 to 1 and get the largest y(t) */
    /* y(t) = 3*(1-t)^2*t*a*sin(alpha)+3*(1-t)*t^2*sin(beta)*b */
    /* which is proportional to (1-t)*t*((1-t)*a*sin(alpha)+t*sin(beta)*b) */
    bigreal t,y,ymax;
    ymax = 0;
    for (int i=0; i<99; ++i) {
        t = (i+1)/100.0;
        y = fabs((1-t)*t*((1-t)*a*sinalpha+t*sinbeta*b));
        if ( y > ymax ) ymax = y;
    }
    return ymax;
}

For quadratic splines one can just take the distance of the control point to the chord (not yet implemented).

Then I have replaced the whole block of else if ( pointtype==pt_tangent ) {...} by:

        } else if ( pointtype==pt_tangent ) {
    /* Added by Linus Romer: Check if previous or next spline is more linear and make it linear */
    if ( sp->next!=NULL && sp->prev!=NULL ) {
        bigreal prevlinearity = Linearity(sp->prev);
        bigreal nextlinearity = Linearity(sp->next);
        if ( prevlinearity>=0 && nextlinearity>=0 ) {
            if ( nextlinearity >= prevlinearity ) { /* make prev linear */
                sp->prev->from->nextcp = sp->prev->from->me;
                sp->prevcp = sp->me;
                //sp->noprevcp = true;
                //sp->prev->from->nonextcp = true;
                unitnext = NormVec(BPSub(sp->me,sp->prev->from->me));
                sp->nextcp = BPAdd(sp->me, BPScale(unitnext,fabs(BPDot(BPSub(sp->nextcp, sp->me), unitnext))));
            } else { /* make next linear */
                sp->next->to->prevcp = sp->next->to->me;
                sp->nextcp = sp->me;
                //sp->nonextcp = true;
                //sp->next->to->noprevcp = true;
                unitprev = NormVec(BPSub(sp->me,sp->next->to->me));
                sp->prevcp = BPAdd(sp->me, BPScale(unitprev,fabs(BPDot(BPSub(sp->prevcp, sp->me), unitprev))));
            }
        } /* else do nothing - this would not make any sense */
    }

And it seems to work okay as far as I have tested it.

The replacment for if( pointtype==pt_curve ) makedflt = true; now looks like the following:

       if( pointtype==pt_curve ) { /* added by Linus Romer */
        if ( oldpointtype==pt_corner ) { 
            makedflt = false; /* just do nothing (later )*/ 
            if ( sp->prev!=NULL && sp->next!=NULL) {
                if ( prevlen!=0 && nextlen!=0 ) { /* take the average direction */
                    sp->nextcp = BPAdd(sp->me, BPScale(BPSub(unitnext, unitprev), .5*nextlen));
                    sp->prevcp = BPSub(sp->me, BPScale(BPSub(unitnext, unitprev), .5*prevlen));
                }
                else if ( prevlen!=0 ) /* and therefore nextlen==0 */
                    sp->nextcp = BPSub(sp->me, BPScale(NormVec(unitprev), .39*BPNorm(BPSub(sp->next->to->me, sp->me))));
                else if ( nextlen!=0 ) /* and therefore prevlen==0 */
                    sp->prevcp = BPSub(sp->me, BPScale(NormVec(unitnext), .39*BPNorm(BPSub(sp->prev->from->me, sp->me))));
                else makedflt = true;
            }
        } 
        else makedflt = true; /* original behaviour */
    }   

At the moment, the .39 is hardwired in my code. I think, one should move #define NICE_PROPORTION .39 from splineutil.c to splineutil.h.

[skef]

For the linearity measure: After you normalize instead of sampling why not just calculate the y-extrema and pick the one that's farthest away? Works for standard curves and ones with a loop, not sure about cusps but not sure that matters much.

One would use the existing 1-D SplineFindExtrema for this if one were so inclined. Might be better to nromalize to a length of 100 (or whatever) instead of one to better match the floating point heuristics in the code.

[skef]

(If you want a code "model" SplineSolveForUTanVec() in utanvec.c does something quite similar.)

[linusromer]

You are probably right about cusps: The probability for random handles to generate cusps is 0 and I have never provoked cusps when drawing glyphs and the same is probably true for all users of FontForge.

I had a look at SplineSolveForUTanVec() and I hope I have understood everything correctly:

• The zero checks ZCHECK are not necessary in my case, I presume?
• I think #define ROTY(p, ut) (-(ut).y*(p).x + (ut).x*(p).y) can be replaced by BPCross()

My implementation for cubic splines now looks as follows:

static bigreal Linearity(Spline *s) { /* added by Linus Romer */
    if ( s->islinear ) return 0;  
    BasePoint ftunit = BPSub(s->to->me, s->from->me);
    bigreal ftlen = BPNorm(ftunit);
    if ( ftlen==0 ) return -1; /* flag for error, no norming possible */
    ftunit = BPScale(ftunit, 1/ftlen);
    BasePoint nextcpscaled = BPScale(BPSub(s->from->nextcp, s->from->me),100./ftlen);
    BasePoint prevcpscaled = BPScale(BPSub(s->to->prevcp, s->from->me),100./ftlen);
    /* we just look a the bezier coefficients in y dimension: */
    Spline1D ys1d;
    ys1d.a = 0;
    ys1d.b = BPCross(ftunit, nextcpscaled); /* rotate it to horizontal*/
    ys1d.c = BPCross(ftunit, prevcpscaled); /* rotate it to horizontal*/
    ys1d.d = 0;
    extended te1, te2;
    SplineFindExtrema(&ys1d, &te1, &te2);
    if ( te1==-1 && te2==-1 ) return 0;
    return fmax(fabs(te1*(1-te1)*((1-te1)*ys1d.b+te1*ys1d.c)), fabs(te2*(1-te2)*((1-te2)*ys1d.b+te2*ys1d.c)));
}

I have compiled it and the tests showed that ys1d.b and ys1d.c are the normed y coordinates of the handles as I have intended. But then I am confused because SplineFindExtrema() produces only te1==-1 and te2==-1. Did I misunderstand SplineFindExtrema()? ys1d should be filled with the y coordinates of the spline that is being examined, true?

[skef]

The zero checks ZCHECK are not necessary in my case, I presume?

I wouldn't think they are.

I think #define ROTY(p, ut) (-(ut).y(p).x + (ut).x(p).y) can be replaced by BPCross()

If it's the same. I wrote this fairly early in the expand stroke replacement project so it could easily be redundant code.

The problem here is the Spline1D values. You can Ctrl-F in https://fontforge.org/docs/techref/splinefont.html for more details but basically a, b, c, and d need to be the expanded coefficients of the cubic equation, not the coordinates. So while you don't need the ZCHECKs you do need some of those calculations.

This is what SplineRefigure3 does in the general case, as in:

        xsp->c = 3*(from->nextcp.x-from->me.x);
        ysp->c = 3*(from->nextcp.y-from->me.y);
        xsp->b = 3*(to->prevcp.x-from->nextcp.x)-xsp->c;
        ysp->b = 3*(to->prevcp.y-from->nextcp.y)-ysp->c;
        xsp->a = to->me.x-from->me.x-xsp->c-xsp->b;
        ysp->a = to->me.y-from->me.y-ysp->c-ysp->b;

(That's the [1,3,3,(1)] of the basis function).

All SplineFindExtrema is going to do is find the roots of the derivative, with a bit of (probably questionable) fiddling around with the lower bits of the floating point calculations.

[skef]

Oh, and you'll also need to strip out any single -1s as those could coincidentally wind up having the maximum Y magnitude.

[linusromer]

Thanks for your hints: So a, b, c, d of Spline1D are the coefficients of the Bézier polynomial in t. (Just stating for others that encounter the same problem: The documentation is not very clear when it says «The Spline1D structures give the equations for the x and y coordinates respectively»)

Yes, I have forgot to check if only te2==-1, I have added it.

I have also added a treatment for quadratic splines. The measures of linearity are the same for quadratic and cubic splines (therefore the factors .5 and 3) as they could theoretically be mixed.

static bigreal Linearity(Spline *s) {
    if ( s->islinear ) return 0;  
    BasePoint ftunit = BPSub(s->to->me, s->from->me);
    bigreal ftlen = BPNorm(ftunit);
    if ( ftlen==0 ) return -1; /* flag for error, no norming possible */
    ftunit = BPScale(ftunit, 1/ftlen);
    BasePoint nextcpscaled = BPScale(BPSub(s->from->nextcp, s->from->me),100./ftlen);
    if ( s->order2 ) return .5*fabs(BPCross(ftunit, nextcpscaled));
    BasePoint prevcpscaled = BPScale(BPSub(s->to->prevcp, s->from->me),100./ftlen);
    /* we just look a the bezier coefficients of the polynomial in t */
    /* in y dimension (divided by 3): */
    Spline1D ys1d;
    ys1d.d = 0;
    ys1d.c = BPCross(ftunit, nextcpscaled); /* rotate it to horizontal*/
    ys1d.b = BPCross(ftunit, prevcpscaled)-2*ys1d.c;
    ys1d.a = -ys1d.b-ys1d.c;
    extended te1, te2;
    SplineFindExtrema(&ys1d, &te1, &te2);
    if ( te1==-1 && te2==-1 ) return 0;
    if ( te2==-1 ) return 3*fabs(((ys1d.a*te1+ys1d.b)*te1+ys1d.c)*te1);
    return 3*fmax(fabs(((ys1d.a*te1+ys1d.b)*te1+ys1d.c)*te1), 
    fabs(((ys1d.a*te2+ys1d.b)*te2+ys1d.c)*te2));
}

But then I am still wondering what to do when converting a node of a quadratic spline to a tangent. How would one treat a situation like the lower left (the solution at the moment depicted to the lower mid is probably suboptimal, the lower right could be a better solution)?

And while testing I have found out another suboptimal issue: During the conversion of a node that joins two splines, the function SPChangePointType() is called twice. Hence, every calculation is done twice. This is not a tragedy but might confuse other testers as well.

[linusromer]

@skef Could you give me a hint for the manipulation of the control point of a quadratic spline: Do I have to manipulate from->nextcp and to->prevcp simultaneously or is it sufficient to change one of them? And how is a line then defined?

[skef]

The question of handling of nextcp and prevcp for quadratic splines is a bit murky. I always try to sync them up to be safe.

Quadratic lines are a special case -- you set each control point to the same value as the associated on-curve point, so that's the only time they should differ for a quadratic.

[skef]

But then I am still wondering what to do when converting a node of a quadratic spline to a tangent. How would one treat a situation like the lower left (the solution at the moment depicted to the lower mid is probably suboptimal, the lower right could be a better solution)?

"Ties" are one argument for being able to influence which side of a tangent conversion becomes a line by selecting a control point. a) You don't see that as a good approach and b) we've already decided on a linearity measure for the general case. I think if you follow the logic of those decisions through, the result is that ties shouldn't matter and can be broken arbitrarily: If the user favors one direction over another they should manipulate the control points before converting to make the favored direction clear to the algorithm.

[linusromer]

Thanks for your hint. Unfortunately, converting a quadratic spline to a line does not work yet for me. This is how I try to make sp->prev linear (and sp->next starting in the correct direction):

                                sp->prev->islinear = true;
                sp->prev->from->nextcp = sp->prev->from->me;
                sp->prevcp = sp->me;
                if ( sp->prev->order2 ) {
                    BasePoint inter;
                    if ( IntersectLines(&inter,&sp->me,&sp->prev->from->me,&sp->nextcp,&sp->next->to->me) ) {
                        sp->nextcp = inter;
                        sp->next->to->prevcp = inter;
                    } /* else just leave things as they are */
                }

And it looks like this:

[iterumllc]

Whenever you change the position of base or control points of a spline you have to call one of the SplineRefigures (usually just SplineRefigure() afterwards. I'm guessing that's what is missing.