Skip to content

MATLAB functions for .m4a & .wav #1

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
jghezzi27 wants to merge 1 commit into
base: main
Choose a base branch
from
Open

MATLAB functions for .m4a & .wav #1

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
95 changes: 95 additions & 0 deletions BPMCalculate.m
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -0,0 +1,95 @@
function [BPM1, BPM2, BPM3] = BPMCalculate(data, fs, sSize, threshold, setting)
%BPMCalculate reads in a song in a .wav file and calculates its BPM
% Detailed explanation goes here

%% Process Data
if setting == 1
ampValues = ampData(data,fs,sSize);

else
for i = 1:(length(data) / sSize)
% Much of the antics here are due to 1 indexing...
sample = data([(((i - 1) * sSize) + 1):(i * sSize)], 1);
sample = abs(sample);
ampValues(:,i) = sample;
end
end

sample2Time = sSize / fs;

%% Find Peaks
sumVec = sum(ampValues, 1);

sumChangeVec = zeros(length(sumVec),1);
for i = 2:length(sumVec)
sumChangeVec(i) = sumVec(i) - sumVec(i - 1);
end

peaks = zeros(length(sumChangeVec));
for i = 1:length(sumChangeVec)
if sumChangeVec(i) > threshold
peaks(i) = 1;
end
end

peakIndices = find(peaks);
if isempty(peakIndices) || length(peakIndices) == 1 || length(peakIndices) == 2
BPM1 = nan;
BPM2 = nan;
BPM3 = nan;
return;
end

extraDistance = 0;
nextIndex = 1;
for i = 2:length(peakIndices)
% Special case where distance between is 1
if peakIndices(i) - peakIndices(i-1) == 1
extraDistance = extraDistance + 1;
continue
end

% Special case where extra distance needs to be added on
if extraDistance ~= 0
% Add on extra distance from previous 1s
peakIndices(i) = peakIndices(i) + extraDistance;
peakSampleDist(nextIndex) = peakIndices(i) - peakIndices(i-1);

% Remove that extra distance so it doesn't affect future points
peakIndices(i) = peakIndices(i) - extraDistance;
extraDistance = 0;
nextIndex = nextIndex + 1;
continue
end

peakSampleDist(nextIndex) = peakIndices(i) - peakIndices(i-1);
nextIndex = nextIndex + 1;
end

%% Convert peakSampleDist to peakBPMs
peakTimeDist = peakSampleDist .* sample2Time;

peakBPMs = zeros(length(peakTimeDist),1);
for i = 1:length(peakTimeDist)
peakBPMs(i) = (peakTimeDist(i) ^ (-1)) * 60;
end

% Convert BPMs to BPMs between 75 and 150
i = 1;
while i < length(peakBPMs) + 1
if peakBPMs(i) > 150
peakBPMs(i) = peakBPMs(i) / 2;

elseif peakBPMs(i) < 75
peakBPMs(i) = peakBPMs(i) * 2;

else
i = i + 1;
end
end

BPM1 = mean(peakBPMs);
BPM2 = median(peakBPMs);
BPM3 = mode(peakBPMs);

end
14 changes: 14 additions & 0 deletions ampData.m
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -0,0 +1,14 @@
function [ampData] = ampData(data, fs, sSize)
%ampData takes in .wav or .m4a file data and outputs amplitude data for samples
%with size sSize
% Detailed explanation goes here

%% Make a fft for each time sample and store in an array
for i = 1:(length(data) / sSize)
% Much of the antics here are due to 1 indexing...
sample = data([(((i - 1) * sSize) + 1):(i * sSize)], 1);
[~, singleSided] = fftProcess(sample, length(sample), fs);
ampData(:,i) = singleSided;
end

end
93 changes: 93 additions & 0 deletions bpmMain.m
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -0,0 +1,93 @@
clear; close all;
% Optimal Setting: Threshold: 0.0920 Resolution: 960

% %% Resolution Test 1
% [data, fs] = audioread("102ByAndBy.m4a");
% actBPM = 102;
%
% resolutions = [1080, 1200, 1320, 1440, 1560, 1680, 1800, 2040, 2160, 2280, 2400];
% for j = 1:length(resolutions)
%
% Setting 1
% ampValues = ampData(data,fs,sSize);
%
% Setting 2
% ampValues = zeros(resolutions(j), floor(length(data) / resolutions(j)));
% for i = 1:(length(data) / resolutions(j))
% % Much of the antics here are due to 1 indexing...
% sample = data([(((i - 1) * resolutions(j)) + 1):(i * resolutions(j))], 1);
% sample = abs(sample);
% ampValues(:,i) = sample;
% end
%
% % Data Processing
% sumVec = sum(ampValues, 1);
%
% sumChangeVec = zeros(length(sumVec),1);
% for i = 2:length(sumVec)
% sumChangeVec(i) = sumVec(i) - sumVec(i - 1);
% end
%
% sample2Time = resolutions(j) / fs;
% sampleTimeVec = [1:size(ampValues, 2)] .* sample2Time;
%
% % This is solely for testing and not part of the actual algorithm
% sampleBPMVec = sampleTimeVec .* actBPM ./ 60;
%
% figure();
% plot(sampleBPMVec, sumChangeVec);
% axis([2.5, 12, 0, 20]);
% end

