def SNR_kelly(spiketrain):
    """
    returns the SNR of the waveforms of spiketrains, as computed in
    Kelly et al (2007):
    * compute the mean waveform
    * define the signal as the peak-to-through of such mean waveform
    * define the noise as double the std.dev. of the values of all waveforms,
      each normalised by subtracting the mean waveform

    Parameters:
    -----------
    spiketrain : SpikeTrain
        spike train loaded with rgio (has attribute "vaweforms")

    Returns:
    --------
    snr: float
        The SNR of the input spike train
    """
    mean_waveform = spiketrain.waveforms.mean(axis=0)
    signal = mean_waveform.max() - mean_waveform.min()
    SD = (spiketrain.waveforms - mean_waveform).std()
    return signal / (2. * SD)