1// Package align finds how two readings of one ground line up.2//3// Two captures of one map draw the same places at different scales and offsets:4// a scan of the official map and a rasterized in-game rendering do not share a5// pixel grid, and neither says how to reach the other. What they do share is6// places -- fast travel points, shops, landmarks -- pinned by name in both.7// Every name both readings pin, unambiguously, is a calibration point, and8// enough of them determine the affine transformation between the two pixel9// spaces deterministically: the same anchors always align the same way, with a10// residual report that says how well.11//12// The fit refuses rather than guesses. It stands only on names that are13// unambiguous in both readings, needs enough of them to mean something, trims14// its worst tenth twice to shed the places one side simply pinned sloppily, and15// declines to answer at all when what remains still does not close. A bad16// alignment drawn confidently is worse than none: it moves every feature it17// touches to somewhere plausible and wrong.18package align1920import (21 "fmt"22 "math"23 "regexp"24 "sort"25 "strings"26)2728// Anchor is one named point in a reading's own world pixels.29type Anchor struct {30 Title string31 X float6432 Y float6433}3435// Affine carries base = A·donor + t, spelled per axis: x' = AX·x + BX·y + CX,36// and y' likewise. The cross terms carry any small rotation or shear between37// the renderings; north-up readings leave them near zero.38type Affine struct {39 AX, BX, CX float6440 AY, BY, CY float6441}4243// Apply sends a point from the fitted reading's space into the base space.44func (a Affine) Apply(x, y float64) (float64, float64) {45 return a.AX*x + a.BX*y + a.CX, a.AY*x + a.BY*y + a.CY46}4748// Invert solves the transformation the other way. The false result is a49// degenerate transformation, which nothing may be inverted through.50func (a Affine) Invert() (Affine, bool) {51 det := a.AX*a.BY - a.BX*a.AY52 if math.Abs(det) < 1e-12 {53 return Affine{}, false54 }55 inv := Affine{56 AX: a.BY / det, BX: -a.BX / det,57 AY: -a.AY / det, BY: a.AX / det,58 }59 inv.CX = -(inv.AX*a.CX + inv.BX*a.CY)60 inv.CY = -(inv.AY*a.CX + inv.BY*a.CY)61 return inv, true62}6364// Scale reports how many base pixels one fitted pixel spans, averaged over the65// axes.66func (a Affine) Scale() float64 {67 return (math.Hypot(a.AX, a.AY) + math.Hypot(a.BX, a.BY)) / 268}6970// Report says what the fit stood on and how well it closed, in base pixels. Its71// string form is what a merge ledger records as the alignment, so two builds72// that aligned the same way say so in the same words.73type Report struct {74 Anchors int75 MedianPx float6476 P90Px float6477 WorstPx float6478}7980func (r Report) String() string {81 return fmt.Sprintf("%d anchors, median %.1fpx, p90 %.1fpx, worst %.1fpx",82 r.Anchors, r.MedianPx, r.P90Px, r.WorstPx)83}8485// The fit's own limits.86const (87 // MinAnchors is how many unambiguous shared names a fit needs before it88 // means anything.89 MinAnchors = 1290 // TrimRounds and TrimKeep are the trimming: two rounds, each keeping the91 // nine tenths that fit best, never dropping below MinAnchors.92 TrimRounds = 293 TrimKeep = 0.994 // MaxMedianPx is how far the middle anchor may sit from where the fit95 // predicts before the fit is refused as one that does not close.96 MaxMedianPx = 9697 // minPivot is the smallest pivot the elimination will divide by. Below it98 // the anchors are degenerate -- collinear, or all in one place -- and no99 // transformation stands.100 minPivot = 1e-12101)102103// Fit finds the affine sending donor pixels onto base pixels from the places104// both name. It is pure: the same anchors always yield the same answer.105func Fit(donor, base []Anchor) (Affine, Report, error) {106 matchedDonor, matchedBase := matchByTitle(donor, base)107 if len(matchedDonor) < MinAnchors {108 return Affine{}, Report{}, fmt.Errorf(109 "only %d unambiguous shared names; %d needed", len(matchedDonor), MinAnchors)110 }111112 keep := make([]int, len(matchedDonor))113 for index := range keep {114 keep[index] = index115 }116 var affine Affine117 for round := 0; ; round++ {118 var ok bool119 affine, ok = solve(matchedDonor, matchedBase, keep)120 if !ok {121 return Affine{}, Report{}, fmt.Errorf("anchors are degenerate; no transformation stands")122 }123 if round == TrimRounds {124 break125 }126 sort.Slice(keep, func(a, b int) bool {127 return residual(affine, matchedDonor[keep[a]], matchedBase[keep[a]]) <128 residual(affine, matchedDonor[keep[b]], matchedBase[keep[b]])129 })130 trimmed := int(float64(len(keep)) * TrimKeep)131 if trimmed < MinAnchors {132 trimmed = min(MinAnchors, len(keep))133 }134 keep = keep[:trimmed]135 }136137 residuals := make([]float64, len(keep))138 for at, index := range keep {139 residuals[at] = residual(affine, matchedDonor[index], matchedBase[index])140 }141 sort.Float64s(residuals)142 report := Report{143 Anchors: len(keep),144 MedianPx: residuals[len(residuals)/2],145 P90Px: residuals[int(float64(len(residuals))*0.9)],146 WorstPx: residuals[len(residuals)-1],147 }148 if report.MedianPx > MaxMedianPx {149 return Affine{}, report, fmt.Errorf("alignment does not close: %s", report)150 }151 return affine, report, nil152}153154func residual(a Affine, donor, base Anchor) float64 {155 x, y := a.Apply(donor.X, donor.Y)156 return math.Hypot(x-base.X, y-base.Y)157}158159// matchByTitle pairs anchors by normalized name, using only names that appear160// exactly once in each reading: a name two features share could pair either161// way, and a guessed pair is a corrupted measurement.162func matchByTitle(donor, base []Anchor) ([]Anchor, []Anchor) {163 donorByName := uniqueByName(donor)164 baseByName := uniqueByName(base)165 names := make([]string, 0, len(donorByName))166 for name := range donorByName {167 if _, shared := baseByName[name]; shared {168 names = append(names, name)169 }170 }171 sort.Strings(names)172 matchedDonor := make([]Anchor, len(names))173 matchedBase := make([]Anchor, len(names))174 for at, name := range names {175 matchedDonor[at] = donorByName[name]176 matchedBase[at] = baseByName[name]177 }178 return matchedDonor, matchedBase179}180181func uniqueByName(anchors []Anchor) map[string]Anchor {182 counts := make(map[string]int, len(anchors))183 held := make(map[string]Anchor, len(anchors))184 for _, anchor := range anchors {185 name := NormalizeTitle(anchor.Title)186 if name == "" {187 continue188 }189 counts[name]++190 held[name] = anchor191 }192 for name, count := range counts {193 if count > 1 {194 delete(held, name)195 }196 }197 return held198}199200var (201 fastTravelSuffix = regexp.MustCompile(`\s*[-–—]?\s*fast travel( point)?$`)202 nonAlphanumeric = regexp.MustCompile(`[^a-z0-9]+`)203)204205// NormalizeTitle folds the spellings apart readings give one place. The fast206// travel suffix goes because one reading titles the place and another titles207// the service at it, and they are the same anchor. Anything comparing features208// across readings compares through this, so a name means the same thing to the209// fit and to whoever consumes it.210func NormalizeTitle(title string) string {211 title = strings.ToLower(title)212 title = fastTravelSuffix.ReplaceAllString(title, "")213 return strings.TrimSpace(nonAlphanumeric.ReplaceAllString(title, " "))214}215216// Tokens is a normalized name's words, as a set.217func Tokens(normalized string) map[string]bool {218 tokens := make(map[string]bool)219 for _, token := range strings.Fields(normalized) {220 tokens[token] = true221 }222 return tokens223}224225// solve fits one axis at a time: three unknowns against every kept anchor, by226// the normal equations, eliminated in place.227func solve(donor, base []Anchor, keep []int) (Affine, bool) {228 xRow, ok := solveAxis(donor, base, keep, func(a Anchor) float64 { return a.X })229 if !ok {230 return Affine{}, false231 }232 yRow, ok := solveAxis(donor, base, keep, func(a Anchor) float64 { return a.Y })233 if !ok {234 return Affine{}, false235 }236 return Affine{237 AX: xRow[0], BX: xRow[1], CX: xRow[2],238 AY: yRow[0], BY: yRow[1], CY: yRow[2],239 }, true240}241242func solveAxis(donor, base []Anchor, keep []int, target func(Anchor) float64) ([3]float64, bool) {243 var sxx, sxy, syy, sx, sy, n float64244 var stx, sty, st float64245 for _, index := range keep {246 x, y := donor[index].X, donor[index].Y247 t := target(base[index])248 sxx += x * x249 sxy += x * y250 syy += y * y251 sx += x252 sy += y253 n++254 stx += t * x255 sty += t * y256 st += t257 }258 matrix := [3][3]float64{{sxx, sxy, sx}, {sxy, syy, sy}, {sx, sy, n}}259 vector := [3]float64{stx, sty, st}260 for i := range 3 {261 pivot := matrix[i][i]262 if math.Abs(pivot) < minPivot {263 return [3]float64{}, false264 }265 for j := i + 1; j < 3; j++ {266 factor := matrix[j][i] / pivot267 for k := range 3 {268 matrix[j][k] -= factor * matrix[i][k]269 }270 vector[j] -= factor * vector[i]271 }272 }273 var out [3]float64274 for i := 2; i >= 0; i-- {275 out[i] = vector[i]276 for k := i + 1; k < 3; k++ {277 out[i] -= matrix[i][k] * out[k]278 }279 out[i] /= matrix[i][i]280 }281 return out, true282}