See <https://en.wikipedia.org/wiki/MacAdam_ellipse>,
    <https://doi.org/10.1364%2FJOSA.32.000247>.
779772
def plot_macadam(
    ellipse_scaling=10,
    plot_filter_positions=False,
    plot_standard_deviations=False,
    plot_rgb_triangle=True,
    plot_mesh=True,
    n=1,
    xy_to_2d=lambda xy: xy,
    axes_labels=("x", "y"),
):
    
    dir_path = os.path.dirname(os.path.realpath(__file__))
    with open(os.path.join(dir_path, "data/macadam1942/table3.yaml")) as f:
        data = yaml.safe_load(f)

    # if plot_filter_positions:
    #     with open(os.path.join(dir_path, 'data/macadam1942/table1.yaml')) as f:
    #         filters_xyz = yaml.safe_load(f)
    #     filters_xyz = {
    #         key: 100 * numpy.array(value) for key, value in filters_xyz.items()
    #         }
    #     for key, xyz in filters_xyz.items():
    #         x, y = xyz100_to_2d(xyz)
    #         plt.plot(x, y, 'xk')
    #         ax.annotate(key, (x, y))

    # collect the ellipse centers and offsets
    centers = []
    offsets = []
    for datak in data:
        # collect ellipse points
        _, _, _, _, delta_y_delta_x, delta_s = numpy.array(datak["data"]).T

        offset = (
            numpy.array([numpy.ones(delta_y_delta_x.shape[0]), delta_y_delta_x])
            / numpy.sqrt(1 + delta_y_delta_x ** 2)
            * delta_s
        )

        if offset.shape[1] < 2:
            continue

        centers.append([datak["x"], datak["y"]])
        offsets.append(numpy.column_stack([+offset, -offset]))

    centers = numpy.array(centers)

    _plot_ellipse_data(
        centers,
        offsets,
        ellipse_scaling=ellipse_scaling,
        xy_to_2d=xy_to_2d,
        plot_mesh=plot_mesh,
        n=n,
        plot_rgb_triangle=plot_rgb_triangle,
    )
    return
330413
function circleOverlap(r1, r2, d) {
        // no overlap
        if (d >= r1 + r2) {
            return 0;
        }

        // completely overlapped
        if (d <= Math.abs(r1 - r2)) {
            return Math.PI * Math.min(r1, r2) * Math.min(r1, r2);
        }

        var w1 = r1 - (d * d - r2 * r2 + r1 * r1) / (2 * d),
            w2 = r2 - (d * d - r1 * r1 + r2 * r2) / (2 * d);
        return circleArea(r1, w1) + circleArea(r2, w2);
    }
221036
def clip_ellipsis(self, x, y, radius, ratio=1.0):
        """"""

        self.gsave()
        self.newpath()
        self.moveto(xscale(x) + nscale(radius), yscale(y))
        self.path_arc(xscale(x), yscale(y), nscale(radius), ratio, 0, 360)
        self.closepath()
        self.__clip_stack.append(self.__clip_box)
        self.clip_sub()
170543
function (x, y, w, h, wrapper) {

			if (wrapper) {
				w = h = 2 * wrapper.r;
			}

			// Center correction, #1682
			if (wrapper && wrapper.isCircle) {
				x -= w / 2;
				y -= h / 2;
			}

			// Return the path
			return [
				'wa', // clockwisearcto
				x, // left
				y, // top
				x + w, // right
				y + h, // bottom
				x + w, // start x
				y + h / 2,     // start y
				x + w, // end x
				y + h / 2,     // end y
				//'x', // finish path
				'e' // close
			];
		}
470102
private void  ellipse(Attributes attributes) throws SVGParseException
   {
      debug("<ellipse>");

      if (currentElement == null)
         throw new SVGParseException("Invalid document. Root element must be <svg>");
      SVG.Ellipse  obj = new SVG.Ellipse();
      obj.document = svgDocument;
      obj.parent = currentElement;
      parseAttributesCore(obj, attributes);
      parseAttributesStyle(obj, attributes);
      parseAttributesTransform(obj, attributes);
      parseAttributesConditional(obj, attributes);
      parseAttributesEllipse(obj, attributes);
      currentElement.addChild(obj);     
   }
