Основная задача.
Имею два АМ модулированных сигнала: 38400Гц(модулирован 160Гц) и 41000Гц (модулирован 166Гц).
На приеме, необходимо вычислять отношение фаз модулирующих колебаний с точностью 0.5-1%.
Для простоты, ЦАП на передатчике и АЦП на приемнике тактируются одной частотой (100кГц), в добавок, начало выборки синхронизировано с начальной фазой модулирующего колебания (если рассматривать только один АМ сигнал).
Для удобства поиска/отладки алгортма я подружил свой макет с MatLab-ом, всю выборку (размер которой могу управлять) вываливаю в PC, далее кручу-верчу ... (только обмануть пока не могу )

Итак, для вычисления фазы огибающей мне необходимо знать фазу несущей, хотя ЦАП передатчика синхронен с АЦП приемника, имеется среда передачи и аналоговый тракт, который скорее всего, с течением времени/температуры будет двигать фазу несущей (еще, в идеале, АЦП и ЦАП будут работать асинхронно). Т.е. мне как минимум надо вычислить фазу несущей, фазу одной из огибающих и соответственно вычислить искомую фазу.
Беда в том, что в своих экспериментах я никак не могу получить нужную точность вычисления фазы несущей, у меня меньше 6% не получается . Значение фазы беру из Фурье-преобразования входной выборки с наложеным окном Ханнинга, АЦП 12 разрядов, ЦАП 4 разряда(на прецизионной 0,1% R2R резистивной матрице). Смущает еще - с ростом количества отчетов в выборке, точность вычисления фазы не меняется (строил гистограммы на разные длины выборки, процесс похож на нормальный но не "обужается" с ростом количества отсчетов).

А какие еще возможны способы вычисления фазы, имеет смысл в моей задаче использовать вейвлет?

Уж не хотите ли Вы сказать, что пытаетесь найти фазу сигналов 38.4кГц и 41.0kГц имея выборку с частотой дискретизации 100kГц? Если да, тогда полученная Вами точность 6% - впечатляет. Если же речь идет о модулируемом сигнале 160 и 166 Гц соответственно, тогда, не понятно почему получается такая низкая точность. Для повышения точности попробуйте увеличить длину выборки, а также постарайтесь чтобы в выборку входило строго целое число периодов гармоники, фазу которой Вы вычисляете.

Сообщение отредактировал defunct - Jan 30 2006, 22:03

6% я получаю из фурье-преобразования что на несущие, что на боковые. Математически, Фн=Ф-Ф0 Фв=Ф+Ф0, соответсвенно меньше 6% я не получу. Длину выборки "бесконечно" увеличивать я не могу, максимум 1000 отсчетов (или 10мСек), причем как уже сказал ранее точность не меняется что с 1000 что с 10000 отсчетами. Периоды у меня "целые", я для этого специальную прогу в матлабе накатал, да и просто в качестве эксперимента, перебрал длительность по одному отсчету на период огибающей.

Хм-м, хотя последовательности в самом генераторе у меня хранятся безразрывные, про кратность количества точек я забыл , завтра проверю.

Сообщение отредактировал 3.14 - Jan 30 2006, 22:24

1. Для оценки фазы модулирующего колебания совсем не обязательно знать несущую с точностью до фазы. В данном случае достаточно реализовать некогерентный АМ-приемник, т.е. некоторую частотно-избирательную систему (например фильтр), выделяющую сигнал на данной частоте, после которого идет амплитудный детектор, а после него задача сводится к оценке фазы низкочастотного гармонического колебания. Эту задачу обычно решают, используя квадратурный корреляционный детектор. При этом интервал интегрирования корреляторов нужно брать максимально возможным, привязывать его к частоте сигнала не обязательно. Если надо, и подвернется свободная минутка, могу модель в Simulink набрасать.

2. Необходимо вычислять именно отношение фаз??? не разность случайно? если разность - то в итоге получится гармонический сигнал с разностной частотой.

3. 4-разрядный ЦАП обеспечивает шумы квантования ~-24дБ, то есть отношение сигнал/шум на входе АЦП порядка 20дБ, причем шумы имеют ярко выраженный дискретный характер - будут сильно влиять на точность.

