Skip to content

Commit

Permalink
options to precompute a contamination model before median filtering
Browse files Browse the repository at this point in the history
  • Loading branch information
@gbrammer
gbrammer committed Sep 14, 2023
1 parent 05eac6a commit 7d26844
Show file tree
Hide file tree
Showing 2 changed files with 58 additions and 15 deletions.
55 changes: 43 additions & 12 deletions grizli/multifit.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -1186,7 +1186,7 @@ def subtract_sep_background(self, mask_threshold=3, mask_iter=3, bw=256, bh=256,
flt.grism.header['BKGBH'] = bh, 'sep background bh'


def subtract_median_filter(self, filter_size=71, filter_central=10, revert=True, filter_footprint=None, subtract_model=False, second_pass_filtering=False, box_filter_sn=3, box_filter_width=3):
def subtract_median_filter(self, filter_size=71, filter_central=10, revert=True, filter_footprint=None, subtract_model=False, second_pass_filtering=False, box_filter_sn=3, box_filter_width=3, put_model_in_median=False, verbose=True, mask_sn_threshold=None, mask_sn_dilate_iters=5):
"""
Remove a median filter calculated along the dispersion axis
"""
Expand All @@ -1210,16 +1210,13 @@ def subtract_median_filter(self, filter_size=71, filter_central=10, revert=True,
msg += f' filter_size={filter_size} filter_central={filter_central}'
msg += f' [{_filter_name}]'

utils.log_comment(utils.LOGFILE, msg, verbose=True)
utils.log_comment(utils.LOGFILE, msg, verbose=verbose)

sci_i = flt.grism.data['SCI']*1

if revert & ('MED' in flt.grism.data):
sci_i += flt.grism.data['MED']

if subtract_model:
sci_i -= flt.model


err_i = flt.grism.data['ERR']
dq_i = flt.grism.data['DQ']

Expand All @@ -1229,6 +1226,12 @@ def subtract_median_filter(self, filter_size=71, filter_central=10, revert=True,
ivar = 1/err_i**2
ivar[~ok] = 0

if subtract_model:
msg = 'subtract_median_filter: subtract model before median'
utils.log_comment(utils.LOGFILE, msg, verbose=verbose)

sci_i -= flt.model

#filter_sci = np.zeros_like(sci_i)
#sh = sci_i.shape

Expand All @@ -1245,16 +1248,18 @@ def subtract_median_filter(self, filter_size=71, filter_central=10, revert=True,
if second_pass_filtering:
if nbutils is None:
# need numba installed
msg = 'subtract_median_filter: `numba` not found, skip second filter pass.'
utils.log_comment(utils.LOGFILE, msg, verbose=True)
msg = 'subtract_median_filter: `numba` not found, '
msg += 'skip second filter pass.'
utils.log_comment(utils.LOGFILE, msg, verbose=verbose)
else:
# run filter again, but mask pixels that show significant residuals (e.g. strong emission lines)
msg = f'subtract_median_filter: rerun filtering masking '
msg += f' S/N>{box_filter_sn} pixels in residual'
utils.log_comment(utils.LOGFILE, msg, verbose=True)
utils.log_comment(utils.LOGFILE, msg, verbose=verbose)

# first do some binning/box filtering to identify significantly detected lines in individual exposures
box_filter_footprint = np.ones((box_filter_width, box_filter_width),
box_filter_footprint = np.ones((box_filter_width,
box_filter_width),
dtype=int)
box_filter_clean = nd.generic_filter(sci_i-filter_sci,
nbutils.nansum,
Expand All @@ -1272,9 +1277,35 @@ def subtract_median_filter(self, filter_size=71, filter_central=10, revert=True,
footprint=filter_footprint)

filter_sci[~np.isfinite(filter_sci)] = 0

if mask_sn_threshold:
flt.grism.header['MEDSNTH'] = (mask_sn_threshold,
'Median mask threshold')
flt.grism.header['MEDSNIT'] = (mask_sn_dilate_iters,
'Median mask threshold dilations')

_msk = np.abs(filter_sci) > mask_sn_threshold*err_i
_msk = nd.binary_dilation(_msk, iterations=mask_sn_dilate_iters)
filter_sci[~_msk] = 0

if put_model_in_median & subtract_model:
msg = 'subtract_median_filter: put model in median and reset '
msg += '`object_dispersers`'
utils.log_comment(utils.LOGFILE, msg, verbose=verbose)

flt.grism.data['MED'] = (filter_sci + flt.model)*ok

flt.grism.data['SCI'] -= filter_sci*ok
flt.grism.data['MED'] = filter_sci*ok
# Reset model and object_dispersers
flt.object_dispersers = OrderedDict()
flt.model *= 0
else:
flt.grism.data['MED'] = filter_sci*ok

flt.grism.data['SCI'] -= flt.grism.data['MED']

flt.grism.header['MEDSMOD'] = subtract_model, 'Model subtracted first'
flt.grism.header['MEDWMOD'] = (put_model_in_median,
'Model included in MED')
flt.grism.header['MEDSIZE'] = filter_size, 'Median filter size'
flt.grism.header['MEDCLIP'] = (filter_central,
'Masked central pixels of the filter')
Expand Down
18 changes: 15 additions & 3 deletions grizli/pipeline/auto_script.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -3201,7 +3201,7 @@ def load_GroupFLT(field_root='j142724+334246', PREP_PATH='../Prep', force_ref=No
return [grp]


def grism_prep(field_root='j142724+334246', PREP_PATH='../Prep', EXTRACT_PATH='../Extractions', ds9=None, refine_niter=3, gris_ref_filters=GRIS_REF_FILTERS, force_ref=None, files=None, split_by_grism=True, refine_poly_order=1, refine_fcontam=0.5, cpu_count=0, mask_mosaic_edges=False, prelim_mag_limit=25, refine_mag_limits=[18, 24], init_coeffs=[1.1, -0.5], grisms_to_process=None, pad=(64, 256), model_kwargs={'compute_size': True}, sep_background_kwargs=None, subtract_median_filter=False, median_filter_size=71, median_filter_central=10, second_pass_filtering=False, box_filter_sn=3, box_filter_width=3, use_jwst_crds=False):
def grism_prep(field_root='j142724+334246', PREP_PATH='../Prep', EXTRACT_PATH='../Extractions', ds9=None, refine_niter=3, gris_ref_filters=GRIS_REF_FILTERS, force_ref=None, files=None, split_by_grism=True, refine_poly_order=1, refine_fcontam=0.5, cpu_count=0, mask_mosaic_edges=False, prelim_mag_limit=25, refine_mag_limits=[18, 24], init_coeffs=[1.1, -0.5], grisms_to_process=None, pad=(64, 256), model_kwargs={'compute_size': True}, sep_background_kwargs=None, subtract_median_filter=False, median_filter_size=71, median_filter_central=10, second_pass_filtering=False, box_filter_sn=3, box_filter_width=3, median_mask_sn_threshold=None, median_mask_dilate=8, prelim_model_for_median=False, use_jwst_crds=False):
"""
Contamination model for grism exposures
"""
Expand Down Expand Up @@ -3236,14 +3236,26 @@ def grism_prep(field_root='j142724+334246', PREP_PATH='../Prep', EXTRACT_PATH='.

if subtract_median_filter:

################
# Compute preliminary model before median?
if prelim_model_for_median:
grp.compute_full_model(fit_info=None, verbose=True, store=False,
mag_limit=prelim_mag_limit,
coeffs=init_coeffs,
cpu_count=cpu_count,
model_kwargs=model_kwargs)

# Subtract a median along the dispersion direction
# and don't do any contamination modeling
refine_niter = 0
grp.subtract_median_filter(filter_size=median_filter_size,
filter_central=median_filter_central,
second_pass_filtering=second_pass_filtering,
box_filter_sn=box_filter_sn,
box_filter_width=box_filter_width)
box_filter_width=box_filter_width,
mask_sn_threshold=median_mask_sn_threshold,
mask_sn_dilate_iters=median_mask_dilate,
subtract_model=prelim_model_for_median,
put_model_in_median=prelim_model_for_median)

else:
################
Expand Down

0 comments on commit 7d26844

Please sign in to comment.