MulensModel.uniformcausticsampling module¶
- class MulensModel.uniformcausticsampling.UniformCausticSampling(s, q, n_points=10000)¶
Bases:
object
Uniform sampling of a binary lens caustic. Note that calculations take some time for given (s, q). Keep that in mind, when optimizing your fitting routine.
- Arguments :
- s: float
Separation of the two lens components relative to the Einstein ring size.
- q: float
Mass ratio of the two lens components.
- n_points: int
Number of points used for internal integration. Default value should work fine.
Instead of standard parameters (t_0, u_0, t_E, alpha), here we use four other parameters: two epochs of caustic crossing (t_caustic_in, t_caustic_out) and two curvelinear coordinates of caustic crossing (x_caustic_in, x_caustic_out).
The curvelinear coordinate, x_caustic, is defined so that going from 0 to 1 draws all caustics for given separation and mass ratio. We use 0-1 range, which is a different convention than in the papers cited below (we also use different symbols for epochs of caustic crossing and curvelinear coordinates). For a wide topology (i.e., 2 caustics), there is a value between 0 and 1 (called
x_caustic_sep
) which separates the caustics and a trajectory exists only if x_caustic_in and x_caustic_out correspond to the same caustic, i.e., both are smaller thanx_caustic_sep
or both are larger thanx_caustic_sep
. For a close topology (i.e., 3 caustics), there are two such separating values.For description of the curvelinear coordinates, see:
Cassan A. et al. 2010 A&A 515, 52 “Bayesian analysis of caustic-crossing microlensing events”
In order to visualize the curvelinear coordinates, you can run a code like:
import matplotlib.pyplot as plt import numpy as np sampling = UniformCausticSampling(s=1.1, q=0.3) color = np.linspace(0., 1., 200) points = [sampling.caustic_point(c) for c in color] x = [p.real for p in points] y = [p.imag for p in points] plt.scatter(x, y, c=color) plt.axis('equal') plt.colorbar() plt.show()
This will show an intermediate topology. Change s=1.1 to s=2. to plot a wide topology, or to s=0.7 to plot a close topology.
To be specific, the central caustics are plotted counter-clockwise and x_caustic=0. corresponds to right-hand point where the caustic crosses the X-axis. For a wide topology, the planetary caustic is plotted in a similar way. For a close topology, the lower planetary caustic is plotted counter-clockwise and the upper planetary caustic is symmetric, thus plotted clockwise. For planetary caustics in a close topology, the zero-point of x_caustic values is defined in a very complicated way, however it is a smooth function of s and q.
For more advanced fitting of binary lens events see:
- get_standard_parameters(x_caustic_in, x_caustic_out, t_caustic_in, t_caustic_out)¶
Get standard binary lens parameters (i.e.,
t_0
,u_0
,t_E
,alpha
; seeModelParameters
) based on provided curvelinear parameters.Note that this function quite frequently raises
ValueError
exception. This is because not all (s
,q
,x_caustic_in
, andx_caustic_out
) correspond to real trajectories. The returned values are in conventions used byModel
.- Keywords :
- x_caustic_in: float
Curvelinear coordinate of caustic entrance. Must be in (0, 1) range.
- x_caustic_out: float
Curvelinear coordinate of caustic exit. Must be in (0, 1) range.
- t_caustic_in: float
Epoch of caustic entrance.
- t_caustic_out: float
Epoch of caustic exit.
- Returns :
- parameters: dict
Dictionary with standard binary parameters, i.e, keys are
t_0
,u_0
,t_E
, andalpha
.
- get_x_in_x_out(u_0, alpha)¶
Calculate where given trajectory crosses the caustic.
- Parameters :
- u_0: float
The parameter u_0 of source trajectory, i.e., impact parameter.
- alpha: float
Angle defining the source trajectory.
- Returns :
- x_caustic_points: list of float
Caustic coordinates of points where given trajectory crosses the caustic. The length is 0, 2, 4, or 6. Note that if there are 4 or 6 points, then only some pairs will produce real trajectories.
- get_uniform_sampling(n_points, n_min_for_caustic=10, caustic=None)¶
Sample uniformly (x_caustic_in, x_caustic_out) space according to Jacobian and requirement that x_caustic_in corresponds to caustic entrance, and x_caustic_out corresponds to caustic exit. The Jacobian is defined by Eq. 23 of:
Cassan A. et al. 2010 A&A 515, 52 “Bayesian analysis of caustic-crossing microlensing events”
and the above requirement is defined under Eq. 27 of that paper.
Relative number of points per caustic is not yet specified. Points do not repeat. This function is useful for sampling starting distribution for model fitting. For example sampling see Cassan et al. (2010) bottom panel of Fig. 1.
- Parameters :
- n_points: int
number of points to be returned
- n_min_for_caustic: int
minimum number of points in each caustic
- caustic: int or None
Select which caustic will be sampled. None means all caustics. Can be 1, 2, or 3 but has to be <=
n_caustics
.
- Returns :
- x_caustic_in: np.ndarray
Randomly drawn entrance points.
- x_caustic_out: np.ndarray
Corresponding randomly drawn exit points.
- jacobian(x_caustic_in, x_caustic_out)¶
Evaluates Eq. 23 from Cassan et al. (2010) with condition under Eq. 27.
- Parameters :
- x_caustic_in: float
Point of caustic entrance.
- x_caustic_out: float
Point of caustic exit.
- Returns :
- jacobian: float
Value of Jacobian. Returns 0. if trajectory does not exist.
- check_valid_trajectory(x_caustic_in, x_caustic_out)¶
Check if given (x_caustic_in, x_caustic_out) define an existing trajectory. An obvious case, when they don’t is when both caustic points are on the same fold, but other cases exists.
- Parameters :
- x_caustic_in: float
Coordinate of putative caustic entrance.
- x_caustic_out: float
Coordinate of putative caustic exit.
- Returns :
- check: bool
True if input defines a trajectory, False if it does not.
- caustic_point(x_caustic)¶
Calculate caustic position corresponding to given x_caustic.
- Keywords :
- x_caustic: float
Curvelinear coordinate of the point considered. Has to be in 0-1 range.
- Returns :
- point: numpy.complex128
Caustic point in complex coordinates.
- which_caustic(x_caustic)¶
Indicates on which caustic given point is.
- Keywords :
- x_caustic: float
Curvelinear coordinate to be checked
- Returns :
- i_caustic: int
Number indicating the caustic:
1
- central caustic,2
- planetary caustic; for close configuration it is the lower of the two planetary caustics,3
- upper planetary caustic.
- property n_caustics¶
int
Number of caustics: 1 for resonant topology, 2 for wide topology, or 3 for close topology.
- property s¶
float
separation of the two lens components relative to Einstein ring size
- property q¶
float
Mass ratio.