Erdbeben hören: Raytracing durch einen ganzen Planeten

Eine S-Welle, die die Kern-Mantel-Grenze erreicht, ist am Ende. Der äußere Kern ist flüssig, Flüssigkeit hat keine Scherfestigkeit, und eine Scherwelle braucht etwas, das sie scheren kann — 2891,5 km tief hört die Phase also auf. Was das auf der Gegenseite eines großen Erdbebens hinterlässt, ist eine seismische Schattenzone von der Größe eines Kontinents, und sie gehört zu den Beobachtungen, die uns verrieten, dass die Mitte des Planeten geschmolzen ist — Jahrzehnte, bevor irgendjemand sie modellieren konnte.

Dieser Schatten steht in Tremor auf dem Schirm, sobald Sie ein großes Ereignis auswählen, und wie er dorthin kommt, ist der gesamte Entwurf im Kleinen. Kein Code zeichnet einen Schatten. Nirgends im Renderer existiert ein Schwellwert für die Winkeldistanz. Wir verfolgen Strahlen durch ein Geschwindigkeitsmodell, fragen jeden von ihnen, wo er zur aktuellen Replay-Zeit steckt, und schreiben die Antworten in einen vorab reservierten Puffer.

Strahlen, die am Kern gestorben sind, schreiben NaN. Der Line-Renderer bricht seinen Zug an einem NaN ab und beginnt einen neuen. Die Lücken sind die Schattenzonen.

Das ak135-Laufzeitmodell in 22 Stützstellen

Alles Weitere ruht auf einer Tabelle. ak135 ist ein gängiges eindimensionales Referenzmodell der Erde: zu jeder Tiefe eine Kompressions- und eine Schergeschwindigkeit, gemittelt über den ganzen Planeten — grob genug, um in einem Absatz beschrieben zu werden, und gut genug, dass Seismologen damit bis heute echte Erdbeben lokalisieren. Unsere Fassung sind 22 Zeilen [depthKm, vp, vs], stückweise linear im Radius statt in der Tiefe interpoliert und rund 0,1 km/s am veröffentlichten Modell dran. Die Quelle nennt das Sichtqualität und meint es auch so: ein Zeichenwerkzeug, keine ak135-Laufzeittabelle, gegen die Sie einen Phaseneinsatz stoppen würden.

Der Kniff in der Tabelle ist, wie Diskontinuitäten kodiert sind: Eine Tiefe, die zweimal auftaucht, ist eine Diskontinuität. Der Builder läuft die Stützstellenliste ab, zerlegt sie überall dort in stetige Abschnitte, wo sich eine Tiefe wiederholt, und notiert den Sprung. Sonst kennt nichts in der Engine die Grenzen beim Namen — Moho, die 410er, die 660er, die Kern-Mantel-Grenze und die Grenze zum inneren Kern sind eine einzige Datenstruktur, und eine weitere hinzuzufügen hieße, zwei Zeilen zu ergänzen.

Der flüssige äußere Kern ist in der Engine kein Sonderfall. Er ist eine Null in einer Tabelle.

Diese Null leistet enorm viel Arbeit. Erreicht ein S-Strahl die Kern-Mantel-Grenze, schlägt der Brechungsschritt die Schergeschwindigkeit auf der anderen Seite nach, findet nichts, in das er brechen könnte, und beendet den Strahl — der Tracer erfährt nie, dass er auf einen Kern getroffen ist. Unterhalb der Grenze zum inneren Kern kehrt die Schergeschwindigkeit zurück, weil der innere Kern fest ist; der Tracer würde dort fröhlich weiterlaufen, falls je eine Scherwelle ankäme. Keine tut es, und eine Assertion prüft, dass es nie eine getan hat.

360 Strahlen durch einen Planeten krümmen

Eine kugelsymmetrische Erde schenkt einem eine Dimension: Ein Strahl, der eine Quelle verlässt, verlässt nie die Ebene, die Quelle und Planetenmittelpunkt enthält. Der Tracer ist deshalb zweidimensional, schiebt jeden Strahl durch diese Großkreisebene und legt ihn als Radius, Winkeldistanz und aufsummierte Zeit ab — drei Floats je gespeicherter Probe, in einem flachen Float32Array, eine Probe je sechs Schritte.

