Акустический журнал, 2021, T. 67, № 2, стр. 174-184

ВВЕДЕНИЕ

Оценка углового положения источника в мелководном канале с помощью горизонтальной линейной антенны является одной из важных прикладных задач гидроакустики. Как известно (см., например, [1–3]), ее решение, основанное на методе сканирования, приводит к заметным ошибкам в определении направления на источник, растущим с увеличением угла прихода (отсчитываемого от нормали к апертуре антенны) и числа распространяющихся нормальных волн. Традиционные алгоритмы (типа метода Кейпона, MUSIC, максимума правдоподобия) позволяют получить несмещенные оценки угловых координат (при условии точно известных параметров волновода), однако в процессе локализации они используют трудоемкую процедуру одновременного поиска по глубине, дальности и азимутальному углу, что требует больших вычислительных затрат. В этой связи возникает необходимость разработки устойчивых методов оценивания, позволяющих определить направление на источник без знания глубины его погружения и расстояния до приемной антенной решетки (АР).

Один из способов решения такого рода обратной задачи был предложен в работе [4] (см. также [5], где приведена альтернативная формулировка алгоритма). Соответствующая процедура оценивания, получившая название метода подпространственного пересечения или метода SI (“subspace intersection”), опирается на априорное знание волновых чисел распространяющихся нормальных волн. Однако в изменчивых и всегда не полностью известных условиях морской среды несоответствие (рассогласование) между расчетной моделью канала и реальным акустическим волноводом может приводить к значительному ухудшению работоспособности данного способа локализации.

В последнее десятилетие в теории обработки сигналов антенными решетками интенсивно развивается направление, связанное с определением минимального углового расстояния между источниками, при котором они могут быть корректно локализованы с использованием алгоритмов сверхразрешения. На сегодняшний день существуют три общепринятых критерия, позволяющих определить указанное расстояние. Первый, предложенный Коксом [6], представляет собой необходимое условие одновременного существования двух близко расположенных минимумов целевой функции. Однако такое рассмотрение не является универсальным, поскольку зависит от конкретного вида используемого метода оценивания. В этой связи наибольшее распространение получил критерий, основанный на статистической проверке гипотез о наличии одного или двух сигналов в принятой смеси [7–9], а также критерий Смита [10], в соответствии с которым источники считаются разрешенными, если угловое разнесение между ними превосходит среднеквадратичную ошибку оценивания этого разнесения.

В настоящей работе построен адаптивный алгоритм пониженного ранга NM-RARE (“normal mode based rank reduction”), предназначенный для определения угловых положений источников без знания их пространственных координат в волноводе с неточно известными параметрами. Проводимое рассмотрение основано на наихудшем сценарии приема, учитывающeм отличие ожидаемой реплики от истинной, и позволяющeм минимизировать эффекты рассогласования различной природы [11, 12]. Представлены результаты сравнительного анализа эффективности данного способа оценивания с методом сканирования и методом подпространственного пересечения. На основе критерия Смита найден теоретический предел углового разрешения в зависимости от отношения сигнал/шум и объема входной выборки.

1. ПОСТАНОВКА ЗАДАЧИ. ИСХОДНЫЕ СООТНОШЕНИЯ

Рассмотрим акустический волновод, в котором звуковое поле создается $J$ источниками, излучающими детерминированные узкополосные сигналы ${{s}_{j}}(t)$ $(j = 1, \ldots ,J)$ с одинаковой несущей частотой. Прием осуществляется линейной горизонтальной АР, состоящей из $N$ элементов, расположенных на горизонте ${{z}_{a}}.$ Положение $j$-го источника определяется глубиной его погружения ${{z}_{j}},$ расстоянием ${{r}_{j}}$ до приемной АР и азимутальным углом ${{\phi }_{j}},$ отсчитываемым от нормали к апертуре АР. Геометрия задачи показана на рис. 1. (Начало координат по дальности выбрано в месте установки первого элемента АР.)

В узкополосном приближении поле на входе АР характеризуется N-мерным вектором наблюдения ${{{\mathbf{x}}}_{l}}$:

Здесь $l$ – номер выборочного отсчета, ${{{\mathbf{\theta }}}_{j}} = {{({{r}_{j}},{{z}_{j}})}^{T}}$ (верхний индекс $T$ означает операцию транспонирования), ${\mathbf{g}}({{\phi }_{j}},{{{\mathbf{\theta }}}_{j}})$ = = $[G(\left. {{{{\mathbf{r}}}_{1}},{{z}_{a}}} \right|{{{\mathbf{r}}}_{j}},{{z}_{j}}), \ldots ,G{{(\left. {{{{\mathbf{r}}}_{N}},{{z}_{a}}} \right|{{{\mathbf{r}}}_{j}},{{z}_{j}}]}^{T}}$ – вектор отклика АР при приеме сигнала от$j$-го источника, $\{ G(\left. {{{{\mathbf{r}}}_{n}},{{z}_{a}}} \right|{{{\mathbf{r}}}_{j}},{{z}_{j}})\} _{{n = 1}}^{N}$ – функции Грина среды распространения, ${{{\mathbf{n}}}_{l}}$ – вектор аддитивного белого шума, а $L$ – объем входной выборки. Задача состоит в построении адаптивного алгоритма обработки, позволяющего по принятой выборке $\{ {{{\mathbf{x}}}_{l}}\} _{{l = 1}}^{L}$ оценить угловые положения источников без знания их пространственных координат, и определении наименьшего углового разнесения, при котором источники могут быть корректно разрешены.

При дальнейшем анализе будем считать, что ${{{\mathbf{n}}}_{l}}$ является случайным гауссовым вектором с нулевым средним значением и характеризуется ковариационной матрицей $\left\langle {{{{\mathbf{n}}}_{l}}{\mathbf{n}}_{l}^{ + }} \right\rangle = \sigma _{n}^{2}{\kern 1pt} {{{\mathbf{I}}}_{N}},$ где $\sigma _{n}^{2}\,$ – неизвестный уровень шума, ${{{\mathbf{I}}}_{N}}$ – единичная матрица размерности $N \times N,$ а ${{( \cdot )}^{ + }}$ и $\left\langle \cdot \right\rangle $ означают операции эрмитового сопряжения и статистического усреднения, соответственно.

В рамках волнового подхода функция Грина $G(\left. {{{{\mathbf{r}}}_{n}},{{z}_{a}}} \right|{{{\mathbf{r}}}_{j}},{{z}_{j}})$ может быть представлена в виде суперпозиции конечного числа $M$ распространяющихся нормальных мод. При расположении источников в дальней зоне антенны (когда $\left| {{{{\mathbf{r}}}_{n}} - {{{\mathbf{r}}}_{j}}} \right| \approx {{r}_{j}} - d(n - 1)\sin ({{\phi }_{j}})$, где $d$ – межэлементное расстояние) эта функция записывается следующим образом

Здесь ${{\varphi }_{m}}({{z}_{a}})$ и ${{\varphi }_{m}}({{z}_{j}})$ – собственные функции $m$‑ой моды на глубине расположения приемной АР и $j$-го источника излучения, соответственно, а ${{k}_{m}}$ – горизонтальное волновое число.

С использованием модового описания для вектора отклика АР имеем:

где ${\mathbf{b}}({\mathbf{\theta }})$ – вектор размерности $M \times 1$, компонентами которого являются амплитуды мод, определяемые формулой (2), а ${\mathbf{U}}(\phi )$ – матрица размерности $N \times M$ вида

При $M < N$ ранг матрицы ${\mathbf{U}}(\varphi )$ равен $M$ (при условии, что векторы $\{ {{{\mathbf{u}}}_{m}}(\phi )\} _{{m = 1}}^{M}$ линейно независимы).

С учетом (3) исходный вектор наблюдения (1) может быть переписан следующим образом

Здесь ${\mathbf{\varphi }} = {{({{\phi }_{1}}, \ldots ,{{\phi }_{J}})}^{T}}$ – искомый вектор направлений размерности $J{\kern 1pt} \times 1,$ ${\mathbf{\theta }} = {{({\mathbf{\theta }}_{1}^{T}, \ldots {\mathbf{\theta }}_{J}^{T})}^{T}}$ – вектор размерности $2J{\kern 1pt} \times 1,$ определяющий пространственные положения источников, ${\mathbf{G}}({\mathbf{\varphi }},{\mathbf{\theta }})$ = $ = [{\mathbf{U}}({{\phi }_{1}}){\mathbf{b}}({{{\mathbf{\theta }}}_{1}}), \ldots ,{\mathbf{U}}({{\phi }_{J}}){\mathbf{b}}({{{\mathbf{\theta }}}_{J}})]$ – передаточная матрица канала размерности $N \times J,$ а ${{{\mathbf{s}}}_{l}} = {{[{{s}_{1}}(l), \ldots ,{{s}_{J}}(l)]}^{T}}$ – вектор комплексных огибающих излученных сигналов размерности $J \times 1.$

Одним из наиболее распространенных способов локализации является метод MUSIC [13]. Эта процедура основана на использовании информации, содержащейся в системе собственных векторов $\{ {\mathbf{\hat {\psi }}}\} _{{i = 1}}^{N}$ выборочной матрицы ${{{\mathbf{\hat {\Gamma }}}}_{{\mathbf{x}}}} = ({1 \mathord{\left/ {\vphantom {1 L}} \right. \kern-0em} L})\sum\nolimits_{l = 1}^L {{\mathbf{x}}({{t}_{l}}){{{\mathbf{x}}}^{ + }}({{t}_{l}})} :$

где ${{\hat {\lambda }}_{1}}, \ldots ,{{\hat {\lambda }}_{N}}$ – положительные собственные числа, пронумерованные в порядке убывания, т. е. ${{\hat {\lambda }}_{1}} \geqslant {{\hat {\lambda }}_{2}} \geqslant \cdots \geqslant {{\hat {\lambda }}_{N}}$. Первые $J$ старших собственных векторов формируют сигнальное подпространство ${{{\mathbf{\hat {\Psi }}}}_{{\mathbf{s}}}} = [{{{\mathbf{\hat {\psi }}}}_{1}} \cdots {{{\mathbf{\hat {\psi }}}}_{J}}]$, а $N - J$ оставшихся векторов – шумовое ${{{\mathbf{\hat {\Psi }}}}_{{\mathbf{n}}}} = [{{{\mathbf{\hat {\psi }}}}_{{J + 1}}}, \ldots ,{{{\mathbf{\hat {\psi }}}}_{N}}].$

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

где ${{{\mathbf{\hat {\Pi }}}}_{{\mathbf{n}}}} = {{{\mathbf{\hat {\Psi }}}}_{{\mathbf{n}}}}{\mathbf{\hat {\Psi }}}_{{\mathbf{n}}}^{ + }$ – проекционная матрица на шумовое подпространство. Процедура оценивания сводится к поиску $J$ минимумов целевой функции ${{{\mathbf{g}}}^{ + }}(\phi ,{\mathbf{\theta }}){{{\mathbf{\hat {\Pi }}}}_{{\mathbf{n}}}}{\mathbf{g}}(\phi ,{\mathbf{\theta }})$ в трехмерной области параметров. Последнее требует больших вычислительных затрат.

Ниже нас будет интересовать оценка угловых координат источников в волноводе. Тогда, рассматривая модовый вектор ${\mathbf{b}}({\mathbf{\theta }})$ (зависящий от неинформационного параметра $\,{\mathbf{\theta }}$) как неизвестный и принимая во внимание представление (3), переформулируем критерий MUSIC (4) следующим образом

где $C(\phi ) = {{U}^{ + }}(\phi ){{{\mathbf{\hat {\Pi }}}}_{{\mathbf{n}}}}U(\phi )$ – матрица размерности $M \times M,$ ранг которой равен $rank\{ {\mathbf{C}}\} = $ $ = rank\{ {\mathbf{C}}\} = \min {\kern 1pt} {\kern 1pt} (N - J,M),$ и при совпадает с числом распространяющихся нормальных волн $M.$

Входящий в (5) неизвестный вектор ${\mathbf{b}}$ может быть найден лишь с точностью до произвольного комплексного множителя, поэтому для однозначного определения ${\mathbf{b}}$ можно наложить на него дополнительное ограничение ${{\left\| {\mathbf{b}} \right\|}^{2}} = 1$ (исключающее тривиальное решение ${\mathbf{b}} = 0$). В этом случае минимум квадратичной формы ${{b}^{ + }}C(\phi )\,{\mathbf{b}}$ реализуется при условии совпадения ${\mathbf{b}}$ с собственным вектором, отвечающим наименьшему собственному значению ${{\lambda }_{{\min }}}\{ {\kern 1pt} {\kern 1pt} C(\phi )\} $ матрицы $C(\phi ).$ В результате угловые положения источников могут быть найдены из следующего критерия:

Метод оценивания (7) аналогичен алгоритму пониженного ранга RARE [14, 15], используемому в технике частично калиброванных АР, в котором роль числа подрешеток играет число мод $M,$ а вектор амплитудно-фазовой калибровки соответствующих подрешеток заменяется вектором ${\mathbf{b}}.$ Важно подчеркнуть, что матрица ${\mathbf{C}}$ не содержит информации о глубине источников и их удалении от АР; следовательно, критерий (7) позволяет определить направления на источники путем одномерного поиска $J$ максимумов выходной мощности процессора NM-RARE.

2. ПОСТРОЕНИЕ АДАПТИВНОЙ ПРОЦЕДУРЫ NM–RARE

Приведенный алгоритм (7) предполагает априорное знание матрицы направлений $U(\phi )$. Однако на практике в качестве этой матрицы (вследствие неполной информации о канале распространения) используется некоторая оценочная матрица ${{{\mathbf{U}}}_{0}}(\phi ),$ рассчитываемая для номинальных акустических характеристик волновода. При наличии рассогласования между $U(\phi )$ и ${{{\mathbf{U}}}_{0}}(\phi )$ предложенный способ оценивания нуждается в уточнении.

При построении адаптивной процедуры NM‑RARE, основанной на наихудшем сценарии приема, будем предполагать возможность контролируемого отклонения ожидаемой матрицы ${{{\mathbf{U}}}_{0}}(\phi )$ от истинной ${\mathbf{U}}(\phi )$: норма Фробениуса матрицы рассогласования не должна превышать заданную величину: $\left\| {U(\phi ) - {{U}_{0}}(\phi )} \right\|_{F}^{2} \leqslant \,\varepsilon $, где $\varepsilon $ – положительный параметр регуляризации. Адаптация к неизвестным условиям приема состоит в нахождении робастной матрицы ${\mathbf{U}}(\phi ,\varepsilon )$, удовлетворяющей указанному ограничению, условию нормировки и обеспечивающей минимум целевой функции ${{b}^{ + }}{{U}^{ + }}(\phi ){{{\mathbf{\hat {\Pi }}}}_{n}}U(\phi )b$ для всех возможных значений нормированных векторов ${\mathbf{b}}{\kern 1pt} :$

Решение оптимизационной задачи (8) может быть найдено с помощью метода неопределенных множителей Лагранжа аналогично тому, как это сделано в [16, 17]. В итоге искомые положения источников находятся из условия максимума выходной мощности адаптивного процессора NM–RARE:

где а $Tr\,( \cdot )$ означает след матрицы. При $\varepsilon = 0$ функция $f(\phi ,\varepsilon ) = 1$ и, следовательно, неадаптивный метод NM-RARE является частным случаем предложенного способа оценивания.

Отметим также, что в процессе поиска азимутальных углов параметр регуляризации $\varepsilon $ должен удовлетворять неравенству $\varepsilon < 2M[1 - \sqrt {1 - {{V}_{0}}(\phi )} ],$ поскольку в противном случае целевая функция всюду обращается в нуль.

3. ГРАНИЦА КРАМЕРА–РАО И ПРЕДЕЛ УГЛОВОГО РАЗРЕШЕНИЯ

Обозначим через ${\mathbf{\hat {\varphi }}} = {{({{\hat {\phi }}_{1}}, \cdots {{\hat {\phi }}_{J}})}^{T}}$ оценку вектора угловых координат источников (полученную без использования информации о глубинах их погружения и расстояний до приемной АР). Ковариационная матрица ошибки этого вектора удовлетворяет неравенству

где ${\mathbf{CRB}}\,({\mathbf{\varphi }})$ – нижняя граница Крамера–Рао, определяющая предел точности измерения соответствующих координат. Для рассматриваемого сценария (когда форма огибающих излученных сигналов и векторы модовых амплитуд являются неинформативными параметрами) матрица ${\mathbf{CRB}}\,({\mathbf{\varphi }})$ размерности $J \times J$ дается выражением в котором блоки ${\mathbf{F}} \in {{C}^{{J \times J}}},$ ${\mathbf{{\rm K}}} \in {{C}^{{J \times MJ}}}$ и ${\mathbf{\Sigma }} \in {{C}^{{MJ \times MJ}}}$ равны

Здесь ${\mathbf{\Pi }}_{{\mathbf{G}}}^{ \bot } = {{{\mathbf{I}}}_{N}} - {\mathbf{G}}{{({{{\mathbf{G}}}^{ + }}{\mathbf{G}})}^{{ - 1}}}{{{\mathbf{G}}}^{ + }} \in {{C}^{{N \times N}}},$ ${{{\mathbf{\hat {R}}}}_{s}} = ({1 \mathord{\left/ {\vphantom {1 L}} \right. \kern-0em} L})$ × $ \times \,\,\sum\nolimits_{l = 1}^L {\mathbf{s}} {}_{l}{\mathbf{s}}_{l}^{ + } \in {{C}^{{J \times J}}}$ – выборочная сигнальная матрица, ${{{\mathbf{1}}}_{M}} = {{(1, \cdots ,1)}^{T}}$ – вектор размерности $M \times 1$ с единичными компонентами, а символы $ \circ $ и $ \otimes $ означают произведение Адамара и Кронекера, соответственно. Матрицы ${{{\mathbf{D}}}_{{\mathbf{\varphi }}}} \in {{C}^{{N \times J}}}$ и ${{{\mathbf{D}}}_{{\mathbf{\eta }}}} \in {{C}^{{N \times MJ}}}{\kern 1pt} ,$ входящие в (10), определены соотношениями

где ${{{\mathbf{e}}}_{m}} = {{(\underbrace {0, \ldots ,0,1}_m,0, \ldots 0)}^{T}}$представляет собой $m$‑ый столбец единичной матрицы размерности $M \times M.$ Вывод формулы (10) приведен в Приложении.

В случае, когда глубины погружения источников и их расстояния до приемной АР известны априори (следовательно, известны векторы соответствующих модовых амплитуд), формула (10) упрощается и переходит в классическое выражение, полученное в работе [18] применительно к детерминированной модели излучаемых сигналов:

Заметим, что в отличие от (14) соотношение (10) содержит слагаемое $ - {\mathbf{{\rm K}}}{{{\mathbf{\Sigma }}}^{{ - 1}}}{{{\mathbf{{\rm K}}}}^{T}}.$ Поскольку согласно (11) ${\mathbf{\Sigma }}$ является неотрицательно определенной эрмитовой матрицей, то ${\mathbf{{\rm K}}}{{{\mathbf{\Sigma }}}^{{ - 1}}}{{{\mathbf{{\rm K}}}}^{T}} \geqslant 0$ и, следовательно, справедливо неравенство ${\mathbf{CRB}}\,({\mathbf{\varphi }}) \geqslant {\mathbf{CR}}{{{\mathbf{B}}}_{{det}}}\,({\mathbf{\varphi }}).$ Последнее отражает очевидный факт, в соответствии с которым введение в модель дополнительного неизвестного параметра приводит к потере точности искомой оценки (т.е. увеличению ее дисперсии).

Ниже основное внимание будет уделено решению задачи локализации двух акустических источников с близкими угловыми положениями. В рассматриваемой ситуации матрица Крамера–Рао размерности $2 \times 2$ может быть представлена в виде

где ${\mathbf{\Phi }}({\mathbf{\varphi }}) = ({{2L} \mathord{\left/ {\vphantom {{2L} {\sigma _{n}^{2}}}} \right. \kern-0em} {\sigma _{n}^{2}}})\operatorname{Re} ({\mathbf{F}} - {\mathbf{{\rm K}}}{\kern 1pt} {{{\mathbf{\Sigma }}}^{{ - 1}}}{{{\mathbf{K}}}^{ + }})$ имеет смысл информационной матрицы Фишера. Отметим, что для обращения ${\mathbf{\Phi }}\,({\mathbf{\varphi }})$можно воспользоваться формулой в которой $\mu = {{{{f}_{{12}}}{{f}_{{21}}}} \mathord{\left/ {\vphantom {{{{f}_{{12}}}{{f}_{{21}}}} {{{f}_{{11}}}{{f}_{{22}}}}}} \right. \kern-0em} {{{f}_{{11}}}{{f}_{{22}}}}}.$ В свою очередь, выборочная сигнальная матрица ${{{\mathbf{\hat {R}}}}_{s}}$, входящая в ${\mathbf{\Phi }}\,({\mathbf{\varphi }}),$ в данном случае записывается следующим образом где $\sigma _{{1,2}}^{2} = ({1 \mathord{\left/ {\vphantom {1 L}} \right. \kern-0em} L})\sum\nolimits_{l = 1}^L {{{{\left| {{{s}_{{1,2}}}(l)} \right|}}^{2}}} $ – уровни излучения, а параметр $\rho $ описывает степень корреляции между излученными сигналами.

Нас будет интересовать минимальное угловое расстояние между источниками $\delta = \left| {{{\phi }_{1}} - {{\phi }_{2}}} \right|,$ при котором они могут быть одновременно корректно локализованы. Согласно критерию Смита [10], два источника разрешимы, если $\delta $ превосходит среднеквадратичную дисперсию оценки соответствующего углового расстояния:

В свою очередь, дисперсия $\operatorname{var} {\kern 1pt} {\kern 1pt} (\delta ),$ фигурирующая в (15), удовлетворяет неравенству $\operatorname{var} \,(\delta ) \geqslant CRB(\delta )$, где $CRB(\delta )$ – соответствующая граница Крамера–Рао, равная [19]

а $CRB({{\phi }_{1}}) = {{[{\mathbf{CRB}}{\kern 1pt} ({\mathbf{\varphi }})]}_{{1,1}}},$$CRB({{\phi }_{1}},{{\phi }_{2}}) = {{[{\mathbf{CRB}}{\kern 1pt} ({\mathbf{\varphi }})]}_{{1,2}}},$ $CRB{\kern 1pt} {\kern 1pt} ({{\phi }_{2}}) = {{[{\mathbf{CRB}}{\kern 1pt} {\kern 1pt} ({\mathbf{\varphi }})]}_{{2,2}}}.$

Предел углового разрешения находится из условия равенства обеих частей критерия (15) и замены $\operatorname{var} {\kern 1pt} {\kern 1pt} (\delta )$ на $CRB{\kern 1pt} {\kern 1pt} (\delta ).$ В результате искомая величина $\delta $ является наименьшим положительным корнем уравнения

Отметим, что для заданной геометрии АР найденное таким образом $\delta $ зависит от угла прихода, отношения сигнал/шум, объема входной выборки и степени коррелированности излучаемых сигналов.

4. РЕЗУЛЬТАТЫ СТАТИСТИЧЕСКОГО МОДЕЛИРОВАНИЯ

Приведем результаты статистического моделирования, иллюстрирующие работоспособность предложенного алгоритма (9) и сравним его эффективность с традиционным методом сканирования и методом SI, согласно которому векторы сигнального подпространства ${{{\mathbf{\hat {\Psi }}}}_{{\mathbf{s}}}}$ и векторы модового подпространства ${\mathbf{U}}(\phi )$ при $\phi \in ({{\phi }_{1}}, \ldots ,{{\phi }_{J}})$ становятся линейно зависимыми. В результате искомые положения источников могут быть найдены из условия пересечения соответствующих подпространств [5]:

где ${{\lambda }_{{\min }}}\,\{ \cdot \} $ означает минимальное собственное значение матрицы, стоящей в скобках, а ${\mathbf{P}}_{{\mathbf{U}}}^{ \bot }(\phi ) = {{{\mathbf{I}}}_{N}} - {\mathbf{U}}(\phi )$${{[{{{\mathbf{U}}}^{ + }}(\phi ){\mathbf{U}}(\phi )]}^{{ - 1}}}{{{\mathbf{U}}}^{ + }}(\phi )$ – проекционная матрица размерности $N \times N.$

Подчеркнем, что метод сканирования не учитывает волноводный характер распространения звука, а процедура (17) опирается на априорное знание среды распространения.

В качестве примера рассмотрим мелководный канал глубины $H = 160$ м с характерной летней гидрологией, изображенной на рис. 1б, в котором звуковое поле создается двумя детерминированными источниками одинаковой мощности, излучающими узкополосные сигналы с несущей частотой 250 Гц. В рамках численного эксперимента дно моделировалось жидким поглощающим полупространством с плотностью ${{\rho }_{b}} = 1.8$ г/см³, скоростью звука ${{c}_{b}} = 1750$ м/с и коэффициентом поглощения $\beta = 0.13$ дБ/$\lambda $, а при расчете ожидаемой матрицы ${{{\mathbf{U}}}_{0}}(\phi )$ в качестве номинальных геоакустических параметров дна использовались значения $H = 162.5$ м, ${{\rho }_{b}} = 1.75$ г/см³, ${{c}_{b}} = 1725$ м/с и $\beta = 0.1$ дБ/$\lambda $. Для рассматриваемой несущей частоты и указанных акустических характеристик канала полное число мод $M$ составляло 28.

Предполагалось, что источники находятся на горизонтах 40 и 70 м и удалены соответственно на расстояния 15 и 10 км от горизонтальной приемной антенны, состоящей из $N = 20$ элементов, расположенных через $3$ м на глубине $155$ м. Направление на первый источник ${{\phi }_{1}} = 30^\circ ,$ а угловое положение второго задавалось в виде ${{\phi }_{2}} = {{\phi }_{1}} + \delta ,$ где ${\delta } \in (0.2^\circ ...6^\circ ).$ При вычислениях коэффициент корреляции $\rho $ брался равным 0.4, а входные отношения сигнал/шум, определяемые соотношениями

считались одинаковыми, т. е. $SN{{R}_{1}} = SN{{R}_{2}} = SNR.$

Обратим внимание, что для данных $N,$ $J$ и $M$ условие (6) не выполняется, и для реализации предложенного метода необходимо ограничить число мод, участвующих в процессе локализации (причем максимально допустимое значение используемых мод не должно превышать $N - J = 18$). Ниже при построении матрицы ${{{\mathbf{U}}}_{0}}(\varphi )$ учитывались первые 16 нормальных волн волновода. В частности, при таком выборе $M$ применение алгоритма (9) для $\varepsilon = 0.25$ обеспечивает минимум среднеквадратических ошибок (СКО) оценивания угловых положений источников, отстоящих друг от друга на расстояние $\delta = 6^\circ ,$ как показано на рис. 2. Соответствующие ошибки рассчитывались по формуле

в которой $\hat {\phi }_{j}^{{(q)}}$ – оценка угловой координаты $j$-го источника для $q$-ой реализации вектора наблюдения, $J = 2,$ а общее число независимых реализаций $Q$ бралось равным $1000.$ Выборочная ковариационная матрица формировалась по $L = 100$ временным отсчетам, а входное $SNR$ составляло 0 дБ.

Для используемых методов оценивания на рис. 3а и 3б изображены угловые зависимости нормированной (на максимальное значение) мощности на выходе АР, принимающей сигналы от двух источников, разнесенных друг относительно друга на расстояние $6^\circ $ и $2^\circ $, соответственно. Кривая 1 на рис. 3 отвечает традиционному способу сканирования, кривая 2 – методу SI, а кривая 3 соответствует алгоритму NM-RARE (9). Из рис. 3б видно, что в данном примере лишь адаптивный метод в состоянии различить два источника и одновременно оценить искомые направления.

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

На рис. 4а для рассматриваемых способов обработки представлены результаты расчета ${{P}_{{RES}}}$ в зависимости от величины углового расстояния между источниками. Очевидно, что наилучшие потенциальные возможности демонстрирует предложенный метод локализации (9), позволяющий с вероятностью 0.9 разрешить источники, отстоящие друг относительно друга на величину порядка $1.5^\circ $.

На рис. 4б показано поведение среднеквадратической ошибки измерения угловых координат в зависимости от расстояния между источниками. Из приведенных кривых следует, что применение алгоритма NM-RARE позволяет значительно повысить точность оценивания по сравнению с неадаптивными методами.

В заключение этого раздела приведем результаты расчета предела углового разрешения, являющегося наименьшим положительным корнем уравнения (16). Для рассматриваемой постановки задачи на рис. 5а изображена зависимость величины $\delta $ от входного $SNR$ (при этом выборочная ковариационная матрица формировалась по $L = 100$ временным отсчетам), а на рис. 5б – зависимость $\delta $ от числа выборок $L$ при $SNR = 0$ дБ. Параметром кривых является угловое положение ${{\phi }_{1}}$ фиксированного источника. Из этого рисунка видно, что минимально возможное расстояние между источниками является монотонно убывающей функцией $SNR$ и $L,$ при этом разрешающая способность антенны снижается с ростом угла прихода. Важно подчеркнуть, что полученный предел разрешения является фундаментальной величиной, независящей от вида используемого алгоритма.

ЗАКЛЮЧЕНИЕ

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

На основе критерия Смита определено наименьшее угловое расстояние между источниками, при котором возможна их корректная локализация, в зависимости от входного отношения сигнал/шум и объема входной выборки.

Работа выполнена при финансовой поддержке РНФ (грант № 20-19-00383).