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

Постановка задачи

Пусть имеется комплексная величина

Пусть ее квадратурные компоненты независимо распределены по гауссовскому закону с некоторым средним значением, определяемым символом модуляции, и одинаковой дисперсией , определяемой уровнем шума системы:

Здесь

определяет среднее значение соответствующей квадратуры, а

среднеквадратическое отклонение (СКО).

Аргументу КА соответствует аргумент комплексного числа

Здесь коэффициент определяется четвертью комплексной плоскости.

Задача: найти плотность распределения вероятности случайной величины "аргумент КА". Однако, ввиду влияния коэффициента это сделать непросто, поэтому сначала найдем как распределен соответствующий арктангенс, который можно назвать "фазой ", а затем сделаем "финт ушами". Заметим, что для различения аргумента и фазы введены разные обозначения.

Решение задачи

Первым делом найдем распределение не фазы КА, а тангенса фазы

В нашем случае имеем двумерное отображение:

Здесь случайные величины

преобразуются в величины

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

Тогда плотность вероятности тангенса фазы

Здесь - гауссовские плотности вероятностей для соответствующих квадратур. Была сделана подстановка

Интегрирование по лишней переменной позволяет "маржинализировать" совместную плотность вероятности. Использование системы компьютерной алгебры (СКА) дает следующий ответ:

где

вспомогательная функция-инфлюенция, показывающая вклад средних значений. Если средние значения равны нулю, то плотность вероятности тангенса фазы превращается в стандартное распределение Коши. Методика отыскания плотности для двумерного случая взята из Википедии.

Далее найдем плотность вероятности непосредственно фазы КА, воспользовавшись известной формулой

В нашем случае последовательно имеем:

Заметим, что фаза распределена на интервале

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

Подставим вместо величину

и вместо будем писать просто . Далее введем отношение сигнал-шум

Тогда плотность вероятности фазы будет иметь более простой вид:

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

Наконец, пришло время "финта ушами". Так как из-за влияния арктангенса фаза распределена на интервале длиной , то в найденной плотности вероятности сделаем прямое деление фазы пополам и просто заменим фазу на аргумент, дополнительно перенормировав плотность, а также откоррекрировав SNR специальным коэффициентом, который подобран на основании анализа моделирования квадратур:

Здесь аргумент распределен на полном интервале длиной

Наглость - второе счастье!

import numpy as np
from scipy import special
import matplotlib.pyplot as plt
def F(phi, phi0, snr):
    return (
        np.sqrt(2)
        / np.sqrt(1 + np.tan(0.5 * phi0) ** 2)
        * snr
        * (1 + np.tan(0.5 * phi) * np.tan(0.5 * phi0))
        / np.sqrt(1 + np.tan(0.5 * phi) ** 2)
    )
def W(phi, phi0, snr):
    _F = F(phi, phi0, snr)
    return (
        (1 / 2 / np.pi)
        * np.exp(-2 * snr**2)
        * (1 + np.sqrt(np.pi) * _F * np.exp(_F**2) * special.erf(_F))
    )
# Параметры генерации
size = 8192  # Количество элементов в векторах
mean1, std1 = 0, 1  # Среднее и стандартное отклонение для первого вектора
mean2, std2 = +1, 1  # Среднее и стандартное отклонение для второго вектора
assert std1 == std2
snr = max(mean1, mean2) / std1
phi0 = np.angle(mean1 + 1j * mean2)
# Генерация двух векторов нормальных случайных величин
vector1 = np.random.normal(loc=mean1, scale=std1, size=size)
vector2 = np.random.normal(loc=mean2, scale=std2, size=size)
argument = np.angle(vector1 + 1j * vector2)
# Построение гистограммы
plt.figure(figsize=(8, 6))
plt.hist(
    argument,
    bins=30,
    density=True,
    alpha=0.7,
    color="blue",
    label="Historgam",
    edgecolor="black",
)
plt.xlim([-np.pi, +np.pi])
xmin, xmax = plt.xlim()  # Get the x-axis limits from the histogram
phi = np.linspace(xmin, xmax, 100)
pdf = W(phi, phi0, snr)
plt.plot(phi, pdf, "r", label="PDF")
plt.title("Распределение аргумента комплексной амплитуды")
plt.xlabel("Значение аргумента, рад.")
plt.legend()
plt.tight_layout()
plt.show()