Die Integration ist klassische Strahlentheorie: Mit dem Brechungsindex n = 1/v sagt die Strahlengleichung d/ds(n·d̂) = ∇n, dass ein Strahl sich zum höheren Index hin krümmt — also zum langsameren Material. In einem radialsymmetrischen Medium zeigt der Gradient entlang des Radius, die ganze Krümmung kollabiert damit auf ein Dutzend Zeilen skalarer Arithmetik, ausgeführt einmal je Acht-Kilometer-Schritt.

if (r > 0.5) {
  const dndr = -velocityGradientAt(r, phase) / (v * v);
  const k = (dndr * DS_KM) / r;
  const n = 1 / v;
  const wx = n * dx + k * x;
  const wz = n * dz + k * z;
  const wl = Math.hypot(wx, wz);
  if (wl > 0) {
    dx = wx / wl;
    dz = wz / wl;
  }
}

Geradliniges Schreiten genügt, bis ein Strahl auf eine Diskontinuität trifft — dort zeigt sich die Näherung sofort: Überquert er die Kern-Mantel-Grenze ein paar Kilometer zu spät, stimmt der Brechungswinkel nicht, und der Strahl landet irgendwo, wo er nichts zu suchen hat. Übergänge werden deshalb exakt gelöst. stepToRadius() bestimmt die Wurzel der quadratischen Gleichung, die den Strahl präzise auf die Grenze setzt, der Schritt wird dort abgeschnitten, und das Snelliussche Brechungsgesetz greift in der lokalen Tangential-/Normalenbasis: sin θ₂ = sin θ₁ · v₂/v₁. Übersteigt das eins, gibt es keinen transmittierten Strahl, und die Richtung spiegelt sich — Totalreflexion, und genau so tauchen die am Kern reflektierten Phasen auf, ohne dass jemand danach gefragt hätte. Andernfalls wird die gebrochene Richtung neu aufgebaut und 0,02 km freigerückt, damit der nächste Schritt denselben Übergang nicht noch einmal erkennt.

Dreihundertsechzig Strahlen fächern vom Hypozentrum aus, mit Abstrahlwinkeln zwischen 0,25° und 179,75° und einer auf 1–700 km begrenzten Herdtiefe. Ein voller Fächer kostet 10–20 ms pro Phase — auf dem Haupt-Thread, ohne Web Worker, was nur deshalb vertretbar ist, weil es einmal pro Auswahl passiert und nicht einmal pro Frame.

Die Phasen, die wir nie geschrieben haben

Die Seismologie hat ein reiches Phasenalphabet — P, S, PKP, PKIKP, PcP, ScS — und der verlockende Weg, es zu rendern, wäre, jede einzeln zu behandeln: für eine Kernphase diesen Bogen zeichnen, die direkte P jenseits jenes Winkels ausblenden. Von alldem hat Tremor nichts. Der Tracer kennt Geschwindigkeit, Gradienten, Grenzen und das Snelliussche Gesetz; jede Phase im Bild folgt aus diesen vier Dingen.

Der S-Fächer endet am Kern, weil die Schergeschwindigkeit dort null ist. Der P-Fächer tut etwas Interessanteres: Beim Übertritt in den äußeren Kern fällt die Kompressionsgeschwindigkeit von 13,66 auf 8,0 km/s, die Strahlen brechen also scharf nach unten, die direkten Oberflächenankünfte laufen jenseits von 104° aus, und die durch den Kern gelaufene Energie taucht viel weiter herum als PKP wieder auf. Beide Ränder des P-Schattens ergeben sich aus der Brechung, nicht aus einem if — und genau deshalb muss man sie prüfen.

Also bringt die Engine eine ausführbare Selbstprüfung mit, die in der Entwicklung läuft, im selben Geist wie die ISS-Assertion des Orrery: Eine Physik-Engine, die niemand mit bloßem Auge kontrollieren kann, braucht einen Test, der es kann.