%% Resolution Test 2
resolutions = [1080, 1200, 1320, 1440, 1560, 1680, 1800, 2040, 2160, 2280, 2400];
songList = ["75LeanOnMe", "77Power", "78NoWoman", "82Hopeless", "86Blinding", ...
"98TheMan", "102ByAndBy", "104InMyLife", "107Snow", "114Rainbow", ...
"116IceIceBaby", "120ProveIt", "121Talk", "128September", "137LetItBe", ...
"140RadioNW", "142AllLights"];
numSongs = length(songList);
for i = 1:numSongs
songList(i) = songList(i) + ".m4a";
end

% Threshold Vector
starts = [3, 2, 2, 2, 2, 2, 4, 4, 4, 4, 4];
ends = [8, 8, 8, 12, 10, 10, 12, 15, 12, 12, 12];
numThreshValues = 11;
threshold = zeros(numThreshValues, length(resolutions));
for j = 1:numThreshValues
for k = 1:length(resolutions)
threshold(j, k) = starts(k) + ((ends(k) - starts(k)) * (j - 1) / (numThreshValues - 1));
end
end

meanExpBPM = zeros(numSongs,numThreshValues,length(resolutions));
medExpBPM = meanExpBPM;
modeExpBPM = medExpBPM;

% Actual BPM Vector
actBPM = zeros(numSongs, numThreshValues, length(resolutions));
actBPMplane = zeros(numSongs,numThreshValues);
for j = 1:numThreshValues
actBPMplane(:,j) = [75, 77, 78, 82, 86, 98, 102, 104, 107, 114, 116, 120, 121, 128, 137, 140, 142];
end

for j = 1:length(resolutions)
actBPM(:,:,j) = actBPMplane;
end

for i = 1:numSongs
[data, fs] = audioread(songList(i));
%data = data(1:2000000);
for j = 1:numThreshValues
for k = 1:length(resolutions)
[meanExpBPM(i,j,k), medExpBPM(i,j,k), modeExpBPM(i,j,k)] = ...
BPMCalculate(data, fs, resolutions(k), threshold(j, k), 1);
end
end
end

%% Data Analysis
[meanError, medError, modeError, meanBigError, medBigError, modeBigError, ...
settingTest, averageError, difficulty] = ...
dataAnalysis(meanExpBPM, medExpBPM, modeExpBPM, actBPM);
171 changes: 171 additions & 0 deletions dataAnalysis.m
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -0,0 +1,171 @@
function [meanError, medError, modeError, meanBigError, medBigError, modeBigError, settingTest, averageError, difficulty] = dataAnalysis(meanExpBPM, medExpBPM, modeExpBPM, actBPM)
%dataAnalysis Takes experimental and theoretical BPM data for various
%settings and returns analysis on data
% Detailed explanation goes here
%% Calculate Error
% Row: song
% Column: Threshold
% Plane: Resolution

meanDifferences = meanExpBPM - actBPM;
medDifferences = medExpBPM - actBPM;
modeDifferences = modeExpBPM - actBPM;

meanError = 100 .* abs(meanDifferences) ./ actBPM;
medError = 100 .* abs(medDifferences) ./ actBPM;
modeError = 100 .* abs(modeDifferences) ./ actBPM;

% Adjust Error that is for values twice too big
for i = 1:size(meanError, 1)
for j = 1:size(meanError, 2)
for k = 1:size(meanError, 3)
if meanError(i,j,k) > 80
meanExpBPM(i,j,k) = 0.5 * meanExpBPM(i,j,k);
meanDifferences(i,j,k) = meanExpBPM(i,j,k) - actBPM(i,j,k);
meanError(i,j,k) = 100 .* abs(meanDifferences(i,j,k)) ./ actBPM(i,j,k);
end

if medError(i,j,k) > 80
medExpBPM(i,j,k) = 0.5 * medExpBPM(i,j,k);
medDifferences(i,j,k) = medExpBPM(i,j,k) - actBPM(i,j,k);
medError(i,j,k) = 100 .* abs(medDifferences(i,j,k)) ./ actBPM(i,j,k);
end

if modeError(i,j,k) > 80
modeExpBPM(i,j,k) = 0.5 * modeExpBPM(i,j,k);
modeDifferences(i,j,k) = modeExpBPM(i,j,k) - actBPM(i,j,k);
modeError(i,j,k) = 100 .* abs(modeDifferences(i,j,k)) ./ actBPM(i,j,k);
end
end
end
end

%% Find Where the Function was innaccurate or Failed
% Row: song
% Column: Threshold
% Plane: Resolution

meanBigError = zeros(size(meanError, 1), size(meanError, 2), size(meanError, 3));
medBigError = meanBigError;
modeBigError = meanBigError;
bigErrorPercent = 5;

