import bpy
import math
from mathutils import Vector

# Blender 4.0.2: Flat-Earth / azimuthal-equidistant disk with Orion stars.
# Distances are treated as kilometres (1 Blender unit = 1 km).
IMAGE_PATH = ""
DISK_RADIUS = 20000.0
GRID_Z = 60.0
LINE_RADIUS = 18.0
STAR_HEIGHT = 5000.0
STAR_RADIUS = 65.0

# Clear the current scene.
bpy.ops.object.select_all(action='SELECT')
bpy.ops.object.delete(use_global=False)
for datablocks in (bpy.data.meshes, bpy.data.curves, bpy.data.materials,
                   bpy.data.cameras, bpy.data.lights):
    for datablock in list(datablocks):
        if datablock.users == 0:
            datablocks.remove(datablock)

# Create circular disk mesh.
segments = 256
verts = [(0.0, 0.0, 0.0)]
faces = []
for i in range(segments):
    a = 2.0 * math.pi * i / segments
    verts.append((DISK_RADIUS * math.cos(a), DISK_RADIUS * math.sin(a), 0.0))
for i in range(segments):
    faces.append((0, i + 1, ((i + 1) % segments) + 1))
mesh = bpy.data.meshes.new("Earth disk mesh")
mesh.from_pydata(verts, [], faces)
mesh.update()
disk = bpy.data.objects.new("Earth disk | diameter 40000 km", mesh)
bpy.context.collection.objects.link(disk)

# UV coordinates for AE.jpg; rotate the map 90 degrees counter-clockwise.
uv_layer = mesh.uv_layers.new(name="MapUV")
for poly in mesh.polygons:
    for loop_index in poly.loop_indices:
        vi = mesh.loops[loop_index].vertex_index
        x, y, _ = mesh.vertices[vi].co
        u = 0.5 + x / (2.0 * DISK_RADIUS)
        v = 0.5 + y / (2.0 * DISK_RADIUS)
        uv_layer.data[loop_index].uv = (v, 1.0 - u)

mat = bpy.data.materials.new("Earth map material")
mat.use_nodes = True
nodes = mat.node_tree.nodes
links = mat.node_tree.links
bsdf = nodes.get("Principled BSDF")
if bsdf:
    bsdf.inputs["Roughness"].default_value = 1.0
try:
    image = bpy.data.images.load(IMAGE_PATH, check_existing=True)
    tex = nodes.new("ShaderNodeTexImage")
    tex.image = image
    links.new(tex.outputs["Color"], bsdf.inputs["Base Color"])
except Exception as e:
    print("Could not load map image:", IMAGE_PATH, e)
    mat.diffuse_color = (0.12, 0.35, 0.15, 1.0)
disk.data.materials.append(mat)

# Latitude/longitude grid.
grid_mat = bpy.data.materials.new("Latitude longitude grid material")
grid_mat.diffuse_color = (0.02, 0.02, 0.02, 1.0)
grid_mat.use_nodes = True
grid_bsdf = grid_mat.node_tree.nodes.get("Principled BSDF")
if grid_bsdf:
    grid_bsdf.inputs["Base Color"].default_value = (0.015, 0.015, 0.015, 1.0)
    grid_bsdf.inputs["Roughness"].default_value = 1.0

def add_polyline(name, points, cyclic=False, material=grid_mat):
    curve = bpy.data.curves.new(name + " curve", type='CURVE')
    curve.dimensions = '3D'
    curve.resolution_u = 1
    curve.bevel_depth = LINE_RADIUS
    curve.bevel_resolution = 1
    spline = curve.splines.new('POLY')
    spline.points.add(len(points) - 1)
    for pnt, co in zip(spline.points, points):
        pnt.co = (co[0], co[1], co[2], 1.0)
    spline.use_cyclic_u = cyclic
    obj = bpy.data.objects.new(name, curve)
    bpy.context.collection.objects.link(obj)
    obj.data.materials.append(material)
    return obj

# Latitude circles for an azimuthal-equidistant map centered on the North Pole.
for lat in [80, 60, 40, 20, 0, -20, -40, -60, -80]:
    r = DISK_RADIUS * (90.0 - lat) / 180.0
    pts = []
    for i in range(361):
        a = 2.0 * math.pi * i / 360.0
        pts.append((r * math.cos(a), r * math.sin(a), GRID_Z))
    add_polyline("Latitude %s deg" % lat, pts, cyclic=True)

# Longitude rays every 15 degrees.
for lon in range(0, 360, 15):
    a = math.radians(lon)
    pts = []
    for i in range(201):
        r = DISK_RADIUS * i / 200.0
        pts.append((r * math.cos(a), r * math.sin(a), GRID_Z))
    add_polyline("Longitude %s deg" % lon, pts)

