← Back to branch summary

~ubuntu-branches/ubuntu/jaunty/luatex/jaunty

« back to all changes in this revision

Viewing changes to src/libs/luafontforge/fontforge/fontforge/splineoverlap.c

  • Committer: Bazaar Package Importer
  • Author(s): Norbert Preining
  • Date: 2007-09-24 12:56:11 UTC
  • Revision ID: james.westby@ubuntu.com-20070924125611-a8ge689azbptxvla
Tags: upstream-0.11.2
Import upstream version 0.11.2

Show diffs side-by-side

added added

removed removed

Lines of Context:
src/libs/luafontforge/fontforge/fontforge/splineoverlap.c
 
1
/* Copyright (C) 2000-2007 by George Williams */
 
2
/*
 
3
 * Redistribution and use in source and binary forms, with or without
 
4
 * modification, are permitted provided that the following conditions are met:
 
5
 
 
6
 * Redistributions of source code must retain the above copyright notice, this
 
7
 * list of conditions and the following disclaimer.
 
8
 
 
9
 * Redistributions in binary form must reproduce the above copyright notice,
 
10
 * this list of conditions and the following disclaimer in the documentation
 
11
 * and/or other materials provided with the distribution.
 
12
 
 
13
 * The name of the author may not be used to endorse or promote products
 
14
 * derived from this software without specific prior written permission.
 
15
 
 
16
 * THIS SOFTWARE IS PROVIDED BY THE AUTHOR ``AS IS'' AND ANY EXPRESS OR IMPLIED
 
17
 * WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF
 
18
 * MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO
 
19
 * EVENT SHALL THE AUTHOR BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
 
20
 * SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
 
21
 * PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS;
 
22
 * OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY,
 
23
 * WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR
 
24
 * OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF
 
25
 * ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
 
26
 */
 
27
#include "pfaedit.h"
 
28
#include "splinefont.h"
 
29
#include "edgelist2.h"
 
30
#include <math.h>
 
31
#ifdef HAVE_IEEEFP_H
 
32
# include <ieeefp.h>            /* Solaris defines isnan in ieeefp rather than math.h */
 
33
#endif
 
34
#include <stdarg.h>
 
35
 
 
36
#include <gwidget.h>            /* For PostNotice */
 
37
 
 
38
/* First thing we do is divide each spline into a set of sub-splines each of */
 
39
/*  which is monotonic in both x and y (always increasing or decreasing)     */
 
40
/* Then we compare each monotonic spline with every other one and see if they*/
 
41
/*  intersect.  If they do, split each up into sub-sub-segments and create an*/
 
42
/*  intersection point (note we need to be a little careful if an intersec-  */
 
43
/*  tion happens at an end point. We don't need to create a intersection for */
 
44
/*  two adjacent splines, there isn't a real intersection... but if a third  */
 
45
/*  spline crosses that point (or ends there) then all three (four) splines  */
 
46
/*  need to be joined into an intersection point)                            */
 
47
/* Nasty things happen if splines are coincident. They will almost never be  */
 
48
/*  perfectly coincident and will keep crossing and recrossing as rounding   */
 
49
/*  errors suggest one is before the other. Look for coincident splines and  */
 
50
/*  treat the places they start and stop being coincident as intersections   */
 
51
/*  then when we find needed splines below look for these guys and ignore    */
 
52
/*  recrossings of splines which are close together                          */
 
53
/* Figure out if each monotonic sub-spline is needed or not                  */
 
54
/*  (Note: It was tempting to split the bits up into real splines rather     */
 
55
/*   than keeping them as sub-sections of the original. Unfortunately this   */
 
56
/*   splitting introduced rounding errors which meant that we got more       */
 
57
/*   intersections, which meant that splines could be both needed and un.    */
 
58
/*   so I don't do that until later)                                         */
 
59
/*  if the spline hasn't been tagged yet:                                    */
 
60
/*   does the spline change greater in x or y?                               */
 
61
/*   draw a line parallel to the OTHER axis which hits our spline and doesn't*/
 
62
/*    hit any endpoints (or intersections, which are end points too now)     */
 
63
/*   count the winding number (as we do this we can mark other splines as    */
 
64
/*    needed or not) and figure out if our spline is needed                  */
 
65
/*  So run through the list of intersections                                 */
 
66
/*      At an intersection there should be an even number of needed monos.   */
 
67
/*      Use this as the basis of a new splineset, trace it around until      */
 
68
/*        we get back to the start intersection (should happen)              */
 
69
/*        (Note: We may need to reverse a monotonic sub-spline or two)       */
 
70
/*        As we go, mark each monotonic as having been used                  */
 
71
/*  Keep doing this until all needed exits from all intersections have been  */
 
72
/*      used.                                                                */
 
73
/*  The free up our temporary data structures, merge in any open splinesets  */
 
74
/*      free the old closed splinesets                                       */
 
75
 
 
76
typedef struct mlist {
 
77
    Spline *s;
 
78
    Monotonic *m;                       /* May get slightly munched but will */
 
79
                        /* always have right spline. we fix when we need it */
 
80
    extended t;
 
81
    int isend;
 
82
    BasePoint unit;
 
83
    struct mlist *next;
 
84
} MList;
 
85
 
 
86
typedef struct intersection {
 
87
    MList *monos;
 
88
    BasePoint inter;
 
89
    struct intersection *next;
 
90
} Intersection;
 
91
 
 
92
static char *glyphname=NULL;
 
93
 
 
94
static void SOError(char *format,...) {
 
95
    va_list ap;
 
96
    va_start(ap,format);
 
97
    if ( glyphname==NULL )
 
98
        fprintf(stderr, "Internal Error: " );
 
99
    else
 
100
        fprintf(stderr, "Internal Error in %s: ", glyphname );
 
101
    vfprintf(stderr,format,ap);
 
102
}
 
103
 
 
104
static Monotonic *SplineToMonotonic(Spline *s,extended startt,extended endt,
 
105
        Monotonic *last,int exclude) {
 
106
    Monotonic *m;
 
107
    BasePoint start, end;
 
108
 
 
109
    start.x = ((s->splines[0].a*startt+s->splines[0].b)*startt+s->splines[0].c)*startt
 
110
                + s->splines[0].d;
 
111
    start.y = ((s->splines[1].a*startt+s->splines[1].b)*startt+s->splines[1].c)*startt
 
112
                + s->splines[1].d;
 
113
    end.x = ((s->splines[0].a*endt+s->splines[0].b)*endt+s->splines[0].c)*endt
 
114
                + s->splines[0].d;
 
115
    end.y = ((s->splines[1].a*endt+s->splines[1].b)*endt+s->splines[1].c)*endt
 
116
                + s->splines[1].d;
 
117
    if ( (real) (((start.x+end.x)/2)==start.x || (real) ((start.x+end.x)/2)==end.x) &&
 
118
            (real) (((start.y+end.y)/2)==start.y || (real) ((start.y+end.y)/2)==end.y) ) {
 
119
        /* The distance between the two extrema is so small */
 
120
        /*  as to be unobservable. In other words we'd end up with a zero*/
 
121
        /*  length spline */
 
122
        if ( endt==1.0 && last!=NULL && last->s==s )
 
123
            last->tend = endt;
 
124
return( last );
 
125
    }
 
126
 
 
127
    m = chunkalloc(sizeof(Monotonic));
 
128
    m->s = s;
 
129
    m->tstart = startt;
 
130
    m->tend = endt;
 
131
    m->exclude = exclude;
 
132
 
 
133
    if ( end.x>start.x ) {
 
134
        m->xup = true;
 
135
        m->b.minx = start.x;
 
136
        m->b.maxx = end.x;
 
137
    } else {
 
138
        m->b.minx = end.x;
 
139
        m->b.maxx = start.x;
 
140
    }
 
141
    if ( end.y>start.y ) {
 
142
        m->yup = true;
 
143
        m->b.miny = start.y;
 
144
        m->b.maxy = end.y;
 
145
    } else {
 
146
        m->b.miny = end.y;
 
147
        m->b.maxy = start.y;
 
148
    }
 
149
 
 
150
    if ( last!=NULL ) {
 
151
        last->next = m;
 
152
        last->linked = m;
 
153
        m->prev = last;
 
154
    }
 
155
return( m );
 
156
}
 
157
 
 
158
static int SSIsSelected(SplineSet *spl) {
 
159
    SplinePoint *sp;
 
160
 
 
161
    for ( sp=spl->first; ; ) {
 
162
        if ( sp->selected )
 
163
return( true );
 
164
        if ( sp->next==NULL )
 
165
return( false );
 
166
        sp = sp->next->to;
 
167
        if ( sp==spl->first )
 
168
return( false );
 
169
    }
 
170
}
 
171
 
 
172
static int BpSame(BasePoint *bp1, BasePoint *bp2) {
 
173
    BasePoint mid;
 
174
 
 
175
    mid.x = (bp1->x+bp2->x)/2; mid.y = (bp1->y+bp2->y)/2;
 
176
    if ( (bp1->x==mid.x || bp2->x==mid.x) &&
 
177
            (bp1->y==mid.y || bp2->y==mid.y))
 
178
return( true );
 
179
 
 
180
return( false );
 
181
}
 
182
 
 
183
static int SSRmNullSplines(SplineSet *spl) {
 
184
    Spline *s, *first, *next;
 
185
 
 
186
    first = NULL;
 
187
    for ( s=spl->first->next ; s!=first; s=next ) {
 
188
        next = s->to->next;
 
189
        if ( ((s->splines[0].a>-.01 && s->splines[0].a<.01 &&
 
190
                s->splines[0].b>-.01 && s->splines[0].b<.01 &&
 
191
                s->splines[1].a>-.01 && s->splines[1].a<.01 &&
 
192
                s->splines[1].b>-.01 && s->splines[1].b<.01) ||
 
193
                /* That describes a null spline (a line between the same end-point) */
 
194
               RealNear((s->from->nextcp.x-s->from->me.x)*(s->to->me.y-s->to->prevcp.y)-
 
195
                        (s->from->nextcp.y-s->from->me.y)*(s->to->me.x-s->to->prevcp.x),0)) &&
 
196
                /* And the above describes a point with a spline between it */
 
197
                /*  and itself where the spline covers no area (the two cps */
 
198
                /*  point in the same direction) */
 
199
                BpSame(&s->from->me,&s->to->me)) {
 
200
            if ( next==s )
 
201
return( true );
 
202
            if ( next->from->selected ) s->from->selected = true;
 
203
            s->from->next = next;
 
204
            s->from->nextcp = next->from->nextcp;
 
205
            s->from->nonextcp = next->from->nonextcp;
 
206
            s->from->nextcpdef = next->from->nextcpdef;
 
207
            SplinePointFree(next->from);
 
208
            if ( spl->first==next->from )
 
209
                spl->last = spl->first = s->from;
 
210
            next->from = s->from;
 
211
            SplineFree(s);
 
212
        } else {
 
213
            if ( first==NULL )
 
214
                first = s;
 
215
        }
 
216
    }
 
217
return( false );
 
218
}
 
219
 
 
220
static Monotonic *SSToMContour(SplineSet *spl, Monotonic *start,
 
221
        Monotonic **end, enum overlap_type ot) {
 
222
    extended ts[4];
 
223
    Spline *first, *s;
 
224
    Monotonic *head=NULL, *last=NULL;
 
225
    int cnt, i, selected = false;
 
226
    extended lastt;
 
227
 
 
228
    if ( spl->first->prev==NULL )
 
229
return( start );                /* Open contours have no interior, ignore 'em */
 
230
    if ( spl->first->prev->from==spl->first &&
 
231
            spl->first->noprevcp && spl->first->nonextcp )
 
232
return( start );                /* Let's just remove single points */
 
233
 
 
234
    if ( ot==over_rmselected || ot==over_intersel || ot==over_fisel || ot==over_exclude ) {
 
235
        selected = SSIsSelected(spl);
 
236
        if ( ot==over_rmselected || ot==over_intersel || ot==over_fisel ) {
 
237
            if ( !selected )
 
238
return( start );
 
239
            selected = false;
 
240
        }
 
241
    }
 
242
 
 
243
    /* We blow up on zero length splines. And a zero length contour is nasty */
 
244
    if ( SSRmNullSplines(spl))
 
245
return( start );
 
246
 
 
247
    first = NULL;
 
248
    for ( s=spl->first->next; s!=first; s=s->to->next ) {
 
249
        if ( first==NULL ) first = s;
 
250
        cnt = Spline2DFindExtrema(s,ts);
 
251
        lastt = 0;
 
252
        for ( i=0; i<cnt; ++i ) {
 
253
            last = SplineToMonotonic(s,lastt,ts[i],last,selected);
 
254
            if ( head==NULL ) head = last;
 
255
            lastt=ts[i];
 
256
        }
 
257
        if ( lastt!=1.0 ) {
 
258
            last = SplineToMonotonic(s,lastt,1.0,last,selected);
 
259
            if ( head==NULL ) head = last;
 
260
        }
 
261
    }
 
262
    head->prev = last;
 
263
    last->next = head;
 
264
    if ( start==NULL )
 
265
        start = head;
 
266
    else
 
267
        (*end)->linked = head;
 
268
    *end = last;
 
269
return( start );
 
270
}
 
271
 
 
272
Monotonic *SSsToMContours(SplineSet *spl, enum overlap_type ot) {
 
273
    Monotonic *head=NULL, *last = NULL;
 
274
 
 
275
    while ( spl!=NULL ) {
 
276
        if ( spl->first->prev!=NULL )
 
277
            head = SSToMContour(spl,head,&last,ot);
 
278
        spl = spl->next;
 
279
    }
 
280
return( head );
 
281
}
 
282
 
 
283
static void _AddSpline(Intersection *il,Monotonic *m,extended t,int isend) {
 
284
    MList *ml;
 
285
 
 
286
    for ( ml=il->monos; ml!=NULL; ml=ml->next ) {
 
287
        if ( ml->s==m->s && RealNear( ml->t,t ) && ml->isend==isend )
 
288
return;
 
289
    }
 
290
 
 
291
    ml = chunkalloc(sizeof(MList));
 
292
    ml->next = il->monos;
 
293
    il->monos = ml;
 
294
    ml->s = m->s;
 
295
    ml->m = m;                  /* This may change. We'll fix it up later */
 
296
    ml->t = t;
 
297
    ml->isend = isend;
 
298
    if ( isend ) {
 
299
        if ( m->end!=NULL && m->end!=il )
 
300
            SOError("Resetting end.\n");
 
301
        m->end = il;
 
302
    } else {
 
303
        if ( m->start!=NULL && m->start!=il )
 
304
            SOError("Resetting start.\n");
 
305
        m->start = il;
 
306
    }
 
307
return;
 
308
}
 
