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):
        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
        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)
        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))
        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))
        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)
        x_mins = x_mins[solver_precision <= precision_limit]
        y_mins = y_mins[solver_precision <= precision_limit]
        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)
        return x_mins, y_mins