Elektronik / Elektronik Kaynakları/

Teoriden Gerçeğe: Goertzel Algoritması (Vol. I)

Sponsorlu Bağlantılar

Ülkemizde özellikle ‘kopyala-yapıştır programcılık’ mantığıyla çalışan kişi sayısının artmasından bıkan biri olarak, özellikle bir işin nasıl olması gerektiğine kendimce bir yanıt vermeyi amaçladığım bu yazı dizimde, final döneminde rastladığım, kısa notlar ve örneklerle pekiştirdiğim Goertzel algoritmasına değineceğim.

Kısaca değinmek gerekirse, Dr. Gerald Goertzel’in 1958 yılında yayınladığı Goertzel algoritması sinyalin içerisindeki frekans bileşenlerini çözen dijital sinyal işleme tekniğidir. Temel FFT algoritması tanımlı bir bant aralığındaki frekans bileşenlerini bizlere sunarken Goertzel algoritması sadece belirlenen, spesifik frekans bileşenlerini bizlere sunmaktadır.

Algoritmanın matematiksel ifadesi aşağıdadır:

Goertzel-Algoritmasi-Teorik

Goertzel-Algoritmasi-Teorik-2

a) Örnekleme Oranı

Goertzel algoritmasında, aynen FFT’de olduğu gibi, örnekleme frekansı, belirlenmesi istenen frekansların en büyüğünün en az iki katı olmalıdır.

Goertzel-Algoritmasi-Teorik-ornekleme-frekans-3

b) Blok Boyutu

Goertzel blok boyutu N aynı FFT’de olduğu gibi filtre sonucunda çıkacak frekans çözünürlüğüne etki eder. Blok boyutunu N ile temsil edersek frekans çözünürlüğü aralığı fs/N olacaktır. Örneğin örnekleme frekansı fs=8kHz, filtre uzunluğu N=100 seçilirse sonuçta çıkacak frekans çözünürlüğü 80Hz olacaktır.

Blok boyutunu büyük seçmek, özellikle düşük frekans bileşenlerinin daha iyi kararlılıkla belirlenmesine yol açsa da büyük seçilen blok boyutu Goertzel algoritması hesaplamasını bir o kadar uzatacaktır. Dolayısıyla blok boyutu seçiminde bu iki kriter göz önüne alınarak optimum seçimin yapılması gerekmektedir. Son olarak seçilecek olan N değeri, FFT’nin aksine, 2’nin katları olmak zorunda değildir.

c) Sabitler

Goertzel algoritmasını hesaplamadan önce Goertzel algoritmasının kullanacağı katsayıları hesaplamak işlemler sırasında algoritma çözümüne hız kazandıracaktır.

Goertzel-Algoritmasi-Teorik-4

Yukarıdaki katsayılar her bir tespit edilecek frekans için ayrı ayrı bulunur. Yukarıdaki N sayısı blok boyutunu temsil etmektedir.

d) Temel Goertzel

Temel Goertzel algoritması DFT ve FFT’de olduğu gibi geriye sanal frekans bileşenlerini sunar. Genlik ve faz bilgileri geriye dönen bu sanal bileşenlerden yola çıkılarak hesaplanabilir.

Geortzel algoritması hesabında Q0, Q1 ve Q2’den oluşan 3 değişken sürekli hesaplanmak zorundadır.
Bu değişkenlerin başlangıç değerleri sıfır olmak koşuluyla Goertzel algoritması aşağıdaki gibi yazılabilir.

Goertzel-Algoritmasi-Teorik-5

Yukarıdaki denklemler N kez tekrarlanır.
Daha sonra ortaya çıkan Q0, Q1 ve Q2 sayıları aşağıdaki denklemlerde kullanılır.

Goertzel-Algoritmasi-Teorik-6

Yine bu denklemler belirlenecek frekans sayısı kadar tekrarlanarak her bir frekans için genlik değeri bulunmuş olur.

e) Optimize Edilmiş Goertzel

Optimize edilmiş Goertzel algoritması faz açıları yerine sadece genliklere odaklanır ve temel Goertzel algoritmasından daha az işlem yapılmasını sağlar.

Q0, Q1 ve Q2 katsayıları hesaplandıktan sonra temel Goertzel’den farklı olarak yapılması gereken tek işlem;

Goertzel-Algoritmasi-Teorik-7

Mikrodenetleyiciler ile bu işlem yapılırken genliği büyük olan işleme alınacaksa karekök alma işlemine gerek yoktur. Dolayısıyla buradan denklem aşağıdaki forma dönüşür.

Goertzel-Algoritmasi-Teorik-8

f) FFT ve Goertzel Algoritmasının Karşılaştırılması

Yüksek sayıda örnek için FFT her ne kadar avantajlı gözükse de Goertzel düşük sayıda frekansı çözümlemek için FFT’den daha hızlıdır. Genel anlamda N blok boyutu için FFT N(log2)N Goertzel (5/6)(log2)N kadar işlem yapmaktadır.

g) Matlab Uygulaması

