def power_spectral_density(x, time_step, freq_range = None):
    """
    returns the *single sided* power spectral density of the time trace x which is sampled at intervals time_step
    Args:
        x (array):  timetrace
        time_step (float): sampling interval of x
        freq_range (array or tuple): frequency range in the form [f_min, f_max] to return only the spectrum within this range

    Returns:

    """
    N = len(x)
    P = 2 * np.abs(np.fft.rfft(x))**2 / N * time_step
    F = np.fft.rfftfreq(len(x), time_step)

    if freq_range is not None:
        brange = np.all([F >= freq_range[0], F <= freq_range[1]], axis=0)
        P = P[brange]
        F = F[brange]

    return F, P
returns the *single sided* power spectral density of the time trace x which is sampled at intervals time_step
    Args:
        x (array):  timetrace
        time_step (float): sampling interval of x
        freq_range (array or tuple): frequency range in the form [f_min, f_max] to return only the spectrum within this range

    Returns: