fio-l2-cairo-mixedgeom-split

Rerun

import geopandas as gpd

# Load the data
gdf = gpd.read_file('heritage_sites.geojson')

# Check geometry types and columns
print("Columns:", gdf.columns.tolist())
print("Geometry types count:")
print(gdf.geometry.type.value_counts())
print("\nFirst few rows:")
print(gdf.head())

import geopandas as gpd

def solve():
    # Load the input file
    input_file = 'heritage_sites.geojson'
    output_file = 'heritage.gpkg'
    
    # Egypt National Grid is EPSG:22992 (Egypt Red Belt) 
    # or UTM zone 36N (EPSG:32636) or 35N (EPSG:32635). 
    # However, "Egypt's national grid" usually refers to the Egypt Red Belt or similar.
    # Let's use EPSG:22992. If that's not standard in common libraries, 
    # I will fallback to UTM 36N (EPSG:32636) which is widely used in Egypt.
    # But I'll try 22992 first.
    target_crs = "EPSG:22992"
    
    try:
        gdf = gpd.read_file(input_file)
    except Exception as e:
        print(f"Error reading file: {e}")
        return

    # Ensure all geometries are single-part
    # explode() converts MultiPolygons/MultiLineStrings/MultiPoints into single-part features
    gdf_exploded = gdf.explode(index_parts=False)

    # Define a helper to filter and transform
    def process_layer(geometry_type):
        # Filter by geometry type
        # Note: explode() might result in GeometryCollection if we aren't careful, 
        # but usually we want to filter the resulting single-part geometries.
        
        # Filter based on the geometry type of the exploded features
        layer_gdf = gdf_exploded[gdf_exploded.geometry.type.isin(geometry_type)]
        
        if not layer_gdf.empty:
            # Reproject
            layer_gdf = layer_gdf.to_crs(target_crs)
            return layer_gdf
        return None

    # Create the layers
    # Points: Point
    # Lines: LineString
    # Polygons: Polygon
    
    # Note: We use the intersection of the type list to handle different levels
    points_layer = process_layer(['Point'])
    lines_layer = process_layer(['LineString'])
    polygons_layer = process_layer(['Polygon'])

    # Write to GPKG
    # We use engine='pyogrio' if available for better performance, but standard fiona is safer.
    # We must write each layer to the same file but with different layer names.
    
    import fiona
    from sqlalchemy import create_engine

    # Using geopandas write_to_file with mode='a' to add layers
    # First, if the file exists, remove it to start fresh
    import os
    if os.path.exists(output_file):
        os.remove(output_file)

    if points_layer is not None:
        points_layer.to_file(output_file, layer='points', driver="GPKG")
        print("Layer 'points' written.")
    
    if lines_layer is not None:
        lines_layer.to_file(output_file, layer='lines', driver="GPKG", mode='a')
        print("Layer 'lines' written.")
        
    if polygons_layer is not None:
        polygons_layer.to_file(output_file, layer='polygons', driver="GPKG", mode='a')
        print("Layer 'polygons' written.")

    print(f"Successfully created {output_file}")

if __name__ == "__main__":
    solve()

import geopandas as gpd
import os

def solve():
    # Load the input file
    input_file = 'heritage_sites.geojson'
    output_file = 'heritage.gpkg'
    
    # Egypt's national grid is often referred to as Egypt Red Belt (EPSG:22992).
    # However, many practitioners use UTM 36N (EPSG:32636) for Egypt.
    # I will use EPSG:22992 as requested (Egypt's national grid).
    target_crs = "EPSG:22992"
    
    try:
        gdf = gpd.read_file(input_file)
    except Exception as e:
        print(f"Error reading file: {e}")
        return

    # Ensure all geometries are single-part
    # explode() converts MultiPolygons/MultiLineStrings/MultiPoints into single-part features
    gdf_exploded = gdf.explode(index_parts=False).reset_index(drop=True)

    # Define a helper to filter and transform
    def process_layer(geometry_types):
        # Filter based on the geometry type of the exploded features
        layer_gdf = gdf_exploded[gdf_exploded.geometry.type.isin(geometry_types)].copy()
        
        if not layer_gdf.empty:
            # Reproject to Egypt National Grid
            layer_gdf = layer_gdf.to_crs(target_crs)
            return layer_gdf
        return None

    # Create the layers
    # Points: Point
    # Lines: LineString
    # Polygons: Polygon
    points_layer = process_layer(['Point'])
    lines_layer = process_layer(['LineString'])
    polygons_layer = process_layer(['Polygon'])

    # Write to GPKG
    # First, if the file exists, remove it to start fresh
    if os.path.exists(output_file):
        os.remove(output_file)

    # Use write_to_file. For the first layer, we don't need mode='a'.
    # For subsequent layers, we use mode='a' (append).
    
    first_layer = True

    if points_layer is not None:
        mode = 'w' if first_layer else 'a'
        points_layer.to_file(output_file, layer='points', driver="GPKG", mode=mode)
        print("Layer 'points' written.")
        first_layer = False
    
    if lines_layer is not None:
        mode = 'w' if first_layer else 'a'
        lines_layer.to_file(output_file, layer='lines', driver="GPKG", mode=mode)
        print("Layer 'lines' written.")
        first_layer = False
        
    if polygons_layer is not None:
        mode = 'w' if first_layer else 'a'
        polygons_layer.to_file(output_file, layer='polygons', driver="GPKG", mode=mode)
        print("Layer 'polygons' written.")

    print(f"Successfully created {output_file}")

if __name__ == "__main__":
    solve()

import geopandas as gpd

def verify():
    output_file = 'heritage.gpkg'
    
    try:
        # Read layers individually to verify they exist and have correct properties
        points = gpd.read_file(output_file, layer='points')
        lines = gpd.read_file(output_file, layer='lines')
        polygons = gpd.read_file(output_file, layer='polygons')
        
        print("--- Points Layer ---")
        print(f"Count: {len(points)}")
        print(f"CRS: {points.crs}")
        print(f"Columns: {points.columns.tolist()}")
        print(f"Geometry Type: {points.geometry.type.unique()}")
        
        print("\n--- Lines Layer ---")
        print(f"Count: {len(lines)}")
        print(f"CRS: {lines.crs}")
        print(f"Columns: {lines.columns.tolist()}")
        print(f"Geometry Type: {lines.geometry.type.unique()}")

        print("\n--- Polygons Layer ---")
        print(f"Count: {len(polygons)}")
        print(f"CRS: {polygons.crs}")
        print(f"Columns: {polygons.columns.tolist()}")
        print(f"Geometry Type: {polygons.geometry.type.unique()}")

        # Check for site_id
        if 'site_id' in points.columns and 'site_id' in lines.columns and 'site_id' in polygons.columns:
            print("\nSuccess: 'site_id' present in all layers.")
        else:
            print("\nError: 'site_id' missing in one or more layers.")

    except Exception as e:
        print(f"Verification failed: {e}")

if __name__ == "__main__":
    verify()