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