Rectify 2D spec: new script by debora-pe · Pull Request #2036 · pypeit/PypeIt

debora-pe · 2025-11-22T00:51:28Z

This add a new script that creates a rectified 2D spectral image.

Vet test added here: pypeit/PypeIt-development-suite#383

NOTE: a PR is also been opened here to allow ginga to display this rectified spectral image.

kbwestfall

Thanks for doing this! I have a couple of minor comments.

kbwestfall · 2025-11-24T17:47:57Z

+                    hdu.header['CUNIT2'] = 'Angstrom'
+                    hdu.header['CDELT2'] = dsamp
+                    hdu.header['CRPIX2'] = 0
+                    hdu.header['CRVAL2'] = wave_grid[0]


I don't think this should hold up merging the PR, but I'm curious if there are any Astropy utilities that we can use to build the WCS. In particular, I'm wondering if the CTYPE values follow the right standard. And, we might want to allow people to bin the spectra with log steps in wavelength.

Thanks for pointing this out. I have tried to use the astropy utilities for WCS and it does not really work well, but it may be me not understanding it. The "Angstrom" CUNIT is not really recognized and it ends up being "m" in the header (converting also the CDELT value to m). Also, the standard FITS keyword for wavelength should be "WAVE" but it doesn't look like ds9 reads it very well, so I opted to "LAMBDA" which is the name used also by the AllSlits files from the DEEP2 pipeline. If you have any suggestion I can definitely try other options.

jhennawi

I'm not following why you are interpolating the fluxes and errors given that this is something we try to avoid in pypeit 2d cladding at all costs. The 2D coadding was designed to avoid interpolation, and it already rectifies by construction. So I'm not really seeing why you need to interpolate here.

jhennawi · 2025-11-25T12:55:16Z

