Калибровка уровня записи

Психоакустические измерения приводят к самым точным результатам с калиброванным уровнем входа микрофона. Использовать calibrateMicrophone чтобы совпадать с вашим уровнем записи к чтению метра SPL, можно использовать тональный источник на 1 кГц (такой как онлайновое тональное приложение генератора или сотового телефона) и калиброванный метр SPL. SPL калибровочного тона на 1 кГц должен быть достаточно громким, чтобы доминировать над любым фоновым шумом. Для калибровочного примера с помощью MATLAB в качестве тонального источника на 1 кГц смотрите calibrateMicrophone.

Симулируйте тональную запись и включайте некоторый фоновый шум. Примите чтение метра SPL 83,1 дБ (C-weighted).

FS = 48e3;
t = (1:2*FS)/FS;
s = rng('default');
testTone = 0.46*sin(2*pi*t*1000).' + .1*pinknoise(2*FS);
rng(s)

splMeterReading = 83.1;

Чтобы вычислить калибровочный уровень цепи записи, вызовите calibrateMicrophone и задайте тестовый тон, частоту дискретизации, чтение SPL и взвешивание частоты метра SPL. Чтобы компенсировать возможный фоновый шум и произвести точный калибровочный уровень, совпадайте с установкой взвешивания частоты метра SPL.

calib = calibrateMicrophone(testTone,FS,splMeterReading,"FrequencyWeighting","C-weighting");

Уровни звукового давления (SPL)

Если у вас есть калибровочный фактор для вашей цепи записи, можно сделать акустические измерения. При использовании физического метра вы ограничиваетесь настройками, выбранными в течение времени измерения. С splMeter объект, можно изменить настройки после того, как запись была сделана. Это дает возможность экспериментировать с различным временем и опциями взвешивания частоты.

Загрузите запись механизма и создайте соответствующий временной вектор.

[x,FS] = audioread('Engine-16-44p1-stereo-20sec.wav');
x = x(1:8*FS,1); % use only channel 1 and keep only 8 seconds.
t = (1:size(x,1))/FS;

Создайте splMeter возразите и выберите C-weighting, быстрое время, взвесив и 0,2 вторых интервала для пикового измерения SPL.

spl = splMeter("CalibrationFactor",calib, ...
               "FrequencyWeighting","C-weighting", ...
               "TimeWeighting","Fast", ...
               "TimeInterval",0.2, ...
               "SampleRate",FS);

Постройте SPL и пиковый SPL.

[LCF,~,LCpeak] = spl(x);
plot(t,LCpeak,t,LCF)
legend('LCpeak','LCF','Location',"southeast")
title('SPL Measurement of Engine Noise')
xlabel('Time (seconds)')
ylabel('SPL (dB)')
ylim([70 95])
grid on

Психоакустические метрики

Уровень громкости

Контроль SPL важен для податливости охраны труда. Однако измерения SPL не отражают громкость, как воспринято фактическим прослушивателем. acousticLoudness уровни громкости мер, как воспринято человеческим прослушивателем с нормальным слушанием (никакое нарушение слуха). acousticLoudness функционируйте также показывает, какие диапазоны частот способствуют больше всего перцепционной сенсации громкости.

Используя тот же калибровочный уровень как прежде, и принятие записи свободного поля (значение по умолчанию), стационарная громкость графика.

acousticLoudness(x,FS,calib)

Громкость составляет 23,8 сона, и большая часть шума ниже 3.3 (Шкала коры). Преобразуйте 3.3 Коры в Гц с помощью bark2hz

fprintf("Loudness frequency of 3.3 Bark: %d Hz\n",round(bark2hz(3.3),-1));

Loudness frequency of 3.3 Bark: 330 Hz

acousticLoudness функция возвращает воспринятую громкость в сонах. Чтобы изучить измерение сона, сравните его с SPL (дБ) чтение. Сигнал с громкостью 60 фонов, как воспринимают, является столь громким как тон на 1 кГц, измеренный в SPL на 60 дБ. Преобразование 23,8 сонов к фонам с помощью sone2phon демонстрирует, что восприятие громкости шума механизма является столь громким как тон на 1 кГц, измеренный в SPL на 86 дБ.

fprintf("Equivalent 1 kHz SPL: %d phons\n", round(sone2phon(23.8)));

Equivalent 1 kHz SPL: 86 phons

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

