diff --git a/python/lvmdrp/core/fiberrows.py b/python/lvmdrp/core/fiberrows.py index a309c2ec..bdaebb5b 100644 --- a/python/lvmdrp/core/fiberrows.py +++ b/python/lvmdrp/core/fiberrows.py @@ -870,6 +870,8 @@ def measureArcLines( bg = numpy.ones((self._fibers, nlines), dtype=numpy.float32) * numpy.nan masked = numpy.zeros((self._fibers, nlines), dtype="bool") + stretch_min, stretch_max, stretch_steps = 0.998, 1.002, 40 + spec = self.getSpec(ref_fiber) flux[ref_fiber], cent_wave[ref_fiber], fwhm[ref_fiber], bg[ref_fiber] = spec.fitSepGauss(ref_cent, aperture, fwhm_guess, bg_guess, flux_range, cent_range, fwhm_range, bg_range, axs=axs[ref_fiber][1]) masked[ref_fiber] = numpy.isnan(flux[ref_fiber])|numpy.isnan(cent_wave[ref_fiber])|numpy.isnan(fwhm[ref_fiber]) @@ -899,15 +901,19 @@ def measureArcLines( cc, bhat, mhat = _cross_match_float( ref_spec=last_spec._data, obs_spec=spec._data, - stretch_factors=numpy.linspace(0.99,1.01,20), - shift_range=[-5, 5], + stretch_factors=numpy.linspace(stretch_min, stretch_max, stretch_steps), + shift_range=[-10, 10], normalize_spectra=False, ) + if mhat == stretch_min or mhat == stretch_max: + log.warning(f"boundary of stretch factors: {mhat = } ({stretch_min, stretch_max = })") cent_guess = mhat * last_cent + bhat flux[i], cent_wave[i], fwhm[i], bg[i] = spec.fitSepGauss(cent_guess, aperture, fwhm_guess, bg_guess, flux_range, cent_range, fwhm_range, bg_range, axs=axs_fiber) masked[i] = numpy.isnan(flux[i])|numpy.isnan(cent_wave[i])|numpy.isnan(fwhm[i]) if masked[i].any(): log.warning(f"some lines were not fitted properly in fiber {i}: ") + log.warning(f" guess = {numpy.round(cent_guess, 3)} ({mhat = :.5f}, {bhat = :.5f})") + log.warning(f" mask = {masked[i]}") log.warning(f" flux = {numpy.round(flux[i],3)}") log.warning(f" cent = {numpy.round(cent_wave[i],3)}") log.warning(f" fwhm = {numpy.round(fwhm[i],3)}") @@ -939,13 +945,23 @@ def measureArcLines( cc, bhat, mhat = _cross_match_float( ref_spec=last_spec._data, obs_spec=spec._data, - stretch_factors=numpy.linspace(0.99,1.01,20), - shift_range=[-5, 5], - normalize_spectra=False + stretch_factors=numpy.linspace(stretch_min, stretch_max, stretch_steps), + shift_range=[-10, 10], + normalize_spectra=False, ) + if mhat == stretch_min or mhat == stretch_max: + log.warning(f"boundary of stretch factors: {mhat = } ({stretch_min, stretch_max = })") cent_guess = mhat * last_cent + bhat flux[i], cent_wave[i], fwhm[i], bg[i] = spec.fitSepGauss(cent_guess, aperture, fwhm_guess, bg_guess, flux_range, cent_range, fwhm_range, bg_range, axs=axs_fiber) masked[i] = numpy.isnan(flux[i])|numpy.isnan(cent_wave[i])|numpy.isnan(fwhm[i]) + if masked[i].any(): + log.warning(f"some lines were not fitted properly in fiber {i}: ") + log.warning(f" guess = {numpy.round(cent_guess, 3)} ({mhat = :.5f}, {bhat = :.5f})") + log.warning(f" mask = {masked[i]}") + log.warning(f" flux = {numpy.round(flux[i],3)}") + log.warning(f" cent = {numpy.round(cent_wave[i],3)}") + log.warning(f" fwhm = {numpy.round(fwhm[i],3)}") + log.warning(f" bg = {numpy.round(bg[i],3)}") last_spec = copy(spec) last_cent = copy(cent_wave[i]) diff --git a/python/lvmdrp/core/spectrum1d.py b/python/lvmdrp/core/spectrum1d.py index 82fab291..d05e4df0 100644 --- a/python/lvmdrp/core/spectrum1d.py +++ b/python/lvmdrp/core/spectrum1d.py @@ -290,6 +290,9 @@ def _cross_match_float( best_shift = 0 best_stretch_factor = 1 + # define pixels array centered around middle part of the spectrum + # pixels = numpy.arange(ref_spec.size) - ref_spec.size // 2 + # normalize the peaks to roughly magnitude 1, so that individual very bright # fibers do not dominate the signal if normalize_spectra: @@ -303,6 +306,7 @@ def _cross_match_float( for factor in stretch_factors: # Stretch the first signal stretched_signal1 = zoom(ref_spec_, factor, mode="constant", prefilter=True) + # stretched_signal1 = numpy.interp(pixels * factor, pixels, ref_spec_) # Make the lengths equal len_diff = len(obs_spec_) - len(stretched_signal1) @@ -3496,7 +3500,7 @@ def obtainGaussFluxPeaks(self, pos, sigma, replace_error=1e10, plot=False): if bad_pix is not None and bn.nansum(bad_pix) > 0: error[bad_pix] = replace_error - # pyfits.writeto('B.fits', B.toarray(), overwrite=True) + # pyfits.writeto('B.fits', B.toarray(), overwrite=True) # if plot: # plt.plot(self._data, "ok") # plt.plot(numpy.dot(A * self._error[:, None], out[0]), "-r") diff --git a/python/lvmdrp/functions/rssMethod.py b/python/lvmdrp/functions/rssMethod.py index baf11b95..3f758ded 100644 --- a/python/lvmdrp/functions/rssMethod.py +++ b/python/lvmdrp/functions/rssMethod.py @@ -364,14 +364,16 @@ def determine_wavelength_solution(in_arcs: List[str]|str, out_wave: str, out_lsf # fix cc_max_shift # cross-match spectrum and pixwav map + stretch_min, stretch_max, stretch_steps = 0.95, 1.05, 10000 cc, bhat, mhat = _cross_match_float( ref_spec=pix_spec, obs_spec=arc._data[ref_fiber], - stretch_factors=numpy.linspace(0.8,1.2,10000), + stretch_factors=numpy.linspace(stretch_min, stretch_max, stretch_steps), shift_range=[-cc_max_shift, cc_max_shift], normalize_spectra=False, ) - + if mhat == stretch_min or mhat == stretch_max: + log.warning(f"boundary of stretch factors: {mhat = } ({stretch_min, stretch_max = })") log.info(f"max CC = {cc:.2f} for strech = {mhat:.8f} and shift = {bhat:.8f}") else: mhat, bhat = 1.0, 0.0