resample problem for NIRSpec FS level 3

[stscijgbot-jp]

Issue JP-2083 was created on JIRA by James Muzerolle:

The output of level 3 processing of a set of FS simulated exposures does not look as expected.  The data comprise a set of 9 simulated exposures taken with a 3-point primary nod pattern with spectral subdithers.  For this level 3 test, I only included 3 exposures, the ones exactly at the primary nod positions.  In the s2d product, there are single positive and negative source traces, with a second negative trace cut off at the bottom of the 2D image (see first attachment). There are two additional areas of negative flux near the top, separated by wide swaths of pixels with values of 0.  Looking at the CON extension, the coverage map looks funny - the largest area of overlap has a value of 6, while I would have expected 7 (as 3 nodded exposures went into this); there is only a very small region of pixels with a value of 7 at the extreme right end of the 2D image (see second attachment).  There are also disconnected regions with a value of 0, as indicted in the actual spectrum.  As far as I can tell, the coordinate keywords are all populated correctly; the RA and DEC of the source traces in the level 2b products seem to be consistent among the three exposures.

Note: the level 3 association file I used is also in the central store directory, name "spec3_short_asn.json".  (The other spec3 json file is for the entire set of 9 exposures, whose combined product I haven't included.  That one has an even worse problem with the source traces not overlapping, but I'm still checking that the coordinate keywords are correct.)

[jdavies-st]

The roll angle matters a great deal in the simulation. I.e. each nodded image has to have unique values of RA_REF, DEC_REF and ROLL_REF which will all vary in a quite specific way for each nod.

Correct values for V2_REF, V3_REF, V3I_YANG and VPARITY are also important, which do not change with nods, but are pulled from the SIAF for the given aperture.

All of these keywords need to be simulated correctly and consistently for anything to work.

[stscijgbot-jp]

Comment by Howard Bushouse on JIRA:

James Muzerolle can you remind me what the spatial ordering should be for an observation using 3 primary nods? I think it starts with the source at the "top" of the slit as seen in the frame of our resampled science images and steps downward from there, right? So then the exposures subtracted as background from the first nod would result in two negative traces below the positive trace in the exposure corresponding to the first nod. And if the first nod is used as the "reference" exposure for building the combined output frame during level 3 processing, you would expect to see the 2 negative traces towards the bottom, as shown in your first image. If this were MOS, I'd then say that the additional negative traces near the top could be the negative traces from the next slitlet above, but since this is FS, that can't be the case. So there definitely are problems of some kind.

[stscijgbot-jp]

Comment by James Muzerolle on JIRA:

Howard Bushouse that's correct, this primary dither pattern goes top-middle-bottom.  I would have expected to see 2 negative traces both above and below - if the third dither is shifted to line up with the first, the two negative traces from that should then lie above the overlapping positive traces.

[stscijgbot-jp]

Comment by Howard Bushouse on JIRA:

Ah yes, you're right. You'd have 2 negative traces below the source in the first resampled image by itself, but then once you shift all of the background-subtracted images to line up the positive traces you'd end up with negative traces both above and below. You might only end up seeing 1 negative trace on either side, depending on how big the output frame ends up being in the cross-dispersion (nod) direction.

[stscijgbot-jp]

Comment by James Muzerolle on JIRA:

I did some more investigation by reprocessing the data using different roll angles (the original case used PA=38 degrees), no nod subtraction.  I think there may be a fundamental limitation with the way resample_spec is currently designed to work with RA and DEC coordinates.

I first tried a roll angle of zero, which aligns the slit along roughly constant RA (modulo the distortion).  In that case, the s2d products are almost completely blank with all zeroes except for a very narrow set of pixels at the middle of the 2d cutout (3rd attached figure).  I think this may be happening because the pixel RA values has almost no variation over the 2d input.  Moreover, the range of RA values does not actually encompass the source RA position for any of the input exposures except for the one where the source is near the center of the slit, so the level 3 product is similarly messed up.

I then tried a roll angle of 45 degrees, and get much improved products at both level 2b and 3.  The 4th attachment shows the latter, which includes all 9 of the simulated exposures (3 point primary patter with spectral subdithers).  The source overlap looks good, though there is still the problem of incomplete coverage in the top ~half of the output s2d image, similar to what I showed originally.

Regarding the level 3 combination, I think the problem is that the data WCS are always calculated with respect to the center of the aperture, whether that's a fixed slit or MSA shutter.  The source RA and DEC will only line up exactly with the source trace in the data when the source is placed at that point.  When it's at a different position, particularly in the dispersion direction, there will be an offset.  That offset will be larger the further away the source is from the center of the slit, or the more aligned the slit is with the RA or DEC direction.  For example, the 5th attachment shows the location of the source RA and DEC coordinates (purple and red lines, respectively) in the cal 2d spectrum for dither position 1 in the case of PA=0; DEC is reasonably well-aligned with the source trace, but RA is outside the top edge of the 2d spectrum.  The 6th attachment the same case for PA=45; here both RA and DEC are well-aligned, though slightly separated (the separation between the two is almost non-existent for the dither positions at the center of the slit in the dispersion direction).

