Initial commit: FastMCP server for BLM PLSS queries

Query Section/Township/Range from lat/long coordinates using
the BLM National Cadastral ArcGIS REST service.

Tools:
- get_plss_location: human-readable PLSS description
- get_plss_details: structured data with full metadata

.gitignore

@@ -0,0 +1,40 @@
+# Python
+__pycache__/
+*.py[cod]
+*$py.class
+*.so
+.Python
+build/
+develop-eggs/
+dist/
+downloads/
+eggs/
+.eggs/
+lib/
+lib64/
+parts/
+sdist/
+var/
+wheels/
+*.egg-info/
+.installed.cfg
+*.egg
+
+# Virtual environments
+.venv/
+venv/
+ENV/
+
+# IDE
+.idea/
+.vscode/
+*.swp
+*.swo
+
+# uv
+.python-version
+uv.lock
+
+# OS
+.DS_Store
+Thumbs.db

README.md

@@ -0,0 +1,32 @@
+# mcblmplss
+
+FastMCP server for querying BLM Public Land Survey System (PLSS) data by coordinates.
+
+## Usage
+
+```bash
+# Install and run
+uvx mcblmplss
+
+# Add to Claude Code
+claude mcp add plss "uvx mcblmplss"
+```
+
+## Tools
+
+- **get_plss_location** - Get Section/Township/Range for lat/long coordinates
+- **get_plss_details** - Get full structured PLSS data
+
+## Example
+
+```
+> get_plss_location(latitude=40.0, longitude=-105.0)
+
+Section 10, Township 1N, Range 68W, 6th Meridian
+State: CO
+PLSS ID: CO060010N0680W0SN100
+```
+
+## Coverage
+
+PLSS data is available for 30 states where federal land surveys were conducted (western/midwestern US). Not available for eastern seaboard states or Texas.

pyproject.toml

@@ -0,0 +1,38 @@
+[project]
+name = "mcblmplss"
+version = "2024.12.03"
+description = "FastMCP server for querying BLM Public Land Survey System (PLSS) data by coordinates"
+readme = "README.md"
+requires-python = ">=3.11"
+authors = [
+    {name = "Ryan Malloy", email = "ryan@supported.systems"}
+]
+keywords = ["mcp", "fastmcp", "blm", "plss", "cadastral", "land-survey", "gis"]
+classifiers = [
+    "Development Status :: 4 - Beta",
+    "Intended Audience :: Developers",
+    "Programming Language :: Python :: 3.11",
+    "Programming Language :: Python :: 3.12",
+    "Topic :: Scientific/Engineering :: GIS",
+]
+dependencies = [
+    "fastmcp>=2.13.2",
+    "httpx>=0.28.1",
+]
+
+[project.scripts]
+mcblmplss = "mcblmplss:main"
+
+[build-system]
+requires = ["hatchling"]
+build-backend = "hatchling.build"
+
+[tool.hatch.build.targets.wheel]
+packages = ["src/mcblmplss"]
+
+[tool.ruff]
+line-length = 100
+target-version = "py311"
+
+[tool.ruff.lint]
+select = ["E", "F", "I", "UP"]

src/mcblmplss/__init__.py

@@ -0,0 +1,10 @@
+"""
+mcblmplss - FastMCP server for BLM Public Land Survey System queries.
+
+Query PLSS data (Township, Range, Section) from coordinates using the
+BLM National Cadastral API.
+"""
+
+from mcblmplss.server import main, mcp
+
+__all__ = ["main", "mcp"]

src/mcblmplss/server.py