339175
def add_ellipsoid(self, x0, radii, lcar=None, with_volume=True, holes=None):
        
        if holes is None:
            holes = []

        if holes:
            assert with_volume

        # Add points.
        p = [
            self.add_point(x0, lcar=lcar),
            self.add_point([x0[0] + radii[0], x0[1], x0[2]], lcar=lcar),
            self.add_point([x0[0], x0[1] + radii[1], x0[2]], lcar=lcar),
            self.add_point([x0[0], x0[1], x0[2] + radii[2]], lcar=lcar),
            self.add_point([x0[0] - radii[0], x0[1], x0[2]], lcar=lcar),
            self.add_point([x0[0], x0[1] - radii[1], x0[2]], lcar=lcar),
            self.add_point([x0[0], x0[1], x0[2] - radii[2]], lcar=lcar),
        ]
        # Add skeleton.
        # Alternative for circles:
        # `self.add_circle_arc(a, b, c)`
        c = [
            self.add_ellipse_arc(p[1], p[0], p[6], p[6]),
            self.add_ellipse_arc(p[6], p[0], p[4], p[4]),
            self.add_ellipse_arc(p[4], p[0], p[3], p[3]),
            self.add_ellipse_arc(p[3], p[0], p[1], p[1]),
            self.add_ellipse_arc(p[1], p[0], p[2], p[2]),
            self.add_ellipse_arc(p[2], p[0], p[4], p[4]),
            self.add_ellipse_arc(p[4], p[0], p[5], p[5]),
            self.add_ellipse_arc(p[5], p[0], p[1], p[1]),
            self.add_ellipse_arc(p[6], p[0], p[2], p[2]),
            self.add_ellipse_arc(p[2], p[0], p[3], p[3]),
            self.add_ellipse_arc(p[3], p[0], p[5], p[5]),
            self.add_ellipse_arc(p[5], p[0], p[6], p[6]),
        ]
        # Add surfaces (1/8th of the ball surface).
        ll = [
            # one half
            self.add_line_loop([c[4], c[9], c[3]]),
            self.add_line_loop([c[8], -c[4], c[0]]),
            self.add_line_loop([-c[9], c[5], c[2]]),
            self.add_line_loop([-c[5], -c[8], c[1]]),
            # the other half
            self.add_line_loop([c[7], -c[3], c[10]]),
            self.add_line_loop([c[11], -c[7], -c[0]]),
            self.add_line_loop([-c[10], -c[2], c[6]]),
            self.add_line_loop([-c[1], -c[6], -c[11]]),
        ]
        # Create a surface for each line loop.
        s = [self.add_surface(l) for l in ll]

        # Combine the surfaces to avoid seams
        if self._gmsh_major() == 3:
            s = [self.add_compound_surface(s[:4]), self.add_compound_surface(s[4:])]
        else:
            assert self._gmsh_major() == 4
            # <https://gitlab.onelab.info/gmsh/gmsh/issues/507>
            self.add_raw_code(
                "Compound Surface{{{}}};".format(",".join([surf.id for surf in s[:4]]))
            )
            self.add_raw_code(
                "Compound Surface{{{}}};".format(",".join([surf.id for surf in s[4:]]))
            )

        # Create the surface loop.
        surface_loop = self.add_surface_loop(s)
        # if holes:
        #     # Create an array of surface loops; the first entry is the outer
        #     # surface loop, the following ones are holes.
        #     surface_loop = self.add_array([surface_loop] + holes)
        # Create volume.
        volume = self.add_volume(surface_loop, holes) if with_volume else None

        class Ellipsoid(object):
            dimension = 3

            def __init__(self, x0, radii, surface_loop, volume, lcar=None):
                self.x0 = x0
                self.lcar = lcar
                self.radii = radii
                self.surface_loop = surface_loop
                self.volume = volume
                return

        return Ellipsoid(x0, radii, surface_loop, volume, lcar=lcar)
912299
function(obj) {
                if (obj.normal || (obj.start && obj.end)) { return obj.isParallelTo(this); }
                var theta = this.direction.angleFrom(obj.direction);
                return (Math.abs(theta) <= sylvester.precision || Math.abs(theta - Math.PI) <= sylvester.precision);
              }
