Tomografie van thermosferische dichtheid op basis van Starlink Ephemeris: initieel rapport

Samenvatting

Dit rapport presenteert de eerste succesvolle demonstratie van een tomografische analyse van de thermosferische dichtheid met behulp van orbitale gegevens van een massief Starlink-satellietconstellatie. De thermosferische dichtheid — de neutrale atmosferische dichtheid op hoogten tussen 100 en 1000 km — is essentieel voor de vooruitgang van de wetenschap over de bovenste atmosfeer en ter ondersteuning van ruimtevaarttechnische operaties.

Waar traditionele observatiemethoden vertrouwen op Two-Line Element (TLE) gegevens, presenteert deze studie een nieuwe tomografische analyse met behulp van Starlink Ephemerides: gedetailleerde orbitale informatie die publiekelijk is vrijgegeven door SpaceX. In tegenstelling tot TLE's bieden deze ephemeriden positie- en snelheidsvectoren met een hoge resolutie, waardoor het energieverlies van satellieten door atmosferische sleep (drag) nauwkeurig kan worden geschat.

Er is een methode ontwikkeld om het energieverlies te kwantificeren en een coëfficiënt te bepalen die dit relateert aan de thermosferische dichtheid. De energiedissipatie werd berekend voor ongeveer 1200 satellieten op een hoogte van 482 km en een inclinatie van 53°, wat diende als basis voor de tomografische analyse. Er werd gebruikgemaakt van een sferische harmonische expansie om de dichtheidsverdeling per lengte- en breedtegraad te modelleren, waarbij in dit rapport uitsluitend rekening is gehouden met diurnale (dagelijkse) variaties. Er werd aangenomen dat de dichtheid-hoogteverdeling exponentieel afneemt met een enkele schaalhoogte van 60 km.

De analyse werd uitgevoerd voor de periode van 1 tot 7 september 2025 en legt kenmerken van de thermosferische dichtheid vast die worden voorspeld door het NRLMSIS 2.1-model. De resultaten werden gevalideerd tegen onafhankelijke dichtheidsmetingen van de SWARM-satellieten. De tomografische schattingen vertoonden een hoge consistentie met de SWARM-waarnemingen, bijvoorbeeld in dichtheidsvariaties langs satellietbanen. Kwantitatieve vergelijkingen over 19 gevallen toonden aan dat de geschatte dichtheidswaarden varieerden van 0,6 tot 1,2 keer de SWARM-metingen, met een gemiddelde ratio van 0,95. Deze bevindingen suggereren dat deze methode het potentieel heeft om snelle, hoge-resolutie gegevens van de thermosferische dichtheid te leveren.

1. Inleiding

Thermosferische dichtheid verwijst naar de neutrale atmosferische dichtheid op hoogten van 100 tot 1000 km boven het aardoppervlak. De dichtheid van de bovenste atmosfeer is extreem laag; op een hoogte van 300 km is deze bijvoorbeeld $2 \times 10^{-11}$ keer de dichtheid aan het aardoppervlak. De bovenste atmosfeer bestaat voor meer dan 99% uit een elektrisch neutrale atmosfeer, de "thermosfeer" genoemd, en voor minder dan 1% uit een geïoniseerde atmosfeer, de "ionosfeer". Vanwege het effect van geïoniseerd gas op radiopropagatie is de ionosfeer relatief eenvoudig te observeren, terwijl waarnemingen van de thermosfeer moeilijker zijn. Het ontwikkelen van observatiemethoden voor thermosferische dichtheid is belangrijk voor de verdere vooruitgang van de wetenschap over de bovenste atmosfeer. Bovendien neemt het belang van thermosferische dichtheid vanuit technisch perspectief toe door de recente groei in ruimtevaartontwikkeling.

Een belangrijke meettechniek voor thermosferische dichtheid is het gebruik van atmosferische sleep op satellietbanen (Emmert et al., 2015). Het Two-Line Element (TLE) is het meest gebruikte formaat voor orbitale informatie, terwijl de Vector Covariance Message (VCM) een nieuwer formaat is. Emmert et al. (2021) presenteerden gedetailleerde analytische methoden voor dichtheidsschatting op basis van TLE en VCM. Zij merken op dat satellietbaangegevens nu primair gebaseerd zijn op VCM en dat TLE een product is van VCM, wat suggereert dat het gebruik van VCM analytisch gezien te prefereren is boven TLE. Vanuit het perspectief van gegevensbeschikbaarheid zijn TLE-gegevens publiekelijk beschikbaar voor de meeste satellieten, terwijl VCM-gegevens alleen beschikbaar zijn onder specifieke contracten met databronnen (Emmert et al., 2021).