309
 
 
310
static void AddSpline(Intersection *il,Monotonic *m,extended t) {
 
311
    MList *ml;
 
312
 
 
313
    if ( m->start==il || m->end==il )
 
314
return;
 
315
 
 
316
    for ( ml=il->monos; ml!=NULL; ml=ml->next ) {
 
317
        if ( ml->s==m->s && RealWithin( ml->t,t,.0001 ))
 
318
return;
 
319
    }
 
320
 
 
321
    ml = chunkalloc(sizeof(MList));
 
322
    ml->next = il->monos;
 
323
    il->monos = ml;
 
324
    ml->s = m->s;
 
325
    ml->m = m;                  /* This may change. We'll fix it up later */
 
326
    ml->t = t;
 
327
    ml->isend = true;
 
328
    if ( t-m->tstart < m->tend-t && RealNear(m->tstart,t) ) {
 
329
        if ( m->start!=NULL && m->start!=il )
 
330
            SOError("Resetting start.\n");
 
331
        m->start = il;
 
332
        ml->t = m->tstart;
 
333
        ml->isend = false;
 
334
        _AddSpline(il,m->prev,m->prev->tend,true);
 
335
    } else if ( RealNear(m->tend,t)) {
 
336
        if ( m->end!=NULL && m->end!=il )
 
337
            SOError("Resetting end.\n");
 
338
        m->end = il;
 
339
        ml->t = m->tend;
 
340
        _AddSpline(il,m->next,m->next->tstart,false);
 
341
    } else {
 
342
        /* Ok, if we've got a new intersection on this spline then break up */
 
343
        /*  the monotonic into two bits which end and start at this inter */
 
344
        if ( t<m->tstart || t>m->tend )
 
345
            SOError( "Attempt to subset monotonic rejoin inappropriately: %g should be [%g,%g]\n",
 
346
                    t, m->tstart, m->tend );
 
347
        else {
 
348
            /* It is monotonic, so a subset of it must also be */
 
349
            Monotonic *m2 = chunkalloc(sizeof(Monotonic));
 
350
            BasePoint pt;
 
351
            *m2 = *m;
 
352
            m->next = m2;
 
353
            m2->prev = m;
 
354
            m2->next->prev = m2;
 
355
            m->linked = m2;
 
356
            m->tend = t;
 
357
            m->end = il;
 
358
            m2->start = il;
 
359
            m2->tstart = t;
 
360
            pt.x = ((m->s->splines[0].a*m->tstart+m->s->splines[0].b)*m->tstart+
 
361
                    m->s->splines[0].c)*m->tstart+m->s->splines[0].d;
 
362
            pt.y = ((m->s->splines[1].a*m->tstart+m->s->splines[1].b)*m->tstart+
 
363
                    m->s->splines[1].c)*m->tstart+m->s->splines[1].d;
 
364
            if ( pt.x>il->inter.x ) {
 
365
                m->b.minx = il->inter.x;
 
366
                m->b.maxx = pt.x;
 
367
            } else {
 
368
                m->b.minx = pt.x;
 
369
                m->b.maxx = il->inter.x;
 
370
            }
 
371
            if ( pt.y>il->inter.y ) {
 
372
                m->b.miny = il->inter.y;
 
373
                m->b.maxy = pt.y;
 
374
            } else {
 
375
                m->b.miny = pt.y;
 
376
                m->b.maxy = il->inter.y;
 
377
            }
 
378
            pt.x = ((m2->s->splines[0].a*m2->tend+m2->s->splines[0].b)*m2->tend+
 
379
                    m2->s->splines[0].c)*m2->tend+m2->s->splines[0].d;
 
380
            pt.y = ((m2->s->splines[1].a*m2->tend+m2->s->splines[1].b)*m2->tend+
 
381
                    m2->s->splines[1].c)*m2->tend+m2->s->splines[1].d;
 
382
            if ( pt.x>il->inter.x ) {
 
383
                m2->b.minx = il->inter.x;
 
384
                m2->b.maxx = pt.x;
 
385
            } else {
 
386
                m2->b.minx = pt.x;
 
387
                m2->b.maxx = il->inter.x;
 
388
            }
 
389
            if ( pt.y>il->inter.y ) {
 
390
                m2->b.miny = il->inter.y;
 
391
                m2->b.maxy = pt.y;
 
392
            } else {
 
393
                m2->b.miny = pt.y;
 
394
                m2->b.maxy = il->inter.y;
 
395
            }
 
396
            _AddSpline(il,m2,t,false);
 
397
        }
 
398
    }
 
399
}
 
400
 
 
401
static void SetStartPoint(BasePoint *pt,Monotonic *m) {
 
402
    if ( m->tstart==0 )
 
403
        *pt = m->s->from->me;
 
404
    else if ( m->start!=NULL )
 
405
        *pt = m->start->inter;
 
406
    else {
 
407
        pt->x = ((m->s->splines[0].a*m->tstart+m->s->splines[0].b)*m->tstart +
 
408
                m->s->splines[0].c)*m->tstart + m->s->splines[0].d;
 
409
        pt->y = ((m->s->splines[1].a*m->tstart+m->s->splines[1].b)*m->tstart +
 
410
                m->s->splines[1].c)*m->tstart + m->s->splines[1].d;
 
411
    }
 
412
}
 
413
 
 
414
static void SetEndPoint(BasePoint *pt,Monotonic *m) {
 
415
    if ( m->tend==1.0 )
 
416
        *pt = m->s->to->me;
 
417
    else if ( m->end!=NULL )
 
418
        *pt = m->end->inter;
 
419
    else {
 
420
        pt->x = ((m->s->splines[0].a*m->tend+m->s->splines[0].b)*m->tend +
 
421
                m->s->splines[0].c)*m->tend + m->s->splines[0].d;
 
422
        pt->y = ((m->s->splines[1].a*m->tend+m->s->splines[1].b)*m->tend +
 
423
                m->s->splines[1].c)*m->tend + m->s->splines[1].d;
 
424
    }
 
425
}
 
426
 
 
427
static extended RoundToEndpoints(Monotonic *m,extended t,BasePoint *inter) {
 
428
    BasePoint end;
 
429
    extended bound;
 
430
 
 
431
    if ( t==0 || t==1 ) {
 
432
        if ( t==0 )
 
433
            *inter = m->s->from->me;
 
434
        else
 
435
            *inter = m->s->to->me;
 
436
return( t );
 
437
    }
 
438
 
 
439
    if ( t-m->tstart < m->tend-t ) {
 
440
        bound = m->tstart;
 
441
        SetStartPoint(&end,m);
 
442
    } else {
 
443
        bound = m->tend;
 
444
        SetEndPoint(&end,m);
 
445
    }
 
446
    if ( BpSame(&end,inter) || RealWithin(t,bound,.00001)) {
 
447
        *inter = end;
 
448
return( bound );
 
449
    }
 
450
 
 
451
return( t );
 
452
}
 
453
 
 
454
static extended Grad1(Spline1D *s1, Spline1D *s2,
 
455
        extended t1,extended t2 ) {
 
456
    /* d/dt[12] (m1(t1).x-m2(t2).x)^2 + (m1(t1).y-m2(t2).y)^2 */
 
457
    /* d/dt[12] (m1(t1).x^2 -2m1(t1).x*m2(t2).x + m2(t2).x^2) + (m1(t1).y^2 -2m1(t1).y*m2(t2).y + m2(t2).y^2) */
 
458
    extended val2 = ((s2->a*t2+s2->b)*t2+s2->c)*t2+s2->d;
 
459
 
 
460
return( ((((6*(extended)s1->a*s1->a*t1 +
 
461
            5*2*(extended)s1->a*s1->b)*t1 +
 
462
            4*(s1->b*(extended)s1->b+2*s1->a*(extended)s1->c))*t1 +
 
463
            3*2*(s1->a*(extended)s1->d+s1->b*(extended)s1->c))*t1 +
 
464
            2*(s1->c*(extended)s1->c+2*s1->b*(extended)s1->d))*t1 +
 
465
            2*s1->c*(extended)s1->d  -
 
466
             2*val2 * ((3*s1->a*t1 + 2*s1->b)*t1 + s1->c) );
 
467
}
 
468
 
 
469
static void GradImproveInter(Monotonic *m1, Monotonic *m2,
 
470
        extended *_t1,extended *_t2,BasePoint *inter) {
 
471
    Spline *s1 = m1->s, *s2 = m2->s;
 
472
    extended x1, x2, y1, y2;
 
473
    extended gt1=0, gt2=0, glen=1;
 
474
    extended error, olderr=1e10;
 
475
    extended factor = 4096;
 
476
    extended t1=*_t1, t2=*_t2;
 
477
    extended off, off2, yoff;
 
478
    int cnt=0;
 
479
    /* We want to find (t1,t2) so that (m1(t1)-m2(t2))^2==0 */
 
480
    /* Find the gradiant and move in the reverse direction */
 
481
    /* We know that the current values of (t1,t2) are close to an intersection*/
 
482
    /* so the grad should point correctly */
 
483
    /* d/dt[12] (m1(t1).x-m2(t2).x)^2 + (m1(t1).y-m2(t2).y)^2 */
 
484
    /* d/dt[12] (m1(t1).x^2 -2m1(t1).x*m2(t2).x + m2(t2).x^2) + (m1(t1).y^2 -2m1(t1).y*m2(t2).y + m2(t2).y^2) */
 
485
 
 
486
    forever {
 
487
        x1 = ((s1->splines[0].a*t1 + s1->splines[0].b)*t1 + s1->splines[0].c)*t1 + s1->splines[0].d;
 
488
        x2 = ((s2->splines[0].a*t2 + s2->splines[0].b)*t2 + s2->splines[0].c)*t2 + s2->splines[0].d;
 
489
        y1 = ((s1->splines[1].a*t1 + s1->splines[1].b)*t1 + s1->splines[1].c)*t1 + s1->splines[1].d;
 
490
        y2 = ((s2->splines[1].a*t2 + s2->splines[1].b)*t2 + s2->splines[1].c)*t2 + s2->splines[1].d;
 
491
        error = (x1-x2)*(x1-x2) + (y1-y2)*(y1-y2);
 
492
        if ( error>olderr ) {
 
493
            if ( olderr==1e10 )
 
494
    break;
 
495
            factor *= 2;
 
496
            if ( factor>4096*4096 )
 
497
    break;
 
498
            glen *= 2;
 
499
            t1 += gt1/glen;
 
500
            t2 += gt2/glen;
 
501
    continue;
 
502
        } else
 
503
            factor /= 1.4;
 
504
        if ( error<1e-11 )      /* Error is actually the square of the error */
 
505
    break;                      /* So this isn't as constraining as it looks */
 
506
 
 
507
        gt1 = Grad1(&s1->splines[0],&s2->splines[0],t1,t2) + Grad1(&s1->splines[1],&s2->splines[1],t1,t2);
 
508
        gt2 = Grad1(&s2->splines[0],&s1->splines[0],t2,t1) + Grad1(&s2->splines[1],&s1->splines[1],t2,t1);
 
509
        glen = esqrt(gt1*gt1 + gt2*gt2) * factor;
 
510
        if ( glen==0 )
 
511
    break;
 
512
        *_t1 = t1; *_t2 = t2;
 
513
        t1 -= gt1/glen;
 
514
        t2 -= gt2/glen;
 
515
        if ( isnan(t1) || isnan(t2)) {
 
516
            IError( "Nan in grad" );
 
517
    break;
 
518
        }
 
519
        olderr = error;
 
520
        ++cnt;
 
521
        if ( cnt>1000 )
 
522
    break;
 
523
    }
 
524
#if 0
 
525
    if ( cnt<=1 && error>=1e-11 )
 
526
        fprintf(stderr,"No Improvement\n" );
 
527
    else if ( cnt>1 )
 
528
        fprintf(stderr,"Improvement\n" );
 
529
#endif
 
530
    t1 = *_t1; t2 = *_t2;
 
531
    if ( t1<0 && t1>-.00001 ) *_t1 = t1 = 0;
 
532
    else if ( t1>1 && t1<1.00001 ) *_t1 = t1 = 1.0;
 
533
    if ( t2<0 && t2>-.00001 ) *_t2 = t2 = 0;
 
534
    else if ( t2>1 && t2<1.00001 ) *_t2 = t2 = 1.0;
 
535
    x1 = ((s1->splines[0].a*t1 + s1->splines[0].b)*t1 + s1->splines[0].c)*t1 + s1->splines[0].d;
 
536
    x2 = ((s2->splines[0].a*t2 + s2->splines[0].b)*t2 + s2->splines[0].c)*t2 + s2->splines[0].d;
 
537
    y1 = ((s1->splines[1].a*t1 + s1->splines[1].b)*t1 + s1->splines[1].c)*t1 + s1->splines[1].d;
 
538
    y2 = ((s2->splines[1].a*t2 + s2->splines[1].b)*t2 + s2->splines[1].c)*t2 + s2->splines[1].d;
 
539
    inter->x = (x1+x2)/2; inter->y = (y1+y2)/2;
 
540
 
 
541
    if ( (off=x1-x2)<0 ) off = -off;
 
542
    if ( (yoff=y1-y2)<0 ) yoff = -yoff;
 
543
    off += yoff;
 
544
 
 
545
    if ( t1<.0001 ) {
 
546
        t1 = 0;
 
547
        x1 = s1->splines[0].d;
 
548
        y1 = s1->splines[1].d;
 
549
    } else if ( t1>.9999 ) {
 
550
        t1 = 1.0;
 
551
        x1 = s1->splines[0].a+s1->splines[0].b+s1->splines[0].c+s1->splines[0].d;
 
552
        y1 = s1->splines[1].a+s1->splines[1].b+s1->splines[1].c+s1->splines[1].d;
 
553
    }
 
554
    if ( t2<.0001 ) {
 
555
        t2=0;
 
556
        x2 = s2->splines[0].d;
 
557
        y2 = s2->splines[1].d;
 
558
    } else if ( t2>.9999 ) {
 
559
        t2=1.0;
 
560
        x2 = s2->splines[0].a+s2->splines[0].b+s2->splines[0].c+s2->splines[0].d;
 
561
        y2 = s2->splines[1].a+s2->splines[1].b+s2->splines[1].c+s2->splines[1].d;
 
562
    }
 
563
    if ( (off2=x1-x2)<0 ) off2 = -off2;
 
564
    if ( (yoff=y1-y2)<0 ) yoff = -yoff;
 
565
    off2 += yoff;
 
566
    if ( off2<=off ) {
 
567
        *_t1 = t1; *_t2 = t2;
 
568
        inter->x = (x1+x2)/2; inter->y = (y1+y2)/2;
 
569
    }
 
570
}
 
571
 
 
572
static Intersection *AddIntersection(Intersection *ilist,Monotonic *m1,
 
573
        Monotonic *m2,extended t1,extended t2,BasePoint *inter) {
 
574
    Intersection *il;
 
575
    extended ot1 = t1, ot2 = t2;
 
576
 
 
577
    /* Fixup some rounding errors */
 
578
    GradImproveInter(m1,m2,&t1,&t2,inter);
 
579
    if ( t1<m1->tstart || t1>m1->tend || t2<m2->tstart || t2>m2->tend )
 
580
return( ilist );
 
581
 
 
582
#if 0
 
583
    if (( m1->start==m2->start && m1->start!=NULL && RealNear(t1,m1->tstart) && RealNear(t2,m2->tstart)) ||
 
584
            ( m1->start==m2->end && m1->start!=NULL && RealNear(t1,m1->tstart) && RealNear(t2,m2->tend)) ||
 
585
            ( m1->end==m2->start && m1->end!=NULL && RealNear(t1,m1->tend) && RealNear(t2,m2->tstart)) ||
 
586
            ( m1->end==m2->end && m1->end!=NULL && RealNear(t1,m1->tend) && RealNear(t2,m2->tend)) )
 
587
return( ilist );
 
588
    else
 
589
    if ( RealWithin(t1,1.0,.01) && RealWithin(t2,0.0,.01) && BpSame(&m1->s->to->me,&m2->s->from->me)) {
 
590
        t1 = 1.0;
 
591
        t2 = 0.0;
 
592
    } else if ( RealWithin(t2,1.0,.01) && RealWithin(t1,0.0,.01) && BpSame(&m2->s->to->me,&m1->s->from->me)) {
 
593
        t1 = 0.0;
 
594
        t2 = 1.0;
 
595
    } else {
 
596
#endif
 
597
        t1 = RoundToEndpoints(m1,t1,inter);
 
598
        t2 = RoundToEndpoints(m2,t2,inter);
 
599
        t1 = RoundToEndpoints(m1,t1,inter);     /* Do it twice. rounding t2 can mean we now need to round t1 */
 
600
#if 0
 
601
    }
 
602
#endif
 
603
 
 
604
    if (( m1->s->to == m2->s->from && RealWithin(t1,1.0,.01) && RealWithin(t2,0,.01)) ||
 
605
            ( m1->s->from == m2->s->to && RealWithin(t1,0,.01) && RealWithin(t2,1.0,.01)))
 
606
return( ilist );
 
607
 
 
608
    if (( t1==m1->tstart && m1->start!=NULL &&
 
609
            (inter->x!=m1->start->inter.x || inter->y!=m1->start->inter.y)) ||
 
610
        ( t1==m1->tend && m1->end!=NULL &&
 
611
            (inter->x!=m1->end->inter.x || inter->y!=m1->end->inter.y)))
 
612
        t1 = ot1;
 
613
    if (( t2==m2->tstart && m2->start!=NULL &&
 
614
            (inter->x!=m2->start->inter.x || inter->y!=m2->start->inter.y)) ||
 
615
        ( t2==m2->tend && m2->end!=NULL &&
 
616
            (inter->x!=m2->end->inter.x || inter->y!=m2->end->inter.y)))
 
617
        t2 = ot2;
 
618
 
 
619
    /* The ordinary join of one spline to the next doesn't really count */
 
620
    /*  Or one monotonic sub-spline to the next either */
 
621
    if (( m1->next==m2 && RealNear(t1,m1->tend) && RealNear(t2,m2->tstart)) ||
 
622
            (m2->next==m1 && RealNear(t1,m1->tstart) && RealNear(t2,m2->tend)) )
 
623
return( ilist );
 
624
 
 
625
    if ( RealWithin(m1->tstart,t1,.01) )
 
626
        il = m1->start;
 
627
    else if ( RealWithin(m1->tend,t1,.01) )
 
628
        il = m1->end;
 
629
    else
 
630
        il = NULL;
 
631
    if ( il!=NULL &&
 
632
            ((RealWithin(m2->tstart,t2,.01) && m2->start==il) ||
 
633
             (RealWithin(m2->tend,t2,.01) && m2->end==il)) )
 
634
return( ilist );
 
635
 
 
636
    for ( il = ilist; il!=NULL; il=il->next ) {
 
637
        if ( RealWithin(il->inter.x,inter->x,.01) && RealWithin(il->inter.y,inter->y,.01)) {
 
638
            AddSpline(il,m1,t1);
 
639
            AddSpline(il,m2,t2);
 
640
return( ilist );
 
641
        }
 
642
    }
 
643
 
 
644
    il = chunkalloc(sizeof(Intersection));
 
645
    il->inter = *inter;
 
646
    il->next = ilist;
 
647
    AddSpline(il,m1,t1);
 
648
    AddSpline(il,m2,t2);
 
649
return( il );
 
650
}
 
