code: def find_and_refine_peaks(self, threshold, min_separation=1.0,
                                use_cumul=False):
          
          if use_cumul:
              theMap = self._ts_cumul
          else:
              theMap = self._tsmap
  
          peaks = find_peaks(theMap, threshold, min_separation)
          for peak in peaks:
              o, skydir = fit_error_ellipse(theMap, (peak['ix'], peak['iy']),
                                            dpix=2)
              peak['fit_loc'] = o
              peak['fit_skydir'] = skydir
              if o['fit_success']:
                  skydir = peak['fit_skydir']
              else:
                  skydir = peak['skydir']
  
          return peaks
docs: Run a simple peak-finding algorithm, and fit the peaks to
          paraboloids to extract their positions and error ellipses.
  
          Parameters
          ----------
          threshold : float
              Peak threshold in TS.
  
          min_separation : float
              Radius of region size in degrees.  Sets the minimum allowable
              separation between peaks.
  
          use_cumul : bool
              If true, used the cumulative TS map (i.e., the TS summed
              over the energy bins) instead of the TS Map from the fit
              to and index=2 powerlaw.
  
          Returns
          -------
          peaks    : list
              List of dictionaries containing the location and amplitude of
              each peak.  Output of `~fermipy.sourcefind.find_peaks`