5. Müzik Sinyalleri frekans spektrum analizi

5.1. Spektrum, frekans uzayının parçalara bölünmesi, nota frekansları, ezgi ve tonalite

Bu defterimizden itibaren müzik sinyal işleme konularını ele almaya başlayacağız. Bir önceki defterde sinyal spektrumunu ele almıştık. Burada bıraktığımız yerden başlayıp adım adım önce müzik kayıtları içerisinde ezgi boyutunu, daha sonra (ilerleyen defterlerde) ritim boyutunu analiz etmek için kodlar yazacağız.

Bir önceki defterimizde bir sinyal içerisindeki sinüzoidal bir bileşenin spektrumda (üzerine pencere fonksiyonunun spektrumu kopyalanmış olarak) öne çıktığını ve tepe noktasının pozisyonundan bileşenin frekansını kestirebileceğimizi görmüştük. Zaman-frekans gösterimi olan spektrogramda da harmonik izlerinin kayıttaki müzikal ezgi ile doğrudan bağlantılı olduğunu görmüştük. Aşağıda bir müzik sinyali için spektrogram çizdirip bu konuyu daha detaylı ele almaya başlayalım.

Gerekli kütüphaneleri yükleyelim.

import sys
import os
import numpy as np
import matplotlib.pyplot as plt
import soundfile as sf
from scipy import signal
import urllib.request
from IPython.display import Audio

Bir saksafon kaydını indirip incelemeye başlayalım.

url = 'https://github.com/MTG/sms-tools/raw/master/sounds/sax-phrase-short.wav'
urllib.request.urlretrieve(url,'sax-phrase-short.wav')
muzik_sinyali, ornekleme_fr = sf.read('sax-phrase-short.wav')

fig = plt.figure(figsize=(15,3))
plt.plot(muzik_sinyali)
plt.title('saksafon sinyali örneği');
plt.xlabel('zaman endeksi (n)');
plt.ylabel('x[n]');
Audio(muzik_sinyali, rate=ornekleme_fr)

Şekil 5.1: Örnek müzik sinyali dalga formu

Spektrum/spektrogram çizimleri oluşturmak için daha önceki defterlerde kullanıdığımız fonksiyonları kullanacağız. Aşağıya o fonksiyonları kopyalayalım.

# Genlik spektrumu hesaplayan fonksiyon 
def genlik_spektrumu(x, fft_N):
  '''Verilen x sinyalinin dB skalasında genlik spektrumunu hesaplar ve 
  spektrumun ilk yarısını([0, örnekleme-frekansı/2] aralığını) döndürür

  Parametreler
  ----------
  x : numpy.array
      Giriş sinyali
  fft_N : int
      Spektrum hesabında kullanılan nokta sayısı

  Döndürülenler
  -------
  numpy.array
      Genlik spektrumu (pozitif frekans bölgesi)
  '''
  X = np.fft.fft(x, fft_N)
  genlikX = np.abs(X[:int(fft_N / 2)])
  genlikX[genlikX < np.finfo(float).eps] = np.finfo(float).eps # log operasyonundan önce önlem
  return 20 * np.log10(genlikX)

def spektrogram_cizdir(x, pencereFonksiyonAdi, pencere_genisligi, pencere_kaydirma_miktari, fft_N, ornekleme_fr, eksen_limitleri={}):
  '''Sinyalin küçük kesitlerinin spektrumlarını genlik_spektrumu() fonksiyonu ile hesaplayıp 
  yanyana dizerek spektrogram zaman-frekans temsilini hesaplar ve çizdirir'''
  pencere_fonksiyonu = signal.get_window(pencereFonksiyonAdi, pencere_genisligi)

  # Pencerelerin baslangic noktalarını hesaplayalım
  pencere_baslangiclari = np.arange(0, x.size - pencere_genisligi, pencere_kaydirma_miktari, dtype = int)
  pencereSayisi = pencere_baslangiclari.size

  spktrgrm = np.array([]).reshape(0, int(fft_N / 2)) # spektogram boş bir dizi olarak başlatılıyor
  for k in range(pencereSayisi):
      pencere_baslangici = pencere_baslangiclari[k]
      x_w = x[pencere_baslangici:pencere_baslangici + pencere_genisligi] * pencere_fonksiyonu
      spktrgrm = np.vstack((spktrgrm, genlik_spektrumu(x_w, fft_N))) # k. sinyal kesitinin spektrumunun eklenmesi

  # Spektrogram matrisinin çizdirilmesi
  zamanEkseni = np.arange(spktrgrm.shape[0]) * pencere_kaydirma_miktari / float(ornekleme_fr)
  frekansEkseniHz = np.arange(spktrgrm.shape[1]) * float(ornekleme_fr) / float(fft_N)

  fig = plt.figure(figsize=(15,4))
  plt.pcolormesh(zamanEkseni, frekansEkseniHz, np.transpose(spktrgrm))
  plt.xlim([0, zamanEkseni[-1]])
  plt.ylim([0, frekansEkseniHz[-1]])
  plt.title('Genlik Spektrogramı (dB)')
  plt.ylabel('frekans(Hz)')
  plt.xlabel('zaman(saniye)')
  if 'x' in eksen_limitleri:
    plt.xlim(eksen_limitleri['x'])
  if 'y' in eksen_limitleri:
    plt.ylim(eksen_limitleri['y'])
  plt.show();