Not yet sure how to best deal with this.  Perhaps we could adjust the WCS to shift it to the source frame in each exposure?

[stscijgbot-jp]

Comment by James Muzerolle on JIRA:

Note, all inputs and products for the PA=45 case will shortly be in the central store directory.  Let me know if you want to see the PA=0 case as well.

[stscijgbot-jp]

Comment by James Muzerolle on JIRA:

James Davies [X] RA_REF, DEC_REF, and ROLL_REF were all changed to be consistent with the adopted value of PA_APER.  Though I would definitely appreciate someone double-checking to make sure those are correct.  The other keywords (and the rest of the headers) are taken directly from OTB data generated from an observation with a similar instrument configuration.

[stscijgbot-jp]

Comment by James Muzerolle on JIRA:

Just updated the PA=45 products on central store.  The original versions had grating tilt keyword values that were inconsistent with what were used to generate the simulations, which meant that the wavelength scale was way off.  This shouldn't affect the resample_spec issue.

Also, to add more detail on the pointing keyword determinations, here is the process I used:

For each dither position, values of target RA, DEC, V2, and V3 were taken from an APT pointing file (the first two are essentially arbitrary, as long as they are the same for all cases).  These were then input, along with an aperture PA value, into a script to calculate the *_REF keywords, which used the following pysiaf transforms:

NIRSpec = pysiaf.Siaf('NIRSpec')

aper = NIRSpec['NRS_S200A1_SLIT']

v3pa = apa - aper.V3IdlYAngle

attitude_matrix = rotations.attitude(v2_targ, v3_targ, ra_targ, dec_targ, v3pa)

ra_ref, dec_ref = rotations.pointing(attitude_matrix, aper.V2Ref, aper.V3Ref)

The header keywords are then substituted with these new values:

PA_APER = apa

RA_REF = ra_ref

DEC_REF = dec_ref

ROLL_REF = v3pa

[stscijgbot-jp]

Comment by James Davies [X] on JIRA:

Thanks James. That looks correct.  The pipeline doesn't use PA_APER directly at all, but we do use the others for the final V2/V3 to sky transform.  From your pseudo-code above, it looks like you are computing ROLL_REF via:

 

ROLL_REF = PA_APER - V3I_YANG

 

? With the result wrapped at 360 degrees.

Btw, this is a really hard problem to crack, with lots of ways to make mistakes. We spent a lot of effort 4 years ago getting the mirage simulator to be consistent with the jwst pipeline. And there were a bunch of pitfalls. NIRSpec is even more complicated in its optics.

[stscijgbot-jp]

Comment by James Davies [X] on JIRA:

I'm first looking at the simulated data to see if there's any problems with the consistency of the WCS, pointing and aperture keywords.

I'm comparing this dataset with data coming out of OTB, and run through the pipeline, the above relation between ROLL_REF and PA_APER mention in the comment above holds for our test data suites.

If we compute ROLL_REF using the code that DMS processing uses:

In [1]: from jwst.lib.set_telescope_pointing import compute_local_roll

In [2]: from jwst import datamodels

In [3]: from glob import glob

In [4]: files = glob("*_rate.fits")

In [12]: for file in files:
    ...:     with datamodels.open(file) as model:
    ...:         pa_v3 = model.meta.pointing.pa_v3
    ...:         ra_ref = model.meta.wcsinfo.ra_ref
    ...:         dec_ref = model.meta.wcsinfo.dec_ref
    ...:         v2_ref = model.meta.wcsinfo.v2_ref
    ...:         v3_ref = model.meta.wcsinfo.v3_ref
    ...:         model.meta.wcsinfo.roll_ref = compute_local_roll(pa_v3, ra_ref, dec_ref, v2_ref, v3_ref)
    ...:         model.save(model.meta.filename)

it is pretty close and comes out to ~38.218 deg, with slight variation as there's slight variation in the pointing due to the primary dither pattern.

Note here that ROLL_REF is not the slit aperture position angle. And one can't pick arbitrary ones without redoing the simulation, as it means the different nods won't line up.

Another check:

In [28]: c = datamodels.open("spec3_short_asn.json")