Goertzel algoritması için her ne kadar Matlab kendi fonksiyonunu barındırsa da yukarıdaki denklemleri kullanarak yazdığımız .m kodlarını aşağıda görebilirsiniz. Aşağıdaki kodlardan da görülebileceği üzere DTMF keypad tuşlarından 3 örneklenmiştir.

% Yazar: Fırat DEVECİ
% Konu : Goertzel Algoritması

clear all;
Fs   = 4000;           % Örnekleme frekansı
n=0:1/Fs:0.1;           % Örnek
Goertzel_N=length(n);   % Goertzel örnek sayısı


freq=[697, 770, 852, 941, 1209, 1336, 1477, 1633];  % Bulmak istediğimiz frekans bileşenleri
Bins=length(freq);                                  % Kaç adet frekans bulacağımız bilgisi
samples=sin(2*pi*697*n)+sin(2*pi*1477*n);           % Örneklediğimiz sinyal
sound(samples);                                     % Ses çalınıyor

for i=1:Bins                            % Goertzel katsayıları hesaplanıyor
    k(i)=(0.5+(Goertzel_N*freq(i))/Fs);
    coeffs(i)=2*cos(2*pi*k(i)/Goertzel_N)
end

prev1=[zeros(Bins)];
prev2=[zeros(Bins)];

for i=1:Goertzel_N                      % Hangi frekansın baskın olduğu bulunuyor
    for j=1:Bins                        % Goertzel algoritması
        val=coeffs(j)*prev1(j)-prev2(j)+samples(i);
        prev2(j)=prev1(j);
        prev1(j)=val;
    end
end

for i=1:Bins                            % Genlikler bulunuyor
    magnitude(i)=sqrt((prev1(i)*prev1(i)+prev2(i)^2)-(coeffs(i)*prev1(i)*prev2(i)));
end

freq_indices = round(freq/Fs*Goertzel_N)+1;
dft_data = goertzel(samples,freq_indices);

subplot(3,1,1);
plot(n,samples);
title('Sinyal Şekli');
subplot(3,1,2);
stem(freq,abs(magnitude));              % Ekrana basılıyor
set(gca,'xtick',freq); 
xlabel('Hz');
ylabel('Magnitude');
title('Kendi Goertzel Algoritmamızın Yanıtı');
subplot(3,1,3);
stem(freq,abs(dft_data));
set(gca,'xtick',freq); 
xlabel('Hz');
ylabel('Magnitude');
title('Matlab Goertzel Algoritmasının Yanıtı');

Yukarıdaki kod işletildiğinde ekrana aşağıdaki şekil (üzerine tıklarsanız büyür) ortaya çıkmaktadır. Görüldüğü üzere kendi yazdığımız optimum Goertzel algoritması oldukça sağlıklı sonuç vermektedir.

Goertzel-Matlab

h) ANSI C Uygulaması

Yukarıda Matlab için yaptığımız uygulamayı ANSI C ortamına taşırsak aşağıdaki kod öbeğine ulaşmış oluruz.

/*
Yazar: Fırat DEVECİ
Konu : Goerzel Algoritması
*/

#include "stdio.h"
#include "stdlib.h"
#include "math.h"

#define Fs              4000            // Örnekleme frekansımız, izleyeceğimiz en yüksek frekans değerinin en az iki katı olmalı
#define BINS            8               // 8 temel frekans izlenecek
#define Pi              3.141592654     // Pi sayısı tanımlanıyor
#define GOERTZEL_N      128             // 128 örnek alınacak
#define Level_Segment   8

const int   freqs    [BINS] = {697, 770, 852, 941, 1209, 1336, 1477, 1633};
float       samples  [GOERTZEL_N];

float       coeffs   [BINS];
float       prev1    [BINS];
float       prev2    [BINS];
float       magnitude[BINS];
int         levels   [BINS];
float       n        [GOERTZEL_N];

unsigned char sample_counter    = 0;
unsigned char sample_complete   = 0;

// Goertzel için katsayılar hesaplanıyor

void coefficient(void)
{
    unsigned int k, i;
    for(i=0;i < BINS;i++)
    {
        k           =(unsigned int)(0.5 + (float)GOERTZEL_N * freqs[i] / Fs);
        coeffs[i]   = 2.0 * cos(2.0 * Pi * ((float)k/GOERTZEL_N));
    }

}

int goertzel(void)
{
    int big =0,
        i   =0;
    float val = 0,
          max = 0;

    // Geçici dizilerin içerisi temizleniyor
    for(i=0;i < BINS;i++)
    {
        prev2[i]=0;
        prev1[i]=0;
    }

    // GOERTZEL Algoritması işletiliyor
    for(sample_counter=0;sample_counter < GOERTZEL_N; sample_counter++)
    {
        for(i=0;i < BINS;i++)
        {
            val     = coeffs[i] * prev1[i] - prev2[i] + (float)samples[sample_counter];
            prev2[i]= prev1[i];
            prev1[i]= val;
        }
    }

    // Genlikler belirleniyor
    for(i=0;i < BINS;i++)
    {
        magnitude[i] = (prev1[i] * prev1[i]) + (prev2[i] * prev2[i]) - (coeffs[i] * prev1[i] * prev2[i]);
        if(magnitude[i] > max)
        {
            max=magnitude[i];
            big=i;
        }
    }

    for(i=0;i < BINS;i++)
    {
        levels[i]=(int)(Level_Segment*(float)(magnitude[i]/max));
    }

    return big;
}

