UnrealPrototyping/Tools/Terrain/internal/template/coast.go

package template

import (
	"math"

	"salty/terrain/internal/dt"
	"salty/terrain/internal/field"
	"salty/terrain/internal/noise"
	"salty/terrain/internal/world"
)

// A painted coastline is a drawn line, and a coastline is not a drawn line.
//
// This is the Richardson paradox with a brush in it. An author draws a shore as a smooth curve, because that
// is what a hand and a bezier tool produce; a real coast has bays inside bays inside bays and the length you
// measure depends on the ruler you measure it with. Projected straight, the painting's own smoothness
// survives all the way to the heightmap, and the result reads as exactly what it is - a shape somebody drew -
// however good the erosion downstream is. `internal/coast` does this job on the square canvas, where the
// outline is noise to begin with; the painted planet had a manifest key for it, `coast_jitter_px`, which
// until now nothing anywhere read.
//
// **It is a mask on the waterline, not a warp of the painting.** That distinction was measured rather than
// reasoned. Displacing the point each cell asks the painting about - a domain warp - was tried first and it
// cannot cut a bay: a smooth warp of a smooth boundary is another smooth boundary, just wigglier, and at an
// amplitude large enough to fold it back on itself it drags every inland class boundary the same distance.
// What produces bays and headlands is thresholding a *signed distance field*: how far is this cell from the
// waterline, add fractal noise to that distance in metres, and ask again which side of zero it is on. Land
// juts out where the noise is positive and the sea reaches in where it is negative, at every scale the
// octaves cover, and nothing away from the shore moves at all.
//
// Two things follow from doing it this way, and both are the reason to:
//
//   - A cell that changes sides needs a class, and the distance transform already knows which one: it
//     returns the nearest seed cell as well as the distance to it, so new land takes the class of the land
//     it grew from and new sea takes the class of the water that came in. Sea eaten out of a shore becomes
//     the surf that was lying against it rather than deep ocean.
//   - Small islands have to survive. An islet thirty pixels across, under a noise field whose wavelength is
//     four hundred, sees very nearly a constant - so it either sits still or vanishes whole, and vanishing
//     whole is how an archipelago disappears between two runs. The amplitude is therefore capped per cell at
//     a fraction of the widest land within reach of it, which is a sliding maximum of the land distance.
//     A continent sees the full amplitude; an islet gets nibbled instead of deleted.

// Pass indices for the coast mask's noise, above the painted uplift path's 20..25 so neither can reshuffle
// the other.
const (
	srcCoastMask = 30
)

// islandGuard is how much of the widest land within reach the mask may eat. Two thirds leaves an islet
// recognisably itself while still giving it a ragged edge; at 1 it can take the whole thing.
const islandGuard = 0.66

// Coast is how the painted waterline is roughened before the painting is projected.
type Coast struct {
	// AmplitudePx is the furthest, in template pixels, that the shoreline may move. Zero switches the whole
	// thing off and the painting is used exactly as drawn.
	AmplitudePx float64

	// WavelengthPx is the coarsest octave: the width of the biggest bay it can cut. Octaves halve from
	// there, so the finest detail is this over 2^(Octaves-1). Bays come out about this wide and up to
	// AmplitudePx deep, so the ratio of the two is what decides whether the coast reads as a rough line or
	// as a fjord coast.
	WavelengthPx float64

	// Octaves and Gain are the fractal structure. A gain near 0.5 makes each scale about as prominent as the
	// last, which is the property a real coastline has and a single wobble does not.
	Octaves int
	Gain    float64

	// Scale is a per-pixel multiplier on AmplitudePx at the raster's own resolution, from the annotation
	// layer's coast_jitter marks. Nil is one everywhere, which is every world before D-57.
	//
	// It is what makes a hand-drawn coastline hold. The roughening exists because a drawn shore is smooth and
	// a real one is not, which is true of a shore nobody thought about and false of one somebody traced off a
	// map on purpose; a zero here pins that stretch exactly as painted while the rest of the world is still
	// roughened. Above one chews harder, which is the same knob pointed the other way - a fjord coast wants
	// more than the planet's own amplitude, not less.
	//
	// **A negative entry means the pixel carries no instruction**, which is not the same as one. A mark is a
	// stroke an author drew along a coastline and it lands on whichever side of the waterline their hand was
	// on; if an unmarked cell took the default, a stroke painted on the land would leave the water beside it
	// free to march inland anyway and the coast would move regardless. So an uninstructed cell takes the
	// instruction from the nearest cell on the other side of the waterline, which the distance transform
	// below has already found for a different reason. Painting either side is then enough, and painting over
	// the line - which is what a brush does - is enough twice over.
	Scale []float32

	Seed int64
}

// Amount reports whether this mask does anything.
func (c Coast) Amount() bool {
	return c.AmplitudePx > 0 && c.Octaves > 0 && c.WavelengthPx > 0 && c.Gain > 0
}

// maxScale is the largest multiplier any mark asks for, and at least 1. Uninstructed entries are negative and
// do not count; an unmarked world has no Scale at all and gets 1.
func (c Coast) maxScale() float64 {
	m := 1.0
	for _, v := range c.Scale {
		if float64(v) > m {
			m = float64(v)
		}
	}
	return m
}