In [29]: for model in c:
    ...:     ra = model.slits[0].meta.target.ra
    ...:     dec = model.slits[0].meta.target.dec
    ...:     lam = (model.slits[0].meta.wcsinfo.waverange_start - model.slits[0].meta.wcsinfo.waverange_end) * 1e6
    ...:     print(model.meta.filename, model.slits[0].meta.wcs.backward_transform(ra, dec, lam))
    ...: 
mod_nirspec_002_CLEAR-PRISM_SLIT_fluxcal11_d02_000_491_cts_NRS1_mod_cal.fits (339.9251513181373, 31.55196512040061)
mod_nirspec_002_CLEAR-PRISM_SLIT_fluxcal11_d05_000_491_cts_NRS1_mod_cal.fits (341.27623165087607, 21.218948380959773)
mod_nirspec_002_CLEAR-PRISM_SLIT_fluxcal11_d08_000_491_cts_NRS1_mod_cal.fits (342.3638673259554, 12.882648091211195)

The WCS suggests a ~10 pixel offset up and down the slit for the nods. Checking these against the images confirms that the WCS correctly describes where the source is in the 2D dispersed slit image.

[stscijgbot-jp]

Comment by James Davies [X] on JIRA:

If I use the ROLL_REF angle computed above, which is almost identical as that in the original simulated images described above, and run spec2 on the images without background subtraction, and then resample only the 3 primary nods, turning off outlier_detection, I get something that looks like the following

!Screen Shot 2021-07-09 at 4.37.38 PM.png|thumbnail!

where the combined _s2d image is at the top, and each of the input _cal nods are shown below. Log stretch.

Note that the whole output frame is offset in the y direction slightly by about ~10 pixels, a known issue, see JP-2027.

So, at this point, it looks like resample_spec is doing roughly what it is being asked to do in terms of putting pixels in the correct location, each nod relative to each other.

[stscijgbot-jp]

Comment by James Muzerolle on JIRA:

James Davies [X] thanks for delving into this!

I disagree with your statement above that we have to use the same roll angle.  The simulations are generated using only the V2/V3 values for each nod position, which are the same regardless of the actual position angle of the slit on the sky.  Since the RA and Dec values per pixel are calculated self-consistently with the above transforms, we should be able to specify any roll angle in the header as long as we re-run assign_wcs.  This is only true for single-slit modes such as FS and IFU; MOS simulations can only use one position angle because it affects which shutters need to be opened.

I think your results using the original ROLL_REF track with what I've seen.  However, there seems to be a problem when the position angle is aligned closely with either the RA or Dec directions, such that there is a tiny range in either value among all pixels in the unrectified spectrum.

[stscijgbot-jp]

Comment by James Davies [X] on JIRA:

I wanted to finally add that for the non-primary nods (1 and 3 are almost the same position as the primary nod 2), I think the WCS keywords look correct. If I run the data through the pipeline and check the WCS object, things look consistent in terms of where the WCS says the source lands on the detector, and the WCS round-trips correctly. I.e., if I do the same exercise as above but for all nods, I get:

In [17]: for file in files:
    ...:     with datamodels.open(file) as model:
    ...:         ra = model.slits[0].meta.target.ra
    ...:         dec = model.slits[0].meta.target.dec
    ...:         lam = (model.slits[0].meta.wcsinfo.waverange_start - model.slits[0].meta.wcsinfo.waverange_end) * 1e6
    ...:         print(model.meta.filename, model.slits[0].meta.wcs.backward_transform(ra, dec, lam))
    ...: 
mod_nirspec_002_CLEAR-PRISM_SLIT_fluxcal11_d01_000_491_cts_NRS1_mod_cal.fits (339.66412384806586, 31.54967569825135)
mod_nirspec_002_CLEAR-PRISM_SLIT_fluxcal11_d02_000_491_cts_NRS1_mod_cal.fits (339.9251513181373, 31.55196512040061)
mod_nirspec_002_CLEAR-PRISM_SLIT_fluxcal11_d03_000_491_cts_NRS1_mod_cal.fits (340.1862689921336, 31.55431788511396)
mod_nirspec_002_CLEAR-PRISM_SLIT_fluxcal11_d04_000_491_cts_NRS1_mod_cal.fits (341.01508345537695, 21.216727642439082)
mod_nirspec_002_CLEAR-PRISM_SLIT_fluxcal11_d05_000_491_cts_NRS1_mod_cal.fits (341.27623165087607, 21.218948380959773)
mod_nirspec_002_CLEAR-PRISM_SLIT_fluxcal11_d06_000_491_cts_NRS1_mod_cal.fits (341.53730368833567, 21.221239397171303)
mod_nirspec_002_CLEAR-PRISM_SLIT_fluxcal11_d07_000_491_cts_NRS1_mod_cal.fits (342.10275885258943, 12.88035547632694)
mod_nirspec_002_CLEAR-PRISM_SLIT_fluxcal11_d08_000_491_cts_NRS1_mod_cal.fits (342.3638673259554, 12.882648091211195)
mod_nirspec_002_CLEAR-PRISM_SLIT_fluxcal11_d09_000_491_cts_NRS1_mod_cal.fits (342.62501083188954, 12.885003198630557)

