Stingray Simulator (stingray.simulator)

Introduction

stingray.simulator provides a framework to simulate light curves with given variability distributions. In time series experiments, understanding the certainty is crucial to interpret the derived results in context of physical models. The simulator module provides tools to assess this uncertainty by simulating time series and spectral data.

Stingray simulator supports multiple methods to carry out these simulation. Light curves can be simulated through power-law spectrum, through a user-defined or pre-defined model, or through impulse responses. The module is designed in a way such that all these methods can be accessed using similar set of commands.

Note

stingray.simulator is currently a work-in-progress, and thus it is likely there will still be API changes in later versions of Stingray. Backwards compatibility support between versions will still be maintained as much as possible, but new features and enhancements are coming in future versions.

Getting started

The examples here assume that the following libraries and modules have been imported:

>>> import numpy as np
>>> from stingray import Lightcurve, sampledata
>>> from stingray.simulator import simulator, models

Creating a Simulator Object

Stingray has a simulator class which can be used to instantiate a simulator object and subsequently, perform simulations. We can pass on arguments to this class class to set the properties of the desired light curve.

The simulator object can be instantiated as:

>>> sim = simulator.Simulator(N=1024, mean=0.5, dt=0.125)

Here, N specifies the bins count of the simulated light curve, mean specifies the mean value, and dt is the time resolution. Additional arguments can be provided to specify the rms of the simulated light curve, or to account for the effect of red noise leakage.

Simulate Method

Stingray provides multiple ways to simulate a light curve. However, all these methods follow a common recipe:

>>> sim = simulator.Simulator(N=1024, mean=0.5, dt=0.125)
>>> lc = sim.simulate(2)

Using Power-Law Spectrum

When only an integer argument (beta) is provided to the simulate method, that integer defines the shape of the power law spectrum. Passing beta as 1 gives a flicker-noise distribution, while a beta of 2 generates a random-walk distribution.

from matplotlib import rcParams
rcParams['font.family'] = 'sans-serif'
rcParams['font.sans-serif'] = ['Tahoma']

import matplotlib.pyplot as plt
from stingray.simulator import simulator

# Instantiate simulator object
sim = simulator.Simulator(N=1024, mean=0.5, dt=0.125)
# Specify beta value
lc = sim.simulate(2)

plt.plot(lc.counts, 'g')
plt.title('Random-walk Distribution Simulation', fontsize='16')
plt.xlabel('Counts', fontsize='14', )
plt.ylabel('Flux', fontsize='14')
plt.show()

(Source code, png, hires.png, pdf)

_images/simulator-1.png

Using User-defined Model

Light curve can also be simulated using a user-defined spectrum, which can be passed on as a numpy array.

from matplotlib import rcParams
rcParams['font.family'] = 'sans-serif'
rcParams['font.sans-serif'] = ['Tahoma']

import matplotlib.pyplot as plt
from stingray.simulator import simulator

# Instantiate simulator object
sim = simulator.Simulator(N=1024, mean=0.5, dt=0.125)
# Define a spectrum
w = np.fft.rfftfreq(sim.N, d=sim.dt)[1:]
spectrum = np.power((1/w),2/2)
# Simulate
lc = sim.simulate(spectrum)

plt.plot(lc.counts, 'g')
plt.title('User-defined Model Simulation', fontsize='16')
plt.xlabel('Counts', fontsize='14')
plt.ylabel('Flux', fontsize='14')
plt.show()

(Source code, png, hires.png, pdf)

_images/simulator-2.png

Using Pre-defined Models

One of the pre-defined spectrum models can be used to simulate a light curve. In this case, model name and model parameters (as list iterable) need to be passed on as function arguments.

Using Impulse Response

In order to simulate a light curve using impulse response, we need the original light curve and impulse response. Stingray provides TransferFunction class which can be used to obtain time and energy averaged impulse response by passing in a 2-D intensity profile as the input. A detailed tutorial on obtaining impulse response is provided here.

Here, for the sake of simplicity, we use a simulated impulse response.

from matplotlib import rcParams
rcParams['font.family'] = 'sans-serif'
rcParams['font.sans-serif'] = ['Tahoma']

import matplotlib.pyplot as plt
from stingray import sampledata
from stingray.simulator import simulator

# Obtain a sample light curve
lc = sampledata.sample_data().counts
# Instantiate simulator object
sim = simulator.Simulator(N=1024, mean=0.5, dt=0.125)
# Obtain an artificial impulse response
ir = sim.relativistic_ir()
# Simulate
lc_new = sim.simulate(lc, ir)

plt.plot(lc_new.counts, 'g')
plt.title('Impulse Response based Simulation', fontsize='16')
plt.xlabel('Counts', fontsize='14')
plt.ylabel('Flux', fontsize='14')
plt.show()