Starlink is een communicatienetwerk bestaande uit 9000 of meer satellieten in een lage aardbaan (LEO), uitgerold door SpaceX. Onderzoek naar de bovenste atmosfeer met behulp van orbitale informatie van Starlink-satellieten is reeds uitgevoerd. Er is met succes de relatieve thermosferische dichtheid geschat met behulp van TLE, waarbij de tijd-hoogtevariatie is aangetoond (Yamamoto en Sori, 2026). Opmerkelijk is dat SpaceX gedetailleerde orbitale informatie van Starlink-satellieten publiekelijk heeft gemaakt, equivalent aan VCM-gegevens, bekend als Starlink Ephemerides. Een studie naar het gebruik van Starlink Ephemerides vond een relatie met de thermosferische dichtheid (Ou et al., 2025). Dit was baanbrekend onderzoek, omdat het een toename in atmosferische dichtheid aantoonde door een afname van de mechanische energie van de satelliet tijdens geomagnetische stormen. Echter, een volledige schatting van de thermosferische dichtheid werd toen nog niet bereikt.

Dit document is een initieel rapport over de tomografische analyse van thermosferische dichtheid met behulp van Starlink Ephemerides. We hebben de thermosferische dichtheid met een goede nauwkeurigheid geschat met uitsluitend publiek beschikbare gegevens en de diurnale en globale distributie aangetoond. De huidige resultaten zijn nog preliminair en hebben een lage resolutie, maar we achten het van groot belang om onze ideeën en data-analyse te delen met de onderzoeksgemeenschap.

2. Starlink Ephemeris en mechanische energie

We gebruiken Starlink Ephemerides om de thermosferische dichtheid te schatten. SpaceX publiceert deze gegevens om botsingen tussen Starlink-satellieten en andere objecten te voorkomen. Tot juli 2025 werden de gegevens verstrekt via de Space-Track website; sindsdien worden ze rechtstreeks door SpaceX gepubliceerd.

Een gepubliceerde dataset biedt voorspellingen van de banen van Starlink-satellieten voor 72 uur vanaf het tijdstip van de data-epoch. Omdat de dataset drie keer per dag wordt vernieuwd, is het interval tussen updates ongeveer acht uur. Een enkel gegevensbestand voor een specifieke satelliet registreert de positie en snelheid van de satelliet per minuut in het Mean Equator Mean Equinox (MEME) coördinatensysteem, geverifieerd tegen J2000.0. De gegevens bevatten ook een covariantiematrix voor zes parameters (drie positie- en drie snelheidcomponenten) per tijdstip. Deze informatie is identiek aan de VCM-gegevens beschreven door Emmert et al. (2021).

Starlink-satellieten zijn gegroepeerd in "shells", bestaande uit satellieten op dezelfde hoogte en inclinatie. Satellieten in een shell voeren manoeuvres uit om hun hoogte te behouden; de frequentie van deze manoeuvres neemt toe naarmate de hoogte afneemt. Gegevens uit september 2025 laten zien dat de shells hoogten hebben variërend van 540-580 km, 450-490 km en 360-370 km.

Methodiek van orbitale analyse

Het gebruik van ephemerisgegevens verschilt aanzienlijk van dat van TLE's. Omdat de ephemeris een eenvoudige tijdreeks is van positie-/snelheidsvectoren, moeten we de "Special Perturbations" (SP) propagator gebruiken — een numerieke integrator van de bewegingsvergelijkingen — om de relatie tussen de datareeksen te bepalen. We gebruiken de open-source software General Mission Analysis Tool (GMAT) van NASA voor deze analyses. Onze berekeningen zijn gebaseerd op:

  • Het zwaartekrachtveld van de aarde met perturbaties tot de 20ste orde of hoger.
  • Puntmassa's van de zon en de maan.
  • Atmosferische sleep gebaseerd op het Jacchia-Roberts dichtheidsmodel.
  • Exclusie van kleine fluctuaties door vaste-aarde of oceaangetijden en stralingsdruk.

Relatie tussen ephemeris en dichtheid

Volgend op het werk van Picone et al. (2005) en Emmert et al. (2021), berekenen we de specifieke mechanische energie $E$ (J/kg) van de satelliet. $E$ is de som van de specifieke kinetische energie $K$ en de specifieke potentiële energie $U$:

$$E = K + U = \frac{1}{2}v^{2} - \frac{\mu }{r}\left[ {1 - \sum\limits{{\text{n} = 2}}^{\infty } {J{\text{n}} \left( {\frac{{R{e} }}{r}} \right)^{\text{n}} P{\text{n}} (\sin \phi )} } \right]$$

waarbij:

  • $v$: snelheidsamplitude
  • $r$: afstand tot het middelpunt van de aarde
  • $\phi$: geocentrische breedtegraad
  • $\mu$: gravitatieparameter van de aarde ($GM$)
  • $R_e$: equatoriale straal van de aarde
  • $J_n$: gravitatie-harmonische coëfficiënten
  • $P_n$: Legendre-polynomen

Atmosferische sleep of satellietmanoeuvres resulteren in respectievelijk een daling of stijging van de baan, wat gemeten kan worden als een afname of toename van $E$. De atmosferische sleepkracht $\vec{FD}$ werkt in de tegenovergestelde richting van de relatieve snelheid $\vec{vr}$ ten opzichte van de atmosfeer:

$$\vec{FD} = - \frac{1}{2}\text{C}{\text d} A\rho v{r} \vec{vr}$$

waarbij $\rho$ de dichtheid is, $A$ de dwarsdoorsnede van de satelliet en $Cd$ de sleepcoëfficiënt. Uitgaande van een atmosfeer die meedraait met de aarde, is de relatieve snelheid: $$\vec{vr} = \vec{v} - \vec{\omega_e} \times \vec{r}$$

Het energieverlies per tijdseenheid (vermogensverlies) wordt dan: $$\frac{dE}{{dt}} = - \frac{1}{2}B\rho v{r} (\vec{vr} \cdot \vec{v} )$$

waarbij we definiëren $B = \text{C}{\text{d}} A/ms$ (het omgekeerde van de ballistische coëfficiënt). Wanneer positie en snelheid bekend zijn, kan het vermogensverlies worden bepaald door de dichtheid $\rho$ en de coëfficiënt $B$.

Analyse van energiefluctuaties en manoeuvres

Door de mechanische energie $E$ te analyseren, kunnen we onderscheid maken tussen verschillende effecten:

  1. Zwaartekrachtveld ($E_g$): Kleine fluctuaties ontstaan door het verschil tussen het geoïde en het ellipsoïde van de aarde.
  2. Manoeuvres ($E_{mp}$): Plotselinge stijgingen in energie duiden op actieve baanbehoudmanoeuvres.
  3. Atmosferische sleep ($E_{ap}$): Een constante afname in energie, vaak gecombineerd met kleine oscillaties langs de baan.

Door manoeuvres te identificeren en te verwijderen uit de data, kunnen we het werkelijke geïntegreerde energieverlies door atmosferische sleep ($E{ao}$) isoleren tussen twee tijdstippen $t1$ en $t_2$. In een analyse van ongeveer 1.200 satellieten in de shell op 482 km hoogte werden in 80% van de data manoeuvres gedetecteerd, wat aantoont dat correctie noodzakelijk is voor een hoge databeschikbaarheid.

De coëfficiënt $B$ kan worden geschat door de helling van het energieverlies in GMAT-simulaties te variëren tot deze nul is.

3. Tomografie-analyse

Het doel is een tomografische analyse van de thermosferische dichtheid. Als $y$ het gemeten energieverlies is, schrijven we: $y = E{a} (\vec{x} ) + e{r}$, waarbij $E{a} (\vec{x})$ de modelestimatie is en $er$ de meetfout.

Het model voor het energieverlies is: $$E{a} (\vec{x}) = - \frac{1}{2}B\int{{t{1} }}^{{t{2} }} {\rho (\vec{x})v{r} (\vec{vr} \cdot \vec{v} )dt}$$

Met een groot aantal $N$ datapunten gebruiken we de kleinste-kwadratenmethode om de parameters $\vec{x}$ te vinden die de som van de kwadraten van de residuen minimaliseren.

Modellering van de dichtheid $\rho$

We nemen aan dat $\rho$ verticaal verdeeld is in hydrostatisch evenwicht rond een referentiehoogte $hr$ met een constante schaalhoogte $H$: $$\rho (h,\phi ,\lambda ) = \exp \left( { - \frac{{h - h{r} }}{H}} \right)f\left( {\phi ,\lambda } \right)$$

De functie $f(\phi, \lambda)$ wordt geparametriseerd via sferische harmonische analyse: $$f(\theta ,\lambda ) = \sum\limits{\text{l} = 0}^{\text{L}{\text{max}} } {\sum\limits{\text{m} = 0}^{\text{l}} {\overline{{P{\text{lm}} }} (\cos \theta )(C{\text{lm}} \cos \text{m}\lambda + S{\text{lm}} \sin \text{m}\lambda )} }$$