+                        ivar = imgrect_dict['sciivar'][:nspec_slit, ispat_slit]
+
+                        # Interpolate onto a common wavelength grid
+                        flux_interp = np.interp(wave_grid, wave_slit, flux,


I'm not following what you are doing here. Why are you interpolating the flux values and the errors. This generates correlated errors. We go to great pains in pypeit to avoid interpolation, and so I'm very concerned about seeing this here. The coadd2d algorithm was designed to not interpolate to maintain the integrity of our errors.

Thanks for pointing this out, @jhennawi ! @debora-pe , are you doing this to put all the slits into a single image? Can you instead figure out the full grid (spatial + spectral) you need to accommodate all the slits first, and then use coadd.compute_coadd2d or coadd.rebin2d to perform the rectification of each slit directly into a section of a single image? If we can't, we should be using pypeit.sampling.Resample instead of the interpolation.

jhennawi · 2025-11-25T13:01:12Z

        interp_dspat (bool, optional):
           Interpolate in the spatial coordinate image to faciliate running
           through core.extract.local_skysub_extract. This can be slow.   Default=True.
+        single_frame (bool, optional):


I don't see why you need a new boolean variable here. If nimgs=1, it is clear that you are only performing rectification and not coadding, since you cannot coadd one image.

jhennawi · 2025-11-25T13:02:17Z

        True, the rejection threshold is set based on the number of images to
        combine; see above.  This value is passed directly to
        `astropy.stats.SigmaClip`_ as its ``sigma`` parameter.
    maxiters : :obj:`int`, optional, default=5


Rather than add a new argument, why not just get rid of this warning message altogether.

debora-pe · 2025-12-02T02:41:36Z

I'm not following why you are interpolating the fluxes and errors given that this is something we try to avoid in pypeit 2d cladding at all costs. The 2D coadding was designed to avoid interpolation, and it already rectifies by construction. So I'm not really seeing why you need to interpolate here.

Hey @jhennawi . I have removed the interpolation and the new argument. What do you think about this now?

tbowers7 · 2025-12-02T15:49:42Z

NOTE: a PR is also been opened here to allow ginga to display this rectified spectral image.

@debora-pe , the mentioned ginga PR seems to have been closed without merging. What is the status of displaying the rectified image in ginga?

debora-pe · 2025-12-02T17:18:45Z

NOTE: a PR is also been opened here to allow ginga to display this rectified spectral image.

@debora-pe , the mentioned ginga PR seems to have been closed without merging. What is the status of displaying the rectified image in ginga?

@tbowers7 Thank you for checking this. Eric issued another PR here with the same changes and that was merged. It should be available in the next release by the end of the year.

tbowers7

This looks helpful for QL, but could you explain how this script is different from using (a modified version of) CoAdd2d, which also does image rectification?

Approving.

tbowers7 · 2025-12-02T18:42:08Z

    dwave_eff = dwave*spec_samp_fact
    if log10:
-        ngrid = np.ceil((np.log10(wave_max) - np.log10(wave_min))/dwave_eff).astype(int)
+        ngrid = np.ceil((np.log10(wave_max) - np.log10(wave_min))/dwave_eff).astype(int) + 1
        loglam_grid = np.log10(wave_min) + dwave_eff*np.arange(ngrid)
        wave_grid = np.power(10.0,loglam_grid)
        loglam_grid_mid = np.log10(wave_grid) + dwave_eff/2.0
        wave_grid_mid = np.power(10.0, loglam_grid_mid)
    else:
-        ngrid = np.ceil((wave_max - wave_min)/dwave_eff).astype(int)
+        ngrid = np.ceil((wave_max - wave_min)/dwave_eff).astype(int) + 1
        wave_grid = wave_min + dwave_eff*np.arange(ngrid)
        wave_grid_mid = wave_grid + dwave_eff/2.0


These wavelength grids are ultimately consumed by pypeit.core.coadd, pypeit.coadd2d, and pypeit.sensfunc.

What is the purpose of lengthening the grids by one in the redward direction?

What is the effect on downstream functions of this lengthening? Alternatively, were the previous results wrong or problematic?

I noticed that when I created the wavelength grid providing the wave min and max as arguments, the grid produced ended always one spectral pixel short, i.e., the last value of the wavelength grid was ~1 spectral pixel smaller than the set wave max. I have noticed also from the docstring that wave_grid should have the shape ngrid + 1 and that made me think that there was a missing +1. I'm planning to run the tests to make sure this didn't cause problems.

tbowers7 · 2025-12-02T18:57:25Z

+                    # Interpolate onto a common wavelength grid
+                    for ispat_slit in range(nspat_vec[islit]):
+                        # Get data for this spatial pixel
+                        flux = imgrect_dict['imgminsky'][:nspec_slit, ispat_slit]
+                        ivar = imgrect_dict['sciivar'][:nspec_slit, ispat_slit]
+                        valid = (wave_grid_mid >= wave_slit.min()) & \
+                                (wave_grid_mid <= wave_slit.max())
+
+                        # get the spectral pixel where the slit should start/stop
+                        wstart = np.where(np.isclose(wave_grid_mid, wave_slit[0]))[0][0]
+                        wend = np.where(np.isclose(wave_grid_mid, wave_slit[-1]))[0][0] + 1
+
+                        # Assign to output arrays
+                        image_rect[wstart:wend, spat_left + ispat_slit] = flux
+                        ivar_rect[wstart:wend, spat_left + ispat_slit] = ivar
+                        mask_rect[:, spat_left + ispat_slit] = np.logical_not(valid)
+                        wave_grid_rect2[wstart:wend, spat_left + ispat_slit] = wave_slit


This isn't strictly interpolating, but I'm not exactly sure what's happening to ensure that there is a 1:1 correspondence between the input arrays (flux, ivar, wave_slit) and the output array slice (image_rect[...], etc.). Could you elaborate on how you ensure this without the interpolation?

I have made some clean up here and hopefully explained better with inline comments what it is happening. Each slits in the output image have exactly the same wavelength grid. Therefore, according to what is the wavelength of each wave_slit I find the pixels where to place the slit in the output image. Basically, wave_grid_mid and wave_slit have exactly the same wavelength grid, except that the wave_slit array has the wavelength values for only that particular slit, while wave_grid_mid has the value for the entire grid.

debora-pe · 2025-12-03T00:09:52Z

This looks helpful for QL, but could you explain how this script is different from using (a modified version of) CoAdd2d, which also does image rectification?

Approving

Yes, the purpose it's exactly for QL, this is why it's a very simple FITS file, viewable easily with DS9 and ginga. CoAdd2d does the rectification, but it does not put all the slits in the same wavelength grid, therefore from the location on the coadded image you cannot, e.g., spot features at similar redshifts.

debora-pe · 2025-12-04T22:45:00Z

+                    # get the spectral pixel where the slit should start/stop
+                    _wave_grid_mid = np.round(wave_grid_mid, 4)
+                    _wave_slit = np.round(imgrect_dict['wave_mid'], 4)
+                    wstart = np.where(_wave_grid_mid == _wave_slit[0])[0][0]


I went back and forward between using == or np.close here, and at the end I decided to run
pypeit_rectify_2dspec PypeIt-development-suite/REDUX_OUT/*/*/*/spec2d*.fits
That failed pretty soon for np.close, instead it was successful for all Dev Suite spec2d with ==

debora-pe · 2025-12-04T22:50:48Z

Run the tests and all looks good. There was a problem with a kcwi dataset and consequent vet test failure because of a missing spec1d, but that was a small glitch in the syncing with RAW_DATA and was all fixed after re-running without changes.
One vet test failure is my usual mysterious failure that happens only on my machine
FAILED ../test_scripts.py::test_identify - AssertionError: Bad auto ID

Test Summary
--------------------------------------------------------
--- PYTEST PYPEIT UNIT TESTS PASSED  280 passed, 467 warnings in 443.84s (0:07:23) ---
--- PYTEST UNIT TESTS PASSED  158 passed, 93407 warnings in 1077.08s (0:17:57) ---
--- PYTEST VET TESTS FAILED  2 failed, 68 passed, 120008 warnings in 2262.42s (0:37:42) ---
--- PYPEIT DEVELOPMENT SUITE FAILED 1/284 TESTS  ---
Failed tests:
    keck_kcwi/small_bh2_4200 pypeit
Skipped tests:
Testing Started at 2025-12-02T14:28:53.924552
Testing Completed at 2025-12-03T12:37:23.106625
Total Time: 22:08:29.182073

./pypeit_test reduce -i keck_kcwi -s small_bh2_4200
Running reduce tests.
Running tests on the following instruments:
    keck_kcwi

Reducing data from keck_kcwi for the following setups:
    small_bh2_4200

 0 passed/ 0 failed/ 0 skipped STARTED keck_kcwi/small_bh2_4200 pypeit
 1 passed/ 0 failed/ 0 skipped PASSED  keck_kcwi/small_bh2_4200 pypeit
--- PYPEIT DEVELOPMENT SUITE PASSED 1/1 TESTS  ---

 pytest vet_tests//test_datacube.py::test_coadd_datacube
====================================================================================================================================== test session starts =======================================================================================================================================
platform darwin -- Python 3.13.5, pytest-8.4.2, pluggy-1.6.0
PyQt6 6.9.1 -- Qt runtime 6.9.2 -- Qt compiled 6.9.0
rootdir: /Users/dpelliccia/PypeIt/PypeIt-development-suite
configfile: pyproject.toml
plugins: astropy-0.11.0, mock-3.15.1, remotedata-0.4.1, hypothesis-6.140.3, filter-subpackage-0.2.0, doctestplus-1.4.0, astropy-header-0.2.2, qt-4.5.0, arraydiff-0.6.1, cov-7.0.0
collected 1 item                                                                                                                                                                                                                                                                                 

vet_tests/test_datacube.py .                                                                                                                                                                                                                                                               [100%]

======================================================================================================================================== warnings summary ========================================================================================================================================
vet_tests/test_datacube.py::test_coadd_datacube
   RuntimeWarning: divide by zero encountered in divide (specobj.py:604)

-- Docs: https://docs.pytest.org/en/stable/how-to/capture-warnings.html
=========================================================================================================================== 1 passed, 1 warning in 4783.60s (1:19:43) ============================================================================================================================

debora-pe added 9 commits November 14, 2025 15:08

the beginning

0ebf10b

almost done

ea01063

change name

7dbd836

small changes

bc7886d

finish it

442cbfa

small

04d3903

Merge branch 'develop' into allslits

3a2c8a1

add docs

c935a46

updated docs

66fe60d

debora-pe mentioned this pull request Nov 22, 2025

Vet test for pypeit_rectify_2dspec pypeit/PypeIt-development-suite#383

Merged

updated docs

e52e7ae

kbwestfall self-requested a review November 24, 2025 16:32

kbwestfall approved these changes Nov 24, 2025

View reviewed changes

kbwestfall requested a review from tbowers7 November 24, 2025 17:49

K's comments

f93a81b

jhennawi reviewed Nov 25, 2025

View reviewed changes

remove interp

9ed40bf

kbwestfall added this to the v2.0.0 milestone Dec 1, 2025

debora-pe added 2 commits December 1, 2025 14:44

Merge branch 'develop' into allslits

ecca929

joe's comments

212336d

tbowers7 approved these changes Dec 2, 2025

View reviewed changes

clean up + improvements

c8b241f

debora-pe added 2 commits December 2, 2025 14:21

updated docs

1ce9d6a

last changes

17b2962

debora-pe commented Dec 4, 2025

View reviewed changes

kbwestfall merged commit ed72963 into develop Dec 5, 2025
29 checks passed

kbwestfall deleted the allslits branch December 5, 2025 21:56

Uh oh!

Conversation

debora-pe commented Nov 22, 2025 • edited by tbowers7 Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Uh oh!

kbwestfall left a comment

Choose a reason for hiding this comment

Uh oh!

Uh oh!

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

jhennawi left a comment

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

debora-pe commented Dec 2, 2025

Uh oh!

tbowers7 commented Dec 2, 2025

Uh oh!

debora-pe commented Dec 2, 2025

Uh oh!

tbowers7 left a comment

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

debora-pe commented Dec 3, 2025

Uh oh!

Choose a reason for hiding this comment

Uh oh!

debora-pe commented Dec 4, 2025

Uh oh!

Uh oh!

Reviewers

Assignees

Labels

Projects

Milestone

Development

Uh oh!

4 participants

debora-pe commented Nov 22, 2025 •

edited by tbowers7

Loading