% ==============================================
% Wave Impedance of Hertzian Dipole
%
% Description:
%   Plots the magnitude and phase of wave impedance Zw
%   as a function of normalized distance r (in wavelengths)
%
% Formula:
%   Zw = 377 * ( (2*pi*r)^2/(1 + (2*pi*r)^2) - j/(2*pi*r*(1 + (2*pi*r)^2)) )
% ==============================================

clear;
clc;
close all;

% --- USER INPUT ---
error_percent = input('Enter the percentage error for 377 Ohms (e.g., 1 for 1%): ');

% --- NORMALIZED DISTANCE RANGE ---
r = linspace(0.1, 1, 1000);  % r from 0.1 to 1 (linear scale)

% --- CALCULATE WAVE IMPEDANCE ---
x = 2 * pi * r;  % Normalized parameter

% Real and imaginary parts
real_part = (x.^2) ./ (1 + x.^2);
imag_part = -1 ./ (x .* (1 + x.^2));

% Complex wave impedance
Zw = 377 * (real_part + 1i * imag_part);

% Magnitude and phase
Zw_mag = abs(Zw);
Zw_phase = angle(Zw) * 180 / pi;  % Convert to degrees

% --- FIND DISTANCE WHERE |Zw| APPROACHES 377 WITHIN ERROR ---
target_impedance = 377;
error_margin = error_percent / 100;
lower_bound = target_impedance * (1 - error_margin);
upper_bound = target_impedance * (1 + error_margin);

% Near-field boundary
r_boundary = 1/(2*pi);

% Find indices where |Zw| is within error bounds AND r >= r_boundary
valid_range = find(r >= r_boundary);
within_error = find(Zw_mag(valid_range) >= lower_bound & Zw_mag(valid_range) <= upper_bound);

if !isempty(within_error)
    r_critical = r(valid_range(within_error(1)));  % First point where error condition is met (r >= boundary)
    fprintf('|Zw| reaches 377 ± %.1f%% at r = %.4f wavelengths\n', ...
            error_percent, r_critical);
    %fprintf('This corresponds to %.4f meters at 1 GHz\n', r_critical * 0.3);
else
    fprintf('|Zw| does not reach 377 ± %.1f%% in the range r = [%.4f, 1]\n', ...
            error_percent, r_boundary);
    r_critical = NaN;
end

% --- PLOT SETUP ---
figure;
set(gcf, 'Position', [100, 100, 900, 800]);

% --- SUBPLOT 1: MAGNITUDE ---
subplot(2, 1, 1);
plot(r, Zw_mag, 'b-', 'LineWidth', 3);
hold on;
grid on;

% Near-field boundary
line([r_boundary r_boundary], ylim, 'Color', 'k', 'LineStyle', '--', 'LineWidth', 2, ...
     'DisplayName', 'Near/Far Field Boundary');

% Error threshold lines
if !isempty(within_error)
    plot([min(r) max(r)], [lower_bound lower_bound], 'g:', 'LineWidth', 2, ...
         'DisplayName', sprintf('%.1f%% Error Bounds', error_percent));
    plot([min(r) max(r)], [upper_bound upper_bound], 'g:', 'LineWidth', 2, ...
         'HandleVisibility', 'off');
    plot([r_critical r_critical], [200 400], 'm--', 'LineWidth', 2, ...
         'DisplayName', sprintf('Critical Distance (r=%.3fλ)', r_critical));
end

% Labels and formatting with larger fonts
xlabel('Normalized Distance r (wavelengths)', 'FontSize', 16, 'FontWeight', 'bold');
ylabel('|Z_w| (Ohms)', 'FontSize', 16, 'FontWeight', 'bold');
title('Magnitude of Wave Impedance', 'FontSize', 18, 'FontWeight', 'bold');
ylim([200 400]);
xlim([0.1 1]);
set(gca, 'FontSize', 14);
legend('FontSize', 14, 'Location', 'southeast');

% Annotations with larger fonts
text(r_boundary + 0.02, 220, sprintf('r = 1/(2π) ≈ %.3fλ', r_boundary), ...
     'FontSize', 14, 'BackgroundColor', 'white', 'FontWeight', 'bold');
%text(0.15, 220, 'Near Field', 'FontSize', 14, 'Color', 'blue', 'FontWeight', 'bold');
%text(0.8, 390, 'Far Field', 'FontSize', 14, 'Color', 'blue', 'FontWeight', 'bold');

