MulensModel.modelparameters module

MulensModel.modelparameters.which_parameters(*args)

Prints information on valid parameter combinations that can be used to define a model or a particular effect. May be called with no arguments (returns information on many types of models) or with one argument referring to a specific model (e.g. PSPL) or effect (e.g. parallax).

Valid arguments: str

Model types: ‘PSPL’, ‘FSPL’, ‘PSBL’, ‘FSBL’

Effects: ‘point lens’, ‘binary lens’, ‘finite source’, ‘parallax’, ‘lens orbital motion’

class MulensModel.modelparameters.ModelParameters(parameters)

Bases: object

A class for the basic microlensing model parameters (t_0, u_0, t_E, rho, s, q, alpha, pi_E). Can handle point lens or binary lens. The pi_E assumes NE coordinates (Parallel, Perpendicular coordinates are not supported).

Arguments :
parameters: dictionary
A dictionary of parameters and their values. See which_parameters() for valid parameter combinations.
Attributes :
parameters: dictionary
A dictionary of parameters and their values. Do not use it to change paramter values, instead use e.g.: model_parameters.u_0 = 0.1 or setattr(model_parameters, 'u_0', 0.1).
Example:
Define a point lens model:
params = ModelParameters({'t_0': 2450000., 'u_0': 0.3, 't_E': 35.})
Then you can print the parameters:
print(params)
t_0

float

The time of minimum projected separation between the source and the lens center of mass.

u_0

float

The minimum projected separation between the source and the lens center of mass.

t_star

float

t_star = rho * t_E = source radius crossing time

“day” is the default unit. Can be set as float or astropy.Quantity, but always returns float in units of days.

t_eff

float

t_eff = u_0 * t_E = effective timescale

“day” is the default unit. Can be set as float or astropy.Quantity, but always returns float in units of days.

t_E

float

The Einstein timescale. “day” is the default unit. Can be set as float or astropy.Quantity, but always returns float in units of days.

rho

float

source size as a fraction of the Einstein radius

alpha

astropy.Quantity

The angle of the source trajectory relative to the binary lens axis (or primary-secondary axis). Measured counterclockwise, i.e., according to convention advocated by Skowron et al. 2011 (ApJ, 738, 87), but shifted by 180 deg. May be set as a float –> assumes “deg” is the default unit. Regardless of input value, returns value in degrees.

q

float

mass ratio of the two lens components. Only 2 bodies allowed.

s

float

separation of the two lens components relative to Einstein ring size

pi_E

list of floats

The microlensing parallax vector. Must be set as a vector/list (i.e. [pi_E_N, pi_E_E]). To get the magnitude of pi_E, use pi_E_mag

pi_E_N

float

The North component of the microlensing parallax vector.

pi_E_E

float

The East component of the microlensing parallax vector.

t_0_par

float

The reference time for the calculation of parallax. If not set explicitly, set t_0_par = t_0.

pi_E_mag

float

The magnitude of the microlensing parallax vector.

ds_dt

astropy.Quantity

Change rate of separation s in 1/year. Can be set as AstroPy.Quantity or as float (1/year is assumed default unit). Regardless of input value, returns value in 1/year.

dalpha_dt

astropy.Quantity

Change rate of angle alpha in deg/year. Can be set as AstroPy.Quantity or as float (deg/year is assumed default unit). Regardless of input value, returns value in deg/year.

t_0_kep

float

The reference time for the calculation of lens orbital motion. If not set explicitly, assumes t_0_kep = t_0.

t_0_1

float

The time of minimum projected separation between the source no. 1 and the lens center of mass.

t_0_2

float

The time of minimum projected separation between the source no. 2 and the lens center of mass.

u_0_1

float

The minimum projected separation between the source no. 1 and the lens center of mass.

u_0_2

float

The minimum projected separation between the source no. 2 and the lens center of mass.

t_star_1

float

t_star_1 = rho_1 * t_E_1 = source no. 1 radius crossing time

“day” is the default unit. Can be set as float or astropy.Quantity, but always returns float in units of days.

t_star_2

float

t_star_2 = rho_2 * t_E_2 = source no. 2 radius crossing time

“day” is the default unit. Can be set as float or astropy.Quantity, but always returns float in units of days.

rho_1

float

source no. 1 size as a fraction of the Einstein radius

rho_2

float

source no. 2 size as a fraction of the Einstein radius

get_s(epoch)

Returns the value of separation s at a given epoch or epochs (if orbital motion parameters are set).

Arguments :
epoch: float, list, np.ndarray
The time(s) at which to calculate s.
Returns :
separation: float or np.ndarray
Value(s) of separation for given epochs.
get_alpha(epoch)

Returns the value of angle alpha at a given epoch or epochs (if orbital motion parameters are set).

Arguments :
epoch: float, list, np.ndarray
The time(s) at which to calculate alpha.
Returns :
separation: astropy.Quantity
Value(s) of angle for given epochs in degrees
gamma_parallel

astropy.Quantity

Parallel component of instantaneous velocity of the secondary relative to the primary in 1/year. It is parallel to the primary-secondary axis. Equals ds_dt/s. Cannot be set.

gamma_perp

astropy.Quantity

Perpendicular component of instantaneous velocity of the secondary relative to the primary. It is perpendicular to the primary-secondary axis. It has sign opposite to dalpha_dt and is in rad/yr, not deg/yr. Cannot be set.

gamma

astropy.Quantity

Instantaneous velocity of the secondary relative to the primary in 1/year. Cannot be set.

is_finite_source()

Checks if model has finite source. For binary source models it checks if either of the sources is finite.

Returns:
is_finite_source: boolean
True if at least one source has finite size.
is_static()

Checks if model is static, i.e., orbital motion parameters are not set.

Returns :
is_static: boolean
True if dalpha_dt or ds_dt are set.
n_lenses

int

number of objects in the lens system

n_sources

int

number of luminous sources; it’s possible to be 1 for xallarap model

source_1_parameters

ModelParameters

Parameters of source 1 in multi-source model.

Do not change returned values. To change parameters of the source 1, simply change the parameters of double source instance.

source_2_parameters

ModelParameters

Parameters of source 2 in multi-source model.

Do not change returned values. To change parameters of the source 1, simply change the parameters of double source instance.

as_dict()

Give parameters as a dict.

Returns :
dictionary: dict
The dictionary of model parameters.