What is sound pressure level (SPL) in pspectrum code

I used the pspectrum function to perform frequency analysis on sound pressure data.

My goal is to remove frequency components below a certain threshold and then convert the remaining data into decibels(SPL).

For sound pressure 'P', Sound Pressure Level (SPL) is defined as:

10 *

​

However, I’m confused about how to compute SPL correctly when using the pspectrum function.

Which of the following values should be used as the input pa_data to pspectrum so that the output can be converted into proper SPL values?

1. pa_data = raw pressure P, then divide the output 'ps' by and finally compute 10.

2. pa_data = P/ (Normalize the pressure by ​ before input,) nd then compute 10 .

3. pa_data = 10 * and treating the output as SPL.

4. Others

Below is the part of the code that follows method #2 for reference. (A-weighting process is included)

for i = 1:num_files

file_name = files(i).name;

file_path = fullfile(file_dir, file_name);

data = readmatrix(file_path);

Pa = data(11:end, 1);

% time vector

t = (0:length(data_float)-1) / sample_rate;

% Normalization for Pa → dB

Pa_data = Pa / 0.00002;

% pspectrum

[ps, freq, time] = pspectrum(Pa_data, sample_rate, 'TimeResolution', 0.1, 'spectrogram');

% dB

ps_dB = 10 * log10(ps + 1e-12); % Avoid log(0)

% A-weighting

f_A_weighting = [10 20 25 31.5 40 50 63 80 100 125 160 200 250 315 400 500 630 800 1000 1250 1600 2000 2500 3150 4000 5000 6300 8000 10000 12500 16000 20000];

A_weights = [-70.4 -50.4 -44.5 -39.4 -34.6 -30.2 -26.2 -22.4 -19.1 -16.1 -13.4 -10.9 -8.6 -6.6 -4.8 -3.2 -1.9 -0.8 0 0.6 1 1.2 1.3 1.2 1 0.5 -0.1 -1.1 -2.5 -4.3 -6.6 -9.3];

interp_weights = interp1(f_A_weighting, A_weights, freq, "linear", "extrap");

a_spectrogram_dB = ps_dB + interp_weights(:);

Answer 1

[William Rose]

@지호,

I assume, based on your code, that your original signal has units of Pa.

Thererefore let us make a seven-second-long simulated signal which has units of Pa. It consists of sinusoids at frequencies 500, 1000, and 2000 Hz. Each sinusoid by itself has an intensity which is steady for 1 second intervals. The intensities (in Pa) correspond to 0,30,60,90,60,30,0 dB SPL. (Seven levels, one level per second.)

% Define constants
p0=2e-5;                % pressure (Pa, rms) corresponding to 0 dB SPL
p30=10^(30/20)*p0;      % pressure (Pa, rms) corresponding to 30 dB SPL
p60=10^(60/20)*p0;      % pressure (Pa, rms) corresponding to 60 dB SPL 
p90=10^(90/20)*p0;      % pressure (Pa, rms) corresponding to 90 dB SPL 
fs=1e4;                 % sampling rate (Hz)
f1=500;                 % sound frequency 1 (Hz)
f2=1000;                % sound frequency 2 (Hz)
f3=2000;                % sound frequency 3 (Hz) 
% Create 1-second-long time vector
t1=(1:fs)/fs;           % time vector, 1 second long
% Create 1-second-long waveform with three frequencies, each of which has amplitude=1.
sin3freq=sin(2*pi*f1*t1)+sin(2*pi*f2*t1)+sin(2*pi*f3*t1);
% Create 1-second-long waveforms with different SPLs
% Use sqrt(2) because dB SPL uses RMS pressure, 
% and amplitude A=sqrt(2)*RMS value.
x0=sqrt(2)*p0*sin3freq;         % one second of sound (Pa corresp. to 0 dB)
x30=sqrt(2)*p30*sin3freq;       % one second of sound (Pa corresp. to 30 dB)
x60=sqrt(2)*p60*sin3freq;       % one second of sound (Pa corresp. to 60 dB)
x90=sqrt(2)*p90*sin3freq;       % one second of sound (Pa corresp. to 90 dB)
% x = sound pressure in Pa, corresponding to 0,30,60,90,60,30,0 dB SPL.
x=[x0,x30,x60,x90,x60,x30,x0];  % seven seconds of sound

Compute the spectrogram.

xNorm=x/p0;  % normalize x by the 0 dB SPL reference pressure
% Compute spectrogram
[ps, freq, time] = pspectrum(xNorm, fs, 'TimeResolution', 0.1, 'spectrogram');
% Convert spectrogram power to dB and make the minimum dB=-20
ps_dB=max(pow2db(ps),-20);

I make the minimum dB=-20, because if I don't, there will be spectrogram power levels of -300 dB or something like that, which is un-helpful when plotting.

Display the spectrogram:

% Display spectrogram
surf(time,freq,ps_dB,EdgeColor="none");
view(0,90); title('Spectrogram (dB_{SPL})')
xlabel('Time (s)'); ylabel('Frequency (Hz)'); zlabel('Power (dB)'); 
cb=colorbar; clim([-20,100]);
cb.Label.String='Power (dB_{SPL})';

●

The plot has the expected colors (corresponding to 0,30,60,90,60,30,0 dB SPL) at the expected frequencies and times. This confirms that we are displaying the intensity in dB SPL, as desired. This spectrogram does not include A-weighting.

[William Rose]

@지호,

My method corresponds to choice number 2, in your list of 4 choices. I use the built-in function pow2dB(), which computes ydB=10*log10(y).