if !isempty(within_error)
    text(r_critical + 0.02, 350, sprintf('r = %.3fλ', r_critical), ...
         'FontSize', 14, 'Color', 'm', 'BackgroundColor', 'white', 'FontWeight', 'bold');
end

% --- SUBPLOT 2: PHASE ---
subplot(2, 1, 2);
plot(r, Zw_phase, 'r-', 'LineWidth', 3);
hold on;
grid on;

% Near-field boundary
line([r_boundary r_boundary], ylim, 'Color', 'k', 'LineStyle', '--', 'LineWidth', 2, ...
     'DisplayName', 'Near/Far Field Boundary');

% Labels and formatting with larger fonts
xlabel('Normalized Distance r (wavelengths)', 'FontSize', 16, 'FontWeight', 'bold');
ylabel('arg(Z_w) (degrees)', 'FontSize', 16, 'FontWeight', 'bold');
title('Phase of Wave Impedance', 'FontSize', 18, 'FontWeight', 'bold');
ylim([-90 10]);
xlim([0.1 1]);
set(gca, 'FontSize', 14);
legend('FontSize', 14, 'Location', 'southeast');

% Annotations with larger fonts
text(r_boundary + 0.02, -80, sprintf('r = 1/(2π) ≈ %.3fλ', r_boundary), ...
     'FontSize', 14, 'BackgroundColor', 'white', 'FontWeight', 'bold');
%text(0.15, -85, 'Phase ≈ -90°', 'FontSize', 14, 'Color', 'red', 'FontWeight', 'bold');
%text(0.8, -5, 'Phase ≈ 0°', 'FontSize', 14, 'Color', 'red', 'FontWeight', 'bold');

% --- DISPLAY RESULTS ---
fprintf('\nWave Impedance Characteristics:\n');
fprintf('Near-field boundary: r = 1/(2π) ≈ %.4f wavelengths\n', r_boundary);
fprintf('Far-field impedance: 377 Ohms (free space impedance)\n');
fprintf('Near-field impedance: capacitive (phase ≈ -90°)\n');

% Show specific values at boundary
boundary_idx = find(r >= r_boundary, 1);
fprintf('\nAt near-field boundary (r = %.4fλ):\n', r_boundary);
fprintf('|Zw| = %.1f Ohms\n', Zw_mag(boundary_idx));
fprintf('arg(Zw) = %.1f degrees\n', Zw_phase(boundary_idx));

% Show values at r = 0.5λ
r_half_idx = find(r >= 0.5, 1);
fprintf('\nAt r = 0.5λ:\n');
fprintf('|Zw| = %.1f Ohms\n', Zw_mag(r_half_idx));
fprintf('arg(Zw) = %.1f degrees\n', Zw_phase(r_half_idx));

% Show far-field values at r = 1λ
far_field_idx = find(r >= 1, 1);
fprintf('\nAt r = 1λ:\n');
fprintf('|Zw| = %.1f Ohms\n', Zw_mag(far_field_idx));
fprintf('arg(Zw) = %.1f degrees\n', Zw_phase(far_field_idx));

% Display error analysis results
fprintf('\nError Analysis (%.1f%% tolerance):\n', error_percent);
if !isempty(within_error)
    fprintf('|Zw| first reaches 377 ± %.1f Ohms at r = %.4fλ (r ≥ 1/(2π))\n', ...
            error_percent*377/100, r_critical);
    fprintf('At this distance:\n');
    fprintf('  |Zw| = %.1f Ohms\n', Zw_mag(valid_range(within_error(1))));
    fprintf('  arg(Zw) = %.1f degrees\n', Zw_phase(valid_range(within_error(1))));

    % Additional info: how close to 377 Ohms
    error_actual = abs(Zw_mag(valid_range(within_error(1)))) - target_impedance;
    fprintf('  Actual error: %.2f Ohms (%.2f%%)\n', abs(error_actual), abs(error_actual)/target_impedance*100);
else
    fprintf('|Zw| never reaches 377 ± %.1f Ohms in the range r = [%.4f, 1]λ\n', ...
            error_percent*377/100, r_boundary);

    % Show the closest value in the far-field region
    [min_error, min_idx] = min(abs(Zw_mag(valid_range) - target_impedance));
    r_closest = r(valid_range(min_idx));
    fprintf('Closest approach in far field: r = %.4fλ, |Zw| = %.1f Ohms (error: %.1f%%)\n', ...
            r_closest, Zw_mag(valid_range(min_idx)), min_error/target_impedance*100);
end