And you can see that nods 2, 5 and 8 are offset by the normal primary amount from each other, but the others are basically the same nod position, and I assume that is a subpixel dither of some sort?

That said, when I try to combine nods 1, 2 and 3 only into a single resampled output, I get garbage. They should be pretty much on top of each other, but they are not.

!nods123.png|thumbnail!

So this is where the bug seems to be. Will have to investigate further.

[stscijgbot-jp]

Comment by James Davies [X] on JIRA:

Hi James,

Yes, you're correct that you don't have to use the same roll angle as this is a single point source in a single slit, but if you decide to change ROLL_REF, then you also have to change RA_REF and DEC_REF.

I.e. at a different roll angle, nodding up and down the slit will be putting the aperture reference point (V2_REF, V3_REF) at different positions on the sky (RA_REF, DEC_REF). So if you have a specific RA_REF, DEC_REF for your nods that are produced by the simulator, then that implies a unique ROLL_REF as well, which can be computed by the compute_local_roll() function above.

Does this make sense?

James

[stscijgbot-jp]

Comment by James Davies [X] on JIRA:

Adding some more information here. As mentioned to above, it seems that when I use just the 3 primary nods in this dataset (2, 5, and 8), I get good results in the resampled output, as shown. Same if I separately combine nods 1, 4, 7 and 3, 6, 9. These all combine quite nicely, separately, and if I blink the 3 resampled outputs, they look almost identical. Which is good! But if I start mixing and matching, then things go haywire.

Further, the output WCS in the resampled combined spectrum is always identical, whether I use only the 3 primary nods which look good, or all 9 nods, which look bad. So this implies that the problem isn't in the creation of the output resampled spatial/spectral WCS, but is somehow a problem with either the input WCS or ... something else?

Since there's a known issue with wavecorr, I turned it off when creating the _cal files, and it didn't have any effect on the combined resampled spectrum. Same problem.

Anyway, more investigation needed.

[stscijgbot-jp]

Comment by James Davies [X] on JIRA:

I think I understand the problem here now.

The 3 sets of primary nods are nodded from each other along the spectral axis in the dispersion direction. The slit mapped onto the sky is looking at a different "line" on the sky in these 3 dispersion-direction nods. I.e. slit in the 3 dispersion-direction nods is 3 parallel lines on sky, with no intersecting RA/Dec pairs. Because of this, depending on the orientation, when one tries to figure out where a pixel with a particular RA/Dec pair maps onto the same RA/Dec in the other image, one gets an offset, as the RA/Dec pair maps to Y pixel position, and one has to chose Y --> RA or Y --> Dec. Depending on orientation and which one you choose, you get an offset up or down the slit.

But the real problem is that between the 3 dispersion-direction nods, there are no RA/Dec pairs that are common. Spatial source-in-slit position, yes, RA/Dec, no.

[stscijgbot-jp]

Comment by James Muzerolle on JIRA:

James Davies [X] exactly!  I had trouble articulating this yesterday. I'm working on an experiment where I adjust RA/DEC_REF for each spectral dither with XOFFSET != 0 to the value it would be at the central position. This removes any "extra" spatial information in those dithers, but we don't care about that since the only point is to improve spectral sampling. I think this can also be generalized to MOS using the source_xpos parameter. I'll let you know what I find next week.

[stscijgbot-jp]

Comment by James Davies [X] on JIRA:

Thanks James. Yeah, we're on the same page now. It took me a while to get there, but I understand now.

The reason this happens is that we have the WCS which maps

x, y --> ra, dec, wavelength

and of course the inverse. So when resampling an input image to the output rectified image, we need to compute the mapping of input x, y pixel coordinates to output x, y pixel coordinates. But in the output frame y maps to both RA and Dec, so you have to pick either RA or Dec to map back to y in the inverse transform.

The code currently uses declination

https://github.com/spacetelescope/jwst/blob/bb0f4af73fafd27eaee1f47a56a16a5693ae88ec/jwst/resample/resample_spec.py#L247-L259

unless the declination coordinate is close to parallel to the slit, and then it picks RA.

So instead, I think we need to use RA/Dec to just compute the spatial offsets of the spatial dithers in the slit, and then stack everything on top of each other. In the baseline pipeline, at least currently, there's no way to increase the sampling of the rectified wavelength only, based on the spectral dithers, but that might be something to add in the future.

[stscijgbot-jp]

Comment by James Muzerolle on JIRA:

I tried the following "kludge": adjusting for the lack of overlapping coordinates among the spectral dithers by calculated the REF coordinates at the dither position instead of the slit center. First, I saved the value of the "x_offset" keyword from the original cal file for each dither. Then I used the coordinate transforms to calculate the V2,V3 and RA,DEC coordinates at that offset (as opposed to the default at slit center), and replaced the RA/DEC_REF keywords with those values. Here's the code I used (looped over all 3 spectral dithers at the same primary dither position, where the input file was the original rate file:

root = file[:-9] outfile = root+'offset_rate.fits' # filename for new rate files with adjusted REF coords calfile = root+'cal.fits' # name of original cal file, needed for wcs info cal = datamodels.open(calfile) slit = cal.slits[0]

xoff = slit.meta.dither.x_offset yoff = slit.meta.dither.y_offset

(instead of 0,0, make the ref point at xoff,0 and convert to RA,DEC)

idl2v2v3 = trmodels.IdealToV2V3(slit.meta.wcsinfo.v3yangle,slit.meta.wcsinfo.v2_ref,slit.meta.wcsinfo.v3_ref,slit.meta.wcsinfo.vparity) v2off, v3off = idl2v2v3(xoff,0.) v2v32world = slit.meta.wcs.get_transform('v2v3','world') raoff, decoff, foo = v2v32world(v2off, v3off, slit.meta.wcsinfo.roll_ref)

 

(set RA_REF, DEC_REF in the rate file to the offset values)

hdul = fits.open(file) hdul['SCI'].header['RA_REF'] = raoff hdul['SCI'].header['DEC_REF'] = decoff hdul.writeto(outfile) hdul.close()

After reprocessing through spec2 and spec3, I get a result in which the source trace overlaps pretty well:

!Screen Shot 2021-08-25 at 3.46.18 PM.png|thumbnail!

I'm not sure if this is necessarily the optimal method to fix the problem, but it appears to do the job, at least for this roll angle. I haven't yet tried the more extreme case of PA=0.

[stscijgbot-jp]

Comment by James Muzerolle on JIRA:

As a side note, the resampled spectra here appear to have a slight tilt in which the source trace is not perfectly aligned along the wavelength axis (compare with the purple line, which is the value of the source RA position). I'm not sure if that has to do with the narrow range of input sky coordinate values? The effect appears in all the s2d products for this particular simulation, both spec2 and spec3.

[stscijgbot-jp]

Comment by Howard Bushouse on JIRA:

James Muzerolle Mihai and I are starting to look into this issue again. First thing we want to do is to get a fresh reprocessing of all the simulated images (starting from level-2b) so that we know all of the latest software has been applied. Are the rate files that you added to the directory mentioned above, back in May 2021, still the latest/greatest inputs to use? Or have you done any fiddling with any of the WCS-related values since then?

[stscijgbot-jp]

Comment by James Muzerolle on JIRA:

Howard Bushouse the rate files on central store should have the correct header information. There are two sets - exposures with the original aperture position angle (~180 deg) and a subsequent set I made with apa = 45.

[stscijgbot-jp]

Comment by Howard Bushouse on JIRA:

I've reprocessed all of the data files with the latest Cal software and ref files. All files are in ███████████████████████████████████████████ While all of the individual s2d (resampled) images of the 9 separate spectra look fine, the combined results created by calwebb_spec3 still look odd, for both sets of data (the two different position angles). I'm starting simple and only combining the 3 exposures with different spatial nod positions (i.e. exposures 2, 5, and 8); no dithers in the spectral direction. I'm seeing results qualitatively similar to what James Davies showed (above) in July 2021, in that I see what appears to be a single, combined spectrum in the output level-3 s2d images, but there are oddities in the details of each.

mod_nirspec_002_CLEAR-PRISM_SLIT_fluxcal11_s00001_nirspec_clear-prism-s200a1-subs200a1_s2d.fits: This is the combination of exposures d02, d05, d08 (just the primary spatial offsets in the cross-dispersion direction). d02 has the spectrum near the top of the 2D cutout (image line 31 of the cal file), d05 has the spectrum near the center (image line 21), and d08 has the spectrum near the bottom (image line 12). Naively, I would expect the output s2d frame to cover the union of all the input images and hence I would expect the spectrum in the s2d image to be nearly centered, with a band of above +-10 pixels above and below center having contributions from all 3 inputs, and then as you go further towards the top and bottom the contributions would decrease to 2 and 1 inputs. But the spectrum in this combined image is centered around image line 19, which is well below center. So I'm not sure why the output frame was defined in this way, such that more spatial territory is included above the spectrum than below.

