Python features.shapes函数代码示例

def test_shapes_blank_mask(basic_image):
    """Mask is blank so results should mask shapes without mask."""
    assert np.array_equal(
        list(shapes(
            basic_image,
            mask=np.ones(basic_image.shape, dtype=rasterio.bool_))
        ),
        list(shapes(basic_image))
    )

def test_shapes_band(pixelated_image, pixelated_image_file):
    """Shapes from a band should match shapes from an array."""
    truth = list(shapes(pixelated_image))

    with rasterio.open(pixelated_image_file) as src:
        band = rasterio.band(src, 1)
        assert truth == list(shapes(band))

        # Mask band should function, but will mask out some results
        assert truth[0] == list(shapes(band, mask=band))[0]

def test_shapes_blank_mask(basic_image):
    """ Mask is blank so results should mask shapes without mask """

    with rasterio.drivers():
        assert numpy.array_equal(
            list(shapes(
                basic_image,
                mask=numpy.ones(basic_image.shape, dtype=rasterio.bool_))
            ),
            list(shapes(basic_image))
        )

def union(inputtiles, parsenames):

    tiles = sutils.tile_parser(inputtiles, parsenames)

    xmin, xmax, ymin, ymax = sutils.get_range(tiles)

    zoom = sutils.get_zoom(tiles)

    # make an array of shape (xrange + 3, yrange + 3)
    burn = sutils.burnXYZs(tiles, xmin, xmax, ymin, ymax, 0)

    nw = mercantile.xy(*mercantile.ul(xmin, ymin, zoom))

    se = mercantile.xy(*mercantile.ul(xmax + 1, ymax + 1, zoom))

    aff = Affine(((se[0] - nw[0]) / float(xmax - xmin + 1)), 0.0, nw[0],
        0.0, -((nw[1] - se[1]) / float(ymax - ymin + 1)), nw[1])

    unprojecter = sutils.Unprojecter()

    unionedTiles = [
        {
            'geometry': unprojecter.unproject(feature),
            'properties': {},
            'type': 'Feature'
        } for feature, shapes in features.shapes(np.asarray(np.flipud(np.rot90(burn)).astype(np.uint8), order='C'), transform=aff) if shapes == 1
    ]

    return unionedTiles

def test_shapes_internal_driver_manager():
    """Access to shapes of labeled features"""
    image = numpy.zeros((20, 20), dtype=rasterio.ubyte)
    image[5:15,5:15] = 127
    shapes = ftrz.shapes(image)
    shape, val = next(shapes)
    assert shape['type'] == 'Polygon'

def polygonize(input_tif_dir, output_shp_dir, countries):

    rm_and_mkdir(output_shp_dir)
    for country in countries:
        shp_filename = country + '.shp'
        output_shp_path = os.path.join(output_shp_dir, shp_filename)
        tif_filename = country + '.tif'
        input_tif_path = os.path.join(input_tif_dir, tif_filename)

        with rasterio.open(input_tif_path) as src:
            band = src.read(1)

        mask = band != 255
        shapes = features.shapes(band, mask=mask, transform=src.transform)
        geomvals = list(shapes)

        geom_val_trios = []
        for idx, geom_val in enumerate(geomvals):
            shapely_geom = shape(geomvals[idx][0])
            shapely_val = geomvals[idx][1]
            geom_val_trio = [shapely_geom, shapely_val, country]
            geom_val_trios.append(geom_val_trio)
        gdf = gpd.GeoDataFrame(geom_val_trios, columns={'geometry', 'val', 'country'})
        gdf.crs = {'init': 'epsg:4326', 'no_defs': True}
        gdf.to_file(output_shp_path)

def main(raster_file, vector_file, driver, mask_value):
    
    with rasterio.drivers():
        
        with rasterio.open(raster_file) as src:
            image = src.read_band(1)
        
        if mask_value is not None:
            mask = image == mask_value
        else:
            mask = None
        
        results = (
            {'properties': {'raster_val': v}, 'geometry': s}
            for i, (s, v) 
            in enumerate(
                shapes(image, mask=mask, transform=src.affine)))

        with fiona.open(
                vector_file, 'w', 
                driver=driver,
                crs=src.crs,
                schema={'properties': [('raster_val', 'int')],
                        'geometry': 'Polygon'}) as dst:
            dst.writerecords(results)
    
    return dst.name

