Implement sonicSimilarity extension with real audio feature analysis

Adds internal/audiofeatures: decodes each scanned song to PCM via
ffmpeg and computes actual measured characteristics — tempo (via
autocorrelation of the frame energy envelope, left unset rather than
guessed when a clip is too short to contain a beat), an approximate
loudness figure, spectral centroid/rolloff and zero-crossing rate (via
a hand-written FFT), and a 12-bin chroma vector used both for a
Krumhansl-Kessler key estimate and as similarity-index input. No
external DSP/ML dependency.

The scanner now runs this analysis after indexing each song (as part
of the already-background scan, not on any request path) and stores
results in the existing audio_features table via the new
database.UpsertAudioFeatures/GetAudioFeatures/ListAllAudioFeatures.

internal/similarity is a small in-process nearest-neighbor index:
z-score normalizes each feature dimension across the analyzed library
and answers Nearest (ranked real distances) and Path (a feature-space
straight-line walk from one track to another, hopping to the closest
not-yet-used real song at each step — never a fabricated intermediate
track).

getSonicSimilarTracks and findSonicPath expose these through the
standard sonicMatch response, converting distance to a normalized
(0,1] similarity score. A song that hasn't been analyzed yet is simply
absent from results rather than faked.

Co-Authored-By: Claude Sonnet 5 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_014r98AtcBbevpscWLsfwVvv

internal/audiofeatures/fft.go

@@ -0,0 +1,65 @@
+package audiofeatures
+
+import "math"
+
+// fft computes the in-place iterative radix-2 Cooley-Tukey FFT of a,
+// whose length must be a power of two. It's the only DSP primitive this
+// package needs (spectral centroid/rolloff and the chroma vector all
+// derive from the resulting magnitude spectrum), so a small
+// hand-written implementation avoids pulling in an external numerics
+// dependency for one function.
+func fft(a []complex128) {
+	n := len(a)
+	if n <= 1 {
+		return
+	}
+
+	// Bit-reversal permutation.
+	for i, j := 1, 0; i < n; i++ {
+		bit := n >> 1
+		for ; j&bit != 0; bit >>= 1 {
+			j ^= bit
+		}
+		j ^= bit
+		if i < j {
+			a[i], a[j] = a[j], a[i]
+		}
+	}
+
+	for length := 2; length <= n; length <<= 1 {
+		ang := -2 * math.Pi / float64(length)
+		wlen := complex(math.Cos(ang), math.Sin(ang))
+		for i := 0; i < n; i += length {
+			w := complex(1, 0)
+			half := length / 2
+			for k := 0; k < half; k++ {
+				u := a[i+k]
+				v := a[i+k+half] * w
+				a[i+k] = u + v
+				a[i+k+half] = u - v
+				w *= wlen
+			}
+		}
+	}
+}
+
+// magnitudeSpectrum returns |FFT(frame)| for bins [0, n/2], given a
+// real-valued frame whose length is a power of two.
+func magnitudeSpectrum(frame []float64) []float64 {
+	n := len(frame)
+	buf := make([]complex128, n)
+	for i, v := range frame {
+		buf[i] = complex(v, 0)
+	}
+	fft(buf)
+	mags := make([]float64, n/2+1)
+	for i := range mags {
+		mags[i] = cmplxAbs(buf[i])
+	}
+	return mags
+}
+
+func cmplxAbs(c complex128) float64 {
+	re, im := real(c), imag(c)
+	return math.Sqrt(re*re + im*im)
+}