651
 
 
652
static extended BoundIterateSplineSolve(Spline1D *sp, extended tmin, extended tmax,
 
653
        extended sought,double err) {
 
654
    extended t = IterateSplineSolve(sp,tmin,tmax,sought,err);
 
655
    if ( t<tmin || t>tmax )
 
656
return( -1 );
 
657
 
 
658
return( t );
 
659
}
 
660
 
 
661
static Intersection *FindMonotonicIntersection(Intersection *ilist,Monotonic *m1,Monotonic *m2) {
 
662
    /* I believe that two monotonic cubics can still intersect in two points */
 
663
    /*  so we can't just check if the splines are on oposite sides of each */
 
664
    /*  other at top and bottom */
 
665
    DBounds b;
 
666
    const double error = .0001;
 
667
    BasePoint pt;
 
668
    extended t1,t2;
 
669
 
 
670
    b.minx = m1->b.minx>m2->b.minx ? m1->b.minx : m2->b.minx;
 
671
    b.maxx = m1->b.maxx<m2->b.maxx ? m1->b.maxx : m2->b.maxx;
 
672
    b.miny = m1->b.miny>m2->b.miny ? m1->b.miny : m2->b.miny;
 
673
    b.maxy = m1->b.maxy<m2->b.maxy ? m1->b.maxy : m2->b.maxy;
 
674
 
 
675
    if ( b.maxy==b.miny && b.minx==b.maxx ) {
 
676
        extended x1,y1, x2,y2;
 
677
        if ( m1->next==m2 || m2->next==m1 )
 
678
return( ilist );                /* Not interesting. Only intersection is at an endpoint */
 
679
        if ( ((m1->start==m2->start || m1->end==m2->start) && m2->start!=NULL) ||
 
680
                ((m1->start==m2->end || m1->end==m2->end ) && m2->end!=NULL ))
 
681
return( ilist );
 
682
        pt.x = b.minx; pt.y = b.miny;
 
683
        if ( m1->b.maxx-m1->b.minx > m1->b.maxy-m1->b.miny )
 
684
            t1 = BoundIterateSplineSolve(&m1->s->splines[0],m1->tstart,m1->tend,b.minx,error);
 
685
        else
 
686
            t1 = BoundIterateSplineSolve(&m1->s->splines[1],m1->tstart,m1->tend,b.miny,error);
 
687
        if ( m2->b.maxx-m2->b.minx > m2->b.maxy-m2->b.miny )
 
688
            t2 = BoundIterateSplineSolve(&m2->s->splines[0],m2->tstart,m2->tend,b.minx,error);
 
689
        else
 
690
            t2 = BoundIterateSplineSolve(&m2->s->splines[1],m2->tstart,m2->tend,b.miny,error);
 
691
        if ( t1!=-1 && t2!=-1 ) {
 
692
            x1 = ((m1->s->splines[0].a*t1+m1->s->splines[0].b)*t1+m1->s->splines[0].c)*t1+m1->s->splines[0].d;
 
693
            y1 = ((m1->s->splines[1].a*t1+m1->s->splines[1].b)*t1+m1->s->splines[1].c)*t1+m1->s->splines[1].d;
 
694
            x2 = ((m2->s->splines[0].a*t2+m2->s->splines[0].b)*t2+m2->s->splines[0].c)*t2+m2->s->splines[0].d;
 
695
            y2 = ((m2->s->splines[1].a*t2+m2->s->splines[1].b)*t2+m2->s->splines[1].c)*t2+m2->s->splines[1].d;
 
696
            if ( x1-x2>-.01 && x1-x2<.01 && y1-y2>-.01 && y1-y2<.01 )
 
697
                ilist = AddIntersection(ilist,m1,m2,t1,t2,&pt);
 
698
        }
 
699
    } else if ( b.maxy==b.miny ) {
 
700
        extended x1,x2;
 
701
        if ( m1->next==m2 || m2->next==m1 )
 
702
return( ilist );                /* Not interesting. Only intersection is at an endpoint */
 
703
        t1 = BoundIterateSplineSolve(&m1->s->splines[1],m1->tstart,m1->tend,b.miny,error);
 
704
        t2 = BoundIterateSplineSolve(&m2->s->splines[1],m2->tstart,m2->tend,b.miny,error);
 
705
        if ( t1!=-1 && t2!=-1 ) {
 
706
            x1 = ((m1->s->splines[0].a*t1+m1->s->splines[0].b)*t1+m1->s->splines[0].c)*t1+m1->s->splines[0].d;
 
707
            x2 = ((m2->s->splines[0].a*t2+m2->s->splines[0].b)*t2+m2->s->splines[0].c)*t2+m2->s->splines[0].d;
 
708
            if ( x1-x2>-.01 && x1-x2<.01 ) {
 
709
                pt.x = (x1+x2)/2; pt.y = b.miny;
 
710
                ilist = AddIntersection(ilist,m1,m2,t1,t2,&pt);
 
711
            }
 
712
        }
 
713
    } else if ( b.maxx==b.minx ) {
 
714
        extended y1,y2;
 
715
        if ( m1->next==m2 || m2->next==m1 )
 
716
return( ilist );                /* Not interesting. Only intersection is at an endpoint */
 
717
        t1 = BoundIterateSplineSolve(&m1->s->splines[0],m1->tstart,m1->tend,b.minx,error);
 
718
        t2 = BoundIterateSplineSolve(&m2->s->splines[0],m2->tstart,m2->tend,b.minx,error);
 
719
        if ( t1!=-1 && t2!=-1 ) {
 
720
            y1 = ((m1->s->splines[1].a*t1+m1->s->splines[1].b)*t1+m1->s->splines[1].c)*t1+m1->s->splines[1].d;
 
721
            y2 = ((m2->s->splines[1].a*t2+m2->s->splines[1].b)*t2+m2->s->splines[1].c)*t2+m2->s->splines[1].d;
 
722
            if ( y1-y2>-.01 && y1-y2<.01 ) {
 
723
                pt.x = b.minx; pt.y = (y1+y2)/2;
 
724
                ilist = AddIntersection(ilist,m1,m2,t1,t2,&pt);
 
725
            }
 
726
        }
 
727
    } else if ( b.maxy-b.miny > b.maxx-b.minx ) {
 
728
        extended diff, y, x1,x2, x1o,x2o;
 
729
        extended t1=0,t2=0, t1o=0,t2o=0 ;
 
730
 
 
731
        diff = (b.maxy-b.miny)/32;
 
732
        y = b.miny;
 
733
        x1o = x2o = 0;
 
734
        while ( y<b.maxy ) {
 
735
            while ( y<b.maxy ) {
 
736
                t1o = BoundIterateSplineSolve(&m1->s->splines[1],m1->tstart,m1->tend,b.miny,error);
 
737
                t2o = BoundIterateSplineSolve(&m2->s->splines[1],m2->tstart,m2->tend,b.miny,error);
 
738
                if ( t1o!=-1 && t2o!=-1 )
 
739
            break;
 
740
                y += diff;
 
741
            }
 
742
            x1o = ((m1->s->splines[0].a*t1o+m1->s->splines[0].b)*t1o+m1->s->splines[0].c)*t1o+m1->s->splines[0].d;
 
743
            x2o = ((m2->s->splines[0].a*t2o+m2->s->splines[0].b)*t2o+m2->s->splines[0].c)*t2o+m2->s->splines[0].d;
 
744
            if ( x1o!=x2o )
 
745
        break;
 
746
            y += diff;
 
747
        }
 
748
#if 0
 
749
        if ( x1o==x2o ) {
 
750
            pt.x = x1o; pt.y = y;
 
751
            ilist = AddIntersection(ilist,m1,m2,t1o,t2o,&pt);
 
752
        }
 
753
#endif
 
754
        for ( y+=diff; ; y += diff ) {
 
755
            /* I used to say y<=b.maxy in the above for statement. */
 
756
            /*  that seemed to get rounding errors on the mac, so we do it */
 
757
            /*  like this now: */
 
758
            if ( y>b.maxy ) {
 
759
                if ( y<b.maxy+diff/4 ) y = b.maxy;
 
760
                else
 
761
        break;
 
762
            }
 
763
            t1 = BoundIterateSplineSolve(&m1->s->splines[1],m1->tstart,m1->tend,y,error);
 
764
            t2 = BoundIterateSplineSolve(&m2->s->splines[1],m2->tstart,m2->tend,y,error);
 
765
            if ( t1==-1 || t2==-1 )
 
766
        continue;
 
767
            x1 = ((m1->s->splines[0].a*t1+m1->s->splines[0].b)*t1+m1->s->splines[0].c)*t1+m1->s->splines[0].d;
 
768
            x2 = ((m2->s->splines[0].a*t2+m2->s->splines[0].b)*t2+m2->s->splines[0].c)*t2+m2->s->splines[0].d;
 
769
#if 0
 
770
            if ( x1==x2 && x1o!=x2o ) {
 
771
                pt.x = x1; pt.y = y;
 
772
                ilist = AddIntersection(ilist,m1,m2,t1,t2,&pt);
 
773
                if ( m1->tend==t1 && m1->s==m1->next->s ) m1 = m1->next;
 
774
                if ( m2->tend==t2 && m2->s==m2->next->s ) m2 = m2->next;
 
775
            } else
 
776
#endif
 
777
            if ( x1o!=x2o && (x1o>x2o) != ( x1>x2 ) ) {
 
778
                /* A cross over has occured. (assume we have a small enough */
 
779
                /*  region that three cross-overs can't have occurred) */
 
780
                /* Use a binary search to track it down */
 
781
                extended ytop, ybot;
 
782
                ytop = y;
 
783
                ybot = y-diff;
 
784
                while ( ytop!=ybot ) {
 
785
                    extended ytest = (ytop+ybot)/2;
 
786
                    extended t1t, t2t;
 
787
                    t1t = BoundIterateSplineSolve(&m1->s->splines[1],m1->tstart,m1->tend,ytest,error);
 
788
                    t2t = BoundIterateSplineSolve(&m2->s->splines[1],m2->tstart,m2->tend,ytest,error);
 
789
                    x1 = ((m1->s->splines[0].a*t1t+m1->s->splines[0].b)*t1t+m1->s->splines[0].c)*t1t+m1->s->splines[0].d;
 
790
                    x2 = ((m2->s->splines[0].a*t2t+m2->s->splines[0].b)*t2t+m2->s->splines[0].c)*t2t+m2->s->splines[0].d;
 
791
                    if ( t1t==-1 || t2t==-1 ) {
 
792
                        SOError( "Can't find something in range.\n" );
 
793
                break;
 
794
                    } else if (( x1-x2<error && x1-x2>-error ) || ytop==ytest || ybot==ytest ) {
 
795
                        pt.y = ytest; pt.x = (x1+x2)/2;
 
796
                        ilist = AddIntersection(ilist,m1,m2,t1t,t2t,&pt);
 
797
                        b.maxy = m1->b.maxy<m2->b.maxy ? m1->b.maxy : m2->b.maxy;
 
798
                break;
 
799
                    } else if ( (x1o>x2o) != ( x1>x2 ) ) {
 
800
                        ytop = ytest;
 
801
                    } else {
 
802
                        ybot = ytest;
 
803
                    }
 
804
                }
 
805
                x1 = x1o; x1o = x2o; x2o = x1;
 
806
            } else {
 
807
                x1o = x1; x2o = x2;
 
808
            }
 
809
        }
 
810
    } else {
 
811
        extended diff, x, y1,y2, y1o,y2o;
 
812
        extended t1=0,t2=0, t1o=0,t2o=0 ;
 
813
 
 
814
        diff = (b.maxx-b.minx)/32;
 
815
        x = b.minx;
 
816
        y1o = y2o = 0;
 
817
        while ( x<b.maxx ) {
 
818
            while ( x<b.maxx ) {
 
819
                t1o = BoundIterateSplineSolve(&m1->s->splines[0],m1->tstart,m1->tend,b.minx,error);
 
820
                t2o = BoundIterateSplineSolve(&m2->s->splines[0],m2->tstart,m2->tend,b.minx,error);
 
821
                if ( t1o!=-1 && t2o!=-1 )
 
822
            break;
 
823
                x += diff;
 
824
            }
 
825
            y1o = ((m1->s->splines[1].a*t1o+m1->s->splines[1].b)*t1o+m1->s->splines[1].c)*t1o+m1->s->splines[1].d;
 
826
            y2o = ((m2->s->splines[1].a*t2o+m2->s->splines[1].b)*t2o+m2->s->splines[1].c)*t2o+m2->s->splines[1].d;
 
827
            if ( y1o!=y2o )
 
828
        break;
 
829
            x += diff;
 
830
        }
 
831
#if 0
 
832
        if ( y1o==y2o ) {
 
833
            pt.y = y1o; pt.x = x;
 
834
            ilist = AddIntersection(ilist,m1,m2,t1o,t2o,&pt);
 
835
        }
 
836
#endif
 
837
        y1 = y2 = 0;
 
838
        for ( x+=diff; ; x += diff ) {
 
839
            if ( x>b.maxx ) {
 
840
                if ( x<b.maxx+diff/4 ) x = b.maxx;
 
841
                else
 
842
        break;
 
843
            }
 
844
            t1 = BoundIterateSplineSolve(&m1->s->splines[0],m1->tstart,m1->tend,x,error);
 
845
            t2 = BoundIterateSplineSolve(&m2->s->splines[0],m2->tstart,m2->tend,x,error);
 
846
            if ( t1==-1 || t2==-1 )
 
847
        continue;
 
848
            y1 = ((m1->s->splines[1].a*t1+m1->s->splines[1].b)*t1+m1->s->splines[1].c)*t1+m1->s->splines[1].d;
 
849
            y2 = ((m2->s->splines[1].a*t2+m2->s->splines[1].b)*t2+m2->s->splines[1].c)*t2+m2->s->splines[1].d;
 
850
#if 0
 
851
            if ( y1==y2 && y1o!=y2o ) {
 
852
                pt.y = y1; pt.x = x;
 
853
                ilist = AddIntersection(ilist,m1,m2,t1,t2,&pt);
 
854
                if ( m1->tend==t1 && m1->s==m1->next->s ) m1 = m1->next;
 
855
                if ( m2->tend==t2 && m2->s==m2->next->s ) m2 = m2->next;
 
856
            } else
 
857
#endif
 
858
            if ( (y1o>y2o) != ( y1>y2 ) ) {
 
859
                /* A cross over has occured. (assume we have a small enough */
 
860
                /*  region that three cross-overs can't have occurred) */
 
861
                /* Use a binary search to track it down */
 
862
                extended xtop, xbot;
 
863
                xtop = x;
 
864
                xbot = x-diff;
 
865
                while ( xtop!=xbot ) {
 
866
                    extended xtest = (xtop+xbot)/2;
 
867
                    extended t1t, t2t;
 
868
                    t1t = BoundIterateSplineSolve(&m1->s->splines[0],m1->tstart,m1->tend,xtest,error);
 
869
                    t2t = BoundIterateSplineSolve(&m2->s->splines[0],m2->tstart,m2->tend,xtest,error);
 
870
                    y1 = ((m1->s->splines[1].a*t1t+m1->s->splines[1].b)*t1t+m1->s->splines[1].c)*t1t+m1->s->splines[1].d;
 
871
                    y2 = ((m2->s->splines[1].a*t2t+m2->s->splines[1].b)*t2t+m2->s->splines[1].c)*t2t+m2->s->splines[1].d;
 
872
                    if ( t1t==-1 || t2t==-1 ) {
 
873
                        SOError( "Can't find something in range.\n" );
 
874
                break;
 
875
                    } else if (( y1-y2<error && y1-y2>-error ) || xtop==xtest || xbot==xtest ) {
 
876
                        pt.x = xtest; pt.y = (y1+y2)/2;
 
877
                        ilist = AddIntersection(ilist,m1,m2,t1t,t2t,&pt);
 
878
                        b.maxx = m1->b.maxx<m2->b.maxx ? m1->b.maxx : m2->b.maxx;
 
879
                break;
 
880
                    } else if ( (y1o>y2o) != ( y1>y2 ) ) {
 
881
                        xtop = xtest;
 
882
                    } else {
 
883
                        xbot = xtest;
 
884
                    }
 
885
                }
 
886
                y1 = y1o; y1o = y2o; y2o = y1;
 
887
            } else {
 
888
                y1o = y1; y2o = y2;
 
889
            }
 
890
        }
 
891
    }
 
892
return( ilist );
 
893
}
 
