AIDA: Advanced Ionospheric Data Assimilation

AIDA Real-Time Data Model Output»

Report Bug · Request Feature

Table of Contents

1. About The Project
2. Getting Started
3. Usage
4. License
5. Contact

About The Project

AIDA is a real-time ionosphere/plasmasphere data assimilation model. AIDA uses measurements from ground- and satellite-based Global Navigation Satellite (GNSS) receivers and ionosondes to produce an improved global ionospheric representation. This package contains the necessary software to read the output files produced by the AIDA system, and produce outputs of electron density, Total Electron Content (TEC), MUF3000, and various ionospheric profile parameters (NmF2, foF2, hmF2, etc.)

The AIDA interpreter is a standalone package to read the AIDA output files. It requires output from the AIDA model from https://spaceweather.bham.ac.uk/output/.

To access AIDA model output and forecast products, it is necessary to create a free account at https://spaceweather.bham.ac.uk/accounts/register/.

(back to top)

Getting Started

Installation (pip + git)

Currently the AIDA interpreter must be installed using git. In the future the model interpreter will be made available on PyPI.

python -m pip install git+https://github.com/breid-phys/aida-ionosphere.git

Installation (git clone)

Currently the AIDA interpreter must be installed using git. In the future the model interpreter will be made available on PyPI.

1. Clone the repo

   git clone https://github.com/breid-phys/aida-ionosphere.git

2. Install aida package

   python -m pip install -e /path/to/aida

Configuring the AIDA API

AIDA requires files from the AIDA data assimilation model to produce model output. These files can be automatically downloaded using an API, which requires some configuration. By default, aida will look for a file called api_config.ini in an OS-dependent location. This file must be edited to include two pieces of information.

Linux/Mac

    ~/.config/aida/api_config.ini

Windows

    %USERPROFILE%\AppData\Local\aida\api_config.ini

This file can be created using the function aida.api.configure_api().

API Token: To be able to automatically download output, the api_config.ini file will need to be edited to include your unique API token, which can be found on the SERENE Website. This will require creating an account.

Output Cache: AIDA will need a location to cache output files in order to produce output, which is configurable in api_config.ini. This path is broken into two parts, the folder and the subfolder. This allows AIDA to sort the output into (and create) subdirectories. The folder gives the path to the top-level directory, which must already exist. The subfolder contains a series of tags which allow AIDA to create subfolders based on the date of the output file.

WARNING: There is no automatic cleanup of cached output files

A blank example file can be found in the project directory, which can be copied to the home directory. The location of this example file can be found by calling aida.api.find_api_config().

(back to top)

Usage

The following python packages will be used in the demonstration.

import aida
import numpy as np
import matplotlib.pyplot as plt
from matplotlib import colormaps

AIDA State Object

All AIDA calculations are performed through the AIDAState class.

# Create an empty AIDAState
Model = aida.AIDAState()

An AIDA State must be populated with data. The AIDA State can download and read files automatically using the AIDAState.fromAPI() method. If the config file is located somewhere other than the default location, it can be passed to aida.AIDAState() as an argument APIconfig.

Model.fromAPI(time=np.datetime64("2025-01-01T13:55:00"),
              model='AIDA',
              latency='rapid',
              forecast=90)

To download the latest output for a given model, pass the argument time='latest'.

Model.fromAPI(time='latest', model='AIDA', latency='ultra')

An AIDA State can also be told to read an AIDA output file using the AIDAState.readFile() method.

# load test output file
Model.readFile("./tests/data/output_3_231201_042500.h5")

The AIDAState class has several attributes:

• AIDAState.Version is the version and operating mode of the AIDA model which produced the file
• AIDAState.Time is the time for which this file provides a description of the ionospheric state (in UNIX Epoch format)
• AIDAState.Metadata contains information about the model used to generate the file
  • AIDAState.Metadata.ForecastStart gives the time of the previous AIDA output used to forecast this output (in UNIX Epoch format)
  • AIDAState.Metadata.NeQuickFlux gives the F10.7 cm flux passed to NEQuick for generating the AIDA background model
  • AIDAState.Metadata.NeQuickVersion is the version of NeQuick used for generating the AIDA background model
  • AIDAState.Metadata.InstrumentsGNSS gives a list of GNSS receivers used for generating the output
  • AIDAState.Metadata.InstrumentsIonosonde gives a list of ionosondes used for generating the output

The background state of an AIDA model output can be accessed using the AIDAState.background() method, which returns a copy of the AIDAState object with all parameters set to match the background state.

Example 1: Electron Density (3D Grid)

