From e5372b18ffdc33bc5c50bffe677fcedf385e7d4e Mon Sep 17 00:00:00 2001 From: Paul Hirst Date: Tue, 5 May 2026 12:57:01 -1000 Subject: [PATCH 01/12] Enh/embedded figures - adds a module in gempy/adlibrary to embed non-FITS files (PDFs etc) in astrodata instances and FITS files (#516) * gempy support for embedding files inside a binary table in an astrodata object / FITS file * capture arc solution PDF into astrodata object --- geminidr/core/primitives_spect.py | 4 + gempy/adlibrary/embedded_files.py | 153 +++++++++++++++++++ gempy/adlibrary/tests/test_embedded_files.py | 110 +++++++++++++ gempy/scripts/extract_files | 1 + gempy/scripts/extract_files.py | 45 ++++++ 5 files changed, 313 insertions(+) create mode 100644 gempy/adlibrary/embedded_files.py create mode 100644 gempy/adlibrary/tests/test_embedded_files.py create mode 120000 gempy/scripts/extract_files create mode 100755 gempy/scripts/extract_files.py diff --git a/geminidr/core/primitives_spect.py b/geminidr/core/primitives_spect.py index 6d4f473158..7273fee0fe 100644 --- a/geminidr/core/primitives_spect.py +++ b/geminidr/core/primitives_spect.py @@ -56,6 +56,7 @@ from gempy.library.matching import KDTreeFitter from gempy.library.nddops import NDStacker from gempy.library.spectral import Spek1D +from gempy.adlibrary.embedded_files import embed_file from gwcs.utils import CoordinateFrameError from recipe_system.utils.decorators import parameter_override, capture_provenance from recipe_system.utils.md5 import md5sum @@ -2934,6 +2935,9 @@ def determineWavelengthSolution(self, adinputs=None, **params): pdf.savefig(fig, bbox_inches='tight') plt.close() + # We could add a parameter to enable this if desired + embed_file(ad, plot_filename, contenttype='pdf') + # Timestamp and update the filename gt.mark_history(ad, primname=self.myself(), keyword=timestamp_key) diff --git a/gempy/adlibrary/embedded_files.py b/gempy/adlibrary/embedded_files.py new file mode 100644 index 0000000000..0da207a8c2 --- /dev/null +++ b/gempy/adlibrary/embedded_files.py @@ -0,0 +1,153 @@ + +from astropy.io import fits +import hashlib +""" +This hndles "embedding" files (eg PDFs) in the ad object so that they can be +stored in the FITS file. We create a binary table attribute (tablename, which +defaults to .FIGURES), which has columns for filename, Content-Type (ie mime +type), size, md5, and content, the latter of which is a varibale length binary +field which contains the actual contents of the file. + +This is useful for storing eg pdf figures generated during reduction inside +the reduced data fits file, which is especially useful when the fits file is +stored in the archive. + +This code is not especially memory efficient - multiple copies of the data +may end up in memory. It is not intended for use with large files and will +likely need refactoring if that use cases arises. +""" + +def embed_file(ad, filename, contenttype, data=None, tablename='FIGURES'): + """ + Embed a file in an astrodata instance. + :param ad: The ad instance to embed the file in + :param filename: Filename of the file to embed + :param contenttype: MIME type of the file + :param data: bytes or bytesarray of the actual binary file contents, or a + file-list object from which to read it from. + If NULL, filename will be opened and read + :param tablename: Name of the binary table attribute to store the file in. + Defaults to 'FIGURES' + :return: None + """ + + # Some convenience values are supported for contenttype: + contenttype = 'application/pdf' if contenttype.upper() == 'PDF' else contenttype + contenttype = 'image/png' if contenttype.upper() == 'PNG' else contenttype + contenttype = 'image/jpeg' if contenttype.upper() in ('JPG', 'JPEG') else contenttype + + # If tablename already exists, read the columns into arrays + thdu = getattr(ad, tablename, None) + if thdu is not None: + filenames = thdu['filename'].tolist() + contenttypes = thdu['Content-Type'].tolist() + sizes = thdu['size'].tolist() + md5s = thdu['md5'].tolist() + contents = thdu['content'].tolist() + else: + # Initialize empty column arrays + filenames = [] + contenttypes = [] + sizes = [] + md5s = [] + contents = [] + + # Append the new file to the arrays + filenames.append(filename) + contenttypes.append(contenttype) + + # If data None, read from the filename. + if data is None: + with open(filename, 'rb') as fp: + content = bytearray(fp.read()) + # If data is a file-like, read from it + elif hasattr(data, 'read') and callable(data.read): + print("file like") + content = bytearray(data.read()) + print(type(content)) + elif isinstance(data, bytes): + print("bytes") + content = bytearray(data) + elif isinstance(data, bytearray): + print("bytearray") + content = data + else: + raise RuntimeError("Unknown data type in embed_file") + + contents.append(content) + sizes.append(len(content)) + + h = hashlib.new('md5') + h.update(bytes(content)) + md5s.append(h.hexdigest()) + + # Make a new table with the new arrays + col1 = fits.Column(name='filename', format='128A', array=filenames) + col2 = fits.Column(name='Content-Type', format='32A', array=contenttypes) + col3 = fits.Column(name='size', format='J', array=sizes) + col4 = fits.Column(name='md5', format='32A', array=md5s) + col5 = fits.Column(name='content', format='PB()', + array=contents) + hdu = fits.BinTableHDU.from_columns([col1, col2, col3, col4, col5]) + setattr(ad, tablename, hdu) + +def list_files(ad, fullinfo=False, tablename='FIGURES'): + """ + List the files embedded in the astrodata instance + :param ad: the astrodata instance with files embedded + :param fullinfo: If False (default) return a list of filenames. + If True, return a list of dicts containing 'filename', 'size', 'Content-Type' and 'md5' keys + :param tablename: Name of the binary table attribute to store the file in. + :return: list or dict as above + """ + + hdu = getattr(ad, tablename) + retary = [] + + for row in hdu: + if fullinfo: + d = {} + for item in ('filename', 'size', 'Content-Type', 'md5'): + d[item] = row[item] + retary.append(d) + else: + retary.append(row['filename']) + return retary + +def extract_files(ad, todisk=True, filename=None, check=True, tablename='FIGURES'): + """ + Extract files embedded in an astrodata instance + :param ad: astrodata instance to extract from + :param todisk: If True (default) write files to disk in current directory, + using the filenames specified when they were embedded. If False, return + a bytes instance with the content of the file. + :param filename: If there are multiple files embedded and you want only + one of them, give the filename here. If there are multiple files embedded + and you do not specify a filename and you specify todisk=False, an exception + will be raised + :param check: If True (default) raise an exception if the size or md5 + does not match the content of the file. + :return: if todisk=True, returns a list of files writted to disk. + If todisk=False, returns the binary file content + """ + hdu = getattr(ad, tablename) + retary = [] + if len(hdu) > 1 and todisk is False and filename is None: + raise RuntimeError("Cannot extract multiple files not to disk") + for row in hdu: + if filename is None or row['filename'] == filename: + if check: + if len(row['content']) != row['size']: + raise RuntimeError(f"File size mismatch for embedded file {row['filename']}") + h = hashlib.new('md5') + h.update(bytes(row['content'])) + if h.hexdigest() != row['md5']: + raise RuntimeError(f"MD5 mismatch for embedded file {row['filename']}") + if todisk: + with open(row['filename'], 'wb') as fp: + fp.write(bytes(row['content'])) + retary.append(row['filename']) + else: + return bytes(row['content']) + + return retary \ No newline at end of file diff --git a/gempy/adlibrary/tests/test_embedded_files.py b/gempy/adlibrary/tests/test_embedded_files.py new file mode 100644 index 0000000000..94dc205171 --- /dev/null +++ b/gempy/adlibrary/tests/test_embedded_files.py @@ -0,0 +1,110 @@ +import os + +import astrodata +from gempy.adlibrary.embedded_files import embed_file, list_files, extract_files + +def test_embed_bytes(): + ad = astrodata.AstroData() + assert hasattr(ad, 'FIGURES') is False + + embed_file(ad, 'test.txt', contenttype='text/plain', data=b'abc') + + assert hasattr(ad, 'FIGURES') is True + assert ad.FIGURES[0]['filename'] == 'test.txt' + assert ad.FIGURES[0]['Content-Type'] == 'text/plain' + assert ad.FIGURES[0]['size'] == 3 + assert len(ad.FIGURES[0]['content']) == ad.FIGURES[0]['size'] + assert ad.FIGURES[0]['md5'] =='900150983cd24fb0d6963f7d28e17f72' + assert bytes(ad.FIGURES[0]['content']) == b'abc' + +def test_tablename(): + ad = astrodata.AstroData() + assert hasattr(ad, 'FOOBAR') is False + + embed_file(ad, 'test.txt', tablename='FOOBAR', + contenttype='text/plain', data=b'abc') + + assert hasattr(ad, 'FOOBAR') is True + assert ad.FOOBAR[0]['filename'] == 'test.txt' + +def test_embed_fp(tmp_path): + print("test embed fp") + # Make a temporary file + filename = tmp_path / 'file.dat' + with open(filename, 'wb') as fp: + fp.write(b"Hello, world") + + ad = astrodata.AstroData() + + # Open it again and embed the fp + with open(filename, 'rb') as fp: + embed_file(ad, 'file.dat', contenttype='text/plain', data=fp) + + assert hasattr(ad, 'FIGURES') is True + assert ad.FIGURES[0]['filename'] == 'file.dat' + assert ad.FIGURES[0]['Content-Type'] == 'text/plain' + assert ad.FIGURES[0]['size'] == len("Hello, world") + assert len(ad.FIGURES[0]['content']) == ad.FIGURES[0]['size'] + assert ad.FIGURES[0]['md5'] == 'bc6e6f16b8a077ef5fbc8d59d0b931b9' + assert bytes(ad.FIGURES[0]['content']) == b'Hello, world' + +def test_embed_from_file(tmp_path): + # Make a temporary file + filename = tmp_path / 'file.dat' + with open(filename, 'wb') as fp: + fp.write(b"Hello, world") + + ad = astrodata.AstroData() + + embed_file(ad, filename, contenttype='text/plain') + + assert hasattr(ad, 'FIGURES') is True + assert ad.FIGURES[0]['filename'] == str(filename) + assert ad.FIGURES[0]['filename'].endswith('file.dat') + assert ad.FIGURES[0]['Content-Type'] == 'text/plain' + assert ad.FIGURES[0]['size'] == len("Hello, world") + assert len(ad.FIGURES[0]['content']) == ad.FIGURES[0]['size'] + assert ad.FIGURES[0]['md5'] == 'bc6e6f16b8a077ef5fbc8d59d0b931b9' + assert bytes(ad.FIGURES[0]['content']) == b'Hello, world' + +def test_list_files(): + ad = astrodata.AstroData() + embed_file(ad, 'a.txt', contenttype='text/plain', data=b'abc') + embed_file(ad, 'b.txt', contenttype='text/plain', data=b'defg') + embed_file(ad, 'c.txt', contenttype='text/plain', data=b'hijkl') + + assert list_files(ad) == ['a.txt', 'b.txt', 'c.txt'] + + lod = list_files(ad, fullinfo=True) + a = lod[0] + b = lod[1] + c = lod[2] + + assert a['filename'] == 'a.txt' + assert a['size'] == 3 + assert a['md5'] == '900150983cd24fb0d6963f7d28e17f72' + assert a['Content-Type'] == 'text/plain' + assert b['filename'] == 'b.txt' + assert b['size'] == 4 + assert c['filename'] == 'c.txt' + assert c['size'] == 5 + +def test_extract_files(tmp_path): + ad = astrodata.AstroData() + embed_file(ad, 'a.txt', contenttype='text/plain', data=b'abc') + embed_file(ad, 'b.txt', contenttype='text/plain', data=b'defg') + embed_file(ad, 'c.txt', contenttype='text/plain', data=b'hijkl') + + os.chdir(tmp_path) + + extract_files(ad) + + ld = os.listdir('.') + assert 'a.txt' in ld + assert 'b.txt' in ld + assert 'c.txt' in ld + + with open('b.txt', 'rb') as fp: + b = fp.read() + + assert b == b'defg' \ No newline at end of file diff --git a/gempy/scripts/extract_files b/gempy/scripts/extract_files new file mode 120000 index 0000000000..a3341e55f6 --- /dev/null +++ b/gempy/scripts/extract_files @@ -0,0 +1 @@ +extract_files.py \ No newline at end of file diff --git a/gempy/scripts/extract_files.py b/gempy/scripts/extract_files.py new file mode 100755 index 0000000000..d5da27bf35 --- /dev/null +++ b/gempy/scripts/extract_files.py @@ -0,0 +1,45 @@ +#!/usr/bin/env python +# +# DRAGONS +# gempy.scripts +# fixheader.py +# ----------------------------------------------------------------------------- +""" +This command will extract files that have been embedded within the FITS file +using the gempy embedded_files module. This is typically used to capture PDF +or other figures generated during data reduction and store them in the +reduced data fits file. +""" +from argparse import ArgumentParser +import textwrap +import sys + +import astrodata +from gempy import __version__ + +from gempy.adlibrary.embedded_files import list_files, extract_files + +def main(args=None): + parser = ArgumentParser( + description=f"Embedded file extractor, v{__version__}", + epilog=textwrap.dedent(__doc__)) + parser.add_argument("-v", "--version", action="version", + version=f"v{__version__}") + parser.add_argument('--list', help="List embedded files, " + "as opposed to extracting them to disk", + action='store_true') + parser.add_argument('filename', help='Fits file to extract from') + + args = parser.parse_args(args) + + ad = astrodata.open(args.filename) + + if args.list: + for d in list_files(ad, fullinfo=True): + print(f"{d['filename']}: {d['Content-Type']}, {d['size']} bytes") + else: + for i in extract_files(ad): + print(f"Exracting: {i}") + +if __name__ == '__main__': + sys.exit(main()) \ No newline at end of file From 5a4f702ae6b0b263e7a60f48ffffe1e1bc6e3971 Mon Sep 17 00:00:00 2001 From: Paul Hirst Date: Tue, 5 May 2026 17:34:26 -1000 Subject: [PATCH 02/12] merge enh/checkProcessedArc42 - includes enh/auto_wavecal_help improvements and wavecal monitoring functionality --- geminidr/core/parameters_spect.py | 2 + geminidr/core/primitives_spect.py | 94 +++++++++++- geminidr/core/tests/test_spect.py | 23 ++- .../gmos/recipes/qa/recipes_ARC_LS_SPECT.py | 37 +++++ .../test_determine_wavelength_solution.py | 4 +- .../test_determine_wavelength_solution.py | 1 + gempy/library/fitting.py | 29 ++++ gempy/library/wavecal.py | 140 +++++++++++------- 8 files changed, 265 insertions(+), 65 deletions(-) diff --git a/geminidr/core/parameters_spect.py b/geminidr/core/parameters_spect.py index a72897af0e..c1d7e86748 100644 --- a/geminidr/core/parameters_spect.py +++ b/geminidr/core/parameters_spect.py @@ -579,6 +579,8 @@ class maskBeyondSlitConfig(config.Config): "Minimum fraction of pixel that must be illuminated to not be masked", float, 0.9, min=0., max=1., inclusiveMax=True) +class monitorWavelengthSolutionConfig(config.Config): + suffix = config.Field("Filename suffix", str, "_wavelengthSolutionMonitoring", optional=True) class normalizeFlatConfig(config.core_1Dfitting_config): suffix = config.Field("Filename suffix", str, "_normalized", optional=True) diff --git a/geminidr/core/primitives_spect.py b/geminidr/core/primitives_spect.py index 7273fee0fe..8e9dc30aef 100644 --- a/geminidr/core/primitives_spect.py +++ b/geminidr/core/primitives_spect.py @@ -16,6 +16,7 @@ from importlib import import_module import matplotlib +import numpy import numpy as np from astropy import units as u from astropy.io.ascii.core import InconsistentTableError @@ -2909,22 +2910,44 @@ def determineWavelengthSolution(self, adinputs=None, **params): fit1d.image = image wavecal.update_wcs_with_solution(ext, fit1d, other, config) else: + bad_solutions = [] + single_line_solutions = [] for ext in ad: - if len(ad) > 1: - log.stdinfo(f"Determining solution for extension {ext.id}") - - input_data, fit1d, acceptable_fit = wavecal.get_automated_fit( + input_data, fit1d = wavecal.get_automated_fit( ext, uiparams, p=self, linelist=linelist, bad_bits=DQ.not_signal, absorption=absorption) wavecal.update_wcs_with_solution(ext, fit1d, input_data, config) - if not acceptable_fit: - log.warning("No acceptable wavelength solution found") - else: + if fit1d.image.size: figures.append(wavecal.create_pdf_plot( input_data, fit1d.points[~fit1d.mask], fit1d.image[~fit1d.mask], title=f"{ad.filename}:{ext.id}", absorption=absorption)) + if fit1d.image.size == 1: + single_line_solutions.append(ext.id) + else: # no line matches so not an acceptable solution + bad_solutions.append(ext.id) + + if bad_solutions: + msg = f"{ad.filename}: failed to find an acceptable wavelength solution" + if len(ad) > 1: + if len(ad) > len(bad_solutions): + msg += "in extensions " + ", ".join(str(i) for i in bad_solutions) + else: + msg += " in any extensions" + if self.mode == 'sq' and len(ad) == len(bad_solutions): + raise RuntimeError(msg) + log.warning(msg) + + if single_line_solutions: + msg = f"{ad.filename}: a single-line solution has been applied" + if len(ad) > 1: + if len(ad) > len(single_line_solutions): + msg += "in extensions " + ", ".join(str(i) for i in single_line_solutions) + else: + msg += " in all extensions" + log.warning(msg) + ad.update_filename(suffix=sfx, strip=True) if figures: @@ -3995,6 +4018,63 @@ def maskBeyondSlit(self, adinputs=None, **params): return adinputs + def monitorWavelengthSolution(self, adinputs=None, **params): + """ + This captures some info from the wavelength solution for the GOA + instrument monitoring system into various MWS_ headers for the instrument + monitoring system to pick up. + + Values Captured include the polynomial coefficients and RMS from the + .WAVECAL table, the central wavelength from both the WCS and central_wavlength() + descriptor and the difference between them. + :return: + """ + log = self.log + log.debug(gt.log_message("primitive", self.myself(), "starting")) + + tocapture = ('c0', 'c1', 'c2', 'c3', 'rms', 'inv_rms', 'fwidth') + + for ad in adinputs: + for ext in ad: + wavecal = getattr(ext, 'WAVECAL', None) + if wavecal is None: + continue + + for row in wavecal: + # Capture values from WAVEVAL table into monitoring headers + if row['name'] in tocapture: + keyword = 'MWS_' + row['name'].upper()[:4] + if row['name'] == 'inv_rms': + keyword = 'MWS_IRMS' + ext.hdr[keyword] = (row['coefficients'], + f'WAVECAL {row['name']} coefficient') + + # Get the model and evaluate it at the mean of its domain + model = am.get_named_submodel(ext.wcs.forward_transform, "WAVE") + if model: + # Central Wavelength check + center = np.mean(model.domain) + model_central_wlen = model(center) + # We want actual_central_wavelength rather than central_wavelength + header_central_wlen = ext.actual_central_wavelength(asNanometers=True) + + ext.hdr['MWS_DCWL'] = (model_central_wlen - header_central_wlen, + 'Delta between model and header central_wavelength') + + # Dispersion check + model_dispersion = model(center+0.5) - model(center-0.5) + header_dispersion = ext.dispersion(asNanometers=True) + ext.hdr['MWS_DDIS'] = (model_dispersion - header_dispersion, + 'Delta between model and header dispersion') + + + # Count the number of arc lines in the wavecal table + ext.hdr['MWS_NUML'] = (numpy.count_nonzero(wavecal['peaks']), + 'Number of non-zero peaks in the WAVECAL table') + + ad.update_filename(suffix=params["suffix"], strip=True) + + return adinputs def normalizeFlat(self, adinputs=None, **params): """ diff --git a/geminidr/core/tests/test_spect.py b/geminidr/core/tests/test_spect.py index 53e1278a26..295cee3d71 100644 --- a/geminidr/core/tests/test_spect.py +++ b/geminidr/core/tests/test_spect.py @@ -44,11 +44,12 @@ from specutils.utils.wcs_utils import air_to_vac import astrodata, gemini_instruments -from astrodata.testing import ad_compare, assert_most_close +from astrodata.testing import ad_compare, assert_most_close, download_from_archive from gempy.library import astromodels as am from gempy.library.config.config import FieldValidationError from geminidr.core import primitives_spect from geminidr.f2.primitives_f2_longslit import F2Longslit +from geminidr.gmos.primitives_gmos_longslit import GMOSLongslit from geminidr.niri.primitives_niri_image import NIRIImage from geminidr.niri.primitives_niri_longslit import NIRILongslit from geminidr.gnirs import primitives_gnirs_longslit @@ -574,6 +575,26 @@ def test_adjust_wavelength_zero_point_controlled(filename, center, shift, np.testing.assert_allclose(waves, new_waves, atol=0.1*dw) +@pytest.mark.preprocessed_data +def test_determine_wavelength_solution_exits_with_no_solution_in_sq(path_to_inputs): + """Primitive should crash in SQ mode without a solution""" + ad = astrodata.open(os.path.join(path_to_inputs, "S20260331S0149_mosaic.fits")) + p = GMOSLongslit([ad]) + p.mode = "sq" # it's the default, but just to be sure + with pytest.raises(RuntimeError): + p.determineWavelengthSolution() + + +@pytest.mark.preprocessed_data +def test_determine_wavelength_solution_continues_with_no_solution_in_qa(path_to_inputs, caplog): + """Primitive should crash in SQ mode without a solution""" + ad = astrodata.open(os.path.join(path_to_inputs, "S20260331S0149_mosaic.fits")) + p = GMOSLongslit([ad]) + p.mode = "qa" + p.determineWavelengthSolution() + assert any(record.levelname == 'WARNING' for record in caplog.records) + + @pytest.mark.preprocessed_data @pytest.mark.parametrize('in_file,instrument', [# GNIRS 111/mm LongBlue, off right edge of detector diff --git a/geminidr/gmos/recipes/qa/recipes_ARC_LS_SPECT.py b/geminidr/gmos/recipes/qa/recipes_ARC_LS_SPECT.py index 74f2893d97..6493547e0c 100644 --- a/geminidr/gmos/recipes/qa/recipes_ARC_LS_SPECT.py +++ b/geminidr/gmos/recipes/qa/recipes_ARC_LS_SPECT.py @@ -19,5 +19,42 @@ def makeProcessedArc(p): p.storeProcessedArc() p.writeOutputs() +def checkProcessedArc(p): + """ + Extracts some values from the .WAVECAL extension for instrument monitoring + :param p: + :return: + """ + + p.monitorWavelengthSolution() + p.writeOutputs() + +def checkArc(p): + """ + Reduce an arc, extract some values from the .WAVECAL extension for + instrument monitoring + :param p: + :return: + """ + # In QA mode, if determineWavelengthSolution cannot find a good arc fit it + # uses a linear model based on the headers, whereas in SQ mode it raises an + # exception. For the purposes of checking arcs, we want everything to + # behave as if it is in SQ mode, so we just set that here + p.mode = 'sq' + + p.prepare(require_wcs=False) + p.addDQ() + p.addVAR(read_noise=True) + p.overscanCorrect() + p.ADUToElectrons() + p.addVAR(poisson_noise=True) + p.mosaicDetectors() + # We force center=None to ensure it always gets the central row. + # This is the default anyway, but this prevents it being overridden. + p.determineWavelengthSolution(center=None) + # We don't need determineDistortion to check Arcs, and it's slow. + # p.determineDistortion() + p.monitorWavelengthSolution() + p.writeOutputs() _default = makeProcessedArc diff --git a/geminidr/gnirs/tests/crossdispersed/test_determine_wavelength_solution.py b/geminidr/gnirs/tests/crossdispersed/test_determine_wavelength_solution.py index 294f449df2..590ab10426 100644 --- a/geminidr/gnirs/tests/crossdispersed/test_determine_wavelength_solution.py +++ b/geminidr/gnirs/tests/crossdispersed/test_determine_wavelength_solution.py @@ -328,9 +328,9 @@ def test_regression_determine_wavelength_solution( np.testing.assert_allclose(wavelength[indices], ref_wavelength[indices], atol=tolerance) print(f"Test passed for {ad.filename}, extension {wcalibrated_ext.id}") - except AssertionError: + except AssertionError as e: failed = True - raise + raise AssertionError(f"Test failed for {ad.filename}, extension {wcalibrated_ext.id}\n{repr(model)}\n{repr(ref_model)}") finally: if write_report: do_report(wcalibrated_ext, ref_ext, failed=failed) diff --git a/geminidr/gnirs/tests/longslit/test_determine_wavelength_solution.py b/geminidr/gnirs/tests/longslit/test_determine_wavelength_solution.py index 1d1faaa423..0c63059602 100644 --- a/geminidr/gnirs/tests/longslit/test_determine_wavelength_solution.py +++ b/geminidr/gnirs/tests/longslit/test_determine_wavelength_solution.py @@ -246,6 +246,7 @@ def test_regression_determine_wavelength_solution( pixel_scale = ad[0].pixel_scale() # arcsec / px p.viewer = geminidr.dormantViewer(p, None) + p.mode = 'qa' p.determineWavelengthSolution(**{**determine_wavelength_solution_parameters, **params}) diff --git a/gempy/library/fitting.py b/gempy/library/fitting.py index 9e7d59d1d4..df08f83fac 100644 --- a/gempy/library/fitting.py +++ b/gempy/library/fitting.py @@ -613,3 +613,32 @@ def translate_params(params): if "grow" in params: new_params["grow"] = params["grow"] return new_params + + @classmethod + def create_with_model(cls, model, image, points): + """ + Create a fit_1D object with an already-defined model and an array + of points. + + Parameters + ---------- + model : callable + A model that + + image : array-like + N-dimensional input array containing the values to be fitted. If + it is a `numpy.ma.MaskedArray`, any masked points are ignored when + fitting. + + points : `~numpy.ndarray`, optional + 1-dimensional input array with the x-values of each 1D slice being + fitted, If not given, the independent variable will be treated as + a sequence of integers starting at zero. + model + """ + fit1d = cls(np.arange(5), function="chebyshev", order=1, niter=0) + fit1d._models = model + fit1d.image = np.asarray(image) + fit1d.points = np.asarray(points) + fit1d.mask = np.zeros_like(fit1d.image, dtype=bool) + return fit1d diff --git a/gempy/library/wavecal.py b/gempy/library/wavecal.py index e02c65bace..2f169a6168 100644 --- a/gempy/library/wavecal.py +++ b/gempy/library/wavecal.py @@ -423,7 +423,7 @@ def create_interactive_inputs(ad, ui_params=None, p=None, linelist=None, bad_bits=0, absorption=False): data = {"x": [], "y": [], "meta": []} for ext in ad: - input_data, fit1d, _ = get_automated_fit( + input_data, fit1d = get_automated_fit( ext, ui_params, p=p, linelist=linelist, bad_bits=bad_bits, absorption=absorption) # peak locations and line wavelengths of matched peaks/lines @@ -486,8 +486,6 @@ class holding parameters for the UI, passed from the primitive's Config fit1d : a fit_1D object containing the wavelength solution, plus an "image" attribute that lists the matched arc line wavelengths - acceptable_fit : bool - whether this fit is likely to be good """ input_data = get_all_input_data( ext, p, ui_params.toDict(), linelist=linelist, bad_bits=bad_bits, @@ -499,14 +497,14 @@ class holding parameters for the UI, passed from the primitive's Config dw = np.diff(init_models[0](np.arange(spectrum.size))).mean() kdsigma = fwidth * abs(dw) k = 1 if kdsigma < 3 else 2 - fit1d, acceptable_fit = find_solution( + fit1d = find_solution( init_models, ui_params.toDict(), peaks=peaks, peak_weights=weights[ui_params.weighting], linelist=input_data["linelist"], fwidth=fwidth, kdsigma=kdsigma, k=k, bounds_setter=input_data["bounds_setter"], filename=ext.filename) input_data["fit"] = fit1d - return input_data, fit1d, acceptable_fit + return input_data, fit1d def create_chebyshev(waves, central_wavelength=None, dispersion=None, @@ -773,9 +771,9 @@ def find_solution(init_models, config, peaks=None, peak_weights=None, Returns ------- - length-2 tuple - A tuple of the best-fit model and a boolean denoting whether or not it - is an acceptable fit. + fit_1D: + the best-fit model or the original model, of no acceptable fit is + found. The "image" attribute will be empty if no fit is found. """ log = logutils.get_logger(__name__) min_lines = [int(x) for x in str(config["debug_min_lines"]).split(',')] @@ -805,57 +803,91 @@ def find_solution(init_models, config, peaks=None, peak_weights=None, min_lines_per_fit=min_lines_per_fit, k=k, bounds_setter=bounds_setter) + wavemin, wavemax = sorted(model(model.domain)) + nlines_expected = sum(wavemin <= wline <= wavemax for wline in arc_lines) + # We perform a regular least-squares fit to all the matches # we've made. This allows a high polynomial order to be # used without the risk of it going off the rails matched = np.where(matches > -1)[0] fit_it = fitting.LinearLSQFitter() - if len(matched) > 1: # need at least 2 lines, right? - m_init = models.Chebyshev1D(degree=min(config["order"], len(matched)-1), - domain=domain) + + # We allow this continue if we only matched 1 line but only 1 line is + # expected + if len(matched) > 1 or (len(matched) * nlines_expected == 1 and len(peaks) <= 2): + m_init = models.Chebyshev1D( + degree=min(config["order"], (len(matched) + 1) // 2), + domain=domain) for p, v in zip(model.param_names, model.parameters): if p in m_init.param_names: setattr(m_init, p, v) - #bounds_setter(m_init) + model_bounds = bounds_setter(m_init) + if len(matched) == 1: + m_init.c1.fixed = True #for i in range(len(matched), m_init.degree + 1): # m_init.fixed[f"c{i}"] = True matched_peaks = peaks[matched] matched_arc_lines = arc_lines[matches[matched]] m_final = fit_it(m_init, matched_peaks, matched_arc_lines) + dw = abs(np.diff(m_final(m_final.domain))[0] / np.diff(m_final.domain)[0]) #for p, l in zip(matched_peaks, matched_arc_lines): # print(f"{p:.2f} => {l:.2f}") - # We're close to the correct solution, perform a KDFit - m_init = models.Chebyshev1D(degree=config["order"], domain=domain) - for p, v in zip(m_final.param_names, m_final.parameters): - setattr(m_init, p, v) - dw = abs(np.diff(m_final(m_final.domain))[0] / np.diff(m_final.domain)[0]) - fit_it = matching.KDTreeFitter(sigma=2 * abs(dw), maxsig=5, - k=k, method='Nelder-Mead') - m_final = fit_it(m_init, peaks, arc_lines, in_weights=peak_weights, - ref_weights=arc_weights) - logit(f'{repr(m_final)} {fit_it.statistic}') - - # And then recalculate the matches - match_radius = 4 * fwidth * abs(m_final.c1) / len_data # 2*fwidth pixels - try: - matched = matching.match_sources(m_final(peaks), arc_lines, - radius=match_radius) - incoords, outcoords = zip(*[(peaks[i], arc_lines[m]) - for i, m in enumerate(matched) if m > -1]) - # Probably incoords and outcoords as defined here should go to - # the interactive fitter, but cull to derive the "best" model - fit1d = fit_1D(outcoords, points=incoords, function="chebyshev", - order=min(m_final.degree, len(incoords)-1), - domain=m_final.domain, - niter=config["niter"], sigma_lower=config["lsigma"], - sigma_upper=config["hsigma"]) - fit1d.image = np.asarray(outcoords) - except ValueError: - log.warning("Line-matching failed") - continue - nmatched = np.sum(~fit1d.mask) - logit(f"{filename} {repr(fit1d.model)} {nmatched} {fit1d.rms}") + if len(matched) > 1: + # We're close to the correct solution, perform a KDFit + m_init = models.Chebyshev1D(m_final.degree, domain=domain) + for p, v in zip(m_final.param_names, m_final.parameters): + setattr(m_init, p, v) + fit_it = matching.KDTreeFitter(sigma=2 * abs(dw), maxsig=5, + k=k, method='Nelder-Mead') + m_final = fit_it(m_init, peaks, arc_lines, in_weights=peak_weights, + ref_weights=arc_weights) + logit(f'{repr(m_final)} {fit_it.statistic}') + + # And then recalculate the matches + match_radius = 4 * fwidth * abs(m_final.c1) / len_data # 2*fwidth pixels + try: + matched = matching.match_sources(m_final(peaks), arc_lines, + radius=match_radius) + incoords, outcoords = zip(*[(peaks[i], arc_lines[m]) + for i, m in enumerate(matched) if m > -1]) + # Probably incoords and outcoords as defined here should go to + # the interactive fitter, but cull to derive the "best" model + fit1d = fit_1D(outcoords, points=incoords, function="chebyshev", + order=min(m_final.degree, (len(incoords) + 1) // 2), + domain=m_final.domain, + niter=config["niter"], sigma_lower=config["lsigma"], + sigma_upper=config["hsigma"]) + fit1d.image = np.asarray(outcoords) + except ValueError: + log.warning("Line-matching failed") + continue + nmatched = np.sum(~fit1d.mask) + logit(f"{filename} {repr(fit1d.model)} {nmatched} {fit1d.rms}") + + # A posteriori check that the fit is within the bounds + # (LinearLSQFitter does not allow bounds so we can't force + # the model to be bounded). We allow the tolerances to be twice + # as large as they were during the piecewise fit, as this seems + # to work well empirically. + try: + for p, v in zip(fit1d.model.param_names, fit1d.model.parameters): + bounds = model_bounds.get(p, (-np.inf, np.inf)) + new_bounds = (np.asarray(bounds) + + np.array([-0.5, 0.5]) * np.diff(bounds)[0]) + assert new_bounds[0] <= v <= new_bounds[1] + except AssertionError: + for p, v in zip(fit1d.model.param_names, fit1d.model.parameters): + bounds = model_bounds.get(p, (-np.inf, np.inf)) + new_bounds = (np.asarray(bounds) + + np.array([-0.5, 0.5]) * np.diff(bounds)[0]) + print(p, v, new_bounds) + continue + else: + # Hack a fit1D object with the model and the single match + fit1d = fit_1D.create_with_model( + m_final, image=matched_arc_lines, points=matched_peaks) + nmatched = 1 # Wavelength solution models need to be monotonic. Make that check. waves = fit1d.evaluate(np.arange(len_data)) @@ -870,7 +902,7 @@ def find_solution(init_models, config, peaks=None, peak_weights=None, # Trial and error suggests this criterion works well if fit1d.rms < 0.8 / config["order"] * fwidth * abs(dw) and nmatched >= min_matches_required: #print("RETURNING", fit1d.model.parameters) - return fit1d, True + return fit1d # This seems to be a reasonably ranking for poor models if nmatched > config["order"] + 1: @@ -883,13 +915,8 @@ def find_solution(init_models, config, peaks=None, peak_weights=None, if best_fit1d is None: # Hack a fit1D object that represents the original model with no fitted lines - best_fit1d = fit_1D(np.arange(5), function="chebyshev", order=1, - niter=0) - best_fit1d._models = init_models[0] - best_fit1d.image = np.array([]) - best_fit1d.points = np.array([]) - best_fit1d.mask = np.array([], dtype=bool) - return best_fit1d, True + best_fit1d = fit_1D.create_with_model(init_models[0], [], []) + return best_fit1d def perform_piecewise_fit(model, peaks, arc_lines, pixel_start, kdsigma, @@ -1168,6 +1195,8 @@ def update_wcs_with_solution(ext, fit1d, input_data, config): if len(incoords) > 1: inv_rms = np.std(m_inverse(m_final(incoords)) - incoords) log.stdinfo(f"Inverse model has rms = {inv_rms:.3f} pixels.") + else: + inv_rms = np.nan m_final.name = "WAVE" # always WAVE, never AWAV m_final.inverse = m_inverse @@ -1182,8 +1211,8 @@ def update_wcs_with_solution(ext, fit1d, input_data, config): # while I refactor tests temptable.add_columns([[1], [m_final.degree], [domain[0]], [domain[1]]], names=("ndim", "degree", "domain_start", "domain_end")) - temptable.add_columns([[rms], [input_data["fwidth"]]], - names=("rms", "fwidth")) + temptable.add_columns([[rms], [inv_rms], [input_data["fwidth"]]], + names=("rms", "inv_rms", "fwidth")) if ext.data.ndim > 1: # TODO: Need to update this from the interactive tool's values direction, location = input_data["location"].split() @@ -1276,14 +1305,15 @@ def create_pdf_plot(input_data, peaks, arc_lines, title="", data = input_data["spectrum"] spacing = 0.01 vert_align = "bottom" - xmin, xmax = input_data["init_models"][0].domain + init_model = input_data["init_models"][0] + xmin, xmax = init_model.domain pixels = np.arange(xmin, xmax + 1) data_max = data.max() spacing *= data_max fig, ax = plt.subplots() ax.plot(pixels, data, 'b-') ax.set_ylim(0, data_max * 1.1) - if len(arc_lines) and np.diff(arc_lines)[0] / np.diff(peaks)[0] < 0: + if len(arc_lines) and init_model(xmin) > init_model(xmax): ax.set_xlim(xmax + 1, xmin - 1) else: ax.set_xlim(xmin - 1, xmax + 1) From 9ca14b8abdfc5fcb487f63c9659c33dffad4e171 Mon Sep 17 00:00:00 2001 From: Paul Hirst Date: Wed, 6 May 2026 09:57:49 -1000 Subject: [PATCH 03/12] Comandeered 4.2.1-goawavecal as a version to identify this work in GOA at least temporarily until Kathleen gets back. Updated changes.rst with notes on added improvements --- astrodata/_version.py | 2 +- doc/DRAGONS/changes.rst | 25 +++++++++++++++++++++++++ 2 files changed, 26 insertions(+), 1 deletion(-) diff --git a/astrodata/_version.py b/astrodata/_version.py index 4e45b93839..6fe608c23f 100644 --- a/astrodata/_version.py +++ b/astrodata/_version.py @@ -8,7 +8,7 @@ API = 4 FEATURE = 2 BUG = 1 -TAG = '' +TAG = 'wavecal' def version(short=False, tag=TAG): diff --git a/doc/DRAGONS/changes.rst b/doc/DRAGONS/changes.rst index 5c387405aa..4b3fa2da72 100644 --- a/doc/DRAGONS/changes.rst +++ b/doc/DRAGONS/changes.rst @@ -8,6 +8,31 @@ Change Logs *********** +4.2.1-goawavecal +================ + +Improvements +------------ + +This internal (used in GOA) release includes the following improvements to +the determineWavelengthSolution primitive: + +* Check the final fit a-posteriori to verify that it is within a factor of two + of the initial fit boundary constraints and rejects it if not +* In Science-Quality mode, if no arc lines are found for all orders being + processed, an exception is raised and processing stops. In Quick-Look mode, + the previous behavior of adding a linear wavelength solution based on the + headers is unchanged +* If only one arc line is found, a linear fit based on the headers is shifted + to match the one arc line, a warning is generated, and that fit is returned as + the wavelength solution + +This release also includes functionality in gempy to embed auxillary files +(intended for eg PDF plot figures) in binary tables in astrodata instances, and +hence FITS files. This is used for example to capture the PDFs of the arc line +identifications into the reduced arc FITS files, and a corresponding script +in gempy to extract them. + 4.2.1 ===== From 64e9fb8e432eb3b8652f2db87030e4239e9956f4 Mon Sep 17 00:00:00 2001 From: Paul Hirst Date: Wed, 6 May 2026 10:08:35 -1000 Subject: [PATCH 04/12] oops, fix wavecal to goawavecal for consistency --- astrodata/_version.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/astrodata/_version.py b/astrodata/_version.py index 6fe608c23f..2752bcd30e 100644 --- a/astrodata/_version.py +++ b/astrodata/_version.py @@ -8,7 +8,7 @@ API = 4 FEATURE = 2 BUG = 1 -TAG = 'wavecal' +TAG = 'goawavecal' def version(short=False, tag=TAG): From 7fd17c89bc06314ea8c98c31c65e1ed9913f9c1a Mon Sep 17 00:00:00 2001 From: Paul Hirst Date: Fri, 22 May 2026 11:12:47 -1000 Subject: [PATCH 05/12] fix potential race condition with multiple processes running CachedFileGetter --- recipe_system/cal_service/file_getter.py | 14 +++++++++----- 1 file changed, 9 insertions(+), 5 deletions(-) diff --git a/recipe_system/cal_service/file_getter.py b/recipe_system/cal_service/file_getter.py index 039065e2d4..9bed303b92 100644 --- a/recipe_system/cal_service/file_getter.py +++ b/recipe_system/cal_service/file_getter.py @@ -125,11 +125,15 @@ def get_request(self, url, filename, calmd5): tuples = [] total_bytes = 0 for dirent in dirents: - statobj = dirent.stat() # To avoid multiple stat calls - tuples.append((statobj.st_atime, - dirent.path, - statobj.st_size)) - total_bytes += statobj.st_size + try: + statobj = dirent.stat() # To avoid multiple stat calls + tuples.append((statobj.st_atime, + dirent.path, + statobj.st_size)) + total_bytes += statobj.st_size + except FileNotFoundError: + # Race condition where another process deleted the file + pass # Sort list by increasing atime (ie oldest first) tuples.sort() while total_bytes > (1E9 * self.cachegbs): From 54ad3c0c16b7934c2092195afeccafeffd001e15 Mon Sep 17 00:00:00 2001 From: chris-simpson Date: Wed, 20 May 2026 14:59:53 -1000 Subject: [PATCH 06/12] GMOSLongslit.addIllumMaskToDQ(): avoid crash when no cross-correlation peaks are found and proceed in the same way as when 2+ are found --- geminidr/gmos/primitives_gmos_longslit.py | 17 ++++++++++++----- 1 file changed, 12 insertions(+), 5 deletions(-) diff --git a/geminidr/gmos/primitives_gmos_longslit.py b/geminidr/gmos/primitives_gmos_longslit.py index 3b22e68f0c..a3ff8723db 100644 --- a/geminidr/gmos/primitives_gmos_longslit.py +++ b/geminidr/gmos/primitives_gmos_longslit.py @@ -242,13 +242,20 @@ def addIllumMaskToDQ(self, adinputs=None, suffix=None, illum_mask=None, plt.show() plt.ion() - yshift = mshift - maxima[0] - if len(maxima) > 1 or abs(yshift) > mshift: - log.warning(f"{ad.filename}: cross-correlation peak is" - " untrustworthy so not adding illumination " - "mask. Please re-run with a specified shift.") + if not maxima: + log.warning(f"{ad.filename}: no peaks found in " + "cross-correlation, so not adding " + "illumination mask. Please re-run with a specified shift.") yshift = None log.stdinfo(slit_location_msg) + else: + yshift = mshift - maxima[0] + if len(maxima) > 1 or abs(yshift) > mshift: + log.warning(f"{ad.filename}: cross-correlation peak is" + " untrustworthy so not adding illumination " + "mask. Please re-run with a specified shift.") + yshift = None + log.stdinfo(slit_location_msg) else: yshift = shift From 4056355328f0ac64f59952d695b65fbc97afd86a Mon Sep 17 00:00:00 2001 From: Paul Hirst Date: Wed, 3 Jun 2026 11:28:36 -1000 Subject: [PATCH 07/12] tweak readnoise estimate from overscan - apply N/N-1 ddof correction, and disregard first and last rows due to outlying values --- geminidr/core/primitives_ccd.py | 36 ++++++++++++++++++++++----------- 1 file changed, 24 insertions(+), 12 deletions(-) diff --git a/geminidr/core/primitives_ccd.py b/geminidr/core/primitives_ccd.py index 3570f36145..2f9e25823f 100644 --- a/geminidr/core/primitives_ccd.py +++ b/geminidr/core/primitives_ccd.py @@ -225,7 +225,7 @@ def subtractOverscan(self, adinputs=None, **params): # need to match asec location and size pixels = np.arange(asec.y1, asec.y2) data = np.mean(ext.data[asec.y1:asec.y2, x1:x2], axis=axis) - stds = np.std(ext.data[asec.y1:asec.y2, x1:x2], axis=axis) + stds = np.std(ext.data[asec.y1:asec.y2, x1:x2], axis=axis, ddof=1) else: if y1 > asec.y1: # Bias on top y1 += nbiascontam @@ -240,17 +240,29 @@ def subtractOverscan(self, adinputs=None, **params): data = np.mean(ext.data[y1:y2, asec.x1:asec.x2], axis=axis) stds = np.std(ext.data[y1:y2, asec.x1:asec.x2], axis=axis) - # We have 1 sample (of x2-x1 values) from each of N - # populations each of which has a different population mean, - # (ie the bias leve for that row) but we assume they all - # have the same standard deviation (ie the read noise) and - # we want to estimate that population standard deviation. - # We've calculated the mean (data) and standard deviation - # (stds) of each of the samples, and want to estimate the - # population standard deviation. Because the means are - # different, we can't just calculate the std of all the - # samples. - overstd = np.sqrt(np.mean(stds * stds)) + # We have (x2-x1) samples from each of N populations each + # of which has a different population mean, (ie the bias + # level for that row) but we assume they all have the same + # standard deviation (ie the read noise) and we want to + # estimate that population standard deviation. + # We've calculated the mean (data) and estimates of the + # population standard deviation[*] (stds) of each of the set + # of samples, and want to estimate the population standard + # deviation. Because the means are different, we can't just + # calculate the std of all the samples. + # [*] [by which I mean to say we set ddof=1 in the call to + # np.stds to apply the (N/N-1) factor to estimate the + # population standard deviation from the calculaiton of the + # sample standard deviation. + # We use the mean of these per-row population variance + # estimates as our estimate of the global population + # variance. This seems to work well, but I'm not sure if it + # is truly an optimal estimator. + # It turns out that the first and last row of the readout + # often have pixels with massively outlying values, so we + # disregard those samples + trimmedstds = stds[1:-1] + overstd = np.sqrt(np.mean(trimmedstds * trimmedstds)) # Readnoise descriptor always returns value in electrons if ext.is_in_adu(): From de54cbd0bbe537667619becd10bcaa6f6c57cbe5 Mon Sep 17 00:00:00 2001 From: Paul Hirst Date: Wed, 3 Jun 2026 13:43:36 -1000 Subject: [PATCH 08/12] make the OVERRDNS header readnoise value always be in electrons --- geminidr/core/primitives_ccd.py | 5 ++++- geminidr/gemini/lookups/keyword_comments.py | 2 +- 2 files changed, 5 insertions(+), 2 deletions(-) diff --git a/geminidr/core/primitives_ccd.py b/geminidr/core/primitives_ccd.py index 2f9e25823f..5ef471d458 100644 --- a/geminidr/core/primitives_ccd.py +++ b/geminidr/core/primitives_ccd.py @@ -321,7 +321,10 @@ def subtractOverscan(self, adinputs=None, **params): ext.hdr.set('OVERSCAN', previous_overscan + bias_level, self.keyword_comments['OVERSCAN']) ext.hdr.set('OVERRMS', sigma, self.keyword_comments['OVERRMS']) - ext.hdr.set('OVERRDNS', overstd, self.keyword_comments['OVERRDNS']) + # By convention, "readnoise" is always in electrons + overrdns = overstd * ext.gain() if ext.is_in_adu() else overstd + ext.hdr.set('OVERRDNS', overrdns, self.keyword_comments['OVERRDNS']) + for desc in ('saturation_level', 'non_linear_level'): with suppress(AttributeError, KeyError): ext.hdr[ad._keyword_for(desc)] -= bias_level diff --git a/geminidr/gemini/lookups/keyword_comments.py b/geminidr/gemini/lookups/keyword_comments.py index 8d6769abf4..34a84ada49 100644 --- a/geminidr/gemini/lookups/keyword_comments.py +++ b/geminidr/gemini/lookups/keyword_comments.py @@ -48,7 +48,7 @@ "OBSTYPE": "Observation type", "ORIGNAME": "Original filename prior to processing", "OVERRMS": "RMS of fit to overscan mean values", - "OVERRDNS": "Readnoise estimate derived from overscan", + "OVERRDNS": "Readnoise estimate from overscan [electron]", "OVERSCAN": "Median of Overscan mean values", "OVERSEC": "Section used for overscan calculation", "PIXSCALE": "Pixel scale [arcsec/pixel]", From df38ae60b0872bcff1f206293f431a5ea988a511 Mon Sep 17 00:00:00 2001 From: chris-simpson Date: Thu, 4 Jun 2026 08:38:52 -1000 Subject: [PATCH 09/12] use link_layouts=True in ratio/residual tabs for bokeh 3.9 --- geminidr/interactive/fit/fit1d.py | 4 ++++ 1 file changed, 4 insertions(+) diff --git a/geminidr/interactive/fit/fit1d.py b/geminidr/interactive/fit/fit1d.py index bac54f01b8..6c95361395 100644 --- a/geminidr/interactive/fit/fit1d.py +++ b/geminidr/interactive/fit/fit1d.py @@ -2291,6 +2291,10 @@ def fit1d_figure( tabs = bm.Tabs( tabs=[], sizing_mode="stretch_width", stylesheets=dragons_styles() ) + try: # needed for bokeh 3.9 + tabs.link_layouts = True + except AttributeError: # bokeh < 3.9 + pass tabs.tabs.append(bm.TabPanel(child=p_resid, title="Residuals")) tabs.tabs.append(bm.TabPanel(child=p_ratios, title="Ratios")) From ae5f3643f32192a0da88d63ac29e0d6390d7a76c Mon Sep 17 00:00:00 2001 From: chris-simpson Date: Thu, 4 Jun 2026 08:49:38 -1000 Subject: [PATCH 10/12] actually pass plot_ratios/residuals to the figure-creation code! --- geminidr/interactive/fit/fit1d.py | 2 ++ 1 file changed, 2 insertions(+) diff --git a/geminidr/interactive/fit/fit1d.py b/geminidr/interactive/fit/fit1d.py index 6c95361395..52e4b34116 100644 --- a/geminidr/interactive/fit/fit1d.py +++ b/geminidr/interactive/fit/fit1d.py @@ -1238,6 +1238,8 @@ def build_figures( ylabel=self.ylabel, model=self.model, enable_user_masking=self.enable_user_masking, + plot_ratios=plot_ratios, + plot_residuals=plot_residuals, ) if self.enable_regions: From 183748841a1c51b4a9f39fc441d5b8144bf9dea9 Mon Sep 17 00:00:00 2001 From: chris-simpson Date: Thu, 4 Jun 2026 08:50:44 -1000 Subject: [PATCH 11/12] do not plot fit-to-data ratio (only residuals) in calculateSensitivity and traceApertures --- geminidr/core/primitives_spect.py | 1 + geminidr/interactive/fit/tracing.py | 1 + 2 files changed, 2 insertions(+) diff --git a/geminidr/core/primitives_spect.py b/geminidr/core/primitives_spect.py index 8e9dc30aef..dea0da32e4 100644 --- a/geminidr/core/primitives_spect.py +++ b/geminidr/core/primitives_spect.py @@ -1219,6 +1219,7 @@ def _get_fit1d_input_data(ext, exptime, spec_table): primitive_name="calculateSensitivity", filename_info=filename_info, help_text=CALCULATE_SENSITIVITY_HELP_TEXT, + plot_ratios=False, ui_params=uiparams) geminidr.interactive.server.interactive_fitter(visualizer) diff --git a/geminidr/interactive/fit/tracing.py b/geminidr/interactive/fit/tracing.py index 69e45e505d..62a65962d8 100644 --- a/geminidr/interactive/fit/tracing.py +++ b/geminidr/interactive/fit/tracing.py @@ -75,6 +75,7 @@ def interactive_trace_apertures(ad, tab_labels, fit1d_params, ui_params: UIParam modal_button_label="Trace apertures", modal_message="Tracing apertures...", ui_params=ui_params, + plot_ratios=False, turbo_tabs=True ) From b9f16a182c898a0869bfaf315ba81d048023839c Mon Sep 17 00:00:00 2001 From: chris-simpson Date: Mon, 8 Jun 2026 10:09:43 -1000 Subject: [PATCH 12/12] changelog updated --- doc/DRAGONS/changes.rst | 18 ++++++++++++++++++ 1 file changed, 18 insertions(+) diff --git a/doc/DRAGONS/changes.rst b/doc/DRAGONS/changes.rst index 4b3fa2da72..003c60ffb4 100644 --- a/doc/DRAGONS/changes.rst +++ b/doc/DRAGONS/changes.rst @@ -8,6 +8,24 @@ Change Logs *********** +4.2.2 +===== + +Bug fixes +--------- +** geminidr.gmos ** + +* Prevent ``addIllumMaskToDQ`` from crashing when the GMOS longslit + bridges cannot be located in the spatial profile, and instead output a + warning. + +Compatibility +------------- +A change to the behavior of Tabs in bokeh v3.9.0 caused the residual/ratio +plot in the interactive tools to be corrupted. This release ensures +compatibility with v3.9.0 as well as earlier versions. + + 4.2.1-goawavecal ================