From 959f40ce060ae189bbee667b9ca9f52b65165ba4 Mon Sep 17 00:00:00 2001
From: rcooke <rcooke.ast@gmail.com>
Date: Wed, 18 Dec 2024 11:21:04 +0000
Subject: [PATCH 1/7] iterate object tracing

---
 pypeit/coadd2d.py              |  1 +
 pypeit/core/findobj_skymask.py | 33 +++++++++++++++++++++++++--------
 pypeit/find_objects.py         |  2 ++
 pypeit/par/pypeitpar.py        |  9 +++++++--
 4 files changed, 35 insertions(+), 10 deletions(-)

diff --git a/pypeit/coadd2d.py b/pypeit/coadd2d.py
index 7864a2e63f..15f251e490 100644
--- a/pypeit/coadd2d.py
+++ b/pypeit/coadd2d.py
@@ -1416,6 +1416,7 @@ def compute_offsets(self, offsets):
                     inmask=inmask[iexp,:,:], fwhm=self.par['reduce']['findobj']['find_fwhm'],
                     trim_edg=self.par['reduce']['findobj']['find_trim_edge'],
                     maxdev=self.par['reduce']['findobj']['find_maxdev'],
+                    numiterfit=self.par['reduce']['findobj']['find_numiterfit'],
                     ncoeff=3, snr_thresh=self.par['reduce']['findobj']['snr_thresh'], nperslit=1,
                     find_min_max=self.par['reduce']['findobj']['find_min_max'],
                     show_trace=self.debug_offsets, show_peaks=self.debug_offsets)