137437
private void guessVerticalMetricsFromGlyphBBox() {
        // Approximate capHeight from height of "H"
        // It's most unlikely that a font misses the PCLT table
        // This also assumes that postscriptnames exists ("H")
        // Should look it up int the cmap (that wouldn't help
        // for charsets without H anyway...)
        // Same for xHeight with the letter "x"
        int localCapHeight = 0;
        int localXHeight = 0;
        int localAscender = 0;
        int localDescender = 0;
        for (int i = 0; i < mtxTab.length; i++) {
            if ("H".equals(mtxTab[i].getName())) {
                localCapHeight = mtxTab[i].getBoundingBox()[3];
            } else if ("x".equals(mtxTab[i].getName())) {
                localXHeight = mtxTab[i].getBoundingBox()[3];
            } else if ("d".equals(mtxTab[i].getName())) {
                localAscender = mtxTab[i].getBoundingBox()[3];
            } else if ("p".equals(mtxTab[i].getName())) {
                localDescender = mtxTab[i].getBoundingBox()[1];
            } else {
                // OpenType Fonts with a version 3.0 "post" table don't have glyph names.
                // Use Unicode indices instead.
                List unicodeIndex = mtxTab[i].getUnicodeIndex();
                if (unicodeIndex.size() > 0) {
                    //Only the first index is used
                    char ch = (char)((Integer)unicodeIndex.get(0)).intValue();
                    if (ch == 'H') {
                        localCapHeight = mtxTab[i].getBoundingBox()[3];
                    } else if (ch == 'x') {
                        localXHeight = mtxTab[i].getBoundingBox()[3];
                    } else if (ch == 'd') {
                        localAscender = mtxTab[i].getBoundingBox()[3];
                    } else if (ch == 'p') {
                        localDescender = mtxTab[i].getBoundingBox()[1];
                    }
                }
            }
        }
        log.debug("Ascender from glyph 'd': " + localAscender 
                + " " + (localAscender));
        log.debug("Descender from glyph 'p': " + localDescender
                + " " + (localDescender));
        if (ascender - descender > upem) {
            log.debug("Replacing specified ascender/descender with derived values to get values"
                    + " which fit in the em box.");
            ascender = localAscender;
            descender = localDescender;
        }

        log.debug("xHeight from glyph 'x': " + localXHeight
                + " " + (localXHeight));
        log.debug("CapHeight from glyph 'H': " + localCapHeight
                + " " + (localCapHeight));
        if (capHeight == 0) {
            capHeight = localCapHeight;
            if (localCapHeight == 0) {
                log.warn("capHeight value could not be determined."
                        + " The font may not work as expected.");
            }
        }
        if (xHeight == 0) {
            xHeight = localXHeight;
            if (xHeight == 0) {
                log.warn("xHeight value could not be determined."
                        + " The font may not work as expected.");
            }
        }
    }
1360922
private boolean addTangents(Grid grid, int rowA, int colA, int rowB, int colB) {
		EllipseRotated_F64 a = grid.get(rowA,colA);
		EllipseRotated_F64 b = grid.get(rowB,colB);

		if( !tangentFinder.process(a,b, A0, A1, A2, A3, B0, B1, B2, B3) ) {
			return false;
		}
		Tangents ta = tangents.get(grid.getIndexOfRegEllipse(rowA,colA));
		Tangents tb = tangents.get(grid.getIndexOfRegEllipse(rowB,colB));

		// Which point is 0 or 3 is not defined and can swap arbitrarily.  To fix this problem
		// 0 will be defined as on the 'positive side' of the line connecting the ellipse centers

		double slopeX = b.center.x-a.center.x;
		double slopeY = b.center.y-a.center.y;

		double dx0 = A0.x-a.center.x;
		double dy0 = A0.y-a.center.y;

		double z = slopeX*dy0 - slopeY*dx0;
		if( z < 0 == (rowA == rowB)) {
			Point2D_F64 tmp = A0; A0 = A3; A3 = tmp;
			tmp = B0; B0 = B3; B3 = tmp;
		}

		// add tangent points from the two lines which do not cross the center line
		if( rowA == rowB ) {
			ta.t[ta.countT++].set(A0);
			ta.b[ta.countB++].set(A3);
			tb.t[tb.countT++].set(B0);
			tb.b[tb.countB++].set(B3);
		} else {
			ta.r[ta.countL++].set(A0);
			ta.l[ta.countR++].set(A3);
			tb.r[tb.countL++].set(B0);
			tb.l[tb.countR++].set(B3);
		}
		return true;
	}
