Skip to content
New issue

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

Sign up for GitHub

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

Already on GitHub? Sign in to your account

JP-3244: Provide auto-centering option for slit-like data containing point sources #7763

Draft
wants to merge 8 commits into
base: main
Choose a base branch
from
Draft

JP-3244: Provide auto-centering option for slit-like data containing point sources #7763

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
8 commits
Select commit Hold shift + click to select a range
d4ac776
introduce slit_autocen
tapastro Jul 25, 2023
e7280d0
push limit changes to p_src
tapastro Jul 25, 2023
8515cdb
amend curve_fit bounds for MIRI behavior
tapastro Jul 25, 2023
a299567
enclose curve_fit in try/except to avoid returning errors
tapastro Jul 25, 2023
60e7271
norm each row before median to avoid spectral slope weighting fit
tapastro Aug 1, 2023
b834ec7
edits to allow extraction from cals
tapastro Aug 4, 2023
779e3e6
fix axis slice for dispersion sum
tapastro Aug 9, 2023
8c7f6c1
use ranges
tapastro Aug 14, 2023
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
96 changes: 86 additions & 10 deletions jwst/extract_1d/extract.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -8,6 +8,7 @@
from typing import Union, Tuple, NamedTuple, List
from astropy.modeling import polynomial
from astropy.io import fits
from scipy.optimize import curve_fit
from gwcs import WCS
from stdatamodels import DataModel
from stdatamodels.properties import ObjectNode
Expand Down Expand Up @@ -88,6 +89,10 @@
EXACT = "exact match"


def gaussian(x, a, b, c):
return c / a * np.exp(-1. * (x - b)**2 / (2. * a * a))


