dd-l2-bangkok-multicrs-audit

Rerun

     1	#!/usr/bin/env python3
     2	"""Audit every layer in bangkok_contractor_delivery.gpkg and write crs_audit.csv."""
     3	
     4	import csv
     5	import fiona
     6	import pyproj
     7	
     8	GPKG = "/work/bangkok_contractor_delivery.gpkg"
     9	OUT = "/work/crs_audit.csv"
    10	
    11	# Heuristic: check if text attributes contain bytes that look like latin1-mojibake
    12	# of UTF-8 encoded Thai.  Thai UTF-8 bytes in the range 0xE0-0xFB when
    13	# re-interpreted as latin1 produce sequences like à¹, à¸, etc.
    14	def detect_encoding(props: dict) -> str:
    15	    """Return 'latin1-mojibake' if any string value shows the pattern, else 'utf-8'."""
    16	    for k, v in props.items():
    17	        if isinstance(v, str):
    18	            # Check for common mojibake patterns: sequences like à¹\x9b, à¸\x99 etc.
    19	            # These are the latin1 rendering of UTF-8 Thai bytes.
    20	            # Quick heuristic: look for 'à' followed by a non-ASCII char in the
    21	            # range 0x80-0xBF (which is what latin1 shows for continuation bytes).
    22	            for i, ch in enumerate(v):
    23	                if ch == 'à' and i + 1 < len(v):
    24	                    nxt = v[i + 1]
    25	                    if 0x80 <= ord(nxt) <= 0xBF:
    26	                        return "latin1-mojibake"
    27	                # Also check for other common lead bytes: â, å, ç, è, é, ê, ë, ì, í, î, ï, ñ, ò, ó, ô, õ, ö, ù, ú, û, ü
    28	                if ch in 'âåçèéêëìíîïñòóôõöùúûü' and i + 1 < len(v):
    29	                    nxt = v[i + 1]
    30	                    if 0x80 <= ord(nxt) <= 0xBF:
    31	                        return "latin1-mojibake"
    32	    return "utf-8"
    33	
    34	
    35	rows = []
    36	
    37	layers = list(fiona.listlayers(GPKG))
    38	print(f"Found layers: {layers}")
    39	
    40	for layer_name in layers:
    41	    with fiona.open(GPKG, layer=layer_name) as src:
    42	        crs = src.crs
    43	        schema = src.schema
    44	        geom_type = schema.get("geometry", "Unknown")
    45	        feature_count = len(src)
    46	
    47	        # Format CRS as EPSG:NNNN
    48	        if crs:
    49	            try:
    50	                crs_obj = pyproj.CRS(crs)
    51	                if crs_obj.to_epsg():
    52	                    declared_crs = f"EPSG:{crs_obj.to_epsg()}"
    53	                else:
    54	                    declared_crs = str(crs)
    55	            except Exception:
    56	                declared_crs = str(crs)
    57	        else:
    58	            declared_crs = "None"
    59	
    60	        # Get a sample coordinate
    61	        sample_x = ""
    62	        sample_y = ""
    63	        encoding_detected = "utf-8"
    64	
    65	        if feature_count > 0:
    66	            feat = next(iter(src))
    67	            geom = feat.get("geometry")
    68	            if geom and geom.get("coordinates"):
    69	                coords = geom["coordinates"]
    70	                # For Polygon/MultiPolygon, drill down to first ring
    71	                if geom["type"] in ("Polygon",):
    72	                    coords = coords[0]
    73	                elif geom["type"] in ("MultiPolygon",):
    74	                    coords = coords[0][0]
    75	                elif geom["type"] in ("MultiLineString",):
    76	                    coords = coords[0]
    77	                elif geom["type"] in ("MultiPoint",):
    78	                    coords = coords[0]
    79	
    80	                if coords and len(coords) > 0:
    81	                    pt = coords[0]
    82	                    if isinstance(pt, (list, tuple)) and len(pt) >= 2:
    83	                        sample_x = pt[0]
    84	                        sample_y = pt[1]
    85	
    86	            # Detect encoding from text attributes
    87	            encoding_detected = detect_encoding(feat.get("properties", {}))
    88	
    89	        rows.append({
    90	            "layer_name": layer_name,
    91	            "declared_crs": declared_crs,
    92	            "geometry_type": geom_type,
    93	            "feature_count": feature_count,
    94	            "sample_x": sample_x,
    95	            "sample_y": sample_y,
    96	            "encoding_detected": encoding_detected,
    97	        })
    98	
    99	        print(f"  {layer_name}: CRS={declared_crs}, geom={geom_type}, count={feature_count}, "
   100	              f"sample=({sample_x}, {sample_y}), encoding={encoding_detected}")
   101	
   102	# Write CSV
   103	fieldnames = ["layer_name", "declared_crs", "geometry_type", "feature_count",
   104	              "sample_x", "sample_y", "encoding_detected"]
   105	with open(OUT, "w", newline="") as f:
   106	    writer = csv.DictWriter(f, fieldnames=fieldnames)
   107	    writer.writeheader()
   108	    writer.writerows(rows)
   109	
   110	print(f"\nWrote {len(rows)} rows to {OUT}")

     1	layer_name,declared_crs,geometry_type,feature_count,sample_x,sample_y,encoding_detected
     2	parcels,EPSG:24047,Polygon,4000,657421.5158199039,1509162.3701859994,latin1-mojibake
     3	roads,EPSG:32647,LineString,5000,656835.6850029442,1509487.0031804747,latin1-mojibake
     4	markets,EPSG:4326,Point,1000,100.4533,13.655000000000001,utf-8