179139
function parseTMXShapes() {
         var i = 0;
         var shapes = []; // add an ellipse shape

         if (this.isEllipse === true) {
           // ellipse coordinates are the center position, so set default to the corresonding radius
           shapes.push(new me.Ellipse(this.width / 2, this.height / 2, this.width, this.height).rotate(this.rotation));
         } // add a polygon
         else if (this.isPolygon === true) {
             shapes.push(new me.Polygon(0, 0, this.points).rotate(this.rotation));
           } // add a polyline
           else if (this.isPolyLine === true) {
               var p = this.points;
               var p1, p2;
               var segments = p.length - 1;

               for (i = 0; i < segments; i++) {
                 // clone the value before, as [i + 1]
                 // is reused later by the next segment
                 p1 = new me.Vector2d(p[i].x, p[i].y);
                 p2 = new me.Vector2d(p[i + 1].x, p[i + 1].y);

                 if (this.rotation !== 0) {
                   p1 = p1.rotate(this.rotation);
                   p2 = p2.rotate(this.rotation);
                 }

                 shapes.push(new me.Line(0, 0, [p1, p2]));
               }
             } // it's a rectangle, returns a polygon object anyway
             else {
                 shapes.push(new me.Polygon(0, 0, [new me.Vector2d(), new me.Vector2d(this.width, 0), new me.Vector2d(this.width, this.height), new me.Vector2d(0, this.height)]).rotate(this.rotation));
               } // Apply isometric projection


         if (this.orientation === "isometric") {
           for (i = 0; i < shapes.length; i++) {
             shapes[i].toIso();
           }
         }

         return shapes;
       }
39
function() {
        return {
            theta: Math.atan2(this.z, this.x),
            phi: Math.asin(this.y / this.length())
        };
    }
1232865
function fresh_seg_or_upward_cusp() {
			// trCurrent has at most 1 upper neighbor
			//	and should also have at least 1, since the high-point trapezoid
			//	has been split off another one, which is now above
			var trUpper = trCurrent.uL || trCurrent.uR;

			// trNewLeft and trNewRight CANNOT have been extended from above
			if ( trUpper.dL && trUpper.dR ) {
				// upward cusp: top forms a triangle

				// ATTENTION: the decision whether trNewLeft or trNewRight is the
				//	triangle trapezoid formed by the two segments has already been taken
				//	when selecting trCurrent as the left or right lower neighbor to trUpper !!

				if ( trCurrent == trUpper.dL ) {
					//	*** Case: FUC_UC_LEFT; prev: ----
					// console.log( "fresh_seg_or_upward_cusp: upward cusp, new seg to the left!" );
					//		  upper
					//   -------*-------
					//		   + \
					//	  NL  +   \
					//		 +	NR \
					//		+		\
					trNewRight.uL	= null;			// setAbove; trNewRight.uR, trNewLeft unchanged
					trUpper.dL		= trNewLeft;	// setBelow; dR: unchanged, NEVER null
				} else {
					//	*** Case: FUC_UC_RIGHT; prev: ----
					// console.log( "fresh_seg_or_upward_cusp: upward cusp, new seg from the right!" );
					//		  upper
					//   -------*-------
					//		   / +
					//		  /   +	 NR
					//		 /	NL +
					//		/		+
					trNewLeft.uR	= null;			// setAbove; trNewLeft.uL, trNewRight unchanged
					trUpper.dR		= trNewRight;	// setBelow; dL: unchanged, NEVER null
				}
			} else {
				//	*** Case: FUC_FS; prev: "splitOffLower"
				// console.log( "fresh_seg_or_upward_cusp: fresh segment, high adjacent segment still missing" );
				//		  upper
				//   -------*-------
				//		   +
				//	  NL  +
				//		 +	NR
				//		+
				trNewRight.uL = null;			// setAbove; trNewLeft unchanged, set by "splitOffLower"
				trNewRight.uR = trUpper;
				trUpper.dR = trNewRight;		// setBelow; trUpper.dL unchanged, set by "splitOffLower"
			}
 		}
30378
function sizeofSimple(glyf, glyfSupport, hinting) {

            if (!glyf.contours || 0 === glyf.contours.length) {
                return 0;
            }

            // fixed header + endPtsOfContours
            var result = 12
                + glyf.contours.length * 2
                + glyfSupport.flags.length;

            glyfSupport.xCoord.forEach(function (x) {
                result += 0 <= x && x <= 0xFF ? 1 : 2;
            });

            glyfSupport.yCoord.forEach(function (y) {
                result += 0 <= y && y <= 0xFF ? 1 : 2;
            });

            return result + (hinting && glyf.instructions ? glyf.instructions.length : 0);
        }
1199931
def create_ellipse(width,height,angle):
    """"""
    angle = angle / 180.0 * np.pi
    thetas = np.linspace(0,2*np.pi,200)
    a = width / 2.0
    b = height / 2.0

    x = a*np.cos(thetas)*np.cos(angle) - b*np.sin(thetas)*np.sin(angle)
    y = a*np.cos(thetas)*np.sin(angle) + b*np.sin(thetas)*np.cos(angle)
    z = np.zeros(thetas.shape)
    return np.vstack((x,y,z)).T
