Локализация источника с использованием обобщенной перекрестной корреляции

Открыть Live Script

Этот пример показывает, как определить положение источника широкополосного сигнала с помощью обобщенной перекрестной корреляции (GCC) и триангуляции. Для простоты этот пример ограничен двумерным сценарием, состоящим из одного источника и двух приёмных массивов. Можно расширить этот подход более чем на два датчика или массивы и на три размерности.

Введение

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

Оцените направление источника из каждого сенсорного массива с помощью алгоритма оценки DOA. Для широкополосных сигналов многие хорошо известные алгоритмы оценки направления прибытия, такие как метод Capon или MUSIC, не могут быть применены, потому что они используют различие фаз между элементами, что делает их подходящими только для узкополосных сигналов. В широкополосном случае вместо информации о фазе можно использовать различие во времени прибытия сигнала среди элементов. Чтобы вычислить различия времени прибытия, этот пример использует обобщенную перекрестную корреляцию с фазой алгоритмом преобразования (GCC-PHAT). Из различий во времени прибытия можно вычислить DOA. (Для другого примера узкополосных алгоритмов оценки DOA, см. Направление высокого разрешения оценки прибытия).
Вычислим исходное положение по триангуляции. Сначала проведите прямые линии из массивов по направлениям прибытия. Затем вычислите пересечение этих двух линий. Это исходное расположение. Локализация источника требует знания положения и ориентации приемных датчиков или сенсорных массивов.

Формула триангуляции

Алгоритм триангуляции основан на простых тригонометрических формулах. Предположим, что массивы датчиков расположены в координатах 2-D (0,0) и (L,0), и неизвестное местоположение источника является (x, y). Из знаний о положении сенсорных решеток и двух направлениях прибытия в массивы , $θ_{1}$ и $θ_{2}$ , можно вычислить (x, y) координаты из

$L = y \tan θ_{1} + y \tan θ_{2}$

который вы можете решить для y

$y = L / (\tan θ_{1} + \tan θ_{2})$

а затем для x

$x = y \tan θ_{1}$

В оставшейся части этого примера показано, как можно использовать функции и системные объекты Phased Array System Toolbox™ для вычисления исходного положения.

Геометрия источника и датчика

Установите два приемных ULA с 4 элементами, выровненных по оси X глобальной системы координат и разнесенных на 50 метров. Центр фазы первого массива является (0,0,0). Центр фазы второго массива является (50,0,0). Источник расположен на (30 100) метрах. Как показано на рисунке, приемный массив усиливает точку в направлении + y. Источник передает в направлении -y.

Задайте базовую линию между массивами датчиков.

L = 50;

Создайте 4-элементный приемник ULA всенаправленных микрофонов. Можно использовать ту же phased.ULA Системные object™ для phased.WidebandCollector и phased.GCCEstimator Системные объекты для обоих массивов.

N = 4;
rxULA = phased.ULA('Element',phased.OmnidirectionalMicrophoneElement,...
    'NumElements',N);

Задайте положение и ориентацию первого массива датчиков. Когда вы создаете ULA, элементы массива автоматически разнесены вдоль оси Y. Необходимо повернуть локальные оси массива на 90 °, чтобы выровнять элементы вдоль оси X глобальной системы координат.

rxpos1 = [0;0;0];
rxvel1 = [0;0;0];
rxax1 = azelaxes(90,0);

Задайте положение и ориентацию второго массива датчиков. Выберите локальные оси второго массива, чтобы выровнять их с локальными осями первого массива.

rxpos2 = [L;0;0];
rxvel2 = [0;0;0];
rxax2 = rxax1;

Укажите источник сигнала как один всенаправленный преобразователь.

srcpos = [30;100;0];
srcvel = [0;0;0];
srcax = azelaxes(-90,0);
srcULA = phased.OmnidirectionalMicrophoneElement;

Задайте форму волны

Выберите исходный сигнал, который будет широкополосным LFM-сигналом. Предположим, что рабочая частота системы составляет 300 кГц, и установите пропускную способность сигнала равной 100 кГц. Примите максимальную рабочую область значений 150 м. Затем можно задать интервал повторения импульса (PRI) и частоту повторения импульса (PRF). Примите коэффициент заполнения 10% и установите ширину импульса. Наконец, используйте скорость звука в подводном канале 1500 м/с.

