From 28eec30f25f57a687c895b8bdbdd2081809190cb Mon Sep 17 00:00:00 2001 From: Stephen Bailey Date: Mon, 18 Dec 2023 16:30:50 -0800 Subject: [PATCH 1/2] more robust traceshifts when half CCD is masked --- py/desispec/trace_shifts.py | 21 +++++++++++++++------ 1 file changed, 15 insertions(+), 6 deletions(-) diff --git a/py/desispec/trace_shifts.py b/py/desispec/trace_shifts.py index 47dc50c4c..b306c5cbd 100644 --- a/py/desispec/trace_shifts.py +++ b/py/desispec/trace_shifts.py @@ -417,8 +417,11 @@ def compute_dy_using_boxcar_extraction(xytraceset, image, fibers, width=7, degyy # resampling on common finer wavelength grid flux, ivar, wave = resample_boxcar_frame(qframe.flux, qframe.ivar, qframe.wave, oversampling=4) - # median flux used as internal spectral reference - mflux=np.median(flux,axis=0) + # boolean mask of fibers with good data + good_fibers = (np.sum(ivar>0, axis=1) > 0) + + # median flux of good fibers used as internal spectral reference + mflux=np.median(flux[good_fibers],axis=0) # measure y shifts wavemin = xytraceset.wavemin @@ -661,24 +664,30 @@ def shift_ycoef_using_external_spectrum(psf, xytraceset, image, fibers, # here we get rid of continuum by applying a median filter continuum_win = 17 continuum_foot = np.abs(np.arange(-continuum_win,continuum_win))>continuum_win /2. + # we only keep emission lines and get rid of continuum for ii in range(flux.shape[0]): flux[ii] = flux[ii] - median_filter(flux[ii], footprint=continuum_foot) - mflux = np.median(flux, axis=0) + + # boolean mask of fibers with good data + good_fibers = (np.sum(ivar>0, axis=1) > 0) + num_good_fibers = np.sum(good_fibers) + # median flux used as internal spectral reference + mflux = np.median(flux[good_fibers], axis=0) # we use data variance and MAD from different spectra # to assign variance to a spectrum (1.48 is MAD factor, # pi/2 is a factor from Stddev[median(N(0,1))] mad_factor = 1.48 - mad = np.maximum(np.median(np.abs(flux - mflux[None, :]), + mad = np.maximum(np.median(np.abs(flux[good_fibers] - mflux[None, :]), axis=0), 1e-100) # I prevent it from being zero to avoid the warning below # The exact value does not matter as we're comparing to actual # median(ivar) mivar = np.minimum( - np.median(ivar, axis=0) , - 1./mad_factor**2 / mad**2) * flux.shape[0] * (2. / np.pi) + np.median(ivar[good_fibers], axis=0) , + 1./mad_factor**2 / mad**2) * num_good_fibers * (2. / np.pi) # finally use use the MAD of the background subtracted spectra to # assign further variance limit # this is sort of "effective" noise in the continuum subtracted spectrum From c8178e50a9c7608d2b26dc03d84040e5ab8d9cb8 Mon Sep 17 00:00:00 2001 From: Stephen Bailey Date: Mon, 18 Dec 2023 16:35:15 -0800 Subject: [PATCH 2/2] plot lightgray traces to be visible on white masked regions --- bin/plot_fiber_traces | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/bin/plot_fiber_traces b/bin/plot_fiber_traces index 6a7c92e3e..ea6b547d9 100755 --- a/bin/plot_fiber_traces +++ b/bin/plot_fiber_traces @@ -95,7 +95,7 @@ for i,fiber in enumerate(fibers) : x = tset.x_vs_wave(fiber,wave) y = tset.y_vs_wave(fiber,wave) color=None - if args.image is not None: color="white" + if args.image is not None: color="lightgray" plt.plot(x,y,color=color) @@ -107,7 +107,7 @@ if lines is not None : xl[i] = tset.x_vs_wave(fiber,line) yl[i] = tset.y_vs_wave(fiber,line) color=None - if args.image is not None: color="white" + if args.image is not None: color="lightgray" plt.plot(xl,yl,color=color) plt.xlabel("xccd")