diff --git a/pypeit/core/findobj_skymask.py b/pypeit/core/findobj_skymask.py
index 6e4698c844..543ea3908b 100644
--- a/pypeit/core/findobj_skymask.py
+++ b/pypeit/core/findobj_skymask.py
@@ -168,7 +168,7 @@ def ech_findobj_ineach_order(
     det='DET01', inmask=None, std_trace=None, ncoeff=5, 
     hand_extract_dict=None,
     box_radius=2.0, fwhm=3.0,
-    use_user_fwhm=False, maxdev=2.0, nperorder=2,
+    use_user_fwhm=False, maxdev=2.0, nperorder=2, numiterfit=10,
     extract_maskwidth=3.0, snr_thresh=10.0,
     specobj_dict=None, trim_edg=(5,5),
     show_peaks=False, show_single_fits=False,
@@ -260,6 +260,8 @@ def ech_findobj_ineach_order(
         maxdev (:obj:`float`, optional):
             Maximum deviation of pixels from polynomial fit to trace
             used to reject bad pixels in trace fitting.
+        numiterfit (:obj:`int`, optional):
+            Number of iterations to use when fitting the object traces.
         nperorder (:obj:`int`, optional):
             Maximum number of objects allowed per order.  If there are more
             detections than this number, the code will select the ``nperorder``
@@ -331,7 +333,7 @@ def ech_findobj_ineach_order(
                 spec_min_max=spec_min_max[:,iord],
                 inmask=inmask_iord,std_trace=std_in, 
                 ncoeff=ncoeff, fwhm=fwhm, use_user_fwhm=use_user_fwhm, maxdev=maxdev,
-                hand_extract_dict=hand_extract_dict,  
+                numiterfit=numiterfit, hand_extract_dict=hand_extract_dict,
                 nperslit=nperorder, extract_maskwidth=extract_maskwidth,
                 snr_thresh=snr_thresh, trim_edg=trim_edg, 
                 boxcar_rad=box_radius/plate_scale_ord[iord],
@@ -1038,7 +1040,7 @@ def ech_objfind(image, ivar, slitmask, slit_left, slit_righ, slit_spat_id, order
                 nabove_min_snr=2, pca_explained_var=99.0, 
                 box_radius=2.0, fwhm=3.0,
                 use_user_fwhm=False, maxdev=2.0, 
-                nperorder=2,
+                nperorder=2, numiterfit=10,
                 extract_maskwidth=3.0, snr_thresh=10.0,
                 specobj_dict=None, trim_edg=(5,5),
                 show_peaks=False, show_fits=False, 
@@ -1181,6 +1183,8 @@ def ech_objfind(image, ivar, slitmask, slit_left, slit_righ, slit_spat_id, order
         maxdev (:obj:`float`, optional):
             Maximum deviation of pixels from polynomial fit to trace
             used to reject bad pixels in trace fitting.
+        numiterfit (:obj:`int`, optional):
+            Number of iterations to perform when building the trace fits.
         hand_extract_dict (:obj:`dict`, optional):
             Dictionary with info on manual extraction; see
             :class:`~pypeit.manual_extract.ManualExtractionObj`.
@@ -1290,6 +1294,7 @@ def ech_objfind(image, ivar, slitmask, slit_left, slit_righ, slit_spat_id, order
         use_user_fwhm=use_user_fwhm,
         nperorder=nperorder,
         maxdev=maxdev,
+        numiterfit=numiterfit,
         box_radius=box_radius,
         objfindQA_filename=objfindQA_filename,
         hand_extract_dict=manual_extract_dict)
@@ -1500,7 +1505,7 @@ def get_fwhm(fwhm_in, nsamp, smash_peakflux, spat_fracpos, flux_smash_smth):
 def objs_in_slit(image, ivar, thismask, slit_left, slit_righ, 
                  inmask=None, fwhm=3.0,
                  sigclip_smash=5.0, use_user_fwhm=False, boxcar_rad=7.,
-                 maxdev=2.0, spec_min_max=None, hand_extract_dict=None, std_trace=None,
+                 maxdev=2.0, numiterfit=10, spec_min_max=None, hand_extract_dict=None, std_trace=None,
                  ncoeff=5, nperslit=None, snr_thresh=10.0, trim_edg=(5,5),
                  extract_maskwidth=4.0, specobj_dict=None, find_min_max=None,
                  show_peaks=False, show_fits=False, show_trace=False,
@@ -1588,6 +1593,10 @@ def objs_in_slit(image, ivar, thismask, slit_left, slit_righ,
         maxdev (:obj:`float`, optional):
             Maximum deviation of pixels from polynomial fit to trace
             used to reject bad pixels in trace fitting.
+        numiterfit (:obj:`int`, optional):
+            Number of iterations to use in the iterative trace fitting.
+            The number of iterations will be reduced if the fit has converged
+            to a value less than 0.1 pixels over the entire trace.
         spec_min_max (:obj:`tuple`, optional):
             2-tuple defining the minimum and maximum pixel in the spectral
             direction with useable data for this slit/order.  If None, the
@@ -1904,8 +1913,9 @@ def objs_in_slit(image, ivar, thismask, slit_left, slit_righ,
     msgs.info("Automatic finding routine found {0:d} objects".format(len(sobjs)))
 
     # Fit the object traces
+    tolerance = 0.1  # Tolerance in pixels before the trace fit has converged
     if len(sobjs) > 0:
-        msgs.info('Fitting the object traces')
+        msgs.info('Fitting the traces')
         # Note the transpose is here to pass in the TRACE_SPAT correctly.
         xinit_fweight = np.copy(sobjs.TRACE_SPAT.T).astype(float)
         spec_mask = (spec_vec >= spec_min_max_out[0]) & (spec_vec <= spec_min_max_out[1])
@@ -1914,9 +1924,16 @@ def objs_in_slit(image, ivar, thismask, slit_left, slit_righ,
                                  trace_bpm=np.invert(trc_inmask), fwhm=fwhm, maxdev=maxdev,
                                  idx=sobjs.NAME, debug=show_fits)[0]
         xinit_gweight = np.copy(xfit_fweight)
-        xfit_gweight = fit_trace(image, xinit_gweight, ncoeff, bpm=np.invert(inmask), maxshift=1.,
-                                 trace_bpm=np.invert(trc_inmask), fwhm=fwhm, maxdev=maxdev,
-                                 weighting='gaussian', idx=sobjs.NAME, debug=show_fits)[0]
+        for nn in range(numiterfit):
+            xfit_gweight = fit_trace(image, xinit_gweight, ncoeff, bpm=np.invert(inmask), maxshift=1.,
+                                     trace_bpm=np.invert(trc_inmask), fwhm=fwhm, maxdev=maxdev,
+                                     weighting='gaussian', idx=sobjs.NAME, debug=show_fits)[0]
+            maxdev = np.max(np.abs(xfit_gweight-xinit_gweight))
+            xinit_gweight = np.copy(xfit_gweight)
+            if maxdev <= tolerance:
+                break
+        if maxdev > tolerance:
+            msgs.warn('Trace fitting did not converge.  Check the fits.')
 
         # assign the final trace
         for iobj in range(nobj_reg):
diff --git a/pypeit/find_objects.py b/pypeit/find_objects.py
index 44ba731508..ba619dc045 100644
--- a/pypeit/find_objects.py
+++ b/pypeit/find_objects.py
@@ -819,6 +819,7 @@ def find_objects_pypeline(self, image, ivar, std_trace=None,
                                 use_user_fwhm=self.par['reduce']['extraction']['use_user_fwhm'],
                                 boxcar_rad=self.par['reduce']['extraction']['boxcar_radius'] / self.get_platescale(),  #pixels
                                 maxdev=self.par['reduce']['findobj']['find_maxdev'],
+                                numiterfit=self.par['reduce']['findobj']['find_numiterfit'],
                                 find_min_max=self.par['reduce']['findobj']['find_min_max'],
                                 extract_maskwidth=self.par['reduce']['skysub']['local_maskwidth'],
                                 qa_title=qa_title, nperslit=maxnumber,
@@ -959,6 +960,7 @@ def find_objects_pypeline(self, image, ivar, std_trace=None,
             use_user_fwhm=self.par['reduce']['extraction']['use_user_fwhm'],
             fof_link = self.par['reduce']['findobj']['fof_link'],
             maxdev=self.par['reduce']['findobj']['find_maxdev'],
+            numiterfit=self.par['reduce']['findobj']['find_numiterfit'],
             nperorder=nperorder,
             max_snr=self.par['reduce']['findobj']['ech_find_max_snr'],
             min_snr=self.par['reduce']['findobj']['ech_find_min_snr'],
diff --git a/pypeit/par/pypeitpar.py b/pypeit/par/pypeitpar.py
index a19680b156..fe79cd4edf 100644
--- a/pypeit/par/pypeitpar.py
+++ b/pypeit/par/pypeitpar.py
@@ -4146,7 +4146,7 @@ class FindObjPar(ParSet):
 
     def __init__(self, trace_npoly=None, snr_thresh=None, find_trim_edge=None,
                  find_maxdev=None, find_extrap_npoly=None, maxnumber_sci=None, maxnumber_std=None,
-                 find_fwhm=None, ech_find_max_snr=None, ech_find_min_snr=None,
+                 find_fwhm=None, ech_find_max_snr=None, ech_find_min_snr=None, find_numiterfit=None,
                  ech_find_nabove_min_snr=None, skip_second_find=None, skip_final_global=None,
                  skip_skysub=None, find_negative=None, find_min_max=None, std_spec1d=None,
                  use_std_trace=None, fof_link = None):
@@ -4206,6 +4206,11 @@ def __init__(self, trace_npoly=None, snr_thresh=None, find_trim_edge=None,
         dtypes['find_maxdev'] = [int, float]
         descr['find_maxdev'] = 'Maximum deviation of pixels from polynomial fit to trace used to reject bad pixels in trace fitting.'
 
+        defaults['find_numiterfit'] = 10
+        dtypes['find_numiterfit'] = int
+        descr['find_numiterfit'] = 'Number of iterations to perform on the trace fitting. If the trace fitting ' \
+                                   'has converged, the number of iterations will be less than this number.'
+
         defaults['find_fwhm'] = 5.0
         dtypes['find_fwhm'] = [int, float]
         descr['find_fwhm'] = 'Indicates roughly the fwhm of objects in pixels for object finding'
@@ -4304,7 +4309,7 @@ def from_dict(cls, cfg):
         # Basic keywords
         parkeys = ['trace_npoly', 'snr_thresh', 'find_trim_edge',
                    'find_extrap_npoly', 'maxnumber_sci', 'maxnumber_std',
-                   'find_maxdev', 'find_fwhm', 'ech_find_max_snr',
+                   'find_maxdev', 'find_numiterfit', 'find_fwhm', 'ech_find_max_snr',
                    'ech_find_min_snr', 'ech_find_nabove_min_snr', 'skip_second_find',
                    'skip_final_global', 'skip_skysub', 'find_negative', 'find_min_max',
                    'std_spec1d', 'use_std_trace', 'fof_link']

From 462171a52b0b7f9dd9d190a3ec03edf7a12f967e Mon Sep 17 00:00:00 2001
From: rcooke <rcooke.ast@gmail.com>
Date: Thu, 19 Dec 2024 13:38:17 +0000
Subject: [PATCH 2/7] updated defaults

---
 pypeit/core/findobj_skymask.py | 6 +++---
 pypeit/par/pypeitpar.py        | 2 +-
 2 files changed, 4 insertions(+), 4 deletions(-)

diff --git a/pypeit/core/findobj_skymask.py b/pypeit/core/findobj_skymask.py
index 543ea3908b..3cb5917153 100644
--- a/pypeit/core/findobj_skymask.py
+++ b/pypeit/core/findobj_skymask.py
@@ -168,7 +168,7 @@ def ech_findobj_ineach_order(
     det='DET01', inmask=None, std_trace=None, ncoeff=5, 
     hand_extract_dict=None,
     box_radius=2.0, fwhm=3.0,
-    use_user_fwhm=False, maxdev=2.0, nperorder=2, numiterfit=10,
+    use_user_fwhm=False, maxdev=2.0, nperorder=2, numiterfit=1,
     extract_maskwidth=3.0, snr_thresh=10.0,
     specobj_dict=None, trim_edg=(5,5),
     show_peaks=False, show_single_fits=False,
@@ -1040,7 +1040,7 @@ def ech_objfind(image, ivar, slitmask, slit_left, slit_righ, slit_spat_id, order
                 nabove_min_snr=2, pca_explained_var=99.0, 
                 box_radius=2.0, fwhm=3.0,
                 use_user_fwhm=False, maxdev=2.0, 
-                nperorder=2, numiterfit=10,
+                nperorder=2, numiterfit=1,
                 extract_maskwidth=3.0, snr_thresh=10.0,
                 specobj_dict=None, trim_edg=(5,5),
                 show_peaks=False, show_fits=False, 
@@ -1505,7 +1505,7 @@ def get_fwhm(fwhm_in, nsamp, smash_peakflux, spat_fracpos, flux_smash_smth):
 def objs_in_slit(image, ivar, thismask, slit_left, slit_righ, 
                  inmask=None, fwhm=3.0,
                  sigclip_smash=5.0, use_user_fwhm=False, boxcar_rad=7.,
-                 maxdev=2.0, numiterfit=10, spec_min_max=None, hand_extract_dict=None, std_trace=None,
+                 maxdev=2.0, numiterfit=1, spec_min_max=None, hand_extract_dict=None, std_trace=None,
                  ncoeff=5, nperslit=None, snr_thresh=10.0, trim_edg=(5,5),
                  extract_maskwidth=4.0, specobj_dict=None, find_min_max=None,
                  show_peaks=False, show_fits=False, show_trace=False,
diff --git a/pypeit/par/pypeitpar.py b/pypeit/par/pypeitpar.py
index fe79cd4edf..fe9e0c34fe 100644
--- a/pypeit/par/pypeitpar.py
+++ b/pypeit/par/pypeitpar.py
@@ -4206,7 +4206,7 @@ def __init__(self, trace_npoly=None, snr_thresh=None, find_trim_edge=None,
         dtypes['find_maxdev'] = [int, float]
         descr['find_maxdev'] = 'Maximum deviation of pixels from polynomial fit to trace used to reject bad pixels in trace fitting.'
 
-        defaults['find_numiterfit'] = 10
+        defaults['find_numiterfit'] = 1
         dtypes['find_numiterfit'] = int
         descr['find_numiterfit'] = 'Number of iterations to perform on the trace fitting. If the trace fitting ' \
                                    'has converged, the number of iterations will be less than this number.'

From 5e67de1db734826b25883f17509bf1d534edd12c Mon Sep 17 00:00:00 2001
From: rcooke <rcooke.ast@gmail.com>
Date: Thu, 19 Dec 2024 13:38:26 +0000
Subject: [PATCH 3/7] updated docs

---
 doc/object_finding.rst     | 27 +++++++++++++++++++++++++++
 doc/releases/1.17.1dev.rst |  2 ++
 2 files changed, 29 insertions(+)

diff --git a/doc/object_finding.rst b/doc/object_finding.rst
index 6dde252401..bdc6bcd9ca 100644
--- a/doc/object_finding.rst
+++ b/doc/object_finding.rst
@@ -138,4 +138,31 @@ vs spatial position vector. The best way to choose these pixels is to run PypeIt
 without it set. Then run :ref:`pypeit_show_2dspec` to view the sky-subtracted
 image and decide which pixels to use for object finding. Then re-run ``PypeIt``.
 
+Object Tracing
+==============
+
+For automatically identified objects (i.e., not manual extractions),
+PypeIt improves the object trace by fitting the spatial position of the peak
+as a function of wavelength. In some situations, the object trace is poorly
+determined by the peak location, and the code will fail to trace the object
+correctly. For example, if the slit edges are not well-defined, the object's
+position relative to the slit edges is also poorly defined, and the object trace
+is difficult to determine. In these cases, the user can attempt to iterate
+the object tracing by using a relatively low order polynomial in combination
+with an iterative fit, as follows
+
+.. code-block:: ini
 
+    [reduce]
+        [[findobj]]
+            find_numiterfit = 10
+            trace_npoly = 4
+
+Note that the default value is typically ``trace_npoly=5``. Also note that the
+object tracing will break if successive iterations yields the same object trace
+to within 0.1 pixels. If you notice a relatively poor object trace, sometimes in
+combination with the object counts being masked, an iterative fit may help to
+resolve your problem. If, on the other hand, your object is relatively faint, you
+may benefit from using the trace of a standard star (this is the default behaviour),
+and you can provide a 1D spectrum of a previously reduced standard star with the
+``std_spec1d`` parameter.
diff --git a/doc/releases/1.17.1dev.rst b/doc/releases/1.17.1dev.rst
index 194bea27ad..c7731f4c5a 100644
--- a/doc/releases/1.17.1dev.rst
+++ b/doc/releases/1.17.1dev.rst
@@ -17,6 +17,8 @@ Functionality/Performance Improvements and Additions
 - The WCS for datacubes now adopts the convention of North
   is up and East is left. In previous version of PypeIt,
   East was right.
+- The object trace can now be iterated over with the `find_numiterfit`
+  parameter.
 
 Instrument-specific Updates
 ---------------------------

From af1812fff75f5c89cf7c921d589de35b7d036ce1 Mon Sep 17 00:00:00 2001
From: rcooke <rcooke.ast@gmail.com>
Date: Fri, 20 Dec 2024 12:04:09 +0000
Subject: [PATCH 4/7] cubic bugfix

---
 doc/object_finding.rst         | 15 +++++++--------
 pypeit/core/findobj_skymask.py | 28 ++++++++++------------------
 pypeit/par/pypeitpar.py        |  5 ++---
 3 files changed, 19 insertions(+), 29 deletions(-)

diff --git a/doc/object_finding.rst b/doc/object_finding.rst
index bdc6bcd9ca..5fabb287a8 100644
--- a/doc/object_finding.rst
+++ b/doc/object_finding.rst
@@ -147,21 +147,20 @@ as a function of wavelength. In some situations, the object trace is poorly
 determined by the peak location, and the code will fail to trace the object
 correctly. For example, if the slit edges are not well-defined, the object's
 position relative to the slit edges is also poorly defined, and the object trace
-is difficult to determine. In these cases, the user can attempt to iterate
-the object tracing by using a relatively low order polynomial in combination
-with an iterative fit, as follows
+is difficult to determine. The default is to perform several iterations (typically 9)
+but for some cases this is insufficient. In these cases, the user can attempt to
+increase the number of iterations to improve the object tracing, in combination with
+a relatively low order polynomial, as follows
 
 .. code-block:: ini
 
     [reduce]
         [[findobj]]
-            find_numiterfit = 10
+            find_numiterfit = 100
             trace_npoly = 4
 
-Note that the default value is typically ``trace_npoly=5``. Also note that the
-object tracing will break if successive iterations yields the same object trace
-to within 0.1 pixels. If you notice a relatively poor object trace, sometimes in
-combination with the object counts being masked, an iterative fit may help to
+Note that the default value is typically ``trace_npoly=5``. If you notice a relatively poor object trace, sometimes in
+combination with the object counts being masked, increasing the number of iterations may help to
 resolve your problem. If, on the other hand, your object is relatively faint, you
 may benefit from using the trace of a standard star (this is the default behaviour),
 and you can provide a 1D spectrum of a previously reduced standard star with the
diff --git a/pypeit/core/findobj_skymask.py b/pypeit/core/findobj_skymask.py
index 3cb5917153..e4ba3f964d 100644
--- a/pypeit/core/findobj_skymask.py
+++ b/pypeit/core/findobj_skymask.py
@@ -168,7 +168,7 @@ def ech_findobj_ineach_order(
     det='DET01', inmask=None, std_trace=None, ncoeff=5, 
     hand_extract_dict=None,
     box_radius=2.0, fwhm=3.0,
-    use_user_fwhm=False, maxdev=2.0, nperorder=2, numiterfit=1,
+    use_user_fwhm=False, maxdev=2.0, nperorder=2, numiterfit=9,
     extract_maskwidth=3.0, snr_thresh=10.0,
     specobj_dict=None, trim_edg=(5,5),
     show_peaks=False, show_single_fits=False,
@@ -1040,7 +1040,7 @@ def ech_objfind(image, ivar, slitmask, slit_left, slit_righ, slit_spat_id, order
                 nabove_min_snr=2, pca_explained_var=99.0, 
                 box_radius=2.0, fwhm=3.0,
                 use_user_fwhm=False, maxdev=2.0, 
-                nperorder=2, numiterfit=1,
+                nperorder=2, numiterfit=9,
                 extract_maskwidth=3.0, snr_thresh=10.0,
                 specobj_dict=None, trim_edg=(5,5),
                 show_peaks=False, show_fits=False, 
@@ -1505,7 +1505,7 @@ def get_fwhm(fwhm_in, nsamp, smash_peakflux, spat_fracpos, flux_smash_smth):
 def objs_in_slit(image, ivar, thismask, slit_left, slit_righ, 
                  inmask=None, fwhm=3.0,
                  sigclip_smash=5.0, use_user_fwhm=False, boxcar_rad=7.,
-                 maxdev=2.0, numiterfit=1, spec_min_max=None, hand_extract_dict=None, std_trace=None,
+                 maxdev=2.0, numiterfit=9, spec_min_max=None, hand_extract_dict=None, std_trace=None,
                  ncoeff=5, nperslit=None, snr_thresh=10.0, trim_edg=(5,5),
                  extract_maskwidth=4.0, specobj_dict=None, find_min_max=None,
                  show_peaks=False, show_fits=False, show_trace=False,
@@ -1595,8 +1595,6 @@ def objs_in_slit(image, ivar, thismask, slit_left, slit_righ,
             used to reject bad pixels in trace fitting.
         numiterfit (:obj:`int`, optional):
             Number of iterations to use in the iterative trace fitting.
-            The number of iterations will be reduced if the fit has converged
-            to a value less than 0.1 pixels over the entire trace.
         spec_min_max (:obj:`tuple`, optional):
             2-tuple defining the minimum and maximum pixel in the spectral
             direction with useable data for this slit/order.  If None, the
@@ -1676,11 +1674,12 @@ def objs_in_slit(image, ivar, thismask, slit_left, slit_righ,
         detected.
     """
 
-    #debug_all=True
+    debug_all = False
     if debug_all:
-        show_peaks=True
+        show_peaks = True
         show_fits = True
         show_trace = True
+    show_fits=True
 
     if specobj_dict is None:
         specobj_dict = dict(SLITID=999, DET='DET01', OBJTYPE='unknown', PYPELINE='MultiSlit')
@@ -1920,20 +1919,13 @@ def objs_in_slit(image, ivar, thismask, slit_left, slit_righ,
         xinit_fweight = np.copy(sobjs.TRACE_SPAT.T).astype(float)
         spec_mask = (spec_vec >= spec_min_max_out[0]) & (spec_vec <= spec_min_max_out[1])
         trc_inmask = np.outer(spec_mask, np.ones(len(sobjs), dtype=bool))
-        xfit_fweight = fit_trace(image, xinit_fweight, ncoeff, bpm=np.invert(inmask), maxshift=1.,
+        xfit_fweight = fit_trace(image, xinit_fweight, ncoeff, bpm=np.invert(inmask), maxshift=1., niter=numiterfit,
                                  trace_bpm=np.invert(trc_inmask), fwhm=fwhm, maxdev=maxdev,
                                  idx=sobjs.NAME, debug=show_fits)[0]
         xinit_gweight = np.copy(xfit_fweight)
-        for nn in range(numiterfit):
-            xfit_gweight = fit_trace(image, xinit_gweight, ncoeff, bpm=np.invert(inmask), maxshift=1.,
-                                     trace_bpm=np.invert(trc_inmask), fwhm=fwhm, maxdev=maxdev,
-                                     weighting='gaussian', idx=sobjs.NAME, debug=show_fits)[0]
-            maxdev = np.max(np.abs(xfit_gweight-xinit_gweight))
-            xinit_gweight = np.copy(xfit_gweight)
-            if maxdev <= tolerance:
-                break
-        if maxdev > tolerance:
-            msgs.warn('Trace fitting did not converge.  Check the fits.')
+        xfit_gweight = fit_trace(image, xinit_gweight, ncoeff, bpm=np.invert(inmask), maxshift=1., niter=numiterfit,
+                                 trace_bpm=np.invert(trc_inmask), fwhm=fwhm, maxdev=maxdev,
+                                 weighting='gaussian', idx=sobjs.NAME, debug=show_fits)[0]
 
         # assign the final trace
         for iobj in range(nobj_reg):
diff --git a/pypeit/par/pypeitpar.py b/pypeit/par/pypeitpar.py
index fe9e0c34fe..07c3c0cdc1 100644
--- a/pypeit/par/pypeitpar.py
+++ b/pypeit/par/pypeitpar.py
@@ -4206,10 +4206,9 @@ def __init__(self, trace_npoly=None, snr_thresh=None, find_trim_edge=None,
         dtypes['find_maxdev'] = [int, float]
         descr['find_maxdev'] = 'Maximum deviation of pixels from polynomial fit to trace used to reject bad pixels in trace fitting.'
 
-        defaults['find_numiterfit'] = 1
+        defaults['find_numiterfit'] = 9
         dtypes['find_numiterfit'] = int
-        descr['find_numiterfit'] = 'Number of iterations to perform on the trace fitting. If the trace fitting ' \
-                                   'has converged, the number of iterations will be less than this number.'
+        descr['find_numiterfit'] = 'Number of iterations to perform on the trace fitting.'
 
         defaults['find_fwhm'] = 5.0
         dtypes['find_fwhm'] = [int, float]

From 9a18ef180180cfcf89aca277fe1ce364d1532da6 Mon Sep 17 00:00:00 2001
From: rcooke <rcooke.ast@gmail.com>
Date: Fri, 20 Dec 2024 12:05:23 +0000
Subject: [PATCH 5/7] trace iteration as free parameter

---
 pypeit/core/findobj_skymask.py | 1 -
 1 file changed, 1 deletion(-)

diff --git a/pypeit/core/findobj_skymask.py b/pypeit/core/findobj_skymask.py
index e4ba3f964d..44a0301783 100644
--- a/pypeit/core/findobj_skymask.py
+++ b/pypeit/core/findobj_skymask.py
@@ -1679,7 +1679,6 @@ def objs_in_slit(image, ivar, thismask, slit_left, slit_righ,
         show_peaks = True
         show_fits = True
         show_trace = True
-    show_fits=True
 
     if specobj_dict is None:
         specobj_dict = dict(SLITID=999, DET='DET01', OBJTYPE='unknown', PYPELINE='MultiSlit')

From 9bb74d17dd8fca8c358367562e11ca8d62fdc0a5 Mon Sep 17 00:00:00 2001
From: rcooke <rcooke.ast@gmail.com>
Date: Mon, 6 Jan 2025 11:47:42 +0000
Subject: [PATCH 6/7] rm unused tolerance

---
 pypeit/core/findobj_skymask.py | 1 -
 1 file changed, 1 deletion(-)

diff --git a/pypeit/core/findobj_skymask.py b/pypeit/core/findobj_skymask.py
index 44a0301783..b33ddfc577 100644
--- a/pypeit/core/findobj_skymask.py
+++ b/pypeit/core/findobj_skymask.py
@@ -1911,7 +1911,6 @@ def objs_in_slit(image, ivar, thismask, slit_left, slit_righ,
     msgs.info("Automatic finding routine found {0:d} objects".format(len(sobjs)))
 
     # Fit the object traces
-    tolerance = 0.1  # Tolerance in pixels before the trace fit has converged
     if len(sobjs) > 0:
         msgs.info('Fitting the traces')
         # Note the transpose is here to pass in the TRACE_SPAT correctly.

From 1b0bf79b0f0d7fc7a174d37a99c46a6c896fe24a Mon Sep 17 00:00:00 2001
From: rcooke <rcooke.ast@gmail.com>
Date: Mon, 6 Jan 2025 11:51:56 +0000
Subject: [PATCH 7/7] doc update to mbar

---
 pypeit/coadd3d.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/pypeit/coadd3d.py b/pypeit/coadd3d.py
index 5638cb6f77..10b9cbde6f 100644
--- a/pypeit/coadd3d.py
+++ b/pypeit/coadd3d.py
@@ -261,7 +261,7 @@ def __init__(self, airmass, parangle, pressure, temperature, humidity, cosdec, w
             parangle (:obj:`float`):
                 The parallactic angle of the observations (units=radians, relative to North, towards East is postive)
             pressure (:obj:`float`):
-                The atmospheric pressure during the observations in Pascal. Valid range is from 10kPa - 140 kPa.
+                The atmospheric pressure during the observations in mbar. Valid range is from 100 mbar - 1400 mbar.
             temperature (:obj:`float`):
                 Temperature in degree Celsius. Valid temperate range is -40 to
                 100 degree Celsius.