Python geometry.LineString类代码示例

 def test_line_intersection(self):
     line1 = LineString([(0, 0, 0), (1, 1, 1)])
     line2 = LineString([(0, 1, 1), (1, 0, 0)])
     interxn = line1.intersection(line2)
     self.assertTrue(interxn.has_z)
     self.assertEqual(interxn._ndim, 3)
     self.assertTrue(0.0 <= interxn.z <= 1.0)

 def offset(self, distance):
     self.points.append(self.points[0])
     line = LineString(self.getSequence())
     offset = line.parallel_offset(distance, 'right', join_style=1)
     # return list(offset.coords)
     self.setSequence(list(offset.coords))
     return self

def clipLine(lineGeometry, pt1, pt2):
    """Returns a line cliipped to the extent of two given points"""
    # Assumes pt1, pt2 lie on line
    if lineGeometry is None or lineGeometry.isEmpty() or pt1 is None or pt2 is None:
        return QgsGeometry()
    line = LineString(lineGeometry.geometry())
    d1 = line.project(Point(pt1))
    d2 = line.project(Point(pt2))
    if d1 < d2:
        start = pt1
        ds = d1
        end = pt2
        de = d2
    else:
        start = pt2
        ds = d2
        end = pt1
        de = d1
    clip = []
    clip.append(start)
    for coord in line.coords:
        pt = Point(coord)
        dp = line.project(pt)
        if dp > ds and dp < de:
            clip.append(QgsPointV2(pt.x, pt.y))
    clip.append(end)
    return QgsGeometry.fromPolyline(clip)