def test_shapes(basic_image):
    """ Test creation of shapes from pixel values """

    with Env():
        results = list(shapes(basic_image))

        assert len(results) == 2

        shape, value = results[0]
        assert shape == {
            'coordinates': [
                [(2, 2), (2, 5), (5, 5), (5, 2), (2, 2)]
            ],
            'type': 'Polygon'
        }
        assert value == 1

        shape, value = results[1]
        assert shape == {
            'coordinates': [
                [(0, 0), (0, 10), (10, 10), (10, 0), (0, 0)],
                [(2, 2), (5, 2), (5, 5), (2, 5), (2, 2)]
            ],
            'type': 'Polygon'
        }
        assert value == 0

    def polygonize(self, band_number=1):
        """
        Extract shapes from raster features. This is the inverse operation of rasterizing shapes.
        Uses the `Rasterio <https://mapbox.github.io/rasterio/_modules/rasterio/features.html>'_ library
        for this purpose. The data is loaded into a `geopandas <http://geopandas.org/user.html>`_ GeoDataFrame.
        GeoDataFrame data structures are pandas DataFrames with added functionality, containing a ``geometry``
        column for the `Shapely <http://toblerity.org/shapely/shapely.geometry.html>`_ geometries.
        The raster data should be loaded in the layer before calling this method.

        :param int band_number: The index of the raster band which is to be used as input for extracting \
        gemetrical shapes.

        :returns: geopandas.GeoDataFrame

        """
        raster_data = self.read(band_number)
        mask = raster_data != self.raster_reader.nodata
        T0 = self.raster_reader.affine
        shapes = features.shapes(raster_data, mask=mask, transform=T0)
        df = GeoDataFrame.from_records(shapes, columns=['geometry', 'value'])
        # convert the geometry dictionary from a dictionary format like {'coordinates': [[(-73.5, 83.5),
        #   (-73.5, 83.0),
        #   (-68.0, 83.0),
        #   (-68.0, 83.5),
        #   (-73.5, 83.5)]],
        # 'type': 'Polygon'}
        # to a proper shapely polygon format
        df.geometry = df.geometry.apply(lambda row: Polygon(row['coordinates'][0]))
        df.crs = self.raster_reader.crs
        return df

def test_shapes_connectivity_rook(diagonal_image):
    """
    Diagonals are not connected, so there will be 1 feature per pixel plus
    background.
    """
    with rasterio.Env():
        assert len(list(shapes(diagonal_image, connectivity=4))) == 12

def test_shapes_connectivity_queen(diagonal_image):
    """
    Diagonals are connected, so there will be 1 feature for all pixels plus
    background.
    """
    with rasterio.Env():
        assert len(list(shapes(diagonal_image, connectivity=8))) == 2

def raster_shape(raster_path):
    with rasterio.open(raster_path) as src:

        # read the first band and create a binary mask
        arr = src.read(1)
        ndv = src.nodata
        binarray = (arr == ndv).astype('uint8')

        # extract shapes from raster
        shapes = features.shapes(binarray, transform=src.transform)

        # create geojson feature collection
        fc = {
            'type': 'FeatureCollection',
            'features': []}
        for geom, val in shapes:
            if val == 0:  # not nodata, i.e. valid data
                feature = {
                    'type': 'Feature',
                    'properties': {'name': raster_path},
                    'geometry': geom}
                fc['features'].append(feature)

        # Write to file
        with NamedTemporaryFile(suffix=".geojson", delete=False) as temp:
            temp.file.write(json.dumps(fc))

        return temp.name

def test_rasterize_geometries_symmetric():
    """Make sure that rasterize is symmetric with shapes."""
    transform = (1.0, 0.0, 0.0, 0.0, -1.0, 0.0)
    truth = np.zeros(DEFAULT_SHAPE, dtype=rasterio.ubyte)
    truth[2:5, 2:5] = 1
    s = shapes(truth, transform=transform)
    result = rasterize(s, out_shape=DEFAULT_SHAPE, transform=transform)
    assert np.array_equal(result, truth)

def test_shapes_internal_driver_manager():
    """Make sure this works if driver is managed outside this test"""

    image = numpy.zeros((20, 20), dtype=rasterio.ubyte)
    image[5:15, 5:15] = 127
    shapes = ftrz.shapes(image)
    shape, val = next(shapes)
    assert shape['type'] == 'Polygon'

def test_shapes_invalid_mask_dtype(basic_image):
    """A mask that is the wrong dtype should fail."""
    for dtype in ('int8', 'int16', 'int32'):
        with pytest.raises(ValueError):
            next(shapes(
                basic_image,
                mask=np.ones(basic_image.shape, dtype=dtype)
            ))