СПОСОБ ИЗМЕРЕНИЯ ЧАСТОТЫ СИНУСОИДАЛЬНЫХ СИГНАЛОВ И УСТРОЙСТВО ДЛЯ ЕГО РЕАЛИЗАЦИИ
METHOD OF MEASURING SINE-WAVE SIGNALS AND DEVICE TO THIS END

Изобретение относится к технике связи и может использоваться для измерения частоты синусоидальных сигналов в информационно-измерительных устройствах.
Задача измерения частоты синусоидальных сигналов постоянно привлекает внимание исследователей в связи с ее фундаментальной значимостью при оценке параметров периодических сигналов, смешанных с шумом, для извлечения информации о доплеровском сдвиге в задачах связи, навигации и радиолокации [1-3].
Известен ряд способов измерения частоты синусоидальных сигналов [4-7], однако они не являются статистически оптимальными с точки зрения измерения частоты синусоидальных сигналов на фоне шумов [8, 9].
Известен ряд способов измерения частоты синусоидальных сигналов [10, 11], основанных на использовании преобразования Фурье, при которых за оценку частоты синусоидального сигнала принимается аргумент максимума спектральной плотности мощности (СПМ) смеси сигнала с шумом. Указанные способы измерения частоты являются статистически оптимальными с точки зрения измерения частоты синусоидальных сигналов на фоне шумов [8, 9]. Однако данные способы требуют больших вычислительных затрат, связанных с осуществлением одного или нескольких преобразований Фурье, что не позволяет получать оценку частоты в масштабе времени, близком к реальному.
Известен способ измерения частоты синусоидальных сигналов [12-14], основанный на использовании разностно-фазовой статистики сигнала, который также является статистически оптимальным с точки зрения измерения частоты синусоидальных сигналов на фоне шумов [12, 13], достигающий при больших объемах выборки N нижней границы Крамера-Рао. Способ [14] является существенно более вычислительно эффективным, чем способы [8, 10, 11], так как не требует осуществления преобразований Фурье, наиболее близок к предлагаемому и поэтому принят за прототип.
Согласно этому способу:
1. Получают текущую фазу сигнала как аргумент комплексного числа, в качестве действительной части которого используют синфазные отсчеты, а в качестве мнимой - квадратурные отсчеты, преобразованных в цифровую форму, отфильтрованных составляющих квадратурного разложения сигнала.
k
2. Получают разности ψ смежных текущих фаз сигнала.
1
2
T
3. Формируют двухкомпонентный вектор состояния a(k)=[a(k),а(k)], компоненты которого задают рекуррентными выражениями
4. Определяют оценку частоты сигнала  в соответствии с формулой:
s
где f - частота выборки отсчетов сигнала.
Устройство-прототип [14] содержит подключенные к его входу, через две параллельные цепочки последовательно соединенных перемножителя, фильтра нижних частот (ФНЧ) и аналого-цифрового преобразователя (АЦП), соответствующие входы постоянного запоминающего устройства (ПЗУ), при этом вторые входы перемножителей связаны с общим источником синусоидального сигнала, у первого перемножителя - непосредственно, а у второго перемножителя - через фазовращатель, выход ПЗУ соединен с входом устройства вычитания, между вычитающим входом которого и выходом ПЗУ включен элемент памяти, а тактирующие входы обоих АЦП подключены к общему тактовому генератору (ТГ), к выходу устройства вычитания подключены первое и второе арифметические устройства (АУ), между выходами и вторыми входами которых включены соответственно второй и третий элементы памяти, третьи входы обоих АУ через счетчик соединены с тактовым генератором, между выходами первого и второго АУ включено второе устройство вычитания, выходной сигнал которого пропорционален измеряемой частоте сигнала, при этом первое АУ функционирует в соответствии с выражением
а второе АУ - в соответствии с другим выражением
k
1
1
2
2
k
где ψ - разность смежных текущих фаз сигнала с выхода ПЗУ, а(k) и a(k-1) - текущее и предшествующее значения соответственно первой компоненты вектора состояния с выхода первого АУ, a(k) и a(k-1) - текущее и предшествующее значения второй компоненты вектора состояния с выхода второго АУ, k - индексная переменная k=1, 2, …, формируемая счетчиком по тактам тактового генератора, означающая порядковый номер текущей разности фаз ψ.
Однако, так же как и аналоги [8, 10-13], статистически оптимальный, с точки зрения измерения частоты синусоидальных сигналов на фоне шумов [12, 13], способ-прототип [14] при измерении частоты синусоидальных сигналов от подвижных источников, например при обработке сигналов РЛС с синтезированной апертурой [15], или от объектов дальней космической связи [16] типа Mars Pathfinder на длительных интервалах наблюдения будет вносить погрешность измерения, обусловленную появлением квадратической составляющей фазы сигнала, которая вызвана динамикой объекта наблюдения.
Техническим результатом изобретения является повышение точности измерения частоты синусоидальных сигналов на фоне шумов за счет синтеза статистически оптимального, вычислительно эффективного оценивания в условиях динамического наблюдения с обеспечением возможности измерения скорости изменения частоты сигнала.
i
i
i
Технический результат достигнут тем, что в способе измерения частоты синусоидального сигнала, включающем получение текущей фазы сигнала как аргумента комплексного числа, в качестве действительной части которого используют синфазный отсчет, а в качестве мнимой - квадратурный отсчет, преобразованных в цифровую форму, отфильтрованных составляющих квадратурного разложения сигнала, получение разностей ψ смежных текущих фаз сигнала, где i=1, 2, … - индексная переменная, означающая порядковый номер текущей разности фаз ψ, согласно изобретению накапливают блок данных длиной К=N-1 из последовательных текущих разностей фазы сигнала ψ, i=1, 2, …, K по выборке из N пар синфазных и квадратурных отсчетов сигнала, формируют весовую функцию в соответствии со следующим выражением
формируют другую весовую функцию в соответствии с другим выражением
оценку частоты сигнала определяют в соответствии с формулой
s
где f - частота выборки отсчетов сигнала, а оценку скорости изменения частоты сигнала определяют в соответствии с формулой
Другим техническим результатом изобретения является повышение точности измерения приращения фазы синусоидальных сигналов при условиях динамического наблюдения.
Технический результат достигнут тем, что в способе измерения частоты синусоидальных сигналов, согласно изобретению, дополнительно определяют функцию приращения фазы сигнала в соответствии со следующим выражением
Другим техническим результатом изобретения является устройство для реализации способа измерения частоты синусоидальных сигналов, повышающее точность измерения частоты, за счет использования статистически оптимального, вычислительно эффективного оценивания в условиях динамического наблюдения с одновременным обеспечением возможности измерения скорости изменения частоты сигнала.
i
Технический результат достигнут тем, что в устройстве для реализации способа измерения частоты синусоидальных сигналов, к входу которого через две параллельные цепочки последовательно соединенных перемножителя, фильтра нижних частот (ФНЧ) и аналого-цифрового преобразователя (АЦП) подключены соответствующие входы постоянного запоминающего устройства (ПЗУ), при этом вторые входы перемножителей связаны с общим источником синусоидального сигнала, у первого перемножителя - непосредственно, а у второго перемножителя - через фазовращатель, выход ПЗУ, формирующего текущую фазу сигнала, соединен со входом устройства вычитания, между вычитающим входом которого и выходом ПЗУ включен элемент памяти, а тактирующие входы обоих АЦП подключены к общему тактовому генератору (ТГ), соединенному также со счетчиком, формирующим на выходе индексную переменную i=1, 2, …, K, к выходу устройства вычитания подключены первое и второе арифметические устройства (АУ), третьи входы обоих АУ через счетчик соединены с тактовым генератором, при этом на выходе устройства вычитания формируют последовательные текущие разности ψ, i=1, 2, …, K фазы сигнала по выборке из N пар синфазных и квадратурных отсчетов, согласно изобретению второй вход первого АУ подключен к выходу второго ПЗУ, четвертый вход первого АУ подключен к выходу третьего ПЗУ, а второй вход второго АУ соединен с выходом третьего ПЗУ, при этом первое АУ, выходной сигнал которого пропорционален измеряемой частоте сигнала, функционирует в соответствии с выражением
s
где f - частота ТГ, второе АУ, выходной сигнал которого пропорционален скорости изменения частоты сигнала, функционирует в соответствии с другим выражением
второе ПЗУ формирует весовую функцию из K коэффициентов согласно следующей формуле
а третье ПЗУ формирует другую весовую функцию из K коэффициентов, согласно другой формуле
при этом входы второго и третьего ПЗУ подключены к выходу счетчика, а между входом сброса счетчика и ТГ включено устройство управления
Другим техническим результатом изобретения является устройство для реализации способа, повышающее точность измерения приращения фазы синусоидальных сигналов при условиях динамического наблюдения.
Технический результат достигнут тем, что в устройстве для реализации способа измерения частоты синусоидальных сигналов, согласно изобретению, к выходам первого и второго АУ дополнительно подключено третье АУ, выходной сигнал которого пропорционален функции приращения фазы сигнала, при этом третье АУ функционирует в соответствии со следующим выражением
На фиг.1 представлена структурная схема устройства, в котором реализуется предложенный способ измерения частоты синусоидальных сигналов.
На фиг.2 представлена структурная схема другого варианта устройства, в котором реализуется предложенный способ измерения частоты синусоидальных сигналов.
Согласно предлагаемому способу измерения частоты синусоидальных сигналов:
1. Получают текущую фазу сигнала как аргумент комплексного числа, в качестве действительной части которого используют синфазный отсчет, а в качестве мнимой - квадратурный отсчет, преобразованных в цифровую форму, отфильтрованных составляющих квадратурного разложения сигнала.
k
2. Получают разности ψ смежных текущих фаз сигнала.
i
3. Накапливают блок данных длиной K=N-1 из последовательных текущих разностей фазы сигнала ψ, i=1, 2, …,K по выборке из N пар синфазных и квадратурных отсчетов сигнала.
4. Формируют весовую функцию в соответствии со следующим выражением
5. Формируют другую весовую функцию, в соответствии с другим выражением
6. Определяют оценку частоты сигнала в соответствии с формулой
s
где f - частота выборки отсчетов сигнала.
7. Определяют оценку скорости изменения частоты сигнала в соответствии с формулой
.
Полученная, таким образом, оценка частоты сигнала является наилучшей в среднеквадратическом смысле оценкой в динамической модели наблюдения сигнала.
Покажем это.
При наличии движения сигнала с ускорением его полная фаза в дискретном времени или в дискретной модели наблюдения определяется следующим образом:
0
0
I
I
I
I
mn
где φ - начальная фаза сигнала, ω - его циклическая частота,  - производная циклической частоты сигнала, обусловленная, например, ускорением, а для ошибок измерений, определенных белым шумовым процессом ξ(n), предполагается нулевое среднее и отсутствие корреляции между отсчетами E{ξ(n)}=0, E{ξ(n)ξ{m)}=qδ.
Пусть, так же как и в прототипе, получены разностно-фазовые измерения
0
Требуется получить оценки параметров частоты ω и ее производной . При этом синтезируемые оценки должны, с одной стороны, быть оптимальными (в уточняемом далее смысле), а с другой, должны быть вычислительно эффективными.
Введя в (2) обозначения
0
легко видеть, что искомые величины ω,  однозначно выражаются через введенные параметры
а модель разностно-фазовых измерений можно тождественно переписать как
I
I
где введена разностная случайная последовательность ζ(n)=ξ(n)-ξ(n-1).
Переходя к матричным обозначениям, систему (5) запишем в виде
ς
где С - ковариационная матрица случайной последовательности ζ(n) и введены векторы
а векторы ψ,ς имеют компоненты ψ(n), ζ(n) соответственно.
Как показано в [12-13], ковариационная матрица ошибок разностно-фазовых измерений может быть записана в виде
где элементы нормированной обратной ковариационной матрицы, используемой далее как весовая матрица, задаются соотношением
В условиях (6-8) оценка Гаусса-Маркова (ОГМ) [17] для вектора θ в матричном виде может быть записана как
где можно показать, что при K>2 обратная матрица в (9) существует.
1
2
Если в практической реализации требуется ограничиться элементарными арифметическими операциями, как это сделано в [12-13], то представляет интерес декомпозиция системы (9) на два независимых уравнения относительно θ и θ соответственно. Это может быть относительно легко достигнуто вследствие специального вида матрицы, подлежащей обращению в (9), а именно
где введены величины
,
значения которых вычисляются непосредственно по определению. Соответственно, оценка (9) принимает вид
тогда для отдельных компонент оцениваемых компонент соответствующие оценки могут быть записаны в виде скалярных произведений или взвешенных сумм исходных разностно-фазовых измерений
с весовыми функциями, которые задаются векторами
где
1
2
1
2
i
После тождественных преобразований векторы весовых функций w и w выражений (13)-(14) можно переписать как скалярные функции w(i), w(i) аргумента i, отражающего номер отсчета текущей разности фаз ψ, для блока данных разности фаз длиной K, i=1, 2, …, K, в следующем виде
0
0
s
s
С учетом того, что f=ω/2π и с учетом того, что расстояние во временной области между отсчетами сигнала равно 1/f, где f - частота выборки сигнала, а также с учетом соотношений (4) и (12) выражения для оценки частоты сигнала  и для оценки скорости изменения частоты сигнала  имеют следующий вид
Синтезированные, таким образом, оценки Гаусса-Маркова параметров сигнала являются наилучшими в среднеквадратическом смысле оценками в динамической модели наблюдения сигнала.
Оценки параметров (17)-(18), в отличие от оценок частоты сигнала [12-14], являются статистически эффективными в условиях наличия динамики в изменении частоты и не вносят динамических составляющих погрешности в процессе измерений. В то же время они являются столь же вычислительно эффективными, что и оценки [12-14], поскольку порядок требуемых арифметических операций совпадает.
Предложенный способ измерений позволяет также посредством прямого интегрирования в дискретном времени (суммирования) величин  и  получать оценку функции приращения фазы сигнала в соответствии со следующим выражением
что и используется во втором варианте предлагаемого способа.
Оценка функции приращения фазы сигнала нужна, в частности, в приложениях, связанных с обработкой речевых сигналов [18-19], когда необходимо восстановление однозначной функции фазы (в зарубежной литературе для обозначения такой операции над сигналом используется термин «unwrapping»), являющейся в общем случае разрывной 2π-периодической функцией.
Известен ряд способов для оценки функции приращения фазы сигнала [20-22], однако они не дают однозначного восстановления функции фазы в условиях дребезга скачка фазы на 2π при малых отношениях сигнал/шум и при наличии большой динамики в частоте сигнала [15, 18-22]. Поскольку предложенный способ изначально синтезирован для работы в динамической модели наблюдения, он может обеспечить восстановление однозначной функции фазы посредством прямого интегрирования в дискретном времени величин  и  с использованием выражения (19) с исключением влияния дребезга скачка фазы на 2π.
Согласно второму варианту предлагаемого способа измерения частоты синусоидального сигнала:
1. Получают текущую фазу сигнала как аргумент комплексного числа, в качестве действительной части которого используют синфазный отсчет, а в качестве мнимой - квадратурный отсчет, преобразованных в цифровую форму, отфильтрованных составляющих квадратурного разложения сигнала.
i
2. Получают разности ψ смежных текущих фаз сигнала.
i
3. Накапливают блок данных длиной K=N-1 из последовательных текущих разностей фазы сигнала ψ, i=1, 2, …, K по выборке из N пар синфазных и квадратурных отсчетов сигнала.
4. Формируют весовую функцию в соответствии со следующим выражением
5. Формируют другую весовую функцию, в соответствии с другим выражением
6. Определяют оценку частоты сигнала в соответствии с формулой
s
где f - частота выборки отсчетов сигнала.
7. Определяют оценку скорости изменения частоты сигнала в соответствии с формулой
8. Определяют функцию приращения фазы сигнала в соответствии со следующим выражением
Устройство, реализующее предложенный способ измерения частоты синусоидальных сигналов (смотри фиг.1) содержит подключенные к его входу две параллельные цепочки из последовательно соединенных перемножителя 1, ФНЧ 2, АЦП 3 и, соответственно, из перемножителя 4, ФНЧ 5, АЦП 6. К выходам АЦП 3, 6 подключены соответствующие входы ПЗУ 7, при этом вторые входы перемножителей 1, 4 связаны с общим источником синусоидального сигнала 8, у перемножителя 1 - непосредственно, а у перемножителя 4 - через фазовращатель 9. Выход ПЗУ 7 соединен с входом устройства вычитания 10, между вычитающим входом которого и выходом ПЗУ 7 включен элемент памяти 11. Тактирующие входы обоих АЦП 3, 6 подключены к общему ТГ 12. ТГ 12 соединен также со счетчиком 13, формирующим на выходе индексную переменную i=1, 2, …, К. К выходу устройства вычитания 10 подключены первое арифметическое устройство (АУ) 14 и второе АУ 15. Третьи входы обоих АУ 14 и 15, через счетчик 13, соединены с ТГ 12. Согласно изобретению второй вход первого АУ 14 подключен к выходу второго ПЗУ 16, четвертый вход первого АУ подключен к выходу третьего ПЗУ 17. Входы второго ПЗУ 16 и третьего ПЗУ 17 подключены к выходу счетчика 13, а между входом сброса счетчика 13 и ТГ 12 включено устройство управления 18.
В другом варианте описанного выше устройства (смотри фиг.2) к выходам АУ 14 и 15 подключено третье АУ 19.
Первый вариант предлагаемого устройства (смотри фиг.1), реализующего способ измерения частоты синусоидального сигнала, работает следующим образом.
i
Сигнал в смеси с шумом с входа устройства подается на первые входы перемножителей 1, 4. В продуктах перемножения входного сигнала с сигналом генератора 8, частота которого равна центральной частоте диапазона входных частот сигнала, на выходе перемножителя 1 содержится синфазная составляющая входного сигнала. В продуктах перемножения входного сигнала с сигналом генератора 8, развернутым по фазе на 90° фазовращателем 9, на выходе перемножителя 4 содержится квадратурная составляющая входного сигнала. Каждая из указанных составляющих в полосе низких частот, равной полуширине диапазона частот входных сигналов, фильтруется в ФНЧ 2 и 5 соответственно. Далее АЦП 3 преобразует синфазную составляющую входного сигнала в цифровую форму (в синфазные отсчеты), а АЦП 6 преобразует квадратурную составляющую входного сигнала в цифровую форму (в квадратурные отсчеты). В ПЗУ 7 получают текущую фазу сигнала как аргумент комплексного числа, в качестве действительной части которого используют синфазный отсчет, а в качестве мнимой - квадратурный отсчет. Получают разности смежных текущих фаз сигнала ψ, как результат вычитания текущей фазы с выхода ПЗУ 7 в устройстве вычитания 10 с предыдущим значением фазы, задержанным в элементе памяти 11.
i
1
2
В первом АУ 14 на основе разности ψ смежных текущих фаз сигнала, полученной с выхода устройства вычитания 10, на основе индексной переменной i, поступающей на третий вход АУ 13 из счетчика 13, запускаемого генератором ТГ 12, а также на основе коэффициентов весовых функций w(i) и w(i), получаемых соответственно из второго ПЗУ 16 и сформированных согласно выражению
и из третьего ПЗУ 17 и сформированных согласно другому выражению
определяют оценку частоты сигнала в соответствии с формулой
Накопление суммы в АУ 14, согласно последней формулы, осуществляется по мере роста индексной переменной i и по достижении индексной переменной величины К оценка частоты сигнала  из АУ 14 подается на выход предлагаемого устройства.
i
2
Во втором АУ 15 на основе разности ψ смежных текущих фаз сигнала, полученной с выхода устройства вычитания 10, на основе индексной переменной i, поступающей на третий вход АУ 13 из счетчика 13, запускаемого генератором ТГ 12, а также на основе коэффициентов весовой функции w(i) 5 получаемых из третьего ПЗУ 17 и сформированных согласно выражению
определяют оценку скорости изменения частоты сигнала в соответствии с формулой
Накопление суммы в АУ 15, согласно последней формуле, осуществляется по мере роста индексной переменной i и по достижении индексной переменной величины К оценка скорости изменения частоты сигнала  из АУ 15 также подается на выход предлагаемого устройства.
При этом одновременно с выдачей оценок параметров сигнала на выход устройства осуществляется сброс содержимого промежуточных регистров обоих АУ 14 и 15, а также сброс счетчика 13 устройством управления 18, которое также накапливает и контролирует состояние индексной переменной i по тактам ТГ 12. Сброс счетчика 13 осуществляют через вход сброса счетчика 13, связанного с выходом устройства управления 18, по достижении индексной переменной величины К.
Работа второго варианта предлагаемого устройства (фиг.2), реализующего способ определения координат радиопередатчиков с использованием ГНС, отличается от работы первого тем, что по результатам полученных оценок в АУ 14 (оценка частоты сигнала ) и 15 (оценка скорости изменения частоты сигнала ) на выходе третьего АУ 19 определяют функцию приращения фазы сигнала в соответствии со следующим выражением
Полученные результаты могут быть использованы, в частности, для повышения точности измерения частоты и других параметров сигнала в условиях нелинейной динамики движения источника сигнала.
ИСТОЧНИКИ ИНФОРМАЦИИ
1. Финкельштейн М.И. Основы радиолокации. М.: Радио и связь, 1983, - 454 с.
2. Бортовые устройства спутниковой радионавигации. - Под ред. B.C.Шебшаевича, М.: Транспорт, 1988, - 201 с.
3. Спилкер Дж. Цифровая спутниковая связь. М.: Связь, 1979, - 592 с.
4. Ермолов Р.С. Цифровые частотомеры. Л.: Энергия, 1973, с.13-17.
5. А.с. №541123. Толстых B.C. Способ измерения частоты. М., ЦНИИПИ, 1977, - 2 с.
6. А.с. №1659893. Золотков К.Л. и др. Способ измерения частоты синусоидального напряжения. М., ВНИИПИ, 1986, - 3 с.
7. А.с. №1798717. Минц М.Я. и др. Способ измерения частоты синусоидальных сигналов. М., ВНИИПИ, 1987, - 9 с.
8. Rife D. С., Boorstyn R.R., Single-Tone Parameter Estimation from Discrete Time Observations. // IEEE Trans. on Information Theory, Vol.IT-20, No.5, Sept. 1974.
9. Kay S.М., Fundamentals of Statistical Signal Processing: Estimation Theory, NJ, Prentice-Hall, 1993.
10. А.с. №541123. Гольдштейн И.Н. Цифровой способ измерения частоты и фазы гармонического сигнала. М., ЦНИИПИ, 1974, - 3 с.
11. А.с. №1352390. Кузьменков В.Ю. Способ определения частоты М., ВНИИПИ, 1987, - 4 с.
12. Кау S.М., Statistically/Computationally Efficient Frequency: Estimation. //Proceedings of the ICASSP'88, 1988, p.2292-2295 - аналог.
13. Kay S.М., A Fast and Accurate Single Frequency: Estimator. // IEEE Trans. on ASSP, Vol.37, No.12, December, 1989.
14. Подчиненко Н.Е., Скрипкин А.А., Щербачев В.А., Патент на изобретение РФ №2183839, Способ измерения частоты синусоидальных сигналов и устройство для его реализации, М., ФИПС, 2002 г. - прототип.
15. Cumming I.G., Wong F.H., Digital Processing of Synthetic Aperture Radar Data, Boston-London, Artech House.
16. Harcke, L., and G. Wood, Laboratory and Flight Performance of the Mars Pathfinder (15,16) Convolutionally Encoded Telemetry Link. TDA PR 42-129, January-March 1997, pp.1-11, May 15, 1997.
17. Альберт А. Регрессия, псевдоинверсия и рекуррентное оценивание. М.: Наука, 1977, - 224 с.
18. Применение цифровой обработки сигналов. Под ред. Э. Оппенгейма, М.: Мир, 1980, 552 с.
19. L.R.Rabiner and R.U.Schafer, Digital Processing of Speech Signals, Prentice-Hall.
20. Tribolet I. М., A new phase unwrapping algorithm. // IEEE Trans. Acoust., Speech, Signal Processing, vol. ASSP-26, pp.170-177, 1977.
21. Steiglitz К., Dickinson В., Phase Unwrapping by Factorization, IEEE Trans. Acoust., Speech, Signal Processing, vol. ASSP-30, no.6, pp.984-991, 1982.
22. Loeffler C.М., Leonard R.E., Phase Unwrapping мiа Median Filtering, Proceedings of the ICASSP'84, pp.46.8.1-46.8,3, 1984.