4. Имея сигналы 38.4кГц и 41.0kГц с частотой дискретизации 100kГц теоретически возможно получить бесконечную точность оценки фазы при отсутствии шумов и квантования, ибо частота дискретизации не кратна частотам сигналов.

2 Pathfinder
1) Оба АМ сигнала проходят через УЗ передатчик и приемник АЧХ которых носит резонансный характер, в результате сдвигаются как фазы боковых состовляющих так и амплитуда и в результате детектирования получу биение и как выковыривать от туда фазу не особо представляю, к тому же узкий фильтр затянет/оттянет время измерения, а оно у м еня должно быть в идеале 5мС, но можно и 10-20мС. Еще в этом способе я нахватаю помех от соседнего АМ. А насчет модельки Simulink, если сделаете буду очень признателен.

2) Пока не особо важно, если нет необходимой точности для последующих манипуляций.

3) От части он режется УЗ датчиками, я пробовал считать фурье и из 100000 отсчетов, влияние шума квантования ЦАП должно снизиться, но это не привело к ощутимым улучшениям.

Боюсь маловато такого времени измерения - 5мс, хотя бы один период надо целиком. Вот моделька - с ходу получается погрешность ~0.3 градуса при переходном процессе 30мс. Наверняка можно выжать и лучше если пооптимизировать параметры. Этот метод еще вроде называют методом Ньютона-Рафсона.

Спасибо за модедьку.

Можно попробывать так. Если я правильно понял вашу задачу.

Разрабатываете два бандпас ФИР фильтра на вашу частоту. Один делает косинус преобразование, другой синус преобразовние. На выходе у вас будет сигнал I and Q. Вот вам и угол. Задержа фильтра известна, ее можно учесть. Угол между соседними векторами и есть ваша частота. Вом она тоже известна. Так что вычислисть начальную фазу не составит большого труда. Можно потом еще среднее считать для увеличения точности.

Тоже самое но другим путем. Ставится один фильтр который выделит вашу частоту уберет шум из сигнала (шума должно быть как-можно меньше ). Зная частоту сигнала и частоту сэмплирования подготавливаются две таблицы косинуса и синуса на периуд. Сигнал после фильтра умножается на две эти таблицы и суммируется. Получаем два числа -> вектор. Угол начальная фаза.

Спасибо, интересные предложения, надо осмыслить

to etep
Уточните, пожалуйста, что вы имеете в виду под

Это просто классическое квадратурное детектирование путем домножения на sin(wt) и cos(wt) с последующим low pass FIR filter или что то другое?

Извените что долго не отвечал, был очень занят.
Делается это так.

Разрабатывается lowpass FIR с частатой среза Fc.
Имеем его коэффициенты K(n). Далее из этого фильтра получаем пару Hilbert FIR фильтров. Для этого применяем quadrature bandpass transform.

Ksin(n) = 2*K(n)*sin(2*pi*n*F/Fs)
Kcos(n) = 2*K(n)*cos(2*pi*n*F/Fs)

Где, F - это центральная частота и от нее +-Fc.
Fs - это частота семплирования.
Т.е. получается два Bandpass FIR фильтра с одинаковыми (почти) характеристиками.
Между двумя выходами имеем 90 градусов. Ввиду того что не всегда можно получить одинаковые характеристики, то можно к примеру сделать усреднение углов и уже среднее брать как более точную величину. Ну а дальше вы можете делать все что вам захочится. Задержка фильтра известа, фаза ленейна от частоты. И узнать в какой фазе пришел сигнал не составит большого труда.
Минусы в этой системе это два фильтра придется реализовать и тайм-летенси может быть большим.


Другой вариант на IIR фильтрах.
Разрабатываем lowpass IIR с частатой среза Fc.
Переходим в z-плоскость. На ней имеем комплексносопряженные нули и полюса.
Зная что в этой плоскости:
угол 0 соотвествует частоте 0
90 градусов (против часовой стрелки) Fs/4
180 - Fs/2
270 - 3*Fs/4
Где Fs - частота семплирования.