Second, when you look at the WHT and CON images, the results are even more peculiar. For example, the CON image assigns values of 1, 2, and 4 (powers of 2) for each of the input images. So I would expect that all of the pixels near the center of the combined spectral trace would have values of 1+2+4=7, but in fact they all have values of only 6 (contributions from only images 2 and 3). There are very few places in the entire CON image that include the 1 value, i.e. input image 1 contributed almost nothing to the combined product. Why? The WHT image shows the same behavior (most pixels have values ~1 and ~2, indicating contributions from only 1 or 2 of the inputs, with almost no pixels having a value of ~3).

In addition to these anomalies, the whole combined spectrum just looks very "choppy", especially along the wings of the trace (a few image lines above and below the trace center). So resampling was not done well in these regions.

mod_pa45_nirspec_002_CLEAR-PRISM_SLIT_fluxcal11_s00001_nirspec_clear-prism-s200a1-subs200a1_s2d.fits: This is the same combination of exposures d02, d05, and d08, but the ones with a simulated PA of ~45 degrees. The combined result is qualitatively similar to the first one discussed above, in that the combined spectral trace is again offset towards the bottom of the combined image (instead of centered), and the WHT and CON arrays show the same behavior of only having contributions from 2 of the 3 input images in most places. The overall resampled spectrum looks somewhat cleaner than the one above (no choppiness along the wings of the trace), but still has the issues of not using all the inputs.

Like James and James have done in the past, I verified that the WCS info in the input CAL images looks reasonable, in that pixels in the center of the spectral trace produce RA/Dec values consistent with the target coords listed in the header and RA/Dec is pretty much constant along image rows (dispersion direction). The WCS in the output combined s2d images looks similarly consistent. So the big question remains as to why the 3 inputs don't seem to be getting combined properly.

[stscijgbot-jp]

Comment by Howard Bushouse on JIRA:

Mihai Cara here's a list and explanation of what's in the directory containing the test files.

mod_nirspec_002_CLEAR-PRISM_SLIT_fluxcal11_d0[123456789]_000_491_cts_NRS1_mod_rate.fits: These are the simulated countrate images used as input to the calwebb_spec2 pipeline. These only need to be used if you want to try changing things related to the WCS, because changes would then require reprocessing through calwebb_spec2 to compute a new gwcs object. The 9 images represent 9 different dither positions of the source within the spectrograph slit: 1-3 have the source near the top of the slit (in the cross-dispersion direction) and shift the source to 3 different positions along the dispersion direction; 4-6 have the source near the center of the slit, and 7-9 have the source near the bottom of the slit.

mod_pa45_nirspec_002_CLEAR-PRISM_SLIT_fluxcal11_d0[123456789]_000_491_cts_NRS1_mod_rate.fits: Same as the 9 images above, but at a different position angle on the sky.

mod_nirspec_002_CLEAR-PRISM_SLIT_fluxcal11_d0[123456789]_000_491_cts_NRS1_mod_cal.fits: The 9 outputs from calwebb_spec2 for the first set of images. The cal files are unrectified (no resampling).

mod_pa45_nirspec_002_CLEAR-PRISM_SLIT_fluxcal11_d0[123456789]_000_491_cts_NRS1_mod_cal.fits: The 9 outputs from calwebb_spec2 for the second set of images.

mod_nirspec_002_CLEAR-PRISM_SLIT_fluxcal11_d0[123456789]_000_491_cts_NRS1_mod_s2d.fits: mod_pa45_nirspec_002_CLEAR-PRISM_SLIT_fluxcal11_d0[123456789]_000_491_cts_NRS1_mod_s2d.fits: The 9 resampled images from calwebb_spec2 for the first and second sets of images. Each s2d file is the resampled result of a single cal file.

spec3_short_asn.json: The ASN file used to run calwebb_spec3 on the first set of images. Contains only the 3 images at the primary nod positions (no images with shifts in the spectral direction). spec3_short_pa45_asn.json: The ASN file used to run calwebb_spec3 on the second set of images.

To run calwebb_spec3 from the command line: "strun calwebb_spec3 spec3_short_asn.json"

mod_nirspec_002_CLEAR-PRISM_SLIT_fluxcal11_s00001_nirspec_clear-prism-s200a1-subs200a1_s2d.fits: mod_pa45_nirspec_002_CLEAR-PRISM_SLIT_fluxcal11_s00001_nirspec_clear-prism-s200a1-subs200a1_s2d.fits: The combined, resampled images for the first and second sets of 3 images.

[stscijgbot-jp]

Comment by Mihai Cara on JIRA:

Howard Bushouse James Muzerolle I uploaded "resample_s2d.png" that shows how the spectrum for "spec3_short_asn.json" now looks.

I also forgot to report that there is a PR that hopefully will fix this issue: #6747

