Create 3D Delaunay triangulation of input points.
967405
def delaunay3D(dataset, alpha=0, tol=None, boundary=True):
        """"""
        deln = vtk.vtkDelaunay3D()
        deln.SetInputData(dataset)
        deln.SetAlpha(alpha)
        if tol:
            deln.SetTolerance(tol)
        deln.SetBoundingTriangulation(boundary)
        deln.Update()
        return deln.GetOutput()
589237
private void init(DelaunayTriangulation delaunay, int xCellCount, int yCellCount, BoundingBox region) {
        indexDelaunay = delaunay;
        indexRegion = region;
        x_size = region.getWidth() / yCellCount;
        y_size = region.getHeight() / xCellCount;

        // The grid will hold a trinagle for each cell, so a point (x,y) will lie
        // in the cell representing the grid partition of region to a
        //  xCellCount on yCellCount grid
        grid = new Triangle[xCellCount][yCellCount];

        Triangle colStartTriangle = indexDelaunay.find(middleOfCell(0, 0));
        updateCellValues(0, 0, xCellCount - 1, yCellCount - 1, colStartTriangle);
    }
1137855
def delaunay_2d(self, tol=1e-05, alpha=0.0, offset=1.0, bound=False, inplace=False):
        """"""
        alg = vtk.vtkDelaunay2D()
        alg.SetProjectionPlaneMode(vtk.VTK_BEST_FITTING_PLANE)
        alg.SetInputDataObject(self)
        alg.SetTolerance(tol)
        alg.SetAlpha(alpha)
        alg.SetOffset(offset)
        alg.SetBoundingTriangulation(bound)
        alg.Update()

        mesh = _get_output(alg)
        if inplace:
            self.overwrite(mesh)
        else:
            return mesh
93702
def MeshPoints(*inputobj, **options):
    
    r = options.pop("r", 5)
    c = options.pop("c", "gray")
    alpha = options.pop("alpha", 1)

    mesh, u = _inputsort(inputobj)
    plist = mesh.coordinates()
    if u:
        u_values = np.array([u(p) for p in plist])
    if len(plist[0]) == 2:  # coords are 2d.. not good..
        plist = np.insert(plist, 2, 0, axis=1)  # make it 3d
    if len(plist[0]) == 1:  # coords are 1d.. not good..
        plist = np.insert(plist, 1, 0, axis=1)  # make it 3d
        plist = np.insert(plist, 2, 0, axis=1)

    actor = shapes.Points(plist, r=r, c=c, alpha=alpha)

    actor.mesh = mesh
    if u:
        actor.u = u
        if len(u_values.shape) == 2:
            if u_values.shape[1] in [2, 3]:  # u_values is 2D or 3D
                actor.u_values = u_values
                dispsizes = utils.mag(u_values)
        else:  # u_values is 1D
            dispsizes = u_values
        actor.addPointScalars(dispsizes, "u_values")
    return actor
520439
def euclidean3d(v1, v2):
    """"""
    if not len(v1) == 3 and len(v2) == 3:
        print("Vectors are not in 3D space. Returning None.")
        return None
    return np.sqrt((v1[0] - v2[0]) ** 2 + (v1[1] - v2[1]) ** 2 + (v1[2] - v2[2]) ** 2)
1282115
def _contour(f, vertexlabels=None, contourfunc=None, **kwargs):
    ''''''

    if contourfunc is None:
        contourfunc = plt.tricontour
    if vertexlabels is None:
        vertexlabels = ('1','2','3')
    x = np.linspace(0, 1, 100)
    y = np.linspace(0, np.sqrt(3.0)/2.0, 100)
    points2d = np.transpose([np.tile(x, len(y)), np.repeat(y, len(x))])
    points3d = barycentric(points2d)
    valid = (points3d.sum(axis=1) == 1.0) & ((0.0 <= points3d).all(axis=1))
    points2d = points2d[np.where(valid),:][0]
    points3d = points3d[np.where(valid),:][0]
    z = f(points3d)
    contourfunc(points2d[:,0], points2d[:,1], z, **kwargs)
    _draw_axes(vertexlabels)
    return plt.gcf()
990455
def uniformly_refine_triangulation(self, faces=False, trisect=False):
        

        if faces:
            x_v1, y_v1 = self._add_face_centroids()

        else:
            if not trisect:
                x_v1, y_v1 = self._add_midpoints()
            else:
                x_v1, y_v1 = self._add_tripoints(ratio=0.333333)

        return x_v1, y_v1
