import geopandas as gpd

# --- Input shapefile ---
shapefile_path = "ne_10m_populated_places.shp"
output_geojson = "world_cities_europe_clean.geojson"

# --- Load cities shapefile ---
cities = gpd.read_file(shapefile_path)

# --- Filter Europe + Nordics ---
europe_iso3 = [
    'AUT','BEL','CHE','DEU','DNK','ESP','FIN','FRA',
    'GBR','ITA','NOR','NLD','SWE','ISL','IRL','LUX'
]
cities = cities[cities['ADM0_A3'].isin(europe_iso3)]

# --- Ensure columns ---
if 'POP_MAX' not in cities.columns:
    cities['POP_MAX'] = 0
else:
    cities['POP_MAX'] = cities['POP_MAX'].fillna(0).astype(int)

cities = cities[['NAME','POP_MAX','geometry']]

# --- Remove empty geometries ---
cities = cities[~cities['geometry'].is_empty]

# --- Save to GeoJSON ---
cities.to_file(output_geojson, driver="GeoJSON")
print(f"[+] Saved {len(cities)} cities to {output_geojson}")