int main()
{
    int i,buyuk;
    float a;
    coefficient();


    for(i=0;i < GOERTZEL_N;i++)
        samples[i]=sin(2*Pi*697*i/4000)+sin(2*Pi*1477*i/4000);

    buyuk=goertzel();

    printf("Yazar: Firat DEVECI\n");
    printf("Konu : Goertzel Algoritmasi\n\n");
    printf("En Buyuk Genlik: %d\n\n",buyuk);
    printf("Frekanslar      Katsayilar      Goertzel Sonucu\n");

    for(i=0;i < BINS;i++)
    {
        if(i==buyuk)
            printf("%d\tHz -    %f\t -  %d*\n",freqs[i],coeffs[i],levels[i]);
        else
            printf("%d\tHz -    %f\t -  %d\n",freqs[i],coeffs[i],levels[i]);

    }

        getch();
}

Yukarıdaki kodlar herhangi bir C (Code Blocks vb.) derleyicide derlendiğinde aşağıdaki çıktı ekranı ile karşılaşılır.

Goertzel-ANSI-C

Sonuç

Goertzel algoritması yukarıdaki örneklerden de görüleceği üzere anlaşılması kolay, FFT'ye göre oldukça hızlı ve özellikle düşük güçteki mikrodenetleyiciler ile kolaylıkla kullanılabilecek yapıdadır. Özellikle DTMF kod çözümleme uygulamalarında tercih edilmesinin altında bu yatmaktadır.

Gelecek bölüm: Goertzel algoritmasının mikrodenetleyici ile gerçeklenmesi.

  • Gürkan

    Saygı.

  • fatih

    Güzel çalışma olmuş, elinize sağlık. İyi çalışmalar, başarılar.

  • metin

    Selam Fırat;
    Görüyorum ki olayları gerçekten dibine kadar araştırıp öğreniyorsun. Şu anda Türkiye’de gördüğüm sayılı mühendislerdensin. Tebrik ediyorum.
    Peki bu sağlıklı öğrenmeyi iyi bir ingilizce mi? Sağlam bir altyapı mı? Danışabildiğin kişilerin sağlam ve çoğunlukta olması mı? Okulunuzda verilen eğitim mi sağlıyor?
    Selamlar, Kolay gelsin…

  • FxDev

    İyi günler Metin,
    Aslında bir çok konuyu dibine kadar araştırıp maalesef öğrenemiyorum, buna zamanım elvermiyor. Fakat ilgi alanımla ilgili herşeyi elimden geldiğince senin tabirinle ‘dibine kadar’ öğrenmeye çabalıyorum. Çok iyi olmasa da, normal düzeyde bir İngilizce bilgisi her zaman yardımcı oluyor bana. Danışabileceğim kişileri de özenle seçiyorum, çünkü genelde herkes her şeyi bilen ama hiç bir şey bilmeyenlerden. Okulumuzdaki eğitime gelince siteden de takip etmişsen ben memnun değilim. Bir kaç ders dışında ezberciliğe iten bir yapı var. Bunu zaten okul bittiğinde derin bir yazıyla yazmayı düşünüyorum.

    Sana da kolay gelsin, iyi çalışmalar.

  • Muhittin KAPLAN

    Fırat tebrik ve takdir ediyorum.

  • FxDev

    Teşekkürler hocam ;)

  • Erdal

    Merhaba güzel bir çalışma olmuş.
    Bir önerim olacak. Sayfada çok bilinen sosyal ağlarla paylaşım butonları var bunların yanında birçok web sitesinde bulunan “print” butonunun da bulunmasının faydalı olacağını düşünüyorum.

  • FxDev

    Öncelikle teşekkürler.
    Bu bir süredir benim de aklımda, print ve pdf’e çevirme eklentilerini en kısa zamanda koymaya çalışacağım.

  • Kırıkdr

    Çok teşekkürler çok güzel bir yazı olmus Fırat Hocam…

  • ramazan

    merhaba Fırat hocam c# ile microfondan sesi alıp dtmf tonlarını bulmak istiyorum. 8 khz lik sinyalden 512 lik bloklar alıp sinyayalin grafiğini çiziyoruz.

    sizin vermiş olduğunuz koddaki #define GOERTZEL_N 128 değeri yerine GOERTZEL_N =512 yapıyorum.
    samples dizisi yerine de 512 lik bloğu veriyorum. ancak 3 tuşunu hiç bulamıyorum.bununla ilgili bir fikriniz varmı?

    birde bütün tuşları bulmak için nasıl bir yol izlememi tavsiye edersiniz.

  • ayşe

    Çok güzel anlatım teşekkürler öncelikle..
    Goertzel fonksiyonunda hata alıyorum verdiğiniz kodlarda, ne yapabilirim hocam?