demo_3d_transfer_function.m

function demo_3d_transfer_function(options)
% DEMO_3D_TRANSFER_FUNCTION 3D visualization of transfer function gain and phase
%   Two separate viewports with synchronized camera rotation
%   Each axis rotates around its own fixed origin
%
%   demo_3d_transfer_function()                  - Run the demo
%   demo_3d_transfer_function('SaveVideo', true) - Save as MP4

    arguments
        options.SaveVideo (1,1) logical = false
        options.OutputFile string = "3d_bode_visualization.mp4"
        options.FrameRate (1,1) double = 30
    end

    % Transfer function parameters (from demo_bode_interactive)
    sigma_pole = 0.25;
    omega_pole = 1.0;
    wn = sqrt(sigma_pole^2 + omega_pole^2);
    zeta = sigma_pole / wn;
    K = wn^2;

    % Data ranges for s-plane
    sigma_min = -2; sigma_max = 0.5;
    omega_min = -3; omega_max = 3;

    % Compute surfaces
    resolution = 60;
    [Sigma, Omega] = meshgrid(...
        linspace(sigma_min, sigma_max, resolution), ...
        linspace(omega_min, omega_max, resolution));

    S = Sigma + 1i * Omega;
    G = K ./ (S.^2 + 2*zeta*wn*S + wn^2);

    % Magnitude in dB and phase
    Mag_dB = 20 * log10(abs(G));
    Phase = angle(G);

    % Clip magnitude
    minMag_dB = -40; maxMag_dB = 40;
    Mag_dB = max(min(Mag_dB, maxMag_dB), minMag_dB);

    % Bode data (σ=0)
    omega_bode = linspace(omega_min, omega_max, 200);
    s_bode = 1i * omega_bode;
    G_bode = K ./ (s_bode.^2 + 2*zeta*wn*s_bode + wn^2);
    Mag_bode_dB = max(min(20*log10(abs(G_bode)), maxMag_dB), minMag_dB);
    Phase_bode = angle(G_bode);

    % Colors (3b1b style)
    BG_COLOR = [0.1, 0.1, 0.12];
    SURFACE_GAIN = [0.2, 0.5, 0.9];
    SURFACE_PHASE = [0.9, 0.4, 0.6];
    BODE_GAIN = [0, 1, 1];
    BODE_PHASE = [1, 1, 0];
    AXIS_COLOR = [0.5, 0.5, 0.6];
    TEXT_COLOR = [0.9, 0.9, 0.9];

    % Create figure with two subplots
    fig = figure('Color', BG_COLOR, 'Position', [50, 100, 1400, 600], ...
        'Name', '3D Transfer Function - Bode View');

    %% LEFT: GAIN plot
    ax1 = subplot(1, 2, 1);
    hold(ax1, 'on');
    set(ax1, 'Color', BG_COLOR, 'XColor', AXIS_COLOR, 'YColor', AXIS_COLOR, 'ZColor', AXIS_COLOR);

    % Surface
    surf(ax1, Sigma, Omega, Mag_dB, 'FaceColor', SURFACE_GAIN, 'FaceAlpha', 0.7, 'EdgeColor', 'none');

    % Bode curve at σ=0
    plot3(ax1, zeros(size(omega_bode)), omega_bode, Mag_bode_dB, 'Color', BODE_GAIN, 'LineWidth', 3);

    % Axes setup
    xlabel(ax1, '\sigma', 'Color', TEXT_COLOR, 'FontSize', 14);
    ylabel(ax1, '\omega', 'Color', TEXT_COLOR, 'FontSize', 14);
    zlabel(ax1, 'dB', 'Color', TEXT_COLOR, 'FontSize', 14);
    title(ax1, 'GAIN', 'Color', TEXT_COLOR, 'FontSize', 18, 'FontWeight', 'bold');

    xlim(ax1, [sigma_min, sigma_max]);
    ylim(ax1, [omega_min, omega_max]);
    zlim(ax1, [minMag_dB, maxMag_dB]);

    grid(ax1, 'on');
    ax1.GridColor = [0.3, 0.3, 0.4];
    ax1.GridAlpha = 0.5;
    axis(ax1, 'vis3d');
    view(ax1, -45, 30);

    %% RIGHT: PHASE plot
    ax2 = subplot(1, 2, 2);
    hold(ax2, 'on');
    set(ax2, 'Color', BG_COLOR, 'XColor', AXIS_COLOR, 'YColor', AXIS_COLOR, 'ZColor', AXIS_COLOR);

    % Surface
    surf(ax2, Sigma, Omega, Phase, 'FaceColor', SURFACE_PHASE, 'FaceAlpha', 0.7, 'EdgeColor', 'none');

    % Bode curve at σ=0
    plot3(ax2, zeros(size(omega_bode)), omega_bode, Phase_bode, 'Color', BODE_PHASE, 'LineWidth', 3);

    % Axes setup
    xlabel(ax2, '\sigma', 'Color', TEXT_COLOR, 'FontSize', 14);
    ylabel(ax2, '\omega', 'Color', TEXT_COLOR, 'FontSize', 14);
    zlabel(ax2, 'rad', 'Color', TEXT_COLOR, 'FontSize', 14);
    title(ax2, 'PHASE', 'Color', TEXT_COLOR, 'FontSize', 18, 'FontWeight', 'bold');

    xlim(ax2, [sigma_min, sigma_max]);
    ylim(ax2, [omega_min, omega_max]);
    zlim(ax2, [-pi, pi]);

    grid(ax2, 'on');
    ax2.GridColor = [0.3, 0.3, 0.4];
    ax2.GridAlpha = 0.5;
    axis(ax2, 'vis3d');
    view(ax2, -45, 30);

    drawnow;

    %% Setup video recording
    if options.SaveVideo
        fprintf('Recording video to %s...\n', options.OutputFile);
        v = VideoWriter(options.OutputFile, 'MPEG-4');
        v.FrameRate = options.FrameRate;
        v.Quality = 95;
        open(v);
    end

    fps = options.FrameRate;

    %% Hold initial 3D view for 1 second
    holdFrames = fps * 1;
    for i = 1:holdFrames
        if options.SaveVideo
            writeVideo(v, getframe(fig));
        else
            pause(1/fps);
        end
    end

    %% Animate synchronized rotation to 2D Bode view
    fprintf('Rotating to 2D Bode view...\n');

    startAz = -45; startEl = 30;
    endAz = 90; endEl = 0;
    duration = 3; % seconds
    numFrames = duration * fps;

    for i = 1:numFrames
        t = i / numFrames;
        % Smooth interpolation (ease in-out)
        t_smooth = (1 - cos(t * pi)) / 2;

        az = startAz + (endAz - startAz) * t_smooth;
        el = startEl + (endEl - startEl) * t_smooth;

        view(ax1, az, el);
        view(ax2, az, el);

        drawnow;

        if options.SaveVideo
            writeVideo(v, getframe(fig));
        else
            pause(1/fps);
        end
    end

    %% Hold final 2D Bode view for 2 seconds
    view(ax1, endAz, endEl);
    view(ax2, endAz, endEl);
    drawnow;

    holdFrames = fps * 2;
    for i = 1:holdFrames
        if options.SaveVideo
            writeVideo(v, getframe(fig));
        else
            pause(1/fps);
        end
    end

    %% Finish up
    if options.SaveVideo
        close(v);
        fprintf('Video saved to: %s\n', options.OutputFile);
    end

    fprintf('G(s) = %.2f / (s² + %.2fs + %.2f)\n', K, 2*zeta*wn, wn^2);
    fprintf('Poles: s = %.2f ± j%.2f\n', -sigma_pole, omega_pole);
end