T60/RT60, C50/C80, D in Python

gistfile1.py

import librosa
from scipy.signal import butter, lfilter

SR = 8000
NFFT = 256  # change ~proportionally to SR


def butter_bandpass(lowcut, highcut, fs, order=5):
    nyq = 0.5 * fs
    low = lowcut / nyq
    high = highcut / nyq
    b, a = butter(order, [low, high], btype='band')
    return b, a


def butter_bandpass_filter(data, lowcut, highcut, fs, order=5):
    b, a = butter_bandpass(lowcut, highcut, fs, order=order)
    y = lfilter(b, a, data)
    return y


def _rms(sig):
    return librosa.feature.rms(S=librosa.core.stft(
        butter_bandpass_filter(sig, 80, 4000, SR), NFFT))[0]

def RT(sig, n=10):
    """https://de.wikipedia.org/wiki/Nachhallzeit"""
    rms = pd.Series(_rms(sig))
    t = len(sig)/len(rms)
    idxmax = rms.idxmax()
    rmsmax = rms[idxmax]
    rms = rms[idxmax:]
    rtN = rms[rms <= 10 ** (np.log10(rmsmax) - n/20)].index[0]
    return int(t * (rtN - idxmax))

def C(sig, n=50):
    """https://de.wikipedia.org/wiki/Klarheitsma%C3%9F / https://de.wikipedia.org/wiki/Transparenz_(Akustik)"""
    rms = _rms(sig)
    start = max(0, rms.argmax() - 1)
    n_ms = int(len(rms) / len(sig) * n * SR / 1000)
    return 10 * np.log10(rms[start:start+n_ms].sum() / rms[start+n_ms:].sum())

def D(sig):
    """Center time (Schwerpunktzeit)"""
    rms = _rms(sig)
    return (rms * np.arange(0, len(rms))).sum() / rms.sum() * len(sig) / len(rms)