894
 
 
895
#if 0
 
896
static void SubsetSpline1(Spline1D *subset,Spline1D *orig,
 
897
        extended tstart,extended tend,extended d) {
 
898
    extended f = tend-tstart;
 
899
 
 
900
    subset->a = orig->a*f*f*f;
 
901
    subset->b = (orig->b + 3*orig->a*tstart) *f*f;
 
902
    subset->c = (orig->c + (2*orig->b + 3*orig->a*tstart)*tstart) * f;
 
903
    subset->d = d;
 
904
}
 
905
#endif
 
906
 
 
907
static extended SplineContainsPoint(Monotonic *m,BasePoint *pt) {
 
908
    int which, nw;
 
909
    extended t;
 
910
    BasePoint slope;
 
911
    const double error = .0001;
 
912
 
 
913
    which = ( m->b.maxx-m->b.minx > m->b.maxy-m->b.miny )? 0 : 1;
 
914
    nw = !which;
 
915
    t = BoundIterateSplineSolve(&m->s->splines[which],m->tstart,m->tend,(&pt->x)[which],error);
 
916
    if ( (slope.x = (3*m->s->splines[0].a*t+2*m->s->splines[0].b)*t+m->s->splines[0].c)<0 )
 
917
        slope.x = -slope.x;
 
918
    if ( (slope.y = (3*m->s->splines[1].a*t+2*m->s->splines[1].b)*t+m->s->splines[1].c)<0 )
 
919
        slope.y = -slope.y;
 
920
    if ( t==-1 || (slope.y>slope.x)!=which ) {
 
921
        nw = which;
 
922
        which = 1-which;
 
923
        t = BoundIterateSplineSolve(&m->s->splines[which],m->tstart,m->tend,(&pt->x)[which],error);
 
924
    }
 
925
    if ( t!=-1 && RealWithin((&pt->x)[nw],
 
926
           ((m->s->splines[nw].a*t+m->s->splines[nw].b)*t +
 
927
                m->s->splines[nw].c)*t + m->s->splines[nw].d,.1 ))
 
928
return( t );
 
929
 
 
930
return( -1 );
 
931
}
 
932
 
 
933
/* If two splines are coincident, then pretend they intersect at both */
 
934
/*  end-points and nowhere else */
 
935
static int CoincidentIntersect(Monotonic *m1,Monotonic *m2,BasePoint *pts,
 
936
        extended *t1s,extended *t2s) {
 
937
    const double error = .0001;
 
938
    int cnt=0;
 
939
    extended t, t2, diff;
 
940
 
 
941
    if ( m1==m2 || m1->next==m2 || m1->prev==m2 )
 
942
return( false );                /* Can't be coincident. Adjacent */
 
943
    /* Actually adjacent splines can double back on themselves */
 
944
 
 
945
    if ( (m1->xup==m2->xup && m1->yup==m2->yup) ||
 
946
            ((m1->xup!=m2->xup || (m1->b.minx==m1->b.maxx && m2->b.minx==m2->b.maxx)) ||
 
947
             (m1->yup!=m2->yup || (m1->b.miny==m1->b.maxy && m2->b.miny==m2->b.maxy))))
 
948
        /* A match is possible */;
 
949
    else
 
950
return( false );
 
951
 
 
952
    SetStartPoint(&pts[cnt],m1);
 
953
    t1s[cnt] = m1->tstart;
 
954
    if ( (t2s[cnt] = SplineContainsPoint(m2,&pts[cnt]))!=-1 )
 
955
        ++cnt;
 
956
 
 
957
    SetEndPoint(&pts[cnt],m1);
 
958
    t1s[cnt] = m1->tend;
 
959
    if ( (t2s[cnt] = SplineContainsPoint(m2,&pts[cnt]))!=-1 )
 
960
        ++cnt;
 
961
 
 
962
    if ( cnt!=2 ) {
 
963
        SetStartPoint(&pts[cnt],m2);
 
964
        t2s[cnt] = m2->tstart;
 
965
        if ( (t1s[cnt] = SplineContainsPoint(m1,&pts[cnt]))!=-1 )
 
966
            ++cnt;
 
967
    }
 
968
 
 
969
    if ( cnt!=2 ) {
 
970
        SetEndPoint(&pts[cnt],m2);
 
971
        t2s[cnt] = m2->tend;
 
972
        if ( (t1s[cnt] = SplineContainsPoint(m1,&pts[cnt]))!=-1 )
 
973
            ++cnt;
 
974
    }
 
975
 
 
976
    if ( cnt!=2 )
 
977
return( false );
 
978
 
 
979
    if ( RealWithin(t1s[0],t1s[1],.01) )
 
980
return( false );
 
981
 
 
982
    /* Ok, if we've gotten this far we know that two of the end points are  */
 
983
    /*  on both splines.                                                    */
 
984
#if 0
 
985
    /* Interesting. The following algorithem does not produce a unique      */
 
986
    /*  representation of the subset spline.  Must test each point          */
 
987
    /* Then subset each spline so that [0,1] runs along the range we just   */
 
988
    /*  found for it. Compare these two subsets, if they are almost equal   */
 
989
    /*  then assume they are coincident, and return the two endpoints as    */
 
990
    /*  the intersection points.                                            */
 
991
    SubsetSpline1(&subs1,&m1->s->splines[0],t1s[0],t1s[1],pts[0].x);
 
992
    SubsetSpline1(&subs2,&m2->s->splines[0],t2s[0],t2s[1],pts[0].x);
 
993
    if ( !RealWithin(subs1.a,subs2.a,.1) || !RealWithin(subs1.b,subs2.b,.1) ||
 
994
            !RealWithin(subs1.c,subs2.c,.2) /* Don't need to check d, it's always pts[0].x */)
 
995
return( false );
 
996
 
 
997
    SubsetSpline1(&subs1,&m1->s->splines[1],t1s[0],t1s[1],pts[0].y);
 
998
    SubsetSpline1(&subs2,&m2->s->splines[1],t2s[0],t2s[1],pts[0].y);
 
999
    if ( !RealWithin(subs1.a,subs2.a,.1) || !RealWithin(subs1.b,subs2.b,.1) ||
 
1000
            !RealWithin(subs1.c,subs2.c,.2) /* Don't need to check d, it's always pts[0].x */)
 
1001
return( false );
 
1002
 
 
1003
    t1s[2] = t2s[2] = -1;
 
1004
return( true );
 
1005
#else
 
1006
    t1s[2] = t2s[2] = -1;
 
1007
    if ( !m1->s->knownlinear || !m1->s->knownlinear ) {
 
1008
        if ( t1s[1]<t1s[0] ) {
 
1009
            extended temp = t1s[1]; t1s[1] = t1s[0]; t1s[0] = temp;
 
1010
            temp = t2s[1]; t2s[1] = t2s[0]; t2s[0] = temp;
 
1011
        }
 
1012
        diff = (t1s[1]-t1s[0])/16;
 
1013
        for ( t=t1s[0]+diff; t<t1s[1]-diff/4; t += diff ) {
 
1014
            BasePoint here, slope;
 
1015
            here.x = ((m1->s->splines[0].a*t+m1->s->splines[0].b)*t+m1->s->splines[0].c)*t+m1->s->splines[0].d;
 
1016
            here.y = ((m1->s->splines[1].a*t+m1->s->splines[1].b)*t+m1->s->splines[1].c)*t+m1->s->splines[1].d;
 
1017
            if ( (slope.x = (3*m1->s->splines[0].a*t+2*m1->s->splines[0].b)*t+m1->s->splines[0].c)<0 )
 
1018
                slope.x = -slope.x;
 
1019
            if ( (slope.y = (3*m1->s->splines[1].a*t+2*m1->s->splines[1].b)*t+m1->s->splines[1].c)<0 )
 
1020
                slope.y = -slope.y;
 
1021
            if ( slope.y>slope.x ) {
 
1022
                t2 = BoundIterateSplineSolve(&m2->s->splines[1],t2s[0],t2s[1],here.y,error);
 
1023
                if ( t2==-1 || !RealWithin(here.x,((m2->s->splines[0].a*t2+m2->s->splines[0].b)*t2+m2->s->splines[0].c)*t2+m2->s->splines[0].d,.1))
 
1024
return( false );
 
1025
            } else {
 
1026
                t2 = BoundIterateSplineSolve(&m2->s->splines[0],t2s[0],t2s[1],here.x,error);
 
1027
                if ( t2==-1 || !RealWithin(here.y,((m2->s->splines[1].a*t2+m2->s->splines[1].b)*t2+m2->s->splines[1].c)*t2+m2->s->splines[1].d,.1))
 
1028
return( false );
 
1029
            }
 
1030
        }
 
1031
    }
 
1032
 
 
1033
return( true );
 
1034
#endif
 
1035
}
 
1036
 
 
1037
static void FigureProperMonotonicsAtIntersections(Intersection *ilist) {
 
1038
    MList *ml, *ml2, *mlnext, *prev, *p2;
 
1039
 
 
1040
    while ( ilist!=NULL ) {
 
1041
        for ( ml=ilist->monos; ml!=NULL; ml=ml->next ) {
 
1042
            if ( (ml->t==ml->m->tstart && !ml->isend) ||
 
1043
                    (ml->t==ml->m->tend && ml->isend))
 
1044
                /* It's right */;
 
1045
            else if ( ml->t>ml->m->tstart ) {
 
1046
                while ( ml->t>ml->m->tend ) {
 
1047
                    ml->m = ml->m->next;
 
1048
                    if ( ml->m->s!=ml->s ) {
 
1049
                        SOError("we could not find a matching monotonic\n" );
 
1050
                break;
 
1051
                    }
 
1052
                }
 
1053
            } else {
 
1054
                while ( ml->t<ml->m->tstart ) {
 
1055
                    ml->m = ml->m->prev;
 
1056
                    if ( ml->m->s!=ml->s ) {
 
1057
                        SOError( "we could not find a matching monotonic\n" );
 
1058
                break;
 
1059
                    }
 
1060
                }
 
1061
            }
 
1062
            if ( ml->t==ml->m->tstart && ml->isend )
 
1063
                ml->m = ml->m->prev;
 
1064
            else if ( ml->t==ml->m->tend && !ml->isend )
 
1065
                ml->m = ml->m->next;
 
1066
            if ( ml->t!=ml->m->tstart && ml->t!=ml->m->tend )
 
1067
                SOError( "we could not find a matching monotonic time\n" );
 
1068
        }
 
1069
        for ( prev=NULL, ml=ilist->monos; ml!=NULL; ml = mlnext ) {
 
1070
            mlnext = ml->next;
 
1071
            if ( ml->m->start==ml->m->end ) {
 
1072
                for ( p2 = ml, ml2=ml->next; ml2!=NULL; p2=ml2, ml2 = ml2->next ) {
 
1073
                    if ( ml2->m==ml->m )
 
1074
                break;
 
1075
                }
 
1076
                if ( ml2!=NULL ) {
 
1077
                    if ( ml2==mlnext ) mlnext = ml2->next;
 
1078
                    p2->next = ml2->next;
 
1079
                    chunkfree(ml2,sizeof(*ml2));
 
1080
                }
 
1081
                if ( prev==NULL )
 
1082
                    ilist->monos = mlnext;
 
1083
                else
 
1084
                    prev->next = mlnext;
 
1085
                chunkfree(ml,sizeof(*ml));
 
1086
            }
 
1087
        }
 
1088
#if 0
 
1089
        for ( ml=ilist->monos; ml!=NULL; ml=ml->next ) {
 
1090
            Monotonic *search;
 
1091
            MList *ml2;
 
1092
            extended t;
 
1093
            if ( ml->m->start == ilist ) {
 
1094
                search = ml->m->prev;
 
1095
                t = ( ml->m->tstart==0 ) ? 1.0 : ml->m->tstart;
 
1096
            } else {
 
1097
                search = ml->m->next;
 
1098
                t = ( ml->m->tend==1.0 ) ? 0.0 : ml->m->tend;
 
1099
            }
 
1100
            for ( ml2=ilist->monos; ml2!=NULL && ml2->m!=search; ml2=ml2->next );
 
1101
            if ( ml2==NULL ) {
 
1102
                ml2 = chunkalloc(sizeof(MList));
 
1103
                ml2->m = search;
 
1104
                ml2->s = search->s;
 
1105
                ml2->t = t;
 
1106
                ml2->next = ml->next;
 
1107
                ml->next = ml2;
 
1108
                ml = ml2;
 
1109
            }
 
1110
        }
 
1111
#endif
 
1112
        ilist = ilist->next;
 
1113
    }
 
1114
}
 
1115
 
 
1116
static void Validate(Monotonic *ms, Intersection *ilist) {
 
1117
    MList *ml;
 
1118
    int mcnt;
 
1119
 
 
1120
    while ( ilist!=NULL ) {
 
1121
        for ( mcnt=0, ml=ilist->monos; ml!=NULL; ml=ml->next ) {
 
1122
            if ( ml->m->isneeded ) ++mcnt;
 
1123
            if ( ml->m->start!=ilist && ml->m->end!=ilist )
 
1124
                SOError( "Intersection (%g,%g) not on a monotonic which should contain it.\n",
 
1125
                        ilist->inter.x, ilist->inter.y );
 
1126
        }
 
1127
        if ( mcnt&1 )
 
1128
            SOError( "Odd number of needed monotonic sections at intersection. (%g,%g)\n",
 
1129
                    ilist->inter.x,ilist->inter.y );
 
1130
        ilist = ilist->next;
 
1131
    }
 
1132
 
 
1133
    while ( ms!=NULL ) {
 
1134
        if ( ms->prev->end!=ms->start )
 
1135
            SOError( "Mismatched intersection.\n");
 
1136
        ms = ms->linked;
 
1137
    }
 
1138
}
 
1139
 
 
1140
static Intersection *FindIntersections(Monotonic *ms, enum overlap_type ot) {
 
1141
    Monotonic *m1, *m2;
 
1142
    BasePoint pts[9];
 
1143
    extended t1s[10], t2s[10];
 
1144
    Intersection *ilist=NULL;
 
1145
    int i;
 
1146
    int wasc;
 
1147
 
 
1148
    for ( m1=ms; m1!=NULL; m1=m1->linked ) {
 
1149
        for ( m2=m1->linked; m2!=NULL; m2=m2->linked ) {
 
1150
            if ( m2->b.minx > m1->b.maxx ||
 
1151
                    m2->b.maxx < m1->b.minx ||
 
1152
                    m2->b.miny > m1->b.maxy ||
 
1153
                    m2->b.maxy < m1->b.miny )
 
1154
        continue;               /* Can't intersect */;
 
1155
            wasc = CoincidentIntersect(m1,m2,pts,t1s,t2s);
 
1156
            if ( wasc || m1->s->knownlinear || m2->s->knownlinear ||
 
1157
                    (m1->s->splines[0].a==0 && m1->s->splines[1].a==0 &&
 
1158
                     m2->s->splines[0].a==0 && m2->s->splines[1].a==0 )) {
 
1159
                if ( !wasc && SplinesIntersect(m1->s,m2->s,pts,t1s,t2s)<=0 )
 
1160
        continue;
 
1161
                for ( i=0; i<4 && t1s[i]!=-1; ++i ) {
 
1162
                    if ( t1s[i]>=m1->tstart && t1s[i]<=m1->tend &&
 
1163
                            t2s[i]>=m2->tstart && t2s[i]<=m2->tend ) {
 
1164
                        ilist = AddIntersection(ilist,m1,m2,t1s[i],t2s[i],&pts[i]);
 
1165
                    }
 
1166
                }
 
1167
        continue;
 
1168
            }
 
1169
            ilist = FindMonotonicIntersection(ilist,m1,m2);
 
1170
        }
 
1171
    }
 
1172
 
 
1173
    FigureProperMonotonicsAtIntersections(ilist);
 
1174
 
 
1175
    /* Now suppose we have a contour which intersects nothing? */
 
1176
    /* with no intersections we lose track of it and it will vanish */
 
1177
    /* That's not a good idea. Make sure each contour has at least one inter */
 
1178
    if ( ot!=over_findinter && ot!=over_fisel ) {
 
1179
        for ( m1=ms; m1!=NULL; m1=m2->linked ) {
 
1180
            if ( m1->start==NULL && m1->end==NULL ) {
 
1181
                Intersection *il;
 
1182
                il = chunkalloc(sizeof(Intersection));
 
1183
                il->inter = m1->s->from->me;
 
1184
                il->next = ilist;
 
1185
                AddSpline(il,m1,0);
 
1186
                AddSpline(il,m1->prev,1.0);
 
1187
                ilist = il;
 
1188
            }
 
1189
            for ( m2=m1; m2->linked==m2->next; m2=m2->linked );
 
1190
        }
 
1191
    }
 
1192
 
 
1193
return( ilist );
 
1194
}
 