• Die Laufzeit wächst entlang jedes Strahls monoton. (Ein Strahl, der in der Zeit rückwärts läuft, heißt: Der Integrator ist kaputt.)
• Die direkte P kommt bei 60° nahe 600 s an, innerhalb von ±60 s; bei 90° nahe 780 s innerhalb von ±80 s.
• Die direkte S kommt bei 60° nahe 1090 s an, innerhalb von ±110 s.
• Die letzte Oberflächenankunft der direkten P fällt in den Bereich 90°–118° — weit genug, um ein Modell visueller Güte zu tolerieren, eng genug, um eine kaputte Brechung zu erwischen.
• PKP-Ankünfte existieren jenseits von 130°.
• Kein S-Strahl wird je mehr als 5 km unterhalb der Kern-Mantel-Grenze gefunden.

Eine Isochrone pro Frame, null Allokation

Ein verfolgter Strahl ist ein Weg durch Raum und Zeit; eine Wellenfront ist ein Schnitt hindurch. sampleIsochron() sucht in den Proben jedes Strahls binär nach dem Intervall, das die aktuelle Replay-Zeit enthält, und interpoliert daraus eine Position — ein Radius/Distanz-Paar je Strahl, in ein vorab reserviertes Array. Strahlen, die diese Zeit nie erreicht haben oder am Kern gestorben sind, schreiben (NaN, NaN).

Die Fronten werden als LineSegments gezeichnet und nicht als Polyline, und diese Entscheidung ist tragend. Eine Polyline würde die letzte lebende Probe vor einer Lücke mit der ersten danach verbinden und eine Sehne quer über den Schatten schweißen. Segmente tun das nicht: Ein NaN beendet einen Zug, das nächste gültige Paar beginnt einen neuen, und der Schatten schneidet sich selbst aus der Geometrie heraus. Das Glühband besteht aus demselben Grund aus unabhängigen Quads — ein gemeinsamer Strip würde die Lücke wieder zuschweißen, die die Linien gerade gelassen haben.

Nichts davon wird pro Frame neu gebaut: 720 Probenplätze, 2880 Linienvertices und 1440 Bandquads werden einmal angelegt, und die Arbeit erledigt setDrawRange. Jede Front wird zweimal gezeichnet, gespiegelt auf ±delta, denn eine Wellenfront ist in einer kugelsymmetrischen Erde eine Rotationsfläche um die Quellachse, und der Schnitt legt zwei Halbebenen davon frei. Fronten halten 120 s simulierter Zeit und verblassen mit dem Alter — P von fast weiß über Glut, S von blassem Mint über Petrol —, damit sich ein dichter Replay als Abfolge liest und nicht als Schliere.

Eine aufgeschnittene Erde aus Geometrie, nicht aus Clipping

Um eine Wellenfront durch das Innere laufen zu sehen, muss man den Planeten aufschneiden, und der naheliegende Weg — Clip-Planes — ist für die Schalen das falsche Werkzeug: Fragment-Arbeit auf jedem Pixel jeder Schicht, und am Ende hohle Flächen. Stattdessen werden die Schalen unvollständig gebaut. Eine SphereGeometry mit phiStart = π und phiLength = 1.5π ist eine 270°-Schale mit fehlendem Quadranten, ausgerichtet auf die Standardkamera. Das Loch steckt im Mesh und kostet daher nichts.

Eine Schicht kann diesen Trick nicht nutzen. Küstenlinien und Gradnetz müssen sich drehen — das ist ja ihr Sinn —, sie sind deshalb das Einzige in der Szene, das tatsächlich auf der GPU geclippt wird, durch zwei Ebenen mit aktiviertem clipIntersection. Die Küstenliniendaten werden aus dem Modul dekodiert, das wir für das Orrery gepackt haben, und mit rund 2° je Unterteilung entlang von Großkreisen neu geslerpt, damit lange Sehnen nie unter die Oberfläche durchhängen und in der Kugel verschwinden.

