Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

update the timeouts and repair docs #141

Merged
merged 5 commits into from
Nov 1, 2023
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension


Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
4 changes: 2 additions & 2 deletions config/auxtel.ini
Original file line number Diff line number Diff line change
Expand Up @@ -16,9 +16,9 @@ SPECTRACTOR_DECONVOLUTION_FFM = True
# value of sigma clip parameter for the spectractor deconvolution process PSF2D and FFM
SPECTRACTOR_DECONVOLUTION_SIGMA_CLIP = 100
# maximum time per gradient descent iteration before TimeoutError in seconds
SPECTRACTOR_FIT_TIMEOUT_PER_ITER = 600
SPECTRACTOR_FIT_TIMEOUT_PER_ITER = 1200
# maximum time per gradient descent before TimeoutError in seconds
SPECTRACTOR_FIT_TIMEOUT = 3600
SPECTRACTOR_FIT_TIMEOUT = 7200

[instrument]
# instrument name
Expand Down
3 changes: 2 additions & 1 deletion docs/requirements.txt
Original file line number Diff line number Diff line change
Expand Up @@ -7,4 +7,5 @@ coloredlogs
recommonmark
mistune==2.0.3
docutils<1.0,>=0.16
sphinx-mdinclude
sphinx-mdinclude
sphinx-rtd-theme
2 changes: 1 addition & 1 deletion setup.py
Original file line number Diff line number Diff line change
Expand Up @@ -8,7 +8,7 @@
if os.getenv('CONDA_PREFIX') is None:
reqs = open('requirements.txt', 'r').read().strip().splitlines()

if os.getenv('READTHEDOCS'):
if os.getenv('READTHEDOCS') and 'mpi4py' in reqs:
reqs.remove('mpi4py')

with open('README.md') as file:
Expand Down
18 changes: 11 additions & 7 deletions spectractor/extractor/chromaticpsf.py
Original file line number Diff line number Diff line change
Expand Up @@ -2406,32 +2406,36 @@ def amplitude_derivatives(self):

