MulensModel.model module

class MulensModel.model.Model(parameters=None, coords=None, ra=None, dec=None, ephemerides_file=None)

Bases: object

A Model for a microlensing event with the specified parameters.

Arguments :

parameters: dictionary, ModelParameters

coords: str, astropy.SkyCoords, MulensModel.Coordinates, optional

Sky coordinates of the event. If type is str, then it is assumed that the units are hour angle and degrees for RA and Dec, respectively.

ra, dec: str, optional

Sky coordinates of the event.

ephemerides_file: str, optional

Specify name of the file with satellite ephemerides. See MulensData for more details. Note that if you provide file name here, then it will affect all calculations for this model. In most cases, you want to combine ground-based and satellite data and in those cases set ephemerides_file for specific MulensData instance to pass satellite information.

Attributes :
ephemerides_file: str

Name of file with satellite ephemerides.

Default values for parallax are all True. Use parallax() to turn different parallax effects ON/OFF. If using satellite parallax, you may also specify an ephemerides_file (see MulensData).

Note that you can print an instance of Model, which shows you parameters in a nice way, e.g.,

model = Model(parameters={'t_0': 2456789.0, ....})
print(model)

This will provide information on parameter values, coordinates, methods used for magnification calculations, and limb-darkening coefficients.

plot_magnification(times=None, t_range=None, t_start=None, t_stop=None, dt=None, n_epochs=None, subtract_2450000=False, subtract_2460000=False, satellite_skycoord=None, gamma=None, source_flux_ratio=None, **kwargs)

Plot the model magnification curve.

Keywords :

see plot_lc()

gamma:

see get_magnification()

satellite_skycoord:

see plot_trajectory()

source_flux_ratio: float

If the model has two sources, source_flux_ratio is the ratio of source_flux_2 / source_flux_1

**kwargs:

any arguments accepted by matplotlib.pyplot.plot().

get_lc(times=None, t_range=None, t_start=None, t_stop=None, dt=None, n_epochs=None, source_flux=None, blend_flux=None, source_flux_ratio=None, gamma=None, bandpass=None)

Calculate model light curve in magnitudes.

Keywords :
times: [float, list, numpy.ndarray]

a list of times at which to plot the magnifications

t_range, t_start, t_stop, dt, n_epochs: see set_times()

source_flux: float or list

Explicitly specify the source flux(es) in a system where flux = 1 corresponds to MulensModel.utils.MAG_ZEROPOINT (= 22 mag). If the model has n_source > 1, source_flux may be specified as a list: one value for each source. Alternatively, if source_flux is specified as a float, source_flux_ratio should also be specified. Then, source_flux is taken to be the flux of the first source, and the other source fluxes are derived using source_flux_ratio.

blend_flux: float

Explicitly specify the blend flux in a system where flux = 1 corresponds to MulensModel.utils.MAG_ZEROPOINT (= 22 mag).

source_flux_ratio: float, Optional

If the model has two sources, source_flux_ratio is the ratio of source_flux_2 / source_flux_1.

gamma, bandpass:

see get_magnification()

Returns :
magnitudes: numpy.ndarray

Magnitude values for each epoch.

plot_lc(times=None, t_range=None, t_start=None, t_stop=None, dt=None, n_epochs=None, source_flux=None, blend_flux=None, source_flux_ratio=None, gamma=None, bandpass=None, subtract_2450000=False, subtract_2460000=False, phot_fmt='mag', **kwargs)

Plot the model light curve in magnitudes.

Keywords :
times: [float, list, numpy.ndarray]

a list of times at which to plot the magnifications

t_range, t_start, t_stop, dt, n_epochs:

see set_times()

source_flux, blend_flux, source_flux_ratio:

see get_lc()

gamma, bandpass:

see get_magnification()

subtract_2450000, subtract_2460000: boolean, optional

If True, subtracts 2450000 or 2460000 from the time axis to get more human-scale numbers. If using, make sure to also set the same settings for all other plotting calls (e.g. plot_data())

phot_fmt: str

Specifies whether the photometry is plotted in magnitude or flux space. Accepts either ‘mag’ or ‘flux’. Default = ‘mag’.

**kwargs:

any arguments accepted by matplotlib.pyplot.plot().

plot_caustics(n_points=5000, epoch=None, **kwargs)

Plot the caustic structure. See MulensModel.caustics.Caustics.plot() for binary lenses. For a single lens it just marks (0, 0) point and the first two parameters are ignored.

Additional parameters :
n_points: int, optional

The number of points to calculate along the caustic. Defaults to 5000.

epoch: float, optional

Epoch for which separation s will be used. Important for models with orbital motion. Defaults to t_0_kep, which defaults to t_0.

**kwargs:

keywords accepted by matplotlib.pyplot.scatter()

update_caustics(epoch=None)

Updates caustics property for given epoch.

Parameters :
epoch: float

For orbital motion models, epoch for which separation s is calculated to calculate caustics. Defaults to t_0_kep, which defaults to t_0.

