package overlay

import (
	"math"
	"sort"
)

// Roads: the least-cost paths between the settlements that were just placed.
//
// A road is the one mark whose shape is not a judgement at all. Given where two towns are, the line between
// them is whatever the ground allows - up the valley, round the spur, across the saddle - and that is a
// shortest-path problem with a cost function, not a drawing. It is also the single most tedious thing to
// paint by hand, because getting it right means reading a heightmap pixel by pixel.
//
// Three decisions worth stating:
//
//   - **Water is impassable, so roads never swim.** Each landmass gets its own network. A bridge or a ferry
//     is a deliberate act and belongs to the author, and a generator that guessed at them would put a
//     motorway across a strait it has no idea is thirty kilometres wide.
//   - **A minimum spanning tree, not every pair.** Joining all pairs gives a cobweb; the tree gives exactly
//     enough road to reach everywhere, which is both what a road network minimally is and the thing an
//     author can most easily add to. Edges are weighted by path *cost*, not by straight-line distance, so
//     two towns either side of a range are correctly further apart than the map says.
//   - **It runs on a coarsened grid.** A road at the overlay's full resolution would be a Dijkstra over
//     twenty-nine million cells per settlement. The cost surface is smooth at the scale a road cares about,
//     so it is pooled to a few hundred cells across, solved there, and the resulting polyline is stamped
//     back at full resolution with the mark's real width.

// roadGrid is the coarsened cost surface the paths are solved on.
type roadGrid struct {
	w, h  int
	step  int       // overlay pixels per coarse cell
	cost  []float32 // per coarse cell, +Inf where impassable
	scale float64   // overlay pixels per coarse cell, as a float
}

func buildRoadGrid(d *genData, maxSlopeDeg float64) *roadGrid {
	in := d.in
	// About six hundred cells around the world: fine enough that a coarse cell is well under a kilometre on
	// any world this tool makes, coarse enough that fifty Dijkstras are a second's work.
	step := int(math.Max(1, math.Round(float64(in.W)/600)))
	gw := (in.W + step - 1) / step
	gh := (in.H + step - 1) / step
	g := &roadGrid{w: gw, h: gh, step: step, scale: float64(step), cost: make([]float32, gw*gh)}

	inf := float32(math.Inf(1))
	for gy := 0; gy < gh; gy++ {
		for gx := 0; gx < gw; gx++ {
			// Pool the block: any sea in it makes the cell water, because a road that clips a bay is a road
			// in the sea. The slope taken is the worst in the block, for the same reason.
			var worst float64
			wet := false
			for y := gy * step; y < (gy+1)*step && y < in.H; y++ {
				for x := gx * step; x < (gx+1)*step && x < in.W; x++ {
					i := y*in.W + x
					if in.Sea[i] {
						wet = true
						break
					}
					if s := float64(d.slopeDeg[i]); s > worst {
						worst = s
					}
				}
				if wet {
					break
				}
			}
			gi := gy*gw + gx
			switch {
			case wet:
				g.cost[gi] = inf
			case worst > maxSlopeDeg:
				g.cost[gi] = inf
			default:
				// Slope is what a road pays for. Quadratic rather than linear so that a route prefers a long
				// gentle way round to a short steep one, which is what a real road does.
				t := worst / math.Max(maxSlopeDeg, 1e-6)
				g.cost[gi] = float32(1 + 12*t*t)
			}
		}
	}
	return g
}

func (g *roadGrid) idx(x, y int) int { return y*g.w + x }