Bir önceki defterde dar-bant spektrogramında harmonikler için yatay çizgiler gözlediğimizi ve ilk harmoniğin izinin ezginin frekansının zamanla değişim eğrisi olduğunu da görmüştük. Bu örnek için bu gözlemi tekrar yapalım.

# Dar-bant spektrumu çizdirelim (zamanda kullandığımız pencere genişliği çok küçük olmamalı ki harmonik izlerini gözleyebilelim) 
eksen_limitleri = {'y':(0, 2000)} # Hz cinsinden eksen sınırları
spektrogram_cizdir(muzik_sinyali, pencereFonksiyonAdi='blackman', pencere_genisligi=2048, pencere_kaydirma_miktari=1024, fft_N=2048, ornekleme_fr=ornekleme_fr, eksen_limitleri=eksen_limitleri)

Şekil 5.2: Yakınlaştırılmış dar-bant genlik spektrogramı

Kaydı dinlediğinizde saksafon ezgisinin 6 notadan oluştuğunu işiteceksiniz. Spektrogramımızda 500Hz çevresinde gözlediğimiz kesikli yatay çizgiler müzikal ezgimizin frekansının zamanla değişimini gösteren bir çizim olarak düşünülebilir. Örneği kaydın ortalarına (1.5 saniye) denk gelen uzun notanın bir kısmı için genlik spektrogramını çizdirip inceleyelim. O bölgede 600 Hz civarında bir yatay çizgi gözlüyoruz. Bu bölgenin spektrumunu incelemek istiyoruz.

# 1.5 saniye çevresinde 0.5 saniyelik bir kesit alalım
baslangic_endeks = int((1.5 - 0.25) * ornekleme_fr)
bitis_endeks = int((1.5 + 0.25) * ornekleme_fr)
sinyal_kesidi = muzik_sinyali[baslangic_endeks:bitis_endeks]
sinyal_kesidi = sinyal_kesidi * signal.get_window('hann', sinyal_kesidi.size)

fft_N = ornekleme_fr # fft_N tercihimizi çizimimizde endeksler Hz cinsinden frekansa denk gelecek şekilde yaptık
genlik_spektrum = genlik_spektrumu(sinyal_kesidi, fft_N)

fig = plt.figure(figsize=(16,3))
plt.subplot(1,3,1)
plt.plot(sinyal_kesidi)
plt.title('Zaman uzayındaki pencerelenmiş sinyal kesidi');
plt.xlabel('Zaman (endeksi n)')
plt.subplot(1,3,2)
plt.title('Genlik spektrumu');
plt.plot(genlik_spektrum)
plt.xlabel('Frekans(Hz)')
plt.ylabel('Genlik(dB)')
plt.subplot(1,3,3)
plt.plot(genlik_spektrum[:2000]) # yakınlaştırılmış çizim için ilk 2000Hz'lik spektrumu alalım
plt.xlabel('Frekans(Hz)')
plt.ylabel('Genlik(dB)')
plt.title('Genlik spektrumu (yakınlaştırılmış)');

Şekil 5.3: Pencerelenmiş sinyal dalga formu ve genlik spektrumu

Orta bölgeden aldığımız kesitin spektrumunda sinyalin temel titreşim frekansını ve harmoniklerin frekansını tepeler olarak gözleyebiliyoruz. İlk tepeyi 600 Hz civarında gözlüyoruz (temel titreşim frekansına karşılık geliyor). Bu frekansın tam katlarında (\(k*600\) Hz) da diğer harmoniklerin tepelerini görebiliyoruz. Harmonik izleri spektrogramda (Şekil 5.2) (yine tam katlarda) yatay izler olarak gözleniyor. Bu durumda spektrumdan yola çıkarak birçok ezgi analizi işlemi gerçekleştirebiliriz. Aşağıda bu yönde adımlar atacağız. Ancak öncelikle frekans uzayı ile ilgili önemli konuyu daha ele almalıyız: insanın diğer birçok algısı gibi frekans algısı da logaritmiktir (Bakınız). Müzikal notaların frekansları frekans uzayı logaritmik bölünerek yerleştirilmiştir. Bir işitme testi planlayıp küçük bir demo yapalım.