plot_trajectory(times=None, t_range=None, t_start=None, t_stop=None, dt=None, n_epochs=None, caustics=False, arrow=True, satellite_skycoord=None, arrow_kwargs=None, **kwargs)

Plot the source trajectory.

Keywords (all optional) :

times, t_range, t_start, t_stop, dt, n_epochs:

May all be used to specify exactly when to plot the source trajectory. See also plot_lc() and set_times().

caustics: boolean

plot the caustic structure in addition to the source trajectory. default=False (off). For finer control of plotting features, e.g. color, use plot_caustics() instead.

arrow: boolean

Show the direction of the source motion. Default is True.

satellite_skycoord: astropy.SkyCoord or MulensModel.satelliteskycoord.SatelliteSkyCoord

Allows the user to specify that the trajectory is calculated for a satellite. If astropy.SkyCoord object is provided, then these are satellite positions for all epochs. See also get_satellite_coords()

arrow_kwargs: dict

Kwargs that are passed to pyplot.arrow(). If no color is given here, then we use one specified in **kwargs and if nothing is there, then we use black. The size of the arrow is determined based on limits of current axis. If those are not adequate, then change the size by specifying width keyword and maybe other as well. Note that arrow_kwargs are of dict type and are different than **kwargs.

**kwargs

Controls plotting features of the trajectory. It’s passed to pyplot.plot().

Note that in order to have equal scaling of both axis (i.e., make circles look circular), you have to call appropriate pyplot command. This can be one of these commands:

pyplot.axis('equal')
pyplot.axis('scaled')
pyplot.axis('square')
pyplot.gca().set_aspect('equal')

They have slightly different behavior.

plot_source(times=None, **kwargs)

Plot source: circles of the radius rho at positions corresponding to source positions at times. When the rho is not defined, then X symbols are plotted.

Parameters:
times: float or np.ndarray

epochs for which source positions will be plotted

**kwargs:

Keyword arguments passed to matplotlib.Circle. Examples: color='red', fill=False, linewidth=3, alpha=0.5. When the rho is not defined, then keyword arguments are passed to matplotlib.plot.

Note that it is likely that with default axis scaling, the circles may be plotted as ellipses. To mitigate it, use: plt.gca().set_aspect('equal') or plt.axis('equal') (the other possible options are 'scaled' or 'square').

get_trajectory(times, satellite_skycoord=None)

Get the trajectory of the source for the given set of times. For multi-source models, one trajectory per source is returned.

Parameters :
times: np.ndarray, list of floats, or float

Epochs for which source positions are requested.

satellite_skycoord: astropy.SkyCoord

Allows the user to specify that the trajectory is calculated for a satellite. If astropy.SkyCoord object is provided, then these are satellite positions for all epochs. See also get_satellite_coords()

