Итак, по многочисленным просьбам в теме о распознавании речи
выделяю сабж в отдельный топик.
Для начала имеет смысл напомнить об общей процедуре классификации, приведенной там же.
Здесь только несколько слов о некоторых моментах практической реализации обсуждаемого. В основном будем придерживаться общей схемы процедуры классификации.
Начнем с главного. Теорию вейвлетов можно почитать в разных источниках. Наиболее понятные мне:
есть на сайте http://www.matlab.ru/wavelet/index.asp (ссылки внизу страницы). Очень толковые материалы здесь: http://www.autex.spb.ru/wavelet/, причем много в pdf.
Обращаю ваше внимание на то, что здесь используем именно вейвлет-пакет.
Ниже идет код функции, реализующей вейвлет-преобразование сигнала, записанного в массив sg. Для чего нужна переменная flag, скажу позже, там станет ясно. Физический смысл используемого преобразования состоит в том, что для разложения исходного сигнала по уровням используется его прореживание во временной области (для простоты скажем, разбиение на четные и нечетные отсчеты) плюс некоторые дополнительные преобразования. В частотной области это (опять же упрощенно) равносильно разбиению спектра на две половины. Применямый здесь вейвлет вряд ли вы где-либо найдете . разве что только встретите моего коллегу, Сергея Агиевича (фамилия такая), с которым мы вместе работали. Он придумал этот вейвлет. Говорит, что построил его на основе слайна! Честно скажу, сам я глубоко эту тему не ковырял. Поэтому извиняйте, если не смогу на вопросы ответить. Но эта штука вполне работает, можете не сомневаться.

Код