Установите параметры формы волны LFM и создайте phased.LinearFMWaveform Системные object™.

fc = 300e3;             % 300 kHz
c = 1500;               % 1500 m/s
dmax = 150;             % 150 m
pri = (2*dmax)/c;
prf = 1/pri;
bw = 100.0e3;           % 100 kHz
fs = 2*bw;
waveform = phased.LinearFMWaveform('SampleRate',fs,'SweepBandwidth',bw,...
    'PRF',prf,'PulseWidth',pri/10);

Сигнал передачи может затем быть сгенерирован как

signal = waveform();

Излучение, распространение и сбор сигналов

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

Установите количество поддиапазонов равным 128.

nfft = 128;

Задайте исходный излучателя и коллекторы сенсорной решётки.

radiator = phased.WidebandRadiator('Sensor',srcULA,...
    'PropagationSpeed',c,'SampleRate',fs,...
    'CarrierFrequency',fc,'NumSubbands',nfft);
collector1 = phased.WidebandCollector('Sensor',rxULA,...
    'PropagationSpeed',c,'SampleRate',fs,...
    'CarrierFrequency',fc,'NumSubbands',nfft);
collector2 = phased.WidebandCollector('Sensor',rxULA,...
    'PropagationSpeed',c,'SampleRate',fs,...
    'CarrierFrequency',fc,'NumSubbands',nfft);

Создайте широкополосные распространители сигналов для путей от источника до двух массивов датчиков.

channel1 = phased.WidebandFreeSpace('PropagationSpeed',c,...
    'SampleRate',fs,'OperatingFrequency',fc,'NumSubbands',nfft);
channel2 = phased.WidebandFreeSpace('PropagationSpeed',c,...
    'SampleRate',fs,'OperatingFrequency',fc,'NumSubbands',nfft);

Определите направления распространения от источника до сенсорных массивов. Направления распространения соответствуют локальной системе координат источника.

[~,ang1t] = rangeangle(rxpos1,srcpos,srcax);
[~,ang2t] = rangeangle(rxpos2,srcpos,srcax);

Излучайте сигнал от источника в направлениях массивов.

sigt = radiator(signal,[ang1t ang2t]);

Затем передайте сигнал в массивы датчиков.

sigp1 = channel1(sigt(:,1),srcpos,rxpos1,srcvel,rxvel1);
sigp2 = channel2(sigt(:,2),srcpos,rxpos2,srcvel,rxvel2);

Вычислите направления прихода распространенного сигнала в массивах датчиков. Поскольку реакция коллектора является функцией от направлений прихода в локальной системе координат сенсорного массива, передайте матрицы локальных осей координат в rangeangle функция.

[~,ang1r] = rangeangle(srcpos,rxpos1,rxax1);
[~,ang2r] = rangeangle(srcpos,rxpos2,rxax2);

Соберите сигнал в массивах датчиков приема.

sigr1 = collector1(sigp1,ang1r);
sigr2 = collector2(sigp2,ang2r);

Оценка и триангуляция GCC

Создайте оценки GCC-PHAT.

doa1 = phased.GCCEstimator('SensorArray',rxULA,'SampleRate',fs,...
    'PropagationSpeed',c);
doa2 = phased.GCCEstimator('SensorArray',rxULA,'SampleRate',fs,...
    'PropagationSpeed',c);

Оцените направления прибытия.

angest1 = doa1(sigr1);
angest2 = doa2(sigr2);

Триангулируйте исходное положение используйте формулы, установленные ранее. Поскольку сценарий ограничен плоскостью x-y, задайте нулевую координату z.

yest = L/(abs(tand(angest1)) + abs(tand(angest2)));
xest = yest*abs(tand(angest1));
zest = 0;
srcpos_est = [xest;yest;zest]

srcpos_est = 3×1

   29.9881
  100.5743
         0

Предполагаемое местоположение источника соответствует истинному местоположению в пределах 30 см.

Сводные данные

Этот пример показал, как выполнить локализацию источника с помощью триангуляции. В частности, пример показал, как моделировать, распространять и обрабатывать широкополосные сигналы. Алгоритм GCC-PHAT используется для оценки направления прихода широкополосного сигнала.

Документация