Делаем простую вещь. Поварачиваем каждый ноль и полюсь на угол = 360 * F/Fs
где F - центральная частота нашего будущего bandpass фильтра.
После поворота получили новые коэффициенты фильтра которые теперь комплексные. То есть у нас комплексный фильтр. Который имеет один вход и два выхода. В реализации это будет выглядить как два bandpass IIR фильтра c полосой 2Fc и между выходами ровно 90 градусов (комплексное число на выходе)

Проблемы: все которые имеет IIR фильтр
Если нужно что бы фаза была линейна от частоты, то дальше я вижу пока два варианта решения.
1. Приминить allpass фильтр и уже им скорректированть груповую задержку.
2. Это просто таблица. Где к примеру угол это адрес памяти по которому лежит число на какой угол нужно сделать поворот. Ну а если частоты или частота фиксирована, то это займет не много памяти.

to etep

Огромное спасибо за обстоятельный ответ.
Понимая, что вы очень заняты, не ожидаю быстрого ответа, но по возможности было бы интересно узнать вашу точку зрения.

Вы великолепно владеете теорией создания и преобразования фильтров.
С математической точки зрения, описанный вами алгоритм на discrete time domain, и моя догадка, сформулированная на continuous time domain, дают абсолютно идентичный результат. В вашем варианте фильтр сразу bandpass filter, в моем вначале frequency shift followed by lowpass filter.

Насколько я знаю большинство hardware Hilbert transformer, используемых например в Agilent приборах, построены по описанной мной схеме, с IF=4 * Fs. В этом случае frequency shift реализуется просто домножением (сменой знака) для cos path (1 0 -1 0 ...) sin path (0 1 0 -1 ...), далее идет обычно CIC filter или сразу low pass filter.

В DDC GC4016 Texas Instruments соотношение между IF и Fs не зажато, но алгоритм тот же, вначале frequency shift потом фильтрация.

Мой вопрос следующий, почему в приборах используется большей частью описанный мною вариант, а не ваш, который дает абсолютно тот же результат? У меня предположение, что так проще реализовать в железе, а вот если делать пост обработку на компьютере, то тогда ваш вариант возможно и предпочтительней, если не требуется zero phase feature.

В моих измерениях zero phase feature очень существенный момент, поэтому мой алгоритм вначале frequency shift потом двух проходная фильтрация (функция filtfilt MAtlab).

Прямо сейчас возник вопрос, а что будет, если использовать в filtfilt предложенный вами bandpass filter? Получится ли правильный результат с zero phase shift? Вечером постараюсь с этим поразбираться.

У меня есть еще один академический вопрос, на который мне не удалось найти ответа в книгах по DSP.
Суть его в следующем. Имеем датчик магнитного поля, предствляющий из себя резонансный контур, согласованный с импедансом приемника. С точки зрения DSP такой датчик есть ни что иное как bandpass filter с центральной частотой равной частоте резонанса. Требуется создать компенсирующий фильтер, который бы делал АЧХ датчика в заданном диапазоне близкой к линейной, при линейной фазе. Видимое мне лобовое решение проблемы deconvolution не совсем хорошо, так как сигнал относительно узкополосный и deconvolution значительно повышает шумы. Известно ли вам более правильное решение такой проблемы с использованием filter, полученного на основе знаний о резонансном контуре?

Сообщение отредактировал nadie - Feb 22 2006, 09:48

eteP, nadie, может будем к терминологии с уважением относиться, все-таки в России живем.

to Pathfinder
Объясните, пожалуйста, в чем ваша претензия в части терминологии. Мне существенно проще использовать английскую терминологию, так как все фундаментальные книги по DSP у меня на английском (возможно и есть рускоязычные издания, но они мне не доступны). IMHO, все специалисты работающие в области RF, DSP английский знают, так как и железо и софт используемый для разработок только англоязычный.

to eteP
Более подробно свои соображения по использованию вашего алгоритма изложу немного позднее, сейчас надо с рабочими делами разобраться.

Не знаю, слышит ли меня товарищ из Ижевска, но дело обстоит так:
метод у Вас правильный, у Вас окно не правильное. Просто посмотрите спектр окна Хэнна (или, что одно и тоже Хэннинга). Все вопросы у Вас исчезнут. Используйте окно Кайзера и в пол-процента Вы уложитесь.