spacepy.radbelt.RBmodel

class spacepy.radbelt.RBmodel(grid='L', NL=91, const_kp=False)[source]

1-D Radial diffusion class

This module contains a class for performing and visualizing 1-D radial diffusion simulations of the radiation belts.

Here is an example using the default settings of the model. Each instance must be initialized with (assuming import radbelt as rb):

>>> rmod = rb.RBmodel()

Next, set the start time, end time, and the size of the timestep:

>>> import datetime
>>> start = datetime.datetime(2003,10,14)
>>> end = datetime.datetime(2003,12,26)
>>> delta = datetime.timedelta(hours=1)
>>> rmod.setup_ticks(start, end, delta, dtype='UTC')

Now, run the model over the entire time range using the evolve method:

>>> rmod.evolve()

Finally, visualize the results:

>>> rmod.plot_summary()

Gaussian_source()

Gaussian source term added to radiation belt model.

add_Lmax(Lmax_model)

add last closed drift shell Lmax

add_Lpp(Lpp_model)

add last closed drift shell Lmax

add_PSD_obs([time, PSD, Lstar, satlist])

add PSD observations

add_PSD_twin([dt, Lt])

add observations from PSD database using the ticks list the arguments are the following:

add_omni([keylist])

add omni data to instance according to the tickrange in ticks

add_source([source, A, mu, sigma])

add source parameters A, mu, and sigma for the Gaussian source function

assimilate([method, inflation])

Assimilates data for the radiation belt model using the Ensemble Kalman Filter.

evolve()

calculate the diffusion in L at constant mu,K coordinates

get_DLL(Lgrid, params[, DLL_model])

Calculate DLL as a simple power law function (alpha*L**Bbta) using alpha/beta values from popular models found in the literature and chosen with the kwarg "DLL_model".

plot([Lmax, Lpp, Kp, Dst, clims, title, values])

Create a summary plot of the RadBelt object distribution function.

plot_obs([Lmax, Lpp, Kp, Dst, clims, title, ...])

Create a summary plot of the observations.

set_lgrid([NL])

Using NL grid points, create grid in L.

setup_ticks(start, end, delta[, dtype])

Add time information to the simulation by specifying a start and end time, timestep, and time type (optional).

Gaussian_source()[source]

Gaussian source term added to radiation belt model. The source term is given by the equation:

S = A exp{-(L-mu)^2/(2*sigma^2)}

with A=10^(-8), mu=5.0, and sigma=0.5 as default values

add_Lmax(Lmax_model)[source]

add last closed drift shell Lmax

add_Lpp(Lpp_model)[source]

add last closed drift shell Lmax

add_PSD_obs(time=None, PSD=None, Lstar=None, satlist=None)[source]

add PSD observations

Parameters:
timeTicktock datetime array

array of observation times

PSDlist of numpy arrays

PSD observational data for each time. Each entry in the list is a numpy array with the observations for the corresponding time

Lstarlist of numpy arrays

Lstar location of each PSD observations. Each entry in the list is a numpy array with the location of the observations for the corresponding time

satlistlist of satellite names
Returns:
outlist of dicts

Information of the observational data, where each entry contains the observations and locations of observations for each time specified in the time array. Each list entry is a dictionary with the following information:

TicksTicktock array

time of observations

Lstarnumpy array

location of observations

PSDnumpy array

PSD observation values

satlist of strings

satellite names

MUscalar value

Mu value for the observations

Kscalar value

K value for the observations

add_PSD_twin(dt=0, Lt=1)[source]

add observations from PSD database using the ticks list the arguments are the following:

dt = observation time delta in seconds Lt = observation space delta

add_omni(keylist=None)[source]

add omni data to instance according to the tickrange in ticks

add_source(source=True, A=1e-08, mu=5.0, sigma=0.5)[source]

add source parameters A, mu, and sigma for the Gaussian source function

assimilate(method='EnKF', inflation=0)[source]

Assimilates data for the radiation belt model using the Ensemble Kalman Filter. The algorithm used is the SVD method presented by Evensen in 2003 (Evensen, G., Ocean dynamics, 53, pp.343–367, 2003). To compensate for model errors, three inflation algorithms are implemented. The inflation methodology is specified by the ‘inflation’ argument, and the options are the following:

inflation == 0: Add model error (perturbation for the ensemble) around model state values only where observations are available (DEFAULT).

inflation == 1: Add model error (perturbation for the ensemble) around observation values only where observations are available.

inflation == 2: Inflate around ensemble average for EnKF.

Prior to assimilation, a set of data values has to be speficied by setting the start and end dates, and time step, using the setup_ticks funcion of the radiation belt model:

>>> import spacepy
>>> import datetime
>>> from spacepy import radbelt
>>> start = datetime.datetime(2002,10,23)
>>> end = datetime.datetime(2002,11,4)
>>> delta = datetime.timedelta(hours=0.5)
>>> rmod.setup_ticks(start, end, delta, dtype='UTC')

Once the dates and time step are specified, the data is added using the add_PSD function:

>>> rmod.add_PSD()

The observations are averaged over the time windows, whose interval is give by the time step.

Once the dates and data are set, the assimiation is performed using the ‘assimilate’ function:

>>> rmod.assimilate(inflation=1)

This function will add the PSDa values, which are the analysis state of the radiation belt using the observations within the dates. To plot the analysis simply use the plot funtion:

>>> rmod.plot(values=rmod.PSDa,clims=[-10,-6],Lmax=False,Kp=False,Dst=False)
evolve()[source]

calculate the diffusion in L at constant mu,K coordinates

get_DLL(Lgrid, params, DLL_model='BA2000')[source]

Calculate DLL as a simple power law function (alpha*L**Bbta) using alpha/beta values from popular models found in the literature and chosen with the kwarg “DLL_model”.

The calculated DLL is returned, as is the alpha and beta values used in the calculationp.

The output DLL is in units of units/day.

plot(Lmax=True, Lpp=False, Kp=True, Dst=True, clims=[0, 10], title=None, values=None)[source]

Create a summary plot of the RadBelt object distribution function. For reference, the last closed drift shell, Dst, and Kp are all included. These can be disabled individually using the corresponding Boolean kwargs.

The clims kwarg can be used to manually set the color bar range. To use, set it equal to a two-element list containing minimum and maximum Log_10 value to plot. Default action is to use [0,10] as the log_10 of the color range. This is good enough for most applications.

The title of the top most plot defaults to ‘Summary Plot’ but can be customized using the title kwarg.

The figure object and all three axis objects (PSD axis, Dst axis, and Kp axis) are all returned to allow the user to further customize the plots as necessary. If any of the plots are excluded, None is returned in their stead.

Examples

>>> rb.plot(Lmax=False, Kp=False, clims=[2,10], title='Good work!')

This command would create the summary plot with a color bar range of 100 to 10^10. The Lmax line and Kp values would be excluded. The title of the topmost plot (phase space density) would be set to ‘Good work!’.

plot_obs(Lmax=True, Lpp=False, Kp=True, Dst=True, clims=[0, 10], title=None, values=None)[source]

Create a summary plot of the observations. For reference, the last closed drift shell, Dst, and Kp are all included. These can be disabled individually using the corresponding boolean kwargs.

The clims kwarg can be used to manually set the color bar range. To use, set it equal to a two-element list containing minimum and maximum Log_10 value to plot. Default action is to use [0,10] as the log_10 of the color range. This is good enough for most applications.

The title of the top most plot defaults to ‘Summary Plot’ but can be customized using the title kwarg.

The figure object and all three axis objects (PSD axis, Dst axis, and Kp axis) are all returned to allow the user to further customize the plots as necessary. If any of the plots are excluded, None is returned in their stead.

Examples

>>> rb.plot_obs(Lmax=False, Kp=False, clims=[2,10], title='Observations Plot')

This command would create the summary plot with a color bar range of 100 to 10^10. The Lmax line and Kp values would be excluded. The title of the topmost plot (phase space density) would be set to ‘Good work!’.

set_lgrid(NL=91)[source]

Using NL grid points, create grid in L. Default number of points is 91 (dL=0.1).

setup_ticks(start, end, delta, dtype='ISO')[source]

Add time information to the simulation by specifying a start and end time, timestep, and time type (optional).

Examples

>>> start = datetime.datetime(2003,10,14)
>>> end = datetime.datetime(2003,12,26)
>>> delta = datetime.timedelta(hours=1)
>>> rmod.setup_ticks(start, end, delta, dtype='UTC')