1 /* 2 Copyright 2008-2026 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("JSXGraph example:", 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("JSXGraph example:", 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 // [mi_x_i, ma_x_i, mi_y_i, ma_y_i, mi_xpy_i, mi_xmy_i, ma_xpy_i, ma_xmy_i].includes(i) || 819 [mi_x_i, ma_x_i, mi_y_i, ma_y_i, mi_xpy_i, mi_xmy_i, ma_xpy_i, ma_xmy_i].indexOf(i) >= 0|| 820 (mi_x_i !== mi_xmy_i && this.signedTriangle(ps[mi_x_i], ps[mi_xmy_i], c) >= -eps) || 821 (mi_xmy_i !== ma_y_i && this.signedTriangle(ps[mi_xmy_i], ps[ma_y_i], c) >= -eps) || 822 (ma_y_i !== ma_xpy_i && this.signedTriangle(ps[ma_y_i], ps[ma_xpy_i], c) >= -eps) || 823 (ma_xpy_i !== ma_x_i && this.signedTriangle(ps[ma_xpy_i], ps[ma_x_i], c) >= -eps) || 824 (ma_x_i !== ma_xmy_i && this.signedTriangle(ps[ma_x_i], ps[ma_xmy_i], c) >= -eps) || 825 (ma_xmy_i !== mi_y_i && this.signedTriangle(ps[ma_xmy_i], ps[mi_y_i], c) >= -eps) || 826 (mi_y_i !== mi_xpy_i && this.signedTriangle(ps[mi_y_i], ps[mi_xpy_i], c) >= -eps) || 827 (mi_xpy_i !== mi_x_i && this.signedTriangle(ps[mi_xpy_i], ps[mi_x_i], c) >= -eps) 828 ) { 829 ps_idx.push({ 830 i: i, 831 c: c 832 }); 833 } 834 } 835 N = ps_idx.length; 836 837 // Find the point with the lowest y value 838 mi_idx = 0; 839 mi_x = ps_idx[0].c[1]; 840 mi_y = ps_idx[0].c[2]; 841 for (i = 1; i < N; i++) { 842 if ((ps_idx[i].c[2] < mi_y) || (ps_idx[i].c[2] === mi_y && ps_idx[i].c[1] < mi_x)) { 843 mi_x = ps_idx[i].c[1]; 844 mi_y = ps_idx[i].c[2]; 845 mi_idx = i; 846 } 847 } 848 ps_idx = Type.swap(ps_idx, mi_idx, 0); 849 850 // Our origin o, i.e. the first point. 851 o = ps_idx[0].c; 852 853 // Sort according to the angle around o. 854 ps_idx.sort(function(a_obj, b_obj) { 855 var a = a_obj.c, 856 b = b_obj.c, 857 v = that.signedTriangle(o, a, b); 858 859 if (v === 0) { 860 // if o, a, b are collinear, the point which is further away 861 // from o is considered greater. 862 return Mat.hypot(a[1] - o[1], a[2] - o[2]) - Mat.hypot(b[1] - o[1], b[2] - o[2]); 863 } 864 865 // if v < 0, a is to the left of [o, b], i.e. angle(a) > angle(b) 866 return -v; 867 }); 868 869 // Do the Graham scan. 870 M = 0; 871 for (i = 0; i < N; i++) { 872 while (M > 1 && this.signedTriangle(stack[M - 2].c, stack[M - 1].c, ps_idx[i].c) <= 0) { 873 // stack[M - 1] is to the left of stack[M - 1], ps[i]: discard it 874 stack.pop(); 875 M--; 876 } 877 stack.push(ps_idx[i]); 878 M++; 879 } 880 881 return stack; 882 }, 883 884 // Original method 885 // GrahamScan: function (points, indices) { 886 // var i, 887 // M = 1, 888 // ps = Expect.each(points, Expect.coordsArray), 889 // N = ps.length; 890 // ps = this.sortVertices(ps); 891 // N = ps.length; 892 // for (i = 2; i < N; i++) { 893 // while (this.signedTriangle(ps[M - 1], ps[M], ps[i]) <= 0) { 894 // if (M > 1) { 895 // M -= 1; 896 // } else if (i === N - 1) { 897 // break; 898 // } 899 // i += 1; 900 // } 901 // M += 1; 902 // ps = Type.swap(ps, M, i); 903 // indices = Type.swap(indices, M, i); 904 // } 905 // return ps.slice(0, M); 906 // }, 907 908 /** 909 * Calculate the complex hull of a point cloud by the Graham scan algorithm. 910 * 911 * @param {Array} points An array containing {@link JXG.Point}, {@link JXG.Coords}, and/or arrays. 912 * @param {Boolean} [returnCoords=false] If true, return an array of coords. Otherwise return a list of pointers 913 * to the input list elements. That is, if the input is a list of {@link JXG.Point} elements, the returned list 914 * will contain the points that form the convex hull. 915 * @returns {Array} List containing the convex hull. Format depends on returnCoords. 916 * @see JXG.Math.Geometry.GrahamScan 917 * 918 * @example 919 * // Static example 920 * var i, hull, 921 * p = []; 922 * 923 * p.push( board.create('point', [4, 0], {withLabel:false }) ); 924 * p.push( board.create('point', [0, 4], {withLabel:false }) ); 925 * p.push( board.create('point', [0, 0], {withLabel:false }) ); 926 * p.push( board.create('point', [1, 1], {withLabel:false }) ); 927 * hull = JXG.Math.Geometry.convexHull(p); 928 * for (i = 0; i < hull.length; i++) { 929 * hull[i].setAttribute({color: 'blue'}); 930 * } 931 * 932 * </pre><div id="JXGdfc76123-81b8-4250-96f9-419253bd95dd" class="jxgbox" style="width: 300px; height: 300px;"></div> 933 * <script type="text/javascript"> 934 * (function() { 935 * var board = JXG.JSXGraph.initBoard('JXGdfc76123-81b8-4250-96f9-419253bd95dd', 936 * {boundingbox: [-8, 8, 8,-8], axis: true, showcopyright: false, shownavigation: false}); 937 * var i, hull, 938 * p = []; 939 * 940 * p.push( board.create('point', [4, 0], {withLabel:false }) ); 941 * p.push( board.create('point', [0, 4], {withLabel:false }) ); 942 * p.push( board.create('point', [0, 0], {withLabel:false }) ); 943 * p.push( board.create('point', [1, 1], {withLabel:false }) ); 944 * hull = JXG.Math.Geometry.convexHull(p); 945 * for (i = 0; i < hull.length; i++) { 946 * hull[i].setAttribute({color: 'blue'}); 947 * } 948 * 949 * })(); 950 * 951 * </script><pre> 952 * 953 * @example 954 * // Dynamic version using returnCoords==true: drag the points 955 * var p = []; 956 * 957 * p.push( board.create('point', [4, 0], {withLabel:false }) ); 958 * p.push( board.create('point', [0, 4], {withLabel:false }) ); 959 * p.push( board.create('point', [0, 0], {withLabel:false }) ); 960 * p.push( board.create('point', [1, 1], {withLabel:false }) ); 961 * 962 * var c = board.create('curve', [[], []], {fillColor: 'yellow', fillOpacity: 0.3}); 963 * c.updateDataArray = function() { 964 * var i, 965 * hull = JXG.Math.Geometry.convexHull(p, true); 966 * 967 * this.dataX = []; 968 * this.dataY = []; 969 * 970 * for (i = 0; i < hull.length; i ++) { 971 * this.dataX.push(hull[i][1]); 972 * this.dataY.push(hull[i][2]); 973 * } 974 * this.dataX.push(hull[0][1]); 975 * this.dataY.push(hull[0][2]); 976 * }; 977 * board.update(); 978 * 979 * </pre><div id="JXG61e51909-da0b-432f-9aa7-9fb0c8bb01c9" class="jxgbox" style="width: 300px; height: 300px;"></div> 980 * <script type="text/javascript"> 981 * (function() { 982 * var board = JXG.JSXGraph.initBoard('JXG61e51909-da0b-432f-9aa7-9fb0c8bb01c9', 983 * {boundingbox: [-8, 8, 8,-8], axis: true, showcopyright: false, shownavigation: false}); 984 * var p = []; 985 * 986 * p.push( board.create('point', [4, 0], {withLabel:false }) ); 987 * p.push( board.create('point', [0, 4], {withLabel:false }) ); 988 * p.push( board.create('point', [0, 0], {withLabel:false }) ); 989 * p.push( board.create('point', [1, 1], {withLabel:false }) ); 990 * 991 * var c = board.create('curve', [[], []], {fillColor: 'yellow', fillOpacity: 0.3}); 992 * c.updateDataArray = function() { 993 * var i, 994 * hull = JXG.Math.Geometry.convexHull(p, true); 995 * 996 * this.dataX = []; 997 * this.dataY = []; 998 * 999 * for (i = 0; i < hull.length; i ++) { 1000 * this.dataX.push(hull[i][1]); 1001 * this.dataY.push(hull[i][2]); 1002 * } 1003 * this.dataX.push(hull[0][1]); 1004 * this.dataY.push(hull[0][2]); 1005 * }; 1006 * board.update(); 1007 * 1008 * 1009 * })(); 1010 * 1011 * </script><pre> 1012 * 1013 */ 1014 convexHull: function(points, returnCoords) { 1015 var i, hull, 1016 res = []; 1017 1018 hull = this.GrahamScan(points); 1019 for (i = 0; i < hull.length; i++) { 1020 if (returnCoords) { 1021 res.push(hull[i].c); 1022 } else { 1023 res.push(points[hull[i].i]); 1024 } 1025 } 1026 return res; 1027 }, 1028 1029 // /** 1030 // * Determine if a polygon or a path element is convex, non-convex or complex which are defined like this: 1031 // * <ul> 1032 // * <li> A polygon is convex if for every pair of points, the line segment connecting them does not intersect 1033 // * an edge of the polygon in one point. 1034 // * A single line segment or a a single point is considered as convex. A necessary condition for a polygon 1035 // * to be convex that the angle sum of its interior angles equals ± 2 π. 1036 // * <li> A polygon is non-convex, if it does not self-intersect, but is not convex. 1037 // * <li> A polygon is complex if its the angle sum is not equal to ± 2 π. 1038 // * That is, there must be self-intersection (contiguous coincident points in the path are not treated as self-intersection). 1039 // * </ul> 1040 // * A path element might be specified as an array of coordinate arrays or {@link JXG.Coords}. 1041 // * 1042 // * @param {Array|Polygon|PolygonalChain} points Polygon or list of coordinates 1043 // * @returns {Number} -1: if complex, 0: if non-convex, 1: if convex 1044 // */ 1045 /** 1046 * Determine if a polygon or a path element is convex: 1047 * <p> 1048 * A polygon is convex if for every pair of points, the line segment connecting them does not intersect 1049 * an edge of the polygon in one point. 1050 * A single line segment, a single point, or the empty set is considered as convex. A necessary condition for a polygon 1051 * to be convex that the angle sum of its interior angles equals ± 2 π. 1052 * <p> 1053 * A path element might be specified as an array of coordinate arrays or {@link JXG.Coords}. 1054 * See the discussion at <a href="https://stackoverflow.com/questions/471962/how-do-i-efficiently-determine-if-a-polygon-is-convex-non-convex-or-complex">stackoverflow</a>. 1055 * 1056 * @param {Array|Polygon|PolygonalChain} points Polygon or list of coordinates 1057 * @returns {Boolean} true if convex 1058 * 1059 * @example 1060 * var pol = board.create('polygon', [ 1061 * [-1, -1], 1062 * [3, -1], 1063 * [4, 2], 1064 * [3, 3], 1065 * [0, 4], 1066 * [-3, 1] 1067 * ], { 1068 * vertices: { 1069 * color: 'blue', 1070 * snapToGrid: true 1071 * } 1072 * }); 1073 * 1074 * console.log("JSXGraph example:", JXG.Math.Geometry.isConvex(pol)); 1075 * // > true 1076 * 1077 * 1078 * 1079 * </pre><div id="JXG9b43cc53-15b4-49be-92cc-2a1dfc06665b" class="jxgbox" style="width: 300px; height: 300px;"></div> 1080 * <script type="text/javascript"> 1081 * (function() { 1082 * var board = JXG.JSXGraph.initBoard('JXG9b43cc53-15b4-49be-92cc-2a1dfc06665b', 1083 * {boundingbox: [-8, 8, 8,-8], axis: true, showcopyright: false, shownavigation: false}); 1084 * var pol = board.create('polygon', [ 1085 * [-1, -1], 1086 * [3, -1], 1087 * [4, 2], 1088 * [3, 3], 1089 * [0, 4], 1090 * [-3, 1] 1091 * ], { 1092 * vertices: { 1093 * color: 'blue', 1094 * snapToGrid: true 1095 * } 1096 * }); 1097 * 1098 * console.log("JSXGraph example:", JXG.Math.Geometry.isConvex(pol)); 1099 * })(); 1100 * 1101 * </script><pre> 1102 * 1103 */ 1104 isConvex: function(points) { 1105 var ps, le, i, 1106 eps = Mat.eps * Mat.eps, 1107 old_x, old_y, old_dir, 1108 new_x, new_y, new_dir, 1109 angle, 1110 orient, 1111 angle_sum = 0.0; 1112 1113 if (Type.isArray(points)) { 1114 ps = Expect.each(points, Expect.coordsArray); 1115 } else if (Type.exists(points.type) && points.type === Const.OBJECT_TYPE_POLYGON) { 1116 ps = Expect.each(points.vertices, Expect.coordsArray); 1117 } 1118 le = ps.length; 1119 if (le === 0) { 1120 // Empty set is convex 1121 return true; 1122 } 1123 if (le < 3) { 1124 // Segments and points are convex 1125 return true; 1126 } 1127 1128 orient = null; 1129 old_x = ps[le - 2][1]; 1130 old_y = ps[le - 2][2]; 1131 new_x = ps[le - 1][1]; 1132 new_y = ps[le - 1][2]; 1133 new_dir = Math.atan2(new_y - old_y, new_x - old_x); 1134 for (i = 0; i < le; i++) { 1135 old_x = new_x; 1136 old_y = new_y; 1137 old_dir = new_dir; 1138 new_x = ps[i][1]; 1139 new_y = ps[i][2]; 1140 if (old_x === new_x && old_y === new_y) { 1141 // Repeated consecutive points are ignored 1142 continue; 1143 } 1144 new_dir = Math.atan2(new_y - old_y, new_x - old_x); 1145 angle = new_dir - old_dir; 1146 if (angle <= -Math.PI) { 1147 angle += 2 * Math.PI; 1148 } else if (angle > Math.PI) { 1149 angle -= 2 * Math.PI; 1150 } 1151 if (orient === null) { 1152 if (angle === 0.0) { 1153 continue; 1154 } 1155 orient = (angle > 0) ? 1 : -1; 1156 } else { 1157 if (orient * angle < -eps) { 1158 return false; 1159 } 1160 } 1161 angle_sum += angle; 1162 } 1163 1164 if ((Math.abs(angle_sum / (2 * Math.PI)) - 1) < eps) { 1165 return true; 1166 } 1167 return false; 1168 }, 1169 1170 /** 1171 * A line can be a segment, a straight, or a ray. So it is not always delimited by point1 and point2 1172 * calcStraight determines the visual start point and end point of the line. A segment is only drawn 1173 * from start to end point, a straight line is drawn until it meets the boards boundaries. 1174 * @param {JXG.Line} el Reference to a line object, that needs calculation of start and end point. 1175 * @param {JXG.Coords} point1 Coordinates of the point where line drawing begins. This value is calculated and 1176 * set by this method. 1177 * @param {JXG.Coords} point2 Coordinates of the point where line drawing ends. This value is calculated and set 1178 * by this method. 1179 * @param {Number} margin Optional margin, to avoid the display of the small sides of lines. 1180 * @returns null 1181 * @see Line 1182 * @see JXG.Line 1183 */ 1184 calcStraight: function (el, point1, point2, margin) { 1185 var takePoint1, 1186 takePoint2, 1187 intersection, 1188 intersect1, 1189 intersect2, 1190 straightFirst, 1191 straightLast, 1192 c, p1, p2; 1193 1194 if (!Type.exists(margin)) { 1195 // Enlarge the drawable region slightly. This hides the small sides 1196 // of thick lines in most cases. 1197 margin = 10; 1198 } 1199 1200 straightFirst = el.evalVisProp('straightfirst'); 1201 straightLast = el.evalVisProp('straightlast'); 1202 1203 // If one of the point is an ideal point in homogeneous coordinates 1204 // drawing of line segments or rays are not possible. 1205 if (Math.abs(point1.scrCoords[0]) < Mat.eps) { 1206 straightFirst = true; 1207 } 1208 if (Math.abs(point2.scrCoords[0]) < Mat.eps) { 1209 straightLast = true; 1210 } 1211 1212 // Do nothing in case of line segments (inside or outside of the board) 1213 if (!straightFirst && !straightLast) { 1214 return; 1215 } 1216 1217 // Compute the stdform of the line in screen coordinates. 1218 c = []; 1219 c[0] = 1220 el.stdform[0] - 1221 (el.stdform[1] * el.board.origin.scrCoords[1]) / el.board.unitX + 1222 (el.stdform[2] * el.board.origin.scrCoords[2]) / el.board.unitY; 1223 c[1] = el.stdform[1] / el.board.unitX; 1224 c[2] = -el.stdform[2] / el.board.unitY; 1225 1226 // If p1=p2 1227 if (isNaN(c[0] + c[1] + c[2])) { 1228 return; 1229 } 1230 1231 takePoint1 = false; 1232 takePoint2 = false; 1233 1234 // Line starts at point1 and point1 is inside the board 1235 takePoint1 = 1236 !straightFirst && 1237 Math.abs(point1.usrCoords[0]) >= Mat.eps && 1238 point1.scrCoords[1] >= 0.0 && 1239 point1.scrCoords[1] <= el.board.canvasWidth && 1240 point1.scrCoords[2] >= 0.0 && 1241 point1.scrCoords[2] <= el.board.canvasHeight; 1242 1243 // Line ends at point2 and point2 is inside the board 1244 takePoint2 = 1245 !straightLast && 1246 Math.abs(point2.usrCoords[0]) >= Mat.eps && 1247 point2.scrCoords[1] >= 0.0 && 1248 point2.scrCoords[1] <= el.board.canvasWidth && 1249 point2.scrCoords[2] >= 0.0 && 1250 point2.scrCoords[2] <= el.board.canvasHeight; 1251 1252 // Intersect the line with the four borders of the board. 1253 intersection = this.meetLineBoard(c, el.board, margin); 1254 intersect1 = intersection[0]; 1255 intersect2 = intersection[1]; 1256 1257 /** 1258 * At this point we have four points: 1259 * point1 and point2 are the first and the second defining point on the line, 1260 * intersect1, intersect2 are the intersections of the line with border around the board. 1261 */ 1262 1263 /* 1264 * Here we handle rays where both defining points are outside of the board. 1265 */ 1266 // If both points are outside and the complete ray is outside we do nothing 1267 if (!takePoint1 && !takePoint2) { 1268 // Ray starting at point 1 1269 if ( 1270 !straightFirst && 1271 straightLast && 1272 !this.isSameDirection(point1, point2, intersect1) && 1273 !this.isSameDirection(point1, point2, intersect2) 1274 ) { 1275 return; 1276 } 1277 1278 // Ray starting at point 2 1279 if ( 1280 straightFirst && 1281 !straightLast && 1282 !this.isSameDirection(point2, point1, intersect1) && 1283 !this.isSameDirection(point2, point1, intersect2) 1284 ) { 1285 return; 1286 } 1287 } 1288 1289 /* 1290 * If at least one of the defining points is outside of the board 1291 * we take intersect1 or intersect2 as one of the end points 1292 * The order is also important for arrows of axes 1293 */ 1294 if (!takePoint1) { 1295 if (!takePoint2) { 1296 // Two border intersection points are used 1297 if (this.isSameDir(point1, point2, intersect1, intersect2)) { 1298 p1 = intersect1; 1299 p2 = intersect2; 1300 } else { 1301 p2 = intersect1; 1302 p1 = intersect2; 1303 } 1304 } else { 1305 // One border intersection points is used 1306 if (this.isSameDir(point1, point2, intersect1, intersect2)) { 1307 p1 = intersect1; 1308 } else { 1309 p1 = intersect2; 1310 } 1311 } 1312 } else { 1313 if (!takePoint2) { 1314 // One border intersection points is used 1315 if (this.isSameDir(point1, point2, intersect1, intersect2)) { 1316 p2 = intersect2; 1317 } else { 1318 p2 = intersect1; 1319 } 1320 } 1321 } 1322 1323 if (p1) { 1324 //point1.setCoordinates(Const.COORDS_BY_USER, p1.usrCoords.slice(1)); 1325 point1.setCoordinates(Const.COORDS_BY_USER, p1.usrCoords); 1326 } 1327 1328 if (p2) { 1329 //point2.setCoordinates(Const.COORDS_BY_USER, p2.usrCoords.slice(1)); 1330 point2.setCoordinates(Const.COORDS_BY_USER, p2.usrCoords); 1331 } 1332 }, 1333 1334 /** 1335 * A line can be a segment, a straight, or a ray. so it is not always delimited by point1 and point2. 1336 * 1337 * This method adjusts the line's delimiting points taking into account its nature, the viewport defined 1338 * by the board. 1339 * 1340 * A segment is delimited by start and end point, a straight line or ray is delimited until it meets the 1341 * boards boundaries. However, if the line has infinite ticks, it will be delimited by the projection of 1342 * the boards vertices onto itself. 1343 * 1344 * @param {JXG.Line} el Reference to a line object, that needs calculation of start and end point. 1345 * @param {JXG.Coords} point1 Coordinates of the point where line drawing begins. This value is calculated and 1346 * set by this method. 1347 * @param {JXG.Coords} point2 Coordinates of the point where line drawing ends. This value is calculated and set 1348 * by this method. 1349 * @see Line 1350 * @see JXG.Line 1351 */ 1352 calcLineDelimitingPoints: function (el, point1, point2) { 1353 var distP1P2, 1354 boundingBox, 1355 lineSlope, 1356 intersect1, 1357 intersect2, 1358 straightFirst, 1359 straightLast, 1360 c, 1361 p1, 1362 p2, 1363 takePoint1 = false, 1364 takePoint2 = false; 1365 1366 straightFirst = el.evalVisProp('straightfirst'); 1367 straightLast = el.evalVisProp('straightlast'); 1368 1369 // If one of the point is an ideal point in homogeneous coordinates 1370 // drawing of line segments or rays are not possible. 1371 if (Math.abs(point1.scrCoords[0]) < Mat.eps) { 1372 straightFirst = true; 1373 } 1374 if (Math.abs(point2.scrCoords[0]) < Mat.eps) { 1375 straightLast = true; 1376 } 1377 1378 // Compute the stdform of the line in screen coordinates. 1379 c = []; 1380 c[0] = 1381 el.stdform[0] - 1382 (el.stdform[1] * el.board.origin.scrCoords[1]) / el.board.unitX + 1383 (el.stdform[2] * el.board.origin.scrCoords[2]) / el.board.unitY; 1384 c[1] = el.stdform[1] / el.board.unitX; 1385 c[2] = -el.stdform[2] / el.board.unitY; 1386 1387 // p1=p2 1388 if (isNaN(c[0] + c[1] + c[2])) { 1389 return; 1390 } 1391 1392 takePoint1 = !straightFirst; 1393 takePoint2 = !straightLast; 1394 // Intersect the board vertices on the line to establish the available visual space for the infinite ticks 1395 // Based on the slope of the line we can optimise and only project the two outer vertices 1396 1397 // boundingBox = [x1, y1, x2, y2] upper left, lower right vertices 1398 boundingBox = el.board.getBoundingBox(); 1399 lineSlope = el.getSlope(); 1400 if (lineSlope >= 0) { 1401 // project vertices (x2,y1) (x1, y2) 1402 intersect1 = this.projectPointToLine( 1403 { coords: { usrCoords: [1, boundingBox[2], boundingBox[1]] } }, 1404 el, 1405 el.board 1406 ); 1407 intersect2 = this.projectPointToLine( 1408 { coords: { usrCoords: [1, boundingBox[0], boundingBox[3]] } }, 1409 el, 1410 el.board 1411 ); 1412 } else { 1413 // project vertices (x1, y1) (x2, y2) 1414 intersect1 = this.projectPointToLine( 1415 { coords: { usrCoords: [1, boundingBox[0], boundingBox[1]] } }, 1416 el, 1417 el.board 1418 ); 1419 intersect2 = this.projectPointToLine( 1420 { coords: { usrCoords: [1, boundingBox[2], boundingBox[3]] } }, 1421 el, 1422 el.board 1423 ); 1424 } 1425 1426 /** 1427 * we have four points: 1428 * point1 and point2 are the first and the second defining point on the line, 1429 * intersect1, intersect2 are the intersections of the line with border around the board. 1430 */ 1431 1432 /* 1433 * Here we handle rays/segments where both defining points are outside of the board. 1434 */ 1435 if (!takePoint1 && !takePoint2) { 1436 // Segment, if segment does not cross the board, do nothing 1437 if (!straightFirst && !straightLast) { 1438 distP1P2 = point1.distance(Const.COORDS_BY_USER, point2); 1439 // if intersect1 not between point1 and point2 1440 if ( 1441 Math.abs( 1442 point1.distance(Const.COORDS_BY_USER, intersect1) + 1443 intersect1.distance(Const.COORDS_BY_USER, point2) - 1444 distP1P2 1445 ) > Mat.eps 1446 ) { 1447 return; 1448 } 1449 // if insersect2 not between point1 and point2 1450 if ( 1451 Math.abs( 1452 point1.distance(Const.COORDS_BY_USER, intersect2) + 1453 intersect2.distance(Const.COORDS_BY_USER, point2) - 1454 distP1P2 1455 ) > Mat.eps 1456 ) { 1457 return; 1458 } 1459 } 1460 1461 // If both points are outside and the complete ray is outside we do nothing 1462 // Ray starting at point 1 1463 if ( 1464 !straightFirst && 1465 straightLast && 1466 !this.isSameDirection(point1, point2, intersect1) && 1467 !this.isSameDirection(point1, point2, intersect2) 1468 ) { 1469 return; 1470 } 1471 1472 // Ray starting at point 2 1473 if ( 1474 straightFirst && 1475 !straightLast && 1476 !this.isSameDirection(point2, point1, intersect1) && 1477 !this.isSameDirection(point2, point1, intersect2) 1478 ) { 1479 return; 1480 } 1481 } 1482 1483 /* 1484 * If at least one of the defining points is outside of the board 1485 * we take intersect1 or intersect2 as one of the end points 1486 * The order is also important for arrows of axes 1487 */ 1488 if (!takePoint1) { 1489 if (!takePoint2) { 1490 // Two border intersection points are used 1491 if (this.isSameDir(point1, point2, intersect1, intersect2)) { 1492 p1 = intersect1; 1493 p2 = intersect2; 1494 } else { 1495 p2 = intersect1; 1496 p1 = intersect2; 1497 } 1498 } else { 1499 // One border intersection points is used 1500 if (this.isSameDir(point1, point2, intersect1, intersect2)) { 1501 p1 = intersect1; 1502 } else { 1503 p1 = intersect2; 1504 } 1505 } 1506 } else { 1507 if (!takePoint2) { 1508 // One border intersection points is used 1509 if (this.isSameDir(point1, point2, intersect1, intersect2)) { 1510 p2 = intersect2; 1511 } else { 1512 p2 = intersect1; 1513 } 1514 } 1515 } 1516 1517 if (p1) { 1518 //point1.setCoordinates(Const.COORDS_BY_USER, p1.usrCoords.slice(1)); 1519 point1.setCoordinates(Const.COORDS_BY_USER, p1.usrCoords); 1520 } 1521 1522 if (p2) { 1523 //point2.setCoordinates(Const.COORDS_BY_USER, p2.usrCoords.slice(1)); 1524 point2.setCoordinates(Const.COORDS_BY_USER, p2.usrCoords); 1525 } 1526 }, 1527 1528 /** 1529 * Calculates the visProp.position corresponding to a given angle. 1530 * @param {number} angle angle in radians. Must be in range (-2pi,2pi). 1531 */ 1532 calcLabelQuadrant: function (angle) { 1533 var q; 1534 if (angle < 0) { 1535 angle += 2 * Math.PI; 1536 } 1537 q = Math.floor((angle + Math.PI / 8) / (Math.PI / 4)) % 8; 1538 return ["rt", "urt", "top", "ulft", "lft", "llft", "bot", "lrt"][q]; 1539 }, 1540 1541 /** 1542 * The vectors <tt>p2-p1</tt> and <tt>i2-i1</tt> are supposed to be collinear. If their cosine is positive 1543 * they point into the same direction otherwise they point in opposite direction. 1544 * @param {JXG.Coords} p1 1545 * @param {JXG.Coords} p2 1546 * @param {JXG.Coords} i1 1547 * @param {JXG.Coords} i2 1548 * @returns {Boolean} True, if <tt>p2-p1</tt> and <tt>i2-i1</tt> point into the same direction 1549 */ 1550 isSameDir: function (p1, p2, i1, i2) { 1551 var dpx = p2.usrCoords[1] - p1.usrCoords[1], 1552 dpy = p2.usrCoords[2] - p1.usrCoords[2], 1553 dix = i2.usrCoords[1] - i1.usrCoords[1], 1554 diy = i2.usrCoords[2] - i1.usrCoords[2]; 1555 1556 if (Math.abs(p2.usrCoords[0]) < Mat.eps) { 1557 dpx = p2.usrCoords[1]; 1558 dpy = p2.usrCoords[2]; 1559 } 1560 1561 if (Math.abs(p1.usrCoords[0]) < Mat.eps) { 1562 dpx = -p1.usrCoords[1]; 1563 dpy = -p1.usrCoords[2]; 1564 } 1565 1566 return dpx * dix + dpy * diy >= 0; 1567 }, 1568 1569 /** 1570 * If you're looking from point "start" towards point "s" and you can see the point "p", return true. 1571 * Otherwise return false. 1572 * @param {JXG.Coords} start The point you're standing on. 1573 * @param {JXG.Coords} p The point in which direction you're looking. 1574 * @param {JXG.Coords} s The point that should be visible. 1575 * @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. 1576 */ 1577 isSameDirection: function (start, p, s) { 1578 var dx, 1579 dy, 1580 sx, 1581 sy, 1582 r = false; 1583 1584 dx = p.usrCoords[1] - start.usrCoords[1]; 1585 dy = p.usrCoords[2] - start.usrCoords[2]; 1586 1587 sx = s.usrCoords[1] - start.usrCoords[1]; 1588 sy = s.usrCoords[2] - start.usrCoords[2]; 1589 1590 if (Math.abs(dx) < Mat.eps) { 1591 dx = 0; 1592 } 1593 1594 if (Math.abs(dy) < Mat.eps) { 1595 dy = 0; 1596 } 1597 1598 if (Math.abs(sx) < Mat.eps) { 1599 sx = 0; 1600 } 1601 1602 if (Math.abs(sy) < Mat.eps) { 1603 sy = 0; 1604 } 1605 1606 if (dx >= 0 && sx >= 0) { 1607 r = (dy >= 0 && sy >= 0) || (dy <= 0 && sy <= 0); 1608 } else if (dx <= 0 && sx <= 0) { 1609 r = (dy >= 0 && sy >= 0) || (dy <= 0 && sy <= 0); 1610 } 1611 1612 return r; 1613 }, 1614 1615 /** 1616 * Determinant of three points in the Euclidean plane. 1617 * Zero, if the points are collinear. Used to determine of a point q is left or 1618 * right to a segment defined by points p1 and p2. 1619 * <p> 1620 * Non-homogeneous version. 1621 * 1622 * @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. 1623 * @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. 1624 * @param {Array|JXG.Point} q Point or its coordinates. Point object or array of length 3. First (homogeneous) coordinate is equal to 1. 1625 * @return {Number} Signed area of the triangle formed by these three points. 1626 * 1627 * @see JXG.Math.Geometry.windingNumber 1628 */ 1629 det3p: function (p1, p2, q) { 1630 var pp1, pp2, qq; 1631 1632 if (Type.isPoint(p1)) { 1633 pp1 = p1.Coords(true); 1634 } else { 1635 pp1 = p1; 1636 } 1637 if (Type.isPoint(p2)) { 1638 pp2 = p2.Coords(true); 1639 } else { 1640 pp2 = p2; 1641 } 1642 if (Type.isPoint(q)) { 1643 qq = q.Coords(true); 1644 } else { 1645 qq = q; 1646 } 1647 1648 return (pp1[1] - qq[1]) * (pp2[2] - qq[2]) - (pp2[1] - qq[1]) * (pp1[2] - qq[2]); 1649 }, 1650 1651 /** 1652 * Winding number of a point in respect to a polygon path. 1653 * 1654 * The point is regarded outside if the winding number is zero, 1655 * inside otherwise. The algorithm tries to find degenerate cases, i.e. 1656 * if the point is on the path. This is regarded as "outside". 1657 * If the point is a vertex of the path, it is regarded as "inside". 1658 * 1659 * Implementation of algorithm 7 from "The point in polygon problem for 1660 * arbitrary polygons" by Kai Hormann and Alexander Agathos, Computational Geometry, 1661 * Volume 20, Issue 3, November 2001, Pages 131-144. 1662 * 1663 * @param {Array} usrCoords Homogenous coordinates of the point 1664 * @param {Array} path Array of points / coords determining a path, i.e. the vertices of the polygon / path. The array elements 1665 * do not have to be full points, but have to have a subobject "coords" or should be of type JXG.Coords. 1666 * @param {Boolean} [doNotClosePath=false] If true the last point of the path is not connected to the first point. 1667 * This is necessary if the path consists of two or more closed subpaths, e.g. if the figure has a hole. 1668 * 1669 * @return {Number} Winding number of the point. The point is 1670 * regarded outside if the winding number is zero, 1671 * inside otherwise. 1672 */ 1673 windingNumber: function (usrCoords, path, doNotClosePath) { 1674 var wn = 0, 1675 le = path.length, 1676 x = usrCoords[1], 1677 y = usrCoords[2], 1678 p0, 1679 p1, 1680 p2, 1681 d, 1682 sign, 1683 i, 1684 off = 0; 1685 1686 if (le === 0) { 1687 return 0; 1688 } 1689 1690 doNotClosePath = doNotClosePath || false; 1691 if (doNotClosePath) { 1692 off = 1; 1693 } 1694 1695 // Infinite points are declared outside 1696 if (isNaN(x) || isNaN(y)) { 1697 return 1; 1698 } 1699 1700 if (Type.exists(path[0].coords)) { 1701 p0 = path[0].coords; 1702 p1 = path[le - 1].coords; 1703 } else { 1704 p0 = path[0]; 1705 p1 = path[le - 1]; 1706 } 1707 // Handle the case if the point is the first vertex of the path, i.e. inside. 1708 if (p0.usrCoords[1] === x && p0.usrCoords[2] === y) { 1709 return 1; 1710 } 1711 1712 for (i = 0; i < le - off; i++) { 1713 // Consider the edge from p1 = path[i] to p2 = path[i+1]isClosedPath 1714 if (Type.exists(path[i].coords)) { 1715 p1 = path[i].coords.usrCoords; 1716 p2 = path[(i + 1) % le].coords.usrCoords; 1717 } else { 1718 p1 = path[i].usrCoords; 1719 p2 = path[(i + 1) % le].usrCoords; 1720 } 1721 1722 // If one of the two points p1, p2 is undefined or infinite, 1723 // move on. 1724 if ( 1725 p1[0] === 0 || 1726 p2[0] === 0 || 1727 isNaN(p1[1]) || 1728 isNaN(p2[1]) || 1729 isNaN(p1[2]) || 1730 isNaN(p2[2]) 1731 ) { 1732 continue; 1733 } 1734 1735 if (p2[2] === y) { 1736 if (p2[1] === x) { 1737 return 1; 1738 } 1739 if (p1[2] === y && p2[1] > x === p1[1] < x) { 1740 return 0; 1741 } 1742 } 1743 1744 if (p1[2] < y !== p2[2] < y) { 1745 // Crossing 1746 sign = 2 * (p2[2] > p1[2] ? 1 : 0) - 1; 1747 if (p1[1] >= x) { 1748 if (p2[1] > x) { 1749 wn += sign; 1750 } else { 1751 d = this.det3p(p1, p2, usrCoords); 1752 if (d === 0) { 1753 // Point is on line, i.e. outside 1754 return 0; 1755 } 1756 if (d > 0 + Mat.eps === p2[2] > p1[2]) { 1757 // Right crossing 1758 wn += sign; 1759 } 1760 } 1761 } else { 1762 if (p2[1] > x) { 1763 d = this.det3p(p1, p2, usrCoords); 1764 if (d > 0 + Mat.eps === p2[2] > p1[2]) { 1765 // Right crossing 1766 wn += sign; 1767 } 1768 } 1769 } 1770 } 1771 } 1772 1773 return wn; 1774 }, 1775 1776 /** 1777 * Decides if a point (x,y) is inside of a path / polygon. 1778 * Does not work correct if the path has hole. In this case, windingNumber is the preferred method. 1779 * Implements W. Randolf Franklin's pnpoly method. 1780 * 1781 * See <a href="https://wrf.ecse.rpi.edu/Research/Short_Notes/pnpoly.html">https://wrf.ecse.rpi.edu/Research/Short_Notes/pnpoly.html</a>. 1782 * 1783 * @param {Number} x_in x-coordinate (screen or user coordinates) 1784 * @param {Number} y_in y-coordinate (screen or user coordinates) 1785 * @param {Array} path Array of points / coords determining a path, i.e. the vertices of the polygon / path. The array elements 1786 * do not have to be full points, but have to have a subobject "coords" or should be of type JXG.Coords. 1787 * @param {Number} [coord_type=JXG.COORDS_BY_SCREEN] Type of coordinates used here. 1788 * Possible values are <b>JXG.COORDS_BY_USER</b> and <b>JXG.COORDS_BY_SCREEN</b>. 1789 * Default value is JXG.COORDS_BY_SCREEN. 1790 * @param {JXG.Board} board Board object 1791 * 1792 * @returns {Boolean} if (x_in, y_in) is inside of the polygon. 1793 * @see JXG.Polygon#hasPoint 1794 * @see JXG.Polygon#pnpoly 1795 * @see JXG.Math.Geometry.windingNumber 1796 * 1797 * @example 1798 * var pol = board.create('polygon', [[-1,2], [2,2], [-1,4]]); 1799 * var p = board.create('point', [4, 3]); 1800 * var txt = board.create('text', [-1, 0.5, function() { 1801 * return 'Point A is inside of the polygon = ' + 1802 * JXG.Math.Geometry.pnpoly(p.X(), p.Y(), pol.vertices, JXG.COORDS_BY_USER, board); 1803 * }]); 1804 * 1805 * </pre><div id="JXG4656ed42-f965-4e35-bb66-c334a4529683" class="jxgbox" style="width: 300px; height: 300px;"></div> 1806 * <script type="text/javascript"> 1807 * (function() { 1808 * var board = JXG.JSXGraph.initBoard('JXG4656ed42-f965-4e35-bb66-c334a4529683', 1809 * {boundingbox: [-2, 5, 5,-2], axis: true, showcopyright: false, shownavigation: false}); 1810 * var pol = board.create('polygon', [[-1,2], [2,2], [-1,4]]); 1811 * var p = board.create('point', [4, 3]); 1812 * var txt = board.create('text', [-1, 0.5, function() { 1813 * return 'Point A is inside of the polygon = ' + JXG.Math.Geometry.pnpoly(p.X(), p.Y(), pol.vertices, JXG.COORDS_BY_USER, board); 1814 * }]); 1815 * 1816 * })(); 1817 * 1818 * </script><pre> 1819 * 1820 */ 1821 pnpoly: function (x_in, y_in, path, coord_type, board) { 1822 var i, j, vi, vj, len, 1823 x, y, crds, 1824 v = path, 1825 isIn = false; 1826 1827 if (coord_type === Const.COORDS_BY_USER) { 1828 crds = new Coords(Const.COORDS_BY_USER, [x_in, y_in], board); 1829 x = crds.scrCoords[1]; 1830 y = crds.scrCoords[2]; 1831 } else { 1832 x = x_in; 1833 y = y_in; 1834 } 1835 1836 len = path.length; 1837 for (i = 0, j = len - 2; i < len - 1; j = i++) { 1838 vi = Type.exists(v[i].coords) ? v[i].coords : v[i]; 1839 vj = Type.exists(v[j].coords) ? v[j].coords : v[j]; 1840 1841 if ( 1842 vi.scrCoords[2] > y !== vj.scrCoords[2] > y && 1843 x < 1844 ((vj.scrCoords[1] - vi.scrCoords[1]) * (y - vi.scrCoords[2])) / 1845 (vj.scrCoords[2] - vi.scrCoords[2]) + 1846 vi.scrCoords[1] 1847 ) { 1848 isIn = !isIn; 1849 } 1850 } 1851 1852 return isIn; 1853 }, 1854 1855 /****************************************/ 1856 /**** INTERSECTIONS ****/ 1857 /****************************************/ 1858 1859 /** 1860 * Generate the function which computes the coordinates of the intersection point. 1861 * Primarily used in {@link JXG.Point.createIntersectionPoint}. 1862 * @param {JXG.Board} board object 1863 * @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. 1864 * i determines the intersection point if two points are available: <ul> 1865 * <li>i==0: use the positive square root,</li> 1866 * <li>i==1: use the negative square root.</li></ul> 1867 * @param {Boolean} alwaysintersect. Flag that determines if segments and arc can have an outer intersection point 1868 * on their defining line or circle. 1869 * @returns {Function} Function returning a {@link JXG.Coords} object that determines 1870 * the intersection point. 1871 * 1872 * @see JXG.Point.createIntersectionPoint 1873 */ 1874 intersectionFunction: function (board, el1, el2, i, j, alwaysintersect) { 1875 var func, 1876 that = this, 1877 el1_isArcType = false, 1878 el2_isArcType = false; 1879 1880 el1_isArcType = 1881 el1.elementClass === Const.OBJECT_CLASS_CURVE && 1882 (el1.type === Const.OBJECT_TYPE_ARC || el1.type === Const.OBJECT_TYPE_SECTOR) 1883 ? true 1884 : false; 1885 el2_isArcType = 1886 el2.elementClass === Const.OBJECT_CLASS_CURVE && 1887 (el2.type === Const.OBJECT_TYPE_ARC || el2.type === Const.OBJECT_TYPE_SECTOR) 1888 ? true 1889 : false; 1890 1891 if ( 1892 (el1.elementClass === Const.OBJECT_CLASS_CURVE || 1893 el2.elementClass === Const.OBJECT_CLASS_CURVE) && 1894 (el1.elementClass === Const.OBJECT_CLASS_CURVE || 1895 el1.elementClass === Const.OBJECT_CLASS_CIRCLE) && 1896 (el2.elementClass === Const.OBJECT_CLASS_CURVE || 1897 el2.elementClass === Const.OBJECT_CLASS_CIRCLE) /*&& 1898 !(el1_isArcType && el2_isArcType)*/ 1899 ) { 1900 // curve - curve 1901 // with the exception that both elements are arc types 1902 /** @ignore */ 1903 func = function () { 1904 return that.meetCurveCurve(el1, el2, i, j, el1.board); 1905 }; 1906 } else if ( 1907 (el1.elementClass === Const.OBJECT_CLASS_CURVE && 1908 !el1_isArcType && 1909 el2.elementClass === Const.OBJECT_CLASS_LINE) || 1910 (el2.elementClass === Const.OBJECT_CLASS_CURVE && 1911 !el2_isArcType && 1912 el1.elementClass === Const.OBJECT_CLASS_LINE) 1913 ) { 1914 // curve - line (this includes intersections between conic sections and lines) 1915 // with the exception that the curve is of arc type 1916 /** @ignore */ 1917 func = function () { 1918 return that.meetCurveLine(el1, el2, i, el1.board, Type.evaluate(alwaysintersect)); 1919 }; 1920 } else if ( 1921 el1.type === Const.OBJECT_TYPE_POLYGON || 1922 el2.type === Const.OBJECT_TYPE_POLYGON 1923 ) { 1924 // polygon - other 1925 // Uses the Greiner-Hormann clipping algorithm 1926 // Not implemented: polygon - point 1927 1928 if (el1.elementClass === Const.OBJECT_CLASS_LINE) { 1929 // line - path 1930 /** @ignore */ 1931 func = function () { 1932 var first1 = el1.evalVisProp('straightfirst'), 1933 last1 = el1.evalVisProp('straightlast'), 1934 first2 = el2.evalVisProp('straightfirst'), 1935 last2 = el2.evalVisProp('straightlast'), 1936 a_not; 1937 1938 a_not = (!Type.evaluate(alwaysintersect) && (!first1 || !last1 || !first2 || !last2)); 1939 return that.meetPolygonLine(el2, el1, i, el1.board, a_not); 1940 }; 1941 } else if (el2.elementClass === Const.OBJECT_CLASS_LINE) { 1942 // path - line 1943 /** @ignore */ 1944 func = function () { 1945 var first1 = el1.evalVisProp('straightfirst'), 1946 last1 = el1.evalVisProp('straightlast'), 1947 first2 = el2.evalVisProp('straightfirst'), 1948 last2 = el2.evalVisProp('straightlast'), 1949 a_not; 1950 1951 a_not = (!Type.evaluate(alwaysintersect) && (!first1 || !last1 || !first2 || !last2)); 1952 return that.meetPolygonLine(el1, el2, i, el1.board, a_not); 1953 }; 1954 } else { 1955 // path - path 1956 /** @ignore */ 1957 func = function () { 1958 return that.meetPathPath(el1, el2, i, el1.board); 1959 }; 1960 } 1961 } else if ( 1962 el1.elementClass === Const.OBJECT_CLASS_LINE && 1963 el2.elementClass === Const.OBJECT_CLASS_LINE 1964 ) { 1965 // line - line, lines may also be segments. 1966 /** @ignore */ 1967 func = function () { 1968 var res, 1969 c, 1970 first1 = el1.evalVisProp('straightfirst'), 1971 last1 = el1.evalVisProp('straightlast'), 1972 first2 = el2.evalVisProp('straightfirst'), 1973 last2 = el2.evalVisProp('straightlast'); 1974 1975 /** 1976 * If one of the lines is a segment or ray and 1977 * the intersection point should disappear if outside 1978 * of the segment or ray we call 1979 * meetSegmentSegment 1980 */ 1981 if ( 1982 !Type.evaluate(alwaysintersect) && 1983 (!first1 || !last1 || !first2 || !last2) 1984 ) { 1985 res = that.meetSegmentSegment( 1986 el1.point1.coords.usrCoords, 1987 el1.point2.coords.usrCoords, 1988 el2.point1.coords.usrCoords, 1989 el2.point2.coords.usrCoords 1990 ); 1991 1992 if ( 1993 (!first1 && res[1] < 0) || 1994 (!last1 && res[1] > 1) || 1995 (!first2 && res[2] < 0) || 1996 (!last2 && res[2] > 1) 1997 ) { 1998 // Non-existent 1999 c = [0, NaN, NaN]; 2000 } else { 2001 c = res[0]; 2002 } 2003 2004 return new Coords(Const.COORDS_BY_USER, c, el1.board); 2005 } 2006 2007 return that.meet(el1.stdform, el2.stdform, i, el1.board); 2008 }; 2009 } else { 2010 // All other combinations of circles and lines, 2011 // Arc types are treated as circles. 2012 /** @ignore */ 2013 func = function () { 2014 var res = that.meet(el1.stdform, el2.stdform, i, el1.board), 2015 has = true, 2016 first, 2017 last, 2018 r; 2019 2020 if (Type.evaluate(alwaysintersect)) { 2021 return res; 2022 } 2023 if (el1.elementClass === Const.OBJECT_CLASS_LINE) { 2024 first = el1.evalVisProp('straightfirst'); 2025 last = el1.evalVisProp('straightlast'); 2026 if (!first || !last) { 2027 r = that.affineRatio(el1.point1.coords, el1.point2.coords, res); 2028 if ((!last && r > 1 + Mat.eps) || (!first && r < 0 - Mat.eps)) { 2029 return new Coords(JXG.COORDS_BY_USER, [0, NaN, NaN], el1.board); 2030 } 2031 } 2032 } 2033 if (el2.elementClass === Const.OBJECT_CLASS_LINE) { 2034 first = el2.evalVisProp('straightfirst'); 2035 last = el2.evalVisProp('straightlast'); 2036 if (!first || !last) { 2037 r = that.affineRatio(el2.point1.coords, el2.point2.coords, res); 2038 if ((!last && r > 1 + Mat.eps) || (!first && r < 0 - Mat.eps)) { 2039 return new Coords(JXG.COORDS_BY_USER, [0, NaN, NaN], el1.board); 2040 } 2041 } 2042 } 2043 if (el1_isArcType) { 2044 has = that.coordsOnArc(el1, res); 2045 if (has && el2_isArcType) { 2046 has = that.coordsOnArc(el2, res); 2047 } 2048 if (!has) { 2049 return new Coords(JXG.COORDS_BY_USER, [0, NaN, NaN], el1.board); 2050 } 2051 } 2052 return res; 2053 }; 2054 } 2055 2056 return func; 2057 }, 2058 2059 otherIntersectionFunction: function (input, others, alwaysintersect, precision) { 2060 var func, board, 2061 el1, el2, 2062 that = this; 2063 2064 el1 = input[0]; 2065 el2 = input[1]; 2066 board = el1.board; 2067 /** @ignore */ 2068 func = function () { 2069 var i, k, c, d, 2070 isClose, 2071 le = others.length, 2072 eps = Type.evaluate(precision); 2073 2074 for (i = le; i >= 0; i--) { 2075 if (el1.elementClass === Const.OBJECT_CLASS_CIRCLE && 2076 [Const.OBJECT_CLASS_CIRCLE, Const.OBJECT_CLASS_LINE].indexOf(el2.elementClass) >= 0) { 2077 // circle, circle|line 2078 c = that.meet(el1.stdform, el2.stdform, i, board); 2079 } else if (el1.elementClass === Const.OBJECT_CLASS_CURVE && 2080 [Const.OBJECT_CLASS_CURVE, Const.OBJECT_CLASS_CIRCLE].indexOf(el2.elementClass) >= 0) { 2081 // curve, circle|curve 2082 c = that.meetCurveCurve(el1, el2, i, 0, board); 2083 } else if (el1.elementClass === Const.OBJECT_CLASS_CURVE && el2.elementClass === Const.OBJECT_CLASS_LINE) { 2084 // curve, line 2085 if (Type.exists(el1.dataX)) { 2086 c = JXG.Math.Geometry.meetCurveLine(el1, el2, i, el1.board, Type.evaluate(alwaysintersect)); 2087 } else { 2088 c = JXG.Math.Geometry.meetCurveLineContinuous(el1, el2, i, el1.board); 2089 } 2090 } 2091 2092 if (c === undefined) { 2093 // Intersection point does not exist 2094 continue; 2095 } 2096 2097 // If the intersection is close to one of the points in other 2098 // we have to search for another intersection point. 2099 isClose = false; 2100 for (k = 0; !isClose && k < le; k++) { 2101 if (Type.exists(c) && Type.exists(c.distance)) { 2102 d = c.distance(JXG.COORDS_BY_USER, others[k].coords); 2103 if (d < eps) { 2104 isClose = true; 2105 } 2106 } 2107 } 2108 if (!isClose) { 2109 // We are done, the intersection is away from any other 2110 // intersection point. 2111 return c; 2112 } 2113 } 2114 // Otherwise we return the last intersection point 2115 return c; 2116 }; 2117 return func; 2118 }, 2119 2120 /** 2121 * Returns true if the coordinates are on the arc element, 2122 * false otherwise. Usually, coords is an intersection 2123 * on the circle line. Now it is decided if coords are on the 2124 * circle restricted to the arc line. 2125 * @param {Arc} arc arc or sector element 2126 * @param {JXG.Coords} coords Coords object of an intersection 2127 * @returns {Boolean} 2128 * @private 2129 */ 2130 coordsOnArc: function (arc, coords) { 2131 var angle = this.rad(arc.radiuspoint, arc.center, coords.usrCoords.slice(1)), 2132 alpha = 0.0, 2133 beta = this.rad(arc.radiuspoint, arc.center, arc.anglepoint), 2134 ev_s = arc.evalVisProp('selection'); 2135 2136 if (arc.evalVisProp('orientation') === 'clockwise') { 2137 angle = 2 * Math.PI - angle; 2138 beta = 2 * Math.PI - beta; 2139 } 2140 2141 if ((ev_s === "minor" && beta > Math.PI) || (ev_s === "major" && beta < Math.PI)) { 2142 alpha = beta; 2143 beta = 2 * Math.PI; 2144 } 2145 if (angle < alpha || angle > beta) { 2146 return false; 2147 } 2148 return true; 2149 }, 2150 2151 /** 2152 * Computes the intersection of a pair of lines, circles or both. 2153 * It uses the internal data array stdform of these elements. 2154 * @param {Array} el1 stdform of the first element (line or circle) 2155 * @param {Array} el2 stdform of the second element (line or circle) 2156 * @param {Number|Function} i Index of the intersection point that should be returned. 2157 * @param board Reference to the board. 2158 * @returns {JXG.Coords} Coordinates of one of the possible two or more intersection points. 2159 * Which point will be returned is determined by i. 2160 */ 2161 meet: function (el1, el2, i, board) { 2162 var result, 2163 eps = Mat.eps; 2164 2165 if (Math.abs(el1[3]) < eps && Math.abs(el2[3]) < eps) { 2166 // line line 2167 result = this.meetLineLine(el1, el2, i, board); 2168 } else if (Math.abs(el1[3]) >= eps && Math.abs(el2[3]) < eps) { 2169 // circle line 2170 result = this.meetLineCircle(el2, el1, i, board); 2171 } else if (Math.abs(el1[3]) < eps && Math.abs(el2[3]) >= eps) { 2172 // line circle 2173 result = this.meetLineCircle(el1, el2, i, board); 2174 } else { 2175 // circle circle 2176 result = this.meetCircleCircle(el1, el2, i, board); 2177 } 2178 2179 return result; 2180 }, 2181 2182 /** 2183 * Intersection of the line with the board 2184 * @param {Array} line stdform of the line in screen coordinates 2185 * @param {JXG.Board} board reference to a board. 2186 * @param {Number} margin optional margin, to avoid the display of the small sides of lines. 2187 * @returns {Array} [intersection coords 1, intersection coords 2] 2188 */ 2189 meetLineBoard: function (line, board, margin) { 2190 // Intersect the line with the four borders of the board. 2191 var s = [], 2192 intersect1, 2193 intersect2, 2194 i, j; 2195 2196 if (!Type.exists(margin)) { 2197 margin = 0; 2198 } 2199 2200 // top 2201 s[0] = Mat.crossProduct(line, [margin, 0, 1]); 2202 // left 2203 s[1] = Mat.crossProduct(line, [margin, 1, 0]); 2204 // bottom 2205 s[2] = Mat.crossProduct(line, [-margin - board.canvasHeight, 0, 1]); 2206 // right 2207 s[3] = Mat.crossProduct(line, [-margin - board.canvasWidth, 1, 0]); 2208 2209 // Normalize the intersections 2210 for (i = 0; i < 4; i++) { 2211 if (Math.abs(s[i][0]) > Mat.eps) { 2212 for (j = 2; j > 0; j--) { 2213 s[i][j] /= s[i][0]; 2214 } 2215 s[i][0] = 1.0; 2216 } 2217 } 2218 2219 // line is parallel to "left", take "top" and "bottom" 2220 if (Math.abs(s[1][0]) < Mat.eps) { 2221 intersect1 = s[0]; // top 2222 intersect2 = s[2]; // bottom 2223 // line is parallel to "top", take "left" and "right" 2224 } else if (Math.abs(s[0][0]) < Mat.eps) { 2225 intersect1 = s[1]; // left 2226 intersect2 = s[3]; // right 2227 // left intersection out of board (above) 2228 } else if (s[1][2] < 0) { 2229 intersect1 = s[0]; // top 2230 2231 // right intersection out of board (below) 2232 if (s[3][2] > board.canvasHeight) { 2233 intersect2 = s[2]; // bottom 2234 } else { 2235 intersect2 = s[3]; // right 2236 } 2237 // left intersection out of board (below) 2238 } else if (s[1][2] > board.canvasHeight) { 2239 intersect1 = s[2]; // bottom 2240 2241 // right intersection out of board (above) 2242 if (s[3][2] < 0) { 2243 intersect2 = s[0]; // top 2244 } else { 2245 intersect2 = s[3]; // right 2246 } 2247 } else { 2248 intersect1 = s[1]; // left 2249 2250 // right intersection out of board (above) 2251 if (s[3][2] < 0) { 2252 intersect2 = s[0]; // top 2253 // right intersection out of board (below) 2254 } else if (s[3][2] > board.canvasHeight) { 2255 intersect2 = s[2]; // bottom 2256 } else { 2257 intersect2 = s[3]; // right 2258 } 2259 } 2260 2261 return [ 2262 new Coords(Const.COORDS_BY_SCREEN, intersect1.slice(1), board), 2263 new Coords(Const.COORDS_BY_SCREEN, intersect2.slice(1), board) 2264 ]; 2265 }, 2266 2267 /** 2268 * Intersection of two lines. 2269 * @param {Array} l1 stdform of the first line 2270 * @param {Array} l2 stdform of the second line 2271 * @param {number} i unused 2272 * @param {JXG.Board} board Reference to the board. 2273 * @returns {JXG.Coords} Coordinates of the intersection point. 2274 */ 2275 meetLineLine: function (l1, l2, i, board) { 2276 var s = isNaN(l1[5] + l2[5]) ? [0, 0, 0] : Mat.crossProduct(l1, l2); 2277 2278 // Make intersection of parallel lines more robust: 2279 if (Math.abs(s[0]) < 1.0e-14) { 2280 s[0] = 0; 2281 } 2282 return new Coords(Const.COORDS_BY_USER, s, board); 2283 }, 2284 2285 /** 2286 * Intersection of line and circle. 2287 * @param {Array} lin stdform of the line 2288 * @param {Array} circ stdform of the circle 2289 * @param {number|function} i number of the returned intersection point. 2290 * i==0: use the positive square root, 2291 * i==1: use the negative square root. 2292 * @param {JXG.Board} board Reference to a board. 2293 * @returns {JXG.Coords} Coordinates of the intersection point 2294 */ 2295 meetLineCircle: function (lin, circ, i, board) { 2296 var a, b, c, d, n, A, B, C, k, t; 2297 2298 // Radius is zero, return center of circle 2299 if (circ[4] < Mat.eps) { 2300 if (Math.abs(Mat.innerProduct([1, circ[6], circ[7]], lin, 3)) < Mat.eps) { 2301 return new Coords(Const.COORDS_BY_USER, circ.slice(6, 8), board); 2302 } 2303 2304 return new Coords(Const.COORDS_BY_USER, [NaN, NaN], board); 2305 } 2306 c = circ[0]; 2307 b = circ.slice(1, 3); 2308 a = circ[3]; 2309 d = lin[0]; 2310 n = lin.slice(1, 3); 2311 2312 // Line is assumed to be normalized. Therefore, nn==1 and we can skip some operations: 2313 /* 2314 var nn = n[0]*n[0]+n[1]*n[1]; 2315 A = a*nn; 2316 B = (b[0]*n[1]-b[1]*n[0])*nn; 2317 C = a*d*d - (b[0]*n[0]+b[1]*n[1])*d + c*nn; 2318 */ 2319 A = a; 2320 B = b[0] * n[1] - b[1] * n[0]; 2321 C = a * d * d - (b[0] * n[0] + b[1] * n[1]) * d + c; 2322 2323 k = B * B - 4 * A * C; 2324 if (k > -Mat.eps * Mat.eps) { 2325 k = Math.sqrt(Math.abs(k)); 2326 t = [(-B + k) / (2 * A), (-B - k) / (2 * A)]; 2327 2328 return Type.evaluate(i) === 0 2329 ? new Coords( 2330 Const.COORDS_BY_USER, 2331 [-t[0] * -n[1] - d * n[0], -t[0] * n[0] - d * n[1]], 2332 board 2333 ) 2334 : new Coords( 2335 Const.COORDS_BY_USER, 2336 [-t[1] * -n[1] - d * n[0], -t[1] * n[0] - d * n[1]], 2337 board 2338 ); 2339 } 2340 2341 return new Coords(Const.COORDS_BY_USER, [0, 0, 0], board); 2342 }, 2343 2344 /** 2345 * Intersection of two circles. 2346 * @param {Array} circ1 stdform of the first circle 2347 * @param {Array} circ2 stdform of the second circle 2348 * @param {number|function} i number of the returned intersection point. 2349 * i==0: use the positive square root, 2350 * i==1: use the negative square root. 2351 * @param {JXG.Board} board Reference to the board. 2352 * @returns {JXG.Coords} Coordinates of the intersection point 2353 */ 2354 meetCircleCircle: function (circ1, circ2, i, board) { 2355 var radicalAxis; 2356 2357 // Radius is zero, return center of circle, if on other circle 2358 if (circ1[4] < Mat.eps) { 2359 if ( 2360 Math.abs(this.distance(circ1.slice(6, 8), circ2.slice(6, 8)) - circ2[5]) < 2361 Mat.eps 2362 ) { 2363 return new Coords(Const.COORDS_BY_USER, circ1.slice(6, 8), board); 2364 } 2365 2366 return new Coords(Const.COORDS_BY_USER, [0, 0, 0], board); 2367 } 2368 2369 // Radius is zero, return center of circle, if on other circle 2370 if (circ2[4] < Mat.eps) { 2371 if ( 2372 Math.abs(this.distance(circ2.slice(6, 8), circ1.slice(6, 8)) - circ1[5]) < 2373 Mat.eps 2374 ) { 2375 return new Coords(Const.COORDS_BY_USER, circ2.slice(6, 8), board); 2376 } 2377 2378 return new Coords(Const.COORDS_BY_USER, [0, 0, 0], board); 2379 } 2380 2381 radicalAxis = [ 2382 circ2[3] * circ1[0] - circ1[3] * circ2[0], 2383 circ2[3] * circ1[1] - circ1[3] * circ2[1], 2384 circ2[3] * circ1[2] - circ1[3] * circ2[2], 2385 0, 2386 1, 2387 Infinity, 2388 Infinity, 2389 Infinity 2390 ]; 2391 radicalAxis = Mat.normalize(radicalAxis); 2392 2393 return this.meetLineCircle(radicalAxis, circ1, i, board); 2394 }, 2395 2396 /** 2397 * Segment-wise search for the nr-th intersection of two curves. 2398 * testSegment is always assumed to be true. 2399 * 2400 * @param {JXG.Curve} c1 Curve, Line or Circle 2401 * @param {JXG.Curve} c2 Curve, Line or Circle 2402 * @param {Number} nr the nr-th intersection point will be returned 2403 * @param {JXG.Board} [board=c1.board] Reference to a board object 2404 * @returns {JXG.Coords} intersection as Coords object 2405 * 2406 * @private 2407 * @see JXG.Math.Geometry.meetCurveCurve 2408 */ 2409 meetCurveCurveDiscrete: function (c1, c2, nr, board) { 2410 var co, 2411 i = Type.evaluate(nr); 2412 2413 if (c1.bezierDegree === 3 || c2.bezierDegree === 3) { 2414 co = this.meetBezierCurveRedBlueSegments(c1, c2, i); 2415 } else { 2416 co = this.meetCurveRedBlueSegments(c1, c2, i); 2417 } 2418 return new Coords(Const.COORDS_BY_USER, co, board); 2419 }, 2420 2421 /** 2422 * Apply Newton-Raphson to search for an intersection of two curves 2423 * in a given range of the first curve. 2424 * 2425 * @param {JXG.Curve} c1 Curve, Line or Circle 2426 * @param {JXG.Curve} c2 Curve, Line or Circle 2427 * @param {Array} range Domain for the search of an intersection. The start value 2428 * for the search is chosen to be inside of that range. 2429 * @param {Boolean} testSegment If true require that t1 and t2 are inside of the allowed bounds. 2430 * @returns {Array} [[z, x, y], t1, t2, t, ||c1[t1]-c2[t2]||**2]. The last entry is set to 2431 * 10000 if the intersection is outside of the given domain (range) for the first curve. 2432 * @private 2433 * @see JXG.Math.Geometry._meetCurveCurveIterative 2434 * @see JXG.Math.Numerics.generalizedDampedNewton 2435 * @see JXG.Math.Geometry.meetCurveCurveCobyla 2436 */ 2437 meetCurveCurveNewton: function (c1, c2, range1, range2, testSegment) { 2438 var t1, t2, 2439 co, r, 2440 inphi = (Math.sqrt(5) - 1) * 0.5, 2441 damp = 0.85, // ( 2442 eps3 = Mat.eps * Mat.eps * Mat.eps, 2443 eps2 = Mat.eps * Mat.eps, 2444 2445 ma1 = c1.maxX(), 2446 mi1 = c1.minX(), 2447 ma2 = c2.maxX(), 2448 mi2 = c2.minX(), 2449 2450 F = function(t, n) { 2451 var f1 = c1.Ft(t[0]), 2452 f2 = c2.Ft(t[1]), 2453 e = f1[1] - f2[1], 2454 f = f1[2] - f2[2]; 2455 2456 return [e, f]; 2457 }, 2458 D = function(t, n) { 2459 var h = Mat.eps, 2460 h2 = 2 * h, 2461 f1_1 = c1.Ft(t[0] - h), 2462 f1_2 = c1.Ft(t[0] + h), 2463 f2_1 = c2.Ft(t[1] - h), 2464 f2_2 = c2.Ft(t[1] + h); 2465 return [ 2466 [ (f1_2[1] - f1_1[1]) / h2, 2467 -(f2_2[1] - f2_1[1]) / h2], 2468 [ (f1_2[2] - f1_1[2]) / h2, 2469 -(f2_2[2] - f2_1[2]) / h2] 2470 ]; 2471 }; 2472 2473 t1 = range1[0] + (range1[1] - range1[0]) * (1 - inphi); 2474 t2 = range2[0] + (range2[1] - range2[0]) * (1 - inphi); 2475 2476 // Use damped Newton 2477 // r = Numerics.generalizedDampedNewtonCurves(c1, c2, t1, t2, damp, eps3); 2478 r = Numerics.generalizedDampedNewton(F, D, 2, [t1, t2], damp, eps3, 40); 2479 // r: [t1, t2, F2] 2480 2481 t1 = r[0][0]; 2482 t2 = r[0][1]; 2483 co = c1.Ft(t1); 2484 2485 if ( 2486 t1 < range1[0] - Mat.eps || t1 > range1[1] + Mat.eps || 2487 t2 < range2[0] - Mat.eps || t2 > range2[1] + Mat.eps || 2488 (testSegment && 2489 (t1 < mi1 - eps2 || t1 > ma1 + eps2 || 2490 t2 < mi2 - eps2 || t2 > ma2 + eps2) 2491 ) 2492 ) { 2493 // Damped-Newton found solution outside of range 2494 return [co, t1, t2, 10000]; 2495 } 2496 // console.log(t1, r[3]) 2497 2498 return [co, t1, t2, r[1]]; 2499 }, 2500 2501 /** 2502 * Return a list of the (at most) first i intersection points of two curves. 2503 * Computed iteratively. 2504 * 2505 * @param {JXG.Curve} c1 Curve, Line or Circle 2506 * @param {JXG.Curve} c2 Curve, Line or Circle 2507 * @param {Number} low Lower bound of the search domain (between [0, 1]) 2508 * @param {Number} up Upper bound of the search domain (between [0, 1]) 2509 * @param {Number} i Return a list of the first i intersection points 2510 * @param {Boolean} testSegment If true require that t1 and t2 are inside of the allowed bounds. 2511 * @returns {Array} List of the first i intersection points, given by the parameter t. 2512 * @private 2513 * @see JXG.Math.Geometry.meetCurveCurveNewton 2514 * @see JXG.Math.Geometry.meetCurveCurve 2515 */ 2516 _meetCurveCurveIterative: function(c1, c2, range1, range2, i, testSegment) { 2517 var ret, 2518 t1,// t2, 2519 low1, low2, up1, up2, 2520 eps = Mat.eps * 100, // Minimum difference between zeros 2521 j1, j2, 2522 steps = 20, 2523 d1, d2, 2524 zeros = []; 2525 2526 low1 = range1[0]; 2527 up1 = range1[1]; 2528 low2 = range2[0]; 2529 up2 = range2[1]; 2530 if (up1 < low1 || up2 < low2) { 2531 return []; 2532 } 2533 2534 // console.log('DO iterative', [low1, up1], [low2, up2]) 2535 2536 d1 = (up1 - low1) / steps; 2537 d2 = (up2 - low2) / steps; 2538 for (j1 = 0; j1 < steps; j1++) { 2539 for (j2 = 0; j2 < steps; j2++) { 2540 2541 ret = this.meetCurveCurveNewton(c1, c2, 2542 [low1 + j1 * d1, low1 + (j1 + 1) * d1], 2543 [low2 + j2 * d2, low2 + (j2 + 1) * d2], 2544 testSegment); 2545 2546 if (ret[3] < Mat.eps) { 2547 t1 = ret[1]; 2548 // t2 = ret[2]; 2549 // console.log("\tFOUND", t1, t2, c1.Ft(t1)[2]) 2550 zeros = zeros.concat([t1]); 2551 zeros = Type.toUniqueArrayFloat(zeros, eps); 2552 // console.log(zeros, i) 2553 if (zeros.length > i) { 2554 return zeros; 2555 } 2556 } 2557 } 2558 } 2559 2560 return zeros; 2561 }, 2562 2563 /** 2564 * Compute an intersection of the curves c1 and c2. 2565 * We want to find values t1, t2 such that 2566 * c1(t1) = c2(t2), i.e. (c1_x(t1) - c2_x(t2), c1_y(t1) - c2_y(t2)) = (0, 0). 2567 * 2568 * Available methods: 2569 * <ul> 2570 * <li> discrete, segment-wise intersections 2571 * <li> generalized damped Newton-Raphson 2572 * </ul> 2573 * 2574 * Segment-wise intersection is more stable, but has problems with tangent points. 2575 * Damped Newton-Raphson converges very rapidly but sometimes behaves chaotic. 2576 * 2577 * @param {JXG.Curve} c1 Curve, Line or Circle 2578 * @param {JXG.Curve} c2 Curve, Line or Circle 2579 * @param {Number|Function} nr the nr-th intersection point will be returned. For backwards compatibility: 2580 * if method='newton' and nr is not an integer, {@link JXG.Math.Numerics.generalizedNewton} is called 2581 * directly with nr as start value (not recommended). 2582 * @param {Number} t2ini not longer used. Must be supplied and is ignored. 2583 * @param {JXG.Board} [board=c1.board] Reference to a board object. 2584 * @param {String} [method] Intersection method, possible values are 'newton' and 'segment'. 2585 * If both curves are given by functions (assumed to be continuous), 'newton' is the default, otherwise 2586 * 'segment' is the default. 2587 * @parame {Boolean} testSegment If true require that the intersection is inside of the allowed bounds for both elements (in _meetCurveCurveIterative) 2588 * @returns {JXG.Coords} intersection point 2589 * 2590 * @see JXG.Math.Geometry.meetCurveCurveDiscrete 2591 * @see JXG.Math.Geometry._meetCurveCurveIterative 2592 */ 2593 meetCurveCurve: function (c1, c2, nr, t2ini, board, method, testSegment) { 2594 var co, 2595 zeros, 2596 mi1, ma1, mi2, ma2, 2597 i = Type.evaluate(nr); 2598 2599 board = board || c1.board; 2600 if (method === 'segment' || Type.exists(c1.dataX) || Type.exists(c2.dataX)) { 2601 // Discrete data points, i.e. x-coordinates of c1 or c2 are given in an array) 2602 return this.meetCurveCurveDiscrete(c1, c2, i, board); 2603 } 2604 2605 // Outdated: 2606 // Backwards compatibility if nr is not a positive integer then 2607 // generalizedNewton is still used. 2608 if (Type.exists(method) && method === 'newton' && i < 0 || parseInt(i) !== i) { 2609 co = Numerics.generalizedNewton(c1, c2, i, t2ini); 2610 return new Coords(Const.COORDS_BY_USER, co, board); 2611 } 2612 2613 // Method 'newton' 2614 mi1 = c1.minX(); 2615 ma1 = c1.maxX(); 2616 mi2 = c2.minX(); 2617 ma2 = c2.maxX(); 2618 2619 // console.time('curvecurve') 2620 zeros = this._meetCurveCurveIterative(c1, c2, [mi1, ma1], [mi2, ma2], i, testSegment); 2621 // console.timeEnd('curvecurve') 2622 2623 if (zeros.length > i) { 2624 co = c1.Ft(zeros[i]); 2625 } else { 2626 return [0, NaN, NaN]; 2627 } 2628 2629 return new Coords(Const.COORDS_BY_USER, co, board); 2630 }, 2631 2632 /** 2633 * Intersection of curve with line, 2634 * Order of input does not matter for el1 and el2. 2635 * From version 0.99.7 on this method calls 2636 * {@link JXG.Math.Geometry.meetCurveLineDiscrete}. 2637 * If higher precision is needed, {@link JXG.Math.Geometry.meetCurveLineContinuous} 2638 * has to be used. 2639 * 2640 * @param {JXG.Curve|JXG.Line} el1 Curve or Line 2641 * @param {JXG.Curve|JXG.Line} el2 Curve or Line 2642 * @param {Number|Function} nr the nr-th intersection point will be returned. 2643 * @param {JXG.Board} [board=el1.board] Reference to a board object. 2644 * @param {Boolean} alwaysIntersect If false just the segment between the two defining points are tested for intersection 2645 * @returns {JXG.Coords} Intersection point. In case no intersection point is detected, 2646 * the ideal point [0,1,0] is returned. 2647 */ 2648 meetCurveLine: function (el1, el2, nr, board, alwaysIntersect) { 2649 var v = [0, NaN, NaN], 2650 cu, 2651 li; 2652 2653 if (!Type.exists(board)) { 2654 board = el1.board; 2655 } 2656 2657 if (el1.elementClass === Const.OBJECT_CLASS_CURVE) { 2658 cu = el1; 2659 li = el2; 2660 } else { 2661 cu = el2; 2662 li = el1; 2663 } 2664 2665 if (Type.exists(cu.dataX)) { 2666 // We use the discrete version if 2667 // the curve is not a parametric curve, e.g. implicit plots 2668 v = this.meetCurveLineDiscrete(cu, li, nr, board, !alwaysIntersect); 2669 } else { 2670 v = this.meetCurveLineContinuous(cu, li, nr, board, !alwaysIntersect); 2671 // v = this.meetCurveCurve(cu, li, nr, 0, board, 'newton', !alwaysIntersect); 2672 } 2673 2674 return v; 2675 }, 2676 2677 /** 2678 * Intersection of line and curve, continuous case. 2679 * Finds the nr-th intersection point 2680 * Uses {@link JXG.Math.Geometry.meetCurveLineDiscrete} as a first approximation. 2681 * A more exact solution is then found with {@link JXG.Math.Numerics.root}. 2682 * 2683 * @param {JXG.Curve} cu Curve 2684 * @param {JXG.Line} li Line 2685 * @param {NumberFunction} nr Will return the nr-th intersection point. 2686 * @param {JXG.Board} board 2687 * @param {Boolean} testSegment Test if intersection has to be inside of the segment or somewhere on the 2688 * line defined by the segment 2689 * @returns {JXG.Coords} Coords object containing the intersection. 2690 */ 2691 meetCurveLineContinuous: function (cu, li, nr, board, testSegment) { 2692 var func0, func1, 2693 t, v, x, y, z, 2694 eps = Mat.eps, 2695 epsLow = Mat.eps, 2696 steps, 2697 delta, 2698 tnew, tmin, fmin, 2699 i, ft; 2700 2701 v = this.meetCurveLineDiscrete(cu, li, nr, board, testSegment); 2702 x = v.usrCoords[1]; 2703 y = v.usrCoords[2]; 2704 2705 func0 = function (t) { 2706 var c1, c2; 2707 2708 if (t > cu.maxX() || t < cu.minX()) { 2709 return Infinity; 2710 } 2711 c1 = cu.X(t) - x; 2712 c2 = cu.Y(t) - y; 2713 return c1 * c1 + c2 * c2; 2714 // return c1 * (cu.X(t + h) - cu.X(t - h)) + c2 * (cu.Y(t + h) - cu.Y(t - h)) / h; 2715 }; 2716 2717 func1 = function (t) { 2718 var v = li.stdform[0] + li.stdform[1] * cu.X(t) + li.stdform[2] * cu.Y(t); 2719 return v * v; 2720 }; 2721 2722 // Find t 2723 steps = 50; 2724 delta = (cu.maxX() - cu.minX()) / steps; 2725 tnew = cu.minX(); 2726 fmin = 0.0001; //eps; 2727 tmin = NaN; 2728 for (i = 0; i < steps; i++) { 2729 t = Numerics.root(func0, [ 2730 Math.max(tnew, cu.minX()), 2731 Math.min(tnew + delta, cu.maxX()) 2732 ]); 2733 ft = Math.abs(func0(t)); 2734 if (ft <= fmin) { 2735 fmin = ft; 2736 tmin = t; 2737 if (fmin < eps) { 2738 break; 2739 } 2740 } 2741 2742 tnew += delta; 2743 } 2744 t = tmin; 2745 // Compute "exact" t 2746 t = Numerics.root(func1, [ 2747 Math.max(t - delta, cu.minX()), 2748 Math.min(t + delta, cu.maxX()) 2749 ]); 2750 2751 ft = func1(t); 2752 // Is the point on the line? 2753 if (isNaN(ft) || Math.abs(ft) > epsLow) { 2754 z = 0.0; //NaN; 2755 } else { 2756 z = 1.0; 2757 } 2758 2759 return new Coords(Const.COORDS_BY_USER, [z, cu.X(t), cu.Y(t)], board); 2760 }, 2761 2762 /** 2763 * Intersection of line and curve, discrete case. 2764 * Segments are treated as lines. 2765 * Finding the nr-th intersection point should work for all nr. 2766 * @param {JXG.Curve} cu 2767 * @param {JXG.Line} li 2768 * @param {Number|Function} nr 2769 * @param {JXG.Board} board 2770 * @param {Boolean} testSegment Test if intersection has to be inside of the segment or somewhere on the 2771 * line defined by the segment 2772 * 2773 * @returns {JXG.Coords} Intersection point. In case no intersection point is detected, 2774 * the ideal point [0,1,0] is returned. 2775 */ 2776 meetCurveLineDiscrete: function (cu, li, nr, board, testSegment) { 2777 var i, j, 2778 n = Type.evaluate(nr), 2779 p1, p2, 2780 p, q, 2781 lip1 = li.point1.coords.usrCoords, 2782 lip2 = li.point2.coords.usrCoords, 2783 d, res, 2784 cnt = 0, 2785 len = cu.numberPoints, 2786 ev_sf = li.evalVisProp('straightfirst'), 2787 ev_sl = li.evalVisProp('straightlast'); 2788 2789 // In case, no intersection will be found we will take this 2790 q = new Coords(Const.COORDS_BY_USER, [0, NaN, NaN], board); 2791 2792 if (lip1[0] === 0.0) { 2793 lip1 = [1, lip2[1] + li.stdform[2], lip2[2] - li.stdform[1]]; 2794 } else if (lip2[0] === 0.0) { 2795 lip2 = [1, lip1[1] + li.stdform[2], lip1[2] - li.stdform[1]]; 2796 } 2797 2798 p2 = cu.points[0].usrCoords; 2799 for (i = 1; i < len; i += cu.bezierDegree) { 2800 p1 = p2.slice(0); 2801 p2 = cu.points[i].usrCoords; 2802 d = this.distance(p1, p2); 2803 2804 // The defining points are not identical 2805 if (d > Mat.eps) { 2806 if (cu.bezierDegree === 3) { 2807 res = this.meetBeziersegmentBeziersegment( 2808 [ 2809 cu.points[i - 1].usrCoords.slice(1), 2810 cu.points[i].usrCoords.slice(1), 2811 cu.points[i + 1].usrCoords.slice(1), 2812 cu.points[i + 2].usrCoords.slice(1) 2813 ], 2814 [lip1.slice(1), lip2.slice(1)], 2815 testSegment 2816 ); 2817 } else { 2818 res = [this.meetSegmentSegment(p1, p2, lip1, lip2)]; 2819 } 2820 2821 for (j = 0; j < res.length; j++) { 2822 p = res[j]; 2823 if (0 <= p[1] && p[1] <= 1) { 2824 if (cnt === n) { 2825 /** 2826 * If the intersection point is not part of the segment, 2827 * this intersection point is set to non-existent. 2828 * This prevents jumping behavior of the intersection points. 2829 * But it may be discussed if it is the desired behavior. 2830 */ 2831 if ( 2832 testSegment && 2833 ((!ev_sf && p[2] < 0) || (!ev_sl && p[2] > 1)) 2834 ) { 2835 return q; // break; 2836 } 2837 2838 q = new Coords(Const.COORDS_BY_USER, p[0], board); 2839 return q; // break; 2840 } 2841 cnt += 1; 2842 } 2843 } 2844 } 2845 } 2846 2847 return q; 2848 }, 2849 2850 /** 2851 * Find the n-th intersection point of two curves named red (first parameter) and blue (second parameter). 2852 * We go through each segment of the red curve and search if there is an intersection with a segment of the blue curve. 2853 * This double loop, i.e. the outer loop runs along the red curve and the inner loop runs along the blue curve, defines 2854 * the n-th intersection point. The segments are either line segments or Bezier curves of degree 3. This depends on 2855 * the property bezierDegree of the curves. 2856 * <p> 2857 * This method works also for transformed curves, since only the already 2858 * transformed points are used. 2859 * 2860 * @param {JXG.Curve} red 2861 * @param {JXG.Curve} blue 2862 * @param {Number|Function} nr 2863 */ 2864 meetCurveRedBlueSegments: function (red, blue, nr) { 2865 var i, 2866 j, 2867 n = Type.evaluate(nr), 2868 red1, 2869 red2, 2870 blue1, 2871 blue2, 2872 m, 2873 minX, 2874 maxX, 2875 iFound = 0, 2876 lenBlue = blue.numberPoints, 2877 lenRed = red.numberPoints; 2878 2879 if (lenBlue <= 1 || lenRed <= 1) { 2880 return [0, NaN, NaN]; 2881 } 2882 2883 for (i = 1; i < lenRed; i++) { 2884 red1 = red.points[i - 1].usrCoords; 2885 red2 = red.points[i].usrCoords; 2886 minX = Math.min(red1[1], red2[1]); 2887 maxX = Math.max(red1[1], red2[1]); 2888 2889 blue2 = blue.points[0].usrCoords; 2890 for (j = 1; j < lenBlue; j++) { 2891 blue1 = blue2; 2892 blue2 = blue.points[j].usrCoords; 2893 if ( 2894 Math.min(blue1[1], blue2[1]) < maxX && 2895 Math.max(blue1[1], blue2[1]) > minX 2896 ) { 2897 m = this.meetSegmentSegment(red1, red2, blue1, blue2); 2898 if ( 2899 m[1] >= 0.0 && m[2] >= 0.0 && 2900 // The two segments meet in the interior or at the start points 2901 ((m[1] < 1.0 && m[2] < 1.0) || 2902 // One of the curve is intersected in the very last point 2903 (i === lenRed - 1 && m[1] === 1.0) || 2904 (j === lenBlue - 1 && m[2] === 1.0)) 2905 ) { 2906 if (iFound === n) { 2907 return m[0]; 2908 } 2909 2910 iFound++; 2911 } 2912 } 2913 } 2914 } 2915 2916 return [0, NaN, NaN]; 2917 }, 2918 2919 /** 2920 * (Virtual) Intersection of two segments. 2921 * @param {Array} p1 First point of segment 1 using normalized homogeneous coordinates [1,x,y] 2922 * @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 2923 * @param {Array} q1 First point of segment 2 using normalized homogeneous coordinates [1,x,y] 2924 * @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 2925 * @returns {Array} [Intersection point, t, u] The first entry contains the homogeneous coordinates 2926 * of the intersection point. The second and third entry give the position of the intersection with respect 2927 * to the definiting parameters. For example, the second entry t is defined by: intersection point = p1 + t * deltaP, where 2928 * deltaP = (p2 - p1) when both parameters are coordinates, and deltaP = p2 if p2 is a point at infinity. 2929 * If the two segments are collinear, [[0,0,0], Infinity, Infinity] is returned. 2930 **/ 2931 meetSegmentSegment: function (p1, p2, q1, q2) { 2932 var t, 2933 u, 2934 i, 2935 d, 2936 li1 = Mat.crossProduct(p1, p2), 2937 li2 = Mat.crossProduct(q1, q2), 2938 c = Mat.crossProduct(li1, li2); 2939 2940 if (Math.abs(c[0]) < Mat.eps) { 2941 return [c, Infinity, Infinity]; 2942 } 2943 2944 // Normalize the intersection coordinates 2945 c[1] /= c[0]; 2946 c[2] /= c[0]; 2947 c[0] /= c[0]; 2948 2949 // Now compute in principle: 2950 // t = dist(c - p1) / dist(p2 - p1) and 2951 // u = dist(c - q1) / dist(q2 - q1) 2952 // However: the points q1, q2, p1, p2 might be ideal points - or in general - the 2953 // coordinates might be not normalized. 2954 // Note that the z-coordinates of p2 and q2 are used to determine whether it should be interpreted 2955 // as a segment coordinate or a direction. 2956 i = Math.abs(p2[1] - p2[0] * p1[1]) < Mat.eps ? 2 : 1; 2957 d = p1[i] / p1[0]; 2958 t = (c[i] - d) / (p2[0] !== 0 ? p2[i] / p2[0] - d : p2[i]); 2959 2960 i = Math.abs(q2[1] - q2[0] * q1[1]) < Mat.eps ? 2 : 1; 2961 d = q1[i] / q1[0]; 2962 u = (c[i] - d) / (q2[0] !== 0 ? q2[i] / q2[0] - d : q2[i]); 2963 2964 return [c, t, u]; 2965 }, 2966 2967 /** 2968 * Find the n-th intersection point of two pathes, usually given by polygons. Uses parts of the 2969 * Greiner-Hormann algorithm in JXG.Math.Clip. 2970 * 2971 * @param {JXG.Circle|JXG.Curve|JXG.Polygon} path1 2972 * @param {JXG.Circle|JXG.Curve|JXG.Polygon} path2 2973 * @param {Number|Function} n 2974 * @param {JXG.Board} board 2975 * 2976 * @returns {JXG.Coords} Intersection point. In case no intersection point is detected, 2977 * the ideal point [0,0,0] is returned. 2978 * 2979 */ 2980 meetPathPath: function (path1, path2, nr, board) { 2981 var S, C, len, intersections, 2982 n = Type.evaluate(nr); 2983 2984 S = JXG.Math.Clip._getPath(path1, board); 2985 len = S.length; 2986 if ( 2987 len > 0 && 2988 this.distance(S[0].coords.usrCoords, S[len - 1].coords.usrCoords, 3) < Mat.eps 2989 ) { 2990 S.pop(); 2991 } 2992 2993 C = JXG.Math.Clip._getPath(path2, board); 2994 len = C.length; 2995 if ( 2996 len > 0 && 2997 this.distance(C[0].coords.usrCoords, C[len - 1].coords.usrCoords, 3) < 2998 Mat.eps * Mat.eps 2999 ) { 3000 C.pop(); 3001 } 3002 3003 // Handle cases where at least one of the paths is empty 3004 if (nr < 0 || JXG.Math.Clip.isEmptyCase(S, C, 'intersection')) { 3005 return new Coords(Const.COORDS_BY_USER, [0, 0, 0], board); 3006 } 3007 3008 JXG.Math.Clip.makeDoublyLinkedList(S); 3009 JXG.Math.Clip.makeDoublyLinkedList(C); 3010 3011 intersections = JXG.Math.Clip.findIntersections(S, C, board)[0]; 3012 if (n < intersections.length) { 3013 return intersections[n].coords; 3014 } 3015 return new Coords(Const.COORDS_BY_USER, [0, 0, 0], board); 3016 }, 3017 3018 /** 3019 * Find the n-th intersection point between a polygon and a line. 3020 * @param {JXG.Polygon} path 3021 * @param {JXG.Line} line 3022 * @param {Number|Function} nr 3023 * @param {JXG.Board} board 3024 * @param {Boolean} alwaysIntersect If false just the segment between the two defining points of the line are tested for intersection. 3025 * 3026 * @returns {JXG.Coords} Intersection point. In case no intersection point is detected, 3027 * the ideal point [0,0,0] is returned. 3028 */ 3029 meetPolygonLine: function (path, line, nr, board, alwaysIntersect) { 3030 var i, 3031 n = Type.evaluate(nr), 3032 res, 3033 border, 3034 crds = [0, 0, 0], 3035 len = path.borders.length, 3036 intersections = []; 3037 3038 for (i = 0; i < len; i++) { 3039 border = path.borders[i]; 3040 res = this.meetSegmentSegment( 3041 border.point1.coords.usrCoords, 3042 border.point2.coords.usrCoords, 3043 line.point1.coords.usrCoords, 3044 line.point2.coords.usrCoords 3045 ); 3046 3047 if ( 3048 (!alwaysIntersect || (res[2] >= 0 && res[2] < 1)) && 3049 res[1] >= 0 && 3050 res[1] < 1 3051 ) { 3052 intersections.push(res[0]); 3053 } 3054 } 3055 3056 if (n >= 0 && n < intersections.length) { 3057 crds = intersections[n]; 3058 } 3059 return new Coords(Const.COORDS_BY_USER, crds, board); 3060 }, 3061 3062 /****************************************/ 3063 /**** BEZIER CURVE ALGORITHMS ****/ 3064 /****************************************/ 3065 3066 /** 3067 * Splits a Bezier curve segment defined by four points into 3068 * two Bezier curve segments. Dissection point is t=1/2. 3069 * @param {Array} curve Array of four coordinate arrays of length 2 defining a 3070 * Bezier curve segment, i.e. [[x0,y0], [x1,y1], [x2,y2], [x3,y3]]. 3071 * @returns {Array} Array consisting of two coordinate arrays for Bezier curves. 3072 */ 3073 _bezierSplit: function (curve) { 3074 var p0, p1, p2, p00, p22, p000; 3075 3076 p0 = [(curve[0][0] + curve[1][0]) * 0.5, (curve[0][1] + curve[1][1]) * 0.5]; 3077 p1 = [(curve[1][0] + curve[2][0]) * 0.5, (curve[1][1] + curve[2][1]) * 0.5]; 3078 p2 = [(curve[2][0] + curve[3][0]) * 0.5, (curve[2][1] + curve[3][1]) * 0.5]; 3079 3080 p00 = [(p0[0] + p1[0]) * 0.5, (p0[1] + p1[1]) * 0.5]; 3081 p22 = [(p1[0] + p2[0]) * 0.5, (p1[1] + p2[1]) * 0.5]; 3082 3083 p000 = [(p00[0] + p22[0]) * 0.5, (p00[1] + p22[1]) * 0.5]; 3084 3085 return [ 3086 [curve[0], p0, p00, p000], 3087 [p000, p22, p2, curve[3]] 3088 ]; 3089 }, 3090 3091 /** 3092 * Computes the bounding box [minX, maxY, maxX, minY] of a Bezier curve segment 3093 * from its control points. 3094 * @param {Array} curve Array of four coordinate arrays of length 2 defining a 3095 * Bezier curve segment, i.e. [[x0,y0], [x1,y1], [x2,y2], [x3,y3]]. 3096 * @returns {Array} Bounding box [minX, maxY, maxX, minY] 3097 */ 3098 _bezierBbox: function (curve) { 3099 var bb = []; 3100 3101 if (curve.length === 4) { 3102 // bezierDegree == 3 3103 bb[0] = Math.min(curve[0][0], curve[1][0], curve[2][0], curve[3][0]); // minX 3104 bb[1] = Math.max(curve[0][1], curve[1][1], curve[2][1], curve[3][1]); // maxY 3105 bb[2] = Math.max(curve[0][0], curve[1][0], curve[2][0], curve[3][0]); // maxX 3106 bb[3] = Math.min(curve[0][1], curve[1][1], curve[2][1], curve[3][1]); // minY 3107 } else { 3108 // bezierDegree == 1 3109 bb[0] = Math.min(curve[0][0], curve[1][0]); // minX 3110 bb[1] = Math.max(curve[0][1], curve[1][1]); // maxY 3111 bb[2] = Math.max(curve[0][0], curve[1][0]); // maxX 3112 bb[3] = Math.min(curve[0][1], curve[1][1]); // minY 3113 } 3114 3115 return bb; 3116 }, 3117 3118 /** 3119 * Decide if two Bezier curve segments overlap by comparing their bounding boxes. 3120 * @param {Array} bb1 Bounding box of the first Bezier curve segment 3121 * @param {Array} bb2 Bounding box of the second Bezier curve segment 3122 * @returns {Boolean} true if the bounding boxes overlap, false otherwise. 3123 */ 3124 _bezierOverlap: function (bb1, bb2) { 3125 return bb1[2] >= bb2[0] && bb1[0] <= bb2[2] && bb1[1] >= bb2[3] && bb1[3] <= bb2[1]; 3126 }, 3127 3128 /** 3129 * Append list of intersection points to a list. 3130 * @private 3131 */ 3132 _bezierListConcat: function (L, Lnew, t1, t2) { 3133 var i, 3134 t2exists = Type.exists(t2), 3135 start = 0, 3136 len = Lnew.length, 3137 le = L.length; 3138 3139 if ( 3140 le > 0 && 3141 len > 0 && 3142 ((L[le - 1][1] === 1 && Lnew[0][1] === 0) || 3143 (t2exists && L[le - 1][2] === 1 && Lnew[0][2] === 0)) 3144 ) { 3145 start = 1; 3146 } 3147 3148 for (i = start; i < len; i++) { 3149 if (t2exists) { 3150 Lnew[i][2] *= 0.5; 3151 Lnew[i][2] += t2; 3152 } 3153 3154 Lnew[i][1] *= 0.5; 3155 Lnew[i][1] += t1; 3156 3157 L.push(Lnew[i]); 3158 } 3159 }, 3160 3161 /** 3162 * Find intersections of two Bezier curve segments by recursive subdivision. 3163 * Below maxlevel determine intersections by intersection line segments. 3164 * @param {Array} red Array of four coordinate arrays of length 2 defining the first 3165 * Bezier curve segment, i.e. [[x0,y0], [x1,y1], [x2,y2], [x3,y3]]. 3166 * @param {Array} blue Array of four coordinate arrays of length 2 defining the second 3167 * Bezier curve segment, i.e. [[x0,y0], [x1,y1], [x2,y2], [x3,y3]]. 3168 * @param {Number} level Recursion level 3169 * @returns {Array} List of intersection points (up to nine). Each intersection point is an 3170 * array of length three (homogeneous coordinates) plus preimages. 3171 */ 3172 _bezierMeetSubdivision: function (red, blue, level) { 3173 var bbb, 3174 bbr, 3175 ar, 3176 b0, 3177 b1, 3178 r0, 3179 r1, 3180 m, 3181 p0, 3182 p1, 3183 q0, 3184 q1, 3185 L = [], 3186 maxLev = 5; // Maximum recursion level 3187 3188 bbr = this._bezierBbox(blue); 3189 bbb = this._bezierBbox(red); 3190 3191 if (!this._bezierOverlap(bbr, bbb)) { 3192 return []; 3193 } 3194 3195 if (level < maxLev) { 3196 ar = this._bezierSplit(red); 3197 r0 = ar[0]; 3198 r1 = ar[1]; 3199 3200 ar = this._bezierSplit(blue); 3201 b0 = ar[0]; 3202 b1 = ar[1]; 3203 3204 this._bezierListConcat( 3205 L, 3206 this._bezierMeetSubdivision(r0, b0, level + 1), 3207 0.0, 3208 0.0 3209 ); 3210 this._bezierListConcat( 3211 L, 3212 this._bezierMeetSubdivision(r0, b1, level + 1), 3213 0, 3214 0.5 3215 ); 3216 this._bezierListConcat( 3217 L, 3218 this._bezierMeetSubdivision(r1, b0, level + 1), 3219 0.5, 3220 0.0 3221 ); 3222 this._bezierListConcat( 3223 L, 3224 this._bezierMeetSubdivision(r1, b1, level + 1), 3225 0.5, 3226 0.5 3227 ); 3228 3229 return L; 3230 } 3231 3232 // Make homogeneous coordinates 3233 q0 = [1].concat(red[0]); 3234 q1 = [1].concat(red[3]); 3235 p0 = [1].concat(blue[0]); 3236 p1 = [1].concat(blue[3]); 3237 3238 m = this.meetSegmentSegment(q0, q1, p0, p1); 3239 3240 if (m[1] >= 0.0 && m[2] >= 0.0 && m[1] <= 1.0 && m[2] <= 1.0) { 3241 return [m]; 3242 } 3243 3244 return []; 3245 }, 3246 3247 /** 3248 * @param {Boolean} testSegment Test if intersection has to be inside of the segment or somewhere on the line defined by the segment 3249 */ 3250 _bezierLineMeetSubdivision: function (red, blue, level, testSegment) { 3251 var bbb, bbr, ar, 3252 r0, r1, 3253 m, 3254 p0, p1, q0, q1, 3255 L = [], 3256 maxLev = 5; // Maximum recursion level 3257 3258 bbb = this._bezierBbox(blue); 3259 bbr = this._bezierBbox(red); 3260 3261 if (testSegment && !this._bezierOverlap(bbr, bbb)) { 3262 return []; 3263 } 3264 3265 if (level < maxLev) { 3266 ar = this._bezierSplit(red); 3267 r0 = ar[0]; 3268 r1 = ar[1]; 3269 3270 this._bezierListConcat( 3271 L, 3272 this._bezierLineMeetSubdivision(r0, blue, level + 1), 3273 0.0 3274 ); 3275 this._bezierListConcat( 3276 L, 3277 this._bezierLineMeetSubdivision(r1, blue, level + 1), 3278 0.5 3279 ); 3280 3281 return L; 3282 } 3283 3284 // Make homogeneous coordinates 3285 q0 = [1].concat(red[0]); 3286 q1 = [1].concat(red[3]); 3287 p0 = [1].concat(blue[0]); 3288 p1 = [1].concat(blue[1]); 3289 3290 m = this.meetSegmentSegment(q0, q1, p0, p1); 3291 3292 if (m[1] >= 0.0 && m[1] <= 1.0) { 3293 if (!testSegment || (m[2] >= 0.0 && m[2] <= 1.0)) { 3294 return [m]; 3295 } 3296 } 3297 3298 return []; 3299 }, 3300 3301 /** 3302 * Find the nr-th intersection point of two Bezier curve segments. 3303 * @param {Array} red Array of four coordinate arrays of length 2 defining the first 3304 * Bezier curve segment, i.e. [[x0,y0], [x1,y1], [x2,y2], [x3,y3]]. 3305 * @param {Array} blue Array of four coordinate arrays of length 2 defining the second 3306 * Bezier curve segment, i.e. [[x0,y0], [x1,y1], [x2,y2], [x3,y3]]. 3307 * @param {Boolean} testSegment Test if intersection has to be inside of the segment or somewhere on the line defined by the segment 3308 * @returns {Array} Array containing the list of all intersection points as homogeneous coordinate arrays plus 3309 * preimages [x,y], t_1, t_2] of the two Bezier curve segments. 3310 * 3311 */ 3312 meetBeziersegmentBeziersegment: function (red, blue, testSegment) { 3313 var L, L2, i; 3314 3315 if (red.length === 4 && blue.length === 4) { 3316 L = this._bezierMeetSubdivision(red, blue, 0); 3317 } else { 3318 L = this._bezierLineMeetSubdivision(red, blue, 0, testSegment); 3319 } 3320 3321 L.sort(function (a, b) { 3322 return (a[1] - b[1]) * 10000000.0 + (a[2] - b[2]); 3323 }); 3324 3325 L2 = []; 3326 for (i = 0; i < L.length; i++) { 3327 // Only push entries different from their predecessor 3328 if (i === 0 || L[i][1] !== L[i - 1][1] || L[i][2] !== L[i - 1][2]) { 3329 L2.push(L[i]); 3330 } 3331 } 3332 return L2; 3333 }, 3334 3335 /** 3336 * Find the nr-th intersection point of two Bezier curves, i.e. curves with bezierDegree == 3. 3337 * @param {JXG.Curve} red Curve with bezierDegree == 3 3338 * @param {JXG.Curve} blue Curve with bezierDegree == 3 3339 * @param {Number|Function} nr The number of the intersection point which should be returned. 3340 * @returns {Array} The homogeneous coordinates of the nr-th intersection point. 3341 */ 3342 meetBezierCurveRedBlueSegments: function (red, blue, nr) { 3343 var p, i, j, k, 3344 n = Type.evaluate(nr), 3345 po, tmp, 3346 redArr, 3347 blueArr, 3348 bbr, 3349 bbb, 3350 intersections, 3351 startRed = 0, 3352 startBlue = 0, 3353 lenBlue, lenRed, 3354 L = []; 3355 3356 if (blue.numberPoints < blue.bezierDegree + 1 || red.numberPoints < red.bezierDegree + 1) { 3357 return [0, NaN, NaN]; 3358 } 3359 if (red.bezierDegree === 1 && blue.bezierDegree === 3) { 3360 tmp = red; 3361 red = blue; 3362 blue = tmp; 3363 } 3364 3365 lenBlue = blue.numberPoints - blue.bezierDegree; 3366 lenRed = red.numberPoints - red.bezierDegree; 3367 3368 // For sectors, we ignore the "legs" 3369 if (red.type === Const.OBJECT_TYPE_SECTOR) { 3370 startRed = 3; 3371 lenRed -= 3; 3372 } 3373 if (blue.type === Const.OBJECT_TYPE_SECTOR) { 3374 startBlue = 3; 3375 lenBlue -= 3; 3376 } 3377 3378 for (i = startRed; i < lenRed; i += red.bezierDegree) { 3379 p = red.points; 3380 redArr = [p[i].usrCoords.slice(1), p[i + 1].usrCoords.slice(1)]; 3381 if (red.bezierDegree === 3) { 3382 redArr[2] = p[i + 2].usrCoords.slice(1); 3383 redArr[3] = p[i + 3].usrCoords.slice(1); 3384 } 3385 3386 bbr = this._bezierBbox(redArr); 3387 3388 for (j = startBlue; j < lenBlue; j += blue.bezierDegree) { 3389 p = blue.points; 3390 blueArr = [p[j].usrCoords.slice(1), p[j + 1].usrCoords.slice(1)]; 3391 if (blue.bezierDegree === 3) { 3392 blueArr[2] = p[j + 2].usrCoords.slice(1); 3393 blueArr[3] = p[j + 3].usrCoords.slice(1); 3394 } 3395 3396 bbb = this._bezierBbox(blueArr); 3397 if (this._bezierOverlap(bbr, bbb)) { 3398 intersections = this.meetBeziersegmentBeziersegment(redArr, blueArr); 3399 if (intersections.length === 0) { 3400 continue; 3401 } 3402 for (k = 0; k < intersections.length; k++) { 3403 po = intersections[k]; 3404 if ( 3405 po[1] < -Mat.eps || 3406 po[1] > 1 + Mat.eps || 3407 po[2] < -Mat.eps || 3408 po[2] > 1 + Mat.eps 3409 ) { 3410 continue; 3411 } 3412 L.push(po); 3413 } 3414 if (L.length > n) { 3415 return L[n][0]; 3416 } 3417 } 3418 } 3419 } 3420 if (L.length > n) { 3421 return L[n][0]; 3422 } 3423 3424 return [0, NaN, NaN]; 3425 }, 3426 3427 bezierSegmentEval: function (t, curve) { 3428 var f, 3429 x, 3430 y, 3431 t1 = 1.0 - t; 3432 3433 x = 0; 3434 y = 0; 3435 3436 f = t1 * t1 * t1; 3437 x += f * curve[0][0]; 3438 y += f * curve[0][1]; 3439 3440 f = 3.0 * t * t1 * t1; 3441 x += f * curve[1][0]; 3442 y += f * curve[1][1]; 3443 3444 f = 3.0 * t * t * t1; 3445 x += f * curve[2][0]; 3446 y += f * curve[2][1]; 3447 3448 f = t * t * t; 3449 x += f * curve[3][0]; 3450 y += f * curve[3][1]; 3451 3452 return [1.0, x, y]; 3453 }, 3454 3455 /** 3456 * Generate the defining points of a 3rd degree bezier curve that approximates 3457 * a circle sector defined by three coordinate points A, B, C, each defined by an array of length three. 3458 * The coordinate arrays are given in homogeneous coordinates. 3459 * @param {Array} A First point 3460 * @param {Array} B Second point (intersection point) 3461 * @param {Array} C Third point 3462 * @param {Boolean} withLegs Flag. If true the legs to the intersection point are part of the curve. 3463 * @param {Number} sgn Wither 1 or -1. Needed for minor and major arcs. In case of doubt, use 1. 3464 */ 3465 bezierArc: function (A, B, C, withLegs, sgn) { 3466 var p1, p2, p3, p4, 3467 r, 3468 phi, beta, delta, 3469 // PI2 = Math.PI * 0.5, 3470 x = B[1], 3471 y = B[2], 3472 z = B[0], 3473 dataX = [], 3474 dataY = [], 3475 co, si, 3476 ax, ay, 3477 bx, by, 3478 k, v, d, 3479 matrix; 3480 3481 r = this.distance(B, A); 3482 3483 // x,y, z is intersection point. Normalize it. 3484 x /= z; 3485 y /= z; 3486 3487 phi = this.rad(A.slice(1), B.slice(1), C.slice(1)); 3488 if (sgn === -1) { 3489 phi = 2 * Math.PI - phi; 3490 } 3491 3492 // Always divide the arc into four Bezier arcs. 3493 // Otherwise, the position of gliders on this arc 3494 // will be wrong. 3495 delta = phi / 4; 3496 3497 3498 p1 = A; 3499 p1[1] /= p1[0]; 3500 p1[2] /= p1[0]; 3501 p1[0] /= p1[0]; 3502 3503 p4 = p1.slice(0); 3504 3505 if (withLegs) { 3506 dataX = [x, x + 0.333 * (p1[1] - x), x + 0.666 * (p1[1] - x), p1[1]]; 3507 dataY = [y, y + 0.333 * (p1[2] - y), y + 0.666 * (p1[2] - y), p1[2]]; 3508 } else { 3509 dataX = [p1[1]]; 3510 dataY = [p1[2]]; 3511 } 3512 3513 while (phi > Mat.eps) { 3514 // if (phi > PI2) { 3515 // beta = PI2; 3516 // phi -= PI2; 3517 // } else { 3518 // beta = phi; 3519 // phi = 0; 3520 // } 3521 if (phi > delta) { 3522 beta = delta; 3523 phi -= delta; 3524 } else { 3525 beta = phi; 3526 phi = 0; 3527 } 3528 3529 co = Math.cos(sgn * beta); 3530 si = Math.sin(sgn * beta); 3531 3532 matrix = [ 3533 [1, 0, 0], 3534 [x * (1 - co) + y * si, co, -si], 3535 [y * (1 - co) - x * si, si, co] 3536 ]; 3537 v = Mat.matVecMult(matrix, p1); 3538 p4 = [v[0] / v[0], v[1] / v[0], v[2] / v[0]]; 3539 3540 ax = p1[1] - x; 3541 ay = p1[2] - y; 3542 bx = p4[1] - x; 3543 by = p4[2] - y; 3544 d = Mat.hypot(ax + bx, ay + by); 3545 3546 if (Math.abs(by - ay) > Mat.eps) { 3547 k = ((((ax + bx) * (r / d - 0.5)) / (by - ay)) * 8) / 3; 3548 } else { 3549 k = ((((ay + by) * (r / d - 0.5)) / (ax - bx)) * 8) / 3; 3550 } 3551 3552 p2 = [1, p1[1] - k * ay, p1[2] + k * ax]; 3553 p3 = [1, p4[1] + k * by, p4[2] - k * bx]; 3554 3555 Type.concat(dataX, [p2[1], p3[1], p4[1]]); 3556 Type.concat(dataY, [p2[2], p3[2], p4[2]]); 3557 p1 = p4.slice(0); 3558 } 3559 3560 if (withLegs) { 3561 Type.concat(dataX, [ 3562 p4[1] + 0.333 * (x - p4[1]), 3563 p4[1] + 0.666 * (x - p4[1]), 3564 x 3565 ]); 3566 Type.concat(dataY, [ 3567 p4[2] + 0.333 * (y - p4[2]), 3568 p4[2] + 0.666 * (y - p4[2]), 3569 y 3570 ]); 3571 } 3572 3573 return [dataX, dataY]; 3574 }, 3575 3576 /****************************************/ 3577 /**** PROJECTIONS ****/ 3578 /****************************************/ 3579 3580 /** 3581 * Calculates the coordinates of the projection of a given point on a given circle. I.o.w. the 3582 * nearest one of the two intersection points of the line through the given point and the circles 3583 * center. 3584 * @param {JXG.Point|JXG.Coords} point Point to project or coords object to project. 3585 * @param {JXG.Circle} circle Circle on that the point is projected. 3586 * @param {JXG.Board} [board=point.board] Reference to the board 3587 * @returns {JXG.Coords} The coordinates of the projection of the given point on the given circle. 3588 */ 3589 projectPointToCircle: function (point, circle, board) { 3590 var dist, 3591 P, 3592 x, 3593 y, 3594 factor, 3595 M = circle.center.coords.usrCoords; 3596 3597 if (!Type.exists(board)) { 3598 board = point.board; 3599 } 3600 3601 // gave us a point 3602 if (Type.isPoint(point)) { 3603 dist = point.coords.distance(Const.COORDS_BY_USER, circle.center.coords); 3604 P = point.coords.usrCoords; 3605 // gave us coords 3606 } else { 3607 dist = point.distance(Const.COORDS_BY_USER, circle.center.coords); 3608 P = point.usrCoords; 3609 } 3610 3611 if (Math.abs(dist) < Mat.eps) { 3612 dist = Mat.eps; 3613 } 3614 3615 factor = circle.Radius() / dist; 3616 x = M[1] + factor * (P[1] - M[1]); 3617 y = M[2] + factor * (P[2] - M[2]); 3618 3619 return new Coords(Const.COORDS_BY_USER, [x, y], board); 3620 }, 3621 3622 /** 3623 * Calculates the coordinates of the orthogonal projection of a given point on a given line. I.o.w. the 3624 * intersection point of the given line and its perpendicular through the given point. 3625 * @param {JXG.Point|JXG.Coords} point Point to project. 3626 * @param {JXG.Line} line Line on that the point is projected. 3627 * @param {JXG.Board} [board=point.board|board=line.board] Reference to a board. 3628 * @returns {JXG.Coords} The coordinates of the projection of the given point on the given line. 3629 */ 3630 projectPointToLine: function (point, line, board) { 3631 var v = [0, line.stdform[1], line.stdform[2]], 3632 coords; 3633 3634 if (!Type.exists(board)) { 3635 if (Type.exists(point.coords)) { 3636 board = point.board; 3637 } else { 3638 board = line.board; 3639 } 3640 } 3641 3642 if (Type.exists(point.coords)) { 3643 coords = point.coords.usrCoords; 3644 } else { 3645 coords = point.usrCoords; 3646 } 3647 3648 v = Mat.crossProduct(v, coords); 3649 return new Coords(Const.COORDS_BY_USER, Mat.crossProduct(v, line.stdform), board); 3650 }, 3651 3652 /** 3653 * Calculates the coordinates of the orthogonal projection of a given coordinate array on a given line 3654 * segment defined by two coordinate arrays. 3655 * @param {Array} p Point to project. 3656 * @param {Array} q1 Start point of the line segment on that the point is projected. 3657 * @param {Array} q2 End point of the line segment on that the point is projected. 3658 * @returns {Array} The coordinates of the projection of the given point on the given segment 3659 * and the factor that determines the projected point as a convex combination of the 3660 * two endpoints q1 and q2 of the segment. 3661 */ 3662 projectCoordsToSegment: function (p, q1, q2) { 3663 var t, 3664 denom, 3665 s = [q2[1] - q1[1], q2[2] - q1[2]], 3666 v = [p[1] - q1[1], p[2] - q1[2]]; 3667 3668 /** 3669 * If the segment has length 0, i.e. is a point, 3670 * the projection is equal to that point. 3671 */ 3672 if (Math.abs(s[0]) < Mat.eps && Math.abs(s[1]) < Mat.eps) { 3673 return [q1, 0]; 3674 } 3675 3676 t = Mat.innerProduct(v, s); 3677 denom = Mat.innerProduct(s, s); 3678 t /= denom; 3679 3680 return [[1, t * s[0] + q1[1], t * s[1] + q1[2]], t]; 3681 }, 3682 3683 /** 3684 * Finds the coordinates of the closest point on a Bezier segment of a 3685 * {@link JXG.Curve} to a given coordinate array. 3686 * @param {Array} pos Point to project in homogeneous coordinates. 3687 * @param {JXG.Curve} curve Curve of type "plot" having Bezier degree 3. 3688 * @param {Number} start Number of the Bezier segment of the curve. 3689 * @returns {Array} The coordinates of the projection of the given point 3690 * on the given Bezier segment and the preimage of the curve which 3691 * determines the closest point. 3692 */ 3693 projectCoordsToBeziersegment: function (pos, curve, start) { 3694 var t0, 3695 /** @ignore */ 3696 minfunc = function (t) { 3697 var z = [1, curve.X(start + t), curve.Y(start + t)]; 3698 3699 z[1] -= pos[1]; 3700 z[2] -= pos[2]; 3701 3702 return z[1] * z[1] + z[2] * z[2]; 3703 }; 3704 3705 t0 = JXG.Math.Numerics.fminbr(minfunc, [0.0, 1.0]); 3706 3707 return [[1, curve.X(t0 + start), curve.Y(t0 + start)], t0]; 3708 }, 3709 3710 /** 3711 * Calculates the coordinates of the projection of a given point on a given curve. 3712 * Uses {@link JXG.Math.Geometry.projectCoordsToCurve}. 3713 * 3714 * @param {JXG.Point} point Point to project. 3715 * @param {JXG.Curve} curve Curve on that the point is projected. 3716 * @param {JXG.Board} [board=point.board] Reference to a board. 3717 * @see JXG.Math.Geometry.projectCoordsToCurve 3718 * @returns {Array} [JXG.Coords, position] The coordinates of the projection of the given 3719 * point on the given graph and the relative position on the curve (real number). 3720 */ 3721 projectPointToCurve: function (point, curve, board) { 3722 if (!Type.exists(board)) { 3723 board = point.board; 3724 } 3725 3726 var x = point.X(), 3727 y = point.Y(), 3728 t = point.position, 3729 result; 3730 3731 if (!Type.exists(t)) { 3732 t = curve.evalVisProp('curvetype') === 'functiongraph' ? x : 0.0; 3733 } 3734 result = this.projectCoordsToCurve(x, y, t, curve, board); 3735 // point.position = result[1]; 3736 3737 return result; 3738 }, 3739 3740 /** 3741 * Calculates the coordinates of the projection of a coordinates pair on a given curve. In case of 3742 * function graphs this is the 3743 * intersection point of the curve and the parallel to y-axis through the given point. 3744 * @param {Number} x coordinate to project. 3745 * @param {Number} y coordinate to project. 3746 * @param {Number} t start value for newtons method 3747 * @param {JXG.Curve} curve Curve on that the point is projected. 3748 * @param {JXG.Board} [board=curve.board] Reference to a board. 3749 * @see JXG.Math.Geometry.projectPointToCurve 3750 * @returns {JXG.Coords} Array containing the coordinates of the projection of the given point on the given curve and 3751 * the position on the curve. 3752 */ 3753 projectCoordsToCurve: function (x, y, t, curve, board) { 3754 var newCoords, newCoordsObj, 3755 i, j, mindist, dist, lbda, 3756 v, coords, d, p1, p2, res, minfunc, 3757 t_new, f_new, f_old, dy, 3758 delta, delta1, delta2, steps, 3759 minX, maxX, minX_glob, maxX_glob, 3760 infty = Number.POSITIVE_INFINITY; 3761 3762 if (!Type.exists(board)) { 3763 board = curve.board; 3764 } 3765 3766 if (curve.evalVisProp('curvetype') === 'plot') { 3767 t = 0; 3768 mindist = infty; 3769 if (curve.numberPoints === 0) { 3770 newCoords = [0, 1, 1]; 3771 } else { 3772 newCoords = [curve.Z(0), curve.X(0), curve.Y(0)]; 3773 } 3774 3775 if (curve.numberPoints > 1) { 3776 v = [1, x, y]; 3777 if (curve.bezierDegree === 3) { 3778 j = 0; 3779 } else { 3780 p1 = [curve.Z(0), curve.X(0), curve.Y(0)]; 3781 } 3782 for (i = 0; i < curve.numberPoints - 1; i++) { 3783 if (curve.bezierDegree === 3) { 3784 res = this.projectCoordsToBeziersegment(v, curve, j); 3785 } else { 3786 p2 = [curve.Z(i + 1), curve.X(i + 1), curve.Y(i + 1)]; 3787 res = this.projectCoordsToSegment(v, p1, p2); 3788 } 3789 lbda = res[1]; 3790 coords = res[0]; 3791 3792 if (0.0 <= lbda && lbda <= 1.0) { 3793 dist = this.distance(coords, v); 3794 d = i + lbda; 3795 } else if (lbda < 0.0) { 3796 coords = p1; 3797 dist = this.distance(p1, v); 3798 d = i; 3799 } else if (lbda > 1.0 && i === curve.numberPoints - 2) { 3800 coords = p2; 3801 dist = this.distance(coords, v); 3802 d = curve.numberPoints - 1; 3803 } 3804 3805 if (dist < mindist) { 3806 mindist = dist; 3807 t = d; 3808 newCoords = coords; 3809 } 3810 3811 if (curve.bezierDegree === 3) { 3812 j++; 3813 i += 2; 3814 } else { 3815 p1 = p2; 3816 } 3817 } 3818 } 3819 3820 newCoordsObj = new Coords(Const.COORDS_BY_USER, newCoords, board); 3821 } else { 3822 // 'parameter', 'polar', 'functiongraph' 3823 3824 minX_glob = curve.minX(); 3825 maxX_glob = curve.maxX(); 3826 minX = minX_glob; 3827 maxX = maxX_glob; 3828 3829 if (curve.evalVisProp('curvetype') === 'functiongraph') { 3830 // Restrict the possible position of t 3831 // to the projection of a circle to the x-axis (= t-axis) 3832 dy = Math.abs(y - curve.Y(x)); 3833 if (!isNaN(dy)) { 3834 minX = x - dy; 3835 maxX = x + dy; 3836 } 3837 } 3838 3839 /** 3840 * @ignore 3841 * Find t such that the Euclidean distance between 3842 * [x, y] and [curve.X(t), curve.Y(t)] 3843 * is minimized. 3844 */ 3845 minfunc = function (t) { 3846 var dx, dy; 3847 3848 if (t < minX_glob || t > maxX_glob) { 3849 return Infinity; 3850 } 3851 dx = x - curve.X(t); 3852 dy = y - curve.Y(t); 3853 return dx * dx + dy * dy; 3854 }; 3855 3856 // Search t which minimizes minfunc(t) 3857 // in discrete steps 3858 f_old = minfunc(t); 3859 steps = 50; 3860 delta = (maxX - minX) / steps; 3861 t_new = minX; 3862 for (i = 0; i < steps; i++) { 3863 f_new = minfunc(t_new); 3864 3865 if (f_new < f_old || f_old === Infinity || isNaN(f_old)) { 3866 t = t_new; 3867 f_old = f_new; 3868 } 3869 3870 t_new += delta; 3871 } 3872 3873 // t = Numerics.root(Numerics.D(minfunc), t); 3874 3875 // Ensure that minfunc is defined on the 3876 // enclosing interval [t-delta1, t+delta2] 3877 delta1 = delta; 3878 for (i = 0; i < 20 && isNaN(minfunc(t - delta1)); i++, delta1 *= 0.5); 3879 if (isNaN(minfunc(t - delta1))) { 3880 delta1 = 0.0; 3881 } 3882 delta2 = delta; 3883 for (i = 0; i < 20 && isNaN(minfunc(t + delta2)); i++, delta2 *= 0.5); 3884 if (isNaN(minfunc(t + delta2))) { 3885 delta2 = 0.0; 3886 } 3887 3888 // Finally, apply mathemetical optimization in the determined interval 3889 t = Numerics.fminbr(minfunc, [ 3890 Math.max(t - delta1, minX), 3891 Math.min(t + delta2, maxX) 3892 ]); 3893 3894 // Distinction between closed and open curves is not necessary. 3895 // If closed, the cyclic projection shift will work anyhow 3896 // if (Math.abs(curve.X(minX) - curve.X(maxX)) < Mat.eps && 3897 // Math.abs(curve.Y(minX) - curve.Y(maxX)) < Mat.eps) { 3898 // // Cyclically 3899 // if (t < minX) {console.log(t) 3900 // t = maxX + t - minX; 3901 // } 3902 // if (t > maxX) { 3903 // t = minX + t - maxX; 3904 // } 3905 // } else { 3906 3907 t = t < minX_glob ? minX_glob : t; 3908 t = t > maxX_glob ? maxX_glob : t; 3909 // } 3910 3911 newCoordsObj = new Coords( 3912 Const.COORDS_BY_USER, 3913 [curve.X(t), curve.Y(t)], 3914 board 3915 ); 3916 } 3917 3918 return [curve.updateTransform(newCoordsObj), t]; 3919 }, 3920 3921 /** 3922 * Calculates the coordinates of the closest orthogonal projection of a given coordinate array onto the 3923 * border of a polygon. 3924 * @param {Array} p Point to project. 3925 * @param {JXG.Polygon} pol Polygon element 3926 * @returns {Array} The coordinates of the closest projection of the given point to the border of the polygon. 3927 */ 3928 projectCoordsToPolygon: function (p, pol) { 3929 var i, 3930 len = pol.vertices.length, 3931 d_best = Infinity, 3932 d, 3933 projection, 3934 proj, 3935 bestprojection; 3936 3937 for (i = 0; i < len - 1; i++) { 3938 projection = JXG.Math.Geometry.projectCoordsToSegment( 3939 p, 3940 pol.vertices[i].coords.usrCoords, 3941 pol.vertices[i + 1].coords.usrCoords 3942 ); 3943 3944 if (0 <= projection[1] && projection[1] <= 1) { 3945 d = JXG.Math.Geometry.distance(projection[0], p, 3); 3946 proj = projection[0]; 3947 } else if (projection[1] < 0) { 3948 d = JXG.Math.Geometry.distance(pol.vertices[i].coords.usrCoords, p, 3); 3949 proj = pol.vertices[i].coords.usrCoords; 3950 } else { 3951 d = JXG.Math.Geometry.distance(pol.vertices[i + 1].coords.usrCoords, p, 3); 3952 proj = pol.vertices[i + 1].coords.usrCoords; 3953 } 3954 if (d < d_best) { 3955 bestprojection = proj.slice(0); 3956 d_best = d; 3957 } 3958 } 3959 return bestprojection; 3960 }, 3961 3962 /** 3963 * Calculates the coordinates of the projection of a given point on a given turtle. A turtle consists of 3964 * one or more curves of curveType 'plot'. Uses {@link JXG.Math.Geometry.projectPointToCurve}. 3965 * @param {JXG.Point} point Point to project. 3966 * @param {JXG.Turtle} turtle on that the point is projected. 3967 * @param {JXG.Board} [board=point.board] Reference to a board. 3968 * @returns {Array} [JXG.Coords, position] Array containing the coordinates of the projection of the given point on the turtle and 3969 * the position on the turtle. 3970 */ 3971 projectPointToTurtle: function (point, turtle, board) { 3972 var newCoords, 3973 t, 3974 x, 3975 y, 3976 i, 3977 dist, 3978 el, 3979 minEl, 3980 res, 3981 newPos, 3982 np = 0, 3983 npmin = 0, 3984 mindist = Number.POSITIVE_INFINITY, 3985 len = turtle.objects.length; 3986 3987 if (!Type.exists(board)) { 3988 board = point.board; 3989 } 3990 3991 // run through all curves of this turtle 3992 for (i = 0; i < len; i++) { 3993 el = turtle.objects[i]; 3994 3995 if (el.elementClass === Const.OBJECT_CLASS_CURVE) { 3996 res = this.projectPointToCurve(point, el); 3997 newCoords = res[0]; 3998 newPos = res[1]; 3999 dist = this.distance(newCoords.usrCoords, point.coords.usrCoords); 4000 4001 if (dist < mindist) { 4002 x = newCoords.usrCoords[1]; 4003 y = newCoords.usrCoords[2]; 4004 t = newPos; 4005 mindist = dist; 4006 minEl = el; 4007 npmin = np; 4008 } 4009 np += el.numberPoints; 4010 } 4011 } 4012 4013 newCoords = new Coords(Const.COORDS_BY_USER, [x, y], board); 4014 // point.position = t + npmin; 4015 // return minEl.updateTransform(newCoords); 4016 return [minEl.updateTransform(newCoords), t + npmin]; 4017 }, 4018 4019 /** 4020 * Trivial projection of a point to another point. 4021 * @param {JXG.Point} point Point to project (not used). 4022 * @param {JXG.Point} dest Point on that the point is projected. 4023 * @returns {JXG.Coords} The coordinates of the projection of the given point on the given circle. 4024 */ 4025 projectPointToPoint: function (point, dest) { 4026 return dest.coords; 4027 }, 4028 4029 /** 4030 * 4031 * @param {JXG.Point|JXG.Coords} point 4032 * @param {JXG.Board} [board] 4033 */ 4034 projectPointToBoard: function (point, board) { 4035 var i, 4036 l, 4037 c, 4038 brd = board || point.board, 4039 // comparison factor, point coord idx, bbox idx, 1st bbox corner x & y idx, 2nd bbox corner x & y idx 4040 config = [ 4041 // left 4042 [1, 1, 0, 0, 3, 0, 1], 4043 // top 4044 [-1, 2, 1, 0, 1, 2, 1], 4045 // right 4046 [-1, 1, 2, 2, 1, 2, 3], 4047 // bottom 4048 [1, 2, 3, 0, 3, 2, 3] 4049 ], 4050 coords = point.coords || point, 4051 bbox = brd.getBoundingBox(); 4052 4053 for (i = 0; i < 4; i++) { 4054 c = config[i]; 4055 if (c[0] * coords.usrCoords[c[1]] < c[0] * bbox[c[2]]) { 4056 // define border 4057 l = Mat.crossProduct( 4058 [1, bbox[c[3]], bbox[c[4]]], 4059 [1, bbox[c[5]], bbox[c[6]]] 4060 ); 4061 l[3] = 0; 4062 l = Mat.normalize(l); 4063 4064 // project point 4065 coords = this.projectPointToLine({ coords: coords }, { stdform: l }, brd); 4066 } 4067 } 4068 4069 return coords; 4070 }, 4071 4072 /** 4073 * Calculates the distance of a point to a line. The point and the line are given by homogeneous 4074 * coordinates. For lines this can be line.stdform. 4075 * @param {Array} point Homogeneous coordinates of a point. 4076 * @param {Array} line Homogeneous coordinates of a line ([C,A,B] where A*x+B*y+C*z=0). 4077 * @returns {Number} Distance of the point to the line. 4078 */ 4079 distPointLine: function (point, line) { 4080 var a = line[1], 4081 b = line[2], 4082 c = line[0], 4083 nom; 4084 4085 if (Math.abs(a) + Math.abs(b) < Mat.eps) { 4086 return Number.POSITIVE_INFINITY; 4087 } 4088 4089 nom = a * point[1] + b * point[2] + c; 4090 a *= a; 4091 b *= b; 4092 4093 return Math.abs(nom) / Math.sqrt(a + b); 4094 }, 4095 4096 /** 4097 * Determine the (Euclidean) distance between a point q and a line segment 4098 * defined by two points p1 and p2. In case p1 equals p2, the distance to this 4099 * point is returned. 4100 * 4101 * @param {Array} q Homogeneous coordinates of q 4102 * @param {Array} p1 Homogeneous coordinates of p1 4103 * @param {Array} p2 Homogeneous coordinates of p2 4104 * @returns {Number} Distance of q to line segment [p1, p2] 4105 */ 4106 distPointSegment: function (q, p1, p2) { 4107 var x, y, dx, dy, 4108 den, lbda, 4109 eps = Mat.eps * Mat.eps, 4110 huge = 1000000; 4111 4112 // Difference q - p1 4113 x = q[1] - p1[1]; 4114 y = q[2] - p1[2]; 4115 x = (x === Infinity) ? huge : (x === -Infinity) ? -huge : x; 4116 y = (y === Infinity) ? huge : (y === -Infinity) ? -huge : y; 4117 4118 // Difference p2 - p1 4119 dx = p2[1] - p1[1]; 4120 dy = p2[2] - p1[2]; 4121 dx = (dx === Infinity) ? huge : (dx === -Infinity) ? -huge : dx; 4122 dy = (dy === Infinity) ? huge : (dy === -Infinity) ? -huge : dy; 4123 4124 // If den==0 then p1 and p2 are identical 4125 // In this case the distance to p1 is returned 4126 den = dx * dx + dy * dy; 4127 if (den > eps) { 4128 lbda = (x * dx + y * dy) / den; 4129 if (lbda < 0.0) { 4130 lbda = 0.0; 4131 } else if (lbda > 1.0) { 4132 lbda = 1.0; 4133 } 4134 x -= lbda * dx; 4135 y -= lbda * dy; 4136 } 4137 4138 return Mat.hypot(x, y); 4139 }, 4140 4141 /* ***************************************/ 4142 /* *** 3D CALCULATIONS ****/ 4143 /* ***************************************/ 4144 4145 /** 4146 * Generate the function which computes the data of the intersection between 4147 * <ul> 4148 * <li> plane3d, plane3d, 4149 * <li> plane3d, sphere3d, 4150 * <li> sphere3d, plane3d, 4151 * <li> sphere3d, sphere3d 4152 * </ul> 4153 * 4154 * @param {JXG.GeometryElement3D} el1 Plane or sphere element 4155 * @param {JXG.GeometryElement3D} el2 Plane or sphere element 4156 * @returns {Array} of functions needed as input to create the intersecting line or circle. 4157 * 4158 */ 4159 intersectionFunction3D: function (view, el1, el2) { 4160 var func, 4161 that = this; 4162 4163 if (el1.type === Const.OBJECT_TYPE_PLANE3D) { 4164 if (el2.type === Const.OBJECT_TYPE_PLANE3D) { 4165 // func = () => view.intersectionPlanePlane(el1, el2)[i]; 4166 func = view.intersectionPlanePlane(el1, el2); 4167 } else if (el2.type === Const.OBJECT_TYPE_SPHERE3D) { 4168 func = that.meetPlaneSphere(el1, el2); 4169 } 4170 } else if (el1.type === Const.OBJECT_TYPE_SPHERE3D) { 4171 if (el2.type === Const.OBJECT_TYPE_PLANE3D) { 4172 func = that.meetPlaneSphere(el2, el1); 4173 } else if (el2.type === Const.OBJECT_TYPE_SPHERE3D) { 4174 func = that.meetSphereSphere(el1, el2); 4175 } 4176 } 4177 4178 return func; 4179 }, 4180 4181 /** 4182 * Intersecting point of three planes in 3D. The planes 4183 * are given in Hesse normal form. 4184 * 4185 * @param {Array} n1 Hesse normal form vector of plane 1 4186 * @param {Number} d1 Hesse normal form right hand side of plane 1 4187 * @param {Array} n2 Hesse normal form vector of plane 2 4188 * @param {Number} d2 Hesse normal form right hand side of plane 2 4189 * @param {Array} n3 Hesse normal form vector of plane 1 4190 * @param {Number} d3 Hesse normal form right hand side of plane 3 4191 * @returns {Array} Coordinates array of length 4 of the intersecting point 4192 */ 4193 meet3Planes: function (n1, d1, n2, d2, n3, d3) { 4194 var p = [1, 0, 0, 0], 4195 n31, n12, n23, 4196 denom, 4197 i; 4198 4199 n31 = Mat.crossProduct(n3.slice(1), n1.slice(1)); 4200 n12 = Mat.crossProduct(n1.slice(1), n2.slice(1)); 4201 n23 = Mat.crossProduct(n2.slice(1), n3.slice(1)); 4202 4203 denom = Mat.innerProduct(n1.slice(1), n23, 3); 4204 for (i = 0; i < 3; i++) { 4205 p[i + 1] = (d1 * n23[i] + d2 * n31[i] + d3 * n12[i]) / denom; 4206 } 4207 4208 return p; 4209 }, 4210 4211 /** 4212 * Direction of intersecting line of two planes in 3D. 4213 * 4214 * @param {Array} v11 First vector spanning plane 1 (homogeneous coordinates) 4215 * @param {Array} v12 Second vector spanning plane 1 (homogeneous coordinates) 4216 * @param {Array} v21 First vector spanning plane 2 (homogeneous coordinates) 4217 * @param {Array} v22 Second vector spanning plane 2 (homogeneous coordinates) 4218 * @returns {Array} Coordinates array of length 4 of the direction (homogeneous coordinates) 4219 */ 4220 meetPlanePlane: function (v11, v12, v21, v22) { 4221 var no1, 4222 no2, 4223 v, w; 4224 4225 v = v11.slice(1); 4226 w = v12.slice(1); 4227 no1 = Mat.crossProduct(v, w); 4228 4229 v = v21.slice(1); 4230 w = v22.slice(1); 4231 no2 = Mat.crossProduct(v, w); 4232 4233 w = Mat.crossProduct(no1, no2); 4234 w.unshift(0); 4235 return w; 4236 }, 4237 4238 meetPlaneSphere: function (el1, el2) { 4239 var dis = function () { 4240 return Mat.innerProduct(el1.normal, el2.center.coords, 4) - el1.d; 4241 }; 4242 4243 return [ 4244 // Center 4245 function() { 4246 return Mat.axpy(-dis(), el1.normal, el2.center.coords); 4247 }, 4248 // Normal 4249 el1.normal, 4250 // Radius 4251 function () { 4252 // Radius (returns NaN if spheres don't touch) 4253 var r = el2.Radius(), 4254 s = dis(); 4255 return Math.sqrt(r * r - s * s); 4256 } 4257 ]; 4258 }, 4259 4260 meetSphereSphere: function (el1, el2) { 4261 var skew = function () { 4262 var dist = el1.center.distance(el2.center), 4263 r1 = el1.Radius(), 4264 r2 = el2.Radius(); 4265 return (r1 - r2) * (r1 + r2) / (dist * dist); 4266 }; 4267 return [ 4268 // Center 4269 function () { 4270 var s = skew(); 4271 return [ 4272 1, 4273 0.5 * ((1 - s) * el1.center.coords[1] + (1 + s) * el2.center.coords[1]), 4274 0.5 * ((1 - s) * el1.center.coords[2] + (1 + s) * el2.center.coords[2]), 4275 0.5 * ((1 - s) * el1.center.coords[3] + (1 + s) * el2.center.coords[3]) 4276 ]; 4277 }, 4278 // Normal 4279 function() { 4280 return Stat.subtract(el2.center.coords, el1.center.coords); 4281 }, 4282 // Radius 4283 function () { 4284 // Radius (returns NaN if spheres don't touch) 4285 var dist = el1.center.distance(el2.center), 4286 r1 = el1.Radius(), 4287 r2 = el2.Radius(), 4288 s = skew(), 4289 rIxnSq = 0.5 * (r1 * r1 + r2 * r2 - 0.5 * dist * dist * (1 + s * s)); 4290 return Math.sqrt(rIxnSq); 4291 } 4292 ]; 4293 }, 4294 4295 /** 4296 * Test if parameters are inside of allowed ranges 4297 * 4298 * @param {Array} params Array of length 1 or 2 4299 * @param {Array} r_u First range 4300 * @param {Array} [r_v] Second range 4301 * @returns Boolean 4302 * @private 4303 */ 4304 _paramsOutOfRange: function(params, r_u, r_v) { 4305 return params[0] < r_u[0] || params[0] > r_u[1] || 4306 (params.length > 1 && (params[1] < r_v[0] || params[1] > r_v[1])); 4307 }, 4308 4309 /** 4310 * Given the 2D screen coordinates of a point, finds the nearest point on the given 4311 * parametric curve or surface, and returns its view-space coordinates. 4312 * @param {Array} p Homogeneous 3D coordinates for which the closest point on the curve point is searched. 4313 * @param {JXG.Curve3D|JXG.Surface3D} target Parametric curve or surface to project to. 4314 * @param {Number} n Dimension of the host element to which the coords are projected. 4315 * @param {Array} params New position of point on the target (i.e. it is a return value), 4316 * modified in place during the search, ending up at the nearest point. 4317 * Usually, point.position is supplied for params. 4318 * 4319 * @returns {Array} Array of length 4 containing the coordinates of the nearest point on the curve or surface. 4320 */ 4321 projectCoordsToParametric: function (p, target, n, params) { 4322 // The variables and parameters for the Cobyla constrained 4323 // minimization algorithm are explained in the Cobyla.js comments 4324 var rhobeg, // initial size of simplex (Cobyla) 4325 rhoend, // finial size of simplex (Cobyla) 4326 iprint = 0, // no console output (Cobyla) 4327 maxfun = 200, // call objective function at most 200 times (Cobyla) 4328 _minFunc, // Objective function for Cobyla 4329 f = Math.random() * 0.01 + 0.5, 4330 r_u, r_v, 4331 m = 2 * n; 4332 4333 // adapt simplex size to parameter range 4334 if (n === 1) { 4335 r_u = [Type.evaluate(target.range[0]), Type.evaluate(target.range[1])]; 4336 rhobeg = 0.1 * (r_u[1] - r_u[0]); 4337 4338 } else if (n === 2) { 4339 r_u = [Type.evaluate(target.range_u[0]), Type.evaluate(target.range_u[1])]; 4340 r_v = [Type.evaluate(target.range_v[0]), Type.evaluate(target.range_v[1])]; 4341 rhobeg = 0.1 * Math.min( 4342 r_u[1] - r_u[0], 4343 r_v[1] - r_v[0] 4344 ); 4345 } 4346 rhoend = rhobeg / 5e6; 4347 4348 // Minimize distance of the new position to the original position 4349 _minFunc = function (n, m, w, con) { 4350 var p_new = [ 4351 1, 4352 target.X.apply(target, w), 4353 target.Y.apply(target, w), 4354 target.Z.apply(target, w) 4355 ], 4356 xDiff = p[1] - p_new[1], 4357 yDiff = p[2] - p_new[2], 4358 zDiff = p[3] - p_new[3]; 4359 4360 if (m >= 2) { 4361 con[0] = w[0] - r_u[0]; 4362 con[1] = -w[0] + r_u[1]; 4363 } 4364 if (m >= 4) { 4365 con[2] = w[1] - r_v[0]; 4366 con[3] = -w[1] + r_v[1]; 4367 } 4368 4369 return xDiff * xDiff + yDiff * yDiff + zDiff * zDiff; 4370 }; 4371 4372 // First optimization without range constraints to give a smooth drag experience on 4373 // cyclic structures. 4374 4375 // Set the start values 4376 if (params.length === 0) { 4377 // If length > 0: take the previous position as start values for the optimization 4378 params[0] = f * (r_u[0] + r_u[1]); 4379 if (n === 2) { 4380 params[1] = f * (r_v[0] + r_v[1]); 4381 } 4382 } else { 4383 params[0] = (params[0] <= r_u[0]) ? r_u[0] + Mat.eps : params[0]; 4384 params[0] = (params[0] >= r_u[1]) ? r_u[1] - Mat.eps : params[0]; 4385 if (n === 2) { 4386 params[1] = (params[1] <= r_v[0]) ? r_v[0] + Mat.eps : params[1]; 4387 params[1] = (params[1] >= r_v[1]) ? r_v[1] - Mat.eps : params[1]; 4388 } 4389 } 4390 4391 Mat.Nlp.FindMinimum(_minFunc, n, m, params, rhobeg, rhoend, iprint, maxfun); 4392 4393 // Update p which is used subsequently in _minFunc 4394 p = [ 4395 1, 4396 target.X.apply(target, params), 4397 target.Y.apply(target, params), 4398 target.Z.apply(target, params) 4399 ]; 4400 4401 // If the optimal params are outside of the range: 4402 // Second optimization to obey the range constraints 4403 4404 if (this._paramsOutOfRange(params, r_u, r_v)) { 4405 // Set the start values again 4406 // params[0] = f * (r_u[0] + r_u[1]); 4407 // if (n === 2) { 4408 // params[1] = f * (r_v[0] + r_v[1]); 4409 // } 4410 params[0] = (params[0] <= r_u[0]) ? r_u[0] + Mat.eps : params[0]; 4411 params[0] = (params[0] >= r_u[1]) ? r_u[1] - Mat.eps : params[0]; 4412 if (n === 2) { 4413 params[1] = (params[1] <= r_v[0]) ? r_v[0] + Mat.eps : params[1]; 4414 params[1] = (params[1] >= r_v[1]) ? r_v[1] - Mat.eps : params[1]; 4415 } 4416 Mat.Nlp.FindMinimum(_minFunc, n, m, params, rhobeg, rhoend, iprint, maxfun); 4417 } 4418 4419 return [1, 4420 target.X.apply(target, params), 4421 target.Y.apply(target, params), 4422 target.Z.apply(target, params) 4423 ]; 4424 }, 4425 4426 /** 4427 * Given a the screen coordinates of a point, finds the point on the 4428 * given parametric curve or surface which is nearest in screen space, 4429 * and returns its view-space coordinates. 4430 * @param {Array} pScr Screen coordinates to project. 4431 * @param {JXG.Plane3D|JXG.Curve3D|JXG.Surface3D} target Plane, parametric curve or surface to project to. 4432 * @param {Array} params Parameters of point on the target, initially specifying the starting point of 4433 * the search. The parameters are modified in place during the search, ending up at the nearest point. 4434 * @returns {Array} Array of length 4 containing the coordinates of the nearest point on the curve or surface. 4435 */ 4436 projectScreenCoordsToParametric: function (pScr, target, params, cyclic) { 4437 // The variables and parameters for the Cobyla constrained 4438 // minimization algorithm are explained in the Cobyla.js comments 4439 var rhobeg, // initial size of simplex (Cobyla) 4440 rhoend, // finial size of simplex (Cobyla) 4441 iprint = 0, // no console output (Cobyla) 4442 maxfun = 200, // call objective function at most 200 times (Cobyla) 4443 dim = params.length, 4444 r_u, r_v, 4445 _minFunc; // objective function (Cobyla) 4446 4447 // Adapt simplex size to parameter range 4448 if (dim === 1) { 4449 r_u = [Type.evaluate(target.range[0]), Type.evaluate(target.range[1])]; 4450 rhobeg = 0.1 * (r_u[1] - r_u[0]); 4451 } else if (dim === 2) { 4452 r_u = [Type.evaluate(target.range_u[0]), Type.evaluate(target.range_u[1])]; 4453 r_v = [Type.evaluate(target.range_v[0]), Type.evaluate(target.range_v[1])]; 4454 4455 rhobeg = 0.1 * Math.min( 4456 r_u[1] - r_u[0], 4457 r_v[1] - r_v[0] 4458 ); 4459 } 4460 4461 rhoend = rhobeg / 5e6; 4462 4463 // Minimize screen distance to cursor 4464 _minFunc = function (n, m, w, con) { 4465 var c3d = [ 4466 1, 4467 target.X.apply(target, w), 4468 target.Y.apply(target, w), 4469 target.Z.apply(target, w) 4470 ], 4471 c2d = target.view.project3DTo2D(c3d), 4472 xDiff = pScr[0] - c2d[1], 4473 yDiff = pScr[1] - c2d[2]; 4474 4475 if (n === 1) { 4476 con[0] = w[0] - r_u[0]; 4477 con[1] = -w[0] + r_u[1]; 4478 } else if (n === 2) { 4479 con[0] = w[0] - r_u[0]; 4480 con[1] = -w[0] + r_u[1]; 4481 con[2] = w[1] - r_v[0]; 4482 con[3] = -w[1] + r_v[1]; 4483 } 4484 4485 return xDiff * xDiff + yDiff * yDiff; 4486 }; 4487 4488 if (cyclic) { 4489 // Cyclic 4490 Mat.Nlp.FindMinimum(_minFunc, dim, 0 /*2 * dim*/, params, rhobeg, rhoend, iprint, maxfun); 4491 params[0] = (params[0] + 20 * r_u[1]) % (r_u[1] - r_u[0]); 4492 if (dim === 2) { 4493 params[1] = (params[1] + 20 * r_v[1]) % (r_v[1] - r_v[0]); 4494 } 4495 } else { 4496 Mat.Nlp.FindMinimum(_minFunc, dim, 2 * dim, params, rhobeg, rhoend, iprint, maxfun); 4497 } 4498 4499 return [1, target.X.apply(target, params), target.Y.apply(target, params), target.Z.apply(target, params)]; 4500 }, 4501 4502 project3DTo3DPlane: function (point, normal, foot) { 4503 // TODO: homogeneous 3D coordinates 4504 var sol = [0, 0, 0], 4505 le, 4506 d1, 4507 d2, 4508 lbda; 4509 4510 foot = foot || [0, 0, 0]; 4511 4512 le = Mat.norm(normal); 4513 d1 = Mat.innerProduct(point, normal, 3); 4514 d2 = Mat.innerProduct(foot, normal, 3); 4515 // (point - lbda * normal / le) * normal / le == foot * normal / le 4516 // => (point * normal - foot * normal) == lbda * le 4517 lbda = (d1 - d2) / le; 4518 sol = Mat.axpy(-lbda, normal, point); 4519 4520 return sol; 4521 }, 4522 4523 getPlaneBounds: function (v1, v2, q, s, e) { 4524 var s1, s2, e1, e2, mat, rhs, sol; 4525 4526 if (v1[2] + v2[0] !== 0) { 4527 mat = [ 4528 [v1[0], v2[0]], 4529 [v1[1], v2[1]] 4530 ]; 4531 rhs = [s - q[0], s - q[1]]; 4532 4533 sol = Numerics.Gauss(mat, rhs); 4534 s1 = sol[0]; 4535 s2 = sol[1]; 4536 4537 rhs = [e - q[0], e - q[1]]; 4538 sol = Numerics.Gauss(mat, rhs); 4539 e1 = sol[0]; 4540 e2 = sol[1]; 4541 return [s1, e1, s2, e2]; 4542 } 4543 return null; 4544 }, 4545 4546 /* ***************************************/ 4547 /* *** Various ****/ 4548 /* ***************************************/ 4549 4550 /** 4551 * Helper function to create curve which displays a Reuleaux polygons. 4552 * @param {Array} points Array of points which should be the vertices of the Reuleaux polygon. Typically, 4553 * these point list is the array vertices of a regular polygon. 4554 * @param {Number} nr Number of vertices 4555 * @returns {Array} An array containing the two functions defining the Reuleaux polygon and the two values 4556 * for the start and the end of the paramtric curve. array may be used as parent array of a 4557 * {@link JXG.Curve}. 4558 * 4559 * @example 4560 * var A = brd.create('point',[-2,-2]); 4561 * var B = brd.create('point',[0,1]); 4562 * var pol = brd.create('regularpolygon',[A,B,3], {withLines:false, fillColor:'none', highlightFillColor:'none', fillOpacity:0.0}); 4563 * var reuleauxTriangle = brd.create('curve', JXG.Math.Geometry.reuleauxPolygon(pol.vertices, 3), 4564 * {strokeWidth:6, strokeColor:'#d66d55', fillColor:'#ad5544', highlightFillColor:'#ad5544'}); 4565 * 4566 * </pre><div class="jxgbox" id="JXG2543a843-46a9-4372-abc1-94d9ad2db7ac" style="width: 300px; height: 300px;"></div> 4567 * <script type="text/javascript"> 4568 * var brd = JXG.JSXGraph.initBoard('JXG2543a843-46a9-4372-abc1-94d9ad2db7ac', {boundingbox: [-5, 5, 5, -5], axis: true, showcopyright:false, shownavigation: false}); 4569 * var A = brd.create('point',[-2,-2]); 4570 * var B = brd.create('point',[0,1]); 4571 * var pol = brd.create('regularpolygon',[A,B,3], {withLines:false, fillColor:'none', highlightFillColor:'none', fillOpacity:0.0}); 4572 * var reuleauxTriangle = brd.create('curve', JXG.Math.Geometry.reuleauxPolygon(pol.vertices, 3), 4573 * {strokeWidth:6, strokeColor:'#d66d55', fillColor:'#ad5544', highlightFillColor:'#ad5544'}); 4574 * </script><pre> 4575 */ 4576 reuleauxPolygon: function (points, nr) { 4577 var beta, 4578 pi2 = Math.PI * 2, 4579 pi2_n = pi2 / nr, 4580 diag = (nr - 1) / 2, 4581 d = 0, 4582 makeFct = function (which, trig) { 4583 return function (t, suspendUpdate) { 4584 var t1 = ((t % pi2) + pi2) % pi2, 4585 j = Math.floor(t1 / pi2_n) % nr; 4586 4587 if (!suspendUpdate) { 4588 d = points[0].Dist(points[diag]); 4589 beta = Mat.Geometry.rad( 4590 [points[0].X() + 1, points[0].Y()], 4591 points[0], 4592 points[diag % nr] 4593 ); 4594 } 4595 4596 if (isNaN(j)) { 4597 return j; 4598 } 4599 4600 t1 = t1 * 0.5 + j * pi2_n * 0.5 + beta; 4601 4602 return points[j][which]() + d * Math[trig](t1); 4603 }; 4604 }; 4605 4606 return [makeFct("X", 'cos'), makeFct("Y", 'sin'), 0, pi2]; 4607 } 4608 4609 } 4610 ); 4611 4612 export default Mat.Geometry; 4613