@@ -0,0 +1,242 @@
+"""
+FastMCP server for querying BLM Public Land Survey System (PLSS) data.
+
+Converts lat/long coordinates to Section, Township, Range using the
+BLM National Cadastral ArcGIS REST service.
+"""
+
+import httpx
+from fastmcp import FastMCP
+from pydantic import BaseModel, Field
+
+# BLM Cadastral MapServer endpoint
+BLM_PLSS_URL = "https://gis.blm.gov/arcgis/rest/services/Cadastral/BLM_Natl_PLSS_CadNSDI/MapServer"
+
+# Layer IDs in the MapServer
+LAYER_TOWNSHIP = 1
+LAYER_SECTION = 2
+
+
+class PLSSLocation(BaseModel):
+    """Public Land Survey System location description."""
+
+    section: str | None = Field(None, description="Section number (1-36)")
+    township: str = Field(..., description="Township number with direction (e.g., '4N')")
+    range: str = Field(..., description="Range number with direction (e.g., '6E')")
+    principal_meridian: str = Field(..., description="Principal meridian name")
+    state: str = Field(..., description="State abbreviation")
+    label: str = Field(..., description="Human-readable PLSS label")
+    plss_id: str = Field(..., description="Unique PLSS identifier")
+
+
+class PLSSQueryResult(BaseModel):
+    """Result from a PLSS coordinate query."""
+
+    latitude: float
+    longitude: float
+    township: PLSSLocation | None = None
+    section: PLSSLocation | None = None
+    error: str | None = None
+
+
+def parse_township_attributes(attrs: dict) -> PLSSLocation:
+    """Parse township attributes from BLM API response."""
+    township_num = attrs.get("TWNSHPNO", "").lstrip("0") or "0"
+    township_dir = attrs.get("TWNSHPDIR", "")
+    range_num = attrs.get("RANGENO", "").lstrip("0") or "0"
+    range_dir = attrs.get("RANGEDIR", "")
+
+    return PLSSLocation(
+        section=None,
+        township=f"{township_num}{township_dir}",
+        range=f"{range_num}{range_dir}",
+        principal_meridian=attrs.get("PRINMER", "Unknown"),
+        state=attrs.get("STATEABBR", ""),
+        label=attrs.get("TWNSHPLAB", ""),
+        plss_id=attrs.get("PLSSID", ""),
+    )
+
+
+def parse_section_attributes(attrs: dict, township: PLSSLocation | None) -> PLSSLocation:
+    """Parse section attributes from BLM API response.
+
+    Note: Section layer uses human-readable field names like "First Division Number"
+    while Township layer uses abbreviated names like "TWNSHPNO".
+    """
+    # Try both field name formats (API returns friendly names for section layer)
+    section_num = (
+        attrs.get("First Division Number", "") or attrs.get("FRSTDIVNO", "")
+    ).lstrip("0") or "0"
+
+    # Get the division ID (also uses friendly name)
+    division_id = attrs.get("First Division Identifier", "") or attrs.get("FRSTDIVID", "")
+
+    if township:
+        return PLSSLocation(
+            section=section_num,
+            township=township.township,
+            range=township.range,
+            principal_meridian=township.principal_meridian,
+            state=township.state,
+            label=f"Section {section_num}, {township.label}",
+            plss_id=division_id,
+        )
+
+    # Fallback: extract state from Township Identifier if no township data
+    township_id = attrs.get("Township Identifier", "") or attrs.get("PLSSID", "")
+    return PLSSLocation(
+        section=section_num,
+        township="Unknown",
+        range="Unknown",
+        principal_meridian="Unknown",
+        state=township_id[:2] if len(township_id) >= 2 else "",
+        label=f"Section {section_num}",
+        plss_id=division_id,
+    )
+
+
+async def query_plss(latitude: float, longitude: float) -> PLSSQueryResult:
+    """
+    Query BLM PLSS API for location at given coordinates.
+
+    Uses the ArcGIS REST identify operation to find township and section
+    features that contain the specified point.
+    """
+    # Build identify request parameters
+    # Note: BLM service accepts WGS84 (EPSG:4326) coordinates directly
+    params = {
+        "f": "json",
+        "geometry": f"{longitude},{latitude}",
+        "geometryType": "esriGeometryPoint",
+        "sr": "4326",  # WGS84
+        "layers": f"all:{LAYER_TOWNSHIP},{LAYER_SECTION}",
+        "tolerance": "1",
+        "mapExtent": f"{longitude-1},{latitude-1},{longitude+1},{latitude+1}",
+        "imageDisplay": "100,100,96",
+        "returnGeometry": "false",
+    }
+
+    async with httpx.AsyncClient(timeout=30.0) as client:
+        response = await client.get(f"{BLM_PLSS_URL}/identify", params=params)
+        response.raise_for_status()
+        data = response.json()
+
+    results = data.get("results", [])
+
+    if not results:
+        return PLSSQueryResult(
+            latitude=latitude,
+            longitude=longitude,
+            error="No PLSS data found for this location. The point may be outside surveyed areas.",
+        )
+
+    # Parse results by layer
+    township_data: PLSSLocation | None = None
+    section_data: PLSSLocation | None = None
+
+    for result in results:
+        layer_id = result.get("layerId")
+        attrs = result.get("attributes", {})
+
+        if layer_id == LAYER_TOWNSHIP and township_data is None:
+            township_data = parse_township_attributes(attrs)
+        elif layer_id == LAYER_SECTION and section_data is None:
+            section_data = parse_section_attributes(attrs, township_data)
+
+    return PLSSQueryResult(
+        latitude=latitude,
+        longitude=longitude,
+        township=township_data,
+        section=section_data,
+    )
+
+
+# Initialize FastMCP server
+mcp = FastMCP(
+    name="mcblmplss",
+    instructions="""
+    Query the U.S. Public Land Survey System (PLSS) using coordinates.
+
+    The PLSS divides land into:
+    - Townships: 6x6 mile squares identified by Township (N/S) and Range (E/W)
+    - Sections: 1 square mile (640 acres), numbered 1-36 within townships
+    - Principal Meridians: Reference lines for the survey system
+
+    This server queries the BLM National Cadastral database for authoritative PLSS data.
+    Coverage: 30 western and midwestern states (not eastern seaboard or Texas).
+    """,
+)
+
+
+@mcp.tool()
+async def get_plss_location(
+    latitude: float = Field(description="Latitude in decimal degrees (WGS84)"),
+    longitude: float = Field(description="Longitude in decimal degrees (WGS84)"),
+) -> str:
+    """
+    Get the Public Land Survey System (PLSS) location for coordinates.
+
+    Returns Section, Township, Range, and Principal Meridian for the given
+    latitude/longitude. Example: "Section 12, Township 4N, Range 6E, 6th PM"
+
+    Note: PLSS coverage is limited to 30 states where federal land surveys
+    were conducted. Eastern states and Texas use different systems.
+    """
+    result = await query_plss(latitude, longitude)
+
+    if result.error:
+        return f"Error: {result.error}"
+
+    if result.section:
+        return (
+            f"Section {result.section.section}, "
+            f"Township {result.section.township}, "
+            f"Range {result.section.range}, "
+            f"{result.section.principal_meridian}\n"
+            f"State: {result.section.state}\n"
+            f"PLSS ID: {result.section.plss_id}"
+        )
+    elif result.township:
+        return (
+            f"Township {result.township.township}, "
+            f"Range {result.township.range}, "
+            f"{result.township.principal_meridian}\n"
+            f"State: {result.township.state}\n"
+            f"(Section data not available for this location)"
+        )
+    else:
+        return "Unable to determine PLSS location for these coordinates."
+
+
+@mcp.tool()
+async def get_plss_details(
+    latitude: float = Field(description="Latitude in decimal degrees (WGS84)"),
+    longitude: float = Field(description="Longitude in decimal degrees (WGS84)"),
+) -> PLSSQueryResult:
+    """
+    Get detailed PLSS information as structured data.
+
+    Returns a PLSSQueryResult object with full township and section details
+    including PLSS IDs, principal meridian, and state information.
+
+    Useful when you need programmatic access to all PLSS fields.
+    """
+    return await query_plss(latitude, longitude)
+
+
+def main():
+    """Entry point for the mcblmplss MCP server."""
+    try:
+        from importlib.metadata import version
+
+        package_version = version("mcblmplss")
+    except Exception:
+        package_version = "dev"
+
+    print(f"🗺️  mcblmplss v{package_version} - BLM PLSS Query Server")
+    print("📍 Query Section/Township/Range from coordinates")
+    mcp.run()
+
+
+if __name__ == "__main__":
+    main()