docs: def image_position_from_source(self, sourcePos_x, sourcePos_y, kwargs_lens, min_distance=0.1, search_window=10,
                                     precision_limit=10**(-10), num_iter_max=100, arrival_time_sort=True,
                                     initial_guess_cut=True, verbose=False, x_center=0, y_center=0, num_random=0):
          """
          finds image position source position and lense model
  
          :param sourcePos_x: source position in units of angle
          :param sourcePos_y: source position in units of angle
          :param kwargs_lens: lens model parameters as keyword arguments
          :param min_distance: minimum separation to consider for two images in units of angle
          :param search_window: window size to be considered by the solver. Will not find image position outside this window
          :param precision_limit: required precision in the lens equation solver (in units of angle in the source plane).
          :param num_iter_max: maximum iteration of lens-source mapping conducted by solver to match the required precision
          :param arrival_time_sort: bool, if True, sorts image position in arrival time (first arrival photon first listed)
          :param initial_guess_cut: bool, if True, cuts initial local minima selected by the grid search based on distance criteria from the source position
          :param verbose: bool, if True, prints some useful information for the user
          :param x_center: float, center of the window to search for point sources
          :param y_center: float, center of the window to search for point sources
          :returns: (exact) angular position of (multiple) images ra_pos, dec_pos in units of angle
          :raises: AttributeError, KeyError
          """
  
          # compute number of pixels to cover the search window with the required min_distance
          numPix = int(round(search_window / min_distance) + 0.5)
          x_grid, y_grid = util.make_grid(numPix, min_distance)
          x_grid += x_center
          y_grid += y_center
          # ray-shoot to find the relative distance to the required source position for each grid point
          x_mapped, y_mapped = self.lensModel.ray_shooting(x_grid, y_grid, kwargs_lens)
          absmapped = util.displaceAbs(x_mapped, y_mapped, sourcePos_x, sourcePos_y)
          # select minima in the grid points and select grid points that do not deviate more than the
          # width of the grid point to a solution of the lens equation
          x_mins, y_mins, delta_map = util.neighborSelect(absmapped, x_grid, y_grid)
          if verbose is True:
              print("There are %s regions identified that could contain a solution of the lens equation" % len(x_mins))
          #mag = np.abs(mag)
          #print(x_mins, y_mins, 'before requirement of min_distance')
          if initial_guess_cut is True:
              mag = np.abs(self.lensModel.magnification(x_mins, y_mins, kwargs_lens))
              mag[mag < 1] = 1
              x_mins = x_mins[delta_map <= min_distance*mag*5]
              y_mins = y_mins[delta_map <= min_distance*mag*5]
              if verbose is True:
                  print("The number of regions that meet the plausibility criteria are %s" % len(x_mins))
          x_mins = np.append(x_mins, np.random.uniform(low=-search_window/2+x_center, high=search_window/2+x_center,
                                                       size=num_random))
          y_mins = np.append(y_mins, np.random.uniform(low=-search_window / 2 + y_center, high=search_window / 2 + y_center,
                                               size=num_random))
          # iterative solving of the lens equation for the selected grid points
          x_mins, y_mins, solver_precision = self._findIterative(x_mins, y_mins, sourcePos_x, sourcePos_y, kwargs_lens,
                                                                 precision_limit, num_iter_max, verbose=verbose,
                                                                 min_distance=min_distance)
          # only select iterative results that match the precision limit
          x_mins = x_mins[solver_precision <= precision_limit]
          y_mins = y_mins[solver_precision <= precision_limit]
          # find redundant solutions within the min_distance criterion
          x_mins, y_mins = image_util.findOverlap(x_mins, y_mins, min_distance)
          if arrival_time_sort is True:
              x_mins, y_mins = self.sort_arrival_times(x_mins, y_mins, kwargs_lens)
          #x_mins, y_mins = lenstronomy_util.coordInImage(x_mins, y_mins, numPix, deltapix)
          return x_mins, y_mins