Обратите внимание, что у вас есть 2 канала, что означает, что вы получаете 2-мерный data
. Ваша текущая версия просто выравнивает массив во время некоторых операций, поэтому в нем слишком много элементов.
Есть 2 способа решения проблемы. Первый - использовать только один из каналов:
def spectral_properties(filename):
fs, data = wavfile.read(filename)
# use the first channel only
if data.ndim > 1:
data = data[:, 0]
spec = np.abs(np.fft.rfft(data))
freq = np.fft.rfftfreq(len(data), d=1/fs)
assert len(spec) == len(freq)
amp = spec / spec.sum()
amp_cumsum = amp.cumsum()
assert len(amp_cumsum) == len(freq)
q25 = freq[len(amp_cumsum[amp_cumsum < 0.25])]
q75 = freq[len(amp_cumsum[amp_cumsum < 0.75])]
return (freq * amp).sum(), q25, q75
avg, q25, q75 = spectral_properties('foobar.wav')
print(avg, q25, q75)
Второй - сохранить каналы и указать функции numpy, вдоль оси которых они должны работать. Это также означает, что вычисление квартилей становится менее тривиальным, так как вам нужно найти их отдельно для каждого канала, но это выглядит почти так же просто, как и раньше, из-за понимания списка Python:
def spectral_properties(filename):
fs, data = wavfile.read(filename)
# determine number of channels
num_channels = data.shape[1]
spec = np.abs(np.fft.rfft(data, axis=0))
freq = np.fft.rfftfreq(len(data), d=1/fs)
assert len(spec) == len(freq)
amp = spec / spec.sum(axis=0)
amp_cumsum = amp.cumsum(axis=0)
assert len(amp_cumsum) == len(freq)
q25 = [freq[len(amp_cumsum[:,j][amp_cumsum[:,j] < 0.25])] for j in range(num_channels)]
q75 = [freq[len(amp_cumsum[:,j][amp_cumsum[:,j] < 0.75])] for j in range(num_channels)]
return (freq[:,np.newaxis] * amp).sum(axis=0), q25, q75
avg, q25, q75 = spectral_properties('foobar.wav')
print(avg, q25, q75)
Обратите внимание, что + 1
было проблематично в ваших исходных выражениях для поиска квартилей. Считайте, что все значения меньше 0.25
, кроме последнего. Таким образом, неравенство будет справедливо для n - 1
элементов. Вы add 1
, поэтому вы получаете n
. Но n
слишком высокий показатель для массива freq
длины n
.
Кроме того, я подозреваю, что вы, возможно, захотите поставить квадрат spec
вместо того, чтобы оставить его в величинах.
Обновление:
Вы также можете использовать searchsorted
, чтобы найти квартили, которые должны быть быстрее и легче для чтения:
q25 = freq[np.searchsorted(amp_cumsum, 0.25)]
q75 = freq[np.searchsorted(amp_cumsum, 0.75)]
И
q25 = [freq[np.searchsorted(amp_cumsum[:,j], 0.25)] for j in range(num_channels)]
q75 = [freq[np.searchsorted(amp_cumsum[:,j], 0.75)] for j in range(num_channels)]