(с) Агиевич С.Н., 1999-2000. void wavelet(complex<double> *sg, int im, int flag) {    // Вейвлет   int p1=log10l((long double)im)/log10l((long double)2);   // p1 - степень двойки                                                            // т.е.  int im=pow(2,p1);   double a = 2*M_PI/(im/2);   double b1;   double y;   double *U7 = new double[(int)(im/2)];   double *beta = new double[(int)(im/2)];   complex<double> *sg1 = new complex<double>[im+1]; // массив входных данных   complex<double> *sg2 = new complex<double>[im+1]; // массив входных данных   complex<double> *sigS = new complex<double>[im+1]; // массив входных данных   complex<double> *vx = new complex<double>[im+1]; // ядро Фурье-преобразования   complex<double> *vyx = new complex<double>[im+1]; // ядро Фурье-преобразования (обратного)   complex<double> *wx = new complex<double>[im+1]; // ядро преобразования   complex<double> *wyx = new complex<double>[im+1]; // ядро преобразования   complex<double> *sigD = new complex<double>[im+1]; // ядро преобразования       for(int i=0; i < im/2+1; i++)       {         b1 = a*i;         vx[i] = complex<double>(cos(b1),-sin(b1));         vyx[i] = complex<double>(cos(b1),sin(b1));       }       //взятие половины      if(flag == 1) //1024        {         for(int i=0; i < im/2; i++)          {           sg1[i] = sg[i*2];           sg2[i] = sg[i*2+1];          }        }        else        {           for(int i=0; i < im/2; i++)          {           sg1[i] = sg[i];           sg2[i] = sg[i+im/2];          }          flag = 1;        }        fft(sg1, p1-1, vyx, im/2);        fft(sg2, p1-1, vyx, im/2);       for(int i=0; i < im/2+1; i++)       {         b1 = a*i/2;     // для вейвлетов         wx[i] = complex<double>(cos(b1),-sin(b1));         wyx[i] = complex<double>(cos(b1),sin(b1));       }        // вспомогательные величины       for(int i=0; i < im/2; i++)       {         y=cos(M_PI*i/(im/2));         U7[i]=1-pow(1-y,4)*(pow(y,2)-28*y+61)/                 (61+479*pow(y,2)+179*pow(y,4)+pow(y,6));         beta[i]= U7[i]/(1+pow(U7[i],2));       }       for(int i=0; i < im/2; i++)       {         sigD[i]=(sg2[i]-wx[i]*sg1[i]*U7[i])*(double)(im/2);         sigS[i]=(sg1[i]+beta[i]*wyx[i]*sigD[i]/(double)(im/2))*(double)(im/2);       }       // получено разрежение в частотной области        // переход во временную область       fft(sigD, p1-1, vx, im/2);       fft(sigS, p1-1, vx, im/2);       // возврат значений sg !!!      for(int i=0;i < im/2;i++)      {          sg[i] = sigS[i];          sg[i + (int)(im/2)] = sigD[i];      }       delete[] vyx; vyx = 0;       delete[] sg1; sg1 = 0;       delete[] vx; vx = 0;       delete[] sg2; sg2 = 0;       delete[] sigS; sigS = 0;       delete[] wyx; wyx = 0;       delete[] sigD; sigD = 0;       delete[] wx; wx = 0;       delete[] U7; U7 = 0;       delete[] beta; beta = 0; } // конец wavelet (с) Агиевич С.Н., 1999-2000.

Я здесь не привожу код алгоритма БПФ fft(...), т.к. он уже был показан, например, здесь.

To be continued...

З.Ы. Просьба ждать терпеливо... Очень занят, понимаете ли.
И напоминаю, что если вам нравится, как мы общаемся, то можете оценить свое настроение по 10-балльной шкале внизу страницы

Теперь посмотрим, как этот вейвлет-пакет выглядит.
Для начала представим, что у нас в массиве wavesig имеется wavesig_dim отсчетов (в моем варианте не более 1024) сигнала.

Формируем вейвлет-пакет с помощью рассмотренных выше функций.

Код

int level = 8; // число уровней разложения    int flag=1; // вспомогательная // уровни 0...N - вейвлет-разложение // издеваться будем над sperm - копией wavesig (чтобы сохранить wavesig)    for(int i = 0; i < wavesig_dim; i++)         sperm[i] = wavesig[i];    int ss;    float en;    complex<double> *sgw = new complex<double>[wavesig_dim];    complex<double> *sgc = new complex<double>[wavesig_dim];    for(int i=0;i < wavesig_dim;i++) //1024    {        sgc[i] = complex<double>(sperm[i]);    }    for(int l = 0;l < level; l++)     {      if(l != 0)       {        for(int sl = 0; sl < pow(2,l-1); sl++)        {         for(int i=0; i < wavesig_dim/pow(2,l-1); i++) //1024          {           ss = wavesig_dim/(int)pow(2,l-1)*sl + i;   //1024           sgw[i] = sgc[ss];          }           // разложение и возврат в тот же массив           wavelet(sgw, wavesig_dim/(int)(pow(2,l-1)), flag);           for(int i=0;i<wavesig_dim/(int)(pow(2,l));i++) //1024           {            sgc[i+sl*wavesig_dim/(int)pow(2,l-1)] = sgw[i];            sgc[i+wavesig_dim/(int)(pow(2,l))+sl*wavesig_dim/(int)pow(2,l-1)] =                 sgw[i+wavesig_dim/(int)(pow(2,l))];           }          }       }       else       {         en=0;         for(int i=0;i<wavesig_dim/(int)(pow(2,l));i++) //1024          {           en += fabs(sperm[i]);          }        //нормирование!!!         for(int i=0;i < wavesig_dim;i++) //1024          {           sperm[i] /=(en/wavesig_dim);// e[0];           sgc[i] /=(en/wavesig_dim);// e[0];          }       }         for(int i=0;i<wavesig_dim;i++)            serm[l*wavesig_dim+i] = real(sgc[i]);     }

Чтобы посмотреть собственными глазами, выводим это на форму с компонентом TPaintBox

Код

void __fastcall TWaveForm::FormPaint(TObject *Sender) { // формирование сетки WFPaintBox->Canvas->Pen->Style = psSolid;    WFPaintBox->Canvas->Pen->Color = clWhite;    WFPaintBox->Canvas->Brush->Color = clBlack; WFPaintBox->Canvas->Brush->Style = bsSolid;    int LeftGran = WFPaintBox->Width/20;    int HighGran = WFPaintBox->Height/20;    int RightGran = WFPaintBox->Width-10;    int LowGran = WFPaintBox->Height-10; WFPaintBox->Canvas->Rectangle(LeftGran, HighGran, RightGran, LowGran);    WFPaintBox->Canvas->Pen->Color = clWhite; WFPaintBox->Canvas->Pen->Style = psSolid;    for(int i = 1; i < level; i++)    {         WFPaintBox->Canvas->MoveTo(LeftGran, HighGran +                ((LowGran-HighGran)/level)*i);         WFPaintBox->Canvas->LineTo(RightGran, HighGran +                ((LowGran-HighGran)/level)*i);         for(int j = 1;j <= pow(2,i);j++)         {            WFPaintBox->Canvas->MoveTo(LeftGran + (RightGran -                LeftGran)/pow(2,i)*j, HighGran + (LowGran -                HighGran)/level*i);            WFPaintBox->Canvas->LineTo(LeftGran + (RightGran -            LeftGran)/pow(2,i)*j, LowGran);         }    } // нанесение осей    WFPaintBox->Canvas->Pen->Color = clGreen; WFPaintBox->Canvas->Pen->Style = psDot;    for(int i = 0; i < level; i++)    {         WFPaintBox->Canvas->MoveTo(LeftGran, HighGran +             0.5*((LowGran-HighGran)/level) + ((LowGran-HighGran)/level)*i);         WFPaintBox->Canvas->LineTo(RightGran, HighGran +             0.5*((LowGran-HighGran)/level) + ((LowGran-HighGran)/level)*i);    } // заполнение данными:    TPoint points[1024];   WFPaintBox->Canvas->Pen->Style = psSolid;    WFPaintBox->Canvas->Pen->Color = clYellow;    float x_offset, y_offset;    for(int l = 0;l < level; l++)    {      float max_wavesig = 0;      int j = 0;      while(j < wavesig_dim)       {        if(fabs(serm[l*wavesig_dim+j]) > max_wavesig)            max_wavesig = fabs(serm[l*wavesig_dim+j]);        j += 1;       }         for(int i=0; i < wavesig_dim /* /pow(2,l) */; i++)         {           x_offset = ((float)RightGran - (float)LeftGran)/(float)wavesig_dim;           points[i].x = LeftGran + (float)(x_offset*i);           y_offset = (LowGran-HighGran)/level/2/max_wavesig;         points[i].y = HighGran + (LowGran-HighGran)/level/2 +                l*(LowGran-HighGran)/level -                  (float)(y_offset*serm[l*wavesig_dim+i]);         }         WFPaintBox->Canvas->Polyline(points,wavesig_dim - 1);    } }

Все еще to be continued...

Podval, очень интересные статьи.
Будет ли продолжение ?

И вот еще вопрос - не пробовал ли ты сравнивать информативность
коэффициентов вейвлет- и Фурье-спектров
в процедурах обработки реальной речи ?

Цитата

Будет ли продолжение ?

Естественно, я хочу продолжить, только есть один гемморой: хочу выложить картинки, а никак. А то без картинок неинтересно.

Цитата

И вот еще вопрос - не пробовал ли ты сравнивать информативность коэффициентов вейвлет- и Фурье-спектров в процедурах обработки реальной речи ?

Что имеется в виду под информативностью?
Я могу сказать, что вейвлетами вытаскиваются мелкие детали окраски тембра, по которым можно однозначно идентифицировать диктора. Фурье-анализом так не сделать.
А что еще интересует?

Под информативностью я подразумевал вот что.
(Sorry, если я хочу слишком много )

Сколько коэффициентов вейвлет-спектра стоит использовать
в процессах идентификации (по сравнанию с количеством коэффициентов
амплитудного спектра Фурье, дающим где-то сходный результат) ?
Интересен был бы ответ даже на уровне "больше-меньше".

При использовании спектра Фурье для обработки речи
дело обстоит так, что характерная информационность
сосредоточена в определенной полосе частот,
либо, скажем, в первых 12 коэффициентах кепстра.
Локализован ли так же вейвлет-спектр,
либо для каждого диктора/фонемы по-разному ?
Т.е., можно ли на практике не заниматься выделением спектральных
составляющих, максимально различных для разных, хм, образов,
а проводить типа просто корреляцию с набором образцов ?

Далее. Точность прямого расчета Фурье-спектра в таких задачах
часто оказывается недостаточной из-за ограничений на длину выборки,
обусловленных нестационарностью сигнала.
Как обстоят дела у вейвлетов с точностью определения характеристик
при ограниченной длине выборки ?

И еще.
Насколько результаты вейвлет-преобразования, как бы это сказать,
дикторо-зависимы, либо, наоборот, дикторо-независимы.
В сравнении с Фурье, конечно.

Извиняюсь за многословность...

Ну что тут сказать?

Цитата

Сколько коэффициентов вейвлет-спектра стоит использовать

- это зависит от объектов, подлежащих распознаванию. Не могу однозначно ответить. Дело в том, что существуют некоторые методики, позволяющие вычислить количество уровней разложения, оптимальное по какому-то критерию.

Цитата

Сколько коэффициентов вейвлет-спектра стоит использовать в процессах идентификации (по сравнанию с количеством коэффициентов амплитудного спектра Фурье, дающим где-то сходный результат) ? Интересен был бы ответ даже на уровне "больше-меньше".

В моих задачах получалось, что меньше. Но дело тут было еще и в точности. Точность однозначно выше.

Цитата

При использовании спектра Фурье для обработки речи дело обстоит так, что характерная информационность сосредоточена в определенной полосе частот, либо, скажем, в первых 12 коэффициентах кепстра. Локализован ли так же вейвлет-спектр, либо для каждого диктора/фонемы по-разному ? Т.е., можно ли на практике не заниматься выделением спектральных составляющих, максимально различных для разных, хм, образов, а проводить типа просто корреляцию с набором образцов ?

Вейвлет более локализован.

Цитата

Далее. Точность прямого расчета Фурье-спектра в таких задачах часто оказывается недостаточной из-за ограничений на длину выборки, обусловленных нестационарностью сигнала. Как обстоят дела у вейвлетов с точностью определения характеристик при ограниченной длине выборки ?

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

Цитата

Насколько результаты вейвлет-преобразования, как бы это сказать, дикторо-зависимы, либо, наоборот, дикторо-независимы.

Опять же зависит от того, как сформулировать задачу. В той задаче, что я описал, целью является опознать диктора. Если же ставить целью опознавание семантики, то тут и задачу надо формулировать немного по-другому. Я уверен, что это возможно, но я этим не занимался и не встречал пока готовых продуктов, связанных именно с обработкой речи. Видимо, еще не очень популярны вейвлеты, все-таки относительно недавно появились. Даже производители программ и некоторых игрушек типа говорящий словарь-переводчик пользуются Фурье-анализом. Просто вейвлеты только начинают приживаться.
Я зато встречал публикации медиков, где с помощью вейвлетов диагностируют язву желудка по вейвлет-картинке распределения потенциалов и еще много чего. А ведь поставить диагноз - это задача, вообще говоря, "пациенто-независимая"
Так что, главное - правильно формализовать задачу.
Надеюсь, удовлетворил? Если что, уточняй.

Цитата

Извиняюсь за многословность...

Don't be upset about it Мы здесь для того и собрались.

Спасибо, podval .
Вот ты пишешь, что вейвлет локализован не только по частоте,
но и во времени. А как полагаешь, не превратится ли
локализация его во времени в некоторую проблему при обработке
широко распространенным (в узких кругах ) способом,
когда входной сигнал делится на пакеты, подвергаемые уже
дальнейшей обработке. В таком случае ход времени
определяется номером обрабатываемого пакета, а тут еще и информация
о локализации внутри пакета... Я понимаю, конечно,
информация лишней не бывает, но все же - что ее, игнорировать ?

Еще вот, принято ли и имеет ли смысл, как это принято в технике
Фурье-преобразования, использовать оконные функции в WT ?

И как ты полагаешь, как следует переформулировать задачу,
чтобы обеспечить возможную при данном подходе дикторо-независимость ?

Я конечно понимаю, что все это похоже на технический шпионаж


Это сообщение отредактировал(а) Crait - 3.4.2003, 20:44

Цитата

А как полагаешь, не превратится ли локализация его во времени в некоторую проблему при обработке широко распространенным (в узких кругах  ) способом, когда входной сигнал делится на пакеты, подвергаемые уже дальнейшей обработке. В таком случае ход времени определяется номером обрабатываемого пакета, а тут еще и информация о локализации внутри пакета...

А разве это проблема? Сначала режешь на пакеты, а потом анализируешь вейвлетами. Или я чего-то не того?

Цитата

принято ли и имеет ли смысл, как это принято в технике Фурье-преобразования, использовать оконные функции в WT ?

Разные типы вейвлетов имеют разную форму и разные свойства. Использовать что-то дополнительно вроде не имеет смысла.

Цитата

И как ты полагаешь, как следует переформулировать задачу, чтобы обеспечить возможную при данном подходе дикторо-независимость ?

Думаю, надо работать с фонемами. Возможно, по тому же принципу, как и с использованием Фурье-анализа. Только вместо Фурье используешь другой аппарат - вейвлеты. Ну, может, какие-то свои тонкости появятся. Не занимался я такой задачей. Тут думать надо.

Вернемся к нашим баранам. Если все сделано правильно, то вейвлет-пакет должен выглядеть примерно так:



Здесь использовано 8 уровней разложения в вейвлет-пакет. Нулевой уровень (на рисунке он верхний) представляет собой ни что иное, как исходный сигнал. Обратите внимание, что основная энергия сигнала концентрируется в основном в низкочастотных коэффициентах (левые "половинки"). И если внимательно вглядывались в приведенный выше код (не для слабонервных ), то можете заметить, что геометрический смысл разложения в вейвлет-пакет очень нагляден: мы прореживаем исходную выборку (в нашем случае раскладываем отсчеты сигнала на четные и нечетные и уводим их в "левую", т.е. низкочастотную, и "правую", т.е. высокочастотную, половинки) и производим над полученными подпоследовательностями дополнительные преобразования. Обращаю ваше внимание на то, что длительности полученных подпоследовательностей действительно "половинятся", если измерять длительность именно в количестве отсчетов, но по оси времени они не меняются! Я хочу этим сказать, что такая с позволения сказать древовидная структура вейвлет-пакета представлена таковой только для наглядности. На самом деле все длительности во времени полученных подпоследовательностей одинаковы. И чтобы их изобразить действительно одинаковыми, надо было бы на картинке их растянуть до того же размера, что и исходныя выборка сигнала. Но тогда мы увидим, что в подпоследовательностях отсчеты идут реже, но это не удивительно, мы их на самом деле проредили.
Надеюсь, не очень задурил
Посмотрите, вот, как может концентрироваться энергия сигнала в других участках "спектра" (надеюсь, не забыли, почему слово спектр взято в кавычки?)

Цитата

А разве это проблема? Сначала режешь на пакеты, а потом анализируешь вейвлетами. Или я чего-то не того?

Я спросил об этом потому, что повышенная точность
представления функции ограниченным рядом вейвлетов (аппроксимация)
обусловлена, как мне показалось, именно тем,
что используется информация не только об амплитуде,
еще и о временнОм положении вейвлета.
Корректно ли будет ее игнорировать и не потеряется ли при этом
преимущество перед рядами Фурье?
Хотя, с другой стороны, игнорируем же мы фазовую информацию,
анализируя амплитудный Фурье-спектр...

Это сообщение отредактировал(а) Crait - 5.4.2003, 21:58

Цитата

временнОм положении вейвлета

- это здесь не причем!
Более того, мы об этом еще не говорили, но вейвлеты можно масштабировать!
Если встретишь такое выражение - кратно-масштабный анализ - поинтересуйся, что это за вещь.

Траблы с картинками устранил. Смотрите, наслаждайтесь, и скоро поедем дальше

Вернемся к нашей задаче классификации. Наша уель сейчас - посмотреть, как формируется эталонное описание (эталон) на основе вейвлет-пакета. Используем те же функции, что рассмотрены ранее. Только теперь, когда раскладываем выборку сигнала в вейвлет-пакет, будем суммировать энергию отдельных коэффициентов. Это даст нам картинку распределения энергии сигнала. Собственно, мы такую картинку уже наблюдали (см. рисунки выше), когда рассматривали, как может концентрироваться энергия сигнала в различных участках вейвлет-спектра. Теперь только осталось просуммировать энергию отдельных коэффициентов. Поехали...
Нижеследующая функция coefplan берет выборку размером all сигнала sg, раскладывает в вейвлет-пакет по уровням (level - их количество) и записывает в массив e энергии вейвлет-коэффициентов.

Код

void coefplan(float *sg, double *e, int level, int all) {    complex<double> *vyx = new complex<double>[1025]; // ядро преобразования    complex<double> *sgw = new complex<double>[1024];    complex<double> *sgc = new complex<double>[1024];    int ss;    int num = 0;    double en, en1;    for(int i=0; i < all; i++)     {      sgc[i] = complex<double>(sg[i],0);     }    int flag = 0;    for(int l = 0;l < level; l++)     {      if(l != 0)       {        for(int sl = 0; sl < pow(2,l-1); sl++)        {         for(int i=0; i < all/pow(2,l-1); i++)          {           ss = all/(int)pow(2,l-1)*sl + i;           sgw[i] = sgc[ss];          }           if (all/(int)(pow(2,l-1))>1)           wavelet(sgw, all/(int)(pow(2,l-1)),flag); // разложение в вейвлет-пакет и возврат туда же           flag = 1;         en=0; en1=0;         for(int i=0;i<all/(int)(pow(2,l));i++)           {            en +=sqrt(norm(sgw[i]));            en1 +=sqrt(norm(sgw[i+all/(int)(pow(2,l))]));            sgc[i+sl*all/(int)pow(2,l-1)] = sgw[i];            sgc[i+all/(int)(pow(2,l))+sl*all/(int)pow(2,l-1)]                                               = sgw[i+all/(int)(pow(2,l))];           }         e[num] = en; // в e[] записывается энергия!         num = num + 1;         e[num] = en1; // в e[] записывается энергия!         num = num + 1;        }       }      else       {         en=0;         for(int i=0;i<all/(int)(pow(2,l));i++)          {           en += dabs(sg[i]);          }           e[num] = en; // в e[] записывается энергия!           num +=1;         //нормирование!!!         for(int i=0;i < all;i++)          {           sg[i] /=(e[0]/all);           sgc[i] /=(e[0]/all);          }          e[0] = 1;       }     }    for(int i=0;i < all;i++)    {       sg[i] = real(sgc[i]);    }    delete[] sgw; sgw=0;    delete[] sgc; sgc=0;    delete[] vyx; vyx=0; }

В принципе, ничего сложного: взяли ту самую картинку, что выше, и получили новую:

Для того, чтобы все имеющиеся эталоны были соизмеримы, надо их еще отнормировать к максимальному значению энергии. Тогда максимальный будет иметь значение, равное единице. В этом случае алгоритм не будет ошибаться на разнице эталонов и реальных сигналов в абсолютных значениях энергетики.

Цитата

//нормирование!!!         for(int i=0;i < wavesig_dim;i++) //1024          {           sperm[i] /=(en/wavesig_dim);// e[0];           sgc[i] /=(en/wavesig_dim);// e[0];          }       }         for(int i=0;i<wavesig_dim;i++)            serm[l*wavesig_dim+i] = real(sgc[i]);     }

Добрый день!
Можно покоментировать этот фрагмент? Откуда берётся массив serm и какова его длинна ? Где будут результаты преобразования? В serm или sperm ?