1195
 
 
1196
static int dcmp(const void *_p1, const void *_p2) {
 
1197
    const extended *dpt1 = _p1, *dpt2 = _p2;
 
1198
    if ( *dpt1>*dpt2 )
 
1199
return( 1 );
 
1200
    else if ( *dpt1<*dpt2 )
 
1201
return( -1 );
 
1202
 
 
1203
return( 0 );
 
1204
}
 
1205
 
 
1206
static extended *FindOrderedEndpoints(Monotonic *ms,int which) {
 
1207
    int cnt;
 
1208
    Monotonic *m;
 
1209
    extended *ends;
 
1210
    int i,j,k;
 
1211
 
 
1212
    for ( m=ms, cnt=0; m!=NULL; m=m->linked, ++cnt );
 
1213
    ends = galloc((2*cnt+1)*sizeof(extended));
 
1214
    for ( m=ms, cnt=0; m!=NULL; m=m->linked, cnt+=2 ) {
 
1215
        if ( m->start!=NULL )
 
1216
            ends[cnt] = (&m->start->inter.x)[which];
 
1217
        else if ( m->tstart==0 )
 
1218
            ends[cnt] = (&m->s->from->me.x)[which];
 
1219
        else
 
1220
            ends[cnt] = ((m->s->splines[which].a*m->tstart+m->s->splines[which].b)*m->tstart+
 
1221
                    m->s->splines[which].c)*m->tstart+m->s->splines[which].d;
 
1222
        if ( m->end!=NULL )
 
1223
            ends[cnt+1] = (&m->end->inter.x)[which];
 
1224
        else if ( m->tend==1.0 )
 
1225
            ends[cnt+1] = (&m->s->to->me.x)[which];
 
1226
        else
 
1227
            ends[cnt+1] = ((m->s->splines[which].a*m->tend+m->s->splines[which].b)*m->tend+
 
1228
                    m->s->splines[which].c)*m->tend+m->s->splines[which].d;
 
1229
    }
 
1230
 
 
1231
    qsort(ends,cnt,sizeof(extended),dcmp);
 
1232
    for ( i=0; i<cnt; ++i ) {
 
1233
        for ( j=i; j<cnt && ends[i]==ends[j]; ++j );
 
1234
        if ( j>i+1 ) {
 
1235
            for ( k=i+1; j<cnt; ends[k++] = ends[j++]);
 
1236
            cnt-= j-k;
 
1237
        }
 
1238
    }
 
1239
    ends[cnt] = 1e10;
 
1240
return( ends );
 
1241
}
 
1242
 
 
1243
static int mcmp(const void *_p1, const void *_p2) {
 
1244
    const Monotonic * const *mpt1 = _p1, * const *mpt2 = _p2;
 
1245
    if ( (*mpt1)->other>(*mpt2)->other )
 
1246
return( 1 );
 
1247
    else if ( (*mpt1)->other<(*mpt2)->other )
 
1248
return( -1 );
 
1249
 
 
1250
return( 0 );
 
1251
}
 
1252
 
 
1253
int MonotonicFindAt(Monotonic *ms,int which, extended test, Monotonic **space ) {
 
1254
    /* Find all monotonic sections which intersect the line (x,y)[which] == test */
 
1255
    /*  find the value of the other coord on that line */
 
1256
    /*  Order them (by the other coord) */
 
1257
    /*  then run along that line figuring out which monotonics are needed */
 
1258
    extended t;
 
1259
    Monotonic *m, *mm;
 
1260
    int i, j, k, cnt;
 
1261
    const double error = .0001;
 
1262
    int nw = !which;
 
1263
 
 
1264
    for ( m=ms, i=0; m!=NULL; m=m->linked ) {
 
1265
        if (( which==0 && test >= m->b.minx && test <= m->b.maxx ) ||
 
1266
                ( which==1 && test >= m->b.miny && test <= m->b.maxy )) {
 
1267
            /* Lines parallel to the direction we are testing just get in the */
 
1268
            /*  way and don't add any useful info */
 
1269
            if ( m->s->knownlinear &&
 
1270
                    (( which==1 && m->s->from->me.y==m->s->to->me.y ) ||
 
1271
                        (which==0 && m->s->from->me.x==m->s->to->me.x)))
 
1272
    continue;
 
1273
            t = BoundIterateSplineSolve(&m->s->splines[which],m->tstart,m->tend,test,error);
 
1274
            if ( t==-1 )
 
1275
    continue;
 
1276
            m->t = t;
 
1277
            if ( t==m->tend ) t -= (m->tend-m->tstart)/100;
 
1278
            else if ( t==m->tstart ) t += (m->tend-m->tstart)/100;
 
1279
            m->other = ((m->s->splines[nw].a*t+m->s->splines[nw].b)*t+
 
1280
                    m->s->splines[nw].c)*t+m->s->splines[nw].d;
 
1281
            space[i++] = m;
 
1282
        }
 
1283
    }
 
1284
    cnt = i;
 
1285
 
 
1286
    /* Things get a little tricky at end-points */
 
1287
    for ( i=0; i<cnt; ++i ) {
 
1288
        m = space[i];
 
1289
        if ( m->t==m->tend ) {
 
1290
            /* Ignore horizontal/vertical lines (as appropriate) */
 
1291
            for ( mm=m->next; mm!=m; mm=mm->next ) {
 
1292
                if ( !mm->s->knownlinear )
 
1293
            break;
 
1294
                if (( which==1 && mm->s->from->me.y!=m->s->to->me.y ) ||
 
1295
                        (which==0 && mm->s->from->me.x!=m->s->to->me.x))
 
1296
            break;
 
1297
            }
 
1298
        } else if ( m->t==m->tstart ) {
 
1299
            for ( mm=m->prev; mm!=m; mm=mm->prev ) {
 
1300
                if ( !mm->s->knownlinear )
 
1301
            break;
 
1302
                if (( which==1 && mm->s->from->me.y!=m->s->to->me.y ) ||
 
1303
                        (which==0 && mm->s->from->me.x!=m->s->to->me.x))
 
1304
            break;
 
1305
            }
 
1306
        } else
 
1307
    break;
 
1308
        /* If the next monotonic continues in the same direction, and we found*/
 
1309
        /*  it too, then don't count both. They represent the same intersect */
 
1310
        /* If they are in oposite directions then they cancel each other out */
 
1311
        /*  and that is correct */
 
1312
        if ( mm!=m &&   /* Should always be true */
 
1313
                (&mm->xup)[which]==(&m->xup)[which] ) {
 
1314
            for ( j=cnt-1; j>=0; --j )
 
1315
                if ( space[j]==mm )
 
1316
            break;
 
1317
            if ( j!=-1 ) {
 
1318
                /* remove mm */
 
1319
                for ( k=j+1; k<cnt; ++k )
 
1320
                    space[k-1] = space[k];
 
1321
                --cnt;
 
1322
                if ( i>j ) --i;
 
1323
            }
 
1324
        }
 
1325
    }
 
1326
 
 
1327
    space[cnt] = NULL; space[cnt+1] = NULL;
 
1328
    qsort(space,cnt,sizeof(Monotonic *),mcmp);
 
1329
return(cnt);
 
1330
}
 
1331
 
 
1332
static void FigureNeeds(Monotonic *ms,int which, extended test, Monotonic **space,
 
1333
        enum overlap_type ot, int ignore_close) {
 
1334
    /* Find all monotonic sections which intersect the line (x,y)[which] == test */
 
1335
    /*  find the value of the other coord on that line */
 
1336
    /*  Order them (by the other coord) */
 
1337
    /*  then run along that line figuring out which monotonics are needed */
 
1338
    int i, j, winding, ew, was_close, close;
 
1339
 
 
1340
    MonotonicFindAt(ms,which,test,space);
 
1341
 
 
1342
#if 0           /* Really slow, and it fixes some problems at the expense of causing others */
 
1343
    for ( i=0; space[i+1]!=NULL; ++i ) {
 
1344
        /* If two splines are very close to each other, we may miss an */
 
1345
        /*  intersection. If that has happened, reorder the splines */
 
1346
        if ( space[i+1]->other - space[i]->other < .1 ) {
 
1347
            extended oi, oi1, ti, ti1, test2;
 
1348
            if ( which==1 ) {
 
1349
                if ( space[i+1]->b.miny > space[i]->b.miny )
 
1350
                    test2 = space[i]->b.miny;
 
1351
                else
 
1352
                    test2 = space[i+1]->b.miny;
 
1353
            } else {
 
1354
                if ( space[i+1]->b.minx > space[i]->b.minx )
 
1355
                    test2 = space[i]->b.minx;
 
1356
                else
 
1357
                    test2 = space[i+1]->b.minx;
 
1358
            }
 
1359
            ti = BoundIterateSplineSolve(&space[i]->s->splines[which],
 
1360
                    space[i]->tstart,space[i]->tend,test2,error);
 
1361
            ti1= BoundIterateSplineSolve(&space[i+1]->s->splines[which],
 
1362
                    space[i+1]->tstart,space[i+1]->tend,test2,error);
 
1363
            oi = ((space[i]->s->splines[nw].a*ti+space[i]->s->splines[nw].b)*ti+
 
1364
                    space[i]->s->splines[nw].c)*ti+space[i]->s->splines[nw].d;
 
1365
            oi1= ((space[i+1]->s->splines[nw].a*ti1+space[i+1]->s->splines[nw].b)*ti1+
 
1366
                    space[i+1]->s->splines[nw].c)*ti1+space[i+1]->s->splines[nw].d;
 
1367
            if ( oi1<oi ) {
 
1368
                m = space[i];
 
1369
                space[i] = space[i+1];
 
1370
                space[i+1] = m;
 
1371
 fprintf( stderr, "Flipped\n" );
 
1372
            }
 
1373
        }
 
1374
    }
 
1375
#endif
 
1376
 
 
1377
    winding = 0; ew = 0; was_close = false;
 
1378
    for ( i=0; space[i]!=NULL; ++i ) {
 
1379
        int needed, unneeded, inverted=false;
 
1380
        Monotonic *m;
 
1381
        int new;
 
1382
        int nwinding;
 
1383
      retry:
 
1384
        needed = false, unneeded = false;
 
1385
        nwinding=winding;
 
1386
        new=ew;
 
1387
        m = space[i];
 
1388
        if ( m->exclude )
 
1389
            new += ( (&m->xup)[which] ? 1 : -1 );
 
1390
        else
 
1391
            nwinding += ( (&m->xup)[which] ? 1 : -1 );
 
1392
        if ( ot==over_remove || ot==over_rmselected ) {
 
1393
            if ( winding==0 || nwinding==0 )
 
1394
                needed = true;
 
1395
            else
 
1396
                unneeded = true;
 
1397
        } else if ( ot==over_intersect || ot==over_intersel ) {
 
1398
            if ( (winding>-2 && winding<2 && nwinding>-2 && nwinding<2) ||
 
1399
                    ((winding<=-2 || winding>=2) && (nwinding<=-2 && nwinding>=2)))
 
1400
                unneeded = true;
 
1401
            else
 
1402
                needed = true;
 
1403
        } else if ( ot == over_exclude ) {
 
1404
            if ( (( winding==0 || nwinding==0 ) && ew==0 && new==0 ) ||
 
1405
                    (winding!=0 && (( ew!=0 && new==0 ) || ( ew==0 && new!=0))) )
 
1406
                needed = true;
 
1407
            else
 
1408
                unneeded = true;
 
1409
        }
 
1410
        if ( space[i+1]!=NULL )
 
1411
            close = space[i+1]->other-space[i]->other < 1;
 
1412
        else
 
1413
            close = false;
 
1414
        if (( !close && !was_close ) || ignore_close ) {
 
1415
            if (( m->isneeded || m->isunneeded ) && m->isneeded!=needed ) {
 
1416
                for ( j=i+1; space[j]!=NULL && space[j]->other-m->other<.5; ++j ) {
 
1417
                    if ( space[j]->start==m->start && space[j]->end==m->end &&
 
1418
                            (space[j]->isneeded == needed ||
 
1419
                             (!space[j]->isneeded && !space[j]->isunneeded))) {
 
1420
                        space[i] = space[j];
 
1421
                        space[j] = m;
 
1422
                        m = space[i];
 
1423
                break;
 
1424
                    } else if ( !inverted && space[j]->other-m->other<.001 &&
 
1425
                            (((&space[j]->xup)[which] == (&m->xup)[which] &&
 
1426
                              (space[j]->isneeded == needed ||
 
1427
                               (!space[j]->isneeded && !space[j]->isunneeded))) ||
 
1428
                             ((&space[j]->xup)[which] != (&m->xup)[which] &&
 
1429
                              (space[j]->isneeded != needed ||
 
1430
                               (!space[j]->isneeded && !space[j]->isunneeded)))) ) {
 
1431
                        space[i] = space[j];
 
1432
                        space[j] = m;
 
1433
                        inverted = true;
 
1434
                  goto retry;
 
1435
                    }
 
1436
                }
 
1437
            }
 
1438
            if ( !m->isneeded && !m->isunneeded ) {
 
1439
                m->isneeded = needed; m->isunneeded = unneeded;
 
1440
                m->when_set = test;             /* Debugging */
 
1441
            } else if ( m->isneeded!=needed || m->isunneeded!=unneeded )
 
1442
                SOError( "monotonic is both needed and unneeded.\n" );
 
1443
        }
 
1444
        winding = nwinding;
 
1445
        ew = new;
 
1446
        was_close = close;
 
1447
    }
 
1448
    if ( winding!=0 )
 
1449
        SOError( "Winding number did not return to 0 when %s=%g\n",
 
1450
                which ? "y" : "x", test );
 
1451
}
 
1452
 
 
1453
struct gaps { extended test, len; int which; };
 
1454
 
 
1455
static int gcmp(const void *_p1, const void *_p2) {
 
1456
    const struct gaps *gpt1 = _p1, *gpt2 = _p2;
 
1457
    if ( gpt1->len > gpt2->len )
 
1458
return( 1 );
 
1459
    else if ( gpt1->len < gpt2->len )
 
1460
return( -1 );
 
1461
 
 
1462
return( 0 );
 
1463
}
 
1464
 
 
1465
static void FindNeeded(Monotonic *ms,enum overlap_type ot) {
 
1466
    extended *ends[2];
 
1467
    Monotonic *m, **space;
 
1468
    extended top, bottom, test;
 
1469
    int t,b,i,j,k,cnt,which;
 
1470
    struct gaps *gaps;
 
1471
    extended min_gap;
 
1472
 
 
1473
    if ( ms==NULL )
 
1474
return;
 
1475
 
 
1476
    ends[0] = FindOrderedEndpoints(ms,0);
 
1477
    ends[1] = FindOrderedEndpoints(ms,1);
 
1478
 
 
1479
    for ( m=ms, cnt=0; m!=NULL; m=m->linked, ++cnt );
 
1480
    space = galloc((cnt+2)*sizeof(Monotonic*));
 
1481
    gaps = galloc(2*cnt*sizeof(struct gaps));
 
1482
 
 
1483
    /* Look for the longest splines without interruptions first. These are */
 
1484
    /* least likely to cause problems and will give us a good basis from which*/
 
1485
    /* to make guesses should rounding errors occur later */
 
1486
    for ( j=k=0; j<2; ++j )
 
1487
        for ( i=0; ends[j][i+1]!=1e10; ++i ) {
 
1488
            gaps[k].which = j;
 
1489
            gaps[k].len = (ends[j][i+1]-ends[j][i]);
 
1490
            gaps[k++].test = (ends[j][i+1]+ends[j][i])/2;
 
1491
        }
 
1492
    qsort(gaps,k,sizeof(struct gaps),gcmp);
 
1493
    min_gap = 1e10;
 
1494
    for ( m=ms; m!=NULL; m=m->linked ) {
 
1495
        if ( m->b.maxx-m->b.minx > m->b.maxy-m->b.miny ) {
 
1496
            if ( min_gap > m->b.maxx-m->b.minx ) min_gap = m->b.maxx-m->b.minx;
 
1497
        } else {
 
1498
            if ( m->b.maxy-m->b.miny==0 )
 
1499
                fprintf( stderr, "Foo\n");
 
1500
            if ( min_gap > m->b.maxy-m->b.miny ) min_gap = m->b.maxy-m->b.miny;
 
1501
        }
 
1502
    }
 
1503
    if ( min_gap<.5 ) min_gap = .5;
 
1504
    for ( i=0; i<k && gaps[i].len>=min_gap; ++i )
 
1505
        FigureNeeds(ms,gaps[i].which,gaps[i].test,space,ot,0);
 
1506
 
 
1507
    for ( m=ms; m!=NULL; m=m->linked ) if ( !m->isneeded && !m->isunneeded ) {
 
1508
        if ( m->b.maxx-m->b.minx > m->b.maxy-m->b.miny ) {
 
1509
            top = m->b.maxx;
 
1510
            bottom = m->b.minx;
 
1511
            which = 0;
 
1512
        } else {
 
1513
            top = m->b.maxy;
 
1514
            bottom = m->b.miny;
 
1515
            which = 1;
 
1516
        }
 
1517
        for ( b=0; ends[which][b]<=bottom; ++b );
 
1518
        for ( t=b; ends[which][t]<top; ++t );
 
1519
        --t;
 
1520
        /* b points to an endpoint which is greater than bottom */
 
1521
        /* t points to an endpoint which is less than top */
 
1522
        test = (top+bottom)/2;
 
1523
        for ( i=b; i<=t; ++i ) {
 
1524
            if ( RealNearish(test,ends[which][i]) ) {
 
1525
                if ( i==b )
 
1526
                    test = (bottom+ends[which][i])/2;
 
1527
                else
 
1528
                    test = (ends[which][i-1]+ends[which][i])/2;
 
1529
        break;
 
1530
            }
 
1531
        }
 
1532
        FigureNeeds(ms,which,test,space,ot,1);
 
1533
    }
 
1534
    free(ends[0]);
 
1535
    free(ends[1]);
 
1536
    free(space);
 
1537
    free(gaps);
 
1538
}
 
