# White noise
ms = modulationSpectrum(rnorm(16000), samplingRate = 16000,
logSpec = FALSE, power = TRUE,
amRes = NULL) # analyze the entire sound, giving a single roughness value
str(ms)
# Harmonic sound
s = soundgen(pitch = 440, amFreq = 100, amDep = 50)
ms = modulationSpectrum(s, samplingRate = 16000, amRes = NULL)
ms[c('roughness', 'amMsFreq', 'amMsPurity')] # a single value for each
ms1 = modulationSpectrum(s, samplingRate = 16000, amRes = 5)
ms1[c('roughness', 'amMsFreq', 'amMsPurity')]
# measured over time (low values of amRes mean more precision, so we analyze
# longer segments and get fewer values per sound)
# Embellish
ms = modulationSpectrum(s, samplingRate = 16000, logMPS = TRUE,
xlab = 'Temporal modulation, Hz', ylab = 'Spectral modulation, 1/kHz',
colorTheme = 'matlab', main = 'Modulation spectrum', lty = 3)
# 1D instead of 2D
modulationSpectrum(s, 16000, msType = '1D')
if (FALSE) {
# A long sound with varying AM and a bit of chaos at the end
s_long = soundgen(sylLen = 3500, pitch = c(250, 320, 280),
amFreq = c(30, 55), amDep = c(20, 60, 40),
jitterDep = c(0, 0, 2))
playme(s_long)
ms = modulationSpectrum(s_long, 16000)
# plot AM over time
plot(x = seq(1, 1500, length.out = length(ms$amMsFreq)), y = ms$amMsFreq,
cex = 10^(ms$amMsPurity/20) * 10, xlab = 'Time, ms', ylab = 'AM frequency, Hz')
# plot roughness over time
spectrogram(s_long, 16000, ylim = c(0, 4),
extraContour = list(ms$roughness / max(ms$roughness) * 4000, col = 'blue'))
# As with spectrograms, there is a tradeoff in time-frequency resolution
s = soundgen(pitch = 500, amFreq = 50, amDep = 100, sylLen = 500,
samplingRate = 44100, plot = TRUE)
# playme(s, samplingRate = 44100)
ms = modulationSpectrum(s, samplingRate = 44100,
windowLength = 50, step = 50, amRes = NULL) # poor temporal resolution
ms = modulationSpectrum(s, samplingRate = 44100,
windowLength = 5, step = 1, amRes = NULL) # poor frequency resolution
ms = modulationSpectrum(s, samplingRate = 44100,
windowLength = 15, step = 3, amRes = NULL) # a reasonable compromise
# Start with an auditory spectrogram instead of STFT
modulationSpectrum(s, 44100, specSource = 'audSpec', xlim = c(-100, 100))
modulationSpectrum(s, 44100, specSource = 'audSpec',
logWarpX = c(10, 2), xlim = c(-500, 500),
audSpec_pars = list(nFilters = 32, filterType = 'gammatone', bandwidth = NULL))
# customize the plot
ms = modulationSpectrum(s, samplingRate = 44100,
windowLength = 15, overlap = 80, amRes = NULL,
kernelSize = 17, # more smoothing
xlim = c(-70, 70), ylim = c(0, 4), # zoom in on the central region
quantiles = c(.25, .5, .8), # customize contour lines
col = rev(rainbow(100)), # alternative palette
logWarpX = c(10, 2), # pseudo-log transform
power = 2) # ^2
# Note the peaks at FM = 2/kHz (from "pitch = 500") and AM = 50 Hz (from
# "amFreq = 50")
# Input can be a wav/mp3 file
ms = modulationSpectrum('~/Downloads/temp/16002_Faking_It_Large_clear.wav')
# Input can be path to folder with audio files. Each file is processed
# separately, and the output can contain an MS per file...
ms1 = modulationSpectrum('~/Downloads/temp', kernelSize = 11,
plot = FALSE, averageMS = FALSE)
ms1$summary
names(ms1$original) # a separate MS per file
# ...or a single MS can be calculated:
ms2 = modulationSpectrum('~/Downloads/temp', kernelSize = 11,
plot = FALSE, averageMS = TRUE)
plotMS(ms2$original)
ms2$summary
# Input can also be a list of waveforms (numeric vectors)
ss = vector('list', 10)
for (i in 1:length(ss)) {
ss[[i]] = soundgen(sylLen = runif(1, 100, 1000), temperature = .4,
pitch = runif(3, 400, 600))
}
# lapply(ss, playme)
# MS of the first sound
ms1 = modulationSpectrum(ss[[1]], samplingRate = 16000, scale = 1)
# average MS of all 10 sounds
ms2 = modulationSpectrum(ss, samplingRate = 16000, scale = 1, averageMS = TRUE, plot = FALSE)
plotMS(ms2$original)
# A sound with ~3 syllables per second and only downsweeps in F0 contour
s = soundgen(nSyl = 8, sylLen = 200, pauseLen = 100, pitch = c(300, 200))
# playme(s)
ms = modulationSpectrum(s, samplingRate = 16000, maxDur = .5,
xlim = c(-25, 25), colorTheme = 'seewave',
power = 2)
# note the asymmetry b/c of downsweeps
# "power = 2" returns squared modulation spectrum - note that this affects
# the roughness measure!
ms$roughness
# compare:
modulationSpectrum(s, samplingRate = 16000, maxDur = .5,
xlim = c(-25, 25), colorTheme = 'seewave',
power = 1)$roughness # much higher roughness
# Plotting with or without log-warping the modulation spectrum:
ms = modulationSpectrum(soundgen(), samplingRate = 16000, plot = TRUE)
ms = modulationSpectrum(soundgen(), samplingRate = 16000,
logWarpX = c(2, 2), plot = TRUE)
# logWarp and kernelSize have no effect on roughness
# because it is calculated before these transforms:
modulationSpectrum(s, samplingRate = 16000, logWarpX = c(1, 10))$roughness
modulationSpectrum(s, samplingRate = 16000, logWarpX = NA)$roughness
modulationSpectrum(s, samplingRate = 16000, kernelSize = 17)$roughness
# Log-transform the spectrogram prior to 2D FFT (affects roughness):
modulationSpectrum(s, samplingRate = 16000, logSpec = FALSE)$roughness
modulationSpectrum(s, samplingRate = 16000, logSpec = TRUE)$roughness
# Use a lognormal weighting function to calculate roughness
# (instead of just % in roughRange)
modulationSpectrum(s, 16000, roughRange = NULL,
roughMean = 75, roughSD = 3)$roughness
modulationSpectrum(s, 16000, roughRange = NULL,
roughMean = 100, roughSD = 12)$roughness
# truncate weights outside roughRange
modulationSpectrum(s, 16000, roughRange = c(30, 150),
roughMean = 100, roughSD = 1000)$roughness # very large SD
modulationSpectrum(s, 16000, roughRange = c(30, 150),
roughMean = NULL)$roughness # same as above b/c SD --> Inf
# Complex modulation spectrum with phase preserved
ms = modulationSpectrum(soundgen(), samplingRate = 16000,
returnComplex = TRUE)
plotMS(abs(ms$complex)) # note the symmetry
# compare:
plotMS(ms$original)
}
Run the code above in your browser using DataLab