Случайные процессы и шумы

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

Рассмотрим способы моделирования и анализа случайных процессов в скрипте Engee:

базовые функции для генерации случайных чисел
функция для визуализации плотности распределения
частотный состав шумовых сигналов
осреднение сигнала с аддитивным белым Гауссовским шумом

# Pkg.add("WAV")

Базовые функции для генерации, визуализации и анализа

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

Создадим матрицу целых чисел в диапазоне от 1 до 100 с равномерным распределением при помощи функции rand:

rand(1:100,(3,4))

3×4 Matrix{Int64}:
 55  38  68  18
 91   6  90  26
 50  54  37   8

Создадим вектор из 10000 чисел с плавающей запятой в диапазоне от 0 до 1, и визуализируем его распределение при помощи функции histogram. У данной функции можно опционально настроить количество "бинов", то есть "ворот" для границ значений случайного вектора. Убедимся, что выход функции rand даёт нам равномерное распределение случайной величины:

uniform = rand(10000);
histogram(uniform, nbins=20)

Теперь зададим вектор из 10000 чисел с нормальным (Гауссовским) распределением (единичная дисперсия, нулевое среднее) при помощи функции randn. Также визуализируем его распределение при помощи функции histogram, но на этот раз пронормируем её таким образом, чтобы сумма всех "бинов" равнялась единице. Таким образом мы получаем вероятность попадания числа в векторе в те или иные границы:

gaussian = randn(10000);
histogram(gaussian, nbins=50, normalize=:probability, color=:cyan4)

Посчитаем основную статистику случайных векторов с равномерным и Гауссовским распределениями:

duni, muni = EngeeDSP.Functions.var(uniform);
dgauss, mgauss = EngeeDSP.Functions.var(gaussian);
println("Равномерное распределение: среднее - ", muni, " дисперсия - ", duni)
println("Нормальное распределение: среднее - ", mgauss, " дисперсия - ", dgauss)

Равномерное распределение: среднее - 0.49967003217802675 дисперсия - 0.0834509615643607
Нормальное распределение: среднее - -0.005309704120389261 дисперсия - 1.0058092035592265

Нормальное распределение из суммы равномерных

Центральная предельная теорема (ЦПТ) — это фундаментальное правило теории вероятностей. Оно говорит, что сумма большого числа независимых случайных величин имеет распределение, которое очень близко к нормальному (похоже на колокол), независимо от того, какими были исходные распределения слагаемых. Главное условие — ни одно из слагаемых не должно доминировать над остальными.

Слагаемых = 8 # @param {type:"slider",min:1,max:8,step:1}
N = Слагаемых;

noisemat = rand(10000,N);
sumsignal = sum(noisemat,dims=2);
histogram(sumsignal, nbins=100, 
                    normalize=:probability, 
                    color=:olive, 
                    leg=false,
                    size=(800,250))

Белый Гауссовский шум

Воспользуемся функцией awgn из состава библиотеки EngeeComms для добавления подобного шума к синусоидальному сигналу. Укажем в качестве входного аргумента функции дополнительнй параметр, определяющий отношение сигнал/шум, в значение 8 дБ:

t = 0:1/1000:1;
sine_wave = cos.(10*2pi.*t);
plot(t,sine_wave, grid=true, leg=false,
        xlim=(0,0.2),
        ylim=(-2,2),
        lw=3)
sine_noise, variance = EngeeComms.Functions.awgn(sine_wave,8);
plot!(t,sine_noise, title = "Сигнал с аддитивным белым гауссовским шумом")

Выделим шумовую составляющую в отдельную переменную noise, и отобразим её отсчёты во времени:

noise = sine_noise .- sine_wave;
plot(t,noise,l=:stem, leg=false,
            xlim = (0,0.1),
            ylim = (-1.5,1.5),
            lw=2,
            m=:c,
            ms=2,
            grid=true,
            title="Аддитивный белый гауссовский шум")

Построим спектр этого шума при помощи встроенной функции pspectrum. Убедимся, что шум действительно "белый":

EngeeDSP.Functions.pspectrum(noise; out=:plot)

([6.769341406980952e-9; 4.1379065778417257e-7; … ; 9.60852708997955e-6; 4.978659968030431e-6;;], [0.0; 0.0007671776931843207; … ; 3.1408254758966088; 3.141592653589793;;], Any[], Plot{Plots.PlotlyJSBackend() n=1})

Также убедимся, что он Гауссовский, то есть имеет нормальное распределение:

histogram(noise, nbins=100, grid=true)

Послушаем пример белого шума с нормальным распределением при помощи функции audioplayer:

include("audioplayer.jl") 
beatuiful_sound_to_relax = randn(Int(1e5));
audioplayer(beatuiful_sound_to_relax,44100)

Осреднение сигнала с белым Гауссовским шумом

Создадим "длинный" сигнал - синусоиду с АБГШ. Затем разобъём её на когерентные (с одинаковой начальной и конечной фазой) отрезки по 200 отсчётов. Убедимся, что они успешно накладываются друг на друга:

tt = 0:1/1000:4-1/1000;
long_clean = cos.(10*2pi.*tt);
long_noisy,_ = EngeeComms.Functions.awgn(long_clean,10);
noisy_matrix = reshape(long_noisy,(200,20));
plot(tt[1:200],noisy_matrix[:,1])
plot!(tt[1:200],noisy_matrix[:,11])

Осреднение = 1 # @param {type:"slider",min:1,max:20,step:1}
M = Осреднение;

avg = sum(noisy_matrix[:,1:M], dims=2) ./ M;

plot(tt[1:200],avg, leg=false,
            xlim=(0,0.2),
            grid=true,
            lw=3,
            title = "Осреднение в  $M  раз")

Случайные процессы и шумы