1539
 
 
1540
static void FindUnitVectors(Intersection *ilist) {
 
1541
    MList *ml;
 
1542
    Intersection *il;
 
1543
    BasePoint u;
 
1544
    double len;
 
1545
 
 
1546
    for ( il=ilist; il!=NULL; il=il->next ) {
 
1547
        for ( ml=il->monos; ml!=NULL; ml=ml->next ) {
 
1548
            if ( ml->m->isneeded ) {
 
1549
                Spline *s = ml->m->s;
 
1550
                double t1, t2;
 
1551
                t1 = ml->t;
 
1552
                if ( ml->isend )
 
1553
                    t2 = ml->t - (ml->t-ml->m->tstart)/20.0;
 
1554
                else
 
1555
                    t2 = ml->t + (ml->m->tend-ml->t)/20.0;
 
1556
                u.x = ((s->splines[0].a*t1 + s->splines[0].b)*t1 + s->splines[0].c)*t1 -
 
1557
                        ((s->splines[0].a*t2 + s->splines[0].b)*t2 + s->splines[0].c)*t2;
 
1558
                u.y = ((s->splines[1].a*t1 + s->splines[1].b)*t1 + s->splines[1].c)*t1 -
 
1559
                        ((s->splines[1].a*t2 + s->splines[1].b)*t2 + s->splines[1].c)*t2;
 
1560
                len = u.x*u.x + u.y*u.y;
 
1561
                if ( len!=0 ) {
 
1562
                    len = sqrt(len);
 
1563
                    u.x /= len;
 
1564
                    u.y /= len;
 
1565
                }
 
1566
                ml->unit = u;
 
1567
            }
 
1568
        }
 
1569
    }
 
1570
}
 
1571
 
 
1572
static void TestForBadDirections(Intersection *ilist) {
 
1573
    /* If we have a glyph with at least two contours one drawn clockwise, */
 
1574
    /*  one counter, and these two intersect, then our algorithm will */
 
1575
    /*  not remove what appears to the user to be an overlap. Warn about */
 
1576
    /*  this. */
 
1577
    /* I think it happens iff all exits from an intersection are needed */
 
1578
    MList *ml, *ml2;
 
1579
    int cnt, ncnt;
 
1580
    Intersection *il;
 
1581
 
 
1582
    /* If we have two splines one going from a->b and the other from b->a */
 
1583
    /*  tracing exactly the same route, then they should cancel each other */
 
1584
    /*  out. But depending on the order we hit them they may both be marked */
 
1585
    /*  needed */  /* OverlapBugs.sfd: asciicircumflex */
 
1586
    for ( il=ilist; il!=NULL; il=il->next ) {
 
1587
        for ( ml=il->monos; ml!=NULL; ml=ml->next ) {
 
1588
            if ( ml->m->isneeded && ml->m->s->knownlinear &&
 
1589
                    ml->m->start!=NULL && ml->m->end!=NULL ) {
 
1590
                for ( ml2 = ml->next; ml2!=NULL; ml2=ml2->next ) {
 
1591
                    if ( ml2->m->isneeded && ml2->m->s->knownlinear &&
 
1592
                            ml2->m->start == ml->m->end &&
 
1593
                            ml2->m->end == ml->m->start ) {
 
1594
                        ml2->m->isneeded = false;
 
1595
                        ml->m->isneeded = false;
 
1596
                        ml2->m->isunneeded = true;
 
1597
                        ml->m->isunneeded = true;
 
1598
                break;
 
1599
                    }
 
1600
                }
 
1601
            }
 
1602
        }
 
1603
    }
 
1604
 
 
1605
    while ( ilist!=NULL ) {
 
1606
        cnt = ncnt = 0;
 
1607
        for ( ml = ilist->monos; ml!=NULL; ml=ml->next ) {
 
1608
            ++cnt;
 
1609
            if ( ml->m->isneeded ) ++ncnt;
 
1610
        }
 
1611
#if 0
 
1612
        if ( cnt>=4 && ncnt==cnt ) {
 
1613
            ff_post_notice(_("Warning"),_("Glyph %.40s contains an overlapped region where two contours with opposite orientations intersect. This will not be removed. In many cases doing Element->Correct Direction before Remove Overlap will improve matters."),
 
1614
                    glyphname );
 
1615
            fprintf( stderr, "%s: Fully needed intersection at (%g,%g)\n",
 
1616
                    glyphname,
 
1617
                    ilist->inter.x, ilist->inter.y );
 
1618
    break;
 
1619
        }
 
1620
#endif
 
1621
        ilist = ilist->next;
 
1622
    }
 
1623
}
 
1624
 
 
1625
static void MonoFigure(Spline *s,extended firstt,extended endt, SplinePoint *first,
 
1626
        SplinePoint *end) {
 
1627
    extended f;
 
1628
    Spline1D temp;
 
1629
 
 
1630
    f = endt - firstt;
 
1631
    /*temp.d = first->me.x;*/
 
1632
    /*temp.a = s->splines[0].a*f*f*f;*/
 
1633
    temp.b = (s->splines[0].b + 3*s->splines[0].a*firstt) *f*f;
 
1634
    temp.c = (s->splines[0].c + 2*s->splines[0].b*firstt + 3*s->splines[0].a*firstt*firstt) * f;
 
1635
    first->nextcp.x = first->me.x + temp.c/3;
 
1636
    end->prevcp.x = first->nextcp.x + (temp.b+temp.c)/3;
 
1637
    if ( temp.c>-.01 && temp.c<.01 ) first->nextcp.x = first->me.x;
 
1638
    if ( (temp.b+temp.c)>-.01 && (temp.b+temp.c)<.01 ) end->prevcp.x = end->me.x;
 
1639
 
 
1640
    temp.b = (s->splines[1].b + 3*s->splines[1].a*firstt) *f*f;
 
1641
    temp.c = (s->splines[1].c + 2*s->splines[1].b*firstt + 3*s->splines[1].a*firstt*firstt) * f;
 
1642
    first->nextcp.y = first->me.y + temp.c/3;
 
1643
    end->prevcp.y = first->nextcp.y + (temp.b+temp.c)/3;
 
1644
    if ( temp.c>-.01 && temp.c<.01 ) first->nextcp.y = first->me.y;
 
1645
    if ( (temp.b+temp.c)>-.01 && (temp.b+temp.c)<.01 ) end->prevcp.y = end->me.y;
 
1646
    first->nonextcp = false; end->noprevcp = false;
 
1647
    SplineMake3(first,end);
 
1648
    if ( SplineIsLinear(first->next)) {
 
1649
        first->nextcp = first->me;
 
1650
        end->prevcp = end->me;
 
1651
        first->nonextcp = end->noprevcp = true;
 
1652
        SplineRefigure(first->next);
 
1653
    }
 
1654
}
 
1655
 
 
1656
static Intersection *MonoFollow(Intersection *curil, Monotonic *m) {
 
1657
    Monotonic *mstart=m;
 
1658
 
 
1659
    if ( m->start==curil ) {
 
1660
        while ( m!=NULL && m->end==NULL ) {
 
1661
            m=m->next;
 
1662
            if ( m==mstart )
 
1663
        break;
 
1664
        }
 
1665
        if ( m==NULL )
 
1666
return( NULL );
 
1667
 
 
1668
return( m->end );
 
1669
    } else {
 
1670
        while ( m!=NULL && m->start==NULL ) {
 
1671
            m=m->prev;
 
1672
            if ( m==mstart )
 
1673
        break;
 
1674
        }
 
1675
        if ( m==NULL )
 
1676
return( NULL );
 
1677
 
 
1678
return( m->start );
 
1679
    }
 
1680
}
 
1681
 
 
1682
static int MonoGoesSomewhereUseful(Intersection *curil, Monotonic *m) {
 
1683
    Intersection *nextil = MonoFollow(curil,m);
 
1684
    MList *ml;
 
1685
    int cnt;
 
1686
 
 
1687
    if ( nextil==NULL )
 
1688
return( false );
 
1689
    cnt = 0;
 
1690
    for ( ml=nextil->monos; ml!=NULL ; ml=ml->next )
 
1691
        if ( ml->m->isneeded )
 
1692
            ++cnt;
 
1693
    if ( cnt>=2 )       /* One for the mono that one in, one for another going out... */
 
1694
return( true );
 
1695
 
 
1696
return( false );
 
1697
}
 
1698
 
 
1699
static MList *FindMLOfM(Intersection *curil,Monotonic *finalm) {
 
1700
    MList *ml;
 
1701
 
 
1702
    for ( ml=curil->monos; ml!=NULL; ml=ml->next ) {
 
1703
        if ( ml->m==finalm )
 
1704
return( ml );
 
1705
    }
 
1706
#if 0
 
1707
        if ( ml->isend ) {
 
1708
            for ( m=ml->m; m!=NULL ; m=m->prev ) {
 
1709
                if ( m==finalm )
 
1710
return( ml );
 
1711
                if ( m->start != NULL )
 
1712
            break;
 
1713
            }
 
1714
        } else {
 
1715
            for ( m=ml->m; m!=NULL; m=m->next ) {
 
1716
                if ( m==finalm )
 
1717
return( ml );
 
1718
                if ( m->end!=NULL )
 
1719
            break;
 
1720
            }
 
1721
        }
 
1722
    }
 
1723
#endif
 
1724
return( NULL );
 
1725
}
 
1726
 
 
1727
static SplinePoint *MonoFollowForward(Intersection **curil, MList *ml,
 
1728
        SplinePoint *last, Monotonic **finalm) {
 
1729
    SplinePoint *mid;
 
1730
    Monotonic *m = ml->m, *mstart;
 
1731
 
 
1732
    forever {
 
1733
        for ( mstart = m; m->s==mstart->s; m=m->next) {
 
1734
            if ( !m->isneeded )
 
1735
                SOError( "Expected needed monotonic.\n" );
 
1736
            m->isneeded = false;                /* Mark as used */
 
1737
            if ( m->end!=NULL )
 
1738
        break;
 
1739
        }
 
1740
        if ( m->s==mstart->s ) {
 
1741
            if ( m->end==NULL ) SOError( "Invariant condition does not hold.\n" );
 
1742
            mid = SplinePointCreate(m->end->inter.x,m->end->inter.y);
 
1743
        } else {
 
1744
            m = m->prev;
 
1745
            mid = SplinePointCreate(m->s->to->me.x,m->s->to->me.y);
 
1746
        }
 
1747
        if ( mstart->tstart==0 && m->tend==1.0 ) {
 
1748
            /* I check for this special case to avoid rounding errors */
 
1749
            last->nextcp = m->s->from->nextcp;
 
1750
            last->nonextcp = m->s->from->nonextcp;
 
1751
            mid->prevcp = m->s->to->prevcp;
 
1752
            mid->noprevcp = m->s->to->noprevcp;
 
1753
            SplineMake3(last,mid);
 
1754
        } else {
 
1755
            MonoFigure(m->s,mstart->tstart,m->tend,last,mid);
 
1756
        }
 
1757
        last = mid;
 
1758
        if ( m->end!=NULL ) {
 
1759
            *curil = m->end;
 
1760
            *finalm = m;
 
1761
return( last );
 
1762
        }
 
1763
        m = m->next;
 
1764
    }
 
1765
}
 
1766
 
 
1767
static SplinePoint *MonoFollowBackward(Intersection **curil, MList *ml,
 
1768
        SplinePoint *last, Monotonic **finalm) {
 
1769
    SplinePoint *mid;
 
1770
    Monotonic *m = ml->m, *mstart;
 
1771
 
 
1772
    forever {
 
1773
        for ( mstart=m; m->s==mstart->s; m=m->prev) {
 
1774
            if ( !m->isneeded )
 
1775
                SOError( "Expected needed monotonic.\n" );
 
1776
            m->isneeded = false;                /* Mark as used */
 
1777
            if ( m->start!=NULL )
 
1778
        break;
 
1779
        }
 
1780
        if ( m->s==mstart->s ) {
 
1781
            if ( m->start==NULL ) SOError( "Invariant condition does not hold.\n" );
 
1782
            mid = SplinePointCreate(m->start->inter.x,m->start->inter.y);
 
1783
        } else {
 
1784
            m = m->next;
 
1785
            mid = SplinePointCreate(m->s->from->me.x,m->s->from->me.y);
 
1786
        }
 
1787
        if ( m->s->knownlinear ) mid->pointtype = pt_corner;
 
1788
        if ( mstart->tend==1.0 && m->tstart==0 ) {
 
1789
            /* I check for this special case to avoid rounding errors */
 
1790
            last->nextcp = m->s->to->prevcp;
 
1791
            last->nonextcp = m->s->to->noprevcp;
 
1792
            mid->prevcp = m->s->from->nextcp;
 
1793
            mid->noprevcp = m->s->from->nonextcp;
 
1794
            SplineMake3(last,mid);
 
1795
        } else {
 
1796
            MonoFigure(m->s,mstart->tend,m->tstart,last,mid);
 
1797
        }
 
1798
        last = mid;
 
1799
        if ( m->start!=NULL ) {
 
1800
            *curil = m->start;
 
1801
            *finalm = m;
 
1802
return( last );
 
1803
        }
 
1804
        m = m->prev;
 
1805
    }
 
1806
}
 
