Landcover Classification

Landcover categories

# Define landcover categories
categories = """\
11|Water
24|Building
31|Bare
42|Forest
52|Shrub/Meadow
71|Grass
"""

Setup

# Import modules
import io
import sys
import subprocess
from pathlib import Path

# Find GRASS Python packages
sys.path.append(
  subprocess.check_output(
    ["grass", "--config", "python_path"],
    text=True
    ).strip()
  )

# Import GRASS packages
import grass.script as gs
import grass.jupyter as gj
from grass.tools import Tools

# Initialize GRASS project
gj.init("hilltop")
tools = Tools(overwrite=True)

# Print projection
projection = tools.g_proj(format="shell", flags="p").keyval
print(projection.get("name"))

NAD83(HARN) / Louisiana South

Set region

# Set region
tools.g_region(vector="parcel", res=0.3)

Classify landcover

# Classify ground
tools.r_mapcalc(expression="landcover = if(canopy, 31)")

# Classify buildings
tools.r_mapcalc(expression="landcover = if(canopy >= 4 && canopy < 10 && ndvi < 0, 24, landcover)")

# Classify trees
tools.r_mapcalc(expression="landcover = if(canopy >= 4 && ndvi >= 0, 42, landcover)")

# Classify shrubs & meadows
tools.r_mapcalc(expression="landcover = if(canopy < 4 && canopy >= 0.5 && ndvi > 0.05, 52, landcover)")

# Classify grasses
tools.r_mapcalc(expression="landcover = if(canopy < 0.5 && ndvi > 0.1, 71, landcover)")

# Classify water
tools.r_mapcalc(expression="landcover = if(isnull(lake), landcover, 11)")

# Set landcover colors
tools.r_colors(map="landcover", color="nlcd")

# Set landcover categories
tools.r_category(map="landcover", separator="pipe", rules=io.StringIO(categories))

# Visualize landcover
m = gj.Map(width=800)
m.d_rast(map="landcover")
m.d_legend(raster="landcover", color="white", at=(75, 95, 1, 3), use=(11,24,31,42,52,71))
m.show()
m.save("images/landcover-01.png")