# PyTransit: fast and easy exoplanet transit modelling in Python.
# Copyright (C) 2010-2019 Hannu Parviainen
#
# This program is free software: you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with this program. If not, see <https://www.gnu.org/licenses/>.
#
# This program is free software: you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with this program. If not, see <https://www.gnu.org/licenses/>.
from typing import Union, Optional, List
from numpy import ones, ndarray, asarray, zeros, unique, atleast_1d, issubdtype, integer, float64, full
from ..orbits.orbits_py import ta_ip_calculate_table
[docs]
class TransitModel(object):
"""Base class for the PyTransit exoplanet transit light curve models.
`TransitModel` defines the data-setup half of the interface that every PyTransit transit
model shares, and declares the evaluation half that the subclasses implement. Using any
model follows the same three steps
1. create the model, giving it any model-specific settings,
2. call `set_data` to tell the model *when* the observations were made, and
3. call `evaluate` to compute the fluxes for a set of physical parameters.
The data setup is separated from the evaluation because the observation times change
rarely while the parameters change constantly. Everything that can be precomputed from
the times alone -- index arrays, supersampling offsets, and so on -- is computed once in
`set_data` and reused by every subsequent `evaluate` call. This is what makes the models
efficient inside an optimiser or a sampler.
Attributes
----------
time : ndarray
Mid-exposure times, set by `set_data`.
lcids : ndarray
Light curve index for each exposure, one per element of `time`.
pbids : ndarray
Passband index for each light curve, one per light curve.
epids : ndarray
Epoch index for each light curve, one per light curve.
nsamples : ndarray
Number of supersamples per exposure, one per light curve.
exptimes : ndarray
Exposure time in days, one per light curve.
npt : int
Number of exposures.
nlc : int
Number of light curves.
npb : int
Number of passbands.
Notes
-----
This class is not meant to be used directly. See :class:`~pytransit.models.roadrunner.rrmodel.RoadRunnerModel`
for the recommended general-purpose model.
"""
[docs]
def __init__(self) -> None:
# Declare the basic arrays
# ------------------------
self.time_id: Optional[int] = None
self.time: Optional[ndarray] = None
self.lcids: Optional[ndarray] = None
self.pbids: Optional[ndarray] = None
self.nsamples: Optional[ndarray] = None
self.exptimes: Optional[ndarray] = None
self.epids: Optional[ndarray] = None
self.nlc: int = 0
self.npt: int = 0
self.npb: int = 0
[docs]
def set_data(self, time: Union[ndarray, List],
lcids: Optional[Union[ndarray, List]] = None,
pbids: Optional[Union[ndarray, List]] = None,
nsamples: Optional[Union[ndarray, List]] = None,
exptimes: Optional[Union[ndarray, List]] = None,
epids: Optional[Union[ndarray, List]] = None) -> bool:
"""Set the data for the transit model.
Parameters
----------
time : array-like
Array of mid-exposure times for which the model will be evaluated.
lcids : array-like, optional
Array of integer light curve indices. Must have the same size as the time array.
pbids : array-like, optional
Array of passband indices, one per light curve. Must satisfy `pbids.size == unique(lcids).size`.
nsamples : int or array-like, optional
Number of samples per exposure. Can either be an integer, in which case all the light curves will have the
same supersampling rate, or an array of integers, in which case each light curve can have a different rate.
exptimes : float or array-like, optional
Exposure times, again either for all the modelled data, or one value per light curve.
epids : array-like, optional
Epoch indices that can be used to link a light curve to a specific zero epoch and period (for TTV calculations).
Returns
-------
bool
Whether the data were set. False if the call was a no-op because `time` is the array the
model already holds and nothing else was given, which lets a subclass skip the work it
derives from the data.
"""
# Time samples
# ------------
if id(time) == self.time_id and lcids is None and pbids is None and nsamples is None and exptimes is None and epids is None:
return False
self.time_id = id(time)
self.time = asarray(time, float64)
self.npt = self.time.size
# Light curve indices
# -------------------
# The light curve a datapoint belongs to.
self.lcids = asarray(lcids) if lcids is not None else zeros(self.npt, 'int')
self.nlc = unique(self.lcids).size
if not issubdtype(self.lcids.dtype, integer):
raise ValueError(f"The light curve indices must be given as integers instead of {self.lcids.dtype}.")
if self.lcids.size != self.npt:
raise ValueError(f"Light curve index array size ({self.lcids.size}) should equal to the number of datapoints ({self.npt}).")
# Passband indices
# ----------------
# The passband a light curve belongs to.
self.pbids = asarray(pbids) if pbids is not None else zeros(self.nlc, 'int')
self.npb = unique(self.pbids).size
if not issubdtype(self.pbids.dtype, integer):
raise ValueError(f"The passband indices must be given as integers instead of {self.pbids.dtype}.")
if self.pbids.size != self.nlc:
raise ValueError(f"Passband index array size ({self.pbids.size}) should equal to the number of ligt curves ({self.nlc}).")
if not (self.pbids.max() == (self.npb-1) and self.pbids.min() == 0):
raise ValueError(f"Passband indices (`pbids`) for {self.npb} unique passbands should be given as integers between 0 and {self.npb - 1}.")
# Epoch indices
# -------------
self.epids = asarray(epids) if epids is not None else zeros(self.nlc, 'int')
# Supersampling
# -------------
# A number of samples and the exposure time for each light curve.
self.nsamples = atleast_1d(nsamples) if nsamples is not None else ones(self.nlc, 'int')
# Float, not int: the exposure times are a duration in days, and an integer array would
# silently truncate an exposure time assigned into it later.
self.exptimes = atleast_1d(exptimes) if exptimes is not None else zeros(self.nlc, 'float64')
# A single value applies to every light curve, as documented. Broadcasting it here rather
# than leaving a length-one array is what makes that true: the models index these per light
# curve, and a length-one array is read past its end for every light curve but the first,
# which in compiled code is a garbage sample count rather than an error.
if self.nsamples.size == 1:
self.nsamples = full(self.nlc, self.nsamples[0])
if self.exptimes.size == 1:
self.exptimes = full(self.nlc, self.exptimes[0])
if self.nsamples.size != self.nlc:
raise ValueError(f"The number of sample counts ({self.nsamples.size}) should be one or "
f"equal to the number of light curves ({self.nlc}).")
if self.exptimes.size != self.nlc:
raise ValueError(f"The number of exposure times ({self.exptimes.size}) should be one or "
f"equal to the number of light curves ({self.nlc}).")
return True
def __call__(self, *nargs, **kwargs):
raise NotImplementedError
# Default evaluation methods
# --------------------------
def evaluate(self, k: Union[float, ndarray], ldc: ndarray, t0: Union[float, ndarray], p: Union[float, ndarray],
a: Union[float, ndarray], i: Union[float, ndarray], e: Union[float, ndarray] = None, w: Union[float, ndarray] = None,
copy: bool = True) -> ndarray:
raise NotImplementedError