1807
 
 
1808
static SplineSet *JoinAContour(Intersection *startil,MList *ml) {
 
1809
    SplineSet *ss = chunkalloc(sizeof(SplineSet));
 
1810
    SplinePoint *last;
 
1811
    Intersection *curil;
 
1812
    int allexclude = ml->m->exclude;
 
1813
    Monotonic *finalm;
 
1814
    MList *lastml;
 
1815
 
 
1816
    ss->first = last = SplinePointCreate(startil->inter.x,startil->inter.y);
 
1817
    curil = startil;
 
1818
    forever {
 
1819
        if ( allexclude && !ml->m->exclude ) allexclude = false;
 
1820
        finalm = NULL;
 
1821
        if ( ml->m->start==curil ) {
 
1822
            last = MonoFollowForward(&curil,ml,last,&finalm);
 
1823
        } else if ( ml->m->end==curil ) {
 
1824
            last = MonoFollowBackward(&curil,ml,last,&finalm);
 
1825
        } else {
 
1826
            SOError( "Couldn't find endpoint (%g,%g).\n",
 
1827
                    curil->inter.x, curil->inter.y );
 
1828
            ml->m->isneeded = false;            /* Prevent infinite loops */
 
1829
            ss->last = last;
 
1830
    break;
 
1831
        }
 
1832
        if ( curil==startil ) {
 
1833
            ss->first->prev = last->prev;
 
1834
            ss->first->prevcp = last->prevcp;
 
1835
            ss->first->noprevcp = last->noprevcp;
 
1836
            last->prev->to = ss->first;
 
1837
            SplinePointFree(last);
 
1838
            ss->last = ss->first;
 
1839
    break;
 
1840
        }
 
1841
        lastml = FindMLOfM(curil,finalm);
 
1842
        if ( lastml==NULL ) {
 
1843
            IError("Could not find finalm");
 
1844
            /* Try to preserve direction */
 
1845
            for ( ml=curil->monos; ml!=NULL && (!ml->m->isneeded || ml->m->end==curil); ml=ml->next );
 
1846
            if ( ml==NULL )
 
1847
                for ( ml=curil->monos; ml!=NULL && !ml->m->isneeded; ml=ml->next );
 
1848
        } else {
 
1849
            int k; MList *bestml; double bestdot;
 
1850
            for ( k=0; k<2; ++k ) {
 
1851
                bestml = NULL; bestdot = -2;
 
1852
                for ( ml=curil->monos; ml!=NULL ; ml=ml->next ) {
 
1853
                    if ( ml->m->isneeded && (ml->m->start==curil || k) ) {
 
1854
                        double dot = lastml->unit.x*ml->unit.x + lastml->unit.y*ml->unit.y;
 
1855
                        if ( dot>bestdot ) {
 
1856
                            bestml = ml;
 
1857
                            bestdot = dot;
 
1858
                        }
 
1859
                    }
 
1860
                }
 
1861
                if ( bestml!=NULL )
 
1862
            break;
 
1863
            }
 
1864
            ml = bestml;
 
1865
        }
 
1866
        if ( ml==NULL ) {
 
1867
            for ( ml=curil->monos; ml!=NULL ; ml=ml->next )
 
1868
                if ( ml->m->isunneeded && ml->m->start==curil &&
 
1869
                        MonoFollow(curil,ml->m)==startil )
 
1870
            break;
 
1871
            if ( ml==NULL )
 
1872
                for ( ml=curil->monos; ml!=NULL ; ml=ml->next )
 
1873
                    if ( ml->m->isunneeded && ml->m->end==curil &&
 
1874
                            MonoFollow(curil,ml->m)==startil )
 
1875
                break;
 
1876
            if ( ml!=NULL ) {
 
1877
                SOError("Closing contour with unneeded path\n" );
 
1878
                ml->m->isneeded = true;
 
1879
            }
 
1880
        }
 
1881
        if ( ml==NULL ) {
 
1882
            SOError( "couldn't find a needed exit from an intersection\n" );
 
1883
            ss->last = last;
 
1884
    break;
 
1885
        }
 
1886
    }
 
1887
    SPLCatagorizePoints(ss);
 
1888
    if ( allexclude && SplinePointListIsClockwise(ss))
 
1889
        SplineSetReverse(ss);
 
1890
return( ss );
 
1891
}
 
1892
 
 
1893
static SplineSet *FindMatchingContour(SplineSet *head,SplineSet *cur) {
 
1894
    SplineSet *test;
 
1895
 
 
1896
    for ( test=head; test!=NULL; test=test->next ) {
 
1897
        if ( test->first->prev==NULL &&
 
1898
                test->first->me.x==cur->last->me.x && test->first->me.y==cur->last->me.y &&
 
1899
                test->last->me.x==cur->first->me.x && test->last->me.y==cur->first->me.y )
 
1900
    break;
 
1901
    }
 
1902
    if ( test==NULL ) {
 
1903
        for ( test=head; test!=NULL; test=test->next ) {
 
1904
            if ( test->first->prev==NULL &&
 
1905
                    test->last->me.x==cur->last->me.x && test->last->me.y==cur->last->me.y &&
 
1906
                    test->first->me.x==cur->first->me.x && test->first->me.y==cur->first->me.y ) {
 
1907
                SplineSetReverse(cur);
 
1908
        break;
 
1909
            }
 
1910
        }
 
1911
    }
 
1912
    if ( test==NULL ) {
 
1913
        for ( test=head; test!=NULL; test=test->next ) {
 
1914
            if ( test->first->prev==NULL &&
 
1915
                    ((test->first->me.x==cur->last->me.x && test->first->me.y==cur->last->me.y) ||
 
1916
                     (test->last->me.x==cur->first->me.x && test->last->me.y==cur->first->me.y )))
 
1917
        break;
 
1918
        }
 
1919
    }
 
1920
    if ( test==NULL ) {
 
1921
        for ( test=head; test!=NULL; test=test->next ) {
 
1922
            if ( test->first->prev==NULL &&
 
1923
                    ((test->last->me.x==cur->last->me.x && test->last->me.y==cur->last->me.y) ||
 
1924
                     (test->first->me.x==cur->first->me.x && test->first->me.y==cur->first->me.y ))) {
 
1925
                SplineSetReverse(cur);
 
1926
        break;
 
1927
            }
 
1928
        }
 
1929
    }
 
1930
return( test );
 
1931
}
 
1932
 
 
1933
static SplineSet *JoinAllNeeded(Intersection *ilist) {
 
1934
    Intersection *il;
 
1935
    SplineSet *head=NULL, *last=NULL, *cur, *test;
 
1936
    MList *ml;
 
1937
 
 
1938
    for ( il=ilist; il!=NULL; il=il->next ) {
 
1939
        /* Try to preserve direction */
 
1940
        forever {
 
1941
            for ( ml=il->monos; ml!=NULL && (!ml->m->isneeded || ml->m->end==il); ml=ml->next );
 
1942
            if ( ml==NULL )
 
1943
                for ( ml=il->monos; ml!=NULL && !ml->m->isneeded; ml=ml->next );
 
1944
            if ( ml==NULL )
 
1945
        break;
 
1946
            if ( !MonoGoesSomewhereUseful(il,ml->m)) {
 
1947
                SOError("Humph. This monotonic leads nowhere.\n" );
 
1948
        /* break; */
 
1949
            }
 
1950
            cur = JoinAContour(il,ml);
 
1951
            if ( head==NULL )
 
1952
                head = cur;
 
1953
            else {
 
1954
                if ( cur->first->prev==NULL ) {
 
1955
                    /* Open contours are errors. See if we had an earlier error */
 
1956
                    /*  to which we can join this */
 
1957
                    test = FindMatchingContour(head,cur);
 
1958
                    if ( test!=NULL ) {
 
1959
                        if ( test->first->me.x==cur->last->me.x && test->first->me.y==cur->last->me.y ) {
 
1960
                            test->first->prev = cur->last->prev;
 
1961
                            cur->last->prev->to = test->first;
 
1962
                            SplinePointFree(cur->last);
 
1963
                            if ( test->last->me.x==cur->first->me.x && test->last->me.y==cur->first->me.y ) {
 
1964
                                test->last->next = cur->first->next;
 
1965
                                cur->first->next->from = test->last;
 
1966
                                SplinePointFree(cur->first);
 
1967
                                test->last = test->first;
 
1968
                            } else
 
1969
                                test->first = cur->first;
 
1970
                        } else {
 
1971
                            if ( test->last->me.x!=cur->first->me.x || test->last->me.y!=cur->first->me.y )
 
1972
                                SOError( "Join failed");
 
1973
                            else {
 
1974
                                test->last->next = cur->first->next;
 
1975
                                cur->first->next->from = test->last;
 
1976
                                SplinePointFree(cur->first);
 
1977
                                test->last = test->first;
 
1978
                            }
 
1979
                        }
 
1980
                        cur->first = cur->last = NULL;
 
1981
                        SplinePointListFree(cur);
 
1982
                        cur=NULL;
 
1983
                    }
 
1984
                }
 
1985
                if ( cur!=NULL )
 
1986
                    last->next = cur;
 
1987
            }
 
1988
            if ( cur!=NULL )
 
1989
                last = cur;
 
1990
        }
 
1991
    }
 
1992
return( head );
 
1993
}
 
1994
 
 
1995
static SplineSet *MergeOpenAndFreeClosed(SplineSet *new,SplineSet *old,
 
1996
        enum overlap_type ot) {
 
1997
    SplineSet *next;
 
1998
 
 
1999
    while ( old!=NULL ) {
 
2000
        next = old->next;
 
2001
        if ( old->first->prev==NULL ||
 
2002
                (( ot==over_rmselected || ot==over_intersel || ot==over_fisel) &&
 
2003
                  !SSIsSelected(old)) ) {
 
2004
            old->next = new;
 
2005
            new = old;
 
2006
        } else {
 
2007
            old->next = NULL;
 
2008
            SplinePointListFree(old);
 
2009
        }
 
2010
        old = next;
 
2011
    }
 
2012
return(new);
 
2013
}
 
2014
 
 
2015
void FreeMonotonics(Monotonic *m) {
 
2016
    Monotonic *next;
 
2017
 
 
2018
    while ( m!=NULL ) {
 
2019
        next = m->linked;
 
2020
        chunkfree(m,sizeof(*m));
 
2021
        m = next;
 
2022
    }
 
2023
}
 
2024
 
 
2025
static void FreeMList(MList *ml) {
 
2026
    MList *next;
 
2027
 
 
2028
    while ( ml!=NULL ) {
 
2029
        next = ml->next;
 
2030
        chunkfree(ml,sizeof(*ml));
 
2031
        ml = next;
 
2032
    }
 
2033
}
 
2034
 
 
2035
static void FreeIntersections(Intersection *ilist) {
 
2036
    Intersection *next;
 
2037
 
 
2038
    while ( ilist!=NULL ) {
 
2039
        next = ilist->next;
 
2040
        FreeMList(ilist->monos);
 
2041
        chunkfree(ilist,sizeof(*ilist));
 
2042
        ilist = next;
 
2043
    }
 
2044
}
 
2045
 
 
2046
static void MonoSplit(Monotonic *m) {
 
2047
    Spline *s = m->s;
 
2048
    SplinePoint *last = s->from;
 
2049
    SplinePoint *final = s->to;
 
2050
    extended lastt = 0;
 
2051
 
 
2052
    last->next = NULL;
 
2053
    final->prev = NULL;
 
2054
    while ( m!=NULL && m->s==s && m->tend<1 ) {
 
2055
        if ( m->end!=NULL ) {
 
2056
            SplinePoint *mid = SplinePointCreate(m->end->inter.x,m->end->inter.y);
 
2057
            if ( m->s->knownlinear ) mid->pointtype = pt_corner;
 
2058
            MonoFigure(s,lastt,m->tend,last,mid);
 
2059
            lastt = m->tend;
 
2060
            last = mid;
 
2061
        }
 
2062
        m = m->linked;
 
2063
    }
 
2064
    MonoFigure(s,lastt,1.0,last,final);
 
2065
    SplineFree(s);
 
2066
}
 
2067
 
 
2068
static void FixupIntersectedSplines(Monotonic *ms) {
 
2069
    /* If all we want is intersections, then the contours are already correct */
 
2070
    /*  all we need to do is run through the Monotonic list and when we find */
 
2071
    /*  an intersection, make sure it has real splines around it */
 
2072
    Monotonic *m;
 
2073
    int needs_split;
 
2074
 
 
2075
    while ( ms!=NULL ) {
 
2076
        needs_split = false;
 
2077
        for ( m=ms; m!=NULL && m->s==ms->s; m=m->linked ) {
 
2078
            if ( (m->tstart!=0 && m->start!=NULL) || (m->tend!=1 && m->end!=NULL))
 
2079
                needs_split = true;
 
2080
        }
 
2081
        if ( needs_split )
 
2082
            MonoSplit(ms);
 
2083
        ms = m;
 
2084
    }
 
2085
}
 
2086
 
 
2087
static int BpClose(BasePoint *here, BasePoint *there, double error) {
 
2088
    extended dx, dy;
 
2089
 
 
2090
    if ( (dx = here->x-there->x)<0 ) dx= -dx;
 
2091
    if ( (dy = here->y-there->y)<0 ) dy= -dy;
 
2092
    if ( dx<error && dy<error )
 
2093
return( true );
 
2094
 
 
2095
return( false );
 
2096
}
 
2097
 
 
2098
static SplineSet *SSRemoveTiny(SplineSet *base) {
 
2099
    DBounds b;
 
2100
    double error;
 
2101
    extended test, dx, dy;
 
2102
    SplineSet *prev = NULL, *head = base, *ssnext;
 
2103
    SplinePoint *sp, *nsp;
 
2104
 
 
2105
    SplineSetQuickBounds(base,&b);
 
2106
    error = b.maxy-b.miny;
 
2107
    test = b.maxx-b.minx;
 
2108
    if ( test>error ) error = test;
 
2109
    if ( (test = b.maxy)<0 ) test = -test;
 
2110
    if ( test>error ) error = test;
 
2111
    if ( (test = b.maxx)<0 ) test = -test;
 
2112
    if ( test>error ) error = test;
 
2113
    error /= 30000;
 
2114
 
 
2115
    while ( base!=NULL ) {
 
2116
        ssnext = base->next;
 
2117
        for ( sp=base->first; ; ) {
 
2118
            if ( sp->next==NULL )
 
2119
        break;
 
2120
            nsp = sp->next->to;
 
2121
            if ( BpClose(&sp->me,&nsp->me,error) ) {
 
2122
                if ( BpClose(&sp->me,&sp->nextcp,2*error) &&
 
2123
                        BpClose(&nsp->me,&nsp->prevcp,2*error)) {
 
2124
                    /* Remove the spline */
 
2125
                    if ( nsp==sp ) {
 
2126
                        /* Only this spline in the contour, so remove the contour */
 
2127
                        base->next = NULL;
 
2128
                        SplinePointListFree(base);
 
2129
                        if ( prev==NULL )
 
2130
                            head = ssnext;
 
2131
                        else
 
2132
                            prev->next = ssnext;
 
2133
                        base = NULL;
 
2134
        break;
 
2135
                    }
 
2136
                    SplineFree(sp->next);
 
2137
                    if ( nsp->nonextcp ) {
 
2138
                        sp->nextcp = sp->me;
 
2139
                        sp->nonextcp = true;
 
2140
                    } else {
 
2141
                        sp->nextcp = nsp->nextcp;
 
2142
                        sp->nonextcp = false;
 
2143
                    }
 
2144
                    sp->nextcpdef = nsp->nextcpdef;
 
2145
                    sp->next = nsp->next;
 
2146
                    if ( nsp->next!=NULL ) {
 
2147
                        nsp->next->from = sp;
 
2148
                        SplineRefigure(sp->next);
 
2149
                    }
 
2150
                    if ( nsp==base->last )
 
2151
                        base->last = sp;
 
2152
                    if ( nsp==base->first )
 
2153
                        base->first = sp;
 
2154
                    SplinePointFree(nsp);
 
2155
                    if ( sp->next==NULL )
 
2156
        break;
 
2157
                    nsp = sp->next->to;
 
2158
                } else {
 
2159
                    /* Leave the spline, but move the two points together */
 
2160
                    BasePoint new;
 
2161
                    new.x = (sp->me.x+nsp->me.x)/2;
 
2162
                    new.y = (sp->me.y+nsp->me.y)/2;
 
2163
                    dx = new.x-sp->me.x; dy = new.y-sp->me.y;
 
2164
                    sp->me = new;
 
2165
                    sp->nextcp.x += dx; sp->nextcp.y += dy;
 
2166
                    sp->prevcp.x += dx; sp->prevcp.y += dy;
 
2167
                    dx = new.x-nsp->me.x; dy = new.y-nsp->me.y;
 
2168
                    nsp->me = new;
 
2169
                    nsp->nextcp.x += dx; nsp->nextcp.y += dy;
 
2170
                    nsp->prevcp.x += dx; nsp->prevcp.y += dy;
 
2171
                    SplineRefigure(sp->next);
 
2172
                    if ( sp->prev ) SplineRefigure(sp->prev);
 
2173
                    if ( nsp->next ) SplineRefigure(nsp->next);
 
2174
                }
 
2175
            }
 
2176
            sp = nsp;
 
2177
            if ( sp==base->first )
 
2178
        break;
 
2179
        }
 
2180
        if ( sp->prev!=NULL && !sp->noprevcp ) {
 
2181
            int refigure = false;
 
2182
            if ( sp->me.x-sp->prevcp.x>-error && sp->me.x-sp->prevcp.x<error ) {
 
2183
                sp->prevcp.x = sp->me.x;
 
2184
                refigure = true;
 
2185
            }
 
2186
            if ( sp->me.y-sp->prevcp.y>-error && sp->me.y-sp->prevcp.y<error ) {
 
2187
                sp->prevcp.y = sp->me.y;
 
2188
                refigure = true;
 
2189
            }
 
2190
            if ( sp->me.x==sp->prevcp.x && sp->me.y==sp->prevcp.y )
 
2191
                sp->noprevcp = true;
 
2192
            if ( refigure )
 
2193
                SplineRefigure(sp->prev);
 
2194
        }
 
2195
        if ( sp->next!=NULL && !sp->nonextcp ) {
 
2196
            int refigure = false;
 
2197
            if ( sp->me.x-sp->nextcp.x>-error && sp->me.x-sp->nextcp.x<error ) {
 
2198
                sp->nextcp.x = sp->me.x;
 
2199
                refigure = true;
 
2200
            }
 
2201
            if ( sp->me.y-sp->nextcp.y>-error && sp->me.y-sp->nextcp.y<error ) {
 
2202
                sp->nextcp.y = sp->me.y;
 
2203
                refigure = true;
 
2204
            }
 
2205
            if ( sp->me.x==sp->nextcp.x && sp->me.y==sp->nextcp.y )
 
2206
                sp->nonextcp = true;
 
2207
            if ( refigure )
 
2208
                SplineRefigure(sp->next);
 
2209
        }
 
2210
        if ( base!=NULL )
 
2211
            prev = base;
 
2212
        base = ssnext;
 
2213
    }
 
2214
 
 
2215
return( head );
 
2216
}
 