For most cases, the AIDAState.calc() method is used to calculate values of interest. AIDAState.calc() accepts glat, glon, and alt arguments, and the shape of the output is determined by the grid argument. If grid = '3D', the glat, glon, and alt inputs will be broadcast together to produce a regular 3D grid. When grid=3D or grid=2D, AIDA can be more efficient when performing the internal spherical harmonic transformations, which is much more efficient for large, regular grids. When calculating arbitrary points, grid=1D must be used.

# create several arrays of coordinates
glat = np.linspace(-90.0, 90, 100)
glon = [0.0, 180.0]
alt = np.linspace(0.0, 2e3, 150)

# calculate AIDA output
Output = Model.calc(lat=glat, lon=glon, alt=alt, grid="3D", collapse_particles=True)

By default, Output is in the form of an xarray.Dataset object, for more information see the xarray documentation. To get the output in the form a python dict, pass as_dict=True to AIDAState.calc().

The output of AIDAState.calc() will have a field Ne containing the electron density, along with fields for all of the AIDA ionospheric profile parameters. In a 3D grid, the Ne field will have a size (glon, glat, alt), and all other parameters will have size (glon, glat).

By default, the output fields will also have a dimension labelled particle, which corresponds to the ensemble size of the AIDA state. When working with model output, the particle dimension will always be of size 1, and so this extra dimension can be suppressed with the argument collapse_particles=True.

# create figure
fig, axs = plt.subplots(len(glon), 1, squeeze=False)
fig.set_size_inches(12, 9)

for i, lon in enumerate(glon):
    # plot latitudinal section for each longitude
    ax = axs[i, 0]
    pcm = ax.pcolor(Output.glat, Output.alt, np.log10(Output.Ne.sel(glon=lon).T))
    ax.set_title(f"AIDA Ne: lon = {lon}")
    fig.colorbar(pcm, ax=ax, label=f"log10 {Output['Ne'].attrs['units']}")
    pcm.set_clim(8, 12.5)
    pcm.set_cmap(colormaps["magma"])

Example 2: Electron Density (2D Grid)

For irregular grids, the argument grid=2D can be used. This requires glat and glon to be the same size.

# create irregular grid
glat = [[45.0, 65.0, 77.0], [30.0, 40.0, 50.0]]
glon = [[0.0, 345.0, -160.0], [-60.0, 61.1, 62.2]]
alt = np.linspace(90.0, 500.0, 150)

Output = Model.calc(
    lat=glat, lon=glon, alt=alt, grid="2D", collapse_particles=True, as_dict=True
)

In this example, the output is in dict format.

# create figure
fig, axs = plt.subplots(1, 2, squeeze=False)
fig.set_size_inches(6, 6)

for i, lon in enumerate(glon):
    # plot profiles and mark layers
    ax = axs[0, i]
    ax.plot(Output["Ne"][i, :].T, Output["alt"])
    ax.scatter(Output["NmF2"][i, :], Output["hmF2"][i, :])
    ax.scatter(Output["NmF1"][i, :], Output["hmF1"][i, :])
    ax.scatter(Output["NmE"][i, :], Output["hmE"][i, :])

Example 3: Maps

AIDA uses a parameterized ionospheric profile, and many of these profile parameters are of interest, in addition to the electron density. AIDAState.calc() will generate these values along with the electron density. If the alt argument is omitted, electron density Ne will not be calculated.

AIDA can also calculate derived quantities, such as the Total Electron Content (TEC) and MUF3000. These are more computationally expensive, and are disabled by default. They can be included by passing the arguments TEC=True and MUF3000=True to AIDAState.calc().

glat = np.linspace(-90.0, 90)
glon = np.linspace(-180.0, 180.0, 70)

Output = Model.calc(
    lat=glat, lon=glon, grid="3D", TEC=True, MUF3000=True, collapse_particles=True
)

In this example, the AIDA model will be compared to the values given in the background model.

BkgModel = Model.background()
BkgOutput = BkgModel.calc(
    lat=glat, lon=glon, grid="3D", TEC=True, MUF3000=True, collapse_particles=True
)

fig, axs = plt.subplots(3, 5, squeeze=False)
fig.set_size_inches(15, 10)