for i = 1:size(meanError, 1)
for j = 1:size(meanError, 2)
for k = 1:size(meanError, 3)
if meanError(i, j, k) > bigErrorPercent || isnan(meanError(i, j, k))
meanBigError(i, j, k) = 1;
end
if medError(i, j, k) > bigErrorPercent || isnan(medError(i, j, k))
medBigError(i, j, k) = 1;
end
if modeError(i, j, k) > bigErrorPercent || isnan(modeError(i, j, k))
modeBigError(i, j, k) = 1;
end
end
end
end

%% Count the NaN Values for Each Resolution & Threshold
% Row: Threshold
% Column: Resolution
% Plane 1: Mean
% Plane 2: Median
% Plane 3: Mode

nanCounts = zeros(size(meanError, 2), size(meanError, 3), 3);
for j = 1:size(meanError, 2)
for k = 1:size(meanError, 3)
nanCount1 = 0;
nanCount2 = 0;
nanCount3 = 0;

for i = 1:size(meanError, 1)
if isnan(meanError(i, j, k))
nanCount1 = nanCount1 + 1;
end
if isnan(medError(i, j, k))
nanCount2 = nanCount2 + 1;
end
if isnan(modeError(i, j, k))
nanCount3 = nanCount3 + 1;
end
end
nanCounts(j, k, 1) = nanCount1;
nanCounts(j, k, 2) = nanCount2;
nanCounts(j, k, 3) = nanCount3;
end
end

%% Find Percent of Time that Song Failed
% Row: Song
% Column 1: Mean Errors
% Column 2: Median Errors
% Column 3: Mode Errors

difficulty = zeros(size(meanError, 1), 3);
for i = 1:size(meanError, 1)
difficulty(i, 1) = 100 * sum(sum(meanBigError(i, :, :))) / (size(meanError, 2) * size(meanError, 3));
difficulty(i, 2) = 100 * sum(sum(medBigError(i, :, :))) / (size(meanError, 2) * size(meanError, 3));
difficulty(i, 3) = 100 * sum(sum(modeBigError(i, :, :))) / (size(meanError, 2) * size(meanError, 3));
end

%% Find Average Error and Test Each Setting
% Row: Threshold
% Column: Resolution
% Plane 1: Mean
% Plane 2: Median
% Plane 3: Mode

settingTest = zeros(size(meanError, 2), size(meanError, 3), 3);
averageError = settingTest;
for i = 1:size(meanError, 2)
for j = 1:size(meanError, 3)
settingTest(i, j, 1) = sum(meanBigError(:, i, j));
settingTest(i, j, 2) = sum(medBigError(:, i, j));
settingTest(i, j, 3) = sum(modeBigError(:, i, j));
averageError(i, j, 1) = sum(meanError(:, i, j), 'omitnan') / (size(meanError, 1) - nanCounts(i, j, 1));
averageError(i, j, 2) = sum(medError(:, i, j), 'omitnan') / (size(medError, 1) - nanCounts(i, j, 2));
averageError(i, j, 3) = sum(modeError(:, i, j), 'omitnan') / (size(modeError, 1) - nanCounts(i, j, 3));
end
end

% %% Find the Optimal Setting for each Song
% % Row: Song
% % Column 1: Threshold
% % Column 2: Resolution
% % Column 3: Percent Error
% % Plane 1: Mean Error
% % Plane 2: Median Error
% % Plane 3: Mode Error
%
% bestSetting = zeros(size(meanError, 1), 3, 3);
% for i = 1:size(meanError, 1)
% bestError = 100;
% for j = 1:size(meanError, 2)
% for k = 1:size(meanError, 3)
% if meanError(i, j, k) < bestError
% bestSetting(i, 1, 1) = j;
% bestSetting(i, 2, 1) = k;
% bestSetting(i, 3, 1) = meanError(i, j, k);
% bestError = meanError(i, j, k);
% end
%
% if medError(i, j, k) < bestError
% bestSetting(i, 1, 2) = j;
% bestSetting(i, 2, 2) = k;
% bestSetting(i, 3, 2) = medError(i, j, k);
% bestError = medError(i, j, k);
% end
%
% if modeError(i, j, k) < bestError
% bestSetting(i, 1, 3) = j;
% bestSetting(i, 2, 3) = k;
% bestSetting(i, 3, 3) = modeError(i, j, k);
% bestError = modeError(i, j, k);
% end
% end
% end
% end
%
end
9 changes: 9 additions & 0 deletions fftProcess.m
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -0,0 +1,9 @@
function [freqVec, singleSided] = fftProcess(data, recordLength, fs)
%fftProcess Takes in sound data and returns fast fourier transformed data
dataFFT = fft(data);
twoSided = abs(dataFFT ./ recordLength); % Two sided spectrum
singleSided = twoSided(1:floor(recordLength/2 + 1)); % Remove negative frequency peaks
singleSided(2:(end - 1)) = 2 .* singleSided(2:(end - 1));

freqVec = fs .* (0:recordLength/2) / recordLength;
end