Die Szene läuft dann auf zwei Bezugssystemen, und das ist das Stück, auf das wir still am stolzesten sind. Das Schnittsystem bewegt sich nie: Keil, Schichten, Grenzlinien und Polarachse sind statisch, die Clip-Planes im Weltraum sind damit Konstanten, und die Grafik der Schnittfläche backt sich in ein indiziertes Mesh und einen Draw Call. Das geografische System — Küstenlinien, Gradnetz, der ganze Katalog — wird so geslerpt, dass das ausgewählte Hypozentrum auf der +Y-Polarachse reitet, jener Kante, die sich beide freigelegten Flächen teilen. Jede Ebene durch die Quellachse ist ein gültiger Großkreisschnitt einer kugelsymmetrischen Erde, und hier sind zufällig zwei davon Wände, die man sehen kann. Der Querschnitt ist ehrlich, egal welches Ereignis Sie wählen — und für Besucher, die reduzierte Bewegung angefordert haben, sofort da, ganz ohne Animation.

Viertausend Erdbeben in einem Mesh

Der Katalog kommt aus zwei USGS-Übersichts-Feeds — der vergangene Tag in jeder Magnitude plus M4.5 und darüber aus dem vergangenen Monat — über einen Proxy, der jedes Ereignis auf ein Positionstupel abmagert und sein Caching selbst verantwortet. Eine Revalidierung auf Routenebene würde einen 503 genauso bereitwillig cachen wie Daten, was beim Orrery einmal einen Offline-Zustand für sechs Stunden festgefroren hat. Hier hängt die TTL von 90 Sekunden am älteren der beiden Feed-Zeitstempel, ein fehlgeschlagener Feed wird also weiter versucht, während der gesunde weiter ausliefert. Veraltet schlägt leer.

Auf dem Schirm ist der gesamte Katalog — bis zu 4000 Ereignisse — ein einziges InstancedMesh aus Oktaedern. Jede Instanz sitzt beim Radius (6371 − depth) / 1000 Szeneneinheiten: in ihrer echten Herdtiefe statt an die Oberfläche geheftet, sodass sich eine Subduktionszone als Platte liest, die sich in den Mantel neigt. Die Skalierung ist logarithmisch in der Magnitude, 0.012 + 0.02 · 1.9^(mag − 4), an beiden Enden begrenzt — nach oben gedeckelt, damit ein M8 eine Perle bleibt und kein Mond wird, nach unten, damit es die Mikrobeben noch gibt. Die Farbe ist das Alter: heißes Orange über Glut in den ersten 24 Stunden, dann Glut zu Schiefergrau über einen Schwanz von 30 Tagen. Die Geschichte der vergangenen Woche ist sichtbar, bevor Sie irgendwo klicken.

Erdbeben-Sonifikation aus einem selbst geschriebenen miniSEED-Decoder

Seismogramme sind öffentlich. Jeder kann rohe Bodenbewegung von den FDSN-dataselect-Diensten ziehen, und die IU-Stationen des Global Seismographic Network zeichnen seit Jahrzehnten durchgehend auf. Das Letzte, was Tremor mit einem Ereignis macht, ist deshalb, es Ihnen vorzuspielen — echte Erdbeben-Sonifikation, kein Synthesizer, der eine imitiert.

Der Browser schafft das aus zwei voneinander unabhängigen Gründen nicht allein. Die Daten kommen als Steim-komprimiertes miniSEED, ein differenzkodiertes Binärformat, von dem kein Browser je gehört hat, und der Dienst schickt keine CORS-Header, eine Seite könnte die Bytes also nicht einmal holen, wenn sie sie lesen könnte. Beide Hälften wohnen auf unserem Server: Ein Proxy holt, ein selbst geschriebener Decoder packt aus, und beim Client kommt ein normalisierter Int16-Puffer an.

Die Stationswahl ist ein kleines Ranking-Problem mit einer physikalischen Nebenbedingung darin. Tremor kennt 31 IU-Breitbandstationen und will vierzig Minuten Aufzeichnung der Vertikalkomponente, aber am nächsten ist nicht am besten: näher als 2° übersteuert das Instrument bei einem großen Ereignis, jenseits von 90° beginnt der Kern die P abzuschatten. Also führt der teleseismische Bereich das Ranking an, die Distanz entscheidet Gleichstände, drei Kandidaten werden der Reihe nach probiert, und eine Aufzeichnung jünger als 45 Minuten wird rundheraus abgelehnt — Archive brauchen Zeit, und der Client nimmt diese Ablehnung vorweg, statt Sie sie entdecken zu lassen.