Please clone that branch and try it out. Thanks!

[stscijgbot-jp]

Comment by Howard Bushouse on JIRA:

The new results, using the updated code in https://github.com/spacetelescope/jwst/pull/6747, look REALLY good. The combined spectrum is does not have all the choppy features that the previous version had, it is ruler flat along image rows, and the WHT and CON images show the expected patterns of 1, 2, or all 3 inputs being used in the regions expected (all 3 along the spectral trace, dropping to only 2 or 1 as you move further away from the trace). And the results are good for both sets using different position angles.

The only minor nit is that if the nod that has the source at the top of the slit is listed first in the ASN file, that's the image that gets used as the reference for setting up the output WCS frame, so that the combined image also has the spectrum near the top, and hence not giving any results for the extended regions of background that appear above the trace from nods 2 and 3. That's fixable by simply editing the ASN file to list the middle image first, in which case it gets used as the reference and the resulting combined image then has the spectrum centered, with equal areas of background available above and below the trace. Is there perhaps some way to create the output WCS frame as the union of all the inputs, rather than using the first image as the reference?

[stscijgbot-jp]

Comment by Howard Bushouse on JIRA:

The results from using only the 3 primary nod positions versus using all 9 dithers (including dithers along the dispersion direction) are qualitatively identical, with only small variations in actual pixel values along the spectrum. So it's hard to tell whether the shifts in wavelength are being handled properly. Would perhaps be useful to have data with some obvious spectral features (sharp emission/absorption lines) to analyze this better. Guess we need to wait for the first real observations of a nice planetary nebula or something!

[stscijgbot-jp]

Comment by Mihai Cara on JIRA:

??... if the nod that has the source at the top of the slit is listed first in the ASN file, that's the image that gets used as the reference for setting up the output WCS frame, so that the combined image also has the spectrum near the top, and hence not giving any results for the extended regions of background that appear above the trace from nods 2 and 3. ... Is there perhaps some way to create the output WCS frame as the union of all the inputs, rather than using the first image as the reference???

That should not be happening so I would consider that a bug. Output WCS should be a union. I will look into it. Thanks for noticing!

With regard to the dithering in the spectral direction...

??So it's hard to tell whether the shifts in wavelength are being handled properly.??

I would say dithers in the spectral direction should be handled "properly", in my understanding(!). What I do worry about is that your expectations might be different from what I understood should happen. For example, dithering in itself would not double the number of sampling wavelengths in the output image just like it would not affect output image size in the imaging case. Sampling frequencies are determined from the reference image and so pixel "size (or length)" stays the same as input. That is, dithering would not result in an improved sampling/resolution. Choosing a finer sampling is complicated in itself because of the tabular nature of input WCS along the spectral axis. But even more complicated is the fact that we also need to reduce drop size if you actually want to benefit from a finer sampling. And that should be done only for one side of the pixel (spectral side) and that is not supported in the drizzle, I believe.

[stscijgbot-jp]

Comment by James Muzerolle on JIRA:

Mihai Cara Howard Bushouse I've finally had a chance to try out the new code. As Howard has already seen, the results for the FS simulations look great! I also processed some MOS simulations, and ran into two (possibly related) issues in the Spec2 pipeline. In one case where I extracted just a single target, the exposure-level s2d product appears to be truncated, with the NRS1 portion missing ~15% at the red end of the spectrum, and the (short, ~160 pixels long in the unrectified spectrum) NRS2 portion missing ~90%. This might be the same problem that Howard mentioned in github. In another case trying to extract all targets, the pipeline went through multiple targets successfully in resample_spec before crashing at one with an error saying "Not enough data to construct output WCS." (full traceback below). I haven't verified, but I suspect this particular spectrum may be one that has a width of a small number of pixels on one of the detectors. 