// dijkstra returns the cost to every reachable coarse cell from a source, and the predecessor chain to walk
// a path back. A binary heap over a few hundred thousand cells; the graph is eight-connected and X wraps.
func (g *roadGrid) dijkstra(src int) (cost []float32, pred []int32) {
	n := g.w * g.h
	cost = make([]float32, n)
	pred = make([]int32, n)
	inf := float32(math.Inf(1))
	for i := range cost {
		cost[i] = inf
		pred[i] = -1
	}
	if math.IsInf(float64(g.cost[src]), 1) {
		return cost, pred
	}
	cost[src] = 0
	h := &costHeap{keys: []float32{0}, items: []int32{int32(src)}}
	for h.Len() > 0 {
		c := int(h.pop())
		cx, cy := c%g.w, c/g.w
		base := cost[c]
		for _, o := range neighbours8 {
			nx := (cx + o[0] + g.w) % g.w
			ny := cy + o[1]
			if ny < 0 || ny >= g.h {
				continue
			}
			n := g.idx(nx, ny)
			cc := g.cost[n]
			if math.IsInf(float64(cc), 1) {
				continue
			}
			// Diagonal steps cost their real length, or the network shows a bias along the axes.
			step := float32(1.0)
			if o[0] != 0 && o[1] != 0 {
				step = float32(math.Sqrt2)
			}
			next := base + cc*step
			if next < cost[n] {
				cost[n] = next
				pred[n] = int32(c)
				h.push(int32(n), next)
			}
		}
	}
	return cost, pred
}

// costHeap is a binary min-heap of coarse cells. Lazy deletion is not needed because a cell is only pushed
// when its cost strictly improves, and a stale entry pops with a cost no better than the settled one.
type costHeap struct {
	keys  []float32
	items []int32
}

func (h *costHeap) Len() int { return len(h.items) }

func (h *costHeap) push(item int32, key float32) {
	h.keys = append(h.keys, key)
	h.items = append(h.items, item)
	i := len(h.items) - 1
	for i > 0 {
		p := (i - 1) / 2
		if h.keys[p] <= h.keys[i] {
			break
		}
		h.keys[p], h.keys[i] = h.keys[i], h.keys[p]
		h.items[p], h.items[i] = h.items[i], h.items[p]
		i = p
	}
}

func (h *costHeap) pop() int32 {
	top := h.items[0]
	last := len(h.items) - 1
	h.keys[0], h.items[0] = h.keys[last], h.items[last]
	h.keys = h.keys[:last]
	h.items = h.items[:last]
	i := 0
	for {
		l := 2*i + 1
		if l >= last {
			break
		}
		if r := l + 1; r < last && h.keys[r] < h.keys[l] {
			l = r
		}
		if h.keys[l] >= h.keys[i] {
			break
		}
		h.keys[l], h.keys[i] = h.keys[i], h.keys[l]
		h.items[l], h.items[i] = h.items[i], h.items[l]
		i = l
	}
	return top
}

// paintRoads traces a spanning tree over the settlements of each landmass and stamps it.
func (l *Legend) paintRoads(m *Mark, d *genData, out *Raster, idx uint8, placed []Placed) (int, int) {
	if len(placed) < 2 {
		return 0, 0
	}
	g := m.Generate
	maxSlope := g.MaxSlopeDeg
	if maxSlope <= 0 {
		maxSlope = 22
	}
	widthM := g.WidthM
	if widthM <= 0 {
		widthM = m.WidthM
	}
	if widthM <= 0 {
		widthM = 8
	}
	// A road eight metres wide is less than one overlay pixel at 12.9 m, and a mark thinner than a pixel is
	// not a mark. It is painted at least one pixel wide and the true width travels in the legend, which is
	// exactly how `width_m` is meant to be read.
	halfPx := int(math.Max(0, math.Round(widthM/d.in.CellM/2)))

	rg := buildRoadGrid(d, maxSlope)

	// Settlements grouped by landmass: a spanning tree per island, never between islands.
	byRegion := map[int][]int{}
	for i, p := range placed {
		if p.Region < 0 {
			continue
		}
		byRegion[p.Region] = append(byRegion[p.Region], i)
	}
	regions := make([]int, 0, len(byRegion))
	for r := range byRegion {
		regions = append(regions, r)
	}
	sort.Ints(regions) // deterministic order, cross-cutting rule 12

	total, pieces := 0, 0
	for _, r := range regions {
		members := byRegion[r]
		if len(members) < 2 {
			continue
		}
		total += l.connectRegion(rg, d, out, idx, placed, members, halfPx, &pieces)
	}
	return total, pieces
}

