Правильный способ добавления шума к сигналу


Во многих областях я обнаружил, что, добавляя шум, мы упоминаем некоторые спецификации, такие как нулевое среднее и дисперсия. Мне нужно добавить AWGN, цветной шум, равномерный шум различной SNR в Дб. Следующий код показывает, как я создавал и добавлял шум. Я знаю о функции awgn(), но это своего рода черный ящик, не зная, как добавляется шум. Итак, может кто-нибудь объяснить правильный способ генерирования и добавления шума. Спасибо

SNR = [-10:5:30]; %in Db
snr = 10 .^ (0.1 .* SNR);

for I = 1:length(snr)
    noise = 1 / sqrt(2) * (randn(1, N) + 1i * randn(1, N));
    u = y + noise .* snr(I);
end
4 10

4 ответа:

Я добавляю еще один ответ, поскольку мне кажется, что Стивен не совсем прав, а предложение Хорхлера заглянуть внутрь функции awgn является хорошим.

Либо MATLAB, либо Octave (в Communications toolbox) имеют функцию awgn, которая добавляет (белый гауссовский) шум для достижения желаемого уровня мощности сигнала к шуму; ниже приводится соответствующая часть кода (из функции Octave):

  if (meas == 1)  % <-- if using signal power to determine appropriate noise power
    p = sum( abs( x(:)) .^ 2) / length(x(:));
    if (strcmp(type,"dB"))
      p = 10 * log10(p);
    endif
  endif

  if (strcmp(type,"linear"))
    np = p / snr;
  else   % <-- in dB
    np = p - snr;
  endif

  y = x + wgn (m, n, np, 1, seed, type, out);

Как вы можете видеть по пути p (мощность входных данных) является однако ответ Стивена, судя по всему, не совсем верен.

Вы можете попросить функцию вычислить полную мощность вашего массива данных и объединить ее с требуемым значением s/n, которое вы предоставляете, чтобы вычислить соответствующий уровень мощности добавленного шума. Вы делаете это, передавая строку "measured" среди дополнительных входных данных, например (см. здесь для документации Octave или здесь для документации MATLAB):

     y = awgn (x, snr, 'measured')

Это приводит в конечном счете к meas=1 и поэтому meas==1 истинно в приведенном выше коде. Функция awgn затем использует переданный ей сигнал для вычисления мощности сигнала, и из этого и желаемого s/n она затем вычисляет соответствующий уровень мощности для добавленного шума.

Как далее поясняется в документации

По умолчанию предполагается, что snr и pwr находятся в dB и dBW соответственно. Это поведение по умолчанию может быть выбрано с типом, установленным в "децибел". В случае, когда тип установлен в "линейный", предполагается pwr быть в Ваттах и snr это соотношение.

Это означает, что вы можете передать отрицательное или 0 дБ значение snr. Результат также будет зависеть от других передаваемых параметров, таких как строка "измерено".

Для случая MATLAB я предлагаю прочитать документацию , она объясняет, как использовать функцию awgn в различных сценариях. Обратите внимание, что реализации в Octave и MATLAB не идентичны, вычисление мощности шума должно быть одинаковым, но может быть различным опции.

А вот соответствующая часть из wgn (названная выше по awgn):

  if (strcmp(type,"dBW"))
    np = 10 ^ (p/10);
  elseif (strcmp(type,"dBm"))
    np = 10 ^((p - 30)/10);
  elseif (strcmp(type,"linear"))
    np = p;
  endif

  if(!isempty(seed))
    randn("state",seed);
  endif

  if (strcmp(out,"complex"))
    y = (sqrt(imp*np/2))*(randn(m,n)+1i*randn(m,n)); % imp=1 assuming impedance is 1 Ohm
  else
    y = (sqrt(imp*np))*randn(m,n);
  endif

Если вы хотите проверить мощность вашего шума (np), функции awgn и awg предполагают следующие соотношения:

  np = var(y,1);        % linear scale
  np = 10*log10(np);    % in dB 

Где var(...,1) - это популяционная дисперсия для шума y.

Большинство ответов здесь забывают, что SNR определяется в децибелах. Поэтому вы не должны столкнуться с ошибкой "деление на 0", потому что вы действительно должны делить на 10^(targetSNR/10), которая никогда не является отрицательной или нулевой для реального targetSNR.

Вы можете использовать randn() для генерации вектора шума 'awgnNoise' нужной длины. Затем, учитывая заданное значение SNR, вычислите мощность исходного сигнала и мощность вектора шума "awgnoise". Получите правильный коэффициент масштабирования амплитуды для вектора шума и просто масштабируйте его.

Следующий код является примером искажения сигнала с белым шумом, предполагая, что входной сигнал равен 1D и имеет реальное значение.
function out_signal = addAWGN(signal, targetSNR)
sigLength = length(signal); % length
awgnNoise = randn(size(signal)); % orignal noise
pwrSig = sqrt(sum(signal.^2))/sigLength; % signal power
pwrNoise = sqrt(sum(awgnNoise.^2))/sigLength; % noise power

scaleFactor = (pwrSig/pwrNoise)/targetSNR; %find scale factor
awgnNoise = scaleFactor*awgnNoise; 
out_signal = signal + awgnNoise; % add noise

Будьте осторожны с фактором sqrt (2), Когда вы имеете дело с сложный сигнал, если вы хотите генерировать реальную и воображаемую часть отдельно.

Эта проблема "не следует делить на 0" может быть легко решена, если вы добавите условие, чтобы проверить, если targetSNR равен 0, и делать это только если он не равен 0. Когда ваш целевой SNR равен 0, это означает, что это чистый шум.

function out_signal = addAWGN(signal, targetSNR)
sigLength = length(signal); % length
awgnNoise = randn(size(signal)); % orignal noise
pwrSig = sqrt(sum(signal.^2))/sigLength; % signal power
pwrNoise = sqrt(sum(awgnNoise.^2))/sigLength; % noise power
if targetSNR ~= 0
   scaleFactor = (pwrSig/pwrNoise)/targetSNR; %find scale factor
   awgnNoise = scaleFactor*awgnNoise; 
   out_signal = signal + awgnNoise; % add noise
else
   out_signal = awgnNoise; % noise only
end