for i, d in enumerate(["NmF2", "hmF2", "TEC", "MUF3000", "foF2"]):

    ax = axs[0, i]
    pcm = ax.pcolor(Output.glon, Output.glat, Output[d].T)
    ax.set_title(f"AIDA {d}")
    fig.colorbar(pcm, ax=ax, label=Output[d].attrs["units"])
    pcm.set_cmap(colormaps["magma"])
    cl = pcm.get_clim()

    ax = axs[1, i]
    pcm = ax.pcolor(BkgOutput.glon, BkgOutput.glat, BkgOutput[d].T)
    ax.set_title(f"Background {d}")
    fig.colorbar(pcm, ax=ax, label=Output[d].attrs["units"])
    pcm.set_cmap(colormaps["magma"])
    pcm.set_clim(cl)

    ax = axs[2, i]
    pcm = ax.pcolor(BkgOutput.glon, BkgOutput.glat, (Output[d] - BkgOutput[d]).T)
    ax.set_title(f"Difference {d}")
    fig.colorbar(pcm, ax=ax, label=Output[d].attrs["units"])
    cl = np.max(np.abs(pcm.get_clim()))
    pcm.set_clim(-cl, cl)
    pcm.set_cmap(colormaps["bwr"])

Example 4: Time Series

It is possible to calculate timeseries outputs from AIDA using tools like pandas and xarray. The following example calculates NmF2 at three locations for a 15 minute period. Note that AIDA operates using a 5-minute internal timestep, so it is not possible to resolve timescales finer than this. AIDA will not re-download cached outputs.

import pandas
import datetime
import xarray

glat = np.array([45.0,40.0,35.0])
glon = np.array([0,0,0])

times = pandas.date_range(start=np.datetime64('2024-10-01'), 
                          end=np.datetime64('2024-10-01T00:15'), 
                          freq=datetime.timedelta(minutes=0.5))

Output = []

for time in times:
    Model.fromAPI(time=time, model='AIDA', latency='rapid')
    Result = Model.calc(lat=glat, lon=glon, collapse_particles=True)
    
    Output.append(Result.expand_dims(time=[time]))

Output = xarray.concat(Output, dim='time')

plt.plot(Output['time'], Output['NmF2'], label=Output['glat'].data)
plt.ylabel(f"NmF2 ({Output['NmF2'].attrs['units']})")
plt.title('NmF2 at three latitudes, longitude=$0^o$')
plt.legend()

Example 5: Electron Density (Advanced Usage)

For some performance-critical applications, it is desireable to calculate only the electron density without additional overhead. This can be achieved with the AIDAState.calcNe() method. Like AIDAState.calc(), it accepts glat, glon, alt, and grid, but outputs only the electron density Ne as a numpy array. AIDAState.calcNe() is often 2 to 5 times faster than calling AIDAState.calc() for grid points with 1e3 to 1e4 unique lat/lon pairs.

Use of AIDAState.calcNe() is not recommended unless required, such as when calculating sample points for sTEC calculations.

def cartsph(x, y, z):
    # convert from ECEF to GEO
    rad = np.sqrt(x**2 + y**2 + z**2)
    lat = np.arctan2(z, np.sqrt(x**2 + y**2))
    lon = np.arctan2(y, x)

    return rad, lat, lon

def sphcart(rad, lat, lon):
    # convert from GEO to ECEF
    x = rad * np.cos(lat) * np.cos(lon)
    y = rad * np.cos(lat) * np.sin(lon)
    z = rad * np.sin(lat)

    return x, y, z

In this example slant Total Electron Content will be calculated between two points.

# Earth's Radius
Re = 6371e3

# receiver and satellite coordinates (ECEF)
rX, rY, rZ = sphcart(Re, np.deg2rad(45.0), np.deg2rad(0.0))
tX, tY, tZ = sphcart(Re+22000e3, np.deg2rad(-30.0), np.deg2rad(10.0))

# create linear sample points
t = np.linspace(0,1,1000)
x = rX + (tX-rX) * t
y = rY + (tY-rY) * t
z = rZ + (tZ-rZ) * t

# find geo coordinates of sample points
rad, glat_r, glon_r = cartsph(x, y, z)
alt = (rad - Re) * 1e-3
glat = np.rad2deg(glat_r)
glon = np.rad2deg(glon_r)

The electron density can be calculated with:

# output is in 1e11 m-3!
Ne = Model.calcNe(lat=glat, lon=glon, alt=alt, grid="1D")

# linear distance (m)
d = np.sqrt((x - rX) ** 2 + (y - rY) ** 2 + (z - rZ) ** 2)

# integrate along ray
sTEC = np.trapz(Ne, d)  # array([3.36102813e+17])

(back to top)

License

The AIDA interpreter is distributed under the MIT License. See LICENSE for more information.

(back to top)

Contact

Benjamin Reid - SERENE - University of Birmingham - b.reid@bham.ac.uk

Project Link: GitLab

(back to top)