Der miniSEED-Decoder ist ungefähr so hardwarenah, wie Web-Code werden kann: ein DataView, Big-Endian-Zugriffe und ein Blockette-1000-Header, der sagt, welche Kodierung gerade vorliegt. Steim packt Samples in 64-Byte-Frames aus sechzehn 32-Bit-Wörtern, deren erstes eine Nibble-Map ist — sechzehn Zwei-Bit-Codes, die beschreiben, wie die anderen fünfzehn gepackt sind. Steim-1 bietet eine, zwei oder vier Differenzen pro Wort; Steim-2 ergänzt Untercodes, die fünf 6-Bit-, sechs 5-Bit- oder sieben 4-Bit-Differenzen in ein Wort quetschen, plus einen 30-Bit-Platz für Ausreißer. Die Wörter eins und zwei von Frame null sind gar keine Differenzen, sondern die Integrationskonstanten X0 und Xn, was den Strom selbstprüfend macht. Beide Kodierungen sind implementiert; wo das Timing zweier Records um mehr als 1,5 Samples auseinanderliegt, zählen wir die Lücke und hängen trotzdem aneinander, denn für die Sonifikation schlägt Kontinuität die Exaktheit.

Womit wir bei dem Teil wären, in dem kein Code steckt. Seismische Energie lebt um 0,05–10 Hz; das Hören beginnt nahe 20 Hz. Die Lehrbuchlösung ist ein Pitch-Shifter, und wir haben genau keinen geschrieben, weil die Web Audio API einen die Sample-Rate eines Puffers deklarieren lässt. Vierzig Samples pro Sekunde, deklariert als 40 × 150 = 6000 Hz, spielen 150-mal schneller und 150-mal höher — innerhalb des zulässigen Puffer-Rate-Fensters von 3000 bis 768.000 Hz — und legen das seismische Band auf 7,5–1500 Hz. Vierzig Minuten werden sechzehn Sekunden. Das ist keine Transformation der Daten; es sind die Daten, nur schnell gelesen.

Der Dynamikumfang ist das härtere Problem. Ein teleseismischer P-Einsatz kann das Hundertfache der Koda dahinter sein, und linear abgespielt hört der Besucher einen Tick und dann Stille. Der Client kompandiert deshalb mit einem Exponenten nahe 0,6 — im Browser, nicht auf dem Server, denn die Aufgabe des Servers sind ehrliche Zählwerte, und die Formung gehört dorthin, wo der Besucher sie wirken hört. Vor dem Transport wird die Aufzeichnung mittelwertbefreit (SNZO etwa liegt rund 1500 Counts neben null) und spitzenwertnormiert. Die Sparkline zeichnet aus demselben Grund den Spitzenwert je Bin statt des RMS: Die Ankunft ist eine Spitze, und Mitteln ist die eine Operation, die sie verstecken würde.

Der Player selbst ist bewusst klein: ein GainNode, 30-ms-Rampen mit einem Schritt Vorlauf geplant, ein Audio-Kontext, der im Klick-Handler entsteht, damit die Autoplay-Policy nichts zu beanstanden hat, und eine Wiedergabe, die stoppt, wenn der Tab verschwindet. Die Stummschaltung gehört nicht Tremor — sie liest den seitenweiten Synth-Mute, den sich alle Spiele im Lab teilen. Schalten Sie ein Spiel stumm, wird auch der Planet still.

In der Ecke des Bildschirms steht ein Chip mit Laufzeiten in Sichtqualität, aus demselben Grund, aus dem beim Orrery Propagation in Sichtqualität steht. Ein eindimensionales Modell ohne Stationskorrekturen und ohne dreidimensionale Struktur setzt eine Wellenfront dorthin, wo das Auge sie erwartet, nicht dorthin, wo ein Instrument sie aufzeichnen würde. Treu fürs Auge, nicht auf den Meter; nur zum Staunen, für nichts sonst. Das offen zu sagen kostet nichts und kauft dem Rest des Bildes seine Glaubwürdigkeit.

Die ehrlichen Lücken: Das Raytracing läuft auf dem Haupt-Thread, und die Seite hat noch keine automatisierte End-to-End-Abdeckung — nur die Assertions zur Entwicklungszeit, die die Physik prüfen und nicht die Knöpfe. Beides steht auf der Liste.