def detect_r_peaks(ecg_signal, sample_rate, time_units=False, volts=False, resolution=None,
                   device="biosignalsplux", plot_result=False):
    """
    -----
    Brief
    -----
    Python implementation of R peak detection algorithm (proposed by Raja Selvaraj).

    -----------
    Description
    -----------
    Pan-Tompkins algorithm is one of the gold-standard algorithms in R-peak detection on ECG due to its low
    computational complexity, which allows for real-time applications, preserving high accuracy values.

    This function allows the detection of these events in ECG signals using the Pan-Tompkins.

    ----------
    Parameters
    ----------
    ecg_signal : list
        List of ECG acquired samples.

    sample_rate : int
        Sampling frequency.

    time_units : boolean
        If True this function will return the R peak position in seconds.

    volts : boolean
        If True, then the conversion of raw units to mV will be done. Resolution needs to be
        specified.

    resolution : int or None
        Selected resolution for data acquisition.

    device : str
        Specification of the device category.

    plot_result : boolean
        If True it will be presented a graphical representation of the R peak position in the ECG
        signal.

    Returns
    -------
    out : R peak position (ndarray), R peak amplitude (ndarray)
        R peak position (sample number or time instant in seconds) and amplitude (raw or mV).

    """

    if volts is True:
        if resolution is not None:
            # ecg_signal = ((ecg_signal / 2 ** resolution) - 0.5) * 3
            ecg_signal = raw_to_phy("ECG", device, ecg_signal, resolution, option="mV")
        else:
            raise RuntimeError("For converting raw units to mV is mandatory the specification of "
                               "acquisition resolution.")

    if time_units is True:
        time = numpy.linspace(0, len(ecg_signal) / sample_rate, len(ecg_signal))
    else:
        time = numpy.linspace(0, len(ecg_signal) - 1, len(ecg_signal))

    # Filtering Step of Pan-Tompkins Algorithm.
    filtered = _ecg_band_pass_filter(ecg_signal, sample_rate)

    # Differentiation Step of Pan-Tompkins Algorithm.
    differentiated = _differentiate(filtered)

    # Rectification Step of Pan-Tompkins Algorithm.
    squared = _squaring(differentiated)

    # Integration Step of Pan-Tompkins Algorithm.
    integrated = _integration(squared, sample_rate)

    rr_buffer, spk1, npk1, threshold = _buffer_ini(integrated, sample_rate)
    probable_peaks, possible_peaks = _detects_peaks(integrated, sample_rate)
    definitive_peaks = _checkup(probable_peaks, integrated, sample_rate, rr_buffer, spk1, npk1,
                                threshold)
    definitive_peaks = list(map(int, definitive_peaks))

    # Rephasing step.
    definitive_peaks_rephase = numpy.array(definitive_peaks) - 30 * (sample_rate / 1000)
    definitive_peaks_rephase = list(map(int, definitive_peaks_rephase))

    if time_units is True:
        peaks = numpy.array(time)[definitive_peaks_rephase]
    else:
        peaks = definitive_peaks_rephase

    amplitudes = numpy.array(ecg_signal)[definitive_peaks_rephase]

    # If plot is invoked by plot_result flag, then a graphical representation of the R peaks is
    # presented to the user.
    if plot_result is True:
        time_int = numpy.array(time[1:])
        integrated = numpy.array(integrated)

        fig = figure(x_axis_label='Time (s)', y_axis_label='Raw Data',
                     **opensignals_kwargs("figure"))
        fig.line(time_int, integrated, **opensignals_kwargs("line"))
        fig.circle(time_int[definitive_peaks], integrated[definitive_peaks], size=30,
                   color="#00893E", legend="Definitive Peaks")
        fig.circle(time_int[probable_peaks], integrated[probable_peaks], size=20, color="#009EE3",
                   legend="Probable Peaks")
        fig.circle(time_int[possible_peaks], integrated[possible_peaks], size=10, color="#302683",
                   legend="Possible Peaks")

        fig2 = figure(x_axis_label='Time (s)', y_axis_label='Raw Data',
                      **opensignals_kwargs("figure"))
        fig2.line(time, ecg_signal, **opensignals_kwargs("line"))
        fig2.circle(time[definitive_peaks_rephase],
                    numpy.array(ecg_signal)[definitive_peaks_rephase],
                    size=30, color=opensignals_color_pallet(), legend="Definitive Peaks")

        opensignals_style([fig, fig2])

        grid_plot = gridplot([[fig], [fig2]], **opensignals_kwargs("gridplot"))
        show(grid_plot)

    return peaks, amplitudes