Implementatieproces

De analyse is onderverdeeld in drie fasen:

  • Fase 1 (Datageneratie): Schatten van energieverlies door atmosferische sleep en de coëfficiënt $B$ voor elke satelliet.
  • Fase 2 (Voorbereiding): Definiëren van analysecondities, vaststellen van de lokale-tijd longitude en berekenen van sferische harmonische basisfuncties. Initiële waarden worden verkregen uit het NRLMSIS 2.1-model.
  • Fase 3 (Optimalisatie): Numerieke minimalisatie van de residuen tussen model en meting met behulp van scipy.optimize.leastsquares. Hierbij wordt $\log{10}\rho$ gebruikt in plaats van $\rho$ om variaties te matigen en te garanderen dat $\rho > 0$.

4. Resultaten en samenvatting

De analyse werd uitgevoerd met gegevens van 1 tot 7 september 2025, gebruikmakend van een shell van ongeveer 1.200 satellieten op 482 km hoogte met een inclinatie van 53°.

Validatie en Vergelijking

De resultaten tonen aan dat de maximale atmosferische dichtheid zich bevindt tussen 200° en 220° longitude en tussen 0° en 20° breedtegraad. De analyse reproduceert ruwweg de diurnale distributie van het NRLMSIS 2.1-model.

Vergelijkingen met publieke gegevens van de SWARM-satellieten (hoogten 440–470 km en 490–520 km) laten het volgende zien:

  • De tomografie volgt nauwgezet de dichtheidswaarden van de SWARM-satellieten.
  • Er is een goede overeenstemming nabij de maxima van de variatie.
  • Nabij de minima ligt de tomografie dichter bij de SWARM-metingen dan het NRLMSIS-model.
  • Ondanks dat alleen satellieten binnen $\pm 53^\circ$ breedtegraad zijn gebruikt, is er een redelijke globale dichtheidsverdeling gegenereerd.

Over 19 geanalyseerde gevallen varieerden de dichtheidswaarden van de tomografie van 0,5 tot 1,5 keer die van de SWARM-satellieten, met een gemiddelde ratio van 0,95.

Conclusies en beperkingen

Dit is naar belief de eerste publicatie die tomografische analyse van thermosferische dichtheid rapporteert op basis van echte data uit Starlink Ephemerides. De methode heeft het potentieel om de spatiotemporele structuur van de thermosfeer met een ongekende resolutie in kaart te brengen.

Er zijn echter beperkingen:

  1. Aanname over manoeuvres: De aanname dat geplande manoeuvres exact worden uitgevoerd ($E{mo} \approx E{mp}$) vereist verder onderzoek.
  2. Filtering: De ephemeriden worden verwerkt door het trackingsysteem van SpaceX, wat mogelijk een filtereffect heeft dat echte observaties mengt met modelvoorspellingen.
  3. Data-hiaten: Door de inclinatie van 53° is er een gebrek aan data op hoge breedtegraden, waardoor analyses met $\text{L}_{\text{max}} > 1$ momenteel falen.

Toekomstig onderzoek zal zich richten op het gebruik van satellieten met hogere inclinaties en de synthese van gegevens uit meerdere shells.

Beschikbaarheid van gegevens

  • Starlink Ephemerides: Verstrekt door SpaceX via starlink.com/public-files/ephemerides/README.md.
  • NRLMSIS 2.1: Beschikbaar via NASA Community Coordinated Modeling Center.
  • SWARM-gegevens: Beschikbaar via thermosphere.tudelft.nl.
  • GMAT software: Beschikbaar via NASA.
  • Software-stack: Ontwikkeld in Python met NumPy, SciPy, Matplotlib, pyshtools, pymap3d en pymsis.

Afkortingen

  • ESA: European Space Agency
  • GMAT: General Mission Analysis Tool
  • GPS: Global Positioning System
  • MEME: Mean Equator Mean Equinox coördinatensysteem
  • MSIS: Mass Spectrometer Incoherent Scatter radar
  • NASA: National Aeronautics and Space Administration
  • NRL: Naval Research Laboratory
  • TLE: Two-Line Element
  • VCM: Vector Covariance Matrix

Referenties

(De originele lijst met wetenschappelijke referenties blijft behouden voor bronverwijzing, inclusief werken van Bate et al., Emmert et al., Liu et al., Ou et al., Picone et al., Shoemaker et al., Vallado, Van den IJssel et al., Wieczorek et al., en Yamamoto & Sori).