Added tests for linear and quadratic tests, show agreement to better …

…than 0.1%. Also added basic tests for Legendre and Chebyshev polynomial fit. Potentially fixed RTD broken API generation.

.gitignore

@@ -133,3 +133,4 @@ iframe_figures/
 # test output
 test/*.json
 test/*.npy
+test/*.png

README.md

@@ -40,3 +40,6 @@ The [quickstart guide](https://rascal.readthedocs.io/en/latest/tutorial/quicksta
 If you are interested in contributing code to the project, thank you! For those unfamiliar with the process of contributing to an open-source project, you may want to read through Github’s own short informational section on how to submit a [contribution](https://opensource.guide/how-to-contribute/#how-to-submit-a-contribution) or send me a message.
 
 Style -- as long as it passes flake8.
+
+## Disclaimer
+We duplicate some of the relevant metadata, but we do not process the raw metadata. Some of the metadata this software creates contain full path to the files in your system, which most likely includes a user name on your machine. Please be advised it is your responsibility to be compliant with the privacy law(s) that you are oblidged to follow, and it is your responsibility to remove any metadata that may reveal personal information and/or provide information that can reveal any computing vulunerability.

docs/source/conf.py

@@ -32,8 +32,7 @@
 # ones.
 extensions = [
     'sphinx.ext.todo', 'sphinx.ext.viewcode', 'sphinx.ext.autodoc',
-    'sphinx.ext.autosectionlabel', 'sphinx.ext.coverage',
-    'sphinx.ext.mathjax'
+    'sphinx.ext.autosectionlabel', 'sphinx.ext.coverage', 'sphinx.ext.mathjax'
 ]
 
 # Add any paths that contain templates here, relative to this directory.
@@ -44,6 +43,25 @@
 # This pattern also affects html_static_path and html_extra_path.
 exclude_patterns = []
 
+# The suffix of source filenames.
+source_suffix = '.rst'
+
+# The encoding of source files.
+source_encoding = 'utf-8-sig'
+
+# Produce API reference automatically for every public and private methods
+apidoc_args = [
+    '--force',  # Overwrite existing files
+    '--no-toc',  # Don't create a table of contents file
+    '--output-dir=.',  # Directory to place all output
+]
+
+autodoc_mock_imports = [
+    'matplotlib', 'plotly', 'pynverse', 'tdqm', 'numpy', 'scipy'
+]
+
+autoclass_content = 'both'
+
 mathjax_path = "https://cdn.mathjax.org/mathjax/latest/MathJax.js?config=TeX-AMS-MML_HTMLorMML"
 
 # -- Options for HTML output -------------------------------------------------

rascal/calibrator.py

@@ -344,9 +344,10 @@ def __init__(self, peaks, spectrum=None):
         self.hough_lines = None
         self.hough_points = None
         self.ht = HoughTransform()
+
         self.set_calibrator_properties()
         self.set_hough_properties()
-        self.set_ransac_properties()
+        self.set_ransac_properties(sample_size=0)
 
     def _import_matplotlib(self):
         '''
@@ -1592,14 +1593,12 @@ def set_known_pairs(self, pix=(), wave=()):
             'and wave has size {}.'.format(pix.size, wave.size))
 
         if not all(
-                isinstance(p, (float, int)) & (not np.isnan(p))
-                for p in pix):
+                isinstance(p, (float, int)) & (not np.isnan(p)) for p in pix):
 
             raise ValueError("All pix elements have to be numeric.")
 
         if not all(
-                isinstance(w, (float, int)) & (not np.isnan(w))
-                for w in wave):
+                isinstance(w, (float, int)) & (not np.isnan(w)) for w in wave):
 
             raise ValueError("All wave elements have to be numeric.")
 

requirements.txt

@@ -3,4 +3,4 @@ scipy
 pynverse
 matplotlib
 tqdm
-pytest
+astropy

test/test_lt_sprat_manual_atlas.py

@@ -1,14 +1,11 @@
 import numpy as np
 from astropy.io import fits
 from scipy.signal import find_peaks
-from matplotlib import pyplot as plt
 import os
 
 from rascal.calibrator import Calibrator
 from rascal import util
 
-plt.ion()
-
 
 def test_sprat_manual_atlas():
 
@@ -37,18 +34,14 @@ def test_sprat_manual_atlas():
     c.use_plotly()
     assert c.which_plotting_library() == 'plotly'
 
-    if os.name != 'nt':
-        c.plot_arc()
+    c.plot_arc(display=False)
 
     c.set_hough_properties(num_slopes=5000,
                            range_tolerance=500.,
                            xbins=100,
                            ybins=100,
                            min_wavelength=3500.,
                            max_wavelength=8000.)
-    c.set_ransac_properties(sample_size=5,
-                            top_n_candidate=5,
-                            filter_close=True)
     # blend: 4829.71, 4844.33
     # blend: 5566.62, 5581.88
     # blend: 6261.212, 6265.302
@@ -71,63 +64,72 @@ def test_sprat_manual_atlas():
                      temperature=temperature,
                      relative_humidity=relative_humidity)
 
+    c.set_ransac_properties(sample_size=5,
+                            top_n_candidate=5,
+                            filter_close=True)
+
     c.do_hough_transform()
 
     # Run the wavelength calibration
     best_p, rms, residual, peak_utilisation = c.fit(max_tries=250)
 
-    if os.name != 'nt':
-
-        # Plot the solution
-        c.plot_fit(best_p,
-                   spectrum,
-                   plot_atlas=True,
-                   log_spectrum=False,
-                   tolerance=5.)
+    # Plot the solution
+    c.plot_fit(best_p,
+               spectrum,
+               plot_atlas=True,
+               log_spectrum=False,
+               display=False,
+               tolerance=5.)
 
     fit_coeff_new, peak_matched, atlas_matched, residual,\
         peak_utilisation = c.match_peaks(best_p)
 
-    if os.name != 'nt':
-
-        c.plot_fit(fit_coeff_new,
-                   spectrum,
-                   plot_atlas=True,
-                   log_spectrum=False,
-                   tolerance=5.)
+    c.plot_fit(fit_coeff_new,
+               spectrum,
+               plot_atlas=True,
+               log_spectrum=False,
+               display=False,
+               tolerance=5.)
 
     # Show the parameter space for searching possible solution
-    if os.name != 'nt':
-
-        c.plot_search_space()
+    c.plot_search_space(display=False)
 
     print("Stdev error: {} A".format(residual.std()))
     print("Peaks utilisation rate: {}%".format(peak_utilisation * 100))
 
     c.use_matplotlib()
     assert c.which_plotting_library() == 'matplotlib'
 
-    c.plot_arc(display=False)
+    c.plot_arc(display=False, savefig=True, filename='test/test_lt_sprat_arc')
+    c.plot_arc(log_spectrum=True,
+               display=False,
+               savefig=True,
+               filename='test/test_lt_sprat_arc_log')
 
     # Plot the solution
     c.plot_fit(best_p,
                spectrum,
                plot_atlas=True,
                log_spectrum=False,
                tolerance=5.,
-               display=False)
+               display=False,
+               savefig=True,
+               filename='test/test_lt_sprat_fit')
 
     # Show the parameter space for searching possible solution
     c.plot_search_space(display=False)
 
-    c.plot_arc()
+    c.plot_arc(display=False)
 
     # Plot the solution
     c.plot_fit(best_p,
                spectrum,
                plot_atlas=True,
-               log_spectrum=False,
-               tolerance=5.)
+               log_spectrum=True,
+               tolerance=5.,
+               savefig=True,
+               display=False,
+               filename='test/test_lt_sprat_fit_log')
 
     # Show the parameter space for searching possible solution
-    c.plot_search_space()
+    c.plot_search_space(display=False)

test/test_polynomial_fit.py

@@ -0,0 +1,116 @@
+import numpy as np
+from rascal.calibrator import Calibrator
+
+peaks = np.sort(np.random.random(31) * 1000.)
+
+# Line list
+wavelengths_linear = 3000. + 5. * peaks
+wavelengths_quadratic = 3000. + 4 * peaks + 1.0e-3 * peaks**2.
+
+elements_linear = ['Linear'] * len(wavelengths_linear)
+elements_quadratic = ['Quadratic'] * len(wavelengths_quadratic)
+
+
+def test_linear_fit():
+
+    # Initialise the calibrator
+    c = Calibrator(peaks)
+    c.set_calibrator_properties(num_pix=1000)
+    c.set_hough_properties(num_slopes=1000,
+                           range_tolerance=500.,
+                           xbins=200,
+                           ybins=200,
+                           min_wavelength=3000.,
+                           max_wavelength=8000.)
+    c.add_user_atlas(elements=elements_linear, wavelengths=wavelengths_linear)
+    c.set_ransac_properties()
+
+    # Run the wavelength calibration
+    best_p, rms, residual, peak_utilisation = c.fit(max_tries=500, fit_deg=1)
+    # Refine solution
+    best_p, x_fit, y_fit, residual, peak_utilisation = c.match_peaks(
+        best_p, refine=False, robust_refit=True)
+
+    assert (best_p[1] > 5. * 0.999) & (best_p[1] < 5. * 1.001)
+    assert (best_p[0] > 3000. * 0.999) & (best_p[0] < 3000. * 1.001)
+
+
+def test_quadratic_fit():
+
+    # Initialise the calibrator
+    c = Calibrator(peaks)
+    c.set_calibrator_properties(num_pix=1000)
+    c.set_hough_properties(num_slopes=1000,
+                           range_tolerance=500.,
+                           xbins=100,
+                           ybins=100,
+                           min_wavelength=3000.,
+                           max_wavelength=8000.)
+    c.add_user_atlas(elements=elements_quadratic,
+                     wavelengths=wavelengths_quadratic)
+    c.set_ransac_properties()
+
+    # Run the wavelength calibration
+    best_p, rms, residual, peak_utilisation = c.fit(max_tries=1000, fit_deg=2)
+    # Refine solution
+    best_p, x_fit, y_fit, residual, peak_utilisation = c.match_peaks(
+        best_p, refine=False, robust_refit=True)
+
+    assert (best_p[2] > 1e-3 * 0.999) & (best_p[2] < 1e-3 * 1.001)
+    assert (best_p[1] > 4. * 0.999) & (best_p[1] < 4. * 1.001)
+    assert (best_p[0] > 3000. * 0.999) & (best_p[0] < 3000. * 1.001)
+
+
+def test_quadratic_fit_legendre():
+
+    # Initialise the calibrator
+    c = Calibrator(peaks)
+    c.set_calibrator_properties(num_pix=1000)
+    c.set_hough_properties(num_slopes=500,
+                           range_tolerance=200.,
+                           xbins=100,
+                           ybins=100,
+                           min_wavelength=3000.,
+                           max_wavelength=8000.)
+    c.add_user_atlas(elements=elements_quadratic,
+                     wavelengths=wavelengths_quadratic)
+    c.set_ransac_properties(sample_size=10)
+
+    # Run the wavelength calibration
+    best_p, rms, residual, peak_utilisation = c.fit(max_tries=1000,
+                                                    fit_tolerance=5.,
+                                                    candidate_tolerance=2.,
+                                                    fit_deg=2,
+                                                    fit_type='legendre')
+
+    # Legendre 2nd order takes the form
+    assert (best_p[2] > 1e-3 * 0.999) & (best_p[2] < 1e-3 * 1.001)
+    assert (best_p[1] > 4. * 0.999) & (best_p[1] < 4. * 1.001)
+    assert (best_p[0] > 3000. * 0.999) & (best_p[0] < 3000. * 1.001)
+
+
+def test_quadratic_fit_chebyshev():
+
+    # Initialise the calibrator
+    c = Calibrator(peaks)
+    c.set_calibrator_properties(num_pix=1000)
+    c.set_hough_properties(num_slopes=500,
+                           range_tolerance=200.,
+                           xbins=100,
+                           ybins=100,
+                           min_wavelength=3000.,
+                           max_wavelength=8000.)
+    c.add_user_atlas(elements=elements_quadratic,
+                     wavelengths=wavelengths_quadratic)
+    c.set_ransac_properties(sample_size=10)
+
+    # Run the wavelength calibration
+    best_p, rms, residual, peak_utilisation = c.fit(max_tries=1000,
+                                                    fit_tolerance=5.,
+                                                    candidate_tolerance=2.,
+                                                    fit_deg=2,
+                                                    fit_type='chebyshev')
+
+    assert (best_p[2] > 1e-3 * 0.999) & (best_p[2] < 1e-3 * 1.001)
+    assert (best_p[1] > 4. * 0.999) & (best_p[1] < 4. * 1.001)
+    assert (best_p[0] > 3000. * 0.999) & (best_p[0] < 3000. * 1.001)

test/test_xe_calibration.py

@@ -17,7 +17,7 @@
     4734.15, 4807.02, 4921.48, 5028.28, 5618.88, 5823.89, 5893.29, 5934.17,
     6182.42, 6318.06, 6472.841, 6595.56, 6668.92, 6728.01, 6827.32, 6976.18,
     7119.60, 7257.9, 7393.8, 7584.68, 7642.02, 7740.31, 7802.65, 7887.40,
-    7967.34, 8057.258
+    7967.34, 8057.258, 10000.
 ]
 elements = ['Xe'] * len(wavelengths)
 
@@ -48,6 +48,12 @@ def run_sprat_calibration(fit_deg):
                            min_wavelength=3500.,
                            max_wavelength=8000.)
     c.add_user_atlas(elements=elements, wavelengths=wavelengths)
+    c.clear_atlas()
+    assert len(c.atlas) == 0
+    c.add_user_atlas(elements=elements, wavelengths=wavelengths)
+    c.remove_atlas_lines_range(9999.)
+    assert len(c.atlas) == len(wavelengths) - 1
+    c.list_atlas()
 
     # Run the wavelength calibration
     best_p, rms, residual, peak_utilisation = c.fit(max_tries=200,