438090
function ( points, segments, phiStart, phiLength ) {

	BufferGeometry.call( this );

	this.type = '';

	this.parameters = {
		points: points,
		segments: segments,
		phiStart: phiStart,
		phiLength: phiLength
	};

	segments = Math.floor( segments ) || 12;
	phiStart = phiStart || 0;
	phiLength = phiLength || Math.PI * 2;

	// clamp phiLength so it's in range of [ 0, 2PI ]

	phiLength = _Math.clamp( phiLength, 0, Math.PI * 2 );


	// buffers

	var indices = [];
	var vertices = [];
	var uvs = [];

	// helper variables

	var base;
	var inverseSegments = 1.0 / segments;
	var vertex = new Vector3();
	var uv = new Vector2();
	var i, j;

	// generate vertices and uvs

	for ( i = 0; i <= segments; i ++ ) {

		var phi = phiStart + i * inverseSegments * phiLength;

		var sin = Math.sin( phi );
		var cos = Math.cos( phi );

		for ( j = 0; j <= ( points.length - 1 ); j ++ ) {

			// vertex

			vertex.x = points[ j ].x * sin;
			vertex.y = points[ j ].y;
			vertex.z = points[ j ].x * cos;

			vertices.push( vertex.x, vertex.y, vertex.z );

			// uv

			uv.x = i / segments;
			uv.y = j / ( points.length - 1 );

			uvs.push( uv.x, uv.y );


		}

	}

	// indices

	for ( i = 0; i < segments; i ++ ) {

		for ( j = 0; j < ( points.length - 1 ); j ++ ) {

			base = j + i * points.length;

			var a = base;
			var b = base + points.length;
			var c = base + points.length + 1;
			var d = base + 1;

			// faces

			indices.push( a, b, d );
			indices.push( b, c, d );

		}

	}

	// build geometry

	this.setIndex( indices );
	this.addAttribute( 'position', new Float32BufferAttribute( vertices, 3 ) );
	this.addAttribute( 'uv', new Float32BufferAttribute( uvs, 2 ) );

	// generate normals

	this.computeVertexNormals();

	// if the geometry is closed, we need to average the normals along the seam.
	// because the corresponding vertices are identical (but still have different UVs).

	if ( phiLength === Math.PI * 2 ) {

		var normals = this.attributes.normal.array;
		var n1 = new Vector3();
		var n2 = new Vector3();
		var n = new Vector3();

		// this is the buffer offset for the last line of vertices

		base = segments * points.length * 3;

		for ( i = 0, j = 0; i < points.length; i ++, j += 3 ) {

			// select the normal of the vertex in the first line

			n1.x = normals[ j + 0 ];
			n1.y = normals[ j + 1 ];
			n1.z = normals[ j + 2 ];

			// select the normal of the vertex in the last line

			n2.x = normals[ base + j + 0 ];
			n2.y = normals[ base + j + 1 ];
			n2.z = normals[ base + j + 2 ];

			// average normals

			n.addVectors( n1, n2 ).normalize();

			// assign the new values to both normals

			normals[ j + 0 ] = normals[ base + j + 0 ] = n.x;
			normals[ j + 1 ] = normals[ base + j + 1 ] = n.y;
			normals[ j + 2 ] = normals[ base + j + 2 ] = n.z;

		}

	}

}
690137
def transform_to_3d(points,normal,z=0):
    """"""
    d = np.cross(normal, (0, 0, 1))
    M = rotation_matrix(d)
    transformed_points = M.dot(points.T).T + z
    return transformed_points
78475
def triangulate(vertices):
    
    n = len(vertices)
    vertices = np.asarray(vertices)
    zmean = vertices[:, 2].mean()
    vertices_2d = vertices[:, :2]
    segments = np.repeat(np.arange(n + 1), 2)[1:-1]
    segments[-2:] = n - 1, 0

    if _TRIANGLE_AVAILABLE:
        vertices_2d, triangles = _triangulate_cpp(vertices_2d, segments)
    else:
        vertices_2d, triangles = _triangulate_python(vertices_2d, segments)

    vertices = np.empty((len(vertices_2d), 3))
    vertices[:, :2] = vertices_2d
    vertices[:, 2] = zmean
    return vertices, triangles
337471
def construct_mesh_for_1d_array(x_values):
    """"""

    if not isinstance(x_values, np.ndarray):
        try: x_values = np.array(x_values)
        except: raise TypeError("Can not convert to numpy.ndarray: {0}".format(x_values))
    else: assert (x_values.ndim == 1) and (x_values.size >= 2)

    x_mid_points = 0.5 * (x_values[1:] + x_values[:-1])

    x_min_point = x_values[0] - (x_mid_points[0] - x_values[0])
    x_max_point = x_values[-1] + (x_values[-1] - x_mid_points[-1])

    x_mesh_values = np.empty((len(x_values)+1,), dtype=float)
    x_mesh_values[0] = x_min_point
    x_mesh_values[-1] = x_max_point
    x_mesh_values[1:-1] = x_mid_points

    return x_mesh_values
