% ============================
% SCRIPT 2: BODE ANALYSIS
% ============================

clear all; close all;

% ----- USER SETTINGS -----
I_amp = 0.5;   % same amplitude used in Script 1
freqs = [0.1 0.2 0.5 1 2 5];

gain_list  = zeros(size(freqs));
phase_list = zeros(size(freqs));

for k = 1:length(freqs)

    f = freqs(k);
    filename = sprintf("voltage%.3.csv", f);  % recommended naming
    filename = strrep(filename, '.', 'p');        % 0.1 → 0p1

    fprintf("\nReading %s...\n", filename);

    % Load CSV
    data = readmatrix(filename);
    t_log = data(:,1);
    V_log = data(:,2);

    % Remove DC
    V_ac = V_log - mean(V_log);

    % Fit sine
    w = 2*pi*f;
    sine_model = @(b,t) b(1)*sin(w*t + b(2));
    beta0 = [0.01, 0];

    beta = nlinfit(t_log, V_ac, sine_model, beta0);

    V_amp = abs(beta(1));
    gain_list(k)  = V_amp / I_amp;
    phase_list(k) = rad2deg(beta(2));

    fprintf("Gain = %.6f V/A\n", gain_list(k));
    fprintf("Phase = %.2f deg\n", phase_list(k));
end

% ----- BODE PLOTS -----
figure;
subplot(2,1,1);
semilogx(freqs, gain_list, '-o');
xlabel('Frequency (Hz)');
ylabel('Gain (V/A)');
title('Bode Gain');
grid on;

subplot(2,1,2);
semilogx(freqs, phase_list, '-o');
xlabel('Frequency (Hz)');
ylabel('Phase Lag (deg)');
title('Bode Phase');
grid on;