// connectRegion solves the paths among one landmass's settlements and stamps its spanning tree.
func (l *Legend) connectRegion(rg *roadGrid, d *genData, out *Raster, idx uint8,
	placed []Placed, members []int, halfPx int, pieces *int) int {

	n := len(members)
	src := make([]int, n)
	for k, pi := range members {
		p := placed[pi]
		gx := (p.X / rg.step) % rg.w
		gy := p.Y / rg.step
		if gy >= rg.h {
			gy = rg.h - 1
		}
		src[k] = rg.idx(gx, gy)
	}

	// One Dijkstra per settlement, kept: the coarse grid is a few hundred thousand cells and a landmass has
	// a handful of towns, so holding the predecessor chains costs a few megabytes and saves solving twice.
	costs := make([][]float32, n)
	preds := make([][]int32, n)
	for k := range members {
		costs[k], preds[k] = rg.dijkstra(src[k])
	}

	// Prim's, on path cost. Unreachable pairs are skipped, so a landmass whose towns are separated by ground
	// too steep for a road comes out as two networks rather than one impossible line.
	inTree := make([]bool, n)
	inTree[0] = true
	painted := 0
	for added := 1; added < n; added++ {
		bestA, bestB := -1, -1
		best := float32(math.Inf(1))
		for a := 0; a < n; a++ {
			if !inTree[a] {
				continue
			}
			for b := 0; b < n; b++ {
				if inTree[b] {
					continue
				}
				if c := costs[a][src[b]]; c < best {
					best, bestA, bestB = c, a, b
				}
			}
		}
		if bestA < 0 || math.IsInf(float64(best), 1) {
			break // nothing else on this landmass is reachable by road
		}
		inTree[bestB] = true
		painted += stampPath(rg, preds[bestA], src[bestA], src[bestB], d, out, idx, halfPx)
		*pieces++
	}
	return painted
}

// stampPath walks the predecessor chain back from dst to src and paints it at full resolution.
func stampPath(rg *roadGrid, pred []int32, src, dst int, d *genData, out *Raster, idx uint8, halfPx int) int {
	var chain []int
	for c := dst; c >= 0; {
		chain = append(chain, c)
		if c == src {
			break
		}
		p := pred[c]
		if p < 0 {
			return 0 // no route; leave the ground unpainted rather than drawing a guess
		}
		c = int(p)
	}
	painted := 0
	for k := 0; k+1 < len(chain); k++ {
		ax, ay := coarseCentre(rg, chain[k])
		bx, by := coarseCentre(rg, chain[k+1])
		painted += stampSegment(out, d, ax, ay, bx, by, halfPx, idx)
	}
	return painted
}

func coarseCentre(rg *roadGrid, c int) (int, int) {
	gx, gy := c%rg.w, c/rg.w
	return gx*rg.step + rg.step/2, gy*rg.step + rg.step/2
}

// stampSegment draws one straight run between two coarse-cell centres, wrapping in X the short way so a road
// crossing the seam is one road rather than a line back across the whole map.
func stampSegment(out *Raster, d *genData, ax, ay, bx, by, halfPx int, idx uint8) int {
	dx := wrapDelta(bx-ax, out.W)
	dy := by - ay
	steps := int(math.Max(math.Abs(float64(dx)), math.Abs(float64(dy))))
	if steps == 0 {
		return stampDisc(out, d, ax, ay, halfPx, idx, true)
	}
	painted := 0
	for s := 0; s <= steps; s++ {
		t := float64(s) / float64(steps)
		x := ax + int(math.Round(float64(dx)*t))
		y := ay + int(math.Round(float64(dy)*t))
		if y < 0 || y >= out.H {
			continue
		}
		painted += stampDisc(out, d, x, y, halfPx, idx, true)
	}
	return painted
}