568522
def t_coords(self, t_coords):
        """"""
        if not isinstance(t_coords, np.ndarray):
            raise TypeError('Texture coordinates must be a numpy array')
        if t_coords.ndim != 2:
            raise AssertionError('Texture coordinates must be a 2-dimensional array')
        if t_coords.shape[0] != self.n_points:
            raise AssertionError('Number of texture coordinates ({}) must match number of points ({})'.format(t_coords.shape[0], self.n_points))
        if t_coords.shape[1] != 2:
            raise AssertionError('Texture coordinates must only have 2 components, not ({})'.format(t_coords.shape[1]))
        if np.min(t_coords) < 0.0 or np.max(t_coords) > 1.0:
            raise AssertionError('Texture coordinates must be within (0, 1) range.')
        # convert the array
        vtkarr = numpy_to_vtk(t_coords)
        vtkarr.SetName('Texture Coordinates')
        self.GetPointData().SetTCoords(vtkarr)
        self.GetPointData().Modified()
        return
644151
float[] closestPointOnDetailEdges(MeshTile tile, Poly poly, float[] pos, boolean onlyBoundary) {
        int ANY_BOUNDARY_EDGE = (DT_DETAIL_EDGE_BOUNDARY << 0) | (DT_DETAIL_EDGE_BOUNDARY << 2)
                | (DT_DETAIL_EDGE_BOUNDARY << 4);
        int ip = poly.index;
        PolyDetail pd = tile.data.detailMeshes[ip];

        float dmin = Float.MAX_VALUE;
        float tmin = 0;
        float[] pmin = null;
        float[] pmax = null;

        for (int i = 0; i < pd.triCount; i++) {
            int ti = (pd.triBase + i) * 4;
            int[] tris = tile.data.detailTris;
            if (onlyBoundary && (tris[ti + 3] & ANY_BOUNDARY_EDGE) == 0) {
                continue;
            }

            float[][] v = new float[3][];
            for (int j = 0; j < 3; ++j) {
                if (tris[ti + j] < poly.vertCount) {
                    int index = poly.verts[tris[ti + j]] * 3;
                    v[j] = new float[] { tile.data.verts[index], tile.data.verts[index + 1],
                            tile.data.verts[index + 2] };
                } else {
                    int index = (pd.vertBase + (tris[ti + j] - poly.vertCount)) * 3;
                    v[j] = new float[] { tile.data.detailVerts[index], tile.data.detailVerts[index + 1],
                            tile.data.detailVerts[index + 2] };
                }
            }

            for (int k = 0, j = 2; k < 3; j = k++) {
                if ((getDetailTriEdgeFlags(tris[3], j) & DT_DETAIL_EDGE_BOUNDARY) == 0
                        && (onlyBoundary || tris[j] < tris[k])) {
                    // Only looking at boundary edges and this is internal, or
                    // this is an inner edge that we will see again or have already seen.
                    continue;
                }

                Tupple2<Float, Float> dt = distancePtSegSqr2D(pos, v[j], v[k]);
                float d = dt.first;
                float t = dt.second;
                if (d < dmin) {
                    dmin = d;
                    tmin = t;
                    pmin = v[j];
                    pmax = v[k];
                }
            }
        }

        return vLerp(pmin, pmax, tmin);
    }
250382
def point3d(value, lon, lat, depth):
    
    return longitude(lon), latitude(lat), positivefloat(depth)
428692
def _project_on_ellipsoid(c, r, locations):
    """"""
    p0 = locations - c  # original locations

    l2 = 1 / np.sum(p0**2 / r**2, axis=1, keepdims=True)
    p = p0 * np.sqrt(l2)  # initial approximation (projection of points towards center of ellipsoid)

    fun = lambda x: np.sum((x.reshape(p0.shape) - p0)**2)              # minimize distance between new and old points
    con = lambda x: np.sum(x.reshape(p0.shape)**2 / r**2, axis=1) - 1  # new points constrained to surface of ellipsoid
    res = sp.optimize.minimize(fun, p, constraints={'type': 'eq', 'fun': con}, method='SLSQP')

    return res['x'].reshape(p0.shape) + c
710856
def isosurface_from_data(data, isolevel, origin, spacing):
    
    spacing = np.array(extent/resolution)
    if isolevel >= 0:
        triangles = marching_cubes(data, isolevel)
    else: # Wrong traingle unwinding roder -- god only knows why
        triangles = marching_cubes(-data, -isolevel)
    faces = []
    verts = []
    for i, t in enumerate(triangles):
       faces.append([i * 3, i * 3 +1, i * 3 + 2])
       verts.extend(t)
       
    faces = np.array(faces)
    verts = origin + spacing/2 + np.array(verts)*spacing
    return verts, faces
1285
def triangulate_quads(quads):
    
    if len(quads) == 0:
        return quads
    quads = np.asanyarray(quads)
    faces = np.vstack((quads[:, [0, 1, 2]],
                       quads[:, [2, 3, 0]]))
    return faces