Merge branch 'develop'

downsample_analysis.m

@@ -0,0 +1,114 @@
+clear all
+close all 
+clc
+
+% Modify the desired frequency for the whole code here
+desired_freq = 4000;
+% Read and Set up coeficient
+[y, fs] = audioread("sound/modulator22.wav");
+decrease_coef = round(fs/desired_freq);
+
+% Modify Signal with downsample
+% Only keeps one sample out of decrease_coef
+down_y = downsample(y,decrease_coef);
+down_fs = fs/decrease_coef;
+audiowrite("sound/down_output.wav", down_y, down_fs);
+% Plot Down Modified 
+plot(0:1/down_fs:(length(down_y)-1)/down_fs,down_y);
+xlabel("Time (s)");
+ylabel("Amplitude");
+title("Downsample Sound Amplitude Over Time");
+
+% Modify Signal with decimate
+% Only keeps one sample out of decrease_coef
+deci_y = decimate(y,decrease_coef);
+deci_fs = fs/decrease_coef;
+audiowrite("sound/deci_output.wav", deci_y, deci_fs);
+% Plot Down Modified 
+figure;
+plot(0:1/deci_fs:(length(deci_y)-1)/deci_fs,deci_y);
+xlabel("Time (s)");
+ylabel("Amplitude");
+title("Decimate Sound Amplitude Over Time");
+
+
+% Filter Signal with FIR filter
+order = 30;
+cutoff_freq = 1000; % in Hz
+fir_b = fir1(order, cutoff_freq/(fs/2));
+% Show Filter
+figure
+freqz(fir_b)
+% Apply filter
+fir_y = filter(fir_b, 1, y);
+audiowrite("sound/fir_output.wav", fir_y, fs);
+% Plot Down Modified 
+figure;
+plot(0:1/fs:(length(fir_y)-1)/fs,fir_y);
+xlabel("Time (s)");
+ylabel("Amplitude");
+title("FIR Filter Sound Amplitude Over Time");
+
+% Filter Signal with Butter filter
+order = 8; 
+cutoff_freq = 1000; % in Hz
+[b, a] = butter(order, cutoff_freq/(fs/2), 'low');
+% Show Filter
+figure
+freqz(b)
+% Apply filter
+butt_y = filter(b, a, y);
+audiowrite("sound/butt_output.wav", butt_y, fs);
+% Plot Down Modified 
+figure;
+plot(0:1/fs:(length(butt_y)-1)/fs,butt_y);
+xlabel("Time (s)");
+ylabel("Amplitude");
+title("Butter Filter Sound Amplitude Over Time");
+
+% Test FIR Filter Stability
+% Find the transfer function
+H = tf(fir_b, 1);
+% Get the poles of the transfer function
+poles = pole(H);
+% Check if all poles are inside the unit circle
+if all(abs(poles) < 1)
+    disp('The FIR filter is stable');
+else
+    disp('The FIR filter is unstable');
+end
+% Test Butter Filter Stability
+% Find the transfer function
+H = tf(b, a);
+% Get the poles of the transfer function
+poles = pole(H);
+% Check if all poles are inside the unit circle
+if all(abs(poles) < 1)
+    disp('The Butter filter is stable');
+else
+    disp('The Butter filter is unstable');
+end
+
+% Modify Signal with downsample
+% Only keeps one sample out of decrease_coef
+down_fir_y = downsample(fir_y,decrease_coef);
+down_fir_fs = fs/decrease_coef;
+audiowrite("sound/down_fir_output.wav", down_fir_y, down_fir_fs);
+% Plot Down Modified 
+figure
+plot(0:1/down_fir_fs:(length(down_fir_y)-1)/down_fir_fs,down_fir_y);
+xlabel("Time (s)");
+ylabel("Amplitude");
+title("Downsample FIR Sound Amplitude Over Time");
+
+% Modify Signal with downsample
+% Only keeps one sample out of decrease_coef
+down_butt_y = downsample(butt_y,decrease_coef);
+down_butt_fs = fs/decrease_coef;
+audiowrite("sound/down_butt_output.wav", down_butt_y, down_butt_fs);
+% Plot Down Modified 
+figure
+plot(0:1/down_butt_fs:(length(down_butt_y)-1)/down_butt_fs,down_butt_y);
+xlabel("Time (s)");
+ylabel("Amplitude");
+title("Downsample Butter Sound Amplitude Over Time");

spectral_analysis.m

@@ -3,32 +3,37 @@ close all
 clc
 
 [y, fs] = audioread("sound/modulator22.wav");
-ranges = [17000, 20000; 29000, 37000; 41000, 46000];
+ranges = [17000, 21000; 30500, 36000; 41500, 46000]; 
 one = y(ranges(1,1):ranges(1,2));
 two = y(ranges(2,1):ranges(2,2));
 three = y(ranges(3,1):ranges(3,2));
 
-word = one;
+word = three;
 
 n = length(word);   
 f = (0:n-1)*(fs/n);
 f1 = 0;%Hz
-f2 = 4000;%Hz
+f2 = 2500;%Hz
 idx = find(f >= f1 & f <= f2); %define the index of the freq range
 f = f(idx);
 y = fft(word, n);% compute DFT of input signal
 power = abs(y).^2/n;
-power = power(idx);
+power = power(idx);%limit the signal to the frequency ROI 
 [val, ind] = max(power);
 
-%lowpass for the formant
-
-Fc = 2000; % define the cutoff frequency of the low-pass filter
-[b, a] = butter(6, Fc/(fs/2), 'low'); % design a 4th-order Butterworth low-pass filter
-Pxx_filt = filter(b, a, power); % apply the filter to the power spectrum
-length(Pxx_filt)
-length(f)
+%lowpass for the formant, moving average
 
+for j = 1:length(idx)
+  tot = 0;
+  for k = j-18:j
+    if k <=0
+      tot = tot + power(1);
+    else
+      tot = tot + power(k);
+    endif
+  endfor
+  power_avg(j) = 1/6*(tot);
+endfor
 
 figure;
 
@@ -39,7 +44,7 @@ ylabel('Amplitude (a.u.)');
 
 subplot(1,2,2) % freq range plot 
 plot(f,10*log10(power/power(ind))); hold on;
-plot(f, 10*log10(Pxx_filt), 'r');
+plot(f, 10*log10(power_avg/power(ind)), 'r');
 xlabel('Frequency (Hz)')
 ylabel('Power (dB)')