// RoughenCoast returns the painting with its waterline displaced by fractal noise. The receiver is not
// modified: a caller that wants both keeps both, which is what the studio's preview does.
//
// It runs at the paint's own resolution rather than the planet's. That is a third of the cells, the mask's
// scales are quoted in template pixels anyway, and the thing being roughened is the painting - so a template
// re-exported at a different size is the one case where the coast moves, and that is already true of every
// other thing the painting decides.
func (r *Raster) RoughenCoast(l *Legend, p world.Planet, c Coast) *Raster {
	if !c.Amount() || r.W == 0 || r.H == 0 {
		return r
	}
	sea := make([]bool, len(l.Classes))
	for i := range l.Classes {
		sea[i] = l.Classes[i].Sea
	}

	isSea := make([]bool, len(r.Class))
	anySea, anyLand := false, false
	for i, cl := range r.Class {
		isSea[i] = sea[cl]
		if isSea[i] {
			anySea = true
		} else {
			anyLand = true
		}
	}
	if !anySea || !anyLand {
		return r // nothing to move: the painting is all one or all the other
	}

	// Signed distance to the waterline, positive on land, in template pixels, plus the index of the nearest
	// cell on the other side - which is where a cell that changes sides gets its class from.
	signed := make([]float32, len(r.Class))
	other := make([]int32, len(r.Class))

	// Seeded on land: for every sea cell, how far to land and which land cell.
	d2, near := dt.Transform(invert(isSea), r.W, r.H, true)
	for i := range signed {
		if isSea[i] {
			signed[i] = -float32(math.Sqrt(float64(d2[i])))
			other[i] = near[i]
		}
	}
	// Seeded on sea: for every land cell, how far to water and which water cell. Released in turn so the
	// two transforms are never both alive - at 29 million pixels each one is a quarter of a gigabyte.
	d2, near = dt.Transform(isSea, r.W, r.H, true)
	for i := range signed {
		if !isSea[i] {
			signed[i] = float32(math.Sqrt(float64(d2[i])))
			other[i] = near[i]
		}
	}
	d2, near = nil, nil

	// How wide the land is near each cell, so an islet cannot be eaten whole. Only land contributes, so a
	// lone islet reports its own half-width and not the open water around it.
	landOnly := field.New(r.W, r.H, 1)
	for i := range signed {
		if signed[i] > 0 {
			landOnly.Data[i] = signed[i]
		}
	}
	// The window is the furthest the shore could move *anywhere*, which is no longer the plain amplitude: a
	// mark asking for more than the planet's own can reach past it, and a guard measured over too small a
	// window would under-report how wide the land is and let an islet inside such a mark be eaten whole -
	// the one failure this guard exists to stop.
	reach := int(c.AmplitudePx*c.maxScale() + 0.5)
	widest := field.SlidingMax(landOnly, reach, true)

	// The noise, on world coordinates so it wraps at the seam and two runs of the same world agree.
	// noise.Lattice wraps modulo its cell count, so the lattice has to be a whole number of cells in the
	// noise period; the wavelength is quoted in template pixels and converts through the paint's own scale.
	metresPerPx := p.CircumferenceM() / float64(r.W)
	cells := int(p.NoisePeriodM/(c.WavelengthPx*metresPerPx) + 0.5)
	if cells < 1 {
		cells = 1
	}
	u, v := noise.WorldUV(r.W, r.H, metresPerPx, 0, 0, p.NoisePeriodM)
	n := noise.FBMAt(u, v, noise.NewSource(c.Seed, srcCoastMask),
		noise.Params{BaseCells: cells, Octaves: c.Octaves, Gain: c.Gain})

	// Stretched to its own full range before it is used, so the amplitude means what it says. An fBm stack
	// is normalised by the sum of its octave amplitudes, which is the value it would take if every octave
	// agreed at once - they never do, so the realised spread is far narrower than 0..1 and a nominal 48 px
	// was moving the shore about ten. The same trap as the massif fabric's threshold, and the same fix:
	// measure the distribution rather than assume it. Here it is one pass of min/max over the whole painting
	// - legitimate because the mask runs once on the whole map and not per region, so there is no second
	// caller to disagree with.
	n.Normalise()

	out := &Raster{W: r.W, H: r.H, Class: make([]uint8, len(r.Class))}
	copy(out.Class, r.Class)

	field.Rows(r.H, func(y0, y1 int) {
		for i := y0 * r.W; i < y1*r.W; i++ {
			amp := c.AmplitudePx
			if c.Scale != nil {
				sc := c.Scale[i]
				if sc < 0 {
					// Uninstructed: take the instruction from the far side of the waterline. See Coast.Scale.
					if j := other[i]; j >= 0 {
						sc = c.Scale[j]
					}
				}
				if sc >= 0 {
					amp *= float64(sc)
				}
			}
			if g := islandGuard * float64(widest.Data[i]); g < amp {
				amp = g
			}
			if amp <= 0 {
				continue
			}
			d := float64(signed[i]) + amp*(float64(n.Data[i])*2-1)
			nowLand := d > 0
			if nowLand == !isSea[i] {
				continue // this cell did not change sides
			}
			// It did. Take the class of the nearest cell on the side it has joined, which the transform
			// already found: land grows out of the land beside it, and water comes in as the water that
			// was lying against the shore.
			if j := other[i]; j >= 0 {
				out.Class[i] = r.Class[j]
			}
		}
	})
	return out
}

func invert(b []bool) []bool {
	out := make([]bool, len(b))
	for i, v := range b {
		out[i] = !v
	}
	return out
}