[sone,spec] = acousticLoudness(x,FS,calib);
barks = 0.1:0.1:24; % Bark scale for ISO 532-1 loudness
hz = bark2hz(barks);
specPhon = sone2phon(spec);
semilogx(hz,specPhon)
title('Specific Loudness')
subtitle(sprintf('Loudness = %.1f phons',sone2phon(sone)))
xlabel('Frequency (Hz)')
ylabel('Loudness (phons/Bark)')
xlim(hz([1,end]))
grid on

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

acousticLoudness(x,FS,calib,'TimeVarying',true,'TimeResolution','high')

Постройте определенную громкость с частотой в Гц (логарифмическая шкала).

[tvsoneHD,tvspecHD,perc] = acousticLoudness(x,FS,calib,'TimeVarying',true,'TimeResolution','high');
tvspec = tvspecHD(1:4:end,:,:); % for standard resolution measurements
spectimeHD = 0:5e-4:5e-4*(size(tvspecHD,1)-1); % time axis for loudness output
clf; % do not reuse the previous subplot
surf(spectimeHD,hz,sone2phon(tvspecHD).','EdgeColor','interp');
set(gca,'View',[0 90],'YScale','log','YLim',hz([1,end]));
title('Specific Loudness (HD)')
zlabel('Specific Loudness (phons/Bark)')
ylabel('Frequency (Hz)')
xlabel('Time (seconds)')
colorbar

Уровень резкости

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

sharp = acousticSharpness(spec)

sharp = 1.1512

Розовый шум имеет резкость 2 acums. Это означает, что шум механизма склоняется к низким частотам.

Постройте изменяющуюся во времени резкость.

acousticSharpness(x,FS,calib,'TimeVarying',true);

Сила колебания

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

Во-первых, посмотрите на то, какие частоты модуляции присутствуют в сигнале.

N = 2^nextpow2(size(x,1));
xa = abs(x); % Use the rectified signal
pspectrum(xa-mean(xa),FS,'FrequencyLimits',[0 80],'FrequencyResolution',1)
title('Modulation Frequencies')

Частота модуляции достигает максимума на уровне 24,9 Гц. Ниже 30 Гц модуляция воспринята доминируя как колебание. Существует второй пик на уровне 49,7 Гц, который находится в области значений шероховатости.

Используйте acousticFluctuation вычислить воспринятую силу колебания. Шум механизма является относительно постоянным в этой записи, таким образом, у нас есть алгоритм, автоматически обнаруживают самую слышимую частоту колебания (fMod).

acousticFluctuation(x,FS,calib)

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

[vacil,spec,fMod] = acousticFluctuation(tvspec,'ModulationFrequency',24.9);
clf; % do not reuse previous subplot
flucHz = bark2hz(0.5:0.5:23.5);
spectime = 0:2e-3:2e-3*(size(spec,1)-1);
surf(spectime,flucHz,spec.','EdgeColor','interp');
set(gca,'View',[0 90],'YScale','log','YLim',flucHz([1,end]));
title('Specific Fluctuation Strength')
zlabel('Specific Fluctuation Strength (vacils/Bark)')
ylabel('Frequency (Hz)')
xlabel('Time (seconds)')
colorbar

Шероховатость

Используйте acousticRoughness функция, чтобы вычислить воспринятую шероховатость сигнала. Позвольте алгоритму автоматически обнаружить самую слышимую частоту модуляции (fMod).

acousticRoughness(x,FS,calib)

Интерпретируйте результаты в Герц вместо Коры. Чтобы уменьшать расчеты, снова используйте ранее вычисленную изменяющуюся во времени определенную громкость. Задайте частоту модуляции.

[asper,specR,fModR] = acousticRoughness(tvspecHD,'ModulationFrequency',49.7);
clf; % do not reuse previous subplot
rougHz = bark2hz(0.5:0.5:23.5);
surf(spectimeHD,rougHz,specR.','EdgeColor','interp');
set(gca,'View',[0 90],'YScale','log','YLim',rougHz([1,end]));
title('Specific Roughness')
zlabel('Specific Roughness (aspers/Bark)')
ylabel('Frequency (Hz)')
xlabel('Time (seconds)')
colorbar

Качество звука

Для полной оценки качества звука объедините предыдущие метрики, чтобы произвести психоакустическую метрику раздражения (заданный Zwicker и Fastl). Отношение следующие:

Количественное описание было разработано с помощью результатов психоакустических экспериментов:

с:

В этом примере резкость была меньше 1.75, таким образом, это не влияющий фактор. Поэтому можно установить ${\mathit{w}}_{\mathit{s}}$ обнулять.

Громкость процентили, ${\mathit{N}}_5$, второе значение, возвращенное третьим выходом acousticLoudness когда "TimeVarying" установлен в true.

N5 = perc(2);

Вычислите среднюю силу колебания, игнорирующую первую секунду сигнала.

f = mean(vacil(501:end,1));

Вычислите среднюю шероховатость, игнорирующую первую секунду сигнала.

r = mean(asper(2001:end,1));

Вычислите психоакустическую метрику раздражения.

pa = N5 * (1 + abs(2.18/(N5^.4)*(.4*f+.6*r)))

pa = 25.8848

Эффект улучшенной звукоизоляции

Измерьте удар улучшенной звукоизоляции на измеренном SPL и воспринятом шуме.

Симуляция Используя графический набор фильтров EQ

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

geq = graphicEQ("Bandwidth","1 octave","SampleRate",FS,"Gains",[-0.5 -1.25 -3.4 -7 -8.25 -8.4 -8 -7 -6.4 -5.6]);
cf = getCenterFrequencies(splMeter("Bandwidth","1 octave"));

Отобразите частотную характеристику graphicEQ объект.

[B,A] = coeffs(geq);
sos = [B;A].';
[H,w] = freqz(sos,2^16,FS);
semilogx(w,db(abs(H)))
title('Frequency Response of Soundproofing Simulation Filter')
ylabel('Relative SPL (dB)')
xlabel('Frequency (Hz)')
xlim(cf([1,end]))
grid on

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

x2 = geq(x);

Сравните SPL с и без звукоизоляции. Снова используйте те же настройки метра SPL, но используйте отфильтрованную запись.

reset(spl)
[LCFnew,~,LCpeaknew] = spl(x2);
plot(t,LCpeak,t,LCF,t,LCpeaknew,t,LCFnew)
legend('LCpeak (original)', 'LCF (original)', ...
       'LCpeak (with soundproofing)', ...
       'LCF (with soundproofing)', ...
       'Location','southeast')
title('SPL Measurement of Engine Noise')
xlabel('Time (seconds)')
ylabel('SPL (dB)')
ylim([70 95])
grid on

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

acousticLoudness(x2,FS,calib)

Громкость, уменьшенная от 23,8 до 16,3 сонов. Однако может быть легче интерпретировать громкость в фонах. Преобразуйте единицы сона к фонам с помощью sone2phon.

fprintf("Loudness without soundproofing:   \t%.1f phons\n",sone2phon(23.8));

Loudness without soundproofing:   	85.7 phons

fprintf("Loudness with added soundproofing:\t%.1f phons\n",sone2phon(16.3));

Loudness with added soundproofing:	80.3 phons

fprintf("Perceived noise reduction:\t\t%.1f phons (dB SPL at 1 kHz)\n",sone2phon(23.8)-sone2phon(16.3));

Perceived noise reduction:		5.5 phons (dB SPL at 1 kHz)

После звукоизоляции, acousticLoudness показывает, что восприятие шума механизма на приблизительно 5,5 дБ более тихо. Человеческое восприятие звука ограничивается в очень низких частотах, где большая часть шума механизма. Звукоизоляция является более эффективной на более высоких частотах.

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

[~,spec2hd,perc2] = acousticLoudness(x2,FS,calib,"TimeVarying",true,"TimeResolution","high");
spec2 = spec2hd(1:4:end,:,:);
shp = acousticSharpness(spec2,'TimeVarying',true);
new_sharp = mean(shp(501:end))

new_sharp = 1.0796

Резкость уменьшилась, потому что звукоизоляция является более эффективной при высокочастотном затухании. Это ниже порога 1,75, таким образом, это проигнорировано для фактора раздражения.

Теперь вычислите среднее значение силы колебания и шероховатости.

vacil2 = acousticFluctuation(spec2);
f2 = mean(vacil2(501:end,1));
asper2 = acousticRoughness(spec2hd);
r2 = mean(asper2(2001:end,1));

Вычислите новый психоакустический фактор раздражения. Это уменьшилось, с 25,9 до 17,7.

N5hp = perc2(2); % N5 with soundproofing
pahp = N5hp * (1 + abs(2.18/(N5hp^.4)*(.4*f2+.6*r2)))

pahp = 17.7364

Ссылки

[1] Zwicker, Эберхард и Хьюго Фэстл. Психоакустика: факты и модели. Издание 22. Springer Science & Business Media, 2013.