922125
function toLatLng()
	{
		$wgs84 = new ReferenceEllipsoid(ReferenceEllipsoid::WGS84_MAJ, ReferenceEllipsoid::WGS84_MIN);

		$UTM_F0	 = 0.9996;
		$a = $wgs84->maj;
		$eSquared = $wgs84->ecc;
		$ePrimeSquared = $eSquared / (1.0 - $eSquared);
		$e1 = (1 - sqrt(1 - $eSquared)) / (1 + sqrt(1 - $eSquared));

		$x = $this->easting - 500000.0;;
		$y = $this->northing;

		$zoneNumber = $this->lngZone;
		$zoneLetter = $this->latZone;

		$longitudeOrigin = ($zoneNumber - 1.0) * 6.0 - 180.0 + 3.0;

		// Correct y for southern hemisphere
		if ((ord($zoneLetter) - ord("N")) < 0) {
			$y -= 10000000.0;
		}

		$m = $y / $UTM_F0;

		$mu = $m
			/ ($a
			* (1.0
				- $eSquared / 4.0
				- 3.0 * $eSquared * $eSquared / 64.0
				- 5.0
				* pow($eSquared, 3.0)
				/ 256.0));

		$phi1Rad = $mu
			+ (3.0 * $e1 / 2.0 - 27.0 * pow($e1, 3.0) / 32.0) * sin(2.0 * $mu)
			+ (21.0 * $e1 * $e1 / 16.0 - 55.0 * pow($e1, 4.0) / 32.0)
			* sin(4.0 * $mu)
			+ (151.0 * pow($e1, 3.0) / 96.0) * sin(6.0 * $mu);

		$n = $a / sqrt(1.0 - $eSquared * sin($phi1Rad) * sin($phi1Rad));

		$t = tan($phi1Rad) * tan($phi1Rad);
		$c = $ePrimeSquared * cos($phi1Rad) * cos($phi1Rad);

		$r = $a
			* (1.0 - $eSquared)
			/ pow(
			1.0 - $eSquared * sin($phi1Rad) * sin($phi1Rad),
			1.5);

		$d = $x / ($n * $UTM_F0);

		$latitude = ($phi1Rad
			- ($n * tan($phi1Rad) / $r)
			* ($d * $d / 2.0
				- (5.0
				+ (3.0 * $t)
				+ (10.0 * $c)
				- (4.0 * $c * $c)
				- (9.0 * $ePrimeSquared))
				* pow($d, 4.0)
				/ 24.0
				+ (61.0
				+ (90.0 * $t)
				+ (298.0 * $c)
				+ (45.0 * $t * $t)
				- (252.0 * $ePrimeSquared)
				- (3.0 * $c * $c))
				* pow($d, 6.0)
				/ 720.0)) * (180.0 / pi());

		$longitude = $longitudeOrigin + (($d
			- (1.0 + 2.0 * $t + $c) * pow($d, 3.0) / 6.0
			+ (5.0
			- (2.0 * $c)
			+ (28.0 * $t)
			- (3.0 * $c * $c)
			+ (8.0 * $ePrimeSquared)
			+ (24.0 * $t * $t))
			* pow($d, 5.0)
			/ 120.0)
			/ cos($phi1Rad)) * (180.0 / pi());

		return new LatLng($latitude, $longitude);
	}
171171
def circle(self, clear_screen=True, x=50, y=50, radius=40, fill_color='black', outline_color='black'):
        

        if clear_screen:
            self.clear()

        x1 = x - radius
        y1 = y - radius
        x2 = x + radius
        y2 = y + radius

        return self.draw.ellipse((x1, y1, x2, y2), fill=fill_color, outline=outline_color)
779792
function(elem, other_elems) {
                    // Pairwise overlap detection
                    var overlapping = function(a, b) {
                        if (Math.abs(2.0*a.offsetLeft + a.offsetWidth - 2.0*b.offsetLeft - b.offsetWidth) < a.offsetWidth + b.offsetWidth) {
                            if (Math.abs(2.0*a.offsetTop + a.offsetHeight - 2.0*b.offsetTop - b.offsetHeight) < a.offsetHeight + b.offsetHeight) {
                                return true;
                            }
                        }
                        return false;
                    };
                    var i = 0;
                    // Check elements for overlap one by one, stop and return false as soon as an overlap is found
                    for(i = 0; i < other_elems.length; i++) {
                        if (overlapping(elem, other_elems[i])) {
                            return true;
                        }
                    }
                    return false;
                }