2217
 
 
2218
static void RemoveNextSP(SplinePoint *psp,SplinePoint *sp,SplinePoint *nsp,
 
2219
        SplineSet *base) {
 
2220
    if ( psp==nsp ) {
 
2221
        SplineFree(psp->next);
 
2222
        psp->next = psp->prev;
 
2223
        psp->next->from = psp;
 
2224
        SplinePointFree(sp);
 
2225
        SplineRefigure(psp->prev);
 
2226
    } else {
 
2227
        psp->next = nsp->next;
 
2228
        psp->next->from = psp;
 
2229
        psp->nextcp = nsp->nextcp;
 
2230
        psp->nonextcp = nsp->nonextcp;
 
2231
        psp->nextcpdef = nsp->nextcpdef;
 
2232
        SplineFree(sp->prev);
 
2233
        SplineFree(sp->next);
 
2234
        SplinePointFree(sp);
 
2235
        SplinePointFree(nsp);
 
2236
        SplineRefigure(psp->next);
 
2237
    }
 
2238
    if ( base->first==sp || base->first==nsp )
 
2239
        base->first = psp;
 
2240
    if ( base->last==sp || base->last==nsp )
 
2241
        base->last = psp;
 
2242
}
 
2243
 
 
2244
static void RemovePrevSP(SplinePoint *psp,SplinePoint *sp,SplinePoint *nsp,
 
2245
        SplineSet *base) {
 
2246
    if ( psp==nsp ) {
 
2247
        SplineFree(nsp->prev);
 
2248
        nsp->prev = nsp->next;
 
2249
        nsp->prev->to = nsp;
 
2250
        SplinePointFree(sp);
 
2251
        SplineRefigure(nsp->next);
 
2252
    } else {
 
2253
        nsp->prev = psp->prev;
 
2254
        nsp->prev->to = nsp;
 
2255
        nsp->prevcp = nsp->me;
 
2256
        nsp->noprevcp = true;
 
2257
        nsp->prevcpdef = psp->prevcpdef;
 
2258
        SplineFree(sp->prev);
 
2259
        SplineFree(sp->next);
 
2260
        SplinePointFree(sp);
 
2261
        SplinePointFree(psp);
 
2262
        SplineRefigure(nsp->prev);
 
2263
    }
 
2264
    if ( base->first==sp || base->first==psp )
 
2265
        base->first = nsp;
 
2266
    if ( base->last==sp || base->last==psp )
 
2267
        base->last = nsp;
 
2268
}
 
2269
 
 
2270
static SplinePoint *SameLine(SplinePoint *psp, SplinePoint *sp, SplinePoint *nsp) {
 
2271
    BasePoint noff, poff;
 
2272
    real nlen, plen, normal;
 
2273
 
 
2274
    noff.x = nsp->me.x-sp->me.x; noff.y = nsp->me.y-sp->me.y;
 
2275
    poff.x = psp->me.x-sp->me.x; poff.y = psp->me.y-sp->me.y;
 
2276
    nlen = esqrt(noff.x*noff.x + noff.y*noff.y);
 
2277
    plen = esqrt(poff.x*poff.x + poff.y*poff.y);
 
2278
    if ( nlen==0 )
 
2279
return( nsp );
 
2280
    if ( plen==0 )
 
2281
return( psp );
 
2282
    normal = (noff.x*poff.y - noff.y*poff.x)/nlen/plen;
 
2283
    if ( normal<-.0001 || normal>.0001 )
 
2284
return( NULL );
 
2285
 
 
2286
    if ( noff.x*poff.x < 0 || noff.y*poff.y < 0 )
 
2287
return( NULL );         /* Same line, but different directions */
 
2288
    if ( (noff.x>0 && noff.x>poff.x) ||
 
2289
            (noff.x<0 && noff.x<poff.x) ||
 
2290
            (noff.y>0 && noff.y>poff.y) ||
 
2291
            (noff.y<0 && noff.y<poff.y))
 
2292
return( nsp );
 
2293
 
 
2294
return( psp );
 
2295
}
 
2296
 
 
2297
static double AdjacentSplinesMatch(Spline *s1,Spline *s2,int s2forward) {
 
2298
    /* Is every point on s2 close to a point on s1 */
 
2299
    double t, tdiff, t1 = -1;
 
2300
    double xoff, yoff;
 
2301
    double t1start, t1end;
 
2302
    extended ts[2];
 
2303
    int i;
 
2304
 
 
2305
    if ( (xoff = s2->to->me.x-s2->from->me.x)<0 ) xoff = -xoff;
 
2306
    if ( (yoff = s2->to->me.y-s2->from->me.y)<0 ) yoff = -yoff;
 
2307
    if ( xoff>yoff )
 
2308
        SplineFindExtrema(&s1->splines[0],&ts[0],&ts[1]);
 
2309
    else
 
2310
        SplineFindExtrema(&s1->splines[1],&ts[0],&ts[1]);
 
2311
    if ( s2forward ) {
 
2312
        t = 0;
 
2313
        tdiff = 1/16.0;
 
2314
        t1end = 1;
 
2315
        for ( i=1; i>=0 && ts[i]==-1; --i );
 
2316
        t1start = i<0 ? 0 : ts[i];
 
2317
    } else {
 
2318
        t = 1;
 
2319
        tdiff = -1/16.0;
 
2320
        t1start = 0;
 
2321
        t1end = ( ts[0]==-1 ) ? 1.0 : ts[0];
 
2322
    }
 
2323
 
 
2324
    for ( ; (s2forward && t<=1) || (!s2forward && t>=0 ); t += tdiff ) {
 
2325
        double x1, y1, xo, yo;
 
2326
        double x = ((s2->splines[0].a*t+s2->splines[0].b)*t+s2->splines[0].c)*t+s2->splines[0].d;
 
2327
        double y = ((s2->splines[1].a*t+s2->splines[1].b)*t+s2->splines[1].c)*t+s2->splines[1].d;
 
2328
        if ( xoff>yoff )
 
2329
            t1 = IterateSplineSolve(&s1->splines[0],t1start,t1end,x,.001);
 
2330
        else
 
2331
            t1 = IterateSplineSolve(&s1->splines[1],t1start,t1end,y,.001);
 
2332
        if ( t1<0 || t1>1 )
 
2333
return( -1 );
 
2334
        x1 = ((s1->splines[0].a*t1+s1->splines[0].b)*t1+s1->splines[0].c)*t1+s1->splines[0].d;
 
2335
        y1 = ((s1->splines[1].a*t1+s1->splines[1].b)*t1+s1->splines[1].c)*t1+s1->splines[1].d;
 
2336
        if ( (xo = (x-x1))<0 ) xo = -xo;
 
2337
        if ( (yo = (y-y1))<0 ) yo = -yo;
 
2338
        if ( xo+yo>.5 )
 
2339
return( -1 );
 
2340
    }
 
2341
return( t1 );
 
2342
}
 
2343
        
 
2344
void SSRemoveBacktracks(SplineSet *ss) {
 
2345
    SplinePoint *sp;
 
2346
 
 
2347
    if ( ss==NULL )
 
2348
return;
 
2349
    for ( sp=ss->first; ; ) {
 
2350
        if ( sp->next!=NULL && sp->prev!=NULL ) {
 
2351
            SplinePoint *nsp = sp->next->to, *psp = sp->prev->from, *isp;
 
2352
            BasePoint ndir, pdir;
 
2353
            double dot, pdot, nlen, plen, t = -1;
 
2354
 
 
2355
            ndir.x = (nsp->me.x - sp->me.x); ndir.y = (nsp->me.y - sp->me.y);
 
2356
            pdir.x = (psp->me.x - sp->me.x); pdir.y = (psp->me.y - sp->me.y);
 
2357
            nlen = ndir.x*ndir.x + ndir.y*ndir.y; plen = pdir.x*pdir.x + pdir.y*pdir.y;
 
2358
            dot = ndir.x*pdir.x + ndir.y*pdir.y;
 
2359
            if ( (pdot = ndir.x*pdir.y - ndir.y*pdir.x)<0 ) pdot = -pdot;
 
2360
            if ( dot>0 && dot>pdot ) {
 
2361
                if ( nlen>plen && (t=AdjacentSplinesMatch(sp->next,sp->prev,false))!=-1 ) {
 
2362
                    isp = SplineBisect(sp->next,t);
 
2363
                    psp->nextcp.x = psp->me.x + (isp->nextcp.x-isp->me.x);
 
2364
                    psp->nextcp.y = psp->me.y + (isp->nextcp.y-isp->me.y);
 
2365
                    psp->nonextcp = isp->nonextcp;
 
2366
                    psp->next = isp->next;
 
2367
                    isp->next->from = psp;
 
2368
                    SplineFree(isp->prev);
 
2369
                    SplineFree(sp->prev);
 
2370
                    SplinePointFree(isp);
 
2371
                    SplinePointFree(sp);
 
2372
                    SplineRefigure(psp->next);
 
2373
                    if ( ss->first==sp )
 
2374
                        ss->first = psp;
 
2375
                    if ( ss->last==sp )
 
2376
                        ss->last = psp;
 
2377
                    sp=psp;
 
2378
                } else if ( nlen<plen && (t=AdjacentSplinesMatch(sp->prev,sp->next,true))!=-1 ) {
 
2379
                    isp = SplineBisect(sp->prev,t);
 
2380
                    nsp->prevcp.x = nsp->me.x + (isp->prevcp.x-isp->me.x);
 
2381
                    nsp->prevcp.y = nsp->me.y + (isp->prevcp.y-isp->me.y);
 
2382
                    nsp->noprevcp = isp->noprevcp;
 
2383
                    nsp->prev = isp->prev;
 
2384
                    isp->prev->to = nsp;
 
2385
                    SplineFree(isp->next);
 
2386
                    SplineFree(sp->next);
 
2387
                    SplinePointFree(isp);
 
2388
                    SplinePointFree(sp);
 
2389
                    SplineRefigure(nsp->prev);
 
2390
                    if ( ss->first==sp )
 
2391
                        ss->first = psp;
 
2392
                    if ( ss->last==sp )
 
2393
                        ss->last = psp;
 
2394
                    sp=psp;
 
2395
                }
 
2396
            }
 
2397
        }
 
2398
        if ( sp->next==NULL )
 
2399
    break;
 
2400
        sp=sp->next->to;
 
2401
        if ( sp==ss->first )
 
2402
    break;
 
2403
    }
 
2404
}
 
2405
 
 
2406
/* If we have a contour with no width, say a line from A to B and then from B */
 
2407
/*  to A, then it will be ambiguous, depending on how we hit the contour, as  */
 
2408
/*  to whether it is needed or not. Which will cause us to complain. Since    */
 
2409
/*  they contain no area, they achieve nothing, so we might as well say they  */
 
2410
/*  overlap themselves and remove them here */
 
2411
static SplineSet *SSRemoveReversals(SplineSet *base) {
 
2412
    SplineSet *head = base, *prev=NULL, *next;
 
2413
    SplinePoint *sp;
 
2414
    int changed;
 
2415
 
 
2416
    while ( base!=NULL ) {
 
2417
        next = base->next;
 
2418
        changed = true;
 
2419
        while ( changed ) {
 
2420
            changed = false;
 
2421
            if ( base->first->next==NULL ||
 
2422
                    (base->first->next->to==base->first &&
 
2423
                     base->first->nextcp.x==base->first->prevcp.x &&
 
2424
                     base->first->nextcp.y==base->first->prevcp.y)) {
 
2425
                /* remove single points */
 
2426
                if ( prev==NULL )
 
2427
                    head = next;
 
2428
                else
 
2429
                    prev->next = next;
 
2430
                base->next = NULL;
 
2431
                SplinePointListFree(base);
 
2432
                base = prev;
 
2433
        break;
 
2434
            }
 
2435
            for ( sp=base->first; ; ) {
 
2436
                if ( sp->next!=NULL && sp->prev!=NULL &&
 
2437
                        sp->nextcp.x==sp->prevcp.x && sp->nextcp.y==sp->prevcp.y ) {
 
2438
                    SplinePoint *nsp = sp->next->to, *psp = sp->prev->from, *isp;
 
2439
                    if ( psp->me.x==nsp->me.x && psp->me.y==nsp->me.y &&
 
2440
                            psp->nextcp.x==nsp->prevcp.x && psp->nextcp.y==nsp->prevcp.y ) {
 
2441
                        /* We wish to remove sp, sp->next, sp->prev & one of nsp/psp */
 
2442
                        RemoveNextSP(psp,sp,nsp,base);
 
2443
                        changed = true;
 
2444
            break;
 
2445
                    } else if ( sp->nonextcp /* which implies sp->noprevcp */ &&
 
2446
                            psp->nonextcp && nsp->noprevcp &&
 
2447
                            (isp = SameLine(psp,sp,nsp))!=NULL ) {
 
2448
                        /* We have a line that backtracks, but doesn't cover */
 
2449
                        /*  the entire spline, so we intersect */
 
2450
                        /* We want to remove sp, the shorter of sp->next, sp->prev */
 
2451
                        /*  and a bit of the other one. Also reomve one of nsp,psp */
 
2452
                        if ( isp==psp ) {
 
2453
                            RemoveNextSP(psp,sp,nsp,base);
 
2454
                            psp->nextcp = psp->me;
 
2455
                            psp->nonextcp = true;
 
2456
                        } else {
 
2457
                            RemovePrevSP(psp,sp,nsp,base);
 
2458
                        }
 
2459
                        changed = true;
 
2460
            break;
 
2461
                    }
 
2462
                }
 
2463
                if ( sp->next==NULL )
 
2464
            break;
 
2465
                sp = sp->next->to;
 
2466
                if ( sp==base->first )
 
2467
            break;
 
2468
            }
 
2469
        }
 
2470
        SSRemoveBacktracks(base);
 
2471
        prev = base;
 
2472
        base = next;
 
2473
    }
 
2474
return( head );
 
2475
}
 
2476
 
 
2477
SplineSet *SplineSetRemoveOverlap(SplineChar *sc, SplineSet *base,enum overlap_type ot) {
 
2478
    Monotonic *ms;
 
2479
    Intersection *ilist;
 
2480
    SplineSet *ret;
 
2481
    SplineSet *order3 = NULL;
 
2482
    int is_o2 = false;
 
2483
    SplineSet *ss;
 
2484
 
 
2485
    for ( ss=base; ss!=NULL; ss=ss->next )
 
2486
        if ( ss->first->next!=NULL ) {
 
2487
            is_o2 = ss->first->next->order2;
 
2488
    break;
 
2489
        }
 
2490
    if ( is_o2 ) {
 
2491
        order3 = SplineSetsPSApprox(base);
 
2492
        SplinePointListsFree(base);
 
2493
        base = order3;
 
2494
    }
 
2495
 
 
2496
    if ( sc!=NULL )
 
2497
        glyphname = sc->name;
 
2498
 
 
2499
    base = SSRemoveTiny(base);
 
2500
    SSRemoveStupidControlPoints(base);
 
2501
    SplineSetsRemoveAnnoyingExtrema(base,.3);
 
2502
    SSOverlapClusterCpAngles(base,.01);
 
2503
    base = SSRemoveReversals(base);
 
2504
    ms = SSsToMContours(base,ot);
 
2505
    ilist = FindIntersections(ms,ot);
 
2506
    Validate(ms,ilist);
 
2507
    if ( ot==over_findinter || ot==over_fisel ) {
 
2508
        FixupIntersectedSplines(ms);
 
2509
        ret = base;
 
2510
    } else {
 
2511
        FindNeeded(ms,ot);
 
2512
        FindUnitVectors(ilist);
 
2513
        if ( ot==over_remove || ot == over_rmselected )
 
2514
            TestForBadDirections(ilist);
 
2515
        ret = JoinAllNeeded(ilist);
 
2516
        ret = MergeOpenAndFreeClosed(ret,base,ot);
 
2517
    }
 
2518
    FreeMonotonics(ms);
 
2519
    FreeIntersections(ilist);
 
2520
    glyphname = NULL;
 
2521
    if ( order3!=NULL ) {
 
2522
        ss = SplineSetsTTFApprox(ret);
 
2523
        SplinePointListsFree(ret);
 
2524
        ret = ss;
 
2525
    }
 
2526
return( ret );
 
2527
}