2022-03-07 12:11:30,073 - stpipe.Spec2Pipeline.resample_spec - INFO - Update S_REGION to POLYGON ICRS 215.052087247 52.914563330 215.052974849 52.914563330 215.052974849 52.924948972 215.052087247 52.924948972 Traceback (most recent call last): File "/Users/muzerol/NIRSpec/Pipeline/Testing/Clone/jwst/jwst/pipeline/calwebb_spec2.py", line 115, in process result = self.process_exposure_product( File "/Users/muzerol/NIRSpec/Pipeline/Testing/Clone/jwst/jwst/pipeline/calwebb_spec2.py", line 265, in process_exposure_product resampled = self.resample_spec(calibrated) File "/Users/muzerol/anaconda3/envs/jwst_test/lib/python3.10/site-packages/stpipe/step.py", line 430, in run step_result = self.process(*args) File "/Users/muzerol/NIRSpec/Pipeline/Testing/Clone/jwst/jwst/resample/resample_spec_step.py", line 67, in process result = self._process_multislit(input_models) File "/Users/muzerol/NIRSpec/Pipeline/Testing/Clone/jwst/jwst/resample/resample_spec_step.py", line 103, in _process_multislit resamp = resample_spec.ResampleSpecData(container, *self.drizpars) File "/Users/muzerol/NIRSpec/Pipeline/Testing/Clone/jwst/jwst/resample/resample_spec.py", line 78, in init* self.output_wcs = self.build_nirspec_output_wcs() File "/Users/muzerol/NIRSpec/Pipeline/Testing/Clone/jwst/jwst/resample/resample_spec.py", line 117, in build_nirspec_output_wcs raise ValueError("Not enough data to construct output WCS.") ValueError: Not enough data to construct output WCS.

[stscijgbot-jp]

Comment by Howard Bushouse on JIRA:

Yes, we've seen those same kinds of problems when trying to run some of our regression tests on the PR code. See the recent comments that've been added to the PR itself at #6747  So Mihai Cara is working on additional fixes and updates.

[stscijgbot-jp]

Comment by Mihai Cara on JIRA:

As Howard Bushouse mentioned previously (and verified by me) this code has numerous issues. I am working on them.

With regard to:

??I haven't verified, but I suspect this particular spectrum may be one that has a width of a small number of pixels on one of the detectors.??

That is an issue indeed and I would like to know how to proceed when a slice has no valid data from which to build a valid WCS.

There is another issue that Howard Bushouse brought up above: this code currently takes the wavelengths from the reference image to serve as the sampling wavelength for the output image. When there is dithering along the spectral axis, how should I handle that situation? Based on Howard's comment is seems that you would like the spectrum range to be extended to include the entire range from all input images. That can be done. My biggest question is about how to handle the overlap regions. For example, Let's say one input image samples spectrum at 1, 2, 3, 4, 5, 6 (arbitrary units) and the second image samples at 0.1, 1.1, 2.1, 3.1, 4.1, 5.1, 6.1.  How should be the output grid created?

A: 0.1, 1, 2, 3, 4, 5, 6, 6.1 ?

B: 0.1, 1, 2, 3, 4, 5, 6?

C: 0.1, 1, 1.1, 2, 2.1, 3, 3.1, ...?

[stscijgbot-jp]

Comment by James Muzerolle on JIRA:

Apologies for the redundant test info, I had not read through the comments on the PR.

Mihai Cara regarding the short spectra issue - is there a minimum number of pixels needed to get a solution? I guess there might not be a unique number since it may depend on the range for a given target. The extract_2d step removes any slit from consideration if/when no part of the spectral trace is projected onto one of the detectors, but otherwise creates a cutout regardless of width as long as it's > 0; if need be we could always increase that threshold. However, in the case I mentioned that caused a crash, there are valid WCS values over an area of 166x27 pixels (corresponding to the dispersion and spatial directions, respectively). The original version of the code did generate a resampled product with similar size.

The output grid should contain the full range of sky & wavelength values encompassed by the ensemble of exposures.  For the spectral sampling, I think previously it was based on the minimum of the input data? However, thinking about it some more, it might be better to adopt a similar approach to cube_build, which uses a characteristic value per disperser as specified in a parameter reference file. Right now those values are based I think on the minimum achievable value over the entire wavelength range of the disperser, but will be optimized once we get data with unresolved/unblended lines later in Commissioning.

[stscijgbot-jp]

Comment by Mihai Cara on JIRA:

??For the spectral sampling, I think previously it was based on the minimum of the input data???

This probably addresses the issue of range but not the sampling itself.

[stscijgbot-jp]

Comment by Mihai Cara on JIRA:

I also would like to ask whether resample step must be able to support drizzle-combining images with different orientations of the slit on the sky (aperture PA?)?

[stscijgbot-jp]

Comment by Mihai Cara on JIRA:

My latest commit should have fixed most issues. Please give it another try.

[stscijgbot-jp]

Comment by Howard Bushouse on JIRA:

Fixed by #6747

[stscijgbot-jp]

Comment by Mihai Cara on JIRA:

... and by #6780

[stscijgbot-jp]

Comment by Howard Bushouse on JIRA:

Temporarily resetting the ticket status back to "In Progress", because we discovered some problems in the code updates in the recent PR's and are again doing developer-level testing before declaring this ready for INS to test. Should be ready to go soon.

[stscijgbot-jp]

Comment by Howard Bushouse on JIRA:

Final fixes and updates in #6780 and #6802

[stscijgbot-jp]

Comment by James Muzerolle on JIRA:

Tested on FS and MOS data (one disperser each), both level 2 and 3. The results look fantastic, thanks for all your hard work!