Tree Mapping with SAMGeo and Segment Anything Model 2 (SAM 2)¶
This notebook shows how to segment trees from aerial imagery with the Segment Anything Model 2 (SAM 2).
Make sure you use GPU runtime for this notebook. For Google Colab, go to Runtime -> Change runtime type and select GPU as the hardware accelerator.
Install dependencies¶
Uncomment and run the following cell to install the required dependencies.
# %pip install segment-geospatial
import leafmap
from samgeo import SamGeo2
Create an interactive map¶
m = leafmap.Map(center=[-22.17615, -51.253043], zoom=18, height="800px")
m.add_basemap("SATELLITE")
m
Download a sample image¶
Pan and zoom the map to select the area of interest. Use the draw tools to draw a polygon or rectangle on the map
bbox = m.user_roi_bounds()
if bbox is None:
bbox = [-51.2565, -22.1777, -51.2512, -22.175]
image = "Image.tif"
leafmap.map_tiles_to_geotiff(
output=image, bbox=bbox, zoom=19, source="Satellite", overwrite=True
)
You can also use your own image. Uncomment and run the following cell to use your own image.
# image = '/path/to/your/own/image.tif'
Display the downloaded image on the map.
m.layers[-1].visible = False
m.add_raster(image, layer_name="Image")
m
Initialize SAM class¶
Set automatic=False to enable the SAM2ImagePredictor.
sam = SamGeo2(
model_id="sam2-hiera-large",
automatic=False,
)
Specify the image to segment.
sam.set_image(image)
Display the map. Use the drawing tools to draw some rectangles around the features you want to extract, such as trees, buildings.
m
Create bounding boxes¶
If no rectangles are drawn, the default bounding boxes will be used as follows:
if m.user_rois is not None:
boxes = m.user_rois
else:
boxes = [
[-51.2546, -22.1771, -51.2541, -22.1767],
[-51.2538, -22.1764, -51.2535, -22.1761],
]
Segment the image¶
Use the predict() method to segment the image with specified bounding boxes. The boxes parameter accepts a list of bounding box coordinates in the format of [[left, bottom, right, top], [left, bottom, right, top], ...], a GeoJSON dictionary, or a file path to a GeoJSON file.
sam.predict(boxes=boxes, point_crs="EPSG:4326", output="mask.tif", dtype="uint8")
Display the result¶
Add the segmented image to the map.
m.add_raster("mask.tif", cmap="viridis", nodata=0, layer_name="Mask")
m
Use an existing vector dataset as box prompts¶
You can also use an existing vector dataset as box prompts. The following example uses an existing dataset of tree bounding boxes from GitHub.
geojson = (
"https://github.com/opengeos/datasets/releases/download/samgeo/tree_boxes.geojson"
)
Display the bounding boxes on the map.
m = leafmap.Map()
m.add_raster(image, layer_name="image")
style = {
"color": "#ffff00",
"weight": 2,
"fillColor": "#7c4185",
"fillOpacity": 0,
}
m.add_vector(
geojson,
style=style,
zoom_to_layer=True,
layer_name="Bounding boxes",
info_mode=None,
)
m
Segment trees with box prompts¶
Segment trees using the bounding boxes from the vector dataset.
output_masks = "mask2.tif"
sam.predict(boxes=geojson, point_crs="EPSG:4326", output=output_masks, dtype="uint8")
Display the segmented masks on the map.
m.add_raster(output_masks, nodata=0, opacity=0.5, layer_name="Tree masks")
Post-processing¶
You can use the region_groups() method to clean up the segmentation results, such as removing small regions, and filling holes. In addition, you can compute geometric properties of the regions, such as area, perimeter, eccentricity, and solidity.
out_image = "tree_masks.tif"
out_vector = "tree_vector.geojson"
array, gdf = sam.region_groups(
output_masks, min_size=200, out_vector=out_vector, out_image=out_image
)
gdf.head()
Display the cleaned masks¶
m = leafmap.Map()
m.add_raster(image, layer_name="Image")
style = {
"color": "#ffff00",
"weight": 2,
"fillColor": "#7c4185",
"fillOpacity": 0,
}
m.add_raster(
out_image, colormap="tab20", nodata=0, opacity=0.7, layer_name="Tree masks"
)
m.add_vector(out_vector, style=style, zoom_to_layer=True, layer_name="Tree vector")
m.add_vector(
geojson,
style={"color": "blue", "fillOpacity": 0},
layer_name="Bounding boxes",
info_mode=None,
)
m.add_layer_manager()
m
Create a split map¶
m = leafmap.Map()
m.add_raster(image, layer_name="Image")
m.split_map(
out_image,
image,
left_label="Tree masks",
right_label="Aerial imagery",
left_args={"colormap": "tab20", "nodata": 0, "opacity": 0.7},
)
m