(Source code, png, hires.png, pdf)

_images/simulator-3.png

Since, the new light curve is produced by the convolution of original light curveand impulse response, its length is truncated by default for ease of analysis. This can be changed, however, by supplying an additional parameter full. However, at times, we do not need to include lag delay portion in the output light curve. This can be done by changing the final function parameter to filtered. For a more detailed analysis on lag-frequency spectrum, follow the notebook here.

Channel Simulation

The simulator class provides the functionality to simulate light curves independently for each channel. This is useful, for example, when dealing with energy dependent impulse responses where we can create a di↵erent simulation channel for each energy range. The module provides options to count, retrieve and delete channels.:

>>> sim = simulator.Simulator(N=1024, mean=0.5, dt=0.125)
>>> sim.simulate_channel('3.5 - 4.5', 2)
>>> sim.count_channels()
1
>>> lc = sim.get_channel('3.5 - 4.5')
>>> sim.delete_channel('3.5 - 4.5')

Alternatively, assume that we have light curves in the simulated energy channels 3.5 - 4.5, 4.5 - 5.5 and 5.5 - 6.5. These channels can be retreived or deleted in single commands.

>>> sim.count_channels()
0
>>> sim.simulate_channel('3.5 - 4.5', 2)
>>> sim.simulate_channel('4.5 - 5.5', 2)
>>> sim.simulate_channel('5.5 - 6.5', 2)
>>> chans = sim.get_channels(['3.5 - 4.5','4.5 - 5.5','5.5 - 6.5'])
>>> sim.delete_channels(['3.5 - 4.5','4.5 - 5.5','5.5 - 6.5'])

Reference/API

class stingray.simulator.simulator.Simulator(dt=1, N=1024, mean=0, rms=1, red_noise=1, random_state=None)[source]

Methods to simulate and visualize light curves.

Parameters:

dt : int, default 1

time resolution of simulated light curve

N : int, default 1024

bins count of simulated light curve

mean : float, default 0

mean value of the simulated light curve

rms : float, default 1

fractional rms of the simulated light curve, actual rms is calculated by mean*rms

red_noise : int, default 1

multiple of real length of light curve, by which to simulate, to avoid red noise leakage

seed : int, default None

seed value for random processes

count_channels()[source]

Return total number of energy channels.

delete_channel(channel)[source]

Delete an energy channel.

delete_channels(channels)[source]

Delete multiple energy channels.

get_all_channels()[source]

Get lightcurves belonging to all channels.

get_channel(channel)[source]

Get lightcurve belonging to the energy channel.

get_channels(channels)[source]

Get multiple light curves belonging to the energy channels.

powerspectrum(lc, seg_size=None)[source]

Make a powerspectrum of the simulated light curve.

Parameters:

lc : lightcurve.Lightcurve object OR

iterable of lightcurve.Lightcurve objects The light curve data to be Fourier-transformed.

Returns:

power : numpy.ndarray

The array of normalized squared absolute values of Fourier amplitudes

static read(filename, format_='pickle')[source]

Imports Simulator object.

Parameters:

filename : str

Name of the Simulator object to be read.

format_ : str

Available option is ‘pickle.’

Returns:

object : Simulator object

relativistic_ir(t1=3, t2=4, t3=10, p1=1, p2=1.4, rise=0.6, decay=0.1)[source]

Construct a realistic impulse response considering the relativistic effects.

Parameters:

t1 : int

primary peak time

t2 : int

secondary peak time

t3 : int

end time

p1 : float

value of primary peak

p2 : float

value of secondary peak

rise : float

slope of rising exponential from primary peak to secondary peak

decay : float

slope of decaying exponential from secondary peak to end time

Returns:

h : numpy.ndarray

Constructed impulse response

simple_ir(start=0, width=1000, intensity=1)[source]

Construct a simple impulse response using start time, width and scaling intensity. To create a delta impulse response, set width to 1.

Parameters:

start : int

start time of impulse response

width : int

width of impulse response

intensity : float

scaling parameter to set the intensity of delayed emission corresponding to direct emission.

Returns:

h : numpy.ndarray

Constructed impulse response

simulate(*args)[source]

Simulate light curve generation using power spectrum or impulse response.

Parameters:

args

See examples below.

Returns:

lightCurve : LightCurve object

Examples

  • x = simulate(beta)
    For generating a light curve using power law spectrum.
simulate_channel(channel, *args)[source]

Simulate a lightcurve and add it to corresponding energy channel.

Parameters:

channel : str

range of energy channel (e.g., 3.5-4.5)

*args

see description of simulate() for details

Returns:

lightCurve : LightCurve object

write(filename, format_='pickle')[source]

Exports Simulator object.

Parameters:

filename : str

Name of the Simulator object to be created.

format_ : str

Available options are ‘pickle’ and ‘hdf5’.