Реализация алгоритма Harmonic Product Spectrum в java

В настоящее время я работаю над программой гитарного тюнера на Java и пытаюсь реализовать алгоритм гармонического спектра продукта, чтобы определить основную частоту.

На данный момент я сделал метод, который снижает частоту дискретизации моего спектра с коэффициентом f.

Теперь я пытаюсь умножить все мои разные спектры с пониженной частотой дискретизации вместе. Я кодирую в java и работаю с массивами. Следовательно, у меня есть массивы для индексов, которые были подвергнуты пониженной дискретизации, и массивы для различных значений моего спектра. Сейчас я пытаюсь привести все разные массивы к одинаковому размеру и упорядочить их значения так, чтобы они соответствовали правильному индексу с пониженной дискретизацией. У меня много проблем с размером и значениями....

Более того, я пытаюсь реализовать этот алгоритм на основе примера, который у меня есть на листе бумаги... Следовательно, я могу реализовать этот алгоритм только с 4 спектрами с пониженной частотой дискретизации, но я сомневаюсь, что этого будет достаточно, когда я действительно буду использовать реальный звуковой сигнал. .

Вот код метода даунсэмплинга:

import org.jfree.ui.RefineryUtilities;

public class SousEchantillonnage {

public static double[] valuesDownSample(double[] tab, int factor){

    int N = tab.length;

    double[] values = new double[N];

    for (int i = 0; i < N; i++){
        values[i] = tab[i];
    }

    int lengthDownSample = N + (facteur - 1) * (N - 1);

    double[] valuesDownSample = new double[lengthDownSample];
    for (int i = 0; i < N; i++){
        valuesDownSample[i] = values[i];
    }
    for (int i = N; i < lengthDownSample; i ++){
        valuesDownSample[i] = 0;
    }

    return valuesDownSample;
}

public static double[] indexDownSample(double[] tab, int factor){

    int N = tab.length;

    double[] indexes = new double[N];

    for (int i = 0; i < N; i++){
        indexes[i] = i;
    }

    int lengthDownSample = N + (factor - 1) * (N - 1);

    double[] indexDownSample = new double [lengthDownSample];
    for (int i = 0; i < lengthDownSample; i++){
        indexDownSample[i] = i / factor;
    }

    return indexDownSample;
}

Этот метод, кажется, работает.

Вот мой метод для алгоритма HPS:

public static double[] hps(double[] tab){

    int N = tab.length;

    int factor2 = 2;
    int factor3 = 3;
    int factor4 = 4;
    int lengthDownSample2 = N/2 + (factor2 - 1) * (N/2 - 1);
    int lengthDownSample3 = N/2 + (factor3 - 1) * (N/2 - 1);
    int lengthDownSample4 = N/2 + (factor4 - 1) * (N/2 - 1);

            // Gives us the spectrogram of the signal tab
    double[] spectrogram = new double[N];                       
    spectrogramme = FFT.calculFFT(tab);

            // We only need the first values of the spectrogram. The other half is the same.
    double[] spectrogramCut = new double[N/2];          
    for (int i = 0; i < N/2; i++){
        spectrogramCut[i] = spectrogram[i];
    }

            // We create the array that contains the values of spectrogramCut that we downsample by a factor 2
       double[] valuesSpect2 = new double [sizeDownSamp2];  
       valuesSpect2 = SousEchantillonnage.valuesDownSample(spectrogramCut, factor2);

          // We create an array of the indexes of spectrogramCut that we downsample by a factor 2
      double[] indexSpect2 = new double[sizeDownSamp2];
      indexSpect2 = SousEchantillonnage.indexDownSample(spectrogramCut, factor2);

         // We create the array that contains the values of spectrogramCut that we downsample by a factor 3
       double[] valuesSpect3 = new double [sizeDownSamp3];  
       valuesSpect3 = SousEchantillonnage.valuesDownSample(spectrogramCut, factor3);

         // We create an array of the indexes of spectrogramCut that we downsample by a factor 3
    double[] indexSpect3 = new double[sizeDownSamp3];
    indexSpect3 = SousEchantillonnage.indexDownSample(spectrogramCut, factor3);;

       // We create the array that contains the values of spectrogramCut that we            downsample by a factor 4
   double[] valuesSpect4 = new double [sizeDownSamp4];  
   valuesSpect4 = SousEchantillonnage.valuesDownSample(spectrogramCut, factor4);

       // We create an array of the indexes of spectrogramCut that we downsample by a factor 4
    double[] indexSpect4 = new double[sizeDownSamp4];
    indexSpect4 = SousEchantillonnage.indexDownSample(spectrogramCut, factor4);

        int sizeIndex = N/2 + 5 * (N/2 - 1); // size of the array that contains all the       indexes of the downsamples

        // We create this array
    double[] indexDowSamp = new double[sizeIndex];
    indexDowSamp[0] = 0;
    indexDowSamp[1] = 0.25;
    indexDowSamp[2] = 0.333;
    indexDowSamp[3] = 0.5;
    indexDowSamp[4] = 0.666;
    indexDowSamp[5] = 0.75;

    int q = sizeIndex / 6;      // quantity of packets of 6 we can do
    int r = sizeIndex%6;        // what we are left with.

    for (int i = 6; i < q * 6; i += 6){
        for (int j = 0; j < 6; j++){
        indexDowSamp[i + j] = indexDowSamp[i + j - 6] + 1;
        }
    }
    for (int i = 0; i < r; i++){
        indexDowSamp[q * 6 + i] = indexDowSamp[q * 6 + i - 6] + 1;
    }

Я застрял здесь. Я хотел бы сделать метод, который умножает два массива разной длины вместе.

По сути, как вы можете видеть, когда я понижаю частоту спектрограммы, я получаю два массива:

  • тот, у которого есть индексы, которые были даунсемплингами
  • другой, который имеет значения после понижения частоты дискретизации.

Я хотел бы создать массив того же размера, что и indexDownSample: valuesDownSample. Например, у нас есть indexDownSample[0] = 0. Я хотел бы иметь для valuesDownSample[0] произведение valuesSpectCut[0] *valuesSpect2[0]*valuesSpect3[0]*valuesSpect4[0], потому что все эти массивы имеют значение, соответствующее индексу 0 (indexSpectCut[0] = 0, indexSpect2[0] = 0 = indexSpect3[0] = indexSpect4[0]).

для indexDownSample[1]=0.25 мы замечаем, что только indexSpect4[1]= indexDownSample[1] = 0.25 у нас будет по умолчанию 0 для valuesDownSample[1].

И так продолжаем до тех пор, пока не заполним массив.

Если все пройдет гладко, в конце мы должны иметь:

  • valuesDownSample, который содержит различные значения продуктов
  • indexDownSample, который содержит другой индекс с пониженной дискретизацией.

Мне просто нужно найти максимальный пик, чтобы найти мою основную частоту.

Моя единственная проблема в том, что я понятия не имею, как сделать умножение.....

Если у кого-то есть идея, буду очень признателен!


person fireangel3000    schedule 05.05.2012    source источник
comment
Было бы очень полезно, если бы вы описали, с чем у вас на самом деле были проблемы, и показали хотя бы часть того, что вы сделали до сих пор. Как бы то ни было, сейчас мы должны угадать, что предложить…   -  person Donal Fellows    schedule 05.05.2012
comment
почему ява?!? попробуй с++, фортран, что-нибудь еще!   -  person L7ColWinters    schedule 05.05.2012
comment
Спасибо вам обоим за ваши ответы. Мой проект должен быть в java. Это не выбор... Это было навязано. В то время как для кода до сих пор это то, что я сделал для понижения дискретизации. Извините, хотя мои аннотации будут на французском языке.....   -  person fireangel3000    schedule 05.05.2012


Ответы (1)


Хорошо, у меня есть фактический ответ, вот как я решил проблему:

Прежде всего, я не называю это понижением частоты дискретизации, а разделением частот БПФ. Для этого я использую только 2 массива, например:

float findex[1000]; //index of frequencies for harmonic product spectrum
float mindex[1000]; //index of magnitudes for harmonic product spectrum
unsigned int max_findex;    //number of elements in findex[] and mindex[]

Это моя инициализация: max_findex = 0; for (i = 0; i ‹ 1000; i++) { mindex[i] = 0,0; }

Я разрабатываю на плате Analog Devices Sharc EZ Development Board в VisualDSP++, поэтому мой код написан на C.

Я использую большой размер БПФ 131072 при низкой частоте дискретизации 12 кГц. Это дает мне относительно высокую точность 0,09 Гц на бин. Поскольку более высокие гармоники будут более точными, я сначала начну с высшего деления:

//first run with max division (highest accuracy)
for (i = 0; i < new_peak; i++) {
    findex[max_findex] = f(peak[i].index) / 9;
    mindex[max_findex] = peak[i].magnitude;
    if (max_findex < 999) max_findex++;
    else xmitUARTmessage("ERROR max_findex\r\n", 100);
}

У меня есть все мои пики из БПФ в виде такой структуры: пик [].index — это номер ячейки в БПФ, а пик [].величина — величина пика. Функция f() возвращает частоту бина.

Затем я бы перешел к разделу 8, затем к 7 и т. д. и к разделу 1 (фактические частоты идут последними). Для каждого разделенного пика я просматриваю свой массив, чтобы увидеть, есть ли у меня уже частота в этой точке. Я использую +/- 0,2, и я, вероятно, мог бы ужесточить это, но вы должны настроить его на точность вашего БПФ.

char found;

for (u2 = 8; u2 > 0; u2--) {
    for (i = 0; i < new_peak; i++) {
        tempf = f(peak[i].index) / u2;
        found = 0;
        for (u = 0; u < max_findex; u++) {
            //try to find existing frequency
            if (tempf < findex[u] + 0.2 && tempf > findex[u] - 0.2) {
                //found existing frequency
                mindex[u] *= peak[i].magnitude;
                found = 1;
                break;
            }
        } //for u
        if (!found) {
            //no frequency was found, add new one
            findex[max_findex] = tempf;
            mindex[max_findex] = peak[i].magnitude;
            if (max_findex < 999) max_findex++;
            else xmitUARTmessage("ERROR max_findex\r\n", 100);
        }
    } //for i
} //for u2

И это все. Теперь я просто распечатываю свои значения и сортирую их по величине в Excel...

for (i = 0; i < max_findex; i++) {
    sprintf(tempstr, "%.2f,%.2f\r\n", findex[i], mindex[i]);
    xmitUARTmessage(tempstr, 100);
}

Вот некоторые результаты (я показываю только верхние 6): 60,53 1705693250000 60 1558419875000 20 555159950 179,99 264981525 7,5 1317353 8,57 1317353

Входной звуковой сигнал, создающий выходной сигнал, был следующим: прямоугольная волна 60 Гц (основная плюс нечетные гармоники) синусоида 181,5 Гц синусоида 302,5 Гц синусоида 423,5 Гц

Я имитировал частоту ударов 60 + 60,5 Гц, при этом основная частота полностью отсутствовала на частоте 60,5. Это сложный случай, и он сработал :)

person Peter    schedule 01.11.2013