def jacobian(self, params, epsilon, model_input=None):
r"""Generic function to compute the Jacobian matrix of a model, linear parameters being fixed (see Notes),
with analytical or numerical derivatives. Analytical derivatives are performed if self.analytical is True.
Let's write :math:`\theta` the non-linear model parameters. If the model is written as:
with analytical or numerical derivatives. Analytical derivatives are performed if `self.analytical` is True.
Let's write :math:`\theta` the non-linear model parameters. If the model is written as:

.. math::

\mathbf{I} = \mathbf{M}(\theta) \cdot \hat{\mathbf{A}}(\theta),

this jacobian function returns:

.. math::

\frac{\partial \mathbf{I}}{\partial \theta} = \frac{\partial \mathbf{M}}{\partial \theta} \cdot \hat{\mathbf{A}}.

Notes
-----
The gradient descent is performed on the non-linear parameters :math:`\theta` (PSF shape and position). Linear parameters :math:`\mathbf{A}`(amplitudes) are computed on the fly.
The gradient descent is performed on the non-linear parameters :math:`\theta` (PSF shape and position). Linear parameters :math:`\mathbf{A}` (amplitudes) are computed on the fly.
Therefore, :math:`\chi^2` is a function of :math:`\theta` only

.. math ::

\chi^2(\theta) = \chi'^2(\theta, \hat{\mathbf{A}}

whose partial derivatives on :math:`\theta` for gradient descent are
whose partial derivatives on :math:`\theta` for gradient descent are:

.. math ::

\frac{\partial \chi^2}{\partial \theta} = \left.\left(\frac{\partial \chi'^2}{\partial \theta} + \frac{\partial \chi'^2}{\partial \mathbf{A}} \frac{\partial \mathbf{A}}{\partial \theta}\right)\right\vert_{\mathbf{A} = \hat \mathbf{A}}
\frac{\partial \chi^2}{\partial \theta} = \left.\left(\frac{\partial \chi'^2}{\partial \theta} + \frac{\partial \chi'^2}{\partial \mathbf{A}} \frac{\partial \mathbf{A}}{\partial \theta}\right)\right\vert_{\mathbf{A} = \hat{\mathbf{A}}}

By definition, :math:`\left.\partial \chi'^2/\partial \mathbf{A}\right\vert_{\mathbf{A} = \hat \mathbf{A}}=0` then :math:`\chi^2` partial derivatives must be performed with fixed :math:`\mathbf{A} = \hat \mathbf{A}`
for gradient descent. `self.amplitude_priors_method` is temporarily switched to "keep" in `self.jacobian()` to use previously computed :math:`\hat \mathbf{A}` solution.
By definition, :math:`\left.\partial \chi'^2/\partial \mathbf{A}\right\vert_{\mathbf{A} = \hat{\mathbf{A}}}=0` then :math:`\chi^2` partial derivatives must be performed with fixed :math:`\mathbf{A} = \hat{\mathbf{A}}`
for gradient descent. `self.amplitude_priors_method` is temporarily switched to "keep" in `self.jacobian()` to use previously computed :math:`\hat{\mathbf{A}}` solution.


Parameters
Expand Down
48 changes: 23 additions & 25 deletions spectractor/extractor/psf.py
Original file line number Diff line number Diff line change
Expand Up @@ -18,10 +18,10 @@ def evaluate_moffat1d_normalisation(gamma, alpha):

.. math ::

A = \frac{\Gamma(alpha)}{\gamma \sqrt{\pi} \Gamma(alpha -1/2)}
A = \frac{\Gamma(\alpha)}{\gamma \sqrt{\pi} \Gamma(\alpha -1/2)}
\quad\text{with}, \alpha > 1/2

Note that this function is defined only for :math:`alpha > 1/2`.
Note that this function is defined only for :math:`\alpha > 1/2`.

Parameters
----------
Expand All @@ -48,10 +48,10 @@ def evaluate_moffat1d_normalisation_dalpha(norm, alpha):

.. math ::

A = \frac{\Gamma(alpha)}{\gamma \sqrt{\pi} \Gamma(alpha -1/2)}
A = \frac{\Gamma(\alpha)}{\gamma \sqrt{\pi} \Gamma(\alpha -1/2)}
\quad\text{with}, \alpha > 1/2

Note that this function is defined only for :math:`alpha > 1/2`.
Note that this function is defined only for :math:`\alpha > 1/2`.

Parameters
----------
Expand Down Expand Up @@ -83,8 +83,8 @@ def evaluate_moffat1d(y, amplitude, y_c, gamma, alpha, norm): # pragma: no cove
f(y) \propto \frac{A}{\left[ 1 +\left(\frac{y-y_c}{\gamma}\right)^2 \right]^\alpha}
\quad\text{with}, \alpha > 1/2

Note that this function is defined only for :math:`alpha > 1/2`. The normalisation factor
:math:`\frac{\Gamma(alpha)}{\gamma \sqrt{\pi} \Gamma(alpha -1/2)}` is not included as special functions
Note that this function is defined only for :math:`\alpha > 1/2`. The normalisation factor
:math:`\frac{\Gamma(\alpha)}{\gamma \sqrt{\pi} \Gamma(\alpha -1/2)}` is not included as special functions
are not supported by numba library.

Parameters
Expand All @@ -100,7 +100,7 @@ def evaluate_moffat1d(y, amplitude, y_c, gamma, alpha, norm): # pragma: no cove
alpha: float
Exponent :math:`\alpha` of the Moffat function.
norm: float
Normalisation :math:`\frac{\Gamma(alpha)}{\gamma \sqrt{\pi} \Gamma(alpha -1/2)}`.
Normalisation :math:`\frac{\Gamma(\alpha)}{\gamma \sqrt{\pi} \Gamma(\alpha -1/2)}`.

Returns
-------
Expand Down Expand Up @@ -156,10 +156,10 @@ def evaluate_moffat1d_jacobian(y, amplitude, y_c, gamma, alpha, norm, dnormda, f

.. math ::

f(y) \propto \frac{A}{\left[ 1 +\left(\frac{y-y_c}{\gamma}\right)^2 \right]^\alpha}\times \frac{\Gamma(alpha)}{\gamma \sqrt{\pi} \Gamma(alpha -1/2)}
f(y) \propto \frac{A}{\left[ 1 +\left(\frac{y-y_c}{\gamma}\right)^2 \right]^\alpha}\times \frac{\Gamma(\alpha)}{\gamma \sqrt{\pi} \Gamma(\alpha -1/2)}
\quad\text{with}, \alpha > 1/2

Note that this function is defined only for :math:`alpha > 1/2`.
Note that this function is defined only for :math:`\alpha > 1/2`.

Parameters
----------
Expand All @@ -174,7 +174,7 @@ def evaluate_moffat1d_jacobian(y, amplitude, y_c, gamma, alpha, norm, dnormda, f
alpha: float
Exponent :math:`\alpha` of the Moffat function.
norm: float
Normalisation :math:`\frac{\Gamma(alpha)}{\gamma \sqrt{\pi} \Gamma(alpha -1/2)}`.
Normalisation :math:`\frac{\Gamma(\alpha)}{\gamma \sqrt{\pi} \Gamma(\alpha -1/2)}`.
dnormda: float
Derivatives of the normalisation with respect to alpha.
fixed: array_like
Expand Down Expand Up @@ -241,8 +241,8 @@ def evaluate_moffatgauss1d(y, amplitude, y_c, gamma, alpha, eta_gauss, sigma, no
\frac{1}{\left[ 1 +\left(\frac{y-y_c}{\gamma}\right)^2 \right]^\alpha}+ \eta e^{-(y-y_c)^2/(2\sigma^2)}\right\rbrace
\quad\text{ and } \quad \eta < 0, \alpha > 1/2

Note that this function is defined only for :math:`alpha > 1/2`. The normalisation factor for the Moffat+Gauss
:math:`\frac{\Gamma(alpha)}{\gamma \sqrt{\pi} \Gamma(alpha -1/2)} + \eta \sqrt{2\pi} \sigma` is not included as special functions
Note that this function is defined only for :math:`\alpha > 1/2`. The normalisation factor for the Moffat+Gauss
:math:`\frac{\Gamma(\alpha)}{\gamma \sqrt{\pi} \Gamma(\alpha -1/2)} + \eta \sqrt{2\pi} \sigma` is not included as special functions
are not supproted by the numba library.

Parameters
Expand All @@ -262,7 +262,7 @@ def evaluate_moffatgauss1d(y, amplitude, y_c, gamma, alpha, eta_gauss, sigma, no
sigma: float
Width :math:`\sigma` of the Gaussian function.
norm_moffat: float
Normalisation :math:`\frac{\Gamma(alpha)}{\gamma \sqrt{\pi} \Gamma(alpha -1/2)}`.
Normalisation :math:`\frac{\Gamma(\alpha)}{\gamma \sqrt{\pi} \Gamma(\alpha -1/2)}`.

Returns
-------
Expand Down Expand Up @@ -326,12 +326,12 @@ def evaluate_moffatgauss1d_jacobian(y, amplitude, y_c, gamma, alpha, eta_gauss,
.. math ::

f(y) \propto A \left\lbrace
\frac{1}{\left[ 1 +\left(\frac{y-y_c}{\gamma}\right)^2 \right]^\alpha} \times \frac{\Gamma(alpha)}{\gamma \sqrt{\pi} \Gamma(alpha -1/2)}
\frac{1}{\left[ 1 +\left(\frac{y-y_c}{\gamma}\right)^2 \right]^\alpha} \times \frac{\Gamma(\alpha)}{\gamma \sqrt{\pi} \Gamma(\alpha -1/2)}
- \eta e^{-(y-y_c)^2/(2\sigma^2)}\right\rbrace
\quad\text{ and } \quad \eta < 0, \alpha > 1/2

Note that this function is defined only for :math:`alpha > 1/2`. The normalisation factor for the Moffat
:math:`\frac{\Gamma(alpha)}{\gamma \sqrt{\pi} \Gamma(alpha -1/2)} + \eta \sqrt{2\pi} \sigma` is not included as special functions
Note that this function is defined only for :math:`\alpha > 1/2`. The normalisation factor for the Moffat
:math:`\frac{\Gamma(\alpha)}{\gamma \sqrt{\pi} \Gamma(\alpha -1/2)} + \eta \sqrt{2\pi} \sigma` is not included as special functions
are not supproted by the numba library.

Parameters
Expand All @@ -351,7 +351,7 @@ def evaluate_moffatgauss1d_jacobian(y, amplitude, y_c, gamma, alpha, eta_gauss,
sigma: float
Width :math:`\sigma` of the Gaussian function.
norm_moffat: float
Normalisation :math:`\frac{\Gamma(alpha)}{\gamma \sqrt{\pi} \Gamma(alpha -1/2)}`.
Normalisation :math:`\frac{\Gamma(\alpha)}{\gamma \sqrt{\pi} \Gamma(\alpha -1/2)}`.
dnormda: float
Derivatives of the normalisation with respect to alpha.
fixed: array_like
Expand Down Expand Up @@ -430,7 +430,7 @@ def evaluate_moffat2d(x, y, amplitude, x_c, y_c, gamma, alpha): # pragma: no co
\quad\text{with}\quad
\int_{-\infty}^{\infty}\int_{-\infty}^{\infty}f(x, y) \mathrm{d}x \mathrm{d}y = A

Note that this function is defined only for :math:`alpha > 1`.
Note that this function is defined only for :math:`\alpha > 1`.

Note that the normalisation of a 2D Moffat function is analytical so it is not expected that
the sum of the output array is equal to :math:`A`, but lower.
Expand Down Expand Up @@ -589,7 +589,7 @@ def evaluate_moffatgauss2d(x, y, amplitude, x_c, y_c, gamma, alpha, eta_gauss, s
\int_{-\infty}^{\infty}\int_{-\infty}^{\infty}f(x, y) \mathrm{d}x \mathrm{d}y = A
\quad\text{and} \quad \eta < 0

Note that this function is defined only for :math:`alpha > 1`.
Note that this function is defined only for :math:`\alpha > 1`.

Parameters
----------
Expand Down Expand Up @@ -672,8 +672,7 @@ def evaluate_gauss1d(y, amplitude, y_c, sigma): # pragma: no cover

.. math ::

f(x, y) = \frac{A}{\sigma \sqrt{2 \pi}\left\lbrace e^{-\left[ \left(x-x_c\right)^2\right]/(2 \sigma^2)}
\right\rbrace
f(y) = \frac{A}{\sigma \sqrt{2 \pi}\left\lbrace e^{-\left[ \left(y-y_c\right)^2\right]/(2 \sigma^2)}\right\rbrace

.. math ::
\quad\text{with}\quad
Expand Down Expand Up @@ -743,8 +742,7 @@ def evaluate_gauss1d_jacobian(y, amplitude, y_c, sigma, fixed): # pragma: no co

.. math ::

f(x, y) = \frac{A}{\sigma \sqrt{2 \pi}\left\lbrace e^{-\left[ \left(x-x_c\right)^2\right]/(2 \sigma^2)}
\right\rbrace
f(y) = \frac{A}{\sigma \sqrt{2 \pi}\left\lbrace e^{-\left[ \left(y-y_c\right)^2\right]/(2 \sigma^2)}\right\rbrace

.. math ::
\quad\text{with}\quad
Expand Down Expand Up @@ -1206,7 +1204,7 @@ def evaluate(self, pixels, values=None):

.. math::

f(y) \propto \frac{A \Gamma(alpha)}{\gamma \sqrt{\pi} \Gamma(alpha -1/2)},
f(y) \propto \frac{A \Gamma(\alpha)}{\gamma \sqrt{\pi} \Gamma(\alpha -1/2)},
\quad \int_{y_{\text{min}}}^{y_{\text{max}}} f(y) \mathrm{d}y = A

Parameters
Expand Down Expand Up @@ -1576,7 +1574,7 @@ def evaluate(self, pixels, values=None):

.. math::

f(y) \propto \frac{A}{ \frac{\Gamma(alpha)}{\gamma \sqrt{\pi} \Gamma(alpha -1/2)}+\eta\sqrt{2\pi}\sigma},
f(y) \propto \frac{A}{ \frac{\Gamma(\alpha)}{\gamma \sqrt{\pi} \Gamma(\alpha -1/2)}+\eta\sqrt{2\pi}\sigma},
\quad \int_{y_{\text{min}}}^{y_{\text{max}}} f(y) \mathrm{d}y = A

Parameters
Expand Down
3 changes: 3 additions & 0 deletions spectractor/extractor/spectroscopy.py
Original file line number Diff line number Diff line change
Expand Up @@ -430,6 +430,7 @@ def build_detected_line_table(self, amplitude_units="", calibration_only=False):
--------

Creation of a mock spectrum with emission and absorption lines

>>> from spectractor.extractor.spectrum import Spectrum, detect_lines
>>> lambdas = np.arange(300,1000,1)
>>> spectrum = 1e4*np.exp(-((lambdas-600)/200)**2)
Expand All @@ -444,12 +445,14 @@ def build_detected_line_table(self, amplitude_units="", calibration_only=False):
>>> fwhm_func = interp1d(lambdas, 0.01 * lambdas)

Detect the lines

>>> lines = Lines([HALPHA, HBETA, O2_1], hydrogen_only=True,
... atmospheric_lines=True, redshift=0, emission_spectrum=True)
>>> global_chisq = detect_lines(lines, lambdas, spectrum, spectrum_err, fwhm_func=fwhm_func)
>>> assert(global_chisq < 1)

Print the result

>>> spec.lines = lines
>>> t = lines.build_detected_line_table()
"""
Expand Down
6 changes: 4 additions & 2 deletions spectractor/extractor/spectrum.py
Original file line number Diff line number Diff line change
Expand Up @@ -814,7 +814,8 @@ def load_spectrum(self, input_file_name, spectrogram_file_name_override=None,
Examples
--------

# Latest Spectractor output format: do not provide a config file (parameters are loaded from file header)
Latest Spectractor output format: do not provide a config file (parameters are loaded from file header)

>>> from spectractor import parameters
>>> s = Spectrum(config="")
>>> s.load_spectrum('tests/data/reduc_20170530_134_spectrum.fits')
Expand All @@ -826,7 +827,8 @@ def load_spectrum(self, input_file_name, spectrogram_file_name_override=None,
>>> assert parameters.CCD_REBIN == s.header["REBIN"]
>>> assert s.parallactic_angle == s.header["PARANGLE"]

# Spectractor output format older than version <=2.3: must give the config file
Spectractor output format older than version <=2.3: must give the config file

>>> parameters.VERBOSE = False
>>> s = Spectrum(config="./config/ctio.ini")
>>> s.load_spectrum('tests/data/reduc_20170605_028_spectrum.fits')
Expand Down
Loading