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: [list, str, astropy.SkyCoords], optional
Sky Coordinates of the event.
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.
caustics: Caustics
Caustics for given model
data_ref: int or MulensData
Reference dataset. If int then gives index of reference dataset in datasets.

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).

parameters

ModelParameters

Model parameters.

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

is_static()

see MulensModel.modelparameters.ModelParameters.is_static()

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.
magnification(time, satellite_skycoord=None, gamma=0.0, flux_ratio_constraint=None, separate=False)

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.

flux_ratio_constraint: MulensData or str, optional

Data to constrain the flux ratio of sources in binary source models. Can be MulensData instance and in that case this dataset is used to find flux ratio via regression and this flux ratio is applied in calculation of effective magnification. If str is provided, then it indicates the bandpass and we use value set by set_source_flux_ratio_for_band().
separate: boolean, optional
For binary source models, return magnification of each source separately. Default is False and 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).
data_magnification

list

A list of magnifications calculated for every dataset time vector.

get_data_magnification(dataset)

Get the model magnification for a dataset.

Parameters :
dataset: MulensData
Dataset with epochs for which magnification will be given. Satellite and limb darkening information is taken into account.
Returns :
magnification_vector: np.ndarray
Values on magnification.
set_source_flux_ratio(ratio)

Sets flux ratio of sources for binary source models. If you also call set_source_flux_ratio_for_band(), then the value set here will be used when: 1) no band is specified, or 2) band is specified but flux ratio for given band was not specified.

Parameters :
ratio: float or None
The ratio of fluxes of source no. 2 to source no. 1, i.e., flux_source_2/flux_source_1. Setting it to None removes the internal information, i.e., flux ratio will be fitted via regression (unless specific value is provided for bandpass).
set_source_flux_ratio_for_band(band, ratio)

Sets flux ratio for binary source models for given band.

Parameters :
band: str
Band for which constraint is given.
ratio: float
ratio of fluxes of source no. 2 to source no. 1, i.e., flux_source_band_2/flux_source_band_1
datasets

list

datasets linked to given model

set_datasets(datasets, data_ref=0)

Set datasets property

Parameters :
datasets: list of MulensData
Datasets to be stored.

data_ref: int or, MulensData, optional

Reference dataset.
coords

see Coordinates

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.
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=0.0, flux_ratio_constraint=None, **kwargs)

Plot the model magnification curve.

Keywords :

see plot_lc()

gamma: see magnification()

satellite_skycoord: see plot_trajectory()

**kwargs – any arguments accepted by matplotlib.pyplot.plot().

plot_lc(times=None, t_range=None, t_start=None, t_stop=None, dt=None, n_epochs=None, data_ref=None, f_source=None, f_blend=None, subtract_2450000=False, subtract_2460000=False, flux_ratio_constraint=None, **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()

data_ref: int or a MulensData object

Reference dataset to scale the model to. See get_ref_fluxes()
f_source, f_blend: float
Explicitly specify the source and blend fluxes in a system where flux = 1 corresponds to MulensModel.utils.MAG_ZEROPOINT (= 22 mag).
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())

flux_ratio_constraint: instance of MulensData or str, optional

Option for binary source models only. Data or bandpass to constrain the flux ratio of sources.

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

Provide data_ref or (f_source, f_blend) if you want to plot in flux units different than last value of data_ref (defaults to the first dataset).

get_ref_fluxes(data_ref=None)

Get source and blending fluxes for the model by finding the best-fit values compared to data_ref.

Parameters :
data_ref: MulensData or int
Reference dataset. If int, corresponds to the index of the dataset in self.datasets. If None, than the first dataset will be used.
Returns :
f_source: np.ndarray
Sources’ flux; normally of size (1). If it is of size (1) for a double source model, then it is a sum of fluxes of both sources.
f_blend: float
blending flux

Determine the reference flux system from the datasets. The data_ref may either be a dataset or the index of a dataset (if Model.set_datasets() was previously called). If data_ref is not set, it will use the first dataset. If you call this without calling set_datasets() first, there will be an exception and that’s on you.

reset_plot_properties()

This function will be deprecated.

Resets internal plotting properties of all attached datasets.

plot_data(data_ref=None, show_errorbars=None, show_bad=None, color_list=None, marker_list=None, size_list=None, label_list=None, alpha_list=None, zorder_list=None, subtract_2450000=False, subtract_2460000=False, **kwargs)

Plot the data scaled to the model.

Keywords (all optional):
data_ref: see get_ref_fluxes()
If data_ref is not specified, uses the first dataset as the reference for flux scale.
show_errorbars: boolean or None
Do you want errorbars to be shown for all datasets? Default is None, which means the option is taken from each dataset plotting properties (for which default is True). If True, then data are plotted using matplotlib.errorbar(). If False, then data are plotted using matplotlib.scatter().
show_bad: boolean or None
Do you want data marked as bad to be shown? Default is None, which means the option is taken from each dataset plotting properties (for which default is False). If bad data are shown, then they are plotted with ‘x’ marker.
subtract_2450000, subtract_2460000: boolean
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_lc()).
**kwargs:
Passed to matplotlib plotting functions. Contrary to previous behavior, **kwargs are no longer remembered.
get_residuals(data_ref=None, type='mag', data=None)

Calculate the residuals from the model for each dataset at once, or just a single dataset.

Note: if residuals are returned in magnitudes, they are transformed to the magnitude system specified by data_ref, so only suitable for plotting.

Keywords :
data_ref: optional
see get_ref_fluxes()
type: str, optional
specify whether the residuals should be returned in magnitudes (‘mag’) or in flux (‘flux’). Default is ‘mag’.
data: MulensData, optional
dataset for which residuals are returned. If specified, then returned lists are single element.
Returns :
residuals: list
each element of the list is a np.array() with the residuals for the corresponding dataset.
errorbars: list
the scaled errorbars for each point. For plotting errorbars for the residuals.
plot_residuals(show_errorbars=None, color_list=None, marker_list=None, size_list=None, label_list=None, alpha_list=None, zorder_list=None, data_ref=None, subtract_2450000=False, subtract_2460000=False, show_bad=None, **kwargs)

Plot the residuals (in magnitudes) of the model.

For explanation of keywords, see doctrings in plot_data(). Note different order of keywords.

plot_trajectory(times=None, t_range=None, t_start=None, t_stop=None, dt=None, n_epochs=None, caustics=False, show_data=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.
show_data: boolean
Mark epochs of data. Use plot_source_for_datasets() for more control of how data epochs are marked.
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().
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_caustics(n_points=5000, epoch=None, **kwargs)

Plot the caustic structure. See MulensModel.caustics.Caustics.plot().

Additional parameters :
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.
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.
plot_source_for_datasets(**kwargs)

Plot source positions for all linked datasets. Colors used for each dataset are the same used for plotting photometry.

Parameters:
**kwargs:
see plot_source()
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]. For binary source models, respectively, smaller and larger of t_0_1/2 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)
set_default_magnification_method(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(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
Which source given methods apply to? Accepts 1, 2, or None (i.e., all sources).
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}}.
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.
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_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
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
bandpasses

list

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

fit

MulensModel.fit.Fit

MulensModel.fit.Fit instance recently used. It gives access to source and blending fluxes.