# Orion's eight principal stars. Coordinates are the user-specified flat-model
# nadir points; every star is placed directly above its listed point at 5000 km.
STAR_SPECS = [
    ("Betelgeuse",       8.92, 120.98, (1.00, 0.18, 0.06)),
    ("Rigel",           -4.77, 116.00, (0.28, 0.55, 1.00)),
    ("Bellatrix",        7.93, 117.10, (0.65, 0.82, 1.00)),
    ("Mintaka",          2.53, 118.16, (0.75, 0.88, 1.00)),
    ("Alnilam",          1.80, 118.73, (0.75, 0.88, 1.00)),
    ("Alnitak",          1.21, 119.34, (0.75, 0.88, 1.00)),
    ("Saiph",           -5.68, 120.29, (0.45, 0.68, 1.00)),
    ("Meissa",           10.88, 118.44, (0.75, 0.88, 1.00)),
]

def ae_xy(latitude, longitude):
    # Match the map/grid convention in the supplied script.
    r = DISK_RADIUS * (90.0 - latitude) / 180.0
    a = math.radians(longitude)
    return (r * math.cos(a), r * math.sin(a))

def make_emission_material(name, color):
    material = bpy.data.materials.new(name)
    material.diffuse_color = (*color, 1.0)
    material.use_nodes = True
    bsdf = material.node_tree.nodes.get("Principled BSDF")
    if bsdf:
        bsdf.inputs["Base Color"].default_value = (*color, 1.0)
        if "Emission Color" in bsdf.inputs:
            bsdf.inputs["Emission Color"].default_value = (*color, 1.0)
            bsdf.inputs["Emission Strength"].default_value = 3.0
        elif "Emission" in bsdf.inputs:
            bsdf.inputs["Emission"].default_value = (*color, 1.0)
            if "Emission Strength" in bsdf.inputs:
                bsdf.inputs["Emission Strength"].default_value = 3.0
    return material

star_objects = []
for star_name, lat, lon, color in STAR_SPECS:
    x, y = ae_xy(lat, lon)
    bpy.ops.mesh.primitive_uv_sphere_add(
        segments=24, ring_count=12, radius=STAR_RADIUS,
        location=(x, y, STAR_HEIGHT))
    star = bpy.context.object
    star.name = "Orion | %s | nadir %.2fN %.2fE | altitude 5000 km" % (
        star_name, lat, lon)
    star.data.materials.append(make_emission_material(
        "Star emission | " + star_name, color))
    bpy.ops.object.shade_smooth()
    star_objects.append(star)

# Camera positions use the same AE map convention as the city coordinates.
# Every camera points directly at Alnilam at altitude 5000 km.
CITIES = [
    ("Manila",       14.5995, 120.9842),
    ("Tokyo",        35.6762, 139.6503),
    ("Perth",       -31.9505, 115.8605),
    ("Delhi",        28.6139,  77.2090),
    ("Irkutsk",      52.2869, 104.3050),
    ("Beijing",      39.9042, 116.4074),
    ("Auckland",    -36.8509, 174.7645),
]

target = next(obj.location.copy() for obj in star_objects if "Alnilam" in obj.name)
camera_objects = []
for city, lat, lon in CITIES:
    x, y = ae_xy(lat, lon)
    cam_data = bpy.data.cameras.new(city + " Observation Camera")
    cam = bpy.data.objects.new(city + " Observation Camera", cam_data)
    bpy.context.collection.objects.link(cam)
    cam.location = (x, y, 100.0)
    cam_data.type = 'PERSP'
    cam_data.lens = 24.0
    cam_data.clip_start = 0.1
    cam_data.clip_end = 100000.0
    direction = target - cam.location
    cam.rotation_euler = direction.to_track_quat('-Z', 'Y').to_euler()
    camera_objects.append(cam)

# Manila is the active camera by default; select another camera in the Outliner
# and set it as the scene camera to inspect another city.
scene = bpy.context.scene
scene.camera = next(cam for cam in camera_objects if cam.name.startswith("Manila"))
scene.frame_start = 1
scene.frame_end = 1
scene.frame_set(1)
scene.render.resolution_x = 1280
scene.render.resolution_y = 720
scene.render.resolution_percentage = 100

# Predictable viewport.
for screen in bpy.data.screens:
    for area in screen.areas:
        if area.type == 'VIEW_3D':
            space = area.spaces.active
            space.clip_start = 0.1
            space.clip_end = 100000.0
            space.region_3d.view_distance = 50000.0
            space.region_3d.view_location = (0.0, 0.0, 0.0)
            space.region_3d.view_perspective = 'ORTHO'
            space.shading.type = 'MATERIAL'

bpy.ops.object.select_all(action='DESELECT')
disk.select_set(True)
bpy.context.view_layer.objects.active = disk

print("Created AE disk/grid, 8 Orion stars at altitude 5000 km, and 7 city cameras all aimed at Alnilam.")
print("Active camera: Manila. Change scene.camera to another city camera for that view.")
print("Map image path:", IMAGE_PATH)