def continuous_magnitude_peak(disc_peak_index, disc_magnitude, exponent=0.2308):
    """Estimate the value and location of a continuous spectrum peak.

    Applies quadratic interpolation to the discrete peak and its two neighbors
    to estimate the location and value of the continuous peak. Weighs the
    peaks using an exponential function to increase the accuracy of
    the estimate.

    Based on the paper The XQIFFT: Increasing the Accuracy of Quadratic
    Interpolation of Spectral Peaks via Exponeltial Magnitude Spectrum
    Weighting by Kurt James Werner.

    Parameters
    ----------
    disc_peak_index : integer
        Index in the `disc_magnitude` array where the magnitude peak
        is located.
    disc_magnitude : array_like
        Magnitude of the discrete spectrum with the peak of interest.
    exponent : float
        Exponent of the weighting function.

    Returns
    -------
    float
        A value between `0` and `len(disc_magnitude)-1` indicating the
        location of the peak in the continuous spectrum.
    float
        Estimated peak magnitude in the continuous spectrum.

    """
    # Weighing function and its inverse
    omega = lambda x: x**exponent
    inv_omega = lambda x: x**(1/exponent)

    # Weigh alpha, beta and gama
    alpha = omega(disc_magnitude[disc_peak_index-1])
    beta = omega(disc_magnitude[disc_peak_index])
    gamma = omega(disc_magnitude[disc_peak_index+1])

    # Estimate continuous peak index (peak location)
    cont_peak_index = (disc_peak_index
                       + (1/2)*(alpha - gamma)/(alpha - 2*beta + gamma))
    # Estimate continuous peak magnitude
    cont_peak = beta - (1/8)*(alpha - gamma)**2/(alpha - 2*beta + gamma)
    # Unweigh peak magnitude
    cont_peak = inv_omega(cont_peak)

    return cont_peak_index, cont_peak