Returns :
trajectories: :py:class:`~MulensModel.trajectory.Trajectory object or a list of them

Single object for single source model, a list otherwise.

set_times(t_range=None, t_start=None, t_stop=None, dt=None, n_epochs=1000)

Return a list of times. If no keywords are specified, default is 1000 epochs from [\(t_0 - 1.5 * t_E\), \(t_0 + 1.5 * t_E\)] range.

For multi-source models, respectively, minimum and maximum of t_0_N values are used.

Parameters (all optional):
t_range: [list, tuple]

A range of times of the form [t_start, t_stop]

t_start, t_stop: float

a start or stop time.

dt: float

the interval spacing between successive points

n_epochs: int

the number of epochs (evenly spaced)

Returns :
times: np.ndarray

Vector of epochs.

set_magnification_methods(methods, source=None)

Sets methods used for magnification calculation. See MagnificationCurve for a list of implemented methods.

Parameters :
methods: list

List that specifies which methods (str) should be used when (float values for Julian dates). Given method will be used for times between the times between which it is on the list, e.g.,

methods = [
    2455746., 'Quadrupole', 2455746.6, 'Hexadecapole',
    2455746.7, 'VBBL', 2455747., 'Hexadecapole',
    2455747.15, 'Quadrupole', 2455748.]
source: int or None, optional

Which source do the given methods apply to? Accepts 1, 2, or None (i.e., all sources). Default is None

get_magnification_methods(source=None)

Gets methods used for magnification calculation. See set_magnification_methods()

Parameters :
source: int or None, optional

Which source do the given methods apply to? Accepts 1, 2, or None (i.e., all sources). Default is None.

property methods

list of methods used for magnification calculation or dict of lists if there are multiple sources.

property default_magnification_method

Stores information on method to be used, when no method is directly specified. See MagnificationCurve for a list of implemented methods.

Parameters:
method: str

Name of the method to be used.

set_magnification_methods_parameters(methods_parameters)

Set additional parameters for magnification calculation methods.

Parameters :
methods_parameters: dict

Dictionary that for method names (keys) returns dictionary in the form of **kwargs that are passed to given method, e.g., {'VBBL': {'accuracy': 0.005}}.

get_magnification_methods_parameters(method)

Get additional parameters for a specific magnification calculation method or methods.

Parameters :
method: str, list

Name of method or a list of the names for which parameters will be returned.

Returns :
method_parameters: dict

see set_magnification_methods_parameters()

set_limb_coeff_gamma(bandpass, coeff)

Store gamma limb darkening coefficient for given band. See also LimbDarkeningCoeffs.

Parameters :
bandpass: str

Bandpass for the coefficient you provide.

coeff: float

Value of the coefficient.

get_limb_coeff_gamma(bandpass)

Get gamma limb darkening coefficient for given band.

Parameters :
bandpass: str

Bandpass for which coefficient will be provided.

Returns :
gamma: float

limb darkening coefficient

set_limb_coeff_u(bandpass, coeff)

Store u limb darkening coefficient for given band. See also MulensModel.limbdarkeningcoeffs.LimbDarkeningCoeffs.

Parameters :
bandpass: str

Bandpass for which coefficient you provide.

coeff: float

Value of the coefficient.

get_limb_coeff_u(bandpass)

Get u limb darkening coefficient for given band.

Parameters :
bandpass: str

Bandpass for which coefficient will be provided.

Returns :
u: float

limb darkening coefficient

parallax(earth_orbital=None, satellite=None, topocentric=None)

Specifies the types of the microlensing parallax that will be included in calculations.

Parameters :

earth_orbital: boolean, optional

Do you want to include the effect of Earth motion about the Sun? Default is False.

satellite: boolean, optional

Do you want to include the effect of difference due to separation between the Earth and satellite? Note that this separation changes over time. Default is False.

topocentric: boolean, optional

Do you want to include the effect of different positions of observatories on the Earth? Default is False. Note that this is significant only for very high magnification events and if high quality datasets are analyzed. Hence, this effect is rarely needed. Not Implemented yet.

get_parallax()

Returns dict that specifies the types of the microlensing parallax that are included in calculations.

Returns :
parallax: dict

For keys 'earth_orbital', 'satellite', and 'topocentric' returns bool.

get_satellite_coords(times)

Get astropy.SkyCoord object that gives satellite positions for given times. see also MulensModel.satelliteskycoord.SatelliteSkyCoord

Parameters :
times: np.ndarray or list

Epochs for which satellite position is requested.

Returns :
satellite_skycoord: astropy.SkyCoord

SkyCoord giving satellite positions. The parameter representation is set to ‘spherical’. If ephemerides_file is not set, returns None.

get_magnification(time, satellite_skycoord=None, gamma=None, bandpass=None, source_flux_ratio=None, separate=None)

Calculate the model magnification for the given time(s).

Parameters :
time: np.ndarray, list of floats, or float

Times for which magnification values are requested.

satellite_skycoord: astropy.coordinates.SkyCoord, optional

SkyCoord object that gives satellite positions. Must be the same length as time parameter. Use only for satellite parallax calculations.

gamma: float, optional

The limb-darkening coefficient in gamma convention. Default is 0 which means no limb darkening effect.

bandpass: str, optional

The bandpass for setting the limb-darkening coefficient. Expects that you have used set_limb_coeff_gamma() or set_limb_coeff_u(). Only ONE of ‘gamma’ or ‘bandpass’ may be specified.

source_flux_ratio: float, list

If the model has 2 sources, source_flux_ratio is the ratio of source_flux_2 / source_flux_1. For N sources, source_flux_ratio a list of the ratios of source_flux_i / source_flux_1 where i = 2, N.

separate: boolean, optional

For multi source models, return magnification of each source separately. Defaults to True if source_flux_ratio is provided and False otherwise (then only effective magnification is returned).

Returns :
magnification: np.ndarray

A vector of calculated magnification values. For binary source models, the effective magnification is returned (unless separate=True).

get_magnification_curve(time, satellite_skycoord, gamma)

Create a MagnificationCurve object for a given set of times.

Parameters :
time: np.ndarray, list of floats, or float

Times for which magnification values are requested.

satellite_skycoord: astropy.coordinates.SkyCoord, optional

SkyCoord object that gives satellite positions. Must be the same length as time parameter. Use only for satellite parallax calculations.

gamma: float, optional

The limb-darkening coefficient in gamma convention. Default is 0 which means no limb darkening effect.

Return:

py:class:~MulensModel.magnificationcurve.MagnificationCurve

get_magnification_curves(time, satellite_skycoord, gamma)

Create a list of MagnificationCurve objects for multiple sources, given a set of times.

Parameters :
time: np.ndarray, list of floats, or float

Times for which magnification values are requested.

satellite_skycoord: astropy.coordinates.SkyCoord, optional

SkyCoord object that gives satellite positions. Must be the same length as time parameter. Use only for satellite parallax calculations.

gamma: float, optional

The limb-darkening coefficient in gamma convention. Default is 0 which means no limb darkening effect.

Return:

list of py:class:~MulensModel.magnificationcurve.MagnificationCurve

property caustics

Caustics

Caustics for given model

property parameters

ModelParameters

Model parameters.

property n_lenses

int

number of objects in the lens system

property n_sources

int

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

is_static()

see MulensModel.modelparameters.ModelParameters.is_static()

property coords

see Coordinates

property bandpasses

list

List of all bandpasses for which limb darkening coefficients are set.