1 /*
  2     Copyright 2008-2024
  3         Matthias Ehmann,
  4         Michael Gerhaeuser,
  5         Carsten Miller,
  6         Bianca Valentin,
  7         Andreas Walter,
  8         Alfred Wassermann,
  9         Peter Wilfahrt
 10 
 11     This file is part of JSXGraph.
 12 
 13     JSXGraph is free software dual licensed under the GNU LGPL or MIT License.
 14 
 15     You can redistribute it and/or modify it under the terms of the
 16 
 17       * GNU Lesser General Public License as published by
 18         the Free Software Foundation, either version 3 of the License, or
 19         (at your option) any later version
 20       OR
 21       * MIT License: https://github.com/jsxgraph/jsxgraph/blob/master/LICENSE.MIT
 22 
 23     JSXGraph is distributed in the hope that it will be useful,
 24     but WITHOUT ANY WARRANTY; without even the implied warranty of
 25     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 26     GNU Lesser General Public License for more details.
 27 
 28     You should have received a copy of the GNU Lesser General Public License and
 29     the MIT License along with JSXGraph. If not, see <https://www.gnu.org/licenses/>
 30     and <https://opensource.org/licenses/MIT/>.
 31  */
 32 
 33 /*global JXG: true, define: true*/
 34 /*jslint nomen: true, plusplus: true*/
 35 
 36 /**
 37  * @fileoverview This file contains the Math.Geometry namespace for calculating algebraic/geometric
 38  * stuff like intersection points, angles, midpoint, and so on.
 39  */
 40 
 41 import JXG from "../jxg.js";
 42 import Const from "../base/constants.js";
 43 import Coords from "../base/coords.js";
 44 import Mat from "./math.js";
 45 import Numerics from "./numerics.js";
 46 import Type from "../utils/type.js";
 47 import Expect from "../utils/expect.js";
 48 
 49 /**
 50  * Math.Geometry namespace definition. This namespace holds geometrical algorithms,
 51  * especially intersection algorithms.
 52  * @name JXG.Math.Geometry
 53  * @namespace
 54  */
 55 Mat.Geometry = {};
 56 
 57 // the splitting is necessary due to the shortcut for the circumcircleMidpoint method to circumcenter.
 58 
 59 JXG.extend(
 60     Mat.Geometry,
 61     /** @lends JXG.Math.Geometry */ {
 62         /* ***************************************/
 63         /* *** GENERAL GEOMETRIC CALCULATIONS ****/
 64         /* ***************************************/
 65 
 66         /**
 67          * Calculates the angle defined by the points A, B, C.
 68          * @param {JXG.Point|Array} A A point  or [x,y] array.
 69          * @param {JXG.Point|Array} B Another point or [x,y] array.
 70          * @param {JXG.Point|Array} C A circle - no, of course the third point or [x,y] array.
 71          * @deprecated Use {@link JXG.Math.Geometry.rad} instead.
 72          * @see #rad
 73          * @see #trueAngle
 74          * @returns {Number} The angle in radian measure.
 75          */
 76         angle: function (A, B, C) {
 77             var u,
 78                 v,
 79                 s,
 80                 t,
 81                 a = [],
 82                 b = [],
 83                 c = [];
 84 
 85             JXG.deprecated("Geometry.angle()", "Geometry.rad()");
 86             if (A.coords) {
 87                 a[0] = A.coords.usrCoords[1];
 88                 a[1] = A.coords.usrCoords[2];
 89             } else {
 90                 a[0] = A[0];
 91                 a[1] = A[1];
 92             }
 93 
 94             if (B.coords) {
 95                 b[0] = B.coords.usrCoords[1];
 96                 b[1] = B.coords.usrCoords[2];
 97             } else {
 98                 b[0] = B[0];
 99                 b[1] = B[1];
100             }
101 
102             if (C.coords) {
103                 c[0] = C.coords.usrCoords[1];
104                 c[1] = C.coords.usrCoords[2];
105             } else {
106                 c[0] = C[0];
107                 c[1] = C[1];
108             }
109 
110             u = a[0] - b[0];
111             v = a[1] - b[1];
112             s = c[0] - b[0];
113             t = c[1] - b[1];
114 
115             return Math.atan2(u * t - v * s, u * s + v * t);
116         },
117 
118         /**
119          * Calculates the angle defined by the three points A, B, C if you're going from A to C around B counterclockwise.
120          * @param {JXG.Point|Array} A Point or [x,y] array
121          * @param {JXG.Point|Array} B Point or [x,y] array
122          * @param {JXG.Point|Array} C Point or [x,y] array
123          * @see #rad
124          * @returns {Number} The angle in degrees.
125          */
126         trueAngle: function (A, B, C) {
127             return this.rad(A, B, C) * 57.295779513082323; // *180.0/Math.PI;
128         },
129 
130         /**
131          * Calculates the internal angle defined by the three points A, B, C if you're going from A to C around B counterclockwise.
132          * @param {JXG.Point|Array} A Point or [x,y] array
133          * @param {JXG.Point|Array} B Point or [x,y] array
134          * @param {JXG.Point|Array} C Point or [x,y] array
135          * @see #trueAngle
136          * @returns {Number} Angle in radians.
137          */
138         rad: function (A, B, C) {
139             var ax, ay, bx, by, cx, cy, phi;
140 
141             if (A.coords) {
142                 ax = A.coords.usrCoords[1];
143                 ay = A.coords.usrCoords[2];
144             } else {
145                 ax = A[0];
146                 ay = A[1];
147             }
148 
149             if (B.coords) {
150                 bx = B.coords.usrCoords[1];
151                 by = B.coords.usrCoords[2];
152             } else {
153                 bx = B[0];
154                 by = B[1];
155             }
156 
157             if (C.coords) {
158                 cx = C.coords.usrCoords[1];
159                 cy = C.coords.usrCoords[2];
160             } else {
161                 cx = C[0];
162                 cy = C[1];
163             }
164 
165             phi = Math.atan2(cy - by, cx - bx) - Math.atan2(ay - by, ax - bx);
166 
167             if (phi < 0) {
168                 phi += 6.2831853071795862;
169             }
170 
171             return phi;
172         },
173 
174         /**
175          * Calculates a point on the bisection line between the three points A, B, C.
176          * As a result, the bisection line is defined by two points:
177          * Parameter B and the point with the coordinates calculated in this function.
178          * Does not work for ideal points.
179          * @param {JXG.Point} A Point
180          * @param {JXG.Point} B Point
181          * @param {JXG.Point} C Point
182          * @param [board=A.board] Reference to the board
183          * @returns {JXG.Coords} Coordinates of the second point defining the bisection.
184          */
185         angleBisector: function (A, B, C, board) {
186             var phiA,
187                 phiC,
188                 phi,
189                 Ac = A.coords.usrCoords,
190                 Bc = B.coords.usrCoords,
191                 Cc = C.coords.usrCoords,
192                 x,
193                 y;
194 
195             if (!Type.exists(board)) {
196                 board = A.board;
197             }
198 
199             // Parallel lines
200             if (Bc[0] === 0) {
201                 return new Coords(
202                     Const.COORDS_BY_USER,
203                     [1, (Ac[1] + Cc[1]) * 0.5, (Ac[2] + Cc[2]) * 0.5],
204                     board
205                 );
206             }
207 
208             // Non-parallel lines
209             x = Ac[1] - Bc[1];
210             y = Ac[2] - Bc[2];
211             phiA = Math.atan2(y, x);
212 
213             x = Cc[1] - Bc[1];
214             y = Cc[2] - Bc[2];
215             phiC = Math.atan2(y, x);
216 
217             phi = (phiA + phiC) * 0.5;
218 
219             if (phiA > phiC) {
220                 phi += Math.PI;
221             }
222 
223             x = Math.cos(phi) + Bc[1];
224             y = Math.sin(phi) + Bc[2];
225 
226             return new Coords(Const.COORDS_BY_USER, [1, x, y], board);
227         },
228 
229         // /**
230         //  * Calculates a point on the m-section line between the three points A, B, C.
231         //  * As a result, the m-section line is defined by two points:
232         //  * Parameter B and the point with the coordinates calculated in this function.
233         //  * The m-section generalizes the bisector to any real number.
234         //  * For example, the trisectors of an angle are simply the 1/3-sector and the 2/3-sector.
235         //  * Does not work for ideal points.
236         //  * @param {JXG.Point} A Point
237         //  * @param {JXG.Point} B Point
238         //  * @param {JXG.Point} C Point
239         //  * @param {Number} m Number
240         //  * @param [board=A.board] Reference to the board
241         //  * @returns {JXG.Coords} Coordinates of the second point defining the bisection.
242         //  */
243         // angleMsector: function (A, B, C, m, board) {
244         //     var phiA, phiC, phi,
245         //         Ac = A.coords.usrCoords,
246         //         Bc = B.coords.usrCoords,
247         //         Cc = C.coords.usrCoords,
248         //         x, y;
249 
250         //     if (!Type.exists(board)) {
251         //         board = A.board;
252         //     }
253 
254         //     // Parallel lines
255         //     if (Bc[0] === 0) {
256         //         return new Coords(Const.COORDS_BY_USER,
257         //             [1, (Ac[1] + Cc[1]) * m, (Ac[2] + Cc[2]) * m], board);
258         //     }
259 
260         //     // Non-parallel lines
261         //     x = Ac[1] - Bc[1];
262         //     y = Ac[2] - Bc[2];
263         //     phiA =  Math.atan2(y, x);
264 
265         //     x = Cc[1] - Bc[1];
266         //     y = Cc[2] - Bc[2];
267         //     phiC =  Math.atan2(y, x);
268 
269         //     phi = phiA + ((phiC - phiA) * m);
270 
271         //     if (phiA - phiC > Math.PI) {
272         //         phi += 2*m*Math.PI;
273         //     }
274 
275         //     x = Math.cos(phi) + Bc[1];
276         //     y = Math.sin(phi) + Bc[2];
277 
278         //     return new Coords(Const.COORDS_BY_USER, [1, x, y], board);
279         // },
280 
281         /**
282          * Reflects the point along the line.
283          * @param {JXG.Line} line Axis of reflection.
284          * @param {JXG.Point} point Point to reflect.
285          * @param [board=point.board] Reference to the board
286          * @returns {JXG.Coords} Coordinates of the reflected point.
287          */
288         reflection: function (line, point, board) {
289             // (v,w) defines the slope of the line
290             var x0,
291                 y0,
292                 x1,
293                 y1,
294                 v,
295                 w,
296                 mu,
297                 pc = point.coords.usrCoords,
298                 p1c = line.point1.coords.usrCoords,
299                 p2c = line.point2.coords.usrCoords;
300 
301             if (!Type.exists(board)) {
302                 board = point.board;
303             }
304 
305             v = p2c[1] - p1c[1];
306             w = p2c[2] - p1c[2];
307 
308             x0 = pc[1] - p1c[1];
309             y0 = pc[2] - p1c[2];
310 
311             mu = (v * y0 - w * x0) / (v * v + w * w);
312 
313             // point + mu*(-y,x) is the perpendicular foot
314             x1 = pc[1] + 2 * mu * w;
315             y1 = pc[2] - 2 * mu * v;
316 
317             return new Coords(Const.COORDS_BY_USER, [x1, y1], board);
318         },
319 
320         /**
321          * Computes the new position of a point which is rotated
322          * around a second point (called rotpoint) by the angle phi.
323          * @param {JXG.Point} rotpoint Center of the rotation
324          * @param {JXG.Point} point point to be rotated
325          * @param {Number} phi rotation angle in arc length
326          * @param {JXG.Board} [board=point.board] Reference to the board
327          * @returns {JXG.Coords} Coordinates of the new position.
328          */
329         rotation: function (rotpoint, point, phi, board) {
330             var x0,
331                 y0,
332                 c,
333                 s,
334                 x1,
335                 y1,
336                 pc = point.coords.usrCoords,
337                 rotpc = rotpoint.coords.usrCoords;
338 
339             if (!Type.exists(board)) {
340                 board = point.board;
341             }
342 
343             x0 = pc[1] - rotpc[1];
344             y0 = pc[2] - rotpc[2];
345 
346             c = Math.cos(phi);
347             s = Math.sin(phi);
348 
349             x1 = x0 * c - y0 * s + rotpc[1];
350             y1 = x0 * s + y0 * c + rotpc[2];
351 
352             return new Coords(Const.COORDS_BY_USER, [x1, y1], board);
353         },
354 
355         /**
356          * Calculates the coordinates of a point on the perpendicular to the given line through
357          * the given point.
358          * @param {JXG.Line} line A line.
359          * @param {JXG.Point} point Point which is projected to the line.
360          * @param {JXG.Board} [board=point.board] Reference to the board
361          * @returns {Array} Array of length two containing coordinates of a point on the perpendicular to the given line
362          *                  through the given point and boolean flag "change".
363          */
364         perpendicular: function (line, point, board) {
365             var x,
366                 y,
367                 change,
368                 c,
369                 z,
370                 A = line.point1.coords.usrCoords,
371                 B = line.point2.coords.usrCoords,
372                 C = point.coords.usrCoords;
373 
374             if (!Type.exists(board)) {
375                 board = point.board;
376             }
377 
378             // special case: point is the first point of the line
379             if (point === line.point1) {
380                 x = A[1] + B[2] - A[2];
381                 y = A[2] - B[1] + A[1];
382                 z = A[0] * B[0];
383 
384                 if (Math.abs(z) < Mat.eps) {
385                     x = B[2];
386                     y = -B[1];
387                 }
388                 c = [z, x, y];
389                 change = true;
390 
391                 // special case: point is the second point of the line
392             } else if (point === line.point2) {
393                 x = B[1] + A[2] - B[2];
394                 y = B[2] - A[1] + B[1];
395                 z = A[0] * B[0];
396 
397                 if (Math.abs(z) < Mat.eps) {
398                     x = A[2];
399                     y = -A[1];
400                 }
401                 c = [z, x, y];
402                 change = false;
403 
404                 // special case: point lies somewhere else on the line
405             } else if (Math.abs(Mat.innerProduct(C, line.stdform, 3)) < Mat.eps) {
406                 x = C[1] + B[2] - C[2];
407                 y = C[2] - B[1] + C[1];
408                 z = B[0];
409 
410                 if (Math.abs(z) < Mat.eps) {
411                     x = B[2];
412                     y = -B[1];
413                 }
414 
415                 change = true;
416                 if (
417                     Math.abs(z) > Mat.eps &&
418                     Math.abs(x - C[1]) < Mat.eps &&
419                     Math.abs(y - C[2]) < Mat.eps
420                 ) {
421                     x = C[1] + A[2] - C[2];
422                     y = C[2] - A[1] + C[1];
423                     change = false;
424                 }
425                 c = [z, x, y];
426 
427                 // general case: point does not lie on the line
428                 // -> calculate the foot of the dropped perpendicular
429             } else {
430                 c = [0, line.stdform[1], line.stdform[2]];
431                 c = Mat.crossProduct(c, C); // perpendicuar to line
432                 c = Mat.crossProduct(c, line.stdform); // intersection of line and perpendicular
433                 change = true;
434             }
435 
436             return [new Coords(Const.COORDS_BY_USER, c, board), change];
437         },
438 
439         /**
440          * @deprecated Please use {@link JXG.Math.Geometry.circumcenter} instead.
441          */
442         circumcenterMidpoint: function () {
443             JXG.deprecated("Geometry.circumcenterMidpoint()", "Geometry.circumcenter()");
444             this.circumcenter.apply(this, arguments);
445         },
446 
447         /**
448          * Calculates the center of the circumcircle of the three given points.
449          * @param {JXG.Point} point1 Point
450          * @param {JXG.Point} point2 Point
451          * @param {JXG.Point} point3 Point
452          * @param {JXG.Board} [board=point1.board] Reference to the board
453          * @returns {JXG.Coords} Coordinates of the center of the circumcircle of the given points.
454          */
455         circumcenter: function (point1, point2, point3, board) {
456             var u,
457                 v,
458                 m1,
459                 m2,
460                 A = point1.coords.usrCoords,
461                 B = point2.coords.usrCoords,
462                 C = point3.coords.usrCoords;
463 
464             if (!Type.exists(board)) {
465                 board = point1.board;
466             }
467 
468             u = [B[0] - A[0], -B[2] + A[2], B[1] - A[1]];
469             v = [(A[0] + B[0]) * 0.5, (A[1] + B[1]) * 0.5, (A[2] + B[2]) * 0.5];
470             m1 = Mat.crossProduct(u, v);
471 
472             u = [C[0] - B[0], -C[2] + B[2], C[1] - B[1]];
473             v = [(B[0] + C[0]) * 0.5, (B[1] + C[1]) * 0.5, (B[2] + C[2]) * 0.5];
474             m2 = Mat.crossProduct(u, v);
475 
476             return new Coords(Const.COORDS_BY_USER, Mat.crossProduct(m1, m2), board);
477         },
478 
479         /**
480          * Calculates the Euclidean distance for two given arrays of the same length.
481          * @param {Array} array1 Array of Number
482          * @param {Array} array2 Array of Number
483          * @param {Number} [n] Length of the arrays. Default is the minimum length of the given arrays.
484          * @returns {Number} Euclidean distance of the given vectors.
485          */
486         distance: function (array1, array2, n) {
487             var i,
488                 sum = 0;
489 
490             if (!n) {
491                 n = Math.min(array1.length, array2.length);
492             }
493 
494             for (i = 0; i < n; i++) {
495                 sum += (array1[i] - array2[i]) * (array1[i] - array2[i]);
496             }
497 
498             return Math.sqrt(sum);
499         },
500 
501         /**
502          * Calculates Euclidean distance for two given arrays of the same length.
503          * If one of the arrays contains a zero in the first coordinate, and the Euclidean distance
504          * is different from zero it is a point at infinity and we return Infinity.
505          * @param {Array} array1 Array containing elements of type number.
506          * @param {Array} array2 Array containing elements of type number.
507          * @param {Number} [n] Length of the arrays. Default is the minimum length of the given arrays.
508          * @returns {Number} Euclidean (affine) distance of the given vectors.
509          */
510         affineDistance: function (array1, array2, n) {
511             var d;
512 
513             d = this.distance(array1, array2, n);
514 
515             if (
516                 d > Mat.eps &&
517                 (Math.abs(array1[0]) < Mat.eps || Math.abs(array2[0]) < Mat.eps)
518             ) {
519                 return Infinity;
520             }
521 
522             return d;
523         },
524 
525         /**
526          * Affine ratio of three collinear points a, b, c: (c - a) / (b - a).
527          * If r > 1 or r < 0 then c is outside of the segment ab.
528          *
529          * @param {Array|JXG.Coords} a
530          * @param {Array|JXG.Coords} b
531          * @param {Array|JXG.Coords} c
532          * @returns {Number} affine ratio (c - a) / (b - a)
533          */
534         affineRatio: function (a, b, c) {
535             var r = 0.0,
536                 dx;
537 
538             if (Type.exists(a.usrCoords)) {
539                 a = a.usrCoords;
540             }
541             if (Type.exists(b.usrCoords)) {
542                 b = b.usrCoords;
543             }
544             if (Type.exists(c.usrCoords)) {
545                 c = c.usrCoords;
546             }
547 
548             dx = b[1] - a[1];
549 
550             if (Math.abs(dx) > Mat.eps) {
551                 r = (c[1] - a[1]) / dx;
552             } else {
553                 r = (c[2] - a[2]) / (b[2] - a[2]);
554             }
555             return r;
556         },
557 
558         /**
559          * Sort vertices counter clockwise starting with the first point.
560          *
561          * @param {Array} p An array containing {@link JXG.Point}, {@link JXG.Coords}, and/or arrays.
562          *
563          * @returns {Array}
564          */
565         sortVertices: function (p) {
566             var ll,
567                 ps = Expect.each(p, Expect.coordsArray),
568                 N = ps.length,
569                 lastPoint = null;
570 
571             // If the last point equals the first point, we take the last point out of the array.
572             // It may be that the several points at the end of the array are equal to the first point.
573             // The polygonal chain is been closed by JSXGraph, but this may also have been done by the user.
574             // Therefore, we use a while lopp to pop the last points.
575             while (
576                 ps[0][0] === ps[N - 1][0] &&
577                 ps[0][1] === ps[N - 1][1] &&
578                 ps[0][2] === ps[N - 1][2]
579             ) {
580                 lastPoint = ps.pop();
581                 N--;
582             }
583             // Find the point with the lowest y value
584             // for (i = 1; i < N; i++) {
585             //     if ((ps[i][2] < ps[0][2]) ||
586             //         // if the current and the lowest point have the same y value, pick the one with
587             //         // the lowest x value.
588             //         (Math.abs(ps[i][2] - ps[0][2]) < Mat.eps && ps[i][1] < ps[0][1])) {
589             //         console.log(i, 0);
590             //         ps = Type.swap(ps, i, 0);
591             //     }
592             // }
593 
594             ll = ps[0];
595             // Sort ps in increasing order of the angle between a point and the first point ll.
596             // If a point is equal to the first point ll, the angle is defined to be -Infinity.
597             // Otherwise, atan2 would return zero, which is a value which also attained by points
598             // on the same horizontal line.
599             ps.sort(function (a, b) {
600                 var rad1 =
601                     a[2] === ll[2] && a[1] === ll[1]
602                         ? -Infinity
603                         : Math.atan2(a[2] - ll[2], a[1] - ll[1]),
604                     rad2 =
605                         b[2] === ll[2] && b[1] === ll[1]
606                             ? -Infinity
607                             : Math.atan2(b[2] - ll[2], b[1] - ll[1]);
608 
609                 return rad1 - rad2;
610             });
611 
612             // If the last point has been taken out of the array, we put it in again.
613             if (lastPoint !== null) {
614                 ps.push(lastPoint);
615             }
616 
617             return ps;
618         },
619 
620         /**
621          * Signed triangle area of the three points given.
622          *
623          * @param {JXG.Point|JXG.Coords|Array} p1
624          * @param {JXG.Point|JXG.Coords|Array} p2
625          * @param {JXG.Point|JXG.Coords|Array} p3
626          *
627          * @returns {Number}
628          */
629         signedTriangle: function (p1, p2, p3) {
630             var A = Expect.coordsArray(p1),
631                 B = Expect.coordsArray(p2),
632                 C = Expect.coordsArray(p3);
633 
634             return 0.5 * ((B[1] - A[1]) * (C[2] - A[2]) - (B[2] - A[2]) * (C[1] - A[1]));
635         },
636 
637         /**
638          * Determine the signed area of a non-selfintersecting polygon.
639          * Surveyor's Formula
640          *
641          * @param {Array} p An array containing {@link JXG.Point}, {@link JXG.Coords}, and/or arrays.
642          * @param {Boolean} [sort=true]
643          *
644          * @returns {Number}
645          */
646         signedPolygon: function (p, sort) {
647             var i,
648                 N,
649                 A = 0,
650                 ps = Expect.each(p, Expect.coordsArray);
651 
652             if (sort === undefined) {
653                 sort = true;
654             }
655 
656             if (!sort) {
657                 ps = this.sortVertices(ps);
658             } else {
659                 // Make sure the polygon is closed. If it is already closed this won't change the sum because the last
660                 // summand will be 0.
661                 ps.unshift(ps[ps.length - 1]);
662             }
663 
664             N = ps.length;
665 
666             for (i = 1; i < N; i++) {
667                 A += ps[i - 1][1] * ps[i][2] - ps[i][1] * ps[i - 1][2];
668             }
669 
670             return 0.5 * A;
671         },
672 
673         /**
674          * Calculate the complex hull of a point cloud.
675          *
676          * @param {Array} points An array containing {@link JXG.Point}, {@link JXG.Coords}, and/or arrays.
677          *
678          * @returns {Array}
679          */
680         GrahamScan: function (points) {
681             var i,
682                 M = 1,
683                 ps = Expect.each(points, Expect.coordsArray),
684                 N = ps.length;
685 
686             ps = this.sortVertices(ps);
687             N = ps.length;
688 
689             for (i = 2; i < N; i++) {
690                 while (this.signedTriangle(ps[M - 1], ps[M], ps[i]) <= 0) {
691                     if (M > 1) {
692                         M -= 1;
693                     } else if (i === N - 1) {
694                         break;
695                     }
696                     i += 1;
697                 }
698 
699                 M += 1;
700                 ps = Type.swap(ps, M, i);
701             }
702 
703             return ps.slice(0, M);
704         },
705 
706         /**
707          * A line can be a segment, a straight, or a ray. So it is not always delimited by point1 and point2
708          * calcStraight determines the visual start point and end point of the line. A segment is only drawn
709          * from start to end point, a straight line is drawn until it meets the boards boundaries.
710          * @param {JXG.Line} el Reference to a line object, that needs calculation of start and end point.
711          * @param {JXG.Coords} point1 Coordinates of the point where line drawing begins. This value is calculated and
712          * set by this method.
713          * @param {JXG.Coords} point2 Coordinates of the point where line drawing ends. This value is calculated and set
714          * by this method.
715          * @param {Number} margin Optional margin, to avoid the display of the small sides of lines.
716          * @returns null
717          * @see Line
718          * @see JXG.Line
719          */
720         calcStraight: function (el, point1, point2, margin) {
721             var takePoint1,
722                 takePoint2,
723                 intersection,
724                 intersect1,
725                 intersect2,
726                 straightFirst,
727                 straightLast,
728                 c, p1, p2;
729 
730             if (!Type.exists(margin)) {
731                 // Enlarge the drawable region slightly. This hides the small sides
732                 // of thick lines in most cases.
733                 margin = 10;
734             }
735 
736             straightFirst = Type.evaluate(el.visProp.straightfirst);
737             straightLast = Type.evaluate(el.visProp.straightlast);
738 
739             // If one of the point is an ideal point in homogeneous coordinates
740             // drawing of line segments or rays are not possible.
741             if (Math.abs(point1.scrCoords[0]) < Mat.eps) {
742                 straightFirst = true;
743             }
744             if (Math.abs(point2.scrCoords[0]) < Mat.eps) {
745                 straightLast = true;
746             }
747 
748             // Do nothing in case of line segments (inside or outside of the board)
749             if (!straightFirst && !straightLast) {
750                 return;
751             }
752 
753             // Compute the stdform of the line in screen coordinates.
754             c = [];
755             c[0] =
756                 el.stdform[0] -
757                 (el.stdform[1] * el.board.origin.scrCoords[1]) / el.board.unitX +
758                 (el.stdform[2] * el.board.origin.scrCoords[2]) / el.board.unitY;
759             c[1] = el.stdform[1] / el.board.unitX;
760             c[2] = -el.stdform[2] / el.board.unitY;
761 
762             // If p1=p2
763             if (isNaN(c[0] + c[1] + c[2])) {
764                 return;
765             }
766 
767             takePoint1 = false;
768             takePoint2 = false;
769 
770             // Line starts at point1 and point1 is inside the board
771             takePoint1 =
772                 !straightFirst &&
773                 Math.abs(point1.usrCoords[0]) >= Mat.eps &&
774                 point1.scrCoords[1] >= 0.0 &&
775                 point1.scrCoords[1] <= el.board.canvasWidth &&
776                 point1.scrCoords[2] >= 0.0 &&
777                 point1.scrCoords[2] <= el.board.canvasHeight;
778 
779             // Line ends at point2 and point2 is inside the board
780             takePoint2 =
781                 !straightLast &&
782                 Math.abs(point2.usrCoords[0]) >= Mat.eps &&
783                 point2.scrCoords[1] >= 0.0 &&
784                 point2.scrCoords[1] <= el.board.canvasWidth &&
785                 point2.scrCoords[2] >= 0.0 &&
786                 point2.scrCoords[2] <= el.board.canvasHeight;
787 
788             // Intersect the line with the four borders of the board.
789             intersection = this.meetLineBoard(c, el.board, margin);
790             intersect1 = intersection[0];
791             intersect2 = intersection[1];
792 
793             /**
794              * At this point we have four points:
795              * point1 and point2 are the first and the second defining point on the line,
796              * intersect1, intersect2 are the intersections of the line with border around the board.
797              */
798 
799             /*
800              * Here we handle rays where both defining points are outside of the board.
801              */
802             // If both points are outside and the complete ray is outside we do nothing
803             if (!takePoint1 && !takePoint2) {
804                 // Ray starting at point 1
805                 if (
806                     !straightFirst &&
807                     straightLast &&
808                     !this.isSameDirection(point1, point2, intersect1) &&
809                     !this.isSameDirection(point1, point2, intersect2)
810                 ) {
811                     return;
812                 }
813 
814                 // Ray starting at point 2
815                 if (
816                     straightFirst &&
817                     !straightLast &&
818                     !this.isSameDirection(point2, point1, intersect1) &&
819                     !this.isSameDirection(point2, point1, intersect2)
820                 ) {
821                     return;
822                 }
823             }
824 
825             /*
826              * If at least one of the defining points is outside of the board
827              * we take intersect1 or intersect2 as one of the end points
828              * The order is also important for arrows of axes
829              */
830             if (!takePoint1) {
831                 if (!takePoint2) {
832                     // Two border intersection points are used
833                     if (this.isSameDir(point1, point2, intersect1, intersect2)) {
834                         p1 = intersect1;
835                         p2 = intersect2;
836                     } else {
837                         p2 = intersect1;
838                         p1 = intersect2;
839                     }
840                 } else {
841                     // One border intersection points is used
842                     if (this.isSameDir(point1, point2, intersect1, intersect2)) {
843                         p1 = intersect1;
844                     } else {
845                         p1 = intersect2;
846                     }
847                 }
848             } else {
849                 if (!takePoint2) {
850                     // One border intersection points is used
851                     if (this.isSameDir(point1, point2, intersect1, intersect2)) {
852                         p2 = intersect2;
853                     } else {
854                         p2 = intersect1;
855                     }
856                 }
857             }
858 
859             if (p1) {
860                 //point1.setCoordinates(Const.COORDS_BY_USER, p1.usrCoords.slice(1));
861                 point1.setCoordinates(Const.COORDS_BY_USER, p1.usrCoords);
862             }
863 
864             if (p2) {
865                 //point2.setCoordinates(Const.COORDS_BY_USER, p2.usrCoords.slice(1));
866                 point2.setCoordinates(Const.COORDS_BY_USER, p2.usrCoords);
867             }
868         },
869 
870         /**
871          * A line can be a segment, a straight, or a ray. so it is not always delimited by point1 and point2.
872          *
873          * This method adjusts the line's delimiting points taking into account its nature, the viewport defined
874          * by the board.
875          *
876          * A segment is delimited by start and end point, a straight line or ray is delimited until it meets the
877          * boards boundaries. However, if the line has infinite ticks, it will be delimited by the projection of
878          * the boards vertices onto itself.
879          *
880          * @param {JXG.Line} el Reference to a line object, that needs calculation of start and end point.
881          * @param {JXG.Coords} point1 Coordinates of the point where line drawing begins. This value is calculated and
882          * set by this method.
883          * @param {JXG.Coords} point2 Coordinates of the point where line drawing ends. This value is calculated and set
884          * by this method.
885          * @see Line
886          * @see JXG.Line
887          */
888         calcLineDelimitingPoints: function (el, point1, point2) {
889             var distP1P2,
890                 boundingBox,
891                 lineSlope,
892                 intersect1,
893                 intersect2,
894                 straightFirst,
895                 straightLast,
896                 c,
897                 p1,
898                 p2,
899                 takePoint1 = false,
900                 takePoint2 = false;
901 
902             straightFirst = Type.evaluate(el.visProp.straightfirst);
903             straightLast = Type.evaluate(el.visProp.straightlast);
904 
905             // If one of the point is an ideal point in homogeneous coordinates
906             // drawing of line segments or rays are not possible.
907             if (Math.abs(point1.scrCoords[0]) < Mat.eps) {
908                 straightFirst = true;
909             }
910             if (Math.abs(point2.scrCoords[0]) < Mat.eps) {
911                 straightLast = true;
912             }
913 
914             // Compute the stdform of the line in screen coordinates.
915             c = [];
916             c[0] =
917                 el.stdform[0] -
918                 (el.stdform[1] * el.board.origin.scrCoords[1]) / el.board.unitX +
919                 (el.stdform[2] * el.board.origin.scrCoords[2]) / el.board.unitY;
920             c[1] = el.stdform[1] / el.board.unitX;
921             c[2] = -el.stdform[2] / el.board.unitY;
922 
923             // p1=p2
924             if (isNaN(c[0] + c[1] + c[2])) {
925                 return;
926             }
927 
928             takePoint1 = !straightFirst;
929             takePoint2 = !straightLast;
930             // Intersect the board vertices on the line to establish the available visual space for the infinite ticks
931             // Based on the slope of the line we can optimise and only project the two outer vertices
932 
933             // boundingBox = [x1, y1, x2, y2] upper left, lower right vertices
934             boundingBox = el.board.getBoundingBox();
935             lineSlope = el.getSlope();
936             if (lineSlope >= 0) {
937                 // project vertices (x2,y1) (x1, y2)
938                 intersect1 = this.projectPointToLine(
939                     { coords: { usrCoords: [1, boundingBox[2], boundingBox[1]] } },
940                     el,
941                     el.board
942                 );
943                 intersect2 = this.projectPointToLine(
944                     { coords: { usrCoords: [1, boundingBox[0], boundingBox[3]] } },
945                     el,
946                     el.board
947                 );
948             } else {
949                 // project vertices (x1, y1) (x2, y2)
950                 intersect1 = this.projectPointToLine(
951                     { coords: { usrCoords: [1, boundingBox[0], boundingBox[1]] } },
952                     el,
953                     el.board
954                 );
955                 intersect2 = this.projectPointToLine(
956                     { coords: { usrCoords: [1, boundingBox[2], boundingBox[3]] } },
957                     el,
958                     el.board
959                 );
960             }
961 
962             /**
963              * we have four points:
964              * point1 and point2 are the first and the second defining point on the line,
965              * intersect1, intersect2 are the intersections of the line with border around the board.
966              */
967 
968             /*
969              * Here we handle rays/segments where both defining points are outside of the board.
970              */
971             if (!takePoint1 && !takePoint2) {
972                 // Segment, if segment does not cross the board, do nothing
973                 if (!straightFirst && !straightLast) {
974                     distP1P2 = point1.distance(Const.COORDS_BY_USER, point2);
975                     // if  intersect1 not between point1 and point2
976                     if (
977                         Math.abs(
978                             point1.distance(Const.COORDS_BY_USER, intersect1) +
979                             intersect1.distance(Const.COORDS_BY_USER, point2) -
980                             distP1P2
981                         ) > Mat.eps
982                     ) {
983                         return;
984                     }
985                     // if insersect2 not between point1 and point2
986                     if (
987                         Math.abs(
988                             point1.distance(Const.COORDS_BY_USER, intersect2) +
989                             intersect2.distance(Const.COORDS_BY_USER, point2) -
990                             distP1P2
991                         ) > Mat.eps
992                     ) {
993                         return;
994                     }
995                 }
996 
997                 // If both points are outside and the complete ray is outside we do nothing
998                 // Ray starting at point 1
999                 if (
1000                     !straightFirst &&
1001                     straightLast &&
1002                     !this.isSameDirection(point1, point2, intersect1) &&
1003                     !this.isSameDirection(point1, point2, intersect2)
1004                 ) {
1005                     return;
1006                 }
1007 
1008                 // Ray starting at point 2
1009                 if (
1010                     straightFirst &&
1011                     !straightLast &&
1012                     !this.isSameDirection(point2, point1, intersect1) &&
1013                     !this.isSameDirection(point2, point1, intersect2)
1014                 ) {
1015                     return;
1016                 }
1017             }
1018 
1019             /*
1020              * If at least one of the defining points is outside of the board
1021              * we take intersect1 or intersect2 as one of the end points
1022              * The order is also important for arrows of axes
1023              */
1024             if (!takePoint1) {
1025                 if (!takePoint2) {
1026                     // Two border intersection points are used
1027                     if (this.isSameDir(point1, point2, intersect1, intersect2)) {
1028                         p1 = intersect1;
1029                         p2 = intersect2;
1030                     } else {
1031                         p2 = intersect1;
1032                         p1 = intersect2;
1033                     }
1034                 } else {
1035                     // One border intersection points is used
1036                     if (this.isSameDir(point1, point2, intersect1, intersect2)) {
1037                         p1 = intersect1;
1038                     } else {
1039                         p1 = intersect2;
1040                     }
1041                 }
1042             } else {
1043                 if (!takePoint2) {
1044                     // One border intersection points is used
1045                     if (this.isSameDir(point1, point2, intersect1, intersect2)) {
1046                         p2 = intersect2;
1047                     } else {
1048                         p2 = intersect1;
1049                     }
1050                 }
1051             }
1052 
1053             if (p1) {
1054                 //point1.setCoordinates(Const.COORDS_BY_USER, p1.usrCoords.slice(1));
1055                 point1.setCoordinates(Const.COORDS_BY_USER, p1.usrCoords);
1056             }
1057 
1058             if (p2) {
1059                 //point2.setCoordinates(Const.COORDS_BY_USER, p2.usrCoords.slice(1));
1060                 point2.setCoordinates(Const.COORDS_BY_USER, p2.usrCoords);
1061             }
1062         },
1063 
1064         /**
1065          * Calculates the visProp.position corresponding to a given angle.
1066          * @param {number} angle angle in radians. Must be in range (-2pi,2pi).
1067          */
1068         calcLabelQuadrant: function (angle) {
1069             var q;
1070             if (angle < 0) {
1071                 angle += 2 * Math.PI;
1072             }
1073             q = Math.floor((angle + Math.PI / 8) / (Math.PI / 4)) % 8;
1074             return ["rt", "urt", "top", "ulft", "lft", "llft", "lrt"][q];
1075         },
1076 
1077         /**
1078          * The vectors <tt>p2-p1</tt> and <tt>i2-i1</tt> are supposed to be collinear. If their cosine is positive
1079          * they point into the same direction otherwise they point in opposite direction.
1080          * @param {JXG.Coords} p1
1081          * @param {JXG.Coords} p2
1082          * @param {JXG.Coords} i1
1083          * @param {JXG.Coords} i2
1084          * @returns {Boolean} True, if <tt>p2-p1</tt> and <tt>i2-i1</tt> point into the same direction
1085          */
1086         isSameDir: function (p1, p2, i1, i2) {
1087             var dpx = p2.usrCoords[1] - p1.usrCoords[1],
1088                 dpy = p2.usrCoords[2] - p1.usrCoords[2],
1089                 dix = i2.usrCoords[1] - i1.usrCoords[1],
1090                 diy = i2.usrCoords[2] - i1.usrCoords[2];
1091 
1092             if (Math.abs(p2.usrCoords[0]) < Mat.eps) {
1093                 dpx = p2.usrCoords[1];
1094                 dpy = p2.usrCoords[2];
1095             }
1096 
1097             if (Math.abs(p1.usrCoords[0]) < Mat.eps) {
1098                 dpx = -p1.usrCoords[1];
1099                 dpy = -p1.usrCoords[2];
1100             }
1101 
1102             return dpx * dix + dpy * diy >= 0;
1103         },
1104 
1105         /**
1106          * If you're looking from point "start" towards point "s" and you can see the point "p", return true.
1107          * Otherwise return false.
1108          * @param {JXG.Coords} start The point you're standing on.
1109          * @param {JXG.Coords} p The point in which direction you're looking.
1110          * @param {JXG.Coords} s The point that should be visible.
1111          * @returns {Boolean} True, if from start the point p is in the same direction as s is, that means s-start = k*(p-start) with k>=0.
1112          */
1113         isSameDirection: function (start, p, s) {
1114             var dx,
1115                 dy,
1116                 sx,
1117                 sy,
1118                 r = false;
1119 
1120             dx = p.usrCoords[1] - start.usrCoords[1];
1121             dy = p.usrCoords[2] - start.usrCoords[2];
1122 
1123             sx = s.usrCoords[1] - start.usrCoords[1];
1124             sy = s.usrCoords[2] - start.usrCoords[2];
1125 
1126             if (Math.abs(dx) < Mat.eps) {
1127                 dx = 0;
1128             }
1129 
1130             if (Math.abs(dy) < Mat.eps) {
1131                 dy = 0;
1132             }
1133 
1134             if (Math.abs(sx) < Mat.eps) {
1135                 sx = 0;
1136             }
1137 
1138             if (Math.abs(sy) < Mat.eps) {
1139                 sy = 0;
1140             }
1141 
1142             if (dx >= 0 && sx >= 0) {
1143                 r = (dy >= 0 && sy >= 0) || (dy <= 0 && sy <= 0);
1144             } else if (dx <= 0 && sx <= 0) {
1145                 r = (dy >= 0 && sy >= 0) || (dy <= 0 && sy <= 0);
1146             }
1147 
1148             return r;
1149         },
1150 
1151         /**
1152          * Determinant of three points in the Euclidean plane.
1153          * Zero, if the points are collinear. Used to determine of a point q is left or
1154          * right to a segment defined by points p1 and p2.
1155          * <p>
1156          * Non-homogeneous version.
1157          *
1158          * @param  {Array|JXG.Point} p1 First point or its coordinates of the segment. Point object or array of length 3. First (homogeneous) coordinate is equal to 1.
1159          * @param  {Array|JXG.Point} p2 Second point or its coordinates of the segment. Point object or array of length 3. First (homogeneous) coordinate is equal to 1.
1160          * @param  {Array|JXG.Point} q Point or its coordinates. Point object or array of length 3. First (homogeneous) coordinate is equal to 1.
1161          * @return {Number} Signed area of the triangle formed by these three points.
1162          *
1163          * @see #windingNumber
1164          */
1165         det3p: function (p1, p2, q) {
1166             var pp1, pp2, qq;
1167 
1168             if (Type.isPoint(p1)) {
1169                 pp1 = p1.Coords(true);
1170             } else {
1171                 pp1 = p1;
1172             }
1173             if (Type.isPoint(p2)) {
1174                 pp2 = p2.Coords(true);
1175             } else {
1176                 pp2 = p2;
1177             }
1178             if (Type.isPoint(q)) {
1179                 qq = q.Coords(true);
1180             } else {
1181                 qq = q;
1182             }
1183 
1184             return (pp1[1] - qq[1]) * (pp2[2] - qq[2]) - (pp2[1] - qq[1]) * (pp1[2] - qq[2]);
1185         },
1186 
1187         /**
1188          * Winding number of a point in respect to a polygon path.
1189          *
1190          * The point is regarded outside if the winding number is zero,
1191          * inside otherwise. The algorithm tries to find degenerate cases, i.e.
1192          * if the point is on the path. This is regarded as "outside".
1193          * If the point is a vertex of the path, it is regarded as "inside".
1194          *
1195          * Implementation of algorithm 7 from "The point in polygon problem for
1196          * arbitrary polygons" by Kai Hormann and Alexander Agathos, Computational Geometry,
1197          * Volume 20, Issue 3, November 2001, Pages 131-144.
1198          *
1199          * @param  {Array} usrCoords Homogenous coordinates of the point
1200          * @param  {Array} path      Array of points / coords determining a path, i.e. the vertices of the polygon / path. The array elements
1201          * do not have to be full points, but have to have a subobject "coords" or should be of type JXG.Coords.
1202          * @param  {Boolean} [doNotClosePath=false] If true the last point of the path is not connected to the first point.
1203          * This is necessary if the path consists of two or more closed subpaths, e.g. if the figure has a hole.
1204          *
1205          * @return {Number}          Winding number of the point. The point is
1206          *                           regarded outside if the winding number is zero,
1207          *                           inside otherwise.
1208          */
1209         windingNumber: function (usrCoords, path, doNotClosePath) {
1210             var wn = 0,
1211                 le = path.length,
1212                 x = usrCoords[1],
1213                 y = usrCoords[2],
1214                 p0,
1215                 p1,
1216                 p2,
1217                 d,
1218                 sign,
1219                 i,
1220                 off = 0;
1221 
1222             if (le === 0) {
1223                 return 0;
1224             }
1225 
1226             doNotClosePath = doNotClosePath || false;
1227             if (doNotClosePath) {
1228                 off = 1;
1229             }
1230 
1231             // Infinite points are declared outside
1232             if (isNaN(x) || isNaN(y)) {
1233                 return 1;
1234             }
1235 
1236             if (Type.exists(path[0].coords)) {
1237                 p0 = path[0].coords;
1238                 p1 = path[le - 1].coords;
1239             } else {
1240                 p0 = path[0];
1241                 p1 = path[le - 1];
1242             }
1243             // Handle the case if the point is the first vertex of the path, i.e. inside.
1244             if (p0.usrCoords[1] === x && p0.usrCoords[2] === y) {
1245                 return 1;
1246             }
1247 
1248             for (i = 0; i < le - off; i++) {
1249                 // Consider the edge from p1 = path[i] to p2 = path[i+1]isClosedPath
1250                 if (Type.exists(path[i].coords)) {
1251                     p1 = path[i].coords.usrCoords;
1252                     p2 = path[(i + 1) % le].coords.usrCoords;
1253                 } else {
1254                     p1 = path[i].usrCoords;
1255                     p2 = path[(i + 1) % le].usrCoords;
1256                 }
1257 
1258                 // If one of the two points p1, p2 is undefined or infinite,
1259                 // move on.
1260                 if (
1261                     p1[0] === 0 ||
1262                     p2[0] === 0 ||
1263                     isNaN(p1[1]) ||
1264                     isNaN(p2[1]) ||
1265                     isNaN(p1[2]) ||
1266                     isNaN(p2[2])
1267                 ) {
1268                     continue;
1269                 }
1270 
1271                 if (p2[2] === y) {
1272                     if (p2[1] === x) {
1273                         return 1;
1274                     }
1275                     if (p1[2] === y && p2[1] > x === p1[1] < x) {
1276                         return 0;
1277                     }
1278                 }
1279 
1280                 if (p1[2] < y !== p2[2] < y) {
1281                     // Crossing
1282                     sign = 2 * (p2[2] > p1[2] ? 1 : 0) - 1;
1283                     if (p1[1] >= x) {
1284                         if (p2[1] > x) {
1285                             wn += sign;
1286                         } else {
1287                             d = this.det3p(p1, p2, usrCoords);
1288                             if (d === 0) {
1289                                 // Point is on line, i.e. outside
1290                                 return 0;
1291                             }
1292                             if (d > 0 + Mat.eps === p2[2] > p1[2]) {
1293                                 // Right crossing
1294                                 wn += sign;
1295                             }
1296                         }
1297                     } else {
1298                         if (p2[1] > x) {
1299                             d = this.det3p(p1, p2, usrCoords);
1300                             if (d > 0 + Mat.eps === p2[2] > p1[2]) {
1301                                 // Right crossing
1302                                 wn += sign;
1303                             }
1304                         }
1305                     }
1306                 }
1307             }
1308 
1309             return wn;
1310         },
1311 
1312         /**
1313          * Decides if a point (x,y) is inside of a path / polygon.
1314          * Does not work correct if the path has hole. In this case, windingNumber is the preferred method.
1315          * Implements W. Randolf Franklin's pnpoly method.
1316          *
1317          * See <a href="https://wrf.ecse.rpi.edu/Research/Short_Notes/pnpoly.html">https://wrf.ecse.rpi.edu/Research/Short_Notes/pnpoly.html</a>.
1318          *
1319          * @param {Number} x_in x-coordinate (screen or user coordinates)
1320          * @param {Number} y_in y-coordinate (screen or user coordinates)
1321          * @param  {Array} path  Array of points / coords determining a path, i.e. the vertices of the polygon / path. The array elements
1322          * do not have to be full points, but have to have a subobject "coords" or should be of type JXG.Coords.
1323          * @param {Number} [coord_type=JXG.COORDS_BY_SCREEN] Type of coordinates used here.
1324          *   Possible values are <b>JXG.COORDS_BY_USER</b> and <b>JXG.COORDS_BY_SCREEN</b>.
1325          *   Default value is JXG.COORDS_BY_SCREEN.
1326          *
1327          * @returns {Boolean} if (x_in, y_in) is inside of the polygon.
1328          * @see JXG.Polygon.hasPoint
1329          * @see JXG.Polygon.pnpoly
1330          * @see #windingNumber
1331          *
1332          * @example
1333          * var pol = board.create('polygon', [[-1,2], [2,2], [-1,4]]);
1334          * var p = board.create('point', [4, 3]);
1335          * var txt = board.create('text', [-1, 0.5, function() {
1336          *   return 'Point A is inside of the polygon = ' +
1337          *     JXG.Math.Geometry.pnpoly(p.X(), p.Y(), JXG.COORDS_BY_USER, pol.vertices);
1338          * }]);
1339          *
1340          * </pre><div id="JXG4656ed42-f965-4e35-bb66-c334a4529683" class="jxgbox" style="width: 300px; height: 300px;"></div>
1341          * <script type="text/javascript">
1342          *     (function() {
1343          *         var board = JXG.JSXGraph.initBoard('JXG4656ed42-f965-4e35-bb66-c334a4529683',
1344          *             {boundingbox: [-2, 5, 5,-2], axis: true, showcopyright: false, shownavigation: false});
1345          *     var pol = board.create('polygon', [[-1,2], [2,2], [-1,4]]);
1346          *     var p = board.create('point', [4, 3]);
1347          *     var txt = board.create('text', [-1, 0.5, function() {
1348          *     		return 'Point A is inside of the polygon = ' + JXG.Math.Geometry.pnpoly(p.X(), p.Y(), JXG.COORDS_BY_USER, pol.vertices);
1349          *     }]);
1350          *
1351          *     })();
1352          *
1353          * </script><pre>
1354          *
1355          */
1356         pnpoly: function (x_in, y_in, path, coord_type) {
1357             var i, j, vi, vj, len,
1358                 x, y, crds,
1359                 v = path,
1360                 isIn = false;
1361 
1362             if (coord_type === Const.COORDS_BY_USER) {
1363                 crds = new Coords(Const.COORDS_BY_USER, [x_in, y_in], this.board);
1364                 x = crds.scrCoords[1];
1365                 y = crds.scrCoords[2];
1366             } else {
1367                 x = x_in;
1368                 y = y_in;
1369             }
1370 
1371             len = path.length;
1372             for (i = 0, j = len - 2; i < len - 1; j = i++) {
1373                 vi = Type.exists(v[i].coords) ? v[i].coords : v[i];
1374                 vj = Type.exists(v[j].coords) ? v[j].coords : v[j];
1375 
1376                 if (
1377                     vi.scrCoords[2] > y !== vj.scrCoords[2] > y &&
1378                     x <
1379                     ((vj.scrCoords[1] - vi.scrCoords[1]) * (y - vi.scrCoords[2])) /
1380                     (vj.scrCoords[2] - vi.scrCoords[2]) +
1381                     vi.scrCoords[1]
1382                 ) {
1383                     isIn = !isIn;
1384                 }
1385             }
1386 
1387             return isIn;
1388         },
1389 
1390         /****************************************/
1391         /****          INTERSECTIONS         ****/
1392         /****************************************/
1393 
1394         /**
1395          * Generate the function which computes the coordinates of the intersection point.
1396          * Primarily used in {@link JXG.Point#createIntersectionPoint}.
1397          * @param {JXG.Board} board object
1398          * @param {JXG.Line,JXG.Circle_JXG.Line,JXG.Circle_Number|Function} el1,el2,i The result will be a intersection point on el1 and el2.
1399          * i determines the intersection point if two points are available: <ul>
1400          *   <li>i==0: use the positive square root,</li>
1401          *   <li>i==1: use the negative square root.</li></ul>
1402          * See further {@link JXG.Point#createIntersectionPoint}.
1403          * @param {Boolean} alwaysintersect. Flag that determines if segments and arc can have an outer intersection point
1404          * on their defining line or circle.
1405          * @returns {Function} Function returning a {@link JXG.Coords} object that determines
1406          * the intersection point.
1407          */
1408         intersectionFunction: function (board, el1, el2, i, j, alwaysintersect) {
1409             var func,
1410                 that = this,
1411                 el1_isArcType = false,
1412                 el2_isArcType = false;
1413 
1414             el1_isArcType =
1415                 el1.elementClass === Const.OBJECT_CLASS_CURVE &&
1416                     (el1.type === Const.OBJECT_TYPE_ARC || el1.type === Const.OBJECT_TYPE_SECTOR)
1417                     ? true
1418                     : false;
1419             el2_isArcType =
1420                 el2.elementClass === Const.OBJECT_CLASS_CURVE &&
1421                     (el2.type === Const.OBJECT_TYPE_ARC || el2.type === Const.OBJECT_TYPE_SECTOR)
1422                     ? true
1423                     : false;
1424 
1425             if (
1426                 (el1.elementClass === Const.OBJECT_CLASS_CURVE ||
1427                     el2.elementClass === Const.OBJECT_CLASS_CURVE) &&
1428                 (el1.elementClass === Const.OBJECT_CLASS_CURVE ||
1429                     el1.elementClass === Const.OBJECT_CLASS_CIRCLE) &&
1430                 (el2.elementClass === Const.OBJECT_CLASS_CURVE ||
1431                     el2.elementClass === Const.OBJECT_CLASS_CIRCLE) /*&&
1432                 !(el1_isArcType && el2_isArcType)*/
1433             ) {
1434                 // curve - curve
1435                 // with the exception that both elements are arc types
1436                 /** @ignore */
1437                 func = function () {
1438                     return that.meetCurveCurve(el1, el2, i, j, el1.board);
1439                 };
1440             } else if (
1441                 (el1.elementClass === Const.OBJECT_CLASS_CURVE &&
1442                     !el1_isArcType &&
1443                     el2.elementClass === Const.OBJECT_CLASS_LINE) ||
1444                 (el2.elementClass === Const.OBJECT_CLASS_CURVE &&
1445                     !el2_isArcType &&
1446                     el1.elementClass === Const.OBJECT_CLASS_LINE)
1447             ) {
1448                 // curve - line (this includes intersections between conic sections and lines)
1449                 // with the exception that the curve is of arc type
1450                 /** @ignore */
1451                 func = function () {
1452                     return that.meetCurveLine(el1, el2, i, el1.board, Type.evaluate(alwaysintersect));
1453                 };
1454             } else if (
1455                 el1.type === Const.OBJECT_TYPE_POLYGON ||
1456                 el2.type === Const.OBJECT_TYPE_POLYGON
1457             ) {
1458                 // polygon - other
1459                 // Uses the Greiner-Hormann clipping algorithm
1460                 // Not implemented: polygon - point
1461 
1462                 if (el1.elementClass === Const.OBJECT_CLASS_LINE) {
1463                     // line - path
1464                     /** @ignore */
1465                     func = function () {
1466                         var first1 = Type.evaluate(el1.visProp.straightfirst),
1467                             last1 = Type.evaluate(el1.visProp.straightlast),
1468                             first2 = Type.evaluate(el2.visProp.straightfirst),
1469                             last2 = Type.evaluate(el2.visProp.straightlast),
1470                             a_not;
1471 
1472                         a_not = (!Type.evaluate(alwaysintersect) && (!first1 || !last1 || !first2 || !last2));
1473                         return that.meetPolygonLine(el2, el1, i, el1.board, a_not);
1474                     };
1475                 } else if (el2.elementClass === Const.OBJECT_CLASS_LINE) {
1476                     // path - line
1477                     func = function () {
1478                         var first1 = Type.evaluate(el1.visProp.straightfirst),
1479                             last1 = Type.evaluate(el1.visProp.straightlast),
1480                             first2 = Type.evaluate(el2.visProp.straightfirst),
1481                             last2 = Type.evaluate(el2.visProp.straightlast),
1482                             a_not;
1483 
1484                         a_not = (!Type.evaluate(alwaysintersect) && (!first1 || !last1 || !first2 || !last2));
1485                         return that.meetPolygonLine(el1, el2, i, el1.board, a_not);
1486                     };
1487                 } else {
1488                     // path - path
1489                     /** @ignore */
1490                     func = function () {
1491                         return that.meetPathPath(el1, el2, i, el1.board);
1492                     };
1493                 }
1494             } else if (
1495                 el1.elementClass === Const.OBJECT_CLASS_LINE &&
1496                 el2.elementClass === Const.OBJECT_CLASS_LINE
1497             ) {
1498                 // line - line, lines may also be segments.
1499                 /** @ignore */
1500                 func = function () {
1501                     var res,
1502                         c,
1503                         first1 = Type.evaluate(el1.visProp.straightfirst),
1504                         last1 = Type.evaluate(el1.visProp.straightlast),
1505                         first2 = Type.evaluate(el2.visProp.straightfirst),
1506                         last2 = Type.evaluate(el2.visProp.straightlast);
1507 
1508                     /**
1509                      * If one of the lines is a segment or ray and
1510                      * the intersection point should disappear if outside
1511                      * of the segment or ray we call
1512                      * meetSegmentSegment
1513                      */
1514                     if (
1515                         !Type.evaluate(alwaysintersect) &&
1516                         (!first1 || !last1 || !first2 || !last2)
1517                     ) {
1518                         res = that.meetSegmentSegment(
1519                             el1.point1.coords.usrCoords,
1520                             el1.point2.coords.usrCoords,
1521                             el2.point1.coords.usrCoords,
1522                             el2.point2.coords.usrCoords
1523                         );
1524 
1525                         if (
1526                             (!first1 && res[1] < 0) ||
1527                             (!last1 && res[1] > 1) ||
1528                             (!first2 && res[2] < 0) ||
1529                             (!last2 && res[2] > 1)
1530                         ) {
1531                             // Non-existent
1532                             c = [0, NaN, NaN];
1533                         } else {
1534                             c = res[0];
1535                         }
1536 
1537                         return new Coords(Const.COORDS_BY_USER, c, el1.board);
1538                     }
1539 
1540                     return that.meet(el1.stdform, el2.stdform, i, el1.board);
1541                 };
1542             } else {
1543                 // All other combinations of circles and lines,
1544                 // Arc types are treated as circles.
1545                 /** @ignore */
1546                 func = function () {
1547                     var res = that.meet(el1.stdform, el2.stdform, i, el1.board),
1548                         has = true,
1549                         first,
1550                         last,
1551                         r;
1552 
1553                     if (Type.evaluate(alwaysintersect)) {
1554                         return res;
1555                     }
1556                     if (el1.elementClass === Const.OBJECT_CLASS_LINE) {
1557                         first = Type.evaluate(el1.visProp.straightfirst);
1558                         last = Type.evaluate(el1.visProp.straightlast);
1559                         if (!first || !last) {
1560                             r = that.affineRatio(el1.point1.coords, el1.point2.coords, res);
1561                             if ((!last && r > 1 + Mat.eps) || (!first && r < 0 - Mat.eps)) {
1562                                 return new Coords(JXG.COORDS_BY_USER, [0, NaN, NaN], el1.board);
1563                             }
1564                         }
1565                     }
1566                     if (el2.elementClass === Const.OBJECT_CLASS_LINE) {
1567                         first = Type.evaluate(el2.visProp.straightfirst);
1568                         last = Type.evaluate(el2.visProp.straightlast);
1569                         if (!first || !last) {
1570                             r = that.affineRatio(el2.point1.coords, el2.point2.coords, res);
1571                             if ((!last && r > 1 + Mat.eps) || (!first && r < 0 - Mat.eps)) {
1572                                 return new Coords(JXG.COORDS_BY_USER, [0, NaN, NaN], el1.board);
1573                             }
1574                         }
1575                     }
1576                     if (el1_isArcType) {
1577                         has = that.coordsOnArc(el1, res);
1578                         if (has && el2_isArcType) {
1579                             has = that.coordsOnArc(el2, res);
1580                         }
1581                         if (!has) {
1582                             return new Coords(JXG.COORDS_BY_USER, [0, NaN, NaN], el1.board);
1583                         }
1584                     }
1585                     return res;
1586                 };
1587             }
1588 
1589             return func;
1590         },
1591 
1592         otherIntersectionFunction: function (input, others, alwaysintersect, precision) {
1593             var func, board,
1594                 el1, el2,
1595                 that = this;
1596 
1597             el1 = input[0];
1598             el2 = input[1];
1599             board = el1.board;
1600             func = function () {
1601                 var i, k, c, d,
1602                     isClose,
1603                     le = others.length,
1604                     eps = Type.evaluate(precision);
1605 
1606                 for (i = le; i >= 0; i--) {
1607                     if (el1.elementClass === Const.OBJECT_CLASS_CIRCLE &&
1608                         [Const.OBJECT_CLASS_CIRCLE, Const.OBJECT_CLASS_LINE].indexOf(el2.elementClass) >= 0) {
1609                         // circle, circle|line
1610                         c = that.meet(el1.stdform, el2.stdform, i, board);
1611                     } else if (el1.elementClass === Const.OBJECT_CLASS_CURVE &&
1612                         [Const.OBJECT_CLASS_CURVE, Const.OBJECT_CLASS_CIRCLE].indexOf(el2.elementClass) >= 0) {
1613                         // curve, circle|curve
1614                         c = that.meetCurveCurve(el1, el2, i, 0, board, 'segment');
1615                     } else if (el1.elementClass === Const.OBJECT_CLASS_CURVE && el2.elementClass === Const.OBJECT_CLASS_LINE) {
1616                         // curve, line
1617                         if (Type.exists(el1.dataX)) {
1618                             c = JXG.Math.Geometry.meetCurveLine(el1, el2, i, el1.board, Type.evaluate(alwaysintersect));
1619                         } else {
1620                             c = JXG.Math.Geometry.meetCurveLineContinuous(el1, el2, i, el1.board);
1621                         }
1622                     }
1623 
1624                     // If the intersection is close to one of the points in other
1625                     // we have to search for another intersection point.
1626                     isClose = false;
1627                     for (k = 0; !isClose && k < le; k++) {
1628                         d = c.distance(JXG.COORDS_BY_USER, others[k].coords);
1629                         if (d < eps) {
1630                             isClose = true;
1631                         }
1632                     }
1633                     if (!isClose) {
1634                         // We are done, the intersection is away from any other
1635                         // intersection point.
1636                         return c;
1637                     }
1638                 }
1639                 // Otherwise we return the last intersection point
1640                 return c;
1641             };
1642             return func;
1643         },
1644 
1645         /**
1646          * Generate the function which computes the data of the intersection.
1647          */
1648         intersectionFunction3D: function (view, el1, el2, i) {
1649             var func,
1650                 that = this;
1651 
1652             if (el1.type === Const.OBJECT_TYPE_PLANE3D) {
1653                 if (el2.type === Const.OBJECT_TYPE_PLANE3D) {
1654                     func = () => view.intersectionPlanePlane(el1, el2)[i];
1655                 } else if (el2.type === Const.OBJECT_TYPE_SPHERE3D) {
1656                     func = that.meetPlaneSphere(el1, el2);
1657                 }
1658             } else if (el1.type === Const.OBJECT_TYPE_SPHERE3D) {
1659                 if (el2.type === Const.OBJECT_TYPE_PLANE3D) {
1660                     func = that.meetPlaneSphere(el2, el1);
1661                 } else if (el2.type === Const.OBJECT_TYPE_SPHERE3D) {
1662                     func = that.meetSphereSphere(el1, el2);
1663                 }
1664             }
1665 
1666             return func;
1667         },
1668 
1669         /**
1670          * Returns true if the coordinates are on the arc element,
1671          * false otherwise. Usually, coords is an intersection
1672          * on the circle line. Now it is decided if coords are on the
1673          * circle restricted to the arc line.
1674          * @param  {Arc} arc arc or sector element
1675          * @param  {JXG.Coords} coords Coords object of an intersection
1676          * @returns {Boolean}
1677          * @private
1678          */
1679         coordsOnArc: function (arc, coords) {
1680             var angle = this.rad(arc.radiuspoint, arc.center, coords.usrCoords.slice(1)),
1681                 alpha = 0.0,
1682                 beta = this.rad(arc.radiuspoint, arc.center, arc.anglepoint),
1683                 ev_s = Type.evaluate(arc.visProp.selection);
1684 
1685             if ((ev_s === "minor" && beta > Math.PI) || (ev_s === "major" && beta < Math.PI)) {
1686                 alpha = beta;
1687                 beta = 2 * Math.PI;
1688             }
1689             if (angle < alpha || angle > beta) {
1690                 return false;
1691             }
1692             return true;
1693         },
1694 
1695         /**
1696          * Computes the intersection of a pair of lines, circles or both.
1697          * It uses the internal data array stdform of these elements.
1698          * @param {Array} el1 stdform of the first element (line or circle)
1699          * @param {Array} el2 stdform of the second element (line or circle)
1700          * @param {Number|Function} i Index of the intersection point that should be returned.
1701          * @param board Reference to the board.
1702          * @returns {JXG.Coords} Coordinates of one of the possible two or more intersection points.
1703          * Which point will be returned is determined by i.
1704          */
1705         meet: function (el1, el2, i, board) {
1706             var result,
1707                 eps = Mat.eps;
1708 
1709             if (Math.abs(el1[3]) < eps && Math.abs(el2[3]) < eps) {
1710                 // line line
1711                 result = this.meetLineLine(el1, el2, i, board);
1712             } else if (Math.abs(el1[3]) >= eps && Math.abs(el2[3]) < eps) {
1713                 // circle line
1714                 result = this.meetLineCircle(el2, el1, i, board);
1715             } else if (Math.abs(el1[3]) < eps && Math.abs(el2[3]) >= eps) {
1716                 // line circle
1717                 result = this.meetLineCircle(el1, el2, i, board);
1718             } else {
1719                 // circle circle
1720                 result = this.meetCircleCircle(el1, el2, i, board);
1721             }
1722 
1723             return result;
1724         },
1725 
1726         /**
1727          * Intersection of the line with the board
1728          * @param  {Array}     line   stdform of the line in screen coordinates
1729          * @param  {JXG.Board} board  reference to a board.
1730          * @param  {Number}    margin optional margin, to avoid the display of the small sides of lines.
1731          * @returns {Array}            [intersection coords 1, intersection coords 2]
1732          */
1733         meetLineBoard: function (line, board, margin) {
1734             // Intersect the line with the four borders of the board.
1735             var s = [],
1736                 intersect1,
1737                 intersect2,
1738                 i, j;
1739 
1740             if (!Type.exists(margin)) {
1741                 margin = 0;
1742             }
1743 
1744             // top
1745             s[0] = Mat.crossProduct(line, [margin, 0, 1]);
1746             // left
1747             s[1] = Mat.crossProduct(line, [margin, 1, 0]);
1748             // bottom
1749             s[2] = Mat.crossProduct(line, [-margin - board.canvasHeight, 0, 1]);
1750             // right
1751             s[3] = Mat.crossProduct(line, [-margin - board.canvasWidth, 1, 0]);
1752 
1753             // Normalize the intersections
1754             for (i = 0; i < 4; i++) {
1755                 if (Math.abs(s[i][0]) > Mat.eps) {
1756                     for (j = 2; j > 0; j--) {
1757                         s[i][j] /= s[i][0];
1758                     }
1759                     s[i][0] = 1.0;
1760                 }
1761             }
1762 
1763             // line is parallel to "left", take "top" and "bottom"
1764             if (Math.abs(s[1][0]) < Mat.eps) {
1765                 intersect1 = s[0]; // top
1766                 intersect2 = s[2]; // bottom
1767                 // line is parallel to "top", take "left" and "right"
1768             } else if (Math.abs(s[0][0]) < Mat.eps) {
1769                 intersect1 = s[1]; // left
1770                 intersect2 = s[3]; // right
1771                 // left intersection out of board (above)
1772             } else if (s[1][2] < 0) {
1773                 intersect1 = s[0]; // top
1774 
1775                 // right intersection out of board (below)
1776                 if (s[3][2] > board.canvasHeight) {
1777                     intersect2 = s[2]; // bottom
1778                 } else {
1779                     intersect2 = s[3]; // right
1780                 }
1781                 // left intersection out of board (below)
1782             } else if (s[1][2] > board.canvasHeight) {
1783                 intersect1 = s[2]; // bottom
1784 
1785                 // right intersection out of board (above)
1786                 if (s[3][2] < 0) {
1787                     intersect2 = s[0]; // top
1788                 } else {
1789                     intersect2 = s[3]; // right
1790                 }
1791             } else {
1792                 intersect1 = s[1]; // left
1793 
1794                 // right intersection out of board (above)
1795                 if (s[3][2] < 0) {
1796                     intersect2 = s[0]; // top
1797                     // right intersection out of board (below)
1798                 } else if (s[3][2] > board.canvasHeight) {
1799                     intersect2 = s[2]; // bottom
1800                 } else {
1801                     intersect2 = s[3]; // right
1802                 }
1803             }
1804 
1805             return [
1806                 new Coords(Const.COORDS_BY_SCREEN, intersect1.slice(1), board),
1807                 new Coords(Const.COORDS_BY_SCREEN, intersect2.slice(1), board)
1808             ];
1809         },
1810 
1811         /**
1812          * Intersection of two lines.
1813          * @param {Array} l1 stdform of the first line
1814          * @param {Array} l2 stdform of the second line
1815          * @param {number} i unused
1816          * @param {JXG.Board} board Reference to the board.
1817          * @returns {JXG.Coords} Coordinates of the intersection point.
1818          */
1819         meetLineLine: function (l1, l2, i, board) {
1820             var s = isNaN(l1[5] + l2[5]) ? [0, 0, 0] : Mat.crossProduct(l1, l2);
1821 
1822             // Make intersection of parallel lines more robust:
1823             if (Math.abs(s[0]) < 1.0e-14) {
1824                 s[0] = 0;
1825             }
1826             return new Coords(Const.COORDS_BY_USER, s, board);
1827         },
1828 
1829         /**
1830          * Intersection of line and circle.
1831          * @param {Array} lin stdform of the line
1832          * @param {Array} circ stdform of the circle
1833          * @param {number|function} i number of the returned intersection point.
1834          *   i==0: use the positive square root,
1835          *   i==1: use the negative square root.
1836          * @param {JXG.Board} board Reference to a board.
1837          * @returns {JXG.Coords} Coordinates of the intersection point
1838          */
1839         meetLineCircle: function (lin, circ, i, board) {
1840             var a, b, c, d, n, A, B, C, k, t;
1841 
1842             // Radius is zero, return center of circle
1843             if (circ[4] < Mat.eps) {
1844                 if (Math.abs(Mat.innerProduct([1, circ[6], circ[7]], lin, 3)) < Mat.eps) {
1845                     return new Coords(Const.COORDS_BY_USER, circ.slice(6, 8), board);
1846                 }
1847 
1848                 return new Coords(Const.COORDS_BY_USER, [NaN, NaN], board);
1849             }
1850             c = circ[0];
1851             b = circ.slice(1, 3);
1852             a = circ[3];
1853             d = lin[0];
1854             n = lin.slice(1, 3);
1855 
1856             // Line is assumed to be normalized. Therefore, nn==1 and we can skip some operations:
1857             /*
1858              var nn = n[0]*n[0]+n[1]*n[1];
1859              A = a*nn;
1860              B = (b[0]*n[1]-b[1]*n[0])*nn;
1861              C = a*d*d - (b[0]*n[0]+b[1]*n[1])*d + c*nn;
1862              */
1863             A = a;
1864             B = b[0] * n[1] - b[1] * n[0];
1865             C = a * d * d - (b[0] * n[0] + b[1] * n[1]) * d + c;
1866 
1867             k = B * B - 4 * A * C;
1868             if (k > -Mat.eps * Mat.eps) {
1869                 k = Math.sqrt(Math.abs(k));
1870                 t = [(-B + k) / (2 * A), (-B - k) / (2 * A)];
1871 
1872                 return Type.evaluate(i) === 0
1873                     ? new Coords(
1874                         Const.COORDS_BY_USER,
1875                         [-t[0] * -n[1] - d * n[0], -t[0] * n[0] - d * n[1]],
1876                         board
1877                     )
1878                     : new Coords(
1879                         Const.COORDS_BY_USER,
1880                         [-t[1] * -n[1] - d * n[0], -t[1] * n[0] - d * n[1]],
1881                         board
1882                     );
1883             }
1884 
1885             return new Coords(Const.COORDS_BY_USER, [0, 0, 0], board);
1886         },
1887 
1888         /**
1889          * Intersection of two circles.
1890          * @param {Array} circ1 stdform of the first circle
1891          * @param {Array} circ2 stdform of the second circle
1892          * @param {number|function} i number of the returned intersection point.
1893          *   i==0: use the positive square root,
1894          *   i==1: use the negative square root.
1895          * @param {JXG.Board} board Reference to the board.
1896          * @returns {JXG.Coords} Coordinates of the intersection point
1897          */
1898         meetCircleCircle: function (circ1, circ2, i, board) {
1899             var radicalAxis;
1900 
1901             // Radius is zero, return center of circle, if on other circle
1902             if (circ1[4] < Mat.eps) {
1903                 if (
1904                     Math.abs(this.distance(circ1.slice(6, 2), circ2.slice(6, 8)) - circ2[4]) <
1905                     Mat.eps
1906                 ) {
1907                     return new Coords(Const.COORDS_BY_USER, circ1.slice(6, 8), board);
1908                 }
1909 
1910                 return new Coords(Const.COORDS_BY_USER, [0, 0, 0], board);
1911             }
1912 
1913             // Radius is zero, return center of circle, if on other circle
1914             if (circ2[4] < Mat.eps) {
1915                 if (
1916                     Math.abs(this.distance(circ2.slice(6, 2), circ1.slice(6, 8)) - circ1[4]) <
1917                     Mat.eps
1918                 ) {
1919                     return new Coords(Const.COORDS_BY_USER, circ2.slice(6, 8), board);
1920                 }
1921 
1922                 return new Coords(Const.COORDS_BY_USER, [0, 0, 0], board);
1923             }
1924 
1925             radicalAxis = [
1926                 circ2[3] * circ1[0] - circ1[3] * circ2[0],
1927                 circ2[3] * circ1[1] - circ1[3] * circ2[1],
1928                 circ2[3] * circ1[2] - circ1[3] * circ2[2],
1929                 0,
1930                 1,
1931                 Infinity,
1932                 Infinity,
1933                 Infinity
1934             ];
1935             radicalAxis = Mat.normalize(radicalAxis);
1936 
1937             return this.meetLineCircle(radicalAxis, circ1, i, board);
1938         },
1939 
1940         /**
1941          * Compute an intersection of the curves c1 and c2.
1942          * We want to find values t1, t2 such that
1943          * c1(t1) = c2(t2), i.e. (c1_x(t1)-c2_x(t2),c1_y(t1)-c2_y(t2)) = (0,0).
1944          *
1945          * Methods: segment-wise intersections (default) or generalized Newton method.
1946          * @param {JXG.Curve} c1 Curve, Line or Circle
1947          * @param {JXG.Curve} c2 Curve, Line or Circle
1948          * @param {Number|Function} nr the nr-th intersection point will be returned.
1949          * @param {Number} t2ini not longer used.
1950          * @param {JXG.Board} [board=c1.board] Reference to a board object.
1951          * @param {String} [method='segment'] Intersection method, possible values are 'newton' and 'segment'.
1952          * @returns {JXG.Coords} intersection point
1953          */
1954         meetCurveCurve: function (c1, c2, nr, t2ini, board, method) {
1955             var co;
1956 
1957             if (Type.exists(method) && method === "newton") {
1958                 co = Numerics.generalizedNewton(c1, c2, Type.evaluate(nr), t2ini);
1959             } else {
1960                 if (c1.bezierDegree === 3 || c2.bezierDegree === 3) {
1961                     co = this.meetBezierCurveRedBlueSegments(c1, c2, nr);
1962                 } else {
1963                     co = this.meetCurveRedBlueSegments(c1, c2, nr);
1964                 }
1965             }
1966 
1967             return new Coords(Const.COORDS_BY_USER, co, board);
1968         },
1969 
1970         /**
1971          * Intersection of curve with line,
1972          * Order of input does not matter for el1 and el2.
1973          * From version 0.99.7 on this method calls
1974          * {@link JXG.Math.Geometry.meetCurveLineDiscrete}.
1975          * If higher precision is needed, {@link JXG.Math.Geometry.meetCurveLineContinuous}
1976          * has to be used.
1977          *
1978          * @param {JXG.Curve|JXG.Line} el1 Curve or Line
1979          * @param {JXG.Curve|JXG.Line} el2 Curve or Line
1980          * @param {Number|Function} nr the nr-th intersection point will be returned.
1981          * @param {JXG.Board} [board=el1.board] Reference to a board object.
1982          * @param {Boolean} alwaysIntersect If false just the segment between the two defining points are tested for intersection
1983          * @returns {JXG.Coords} Intersection point. In case no intersection point is detected,
1984          * the ideal point [0,1,0] is returned.
1985          */
1986         meetCurveLine: function (el1, el2, nr, board, alwaysIntersect) {
1987             var v = [0, NaN, NaN],
1988                 cu,
1989                 li;
1990 
1991             if (!Type.exists(board)) {
1992                 board = el1.board;
1993             }
1994 
1995             if (el1.elementClass === Const.OBJECT_CLASS_CURVE) {
1996                 cu = el1;
1997                 li = el2;
1998             } else {
1999                 cu = el2;
2000                 li = el1;
2001             }
2002 
2003             v = this.meetCurveLineDiscrete(cu, li, nr, board, !alwaysIntersect);
2004 
2005             return v;
2006         },
2007 
2008         /**
2009          * Intersection of line and curve, continuous case.
2010          * Finds the nr-the intersection point
2011          * Uses {@link JXG.Math.Geometry.meetCurveLineDiscrete} as a first approximation.
2012          * A more exact solution is then found with {@link JXG.Math.Numerics.root}.
2013          *
2014          * @param {JXG.Curve} cu Curve
2015          * @param {JXG.Line} li Line
2016          * @param {NumberFunction} nr Will return the nr-th intersection point.
2017          * @param {JXG.Board} board
2018          * @param {Boolean} testSegment Test if intersection has to be inside of the segment or somewhere on the
2019          * line defined by the segment
2020          * @returns {JXG.Coords} Coords object containing the intersection.
2021          */
2022         meetCurveLineContinuous: function (cu, li, nr, board, testSegment) {
2023             var func0, func1,
2024                 t, v, x, y, z,
2025                 eps = Mat.eps,
2026                 epsLow = Mat.eps,
2027                 steps,
2028                 delta,
2029                 tnew, tmin, fmin,
2030                 i, ft;
2031 
2032             v = this.meetCurveLineDiscrete(cu, li, nr, board, testSegment);
2033             x = v.usrCoords[1];
2034             y = v.usrCoords[2];
2035 
2036             func0 = function (t) {
2037                 var c1, c2;
2038 
2039                 if (t > cu.maxX() || t < cu.minX()) {
2040                     return Infinity;
2041                 }
2042                 c1 = cu.X(t) - x;
2043                 c2 = cu.Y(t) - y;
2044                 return c1 * c1 + c2 * c2;
2045                 // return c1 * (cu.X(t + h) - cu.X(t - h)) + c2 * (cu.Y(t + h) - cu.Y(t - h)) / h;
2046             };
2047 
2048             func1 = function (t) {
2049                 var v = li.stdform[0] + li.stdform[1] * cu.X(t) + li.stdform[2] * cu.Y(t);
2050                 return v * v;
2051             };
2052 
2053             // Find t
2054             steps = 50;
2055             delta = (cu.maxX() - cu.minX()) / steps;
2056             tnew = cu.minX();
2057             fmin = 0.0001; //eps;
2058             tmin = NaN;
2059             for (i = 0; i < steps; i++) {
2060                 t = Numerics.root(func0, [
2061                     Math.max(tnew, cu.minX()),
2062                     Math.min(tnew + delta, cu.maxX())
2063                 ]);
2064                 ft = Math.abs(func0(t));
2065                 if (ft <= fmin) {
2066                     fmin = ft;
2067                     tmin = t;
2068                     if (fmin < eps) {
2069                         break;
2070                     }
2071                 }
2072 
2073                 tnew += delta;
2074             }
2075             t = tmin;
2076             // Compute "exact" t
2077             t = Numerics.root(func1, [
2078                 Math.max(t - delta, cu.minX()),
2079                 Math.min(t + delta, cu.maxX())
2080             ]);
2081 
2082             ft = func1(t);
2083             // Is the point on the line?
2084             if (isNaN(ft) || Math.abs(ft) > epsLow) {
2085                 z = 0.0; //NaN;
2086             } else {
2087                 z = 1.0;
2088             }
2089 
2090             return new Coords(Const.COORDS_BY_USER, [z, cu.X(t), cu.Y(t)], board);
2091         },
2092 
2093         /**
2094          * Intersection of line and curve, discrete case.
2095          * Segments are treated as lines.
2096          * Finding the nr-th intersection point should work for all nr.
2097          * @param {JXG.Curve} cu
2098          * @param {JXG.Line} li
2099          * @param {Number|Function} nr
2100          * @param {JXG.Board} board
2101          * @param {Boolean} testSegment Test if intersection has to be inside of the segment or somewhere on the
2102          * line defined by the segment
2103          *
2104          * @returns {JXG.Coords} Intersection point. In case no intersection point is detected,
2105          * the ideal point [0,1,0] is returned.
2106          */
2107         meetCurveLineDiscrete: function (cu, li, nr, board, testSegment) {
2108             var i, j,
2109                 n = Type.evaluate(nr),
2110                 p1, p2,
2111                 p, q,
2112                 lip1 = li.point1.coords.usrCoords,
2113                 lip2 = li.point2.coords.usrCoords,
2114                 d, res,
2115                 cnt = 0,
2116                 len = cu.numberPoints,
2117                 ev_sf = Type.evaluate(li.visProp.straightfirst),
2118                 ev_sl = Type.evaluate(li.visProp.straightlast);
2119 
2120             // In case, no intersection will be found we will take this
2121             q = new Coords(Const.COORDS_BY_USER, [0, NaN, NaN], board);
2122 
2123             if (lip1[0] === 0.0) {
2124                 lip1 = [1, lip2[1] + li.stdform[2], lip2[2] - li.stdform[1]];
2125             } else if (lip2[0] === 0.0) {
2126                 lip2 = [1, lip1[1] + li.stdform[2], lip1[2] - li.stdform[1]];
2127             }
2128 
2129             p2 = cu.points[0].usrCoords;
2130             for (i = 1; i < len; i += cu.bezierDegree) {
2131                 p1 = p2.slice(0);
2132                 p2 = cu.points[i].usrCoords;
2133                 d = this.distance(p1, p2);
2134 
2135                 // The defining points are not identical
2136                 if (d > Mat.eps) {
2137                     if (cu.bezierDegree === 3) {
2138                         res = this.meetBeziersegmentBeziersegment(
2139                             [
2140                                 cu.points[i - 1].usrCoords.slice(1),
2141                                 cu.points[i].usrCoords.slice(1),
2142                                 cu.points[i + 1].usrCoords.slice(1),
2143                                 cu.points[i + 2].usrCoords.slice(1)
2144                             ],
2145                             [lip1.slice(1), lip2.slice(1)],
2146                             testSegment
2147                         );
2148                     } else {
2149                         res = [this.meetSegmentSegment(p1, p2, lip1, lip2)];
2150                     }
2151 
2152                     for (j = 0; j < res.length; j++) {
2153                         p = res[j];
2154                         if (0 <= p[1] && p[1] <= 1) {
2155                             if (cnt === n) {
2156                                 /**
2157                                  * If the intersection point is not part of the segment,
2158                                  * this intersection point is set to non-existent.
2159                                  * This prevents jumping behavior of the intersection points.
2160                                  * But it may be discussed if it is the desired behavior.
2161                                  */
2162                                 if (
2163                                     testSegment &&
2164                                     ((!ev_sf && p[2] < 0) || (!ev_sl && p[2] > 1))
2165                                 ) {
2166                                     return q; // break;
2167                                 }
2168 
2169                                 q = new Coords(Const.COORDS_BY_USER, p[0], board);
2170                                 return q; // break;
2171                             }
2172                             cnt += 1;
2173                         }
2174                     }
2175                 }
2176             }
2177 
2178             return q;
2179         },
2180 
2181         /**
2182          * Find the n-th intersection point of two curves named red (first parameter) and blue (second parameter).
2183          * We go through each segment of the red curve and search if there is an intersection with a segment of the blue curve.
2184          * This double loop, i.e. the outer loop runs along the red curve and the inner loop runs along the blue curve, defines
2185          * the n-th intersection point. The segments are either line segments or Bezier curves of degree 3. This depends on
2186          * the property bezierDegree of the curves.
2187          * <p>
2188          * This method works also for transformed curves, since only the already
2189          * transformed points are used.
2190          *
2191          * @param {JXG.Curve} red
2192          * @param {JXG.Curve} blue
2193          * @param {Number|Function} nr
2194          */
2195         meetCurveRedBlueSegments: function (red, blue, nr) {
2196             var i,
2197                 j,
2198                 n = Type.evaluate(nr),
2199                 red1,
2200                 red2,
2201                 blue1,
2202                 blue2,
2203                 m,
2204                 minX,
2205                 maxX,
2206                 iFound = 0,
2207                 lenBlue = blue.numberPoints, //points.length,
2208                 lenRed = red.numberPoints; //points.length;
2209 
2210             if (lenBlue <= 1 || lenRed <= 1) {
2211                 return [0, NaN, NaN];
2212             }
2213 
2214             for (i = 1; i < lenRed; i++) {
2215                 red1 = red.points[i - 1].usrCoords;
2216                 red2 = red.points[i].usrCoords;
2217                 minX = Math.min(red1[1], red2[1]);
2218                 maxX = Math.max(red1[1], red2[1]);
2219 
2220                 blue2 = blue.points[0].usrCoords;
2221                 for (j = 1; j < lenBlue; j++) {
2222                     blue1 = blue2;
2223                     blue2 = blue.points[j].usrCoords;
2224 
2225                     if (
2226                         Math.min(blue1[1], blue2[1]) < maxX &&
2227                         Math.max(blue1[1], blue2[1]) > minX
2228                     ) {
2229                         m = this.meetSegmentSegment(red1, red2, blue1, blue2);
2230                         if (
2231                             m[1] >= 0.0 &&
2232                             m[2] >= 0.0 &&
2233                             // The two segments meet in the interior or at the start points
2234                             ((m[1] < 1.0 && m[2] < 1.0) ||
2235                                 // One of the curve is intersected in the very last point
2236                                 (i === lenRed - 1 && m[1] === 1.0) ||
2237                                 (j === lenBlue - 1 && m[2] === 1.0))
2238                         ) {
2239                             if (iFound === n) {
2240                                 return m[0];
2241                             }
2242 
2243                             iFound++;
2244                         }
2245                     }
2246                 }
2247             }
2248 
2249             return [0, NaN, NaN];
2250         },
2251 
2252         /**
2253          * (Virtual) Intersection of two segments.
2254          * @param {Array} p1 First point of segment 1 using normalized homogeneous coordinates [1,x,y]
2255          * @param {Array} p2 Second point or direction of segment 1 using normalized homogeneous coordinates [1,x,y] or point at infinity [0,x,y], respectively
2256          * @param {Array} q1 First point of segment 2 using normalized homogeneous coordinates [1,x,y]
2257          * @param {Array} q2 Second point or direction of segment 2 using normalized homogeneous coordinates [1,x,y] or point at infinity [0,x,y], respectively
2258          * @returns {Array} [Intersection point, t, u] The first entry contains the homogeneous coordinates
2259          * of the intersection point. The second and third entry give the position of the intersection with respect
2260          * to the definiting parameters. For example, the second entry t is defined by: intersection point = p1 + t * deltaP, where
2261          * deltaP = (p2 - p1) when both parameters are coordinates, and deltaP = p2 if p2 is a point at infinity.
2262          * If the two segments are collinear, [[0,0,0], Infinity, Infinity] is returned.
2263          **/
2264         meetSegmentSegment: function (p1, p2, q1, q2) {
2265             var t,
2266                 u,
2267                 i,
2268                 d,
2269                 li1 = Mat.crossProduct(p1, p2),
2270                 li2 = Mat.crossProduct(q1, q2),
2271                 c = Mat.crossProduct(li1, li2);
2272 
2273             if (Math.abs(c[0]) < Mat.eps) {
2274                 return [c, Infinity, Infinity];
2275             }
2276 
2277             // Normalize the intersection coordinates
2278             c[1] /= c[0];
2279             c[2] /= c[0];
2280             c[0] /= c[0];
2281 
2282             // Now compute in principle:
2283             //    t = dist(c - p1) / dist(p2 - p1) and
2284             //    u = dist(c - q1) / dist(q2 - q1)
2285             // However: the points q1, q2, p1, p2 might be ideal points - or in general - the
2286             // coordinates might be not normalized.
2287             // Note that the z-coordinates of p2 and q2 are used to determine whether it should be interpreted
2288             // as a segment coordinate or a direction.
2289             i = Math.abs(p2[1] - p2[0] * p1[1]) < Mat.eps ? 2 : 1;
2290             d = p1[i] / p1[0];
2291             t = (c[i] - d) / (p2[0] !== 0 ? p2[i] / p2[0] - d : p2[i]);
2292 
2293             i = Math.abs(q2[1] - q2[0] * q1[1]) < Mat.eps ? 2 : 1;
2294             d = q1[i] / q1[0];
2295             u = (c[i] - d) / (q2[0] !== 0 ? q2[i] / q2[0] - d : q2[i]);
2296 
2297             return [c, t, u];
2298         },
2299 
2300         /**
2301          * Find the n-th intersection point of two pathes, usually given by polygons. Uses parts of the
2302          * Greiner-Hormann algorithm in JXG.Math.Clip.
2303          *
2304          * @param {JXG.Circle|JXG.Curve|JXG.Polygon} path1
2305          * @param {JXG.Circle|JXG.Curve|JXG.Polygon} path2
2306          * @param {Number|Function} n
2307          * @param {JXG.Board} board
2308          *
2309          * @returns {JXG.Coords} Intersection point. In case no intersection point is detected,
2310          * the ideal point [0,0,0] is returned.
2311          *
2312          */
2313         meetPathPath: function (path1, path2, nr, board) {
2314             var S, C, len, intersections,
2315                 n = Type.evaluate(nr);
2316 
2317             S = JXG.Math.Clip._getPath(path1, board);
2318             len = S.length;
2319             if (
2320                 len > 0 &&
2321                 this.distance(S[0].coords.usrCoords, S[len - 1].coords.usrCoords, 3) < Mat.eps
2322             ) {
2323                 S.pop();
2324             }
2325 
2326             C = JXG.Math.Clip._getPath(path2, board);
2327             len = C.length;
2328             if (
2329                 len > 0 &&
2330                 this.distance(C[0].coords.usrCoords, C[len - 1].coords.usrCoords, 3) <
2331                 Mat.eps * Mat.eps
2332             ) {
2333                 C.pop();
2334             }
2335 
2336             // Handle cases where at least one of the paths is empty
2337             if (nr < 0 || JXG.Math.Clip.isEmptyCase(S, C, "intersection")) {
2338                 return new Coords(Const.COORDS_BY_USER, [0, 0, 0], board);
2339             }
2340 
2341             JXG.Math.Clip.makeDoublyLinkedList(S);
2342             JXG.Math.Clip.makeDoublyLinkedList(C);
2343 
2344             intersections = JXG.Math.Clip.findIntersections(S, C, board)[0];
2345             if (n < intersections.length) {
2346                 return intersections[n].coords;
2347             }
2348             return new Coords(Const.COORDS_BY_USER, [0, 0, 0], board);
2349         },
2350 
2351         /**
2352          * Find the n-th intersection point between a polygon and a line.
2353          * @param {JXG.Polygon} path
2354          * @param {JXG.Line} line
2355          * @param {Number|Function} nr
2356          * @param {JXG.Board} board
2357          * @param {Boolean} alwaysIntersect If false just the segment between the two defining points of the line are tested for intersection.
2358          *
2359          * @returns {JXG.Coords} Intersection point. In case no intersection point is detected,
2360          * the ideal point [0,0,0] is returned.
2361          */
2362         meetPolygonLine: function (path, line, nr, board, alwaysIntersect) {
2363             var i,
2364                 n = Type.evaluate(nr),
2365                 res,
2366                 border,
2367                 crds = [0, 0, 0],
2368                 len = path.borders.length,
2369                 intersections = [];
2370 
2371             for (i = 0; i < len; i++) {
2372                 border = path.borders[i];
2373                 res = this.meetSegmentSegment(
2374                     border.point1.coords.usrCoords,
2375                     border.point2.coords.usrCoords,
2376                     line.point1.coords.usrCoords,
2377                     line.point2.coords.usrCoords
2378                 );
2379 
2380                 if (
2381                     (!alwaysIntersect || (res[2] >= 0 && res[2] < 1)) &&
2382                     res[1] >= 0 &&
2383                     res[1] < 1
2384                 ) {
2385                     intersections.push(res[0]);
2386                 }
2387             }
2388 
2389             if (n >= 0 && n < intersections.length) {
2390                 crds = intersections[n];
2391             }
2392             return new Coords(Const.COORDS_BY_USER, crds, board);
2393         },
2394 
2395         meetPlaneSphere: function (el1, el2) {
2396             var dis = function () {
2397                 return (
2398                     el1.normal[0] * el2.center.X()
2399                     + el1.normal[1] * el2.center.Y()
2400                     + el1.normal[2] * el2.center.Z()
2401                     - el1.d
2402                 );
2403             };
2404             return [
2405                 [
2406                     // Center
2407                     function () {
2408                         return el2.center.X() - dis() * el1.normal[0];
2409                     },
2410                     function () {
2411                         return el2.center.Y() - dis() * el1.normal[1];
2412                     },
2413                     function () {
2414                         return el2.center.Z() - dis() * el1.normal[2];
2415                     }
2416                 ],
2417                 [
2418                     // Normal
2419                     () => el1.normal[0],
2420                     () => el1.normal[1],
2421                     () => el1.normal[2]
2422                 ],
2423                 function () {
2424                     // Radius (returns NaN if spheres don't touch)
2425                     var r = el2.Radius(),
2426                         s = dis();
2427                     return Math.sqrt(r * r - s * s);
2428                 }
2429             ];
2430         },
2431 
2432         meetSphereSphere: function (el1, el2) {
2433             var skew = function () {
2434                 var dist = el1.center.distance(el2.center),
2435                     r1 = el1.Radius(),
2436                     r2 = el2.Radius();
2437                 return (r1 - r2) * (r1 + r2) / (dist * dist);
2438             };
2439             return [
2440                 [
2441                     // Center
2442                     function () {
2443                         var s = skew();
2444                         return 0.5 * ((1 - s) * el1.center.X() + (1 + s) * el2.center.X());
2445                     },
2446                     function () {
2447                         var s = skew();
2448                         return 0.5 * ((1 - s) * el1.center.Y() + (1 + s) * el2.center.Y());
2449                     },
2450                     function () {
2451                         var s = skew();
2452                         return 0.5 * ((1 - s) * el1.center.Z() + (1 + s) * el2.center.Z());
2453                     }
2454                 ],
2455                 [
2456                     // Normal
2457                     () => el2.center.X() - el1.center.X(),
2458                     () => el2.center.Y() - el1.center.Y(),
2459                     () => el2.center.Z() - el1.center.Z()
2460                 ],
2461                 function () {
2462                     // Radius (returns NaN if spheres don't touch)
2463                     var dist = el1.center.distance(el2.center),
2464                         r1 = el1.Radius(),
2465                         r2 = el2.Radius(),
2466                         s = skew(),
2467                         rIxnSq = 0.5 * (r1 * r1 + r2 * r2 - 0.5 * dist * dist * (1 + s * s));
2468                     return Math.sqrt(rIxnSq);
2469                 }
2470             ];
2471         },
2472 
2473         /****************************************/
2474         /****   BEZIER CURVE ALGORITHMS      ****/
2475         /****************************************/
2476 
2477         /**
2478          * Splits a Bezier curve segment defined by four points into
2479          * two Bezier curve segments. Dissection point is t=1/2.
2480          * @param {Array} curve Array of four coordinate arrays of length 2 defining a
2481          * Bezier curve segment, i.e. [[x0,y0], [x1,y1], [x2,y2], [x3,y3]].
2482          * @returns {Array} Array consisting of two coordinate arrays for Bezier curves.
2483          */
2484         _bezierSplit: function (curve) {
2485             var p0, p1, p2, p00, p22, p000;
2486 
2487             p0 = [(curve[0][0] + curve[1][0]) * 0.5, (curve[0][1] + curve[1][1]) * 0.5];
2488             p1 = [(curve[1][0] + curve[2][0]) * 0.5, (curve[1][1] + curve[2][1]) * 0.5];
2489             p2 = [(curve[2][0] + curve[3][0]) * 0.5, (curve[2][1] + curve[3][1]) * 0.5];
2490 
2491             p00 = [(p0[0] + p1[0]) * 0.5, (p0[1] + p1[1]) * 0.5];
2492             p22 = [(p1[0] + p2[0]) * 0.5, (p1[1] + p2[1]) * 0.5];
2493 
2494             p000 = [(p00[0] + p22[0]) * 0.5, (p00[1] + p22[1]) * 0.5];
2495 
2496             return [
2497                 [curve[0], p0, p00, p000],
2498                 [p000, p22, p2, curve[3]]
2499             ];
2500         },
2501 
2502         /**
2503          * Computes the bounding box [minX, maxY, maxX, minY] of a Bezier curve segment
2504          * from its control points.
2505          * @param {Array} curve Array of four coordinate arrays of length 2 defining a
2506          * Bezier curve segment, i.e. [[x0,y0], [x1,y1], [x2,y2], [x3,y3]].
2507          * @returns {Array} Bounding box [minX, maxY, maxX, minY]
2508          */
2509         _bezierBbox: function (curve) {
2510             var bb = [];
2511 
2512             if (curve.length === 4) {
2513                 // bezierDegree == 3
2514                 bb[0] = Math.min(curve[0][0], curve[1][0], curve[2][0], curve[3][0]); // minX
2515                 bb[1] = Math.max(curve[0][1], curve[1][1], curve[2][1], curve[3][1]); // maxY
2516                 bb[2] = Math.max(curve[0][0], curve[1][0], curve[2][0], curve[3][0]); // maxX
2517                 bb[3] = Math.min(curve[0][1], curve[1][1], curve[2][1], curve[3][1]); // minY
2518             } else {
2519                 // bezierDegree == 1
2520                 bb[0] = Math.min(curve[0][0], curve[1][0]); // minX
2521                 bb[1] = Math.max(curve[0][1], curve[1][1]); // maxY
2522                 bb[2] = Math.max(curve[0][0], curve[1][0]); // maxX
2523                 bb[3] = Math.min(curve[0][1], curve[1][1]); // minY
2524             }
2525 
2526             return bb;
2527         },
2528 
2529         /**
2530          * Decide if two Bezier curve segments overlap by comparing their bounding boxes.
2531          * @param {Array} bb1 Bounding box of the first Bezier curve segment
2532          * @param {Array} bb2 Bounding box of the second Bezier curve segment
2533          * @returns {Boolean} true if the bounding boxes overlap, false otherwise.
2534          */
2535         _bezierOverlap: function (bb1, bb2) {
2536             return bb1[2] >= bb2[0] && bb1[0] <= bb2[2] && bb1[1] >= bb2[3] && bb1[3] <= bb2[1];
2537         },
2538 
2539         /**
2540          * Append list of intersection points to a list.
2541          * @private
2542          */
2543         _bezierListConcat: function (L, Lnew, t1, t2) {
2544             var i,
2545                 t2exists = Type.exists(t2),
2546                 start = 0,
2547                 len = Lnew.length,
2548                 le = L.length;
2549 
2550             if (
2551                 le > 0 &&
2552                 len > 0 &&
2553                 ((L[le - 1][1] === 1 && Lnew[0][1] === 0) ||
2554                     (t2exists && L[le - 1][2] === 1 && Lnew[0][2] === 0))
2555             ) {
2556                 start = 1;
2557             }
2558 
2559             for (i = start; i < len; i++) {
2560                 if (t2exists) {
2561                     Lnew[i][2] *= 0.5;
2562                     Lnew[i][2] += t2;
2563                 }
2564 
2565                 Lnew[i][1] *= 0.5;
2566                 Lnew[i][1] += t1;
2567 
2568                 L.push(Lnew[i]);
2569             }
2570         },
2571 
2572         /**
2573          * Find intersections of two Bezier curve segments by recursive subdivision.
2574          * Below maxlevel determine intersections by intersection line segments.
2575          * @param {Array} red Array of four coordinate arrays of length 2 defining the first
2576          * Bezier curve segment, i.e. [[x0,y0], [x1,y1], [x2,y2], [x3,y3]].
2577          * @param {Array} blue Array of four coordinate arrays of length 2 defining the second
2578          * Bezier curve segment, i.e. [[x0,y0], [x1,y1], [x2,y2], [x3,y3]].
2579          * @param {Number} level Recursion level
2580          * @returns {Array} List of intersection points (up to nine). Each intersection point is an
2581          * array of length three (homogeneous coordinates) plus preimages.
2582          */
2583         _bezierMeetSubdivision: function (red, blue, level) {
2584             var bbb,
2585                 bbr,
2586                 ar,
2587                 b0,
2588                 b1,
2589                 r0,
2590                 r1,
2591                 m,
2592                 p0,
2593                 p1,
2594                 q0,
2595                 q1,
2596                 L = [],
2597                 maxLev = 5; // Maximum recursion level
2598 
2599             bbr = this._bezierBbox(blue);
2600             bbb = this._bezierBbox(red);
2601 
2602             if (!this._bezierOverlap(bbr, bbb)) {
2603                 return [];
2604             }
2605 
2606             if (level < maxLev) {
2607                 ar = this._bezierSplit(red);
2608                 r0 = ar[0];
2609                 r1 = ar[1];
2610 
2611                 ar = this._bezierSplit(blue);
2612                 b0 = ar[0];
2613                 b1 = ar[1];
2614 
2615                 this._bezierListConcat(
2616                     L,
2617                     this._bezierMeetSubdivision(r0, b0, level + 1),
2618                     0.0,
2619                     0.0
2620                 );
2621                 this._bezierListConcat(
2622                     L,
2623                     this._bezierMeetSubdivision(r0, b1, level + 1),
2624                     0,
2625                     0.5
2626                 );
2627                 this._bezierListConcat(
2628                     L,
2629                     this._bezierMeetSubdivision(r1, b0, level + 1),
2630                     0.5,
2631                     0.0
2632                 );
2633                 this._bezierListConcat(
2634                     L,
2635                     this._bezierMeetSubdivision(r1, b1, level + 1),
2636                     0.5,
2637                     0.5
2638                 );
2639 
2640                 return L;
2641             }
2642 
2643             // Make homogeneous coordinates
2644             q0 = [1].concat(red[0]);
2645             q1 = [1].concat(red[3]);
2646             p0 = [1].concat(blue[0]);
2647             p1 = [1].concat(blue[3]);
2648 
2649             m = this.meetSegmentSegment(q0, q1, p0, p1);
2650 
2651             if (m[1] >= 0.0 && m[2] >= 0.0 && m[1] <= 1.0 && m[2] <= 1.0) {
2652                 return [m];
2653             }
2654 
2655             return [];
2656         },
2657 
2658         /**
2659          * @param {Boolean} testSegment Test if intersection has to be inside of the segment or somewhere on the line defined by the segment
2660          */
2661         _bezierLineMeetSubdivision: function (red, blue, level, testSegment) {
2662             var bbb, bbr, ar,
2663                 r0, r1,
2664                 m,
2665                 p0, p1, q0, q1,
2666                 L = [],
2667                 maxLev = 5; // Maximum recursion level
2668 
2669             bbb = this._bezierBbox(blue);
2670             bbr = this._bezierBbox(red);
2671 
2672             if (testSegment && !this._bezierOverlap(bbr, bbb)) {
2673                 return [];
2674             }
2675 
2676             if (level < maxLev) {
2677                 ar = this._bezierSplit(red);
2678                 r0 = ar[0];
2679                 r1 = ar[1];
2680 
2681                 this._bezierListConcat(
2682                     L,
2683                     this._bezierLineMeetSubdivision(r0, blue, level + 1),
2684                     0.0
2685                 );
2686                 this._bezierListConcat(
2687                     L,
2688                     this._bezierLineMeetSubdivision(r1, blue, level + 1),
2689                     0.5
2690                 );
2691 
2692                 return L;
2693             }
2694 
2695             // Make homogeneous coordinates
2696             q0 = [1].concat(red[0]);
2697             q1 = [1].concat(red[3]);
2698             p0 = [1].concat(blue[0]);
2699             p1 = [1].concat(blue[1]);
2700 
2701             m = this.meetSegmentSegment(q0, q1, p0, p1);
2702 
2703             if (m[1] >= 0.0 && m[1] <= 1.0) {
2704                 if (!testSegment || (m[2] >= 0.0 && m[2] <= 1.0)) {
2705                     return [m];
2706                 }
2707             }
2708 
2709             return [];
2710         },
2711 
2712         /**
2713          * Find the nr-th intersection point of two Bezier curve segments.
2714          * @param {Array} red Array of four coordinate arrays of length 2 defining the first
2715          * Bezier curve segment, i.e. [[x0,y0], [x1,y1], [x2,y2], [x3,y3]].
2716          * @param {Array} blue Array of four coordinate arrays of length 2 defining the second
2717          * Bezier curve segment, i.e. [[x0,y0], [x1,y1], [x2,y2], [x3,y3]].
2718          * @param {Boolean} testSegment Test if intersection has to be inside of the segment or somewhere on the line defined by the segment
2719          * @returns {Array} Array containing the list of all intersection points as homogeneous coordinate arrays plus
2720          * preimages [x,y], t_1, t_2] of the two Bezier curve segments.
2721          *
2722          */
2723         meetBeziersegmentBeziersegment: function (red, blue, testSegment) {
2724             var L, L2, i;
2725 
2726             if (red.length === 4 && blue.length === 4) {
2727                 L = this._bezierMeetSubdivision(red, blue, 0);
2728             } else {
2729                 L = this._bezierLineMeetSubdivision(red, blue, 0, testSegment);
2730             }
2731 
2732             L.sort(function (a, b) {
2733                 return (a[1] - b[1]) * 10000000.0 + (a[2] - b[2]);
2734             });
2735 
2736             L2 = [];
2737             for (i = 0; i < L.length; i++) {
2738                 // Only push entries different from their predecessor
2739                 if (i === 0 || L[i][1] !== L[i - 1][1] || L[i][2] !== L[i - 1][2]) {
2740                     L2.push(L[i]);
2741                 }
2742             }
2743             return L2;
2744         },
2745 
2746         /**
2747          * Find the nr-th intersection point of two Bezier curves, i.e. curves with bezierDegree == 3.
2748          * @param {JXG.Curve} red Curve with bezierDegree == 3
2749          * @param {JXG.Curve} blue Curve with bezierDegree == 3
2750          * @param {Number|Function} nr The number of the intersection point which should be returned.
2751          * @returns {Array} The homogeneous coordinates of the nr-th intersection point.
2752          */
2753         meetBezierCurveRedBlueSegments: function (red, blue, nr) {
2754             var p, i, j, k,
2755                 n = Type.evaluate(nr),
2756                 po, tmp,
2757                 redArr,
2758                 blueArr,
2759                 bbr,
2760                 bbb,
2761                 intersections,
2762                 startRed = 0,
2763                 startBlue = 0,
2764                 lenBlue, lenRed,
2765                 L = [];
2766 
2767             if (blue.numberPoints < blue.bezierDegree + 1 || red.numberPoints < red.bezierDegree + 1) {
2768                 return [0, NaN, NaN];
2769             }
2770             if (red.bezierDegree === 1 && blue.bezierDegree === 3) {
2771                 tmp = red;
2772                 red = blue;
2773                 blue = tmp;
2774             }
2775 
2776             lenBlue = blue.numberPoints - blue.bezierDegree;
2777             lenRed = red.numberPoints - red.bezierDegree;
2778 
2779             // For sectors, we ignore the "legs"
2780             if (red.type === Const.OBJECT_TYPE_SECTOR) {
2781                 startRed = 3;
2782                 lenRed -= 3;
2783             }
2784             if (blue.type === Const.OBJECT_TYPE_SECTOR) {
2785                 startBlue = 3;
2786                 lenBlue -= 3;
2787             }
2788 
2789             for (i = startRed; i < lenRed; i += red.bezierDegree) {
2790                 p = red.points;
2791                 redArr = [p[i].usrCoords.slice(1), p[i + 1].usrCoords.slice(1)];
2792                 if (red.bezierDegree === 3) {
2793                     redArr[2] = p[i + 2].usrCoords.slice(1);
2794                     redArr[3] = p[i + 3].usrCoords.slice(1);
2795                 }
2796 
2797                 bbr = this._bezierBbox(redArr);
2798 
2799                 for (j = startBlue; j < lenBlue; j += blue.bezierDegree) {
2800                     p = blue.points;
2801                     blueArr = [p[j].usrCoords.slice(1), p[j + 1].usrCoords.slice(1)];
2802                     if (blue.bezierDegree === 3) {
2803                         blueArr[2] = p[j + 2].usrCoords.slice(1);
2804                         blueArr[3] = p[j + 3].usrCoords.slice(1);
2805                     }
2806 
2807                     bbb = this._bezierBbox(blueArr);
2808                     if (this._bezierOverlap(bbr, bbb)) {
2809                         intersections = this.meetBeziersegmentBeziersegment(redArr, blueArr);
2810                         if (intersections.length === 0) {
2811                             continue;
2812                         }
2813                         for (k = 0; k < intersections.length; k++) {
2814                             po = intersections[k];
2815                             if (
2816                                 po[1] < -Mat.eps ||
2817                                 po[1] > 1 + Mat.eps ||
2818                                 po[2] < -Mat.eps ||
2819                                 po[2] > 1 + Mat.eps
2820                             ) {
2821                                 continue;
2822                             }
2823                             L.push(po);
2824                         }
2825                         if (L.length > n) {
2826                             return L[n][0];
2827                         }
2828                     }
2829                 }
2830             }
2831             if (L.length > n) {
2832                 return L[n][0];
2833             }
2834 
2835             return [0, NaN, NaN];
2836         },
2837 
2838         bezierSegmentEval: function (t, curve) {
2839             var f,
2840                 x,
2841                 y,
2842                 t1 = 1.0 - t;
2843 
2844             x = 0;
2845             y = 0;
2846 
2847             f = t1 * t1 * t1;
2848             x += f * curve[0][0];
2849             y += f * curve[0][1];
2850 
2851             f = 3.0 * t * t1 * t1;
2852             x += f * curve[1][0];
2853             y += f * curve[1][1];
2854 
2855             f = 3.0 * t * t * t1;
2856             x += f * curve[2][0];
2857             y += f * curve[2][1];
2858 
2859             f = t * t * t;
2860             x += f * curve[3][0];
2861             y += f * curve[3][1];
2862 
2863             return [1.0, x, y];
2864         },
2865 
2866         /**
2867          * Generate the defining points of a 3rd degree bezier curve that approximates
2868          * a circle sector defined by three coordinate points A, B, C, each defined by an array of length three.
2869          * The coordinate arrays are given in homogeneous coordinates.
2870          * @param {Array} A First point
2871          * @param {Array} B Second point (intersection point)
2872          * @param {Array} C Third point
2873          * @param {Boolean} withLegs Flag. If true the legs to the intersection point are part of the curve.
2874          * @param {Number} sgn Wither 1 or -1. Needed for minor and major arcs. In case of doubt, use 1.
2875          */
2876         bezierArc: function (A, B, C, withLegs, sgn) {
2877             var p1, p2, p3, p4,
2878                 r,
2879                 phi, beta, delta,
2880                 // PI2 = Math.PI * 0.5,
2881                 x = B[1],
2882                 y = B[2],
2883                 z = B[0],
2884                 dataX = [],
2885                 dataY = [],
2886                 co, si,
2887                 ax, ay,
2888                 bx, by,
2889                 k, v, d,
2890                 matrix;
2891 
2892             r = this.distance(B, A);
2893 
2894             // x,y, z is intersection point. Normalize it.
2895             x /= z;
2896             y /= z;
2897 
2898             phi = this.rad(A.slice(1), B.slice(1), C.slice(1));
2899             if (sgn === -1) {
2900                 phi = 2 * Math.PI - phi;
2901             }
2902 
2903             // Always divide the arc into four Bezier arcs.
2904             // Otherwise, the position of gliders on this arc
2905             // will be wrong.
2906             delta = phi / 4;
2907 
2908 
2909             p1 = A;
2910             p1[1] /= p1[0];
2911             p1[2] /= p1[0];
2912             p1[0] /= p1[0];
2913 
2914             p4 = p1.slice(0);
2915 
2916             if (withLegs) {
2917                 dataX = [x, x + 0.333 * (p1[1] - x), x + 0.666 * (p1[1] - x), p1[1]];
2918                 dataY = [y, y + 0.333 * (p1[2] - y), y + 0.666 * (p1[2] - y), p1[2]];
2919             } else {
2920                 dataX = [p1[1]];
2921                 dataY = [p1[2]];
2922             }
2923 
2924             while (phi > Mat.eps) {
2925                 // if (phi > PI2) {
2926                 //     beta = PI2;
2927                 //     phi -= PI2;
2928                 // } else {
2929                 //     beta = phi;
2930                 //     phi = 0;
2931                 // }
2932                 if (phi > delta) {
2933                     beta = delta;
2934                     phi -= delta;
2935                 } else {
2936                     beta = phi;
2937                     phi = 0;
2938                 }
2939 
2940                 co = Math.cos(sgn * beta);
2941                 si = Math.sin(sgn * beta);
2942 
2943                 matrix = [
2944                     [1, 0, 0],
2945                     [x * (1 - co) + y * si, co, -si],
2946                     [y * (1 - co) - x * si, si, co]
2947                 ];
2948                 v = Mat.matVecMult(matrix, p1);
2949                 p4 = [v[0] / v[0], v[1] / v[0], v[2] / v[0]];
2950 
2951                 ax = p1[1] - x;
2952                 ay = p1[2] - y;
2953                 bx = p4[1] - x;
2954                 by = p4[2] - y;
2955                 d = Mat.hypot(ax + bx, ay + by);
2956 
2957                 if (Math.abs(by - ay) > Mat.eps) {
2958                     k = ((((ax + bx) * (r / d - 0.5)) / (by - ay)) * 8) / 3;
2959                 } else {
2960                     k = ((((ay + by) * (r / d - 0.5)) / (ax - bx)) * 8) / 3;
2961                 }
2962 
2963                 p2 = [1, p1[1] - k * ay, p1[2] + k * ax];
2964                 p3 = [1, p4[1] + k * by, p4[2] - k * bx];
2965 
2966                 Type.concat(dataX, [p2[1], p3[1], p4[1]]);
2967                 Type.concat(dataY, [p2[2], p3[2], p4[2]]);
2968                 p1 = p4.slice(0);
2969             }
2970 
2971             if (withLegs) {
2972                 Type.concat(dataX, [
2973                     p4[1] + 0.333 * (x - p4[1]),
2974                     p4[1] + 0.666 * (x - p4[1]),
2975                     x
2976                 ]);
2977                 Type.concat(dataY, [
2978                     p4[2] + 0.333 * (y - p4[2]),
2979                     p4[2] + 0.666 * (y - p4[2]),
2980                     y
2981                 ]);
2982             }
2983 
2984             return [dataX, dataY];
2985         },
2986 
2987         /****************************************/
2988         /****           PROJECTIONS          ****/
2989         /****************************************/
2990 
2991         /**
2992          * Calculates the coordinates of the projection of a given point on a given circle. I.o.w. the
2993          * nearest one of the two intersection points of the line through the given point and the circles
2994          * center.
2995          * @param {JXG.Point|JXG.Coords} point Point to project or coords object to project.
2996          * @param {JXG.Circle} circle Circle on that the point is projected.
2997          * @param {JXG.Board} [board=point.board] Reference to the board
2998          * @returns {JXG.Coords} The coordinates of the projection of the given point on the given circle.
2999          */
3000         projectPointToCircle: function (point, circle, board) {
3001             var dist,
3002                 P,
3003                 x,
3004                 y,
3005                 factor,
3006                 M = circle.center.coords.usrCoords;
3007 
3008             if (!Type.exists(board)) {
3009                 board = point.board;
3010             }
3011 
3012             // gave us a point
3013             if (Type.isPoint(point)) {
3014                 dist = point.coords.distance(Const.COORDS_BY_USER, circle.center.coords);
3015                 P = point.coords.usrCoords;
3016                 // gave us coords
3017             } else {
3018                 dist = point.distance(Const.COORDS_BY_USER, circle.center.coords);
3019                 P = point.usrCoords;
3020             }
3021 
3022             if (Math.abs(dist) < Mat.eps) {
3023                 dist = Mat.eps;
3024             }
3025 
3026             factor = circle.Radius() / dist;
3027             x = M[1] + factor * (P[1] - M[1]);
3028             y = M[2] + factor * (P[2] - M[2]);
3029 
3030             return new Coords(Const.COORDS_BY_USER, [x, y], board);
3031         },
3032 
3033         /**
3034          * Calculates the coordinates of the orthogonal projection of a given point on a given line. I.o.w. the
3035          * intersection point of the given line and its perpendicular through the given point.
3036          * @param {JXG.Point|JXG.Coords} point Point to project.
3037          * @param {JXG.Line} line Line on that the point is projected.
3038          * @param {JXG.Board} [board=point.board|board=line.board] Reference to a board.
3039          * @returns {JXG.Coords} The coordinates of the projection of the given point on the given line.
3040          */
3041         projectPointToLine: function (point, line, board) {
3042             var v = [0, line.stdform[1], line.stdform[2]],
3043                 coords;
3044 
3045             if (!Type.exists(board)) {
3046                 if (Type.exists(point.coords)) {
3047                     board = point.board;
3048                 } else {
3049                     board = line.board;
3050                 }
3051             }
3052 
3053             if (Type.exists(point.coords)) {
3054                 coords = point.coords.usrCoords;
3055             } else {
3056                 coords = point.usrCoords;
3057             }
3058 
3059             v = Mat.crossProduct(v, coords);
3060             return new Coords(Const.COORDS_BY_USER, Mat.crossProduct(v, line.stdform), board);
3061         },
3062 
3063         /**
3064          * Calculates the coordinates of the orthogonal projection of a given coordinate array on a given line
3065          * segment defined by two coordinate arrays.
3066          * @param {Array} p Point to project.
3067          * @param {Array} q1 Start point of the line segment on that the point is projected.
3068          * @param {Array} q2 End point of the line segment on that the point is projected.
3069          * @returns {Array} The coordinates of the projection of the given point on the given segment
3070          * and the factor that determines the projected point as a convex combination of the
3071          * two endpoints q1 and q2 of the segment.
3072          */
3073         projectCoordsToSegment: function (p, q1, q2) {
3074             var t,
3075                 denom,
3076                 s = [q2[1] - q1[1], q2[2] - q1[2]],
3077                 v = [p[1] - q1[1], p[2] - q1[2]];
3078 
3079             /**
3080              * If the segment has length 0, i.e. is a point,
3081              * the projection is equal to that point.
3082              */
3083             if (Math.abs(s[0]) < Mat.eps && Math.abs(s[1]) < Mat.eps) {
3084                 return [q1, 0];
3085             }
3086 
3087             t = Mat.innerProduct(v, s);
3088             denom = Mat.innerProduct(s, s);
3089             t /= denom;
3090 
3091             return [[1, t * s[0] + q1[1], t * s[1] + q1[2]], t];
3092         },
3093 
3094         /**
3095          * Finds the coordinates of the closest point on a Bezier segment of a
3096          * {@link JXG.Curve} to a given coordinate array.
3097          * @param {Array} pos Point to project in homogeneous coordinates.
3098          * @param {JXG.Curve} curve Curve of type "plot" having Bezier degree 3.
3099          * @param {Number} start Number of the Bezier segment of the curve.
3100          * @returns {Array} The coordinates of the projection of the given point
3101          * on the given Bezier segment and the preimage of the curve which
3102          * determines the closest point.
3103          */
3104         projectCoordsToBeziersegment: function (pos, curve, start) {
3105             var t0,
3106                 /** @ignore */
3107                 minfunc = function (t) {
3108                     var z = [1, curve.X(start + t), curve.Y(start + t)];
3109 
3110                     z[1] -= pos[1];
3111                     z[2] -= pos[2];
3112 
3113                     return z[1] * z[1] + z[2] * z[2];
3114                 };
3115 
3116             t0 = JXG.Math.Numerics.fminbr(minfunc, [0.0, 1.0]);
3117 
3118             return [[1, curve.X(t0 + start), curve.Y(t0 + start)], t0];
3119         },
3120 
3121         /**
3122          * Calculates the coordinates of the projection of a given point on a given curve.
3123          * Uses {@link JXG.Math.Geometry.projectCoordsToCurve}.
3124          *
3125          * @param {JXG.Point} point Point to project.
3126          * @param {JXG.Curve} curve Curve on that the point is projected.
3127          * @param {JXG.Board} [board=point.board] Reference to a board.
3128          * @see #projectCoordsToCurve
3129          * @returns {Array} [JXG.Coords, position] The coordinates of the projection of the given
3130          * point on the given graph and the relative position on the curve (real number).
3131          */
3132         projectPointToCurve: function (point, curve, board) {
3133             if (!Type.exists(board)) {
3134                 board = point.board;
3135             }
3136 
3137             var x = point.X(),
3138                 y = point.Y(),
3139                 t = point.position,
3140                 result;
3141 
3142             if (!Type.exists(t)) {
3143                 t = Type.evaluate(curve.visProp.curvetype) === 'functiongraph' ? x : 0.0;
3144             }
3145             result = this.projectCoordsToCurve(x, y, t, curve, board);
3146 
3147             // point.position = result[1];
3148 
3149             return result;
3150         },
3151 
3152         /**
3153          * Calculates the coordinates of the projection of a coordinates pair on a given curve. In case of
3154          * function graphs this is the
3155          * intersection point of the curve and the parallel to y-axis through the given point.
3156          * @param {Number} x coordinate to project.
3157          * @param {Number} y coordinate to project.
3158          * @param {Number} t start value for newtons method
3159          * @param {JXG.Curve} curve Curve on that the point is projected.
3160          * @param {JXG.Board} [board=curve.board] Reference to a board.
3161          * @see #projectPointToCurve
3162          * @returns {JXG.Coords} Array containing the coordinates of the projection of the given point on the given curve and
3163          * the position on the curve.
3164          */
3165         projectCoordsToCurve: function (x, y, t, curve, board) {
3166             var newCoords, newCoordsObj,
3167                 i, j, mindist, dist, lbda,
3168                 v, coords, d, p1, p2, res, minfunc,
3169                 t_new, f_new, f_old, dy,
3170                 delta, delta1, delta2, steps, minX, maxX,
3171                 infty = Number.POSITIVE_INFINITY;
3172 
3173             if (!Type.exists(board)) {
3174                 board = curve.board;
3175             }
3176 
3177             if (Type.evaluate(curve.visProp.curvetype) === "plot") {
3178                 t = 0;
3179                 mindist = infty;
3180                 if (curve.numberPoints === 0) {
3181                     newCoords = [0, 1, 1];
3182                 } else {
3183                     newCoords = [curve.Z(0), curve.X(0), curve.Y(0)];
3184                 }
3185 
3186                 if (curve.numberPoints > 1) {
3187                     v = [1, x, y];
3188                     if (curve.bezierDegree === 3) {
3189                         j = 0;
3190                     } else {
3191                         p1 = [curve.Z(0), curve.X(0), curve.Y(0)];
3192                     }
3193                     for (i = 0; i < curve.numberPoints - 1; i++) {
3194                         if (curve.bezierDegree === 3) {
3195                             res = this.projectCoordsToBeziersegment(v, curve, j);
3196                         } else {
3197                             p2 = [curve.Z(i + 1), curve.X(i + 1), curve.Y(i + 1)];
3198                             res = this.projectCoordsToSegment(v, p1, p2);
3199                         }
3200                         lbda = res[1];
3201                         coords = res[0];
3202 
3203                         if (0.0 <= lbda && lbda <= 1.0) {
3204                             dist = this.distance(coords, v);
3205                             d = i + lbda;
3206                         } else if (lbda < 0.0) {
3207                             coords = p1;
3208                             dist = this.distance(p1, v);
3209                             d = i;
3210                         } else if (lbda > 1.0 && i === curve.numberPoints - 2) {
3211                             coords = p2;
3212                             dist = this.distance(coords, v);
3213                             d = curve.numberPoints - 1;
3214                         }
3215 
3216                         if (dist < mindist) {
3217                             mindist = dist;
3218                             t = d;
3219                             newCoords = coords;
3220                         }
3221 
3222                         if (curve.bezierDegree === 3) {
3223                             j++;
3224                             i += 2;
3225                         } else {
3226                             p1 = p2;
3227                         }
3228                     }
3229                 }
3230 
3231                 newCoordsObj = new Coords(Const.COORDS_BY_USER, newCoords, board);
3232             } else {
3233                 // 'parameter', 'polar', 'functiongraph'
3234 
3235                 if (Type.evaluate(curve.visProp.curvetype) === 'functiongraph') {
3236                     dy = Math.abs(y - curve.Y(x));
3237                     if (!isNaN(dy)) {
3238                         minX = x - dy;
3239                         maxX = x + dy;
3240                     } else {
3241                         minX = curve.minX();
3242                         maxX = curve.maxX();
3243                     }
3244                 } else {
3245                     minX = curve.minX();
3246                     maxX = curve.maxX();
3247                 }
3248 
3249                 /** @ignore */
3250                 minfunc = function (t) {
3251                     var dx, dy;
3252                     if (t < curve.minX() || t > curve.maxX()) {
3253                         return Infinity;
3254                     }
3255                     dx = x - curve.X(t);
3256                     dy = y - curve.Y(t);
3257                     return dx * dx + dy * dy;
3258                 };
3259 
3260                 f_old = minfunc(t);
3261                 steps = 50;
3262 
3263                 delta = (maxX - minX) / steps;
3264                 t_new = minX;
3265 
3266                 for (i = 0; i < steps; i++) {
3267                     f_new = minfunc(t_new);
3268 
3269                     if (f_new < f_old || f_old === Infinity || isNaN(f_old)) {
3270                         t = t_new;
3271                         f_old = f_new;
3272                     }
3273 
3274                     t_new += delta;
3275                 }
3276 
3277                 // t = Numerics.root(Numerics.D(minfunc), t);
3278                 // Ensure that minfunc is defined on the
3279                 // enclsoing interval [t-delta1, t+delta2]
3280                 delta1 = delta;
3281                 for (i = 0;
3282                     i < 20 && isNaN(minfunc(t - delta1));
3283                     i++, delta1 *= 0.5);
3284 
3285                 if (isNaN(minfunc(t - delta1))) {
3286                     delta1 = 0.0;
3287                 }
3288                 delta2 = delta;
3289                 for (i = 0;
3290                     i < 20 && isNaN(minfunc(t + delta2));
3291                     i++, delta2 *= 0.5);
3292                 if (isNaN(minfunc(t + delta2))) {
3293                     delta2 = 0.0;
3294                 }
3295 
3296                 t = Numerics.fminbr(minfunc, [
3297                     Math.max(t - delta1, minX),
3298                     Math.min(t + delta2, maxX)
3299                 ]);
3300 
3301                 // Distinction between closed and open curves is not necessary.
3302                 // If closed, the cyclic projection shift will work anyhow
3303                 // if (Math.abs(curve.X(minX) - curve.X(maxX)) < Mat.eps &&
3304                 //     Math.abs(curve.Y(minX) - curve.Y(maxX)) < Mat.eps) {
3305                 //     // Cyclically
3306                 //     if (t < minX) {console.log(t)
3307                 //         t = maxX + t - minX;
3308                 //     }
3309                 //     if (t > maxX) {
3310                 //         t = minX + t - maxX;
3311                 //     }
3312                 // } else {
3313                 t = t < minX ? minX : t;
3314                 t = t > maxX ? maxX : t;
3315                 // }
3316 
3317                 newCoordsObj = new Coords(
3318                     Const.COORDS_BY_USER,
3319                     [curve.X(t), curve.Y(t)],
3320                     board
3321                 );
3322             }
3323 
3324             return [curve.updateTransform(newCoordsObj), t];
3325         },
3326 
3327         /**
3328          * Calculates the coordinates of the closest orthogonal projection of a given coordinate array onto the
3329          * border of a polygon.
3330          * @param {Array} p Point to project.
3331          * @param {JXG.Polygon} pol Polygon element
3332          * @returns {Array} The coordinates of the closest projection of the given point to the border of the polygon.
3333          */
3334         projectCoordsToPolygon: function (p, pol) {
3335             var i,
3336                 len = pol.vertices.length,
3337                 d_best = Infinity,
3338                 d,
3339                 projection,
3340                 proj,
3341                 bestprojection;
3342 
3343             for (i = 0; i < len - 1; i++) {
3344                 projection = JXG.Math.Geometry.projectCoordsToSegment(
3345                     p,
3346                     pol.vertices[i].coords.usrCoords,
3347                     pol.vertices[i + 1].coords.usrCoords
3348                 );
3349 
3350                 if (0 <= projection[1] && projection[1] <= 1) {
3351                     d = JXG.Math.Geometry.distance(projection[0], p, 3);
3352                     proj = projection[0];
3353                 } else if (projection[1] < 0) {
3354                     d = JXG.Math.Geometry.distance(pol.vertices[i].coords.usrCoords, p, 3);
3355                     proj = pol.vertices[i].coords.usrCoords;
3356                 } else {
3357                     d = JXG.Math.Geometry.distance(pol.vertices[i + 1].coords.usrCoords, p, 3);
3358                     proj = pol.vertices[i + 1].coords.usrCoords;
3359                 }
3360                 if (d < d_best) {
3361                     bestprojection = proj.slice(0);
3362                     d_best = d;
3363                 }
3364             }
3365             return bestprojection;
3366         },
3367 
3368         /**
3369          * Calculates the coordinates of the projection of a given point on a given turtle. A turtle consists of
3370          * one or more curves of curveType 'plot'. Uses {@link JXG.Math.Geometry.projectPointToCurve}.
3371          * @param {JXG.Point} point Point to project.
3372          * @param {JXG.Turtle} turtle on that the point is projected.
3373          * @param {JXG.Board} [board=point.board] Reference to a board.
3374          * @returns {Array} [JXG.Coords, position] Array containing the coordinates of the projection of the given point on the turtle and
3375          * the position on the turtle.
3376          */
3377         projectPointToTurtle: function (point, turtle, board) {
3378             var newCoords,
3379                 t,
3380                 x,
3381                 y,
3382                 i,
3383                 dist,
3384                 el,
3385                 minEl,
3386                 res,
3387                 newPos,
3388                 np = 0,
3389                 npmin = 0,
3390                 mindist = Number.POSITIVE_INFINITY,
3391                 len = turtle.objects.length;
3392 
3393             if (!Type.exists(board)) {
3394                 board = point.board;
3395             }
3396 
3397             // run through all curves of this turtle
3398             for (i = 0; i < len; i++) {
3399                 el = turtle.objects[i];
3400 
3401                 if (el.elementClass === Const.OBJECT_CLASS_CURVE) {
3402                     res = this.projectPointToCurve(point, el);
3403                     newCoords = res[0];
3404                     newPos = res[1];
3405                     dist = this.distance(newCoords.usrCoords, point.coords.usrCoords);
3406 
3407                     if (dist < mindist) {
3408                         x = newCoords.usrCoords[1];
3409                         y = newCoords.usrCoords[2];
3410                         t = newPos;
3411                         mindist = dist;
3412                         minEl = el;
3413                         npmin = np;
3414                     }
3415                     np += el.numberPoints;
3416                 }
3417             }
3418 
3419             newCoords = new Coords(Const.COORDS_BY_USER, [x, y], board);
3420             // point.position = t + npmin;
3421             // return minEl.updateTransform(newCoords);
3422             return [minEl.updateTransform(newCoords), t + npmin];
3423         },
3424 
3425         /**
3426          * Trivial projection of a point to another point.
3427          * @param {JXG.Point} point Point to project (not used).
3428          * @param {JXG.Point} dest Point on that the point is projected.
3429          * @returns {JXG.Coords} The coordinates of the projection of the given point on the given circle.
3430          */
3431         projectPointToPoint: function (point, dest) {
3432             return dest.coords;
3433         },
3434 
3435         /**
3436          *
3437          * @param {JXG.Point|JXG.Coords} point
3438          * @param {JXG.Board} [board]
3439          */
3440         projectPointToBoard: function (point, board) {
3441             var i,
3442                 l,
3443                 c,
3444                 brd = board || point.board,
3445                 // comparison factor, point coord idx, bbox idx, 1st bbox corner x & y idx, 2nd bbox corner x & y idx
3446                 config = [
3447                     // left
3448                     [1, 1, 0, 0, 3, 0, 1],
3449                     // top
3450                     [-1, 2, 1, 0, 1, 2, 1],
3451                     // right
3452                     [-1, 1, 2, 2, 1, 2, 3],
3453                     // bottom
3454                     [1, 2, 3, 0, 3, 2, 3]
3455                 ],
3456                 coords = point.coords || point,
3457                 bbox = brd.getBoundingBox();
3458 
3459             for (i = 0; i < 4; i++) {
3460                 c = config[i];
3461                 if (c[0] * coords.usrCoords[c[1]] < c[0] * bbox[c[2]]) {
3462                     // define border
3463                     l = Mat.crossProduct(
3464                         [1, bbox[c[3]], bbox[c[4]]],
3465                         [1, bbox[c[5]], bbox[c[6]]]
3466                     );
3467                     l[3] = 0;
3468                     l = Mat.normalize(l);
3469 
3470                     // project point
3471                     coords = this.projectPointToLine({ coords: coords }, { stdform: l }, brd);
3472                 }
3473             }
3474 
3475             return coords;
3476         },
3477 
3478         /**
3479          * Given a the coordinates of a point, finds the nearest point on the given
3480          * parametric curve or surface, and returns its view-space coordinates.
3481          * @param {Array} pScr Screen coordinates to project.
3482          * @param {JXG.Curve3D|JXG.Surface3D} target Parametric curve or surface to project to.
3483          * @param {Array} params Parameters of point on the target, initially specifying the starting point of
3484          * the search. The parameters are modified in place during the search, ending up at the nearest point.
3485          * @returns {Array} Array of length 4 containing the coordinates of the nearest point on the curve or surface.
3486          */
3487         projectCoordsToParametric: function (p, target, params) {
3488             // The variables and parameters for the Cobyla constrained
3489             // minimization algorithm are explained in the Cobyla.js comments
3490             var rhobeg, // initial size of simplex (Cobyla)
3491                 rhoend, // finial size of simplex (Cobyla)
3492                 iprint = 0, // no console output (Cobyla)
3493                 maxfun = 200, // call objective function at most 200 times (Cobyla)
3494                 dim = params.length,
3495                 _minFunc; // objective function (Cobyla)
3496 
3497             // adapt simplex size to parameter range
3498             if (dim === 1) {
3499                 rhobeg = 0.1 * (target.range[1] - target.range[0]);
3500             } else if (dim === 2) {
3501                 rhobeg = 0.1 * Math.min(
3502                     target.range_u[1] - target.range_u[0],
3503                     target.range_v[1] - target.range_v[0]
3504                 );
3505             }
3506             rhoend = rhobeg / 5e6;
3507 
3508             // minimize screen distance to cursor
3509             _minFunc = function (n, m, w, con) {
3510                 // var xDiff = p[0] - target.X(...w),
3511                 //     yDiff = p[1] - target.Y(...w),
3512                 //     zDiff = p[2] - target.Z(...w);
3513                 var xDiff = p[0] - target.X.apply(null, w),
3514                     yDiff = p[1] - target.Y.apply(null, w),
3515                     zDiff = p[2] - target.Z.apply(null, w);
3516 
3517                 if (n === 1) {
3518                     con[0] = w[0] - target.range[0];
3519                     con[1] = -w[0] + target.range[1];
3520                 } else if (n === 2) {
3521                     con[0] = w[0] - target.range_u[0];
3522                     con[1] = -w[0] + target.range_u[1];
3523                     con[2] = w[1] - target.range_v[0];
3524                     con[3] = -w[1] + target.range_v[1];
3525                 }
3526                 return xDiff * xDiff + yDiff * yDiff + zDiff * zDiff;
3527             };
3528             Mat.Nlp.FindMinimum(_minFunc, dim, 2 * dim, params, rhobeg, rhoend, iprint, maxfun);
3529 
3530             // return [1, target.X(...params), target.Y(...params), target.Z(...params)];
3531             return [1, target.X.apply(null, params), target.Y.apply(null, params), target.Z.apply(null, params)];
3532         },
3533 
3534         /**
3535          * Given a the screen coordinates of a point, finds the point on the
3536          * given parametric curve or surface which is nearest in screen space,
3537          * and returns its view-space coordinates.
3538          * @param {Array} pScr Screen coordinates to project.
3539          * @param {JXG.Curve3D|JXG.Surface3D} target Parametric curve or surface to project to.
3540          * @param {Array} params Parameters of point on the target, initially specifying the starting point of
3541          * the search. The parameters are modified in place during the search, ending up at the nearest point.
3542          * @returns {Array} Array of length 4 containing the coordinates of the nearest point on the curve or surface.
3543          */
3544         projectScreenCoordsToParametric: function (pScr, target, params) {
3545             // The variables and parameters for the Cobyla constrained
3546             // minimization algorithm are explained in the Cobyla.js comments
3547             var rhobeg, // initial size of simplex (Cobyla)
3548                 rhoend, // finial size of simplex (Cobyla)
3549                 iprint = 0, // no console output (Cobyla)
3550                 maxfun = 200, // call objective function at most 200 times (Cobyla)
3551                 dim = params.length,
3552                 _minFunc; // objective function (Cobyla)
3553 
3554             // adapt simplex size to parameter range
3555             if (dim === 1) {
3556                 rhobeg = 0.1 * (target.range[1] - target.range[0]);
3557             } else if (dim === 2) {
3558                 rhobeg = 0.1 * Math.min(
3559                     target.range_u[1] - target.range_u[0],
3560                     target.range_v[1] - target.range_v[0]
3561                 );
3562             }
3563             rhoend = rhobeg / 5e6;
3564 
3565             // minimize screen distance to cursor
3566             _minFunc = function (n, m, w, con) {
3567                 // var c3d = [
3568                 //         1,
3569                 //         target.X(...w),
3570                 //         target.Y(...w),
3571                 //         target.Z(...w)
3572                 //     ],
3573                 var c3d = [
3574                     1,
3575                     target.X.apply(null, w),
3576                     target.Y.apply(null, w),
3577                     target.Z.apply(null, w)
3578                 ],
3579                     c2d = target.view.project3DTo2D(c3d),
3580                     xDiff = pScr[0] - c2d[1],
3581                     yDiff = pScr[1] - c2d[2];
3582 
3583                 if (n === 1) {
3584                     con[0] = w[0] - target.range[0];
3585                     con[1] = -w[0] + target.range[1];
3586                 } else if (n === 2) {
3587                     con[0] = w[0] - target.range_u[0];
3588                     con[1] = -w[0] + target.range_u[1];
3589                     con[2] = w[1] - target.range_v[0];
3590                     con[3] = -w[1] + target.range_v[1];
3591                 }
3592                 return xDiff * xDiff + yDiff * yDiff;
3593             };
3594             Mat.Nlp.FindMinimum(_minFunc, dim, 2 * dim, params, rhobeg, rhoend, iprint, maxfun);
3595 
3596             // return [1, target.X(...params), target.Y(...params), target.Z(...params)];
3597             return [1, target.X.apply(null, params), target.Y.apply(null, params), target.Z.apply(null, params)];
3598         },
3599 
3600         /**
3601          * Calculates the distance of a point to a line. The point and the line are given by homogeneous
3602          * coordinates. For lines this can be line.stdform.
3603          * @param {Array} point Homogeneous coordinates of a point.
3604          * @param {Array} line Homogeneous coordinates of a line ([C,A,B] where A*x+B*y+C*z=0).
3605          * @returns {Number} Distance of the point to the line.
3606          */
3607         distPointLine: function (point, line) {
3608             var a = line[1],
3609                 b = line[2],
3610                 c = line[0],
3611                 nom;
3612 
3613             if (Math.abs(a) + Math.abs(b) < Mat.eps) {
3614                 return Number.POSITIVE_INFINITY;
3615             }
3616 
3617             nom = a * point[1] + b * point[2] + c;
3618             a *= a;
3619             b *= b;
3620 
3621             return Math.abs(nom) / Math.sqrt(a + b);
3622         },
3623 
3624         /**
3625          * Determine the (Euclidean) distance between a point q and a line segment
3626          * defined by two points p1 and p2. In case p1 equals p2, the distance to this
3627          * point is returned.
3628          *
3629          * @param {Array} q Homogeneous coordinates of q
3630          * @param {Array} p1 Homogeneous coordinates of p1
3631          * @param {Array} p2 Homogeneous coordinates of p2
3632          * @returns {Number} Distance of q to line segment [p1, p2]
3633          */
3634         distPointSegment: function (q, p1, p2) {
3635             var x, y, dx, dy,
3636                 den, lbda,
3637                 eps = Mat.eps * Mat.eps,
3638                 huge = 1000000;
3639 
3640             // Difference q - p1
3641             x = q[1] - p1[1];
3642             y = q[2] - p1[2];
3643             x = (x === Infinity) ? huge : (x === -Infinity) ? -huge : x;
3644             y = (y === Infinity) ? huge : (y === -Infinity) ? -huge : y;
3645 
3646             // Difference p2 - p1
3647             dx = p2[1] - p1[1];
3648             dy = p2[2] - p1[2];
3649             dx = (dx === Infinity) ? huge : (dx === -Infinity) ? -huge : dx;
3650             dy = (dy === Infinity) ? huge : (dy === -Infinity) ? -huge : dy;
3651 
3652             // If den==0 then p1 and p2 are identical
3653             // In this case the distance to p1 is returned
3654             den = dx * dx + dy * dy;
3655             if (den > eps) {
3656                 lbda = (x * dx + y * dy) / den;
3657                 if (lbda < 0.0) {
3658                     lbda = 0.0;
3659                 } else if (lbda > 1.0) {
3660                     lbda = 1.0;
3661                 }
3662                 x -= lbda * dx;
3663                 y -= lbda * dy;
3664             }
3665 
3666             return Mat.hypot(x, y);
3667         },
3668 
3669         /**
3670          * Helper function to create curve which displays a Reuleaux polygons.
3671          * @param {Array} points Array of points which should be the vertices of the Reuleaux polygon. Typically,
3672          * these point list is the array vertices of a regular polygon.
3673          * @param {Number} nr Number of vertices
3674          * @returns {Array} An array containing the two functions defining the Reuleaux polygon and the two values
3675          * for the start and the end of the paramtric curve. array may be used as parent array of a
3676          * {@link JXG.Curve}.
3677          *
3678          * @example
3679          * var A = brd.create('point',[-2,-2]);
3680          * var B = brd.create('point',[0,1]);
3681          * var pol = brd.create('regularpolygon',[A,B,3], {withLines:false, fillColor:'none', highlightFillColor:'none', fillOpacity:0.0});
3682          * var reuleauxTriangle = brd.create('curve', JXG.Math.Geometry.reuleauxPolygon(pol.vertices, 3),
3683          *                          {strokeWidth:6, strokeColor:'#d66d55', fillColor:'#ad5544', highlightFillColor:'#ad5544'});
3684          *
3685          * </pre><div class="jxgbox" id="JXG2543a843-46a9-4372-abc1-94d9ad2db7ac" style="width: 300px; height: 300px;"></div>
3686          * <script type="text/javascript">
3687          * var brd = JXG.JSXGraph.initBoard('JXG2543a843-46a9-4372-abc1-94d9ad2db7ac', {boundingbox: [-5, 5, 5, -5], axis: true, showcopyright:false, shownavigation: false});
3688          * var A = brd.create('point',[-2,-2]);
3689          * var B = brd.create('point',[0,1]);
3690          * var pol = brd.create('regularpolygon',[A,B,3], {withLines:false, fillColor:'none', highlightFillColor:'none', fillOpacity:0.0});
3691          * var reuleauxTriangle = brd.create('curve', JXG.Math.Geometry.reuleauxPolygon(pol.vertices, 3),
3692          *                          {strokeWidth:6, strokeColor:'#d66d55', fillColor:'#ad5544', highlightFillColor:'#ad5544'});
3693          * </script><pre>
3694          */
3695         reuleauxPolygon: function (points, nr) {
3696             var beta,
3697                 pi2 = Math.PI * 2,
3698                 pi2_n = pi2 / nr,
3699                 diag = (nr - 1) / 2,
3700                 d = 0,
3701                 makeFct = function (which, trig) {
3702                     return function (t, suspendUpdate) {
3703                         var t1 = ((t % pi2) + pi2) % pi2,
3704                             j = Math.floor(t1 / pi2_n) % nr;
3705 
3706                         if (!suspendUpdate) {
3707                             d = points[0].Dist(points[diag]);
3708                             beta = Mat.Geometry.rad(
3709                                 [points[0].X() + 1, points[0].Y()],
3710                                 points[0],
3711                                 points[diag % nr]
3712                             );
3713                         }
3714 
3715                         if (isNaN(j)) {
3716                             return j;
3717                         }
3718 
3719                         t1 = t1 * 0.5 + j * pi2_n * 0.5 + beta;
3720 
3721                         return points[j][which]() + d * Math[trig](t1);
3722                     };
3723                 };
3724 
3725             return [makeFct("X", "cos"), makeFct("Y", "sin"), 0, pi2];
3726         },
3727 
3728         meet3Planes: function (n1, d1, n2, d2, n3, d3) {
3729             var p = [0, 0, 0],
3730                 n31,
3731                 n12,
3732                 n23,
3733                 denom,
3734                 i;
3735 
3736             n31 = Mat.crossProduct(n3, n1);
3737             n12 = Mat.crossProduct(n1, n2);
3738             n23 = Mat.crossProduct(n2, n3);
3739             denom = Mat.innerProduct(n1, n23, 3);
3740             for (i = 0; i < 3; i++) {
3741                 p[i] = (d1 * n23[i] + d2 * n31[i] + d3 * n12[i]) / denom;
3742             }
3743             return p;
3744         },
3745 
3746         meetPlanePlane: function (v11, v12, v21, v22) {
3747             var i,
3748                 no1,
3749                 no2,
3750                 v = [0, 0, 0],
3751                 w = [0, 0, 0];
3752 
3753             for (i = 0; i < 3; i++) {
3754                 v[i] = Type.evaluate(v11[i]);
3755                 w[i] = Type.evaluate(v12[i]);
3756             }
3757             no1 = Mat.crossProduct(v, w);
3758 
3759             for (i = 0; i < 3; i++) {
3760                 v[i] = Type.evaluate(v21[i]);
3761                 w[i] = Type.evaluate(v22[i]);
3762             }
3763             no2 = Mat.crossProduct(v, w);
3764 
3765             return Mat.crossProduct(no1, no2);
3766         },
3767 
3768         project3DTo3DPlane: function (point, normal, foot) {
3769             // TODO: homogeneous 3D coordinates
3770             var sol = [0, 0, 0],
3771                 le,
3772                 d1,
3773                 d2,
3774                 lbda;
3775 
3776             foot = foot || [0, 0, 0];
3777 
3778             le = Mat.norm(normal);
3779             d1 = Mat.innerProduct(point, normal, 3);
3780             d2 = Mat.innerProduct(foot, normal, 3);
3781             // (point - lbda * normal / le) * normal / le == foot * normal / le
3782             // => (point * normal - foot * normal) ==  lbda * le
3783             lbda = (d1 - d2) / le;
3784             sol = Mat.axpy(-lbda, normal, point);
3785 
3786             return sol;
3787         },
3788 
3789         getPlaneBounds: function (v1, v2, q, s, e) {
3790             var s1, s2, e1, e2, mat, rhs, sol;
3791 
3792             if (v1[2] + v2[0] !== 0) {
3793                 mat = [
3794                     [v1[0], v2[0]],
3795                     [v1[1], v2[1]]
3796                 ];
3797                 rhs = [s - q[0], s - q[1]];
3798 
3799                 sol = Numerics.Gauss(mat, rhs);
3800                 s1 = sol[0];
3801                 s2 = sol[1];
3802 
3803                 rhs = [e - q[0], e - q[1]];
3804                 sol = Numerics.Gauss(mat, rhs);
3805                 e1 = sol[0];
3806                 e2 = sol[1];
3807                 return [s1, e1, s2, e2];
3808             }
3809             return null;
3810         }
3811     }
3812 );
3813 
3814 export default Mat.Geometry;
3815