5.1.1. Müzikal frekans uzayının bölünmesi

İşitsel test: frekans uzayının doğrusal ve logaritmik bölünmesi: Öncelikle ses sentezlemek için bir fonksiyon tanımlayacağız, daha sonra doğrusal ve logaritmik bölünme ile elde edilen frekanslarda ses sentezleyip dinleyeceğiz.

Verilen bir dizi temel titreşim frekans değerine sahip sinüzoidi arka arkaya dizerek ses sentezleyen bir fonksiyon tanımlayalım. Temel titreşim frekans dizimizi doğrusal ve logaritmik oluşturup elde ettiğimiz sesleri dinleyip karşılaştıracağız.

def sinus_dizisi_sentezle(frekans_dizisi, ornekleme_frekansi=16000,sure_saniye=0.75):
  ''' frekans_dizisi'nde verilen frekanslarda sinüzoidler sentezleyen fonksiyon'''
  t = np.arange(0, sure_saniye, 1/ornekleme_frekansi)
  w = signal.tukey(t.size, 0.05) # sinyal kesitleri arası geçişlerin keskin olmaması için her kesiti tukey penceresi ile çarpacağız
  sinyal = np.array([]) # sinyali boş oluşturup sağına sinyal parçaları ekleyeceğiz
  for frekans in frekans_dizisi:
    sinyal = np.concatenate((sinyal, np.cos(2 * np.pi * frekans * t) * w)) # bir sinüzoid oluşturup sinyalin devamına ekleniyor
  return sinyal

Şimdi frekans serilerimizi oluşturalım. Frekans uzayımızı doğrusal ve logaritmik parçalara bölüp frekans değerleri seçerek iki dizi oluşturacağız.

# Frekans bandı sınırları ve kaç parçaya bölüneceğini keyfi olarak seçelim
baslangic_frekansi = 220 # Hz cinsinden başlangıç frekans değeri
bitis_frekansi = 440 # Hz cinsinden bitiş frekans değeri
parca_sayisi = 12

# Doğrusal bölme: 13 adet frekans değeri aralık doğrusal olarak bölünerek elde ediliyor
delta = (bitis_frekansi - baslangic_frekansi) / parca_sayisi
frekanslar_dogrusal_bolunmus = np.arange(baslangic_frekansi, bitis_frekansi+delta/2, delta)
print("Doğrusal eşit bölünme ile seçilen frekanslar (Hz):\n", frekanslar_dogrusal_bolunmus)

# Logaritmik bölme: 13 adet frekans değeri aralık logaritmik olarak bölünerek elde ediliyor
# Frekans uzayını log2 döşümüne tabi tutup eşit parçalara ayırıp değerler seçelim 
delta = (np.log2(bitis_frekansi) - np.log2(baslangic_frekansi)) / parca_sayisi
frekanslar_log2Hz = np.arange(np.log2(baslangic_frekansi), np.log2(bitis_frekansi)+delta/2, delta)
# log2-Hz frekans değerlerimizi tekrar Hz'e döndürelim
frekanslar_log_bolunmus = np.power(2, frekanslar_log2Hz) 
print("\nLogaritmik eşit bölünme ile seçilen frekanslar (Hz):\n", frekanslar_log_bolunmus)

Doğrusal eşit bölünme ile seçilen frekanslar (Hz):
 [220.         238.33333333 256.66666667 275.         293.33333333
 311.66666667 330.         348.33333333 366.66666667 385.
 403.33333333 421.66666667 440.        ]

Logaritmik eşit bölünme ile seçilen frekanslar (Hz):
 [220.         233.08188076 246.94165063 261.6255653  277.18263098
 293.66476792 311.12698372 329.62755691 349.22823143 369.99442271
 391.99543598 415.30469758 440.        ]

Şimdi bu frekans dizileri için ses sentezleyip dinleyelim. Dinlerken kendimize şu soruyu sormalıyız: hangi dizide sesler birbirini takip ederken frekans değişimleri eşit miktarda artmış hissi veriyor? (Birden fazla kez dinleyiniz.)

# Doğrusal bölünme ile elde edilen frekans dizisinin sentezi
ses_dogrusal_bolunme = sinus_dizisi_sentezle(frekanslar_dogrusal_bolunmus)
Audio(ses_dogrusal_bolunme, rate=16000)

# Logaritmik bölünme ile elde edilen frekans dizisinin sentezi
ses_logaritmik_bolunme = sinus_dizisi_sentezle(frekanslar_log_bolunmus)
Audio(ses_logaritmik_bolunme, rate=16000)