def acc2latlng(caseno):
    """Calculate the latitude/longitude and update into mongodb

    :param caseno: case number of the accident
    """
    # get the cnty_rte and mile post
    acc = acc_col.find_one({'caseno':caseno})

    # check exists
    if 'lat' in acc.keys() and 'lng' in acc.keys():
        return

    rte_nbr = acc['rte_nbr']
    cnty_rte = acc['cnty_rte']
    milepost = acc['milepost']
    if cnty_rte.index(rte_nbr) != 2:
        raise Exception("Unknown format")

    # get the road information
    rid = rte_nbr + cnty_rte[0:2]
    r_info = cnty_rte2linestring(rid)
    max_mp = r_info['rlist'][-1]['endmp']
    ls = LineString( r_info['plist'] )
    deg_mp = milepost / max_mp * ls.length # mile to degree
    acc_pos = ls.interpolate( deg_mp )
    lng = acc_pos.coords[0][0]
    lat = acc_pos.coords[0][1]
    #print ls
    #print acc_pos
    # TODO, update the lat lng information
    acc_col.update(
            {'caseno':caseno},
            {'$set':{'lat':lat, 'lng':lng}}
            )

    def find_common_path(self, current_adc_trajectory, last_adc_trajectory):
        current_path_points = current_adc_trajectory.trajectory_point
        last_path_points = last_adc_trajectory.trajectory_point
        current_path = []
        for point in current_path_points:
            current_path.append([point.path_point.x, point.path_point.y])
        last_path = []
        for point in last_path_points:
            last_path.append([point.path_point.x, point.path_point.y])
        if len(current_path) == 0 or len(last_path) == 0:
            return [], []

        current_ls = LineString(current_path)
        last_ls = LineString(last_path)
        current_start_point = Point(current_path[0])

        dist = last_ls.project(current_start_point)
        cut_lines = self.cut(last_ls, dist)
        if len(cut_lines) == 1:
            return [], []
        last_ls = cut_lines[1]
        dist = current_ls.project(Point(last_path[-1]))
        if dist <= current_ls.length:
            current_ls = self.cut(current_ls, dist)[0]
        else:
            dist = last_ls.project(Point(current_path[-1]))
            last_ls = self.cut(last_ls, dist)[0]
        return current_ls.coords, last_ls.coords

    def _add_slope(self, bounds, altitude1, altitude2, point1, point2, bottom=False):
        altitude_diff = altitude2-altitude1
        altitude_middle = (altitude1+altitude2)/2
        altitude_halfdiff = altitude_diff/2
        altitude_base = altitude1
        line = LineString([point1, point2])

        minx, miny, maxx, maxy = bounds
        points_2d = [(minx-100, miny-100), (maxx+100, miny-100), (maxx+100, maxy+100), (minx-100, maxy+100)]
        points_3d = []
        for i, (x, y) in enumerate(points_2d):
            point = Point((x, y))
            pos = line.project(point)
            while pos <= 0 or pos >= line.length-1:
                line = scale(line, xfact=2, yfact=2, zfact=2)
                altitude_diff *= 2
                altitude_halfdiff *= 2
                altitude_base = altitude_middle-altitude_halfdiff
                pos = line.project(point)
            z = ((pos/line.length)*altitude_diff)+altitude_base
            points_3d.append((x, y, z/1000))

        extrude = abs(altitude1-altitude2)/1000+100
        if bottom:
            extrude = -extrude
        self._add_python(
            'last_slope = add_polygon_3d(name=%(name)r, coords=%(coords)r, extrude=%(extrude)f)' % {
                'name': 'tmpslope',
                'coords': tuple(points_3d),
                'extrude': extrude,
            }
        )

    def _add_crossroads(new_edges, current_edges):
        new_new_edges = set()
        new_current_edges = set()
        split_current_edge = [False] * len(current_edges)

        for new_edge in new_edges:
            found_crossroad = False
            for current_edge_index, current_edge in enumerate(current_edges):
                if State._different_starts_and_ends(new_edge, current_edge):
                    new_line_segment = LineString([(new_edge[0].x, new_edge[0].y), (new_edge[1].x, new_edge[1].y)])
                    current_line_segment = LineString(([(current_edge[0].x, current_edge[0].y),
                                                        (current_edge[1].x, current_edge[1].y)]))

                    intersection = new_line_segment.intersection(current_line_segment)
                    if type(intersection) is ShapelyPoint:
                        new_crossroad = Point(intersection.x, intersection.y, Point.TYPE_CROSSROAD)
                        new_current_edges.update(State._split_edge(current_edge, new_crossroad))
                        new_new_edges.update(State._split_edge(new_edge, new_crossroad))
                        split_current_edge[current_edge_index] = True
                        found_crossroad = True
                        break
            if not found_crossroad:
                new_new_edges.add(new_edge)

        for current_edge_index, current_edge in enumerate(current_edges):
            if not split_current_edge[current_edge_index]:
                new_current_edges.add(current_edge)

        if set(current_edges) == new_current_edges and set(new_edges) == new_new_edges:
            return list(current_edges) + list(new_edges)
        else:
            return State._add_crossroads(list(new_new_edges), list(new_current_edges))

def append_to_included_groups(locs, d):
    for group_key in d.keys():
        ref_gene_poly = LineString([(0.0, float(group_key[0])),(0.0, float(group_key[1]))])
        locs_poly = LineString([(0, locs[2]), (0, locs[3])])
        if ref_gene_poly.intersects(locs_poly):
            d[group_key].append(tuple(locs))
    return d

def cut_line(cut_point, line, eps_mult=1e2):
    dist = line.project(cut_point)
    point = line.interpolate(dist)
    eps = line.distance(point) * eps_mult

    coords = list(line.coords)

    if point.coords[0] in coords:
        i = coords.index(point.coords[0])

        if i == 0:
            return LineString(), line
        if i == len(coords) - 1:
            return line, LineString()

        start_segment = LineString(coords[:i + 1])
        end_segment = LineString(coords[i:])

        return start_segment, end_segment


    for i, p in enumerate(coords[:-1]):
        line_segment = LineString([coords[i], coords[i + 1]])
        line_segment_buffer = line_segment.buffer(eps, resolution=1)

        if line_segment_buffer.contains(point):
            start_segment = LineString(coords[:i + 1] + [point])
            end_segment = LineString([point] + coords[i + 1:])

            return start_segment, end_segment

    raise Exception('point not found in line, consider raising eps_mult')