class Aperture(NamedTuple): # When python 3.6 is no longer supported, consider converting to DataClass
xstart: Union[int, float]
xstop: Union[int, float]
Expand Down Expand Up @@ -220,7 +225,8 @@ def get_extract_parameters(
smoothing_length: Union[int, None],
bkg_fit: str,
bkg_order: Union[int, None],
use_source_posn: Union[bool, None]
use_source_posn: Union[bool, None],
slit_autocen: Union[bool, None]
) -> dict:
"""Get extract1d reference file values.

Expand Down Expand Up @@ -283,6 +289,10 @@ def get_extract_parameters(
If None, the value specified in `ref_dict` will be used, or it will
be set to True if not found in `ref_dict`.

slit_autocen : bool
Switch to turn on auto-centering for point source spectral extraction
in slit-like modes. Default is True.

Returns
-------
extract_params : dict
Expand Down Expand Up @@ -310,6 +320,7 @@ def get_extract_parameters(
extract_params['bkg_fit'] = None # because no background sub.
extract_params['bkg_order'] = 0 # because no background sub.
extract_params['subtract_background'] = False
extract_params['slit_autocen'] = slit_autocen

if use_source_posn is None:
extract_params['use_source_posn'] = False
Expand Down Expand Up @@ -386,14 +397,15 @@ def get_extract_parameters(
if use_source_posn is None: # no value set on command line
if use_source_posn_aper is None: # no value set in ref file
# Use a suitable default
if meta.exposure.type in ['MIR_LRS-FIXEDSLIT', 'MIR_MRS', 'NRS_FIXEDSLIT', 'NRS_IFU', 'NRS_MSASPEC']:
if meta.exposure.type in ['MIR_MRS', 'NRS_IFU']:
use_source_posn = True
log.info(f"Turning on source position correction for exp_type = {meta.exposure.type}")
else:
use_source_posn = False
else: # use the value from the ref file
use_source_posn = use_source_posn_aper
extract_params['use_source_posn'] = use_source_posn
extract_params['slit_autocen'] = slit_autocen

extract_params['extract_width'] = aper.get('extract_width')
extract_params['position_correction'] = 0 # default value
Expand Down Expand Up @@ -1056,7 +1068,8 @@ def __init__(
subtract_background: Union[bool, None] = None,
use_source_posn: Union[bool, None] = None,
match: Union[str, None] = None,
ref_file_type: Union[str, None] = None
ref_file_type: Union[str, None] = None,
slit_autocen: Union[bool, None] = None
):
"""
Parameters
Expand Down Expand Up @@ -1198,6 +1211,7 @@ def __init__(
self.use_source_posn = use_source_posn
self.position_correction = position_correction
self.subtract_background = subtract_background
self.slit_autocen = slit_autocen

self.wcs = None # initial value

Expand Down Expand Up @@ -1417,7 +1431,7 @@ def locn_from_wcs(
# Width (height) in the cross-dispersion direction, from the start of the 2-D cutout (or of the full image)
# to the upper limit of the bounding box.
# This may be smaller than the full width of the image, but it's all we need to consider.
xd_width = int(round(bb[1][1])) # must be an int
xd_width = int(round(bb[1][1])+0.0001) # must be an int
middle = int((bb[0][0] + bb[0][1]) / 2.) # Middle of the bounding_box in the dispersion direction.
x = np.empty(xd_width, dtype=np.float64)
x[:] = float(middle)
Expand Down Expand Up @@ -2524,6 +2538,7 @@ def run_extract1d(
subtract_background: Union[bool, None],
use_source_posn: Union[bool, None],
center_xy: Union[float, None],
slit_autocen: Union[bool, None],
ifu_autocen: Union[bool, None],
ifu_rfcorr: Union[bool, None],
was_source_model: bool = False,
Expand Down Expand Up @@ -2578,6 +2593,10 @@ def run_extract1d(
done by the step. Two values, in x,y order, are used for extraction
from IFU cubes. Default is None.

slit_autocen : bool
Switch to turn on auto-centering for point source spectral extraction
in slit-like modes. Default is True.

ifu_autocen : bool
Switch to turn on auto-centering for point source spectral extraction
in IFU mode. Default is False.
Expand Down Expand Up @@ -2632,6 +2651,7 @@ def run_extract1d(
subtract_background,
use_source_posn,
center_xy,
slit_autocen,
ifu_autocen,
ifu_rfcorr,
was_source_model,
Expand Down Expand Up @@ -2694,6 +2714,7 @@ def do_extract1d(
subtract_background: Union[int, None] = None,
use_source_posn: Union[bool, None] = None,
center_xy: Union[int, None] = None,
slit_autocen: Union[bool, None] = None,
ifu_autocen: Union[bool, None] = None,
ifu_rfcorr: Union[bool, None] = None,
was_source_model: bool = False
Expand Down Expand Up @@ -2758,6 +2779,10 @@ def do_extract1d(
Two values, in x,y order, are used for extraction from IFU cubes.
Default is None.

slit_autocen : bool
Switch to turn on auto-centering for point source spectral extraction
in slit-like modes. Default is True.

ifu_autocen : bool
Switch to turn on auto-centering for point source spectral extraction
in IFU mode. Default is False.
Expand Down Expand Up @@ -2869,7 +2894,7 @@ def do_extract1d(
smoothing_length, bkg_fit, bkg_order, use_source_posn,
prev_offset, exp_type, subtract_background, input_model,
output_model, apcorr_ref_model, log_increment,
is_multiple_slits
is_multiple_slits, slit_autocen
)
except ContinueError:
continue
Expand Down Expand Up @@ -2918,7 +2943,7 @@ def do_extract1d(
smoothing_length, bkg_fit, bkg_order, use_source_posn,
prev_offset, exp_type, subtract_background, input_model,
output_model, apcorr_ref_model, log_increment,
is_multiple_slits
is_multiple_slits, slit_autocen
)
except ContinueError:
continue
Expand Down Expand Up @@ -2959,7 +2984,7 @@ def do_extract1d(
smoothing_length, bkg_fit, bkg_order, use_source_posn,
prev_offset, exp_type, subtract_background, input_model,
output_model, apcorr_ref_model, log_increment,
is_multiple_slits
is_multiple_slits, slit_autocen
)
except ContinueError:
continue
Expand Down Expand Up @@ -3451,6 +3476,54 @@ def extract_one_slit(
extract_model.log_extraction_parameters()
extract_model.assign_polynomial_limits()

if extract_model.slit_autocen is True:
good_pixels = np.where(data == data)
y_range = [np.min(good_pixels[0]), np.max(good_pixels[0]) + 1]
x_range = [np.min(good_pixels[1]), np.max(good_pixels[1]) + 1]
extract_model.xstart = x_range[0]
extract_model.xstop = x_range[1]
extract_model.ystart = y_range[0]
extract_model.ystop = y_range[1]

absdata = data.copy()
absdata = absdata[y_range[0]: y_range[1], x_range[0]: x_range[1]]

absdata[absdata < 0] = 0

normed = absdata / np.nanmax(absdata, axis=extract_params['dispaxis'] - 1, keepdims=True)

y_fit = np.nanmedian(normed, axis=int(-1 * extract_params['dispaxis']))
x_fit = np.arange(len(y_fit))
y_peak = np.argmax(y_fit)
# Provide bounds to fit to avoid wasting time in region where a fit is undesirable
try:
fit_params, _ = curve_fit(gaussian, x_fit, y_fit, bounds=(np.array([0, y_peak - 3, 0.5]),
np.array([8, y_peak + 3, 10]))
)
fit_width, fit_center, _ = fit_params
n_sig = 3.
if extract_model.dispaxis == HORIZONTAL:
extract_model.ystart = max(0, fit_center - (n_sig * fit_width)) + y_range[0]
extract_model.ystop = min(data.shape[-2] - 1, fit_center + (n_sig * fit_width)) + y_range[0]
extract_model.p_src = [[create_poly([extract_model.ystart]),
create_poly([extract_model.ystop])]]
log.debug(f"Slit auto-centering fit_center: {fit_center + y_range[0]}"
f"and fit_width: {fit_width}")

else:
extract_model.xstart = max(0, fit_center - (n_sig * fit_width)) + x_range[0]
extract_model.xstop = min(data.shape[-1] - 1, fit_center + (n_sig * fit_width)) + x_range[0]
extract_model.p_src = [[create_poly([extract_model.xstart]),
create_poly([extract_model.xstop])]]

log.debug(f"Slit auto-centering fit_center: {fit_center + x_range[0]}"
f"and fit_width: {fit_width}")

log.info("Overriding cross-dispersion bounds with slit auto-centering.")
except (ValueError, RuntimeError) as e:
log.warning("Auto-centering failed; reverting to default behavior.")
log.debug(f"Auto-centering error: {e}")

# Log the extraction limits being used
if integ < 1:
if extract_model.src_coeff is not None:
Expand Down Expand Up @@ -3599,7 +3672,8 @@ def create_extraction(extract_ref_dict,
output_model,
apcorr_ref_model,
log_increment,
is_multiple_slits):
is_multiple_slits,
slit_autocen):
if slit is None:
meta_source = input_model
else:
Expand Down Expand Up @@ -3648,7 +3722,8 @@ def create_extraction(extract_ref_dict,
# Turn off use_source_posn if the source is not POINT
if source_type != 'POINT':
use_source_posn = False
log.info(f"Setting use_source_posn to False for source type {source_type}")
slit_autocen = False
log.info(f"Setting use_source_posn and slit auto-centering to False for source type {source_type}")

# Turn off use_source_posn if working on non-primary NRS fixed slits
if is_multiple_slits:
Expand All @@ -3674,7 +3749,8 @@ def create_extraction(extract_ref_dict,
smoothing_length,
bkg_fit,
bkg_order,
use_source_posn
use_source_posn,
slit_autocen
)

if subtract_background is not None:
Expand Down
5 changes: 5 additions & 0 deletions jwst/extract_1d/extract_1d_step.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -147,6 +147,7 @@ class Extract1dStep(Step):
use_source_posn = boolean(default=None) # use source coords to center extractions?
center_xy = float_list(min=2, max=2, default=None) # IFU extraction x/y center
apply_apcorr = boolean(default=True) # apply aperture corrections?
slit_autocen = boolean(default=True) # Auto source centering for slit-like point source data.
ifu_autocen = boolean(default=False) # Auto source centering for IFU point source data.
ifu_rfcorr = boolean(default=False) # Apply 1d residual fringe correction
soss_atoca = boolean(default=True) # use ATOCA algorithm
Expand Down Expand Up @@ -258,6 +259,7 @@ def process(self, input):
self.subtract_background,
self.use_source_posn,
self.center_xy,
self.slit_autocen,
self.ifu_autocen,
self.ifu_rfcorr,
was_source_model=was_source_model
Expand Down Expand Up @@ -294,6 +296,7 @@ def process(self, input):
self.subtract_background,
self.use_source_posn,
self.center_xy,
self.slit_autocen,
self.ifu_autocen,
self.ifu_rfcorr,
was_source_model=was_source_model,
Expand Down Expand Up @@ -333,6 +336,7 @@ def process(self, input):
self.subtract_background,
self.use_source_posn,
self.center_xy,
self.slit_autocen,
self.ifu_autocen,
self.ifu_rfcorr,
was_source_model=was_source_model,
Expand Down Expand Up @@ -474,6 +478,7 @@ def process(self, input):
self.subtract_background,
self.use_source_posn,
self.center_xy,
self.slit_autocen,
self.ifu_autocen,
self.ifu_rfcorr,
was_source_model=False,
Expand Down