def point_projects_to_line(point, line):
"""Get the nearest point index on line
"""
coords = list(line.coords)
p = Point(point)
pd = line.project(p)
for i in range(1, len(coords)):
pp = Point(coords[i-1])
cp = Point(coords[i])
prev = line.project(pp)
cur = line.project(cp)
if cur == pd:
return i
if prev == pd:
return i - 1
if prev < pd < cur:
pdist = p.distance(pp)
cdist = p.distance(cp)
return i-1 if pdist <= cdist else i
return None
python类Point()的实例源码
def validate_intersection(self, line_index, point):
"""
Validate the intersection between given line and point.
A spatial tolerance is considered to confirm the validation.
:param line_index: index of specific line
:param point: validated point
:return: coordinates of intersected point on line, otherwise None
"""
line = self.geoms[line_index]
coords = list(line.coords)
buffered_point = Point(point).buffer(self.point_buffer)
touched = None
if line.intersects(buffered_point):
cut = point_projects_to_line(point, line)
touched = coords[cut]
self._update_cut(line_index, cut)
return touched
def get_region(lat, lon, region_path=REGION_PATH):
"""
Return the 2 character, 1 number code (e.g., IP-3) associated with
the physiographic region that encloses the point specified by
*lat* and *lon*. Return `None` if the point is not interior to any
region. Look for the required files at *region_path*.
"""
try:
regions = shpreader.Reader(os.path.join(region_path, 'ConductivityRegions'))
except:
initialize()
regions = shpreader.Reader(os.path.join(region_path, 'ConductivityRegions'))
geos = regions.geometries()
names = [x.attributes['Name'] for x in regions.records()]
point = Point(convert_lon(lon),
lat)
encloses = [name for name, geo in zip(names, geos) if geo.contains(point)]
if encloses:
assert len(encloses) == 1
return encloses[0]
else:
return None
def get_contour_mask(dd, id, dosegridpoints, contour):
"""Get the mask for the contour with respect to the dose plane."""
doselut = dd['lut']
c = matplotlib.path.Path(list(contour))
# def inpolygon(polygon, xp, yp):
# return np.array(
# [Point(x, y).intersects(polygon) for x, y in zip(xp, yp)],
# dtype=np.bool)
# p = Polygon(contour)
# x, y = np.meshgrid(np.array(dd['lut'][0]), np.array(dd['lut'][1]))
# mask = inpolygon(p, x.ravel(), y.ravel())
# return mask.reshape((len(doselut[1]), len(doselut[0])))
grid = c.contains_points(dosegridpoints)
grid = grid.reshape((len(doselut[1]), len(doselut[0])))
return grid
def intersect_point_layer_with_wkt(layer, wkt, buffer_width):
""" Return all features from given layer that intersects with wkt """
assert QgsWKBTypes.geometryType(
int(layer.wkbType())) == QgsWKBTypes.PointGeometry
line = loads(wkt.replace('Z', ' Z'))
if not line.is_valid:
logging.warning('Invalid feature geometry wkt={}'.format(wkt))
return
# square cap style for the buffer -> less points
buf = line.buffer(buffer_width, cap_style=2)
bbox = QgsRectangle(
buf.bounds[0], buf.bounds[1], buf.bounds[2], buf.bounds[3])
# request features inside bounding-box
for feature in layer.getFeatures(QgsFeatureRequest(bbox)):
p = feature.geometry().asPoint()
if buf.contains(Point(p.x(), p.y())):
yield feature
def create_memory_layer(geometry_type, crs, name, custom_properties=None):
assert geometry_type in (QGis.Point, QGis.Line, QGis.Polygon)
layer = QgsVectorLayer(
"{}?crs={}&index=yes".format(
{
QGis.Point: "Point",
QGis.Line: "LineString",
QGis.Polygon: "Polygon"
}[geometry_type],
crs.authid()
), name, "memory")
layer.setCrs(crs)
if custom_properties:
for key in custom_properties:
layer.setCustomProperty(key, custom_properties[key])
return layer
def sampleShapes(shapes, xypoints, attribute):
"""
Get the attribute value at each of the input XY points for sequence of
input shapes (decimal degrees). Slower than sampling grids.
Args:
shapes: Sequence of shapes (decimal degrees) of predictor variable.
xypoints: 2D numpy array of XY points, in decimal degrees.
attribute: String name of attribute to sample in each of the shapes.
Returns:
1D array of attribute values at each of XY points.
"""
samples = []
for x, y in xypoints:
for tshape in shapes:
polygon = shape(tshape['geometry'])
point = Point(x, y)
if polygon.contains(point):
sample = tshape['properties'][attribute]
samples.append(sample)
return np.array(samples)
def getPointByBoundsWithFilterHandler(self,index, polygon, step, bounds):
logger.info('current index %s' % index)
print 'current index %s' % index
x1 = bounds[0]
y1 = bounds[1]
x2 = bounds[2]
y2 = bounds[3]
points = []
print x1,y1 , x2,y2
print (y2 - y1) / float(step)
print (x2 - x1) / float(step)
# TODO ?????????????
for startY in frange2(y1, y2, step):
for startX in frange2(x1, x2, step):
point = str(startY) + ',' + str(startX)
point_wkt = Point(startX, startY)
isIntersects = self.isIntersectsPolygon(polygon, point_wkt)
if isIntersects:
points.append(point)
print "range points size : ", len(points)
return points
# geometry???polygon??
def testTOW(self):
cm = CAP.CAPMessage(TestCAP._get_test_messages()[0])
self.assertEqual('TO.W', cm.get_event_type())
self.assertEqual('http://alerts.weather.gov/cap/wwacapget.php?x=KS1255FCC0A5BC.TornadoWarning.1255FCC0C36CKS.GLDTORGLD.3a1fc090003ef1448f822dfd9b2ddee2', cm.get_event_id())
self.assertEqual('KGLD.TO.W.0021', cm.vtec[-1].event_id)
self.assertEqual(cm, cm.vtec[-1].container)
self.assertEqual(1463877540.0, cm.get_start_time_sec())
self.assertEqual(1463879700.0, cm.get_end_time_sec())
self.assertFalse(cm.is_effective(when=1463877539.9))
self.assertTrue(cm.is_effective(when=1463877540.0))
self.assertTrue(cm.is_effective(when=1463877541.0))
self.assertTrue(cm.is_effective(when=1463879700.0))
self.assertFalse(cm.is_effective(when=1463879700.1))
self.assertTrue(cm.polygon.contains(Point(38.80, -101.45)))
self.assertFalse(cm.polygon.contains(Point(38.90, -101.45)))
self.assertEqual("CAP [ Sun May 22 00:39:00 2016 TO.W /O.NEW.KGLD.TO.W.0021.160522T0039Z-160522T0115Z/ ['020109', '020199'] ]",str(cm))
def red_pixels(self):
pixels = np.zeros(self._red_pixel_sensors)
for i in range(self._red_pixel_sensors):
rotation = self._robot.rotation - self._fov / 2 + i * self._fov / (self._red_pixel_sensors-1)
ray = affinity.rotate(
LineString([np.array(self._robot.pos.coords)[0],
np.array(self._robot.pos.coords)[0] + (self._height*2, 0)]),
rotation, origin=self._robot.pos, use_radians=True)
dists_red = np.min(np.array([self._robot.pos.distance(point) if not point.is_empty else np.inf
for point in (ray.intersection(cube) for cube in self._red_cubes)]))
dists_obs = np.min(np.array([self._robot.pos.distance(point) if not point.is_empty else np.inf
for point in (ray.intersection(cube) for cube in
#filter(lambda o: np.array([self._robot.pos.distance(Point(p)) < dists_red for p in o.exterior.coords]).any(),
self._obstacles)]))
pixels[i] = dists_red < dists_obs
return pixels
def removeIgnoredPointsRects(rects,polyList):
ridxs = list(range(len(rects)))
for ridx in range(len(rects)):
points = rects[ridx]["annopoints"][0]["point"]
pidxs = list(range(len(points)))
for pidx in range(len(points)):
pt = geometry.Point(points[pidx]["x"][0], points[pidx]["y"][0])
bIgnore = False
for poidx in range(len(polyList)):
poly = polyList[poidx]
if (poly.contains(pt)):
bIgnore = True
break
if (bIgnore):
pidxs.remove(pidx)
points = [points[pidx] for pidx in pidxs]
if (len(points) > 0):
rects[ridx]["annopoints"][0]["point"] = points
else:
ridxs.remove(ridx)
rects = [rects[ridx] for ridx in ridxs]
return rects
def removeIgnoredPoints(gtFramesAll,prFramesAll):
imgidxs = []
for imgidx in range(len(gtFramesAll)):
if ("ignore_regions" in gtFramesAll[imgidx].keys() and
len(gtFramesAll[imgidx]["ignore_regions"]) > 0):
regions = gtFramesAll[imgidx]["ignore_regions"]
polyList = []
for ridx in range(len(regions)):
points = regions[ridx]["point"]
pointList = []
for pidx in range(len(points)):
pt = geometry.Point(points[pidx]["x"][0], points[pidx]["y"][0])
pointList += [pt]
poly = geometry.Polygon([[p.x, p.y] for p in pointList])
polyList += [poly]
rects = prFramesAll[imgidx]["annorect"]
prFramesAll[imgidx]["annorect"] = removeIgnoredPointsRects(rects,polyList)
rects = gtFramesAll[imgidx]["annorect"]
gtFramesAll[imgidx]["annorect"] = removeIgnoredPointsRects(rects,polyList)
return gtFramesAll, prFramesAll
def process_item(d, poly_dts):
key = d.key.name
lat = d['lat']
lon = -d['lon']
dt = d['image_datetime']
point = Point(lon, lat)
# # Based on photo's time, find the umbra whose center point's time is the same
# TODO(dek): fix https://b.corp.google.com/issues/64974121
pdt = dt.replace(tzinfo=None)
if pdt not in poly_dts:
print "Point outside eclipse time window:", key, pdt, lat, lon
return key, False
current_poly = poly_dts[pdt][0]
if not current_poly.contains(point):
print "Point outside eclipse time window:", key, pdt, lat,lon
return key, False
# Now we know this photo point/dt is in totality
return key, True
def generate_polygon(points):
from shapely.geometry import Polygon, Point, LineString
p = []
p.append(points[0][1])
for point in points:
p.append(point[0])
p.append(points[-1][1])
for point in points[::-1]:
p.append(point[2])
eclipse_boundary = Polygon(p)
p = []
for point in points:
p.append(point[1])
center_line = LineString(p)
return eclipse_boundary, center_line
def generate_polygon(points):
from shapely.geometry import Polygon, Point, LineString
p = []
p.append(points[0][1])
for point in points:
p.append(point[0])
p.append(points[-1][1])
for point in points[::-1]:
p.append(point[2])
eclipse_boundary = Polygon(p)
p = []
for point in points:
p.append(point[1])
center_line = LineString(p)
return eclipse_boundary, center_line
def process_item(d, fname, polys, poly_dts):
lat = d['lat']
lon = -d['lon']
dt = d['image_datetime']
point = Point(lon, lat)
# # Based on photo's time, find the umbra whose center point's time is the same
# TODO(dek): fix https://b.corp.google.com/issues/64974121
pdt = dt.replace(tzinfo=None)
if pdt not in poly_dts:
print "Point outside eclipse time window:", fname, pdt, lat, lon
return None
current_poly = poly_dts[pdt][0]
if not current_poly.contains(point):
print "Point outside eclipse time window:", fname, pdt, lat,lon
return None
# Now we know this photo point/dt is in totality
# Find all umbra polys that contain this photo
x = []
for j, p in enumerate(polys):
poly, poly_centroid, poly_dt = p
if poly.contains(point):
x.append((j, (poly_centroid.distance(point), fname, lat, lon, dt)))
return x
def genPolygonFromConf(self, polygon, configuration, refpoint=Point()):
"""
recover the workspace polygon using a configuration point
:param configuration:
:return:
"""
if refpoint.is_empty:
xoff = configuration[0]
yoff = configuration[1]
rotpoly = self.rotate(polygon, (configuration[2])/self.scale)
transpoly = self.translate(rotpoly, Point(xoff, yoff))
# print configuration[0], configuration[1]
# print transpoly.centroid.x, transpoly.centroid.y
return transpoly
else:
xoff = configuration[0]
yoff = configuration[1]
rotpoly = self.rotate(polygon, configuration[2]/self.scale)
transpoly = self.translate(rotpoly, Point(xoff, yoff))
return transpoly
def nodes_for_point(self, point, all_nodes):
point = Point(point.x, point.y)
nodes = {}
if self.nodes:
for node in self.nodes:
node = all_nodes[node]
line = LineString([(node.x, node.y), (point.x, point.y)])
if line.length < 5 and not self.clear_geometry_prep.intersects(line):
nodes[node.i] = (None, None)
if not nodes:
nearest_node = min(tuple(all_nodes[node] for node in self.nodes),
key=lambda node: point.distance(node.point))
nodes[nearest_node.i] = (None, None)
else:
nodes = self.fallback_nodes
return nodes
def testNoAdditionalFields(self):
polygons = [Point((0, 0))]
labels = [1]
timing = WorkflowTiming()
info = WorkflowInformation(polygons, labels, timing=timing)
self.assertEqual(len(info), 1)
self.assertEqual(timing, info.timing)
assert_array_equal(labels, info.labels)
assert_array_equal(np.array(polygons, dtype=np.object), info.polygons)
self.assertSetEqual(set(info.fields), {info.DATA_FIELD_POLYGONS, info.DATA_FIELD_LABELS})
object_info = info[0]
self.assertEqual(len(object_info), 2)
self.assertEqual(object_info.polygon, polygons[0])
self.assertEqual(object_info.label, labels[0])
object_info_iter = next(iter(info))
self.assertEqual(len(object_info_iter), 2)
self.assertEqual(object_info_iter.polygon, polygons[0])
self.assertEqual(object_info_iter.label, labels[0])
def testAdditionalFields(self):
polygons = [Point((0, 0)), Point((0, 1))]
labels = [1, 2]
dispatch = [5, 3]
timing = WorkflowTiming()
info = WorkflowInformation(polygons, labels, timing=timing, dispatches=(dispatch, "dispatch"))
self.assertEqual(len(info), 2)
self.assertEqual(timing, info.timing)
assert_array_equal(labels, info.labels)
assert_array_equal(np.array(polygons, dtype=np.object), info.polygons)
assert_array_equal(dispatch, info.dispatches)
self.assertSetEqual(set(info.fields), {info.DATA_FIELD_POLYGONS, info.DATA_FIELD_LABELS, "dispatches"})
for i, object_info in enumerate(info):
self.assertEqual(info[i], object_info)
self.assertEqual(len(object_info), 3)
self.assertEqual(object_info.polygon, polygons[i])
self.assertEqual(object_info.label, labels[i])
self.assertEqual(object_info.dispatch, dispatch[i])
def testErrors(self):
p1, p2 = Point(0, 0), Point(1, 0)
timing = WorkflowTiming()
with self.assertRaises(ValueError):
WorkflowInformation([p1, p2], [1], timing)
with self.assertRaises(ValueError):
WorkflowInformation([p1, p2], [1, 2], timing, others=([2], "other"))
with self.assertRaises(ValueError):
WorkflowInformation([p1, p2], [1, 2], timing, __init__=([1, 2], "__init__s"))
first = WorkflowInformation([p1, p2], [1, 2], timing, others=([2, 1], "other"))
second = WorkflowInformation([p1, p2], [1, 2], timing)
third = WorkflowInformation([p1, p2], [1, 2], timing, others=([2, 1], "oo"))
with self.assertRaises(TypeError):
first._is_compatible(dict())
with self.assertRaises(ValueError):
first._is_compatible(second)
with self.assertRaises(ValueError):
first._is_compatible(third)
def _geom_to_array(geom: BaseGeometry):
if isinstance(geom, geometry.Point):
yield np.array([(np.nan, GeomTypes.POINT)])
yield np.asarray(geom.coords)
elif isinstance(geom, geometry.LineString):
yield np.array([(np.nan, GeomTypes.LINESTRING)])
yield np.asarray(geom.coords)
elif isinstance(geom, geometry.Polygon):
for interior in geom.interiors:
yield np.array([(np.nan, GeomTypes.POLYGON_HOLE)])
yield np.asarray(interior)
yield np.array([(np.nan, GeomTypes.POLYGON_SHELL)])
yield np.asarray(geom.exterior)
elif isinstance(geom, BaseMultipartGeometry):
return chain.from_iterable(map(geom_to_array, geom))
else:
raise TypeError
def add_geom(fig: Figure, geom: BaseGeometry, **kwargs):
"""Add Shapely geom into Bokeh plot.
Args:
fig (Figure):
geom (BaseGeometry):
"""
if isinstance(geom, Point):
fig.circle(*geom.xy, **kwargs)
elif isinstance(geom, LineString):
fig.line(*geom.xy, **kwargs)
elif isinstance(geom, Polygon):
fig.patch(*geom.exterior.xy, **kwargs)
elif isinstance(geom, BaseMultipartGeometry):
for item in geom:
add_geom(fig, item, **kwargs)
else:
raise TypeError('Object geom {geom} no instance of {types}.'.format(
geom=geom, types=BaseGeometry))
def test_linearize():
t_start = 1
t_stop = 6
xy = np.array(([1., 1.], [2., 2.], [3., 3.], [4., 4.], [5., 5.]))
time = np.array([0., 1., 2., 3., 4.])
pos = nept.Position(xy, time)
trajectory = [[0., 0.], [5., 5.], [10., 10.]]
line = LineString(trajectory)
zone_start = Point([1., 1.])
zone_stop = Point([9., 9.])
expand_by = 1
zone = nept.expand_line(zone_start, zone_stop, line, expand_by)
sliced_pos = pos[t_start:t_stop]
linear = sliced_pos.linearize(line, zone)
assert np.allclose(linear.x, np.array([2.82842712, 4.24264069, 5.65685425, 7.07106781]))
assert np.allclose(linear.time, [1, 2, 3, 4])
def test_position_linearize():
times = np.array([1.0, 2.0, 3.0])
data = np.array([[0.0, 0.5],
[0.5, 0.1],
[1.0, 1.2]])
pos = nept.Position(data, times)
line = LineString([(0.0, 0.0), (1.0, 1.0)])
zone_start = Point([1., 1.])
zone_stop = Point([9., 9.])
expand_by = 1
zone = nept.expand_line(zone_start, zone_stop, line, expand_by)
linear = pos.linearize(line, zone)
assert np.allclose(linear.x, np.array([0.35355339, 0.42426407, 1.41421356]))
assert np.allclose(linear.time, np.array([1., 2., 3.]))
def linearize(self, ideal_path, zone):
""" Projects 2D positions into an 'ideal' linear trajectory.
Parameters
----------
ideal_path : shapely.LineString
zone : shapely.Polygon
Returns
-------
pos : nept.Position
1D position.
"""
zpos = []
for point_x, point_y in zip(self.x, self.y):
point = Point([point_x, point_y])
if zone.contains(point):
zpos.append(ideal_path.project(Point(point_x, point_y)))
zpos = np.array(zpos)
return Position(zpos, self.time)
def split(line, point):
distance_on_line = line.project(point)
coords = list(line.coords)
for i, p in enumerate(coords):
pd = line.project(geometry.Point(p))
if pd == distance_on_line:
return [
geometry.LineString(coords[:i+1]),
geometry.LineString(coords[i:])
]
elif distance_on_line < pd:
cp = line.interpolate(distance_on_line)
ls1_coords = coords[:i]
ls1_coords.append(cp.coords[0])
ls2_coords = [cp.coords[0]]
ls2_coords.extend(coords[i:])
return [geometry.LineString(ls1_coords),
geometry.LineString(ls2_coords)]
def is_in_grid(self, shp):
"""Check if arbitrary polygon is within grid bounding box.
Args:
shp (shape):
Returns:
Bool whether shp is in grid box.
Depends on self.grid_box (shapely prep type) being defined
in environment.
"""
if not hasattr(shp, 'geom_type'):
raise Exception("CoreMSR [is_in_grid] : invalid shp given")
if not isinstance(self.grid_box, type(prep(Point(0,0)))):
raise Exception("CoreMSR [is_in_grid] : invalid prep_adm0 "
"found")
return self.grid_box.contains(shp)
def load_nodes(self, numbersuperpoints=None):
coords = {}
supercoords = {}
with open(self.tnxy, 'rb') as f:
with open(self.tnz, 'rb') as fz:
i = 0
while True:
chunkx = f.read(8)
chunky = f.read(8)
chunkz = fz.read(4)
if any(c == '' for c in [chunkx, chunky, chunkz]):
break
else:
i = i + 1
if i <= numbersuperpoints:
supercoords[i] = Point(struct.unpack('>d', chunkx)[0],
struct.unpack('>d', chunky)[0],
struct.unpack('>f', chunkz)[0])
else:
coords[i] = Point(struct.unpack('>d', chunkx)[0],
struct.unpack('>d', chunky)[0],
struct.unpack('>f', chunkz)[0])
return coords, supercoords
def points_along_boundaries(geoseries, distance=1.0):
"""
Generate a shapely MultiPoint of point features along lines at a specified distance
:param geoseries:
:param distance:
:return:
"""
list_points = []
for line3d in iter_line(geoseries):
line = LineString([xy[0:2] for xy in list(line3d.coords)])
current_dist = distance
line_length = line.length
list_points.append(Point(list(line.coords)[0]))
while current_dist < line_length:
list_points.append(line.interpolate(current_dist))
current_dist += distance
list_points.append(Point(list(line.coords)[-1]))
return MultiPoint(list_points)