def getShortestDisLinePoint(pointLon, pointLat, lineCoords, earthR):
	#get shortest distance between point and line
	line = LineString(lineCoords) 
	point = Point(pointLon, pointLat)
	closestCoord = line.interpolate(line.project(point))
	distance = haversine(pointLon, pointLat, closestCoord.x, closestCoord.y, earthR)
	return distance

def identify_first_point_in_polygon(points, ddd=0.5):
    """
    Identify the cycle beginning for the last point
    :param points: array of objects of PointX
    :param ddd: distance to the fist point ddd>0.
    :return -1 if do not close, or i if close in point i.
    """
    if len(points) < 3:
        return -1

    # last point
    lp = Point(points[-1])

    first_point = len(points) - 2

    # Distance between last point and a line segment
    distance_point_line = 2 * ddd

    # Distance is not less than ddd or less than 3 points.
    while not (distance_point_line < ddd and polyline_length(points[first_point:]) > 2):
        # next segment
        first_point -= 1

        if first_point == -1:
            break

        pt1 = points[first_point + 1]
        pt2 = points[first_point]

        # Distance to the first point to the line segment
        line_segment = LineString([pt1, pt2])
        distance_point_line = line_segment.distance(lp)

    return first_point

    def divide_polygon_for_intersection(self, segments):
        """ Generates multiple polygons based on cutting the 
        fracture faces by line segments. 
        """
        R = self.build_rotation_matrix()
        fracture_poly_list = []
        # face_poly_list = []

        for point in self.points:
            rot_point = np.linalg.solve(R, point - self.center)
            fracture_poly_list.append(rot_point[:2])

        seg_rot = []
        for seg in segments:
            p1 = seg[0]
            p2 = seg[1]

            vec = p2 - p1
            p1 = p1 - 100.0 * vec
            p2 = p2 + 100.0 * vec

            p1_rot = np.linalg.solve(R, p1 - self.center)
            p2_rot = np.linalg.solve(R, p2 - self.center)

            line = LineString((p1_rot[:2], p2_rot[:2]))
            dilated = line.buffer(1.0e-10)
            seg_rot.append(dilated)

        fracture_poly = Polygon(fracture_poly_list)

        divided_polygons = fracture_poly.difference(seg_rot[0])

        return (divided_polygons, fracture_poly)

 def test_line_intersection(self):
     line1 = LineString([(0, 0, 0), (1, 1, 1)])
     line2 = LineString([(0, 1, 1), (1, 0, 0)])
     interxn = line1.intersection(line2)
     self.failUnless(interxn.has_z)
     self.failUnless(interxn._ndim == 3)
     self.failUnless(0.0 <= interxn.z <= 1.0)

 def setUp(self):
     self.point = Point(1, 1)
     self.line1 = LineString(([0, 0], [2, 0]))
     self.line2 = LineString(([3, 0], [3, 6]))
     self.multiline = MultiLineString([
         list(self.line1.coords), list(self.line2.coords)
     ])

def lines_fusion(line1, line2):
    """
    Validate each line segment for intersection. Intersection is checked as
    http://en.wikipedia.org/wiki/Line_segment_intersection

    :param line1: main line
    :param line2: line to be fused
    :return same line1 if there is no intersections with line2. else fused lines.
    """
    # TODO use bounding boxes to speed up
    # Cross validation for each line segment in polylines.
    for i in range(len(line1) - 1, 0, -1):
        p1 = line1[i]
        p2 = line1[i - 1]
        # line segment
        l1 = LineString([p1, p2])

        for j in range(len(line2) - 1, 0, -1):
            p3 = line2[j]
            p4 = line2[j - 1]

            l2 = LineString([p3, p4])
            # Check Line intersection
            intersection = l1.intersection(l2)

            # if intersected
            if not intersection.is_empty:
                # take the line 1 from first point until the intersection and line2
                inter_p = (intersection.coords.xy[0][0], intersection.coords.xy[1][0])
                new_line = line2[:j] + [inter_p] + line1[i:]
                return new_line

    # no fusion
    return line1