NawfalMotii79
338 lines
12 KiB
Matlab
338 lines
12 KiB
Matlab
% wg_alumina_slotted_gain.m
|
||
% Alumina-filled slotted waveguide (openEMS) with speed/final toggle
|
||
% + S11 + Far-field + Directivity + (Realized) Gain.
|
||
% Axes: X=b (thickness), Y=a (broad), Z=length (propagation).
|
||
|
||
clear; clc;
|
||
|
||
%% ---------------- Toggles ----------------
|
||
speedMode = true; % set false for full/final run
|
||
forcePatterns = true; % set true to compute patterns even in FAST mode
|
||
|
||
%% ---------------- User Parameters (mm) ----------------
|
||
f0 = 10.5e9; % Hz
|
||
fbw = 0.20; % S11 band ±10%
|
||
eps_r = 9.8; % alumina
|
||
tanD = 1e-4; % alumina loss tangent
|
||
|
||
a = 8.5; % Y (mm)
|
||
b = 3.5; % X (mm)
|
||
twall = 0.6; % wall thickness (mm)
|
||
|
||
% slot rules
|
||
slotWidth_rule = 0.02; % w ~ 0.02*lambda0
|
||
slotLen_rule = 0.47; % L ~ 0.47*lambda_g
|
||
x_off_baseFrac = 0.25; % base offset ~ 0.25*a
|
||
x_off_minFrac = 0.22; % clamp min
|
||
x_off_maxFrac = 0.30; % clamp max
|
||
kaiser_beta = 1.65; % ~ -25 dB-like taper
|
||
|
||
% coax (SMA-ish) defaults
|
||
r_in_default = 0.30; % mm
|
||
r_out_default = 1.50; % mm
|
||
|
||
% speed/final presets
|
||
if speedMode
|
||
Nslots = 9;
|
||
r_in = 0.50; % thicker pin -> coarser mesh
|
||
r_out = 1.80;
|
||
slot_W_override = 0.8; % widen slot to 0.8 mm
|
||
endCrit = 5e-4;
|
||
ratio = 2.2;
|
||
doNF2FF = false;
|
||
else
|
||
Nslots = 32;
|
||
r_in = r_in_default;
|
||
r_out = r_out_default;
|
||
slot_W_override = [];
|
||
endCrit = 1e-4;
|
||
ratio = 1.7;
|
||
doNF2FF = true;
|
||
end
|
||
% allow forcing patterns on in fast mode
|
||
if forcePatterns
|
||
doNF2FF = true;
|
||
end
|
||
|
||
r_os = r_out + 0.40; % outer shell radius (mm)
|
||
|
||
%% ---------------- Derived guide quantities ----------------
|
||
lambda0 = (3e8/f0)*1e3; % mm
|
||
fc10 = (3e8)/(2*(a*1e-3)*sqrt(eps_r));
|
||
assert(f0>fc10, 'f0 below TE10 cutoff with alumina fill.');
|
||
|
||
lambda_g = lambda0 / sqrt(1 - (fc10/f0)^2);
|
||
dz = lambda_g/2;
|
||
slot_L = slotLen_rule * lambda_g;
|
||
slot_W = slotWidth_rule * lambda0;
|
||
if ~isempty(slot_W_override), slot_W = slot_W_override; end
|
||
x_off_base = x_off_baseFrac * a;
|
||
|
||
% longitudinal layout (Z)
|
||
z_first = -((Nslots-1)/2)*dz;
|
||
z_feed = z_first - lambda_g/4;
|
||
z_short = z_feed - lambda_g/4;
|
||
|
||
% margins along Z
|
||
if speedMode
|
||
front_margin = dz;
|
||
back_margin = dz;
|
||
else
|
||
front_margin = 2*dz;
|
||
back_margin = 2*dz;
|
||
end
|
||
|
||
z_min = z_short - back_margin;
|
||
z_max = -z_first + front_margin;
|
||
Lz = (z_max - z_min);
|
||
|
||
fprintf('lambda0=%.2f mm, lambda_g=%.2f mm, dz=%.2f mm\n',lambda0,lambda_g,dz);
|
||
fprintf('z_first=%.2f mm, z_feed=%.2f mm, z_short=%.2f mm\n',z_first,z_feed,z_short);
|
||
|
||
%% ---------------- openEMS / CSX ----------------
|
||
unit = 1e-3; % mm
|
||
FDTD = InitFDTD('EndCriteria', endCrit);
|
||
FDTD = SetGaussExcite(FDTD, f0, (f0*fbw)/2);
|
||
FDTD = SetBoundaryCond(FDTD, {'PML_8','PML_8','PML_8','PML_8','PML_8','PML_8'});
|
||
|
||
CSX = InitCSX();
|
||
|
||
% materials
|
||
CSX = AddMetal(CSX,'PEC');
|
||
|
||
CSX = AddMaterial(CSX,'ALUMINA');
|
||
CSX = SetMaterialProperty(CSX,'ALUMINA','Epsilon',eps_r,'Mue',1, ...
|
||
'Kappa', 2*pi*f0*8.854187817e-12*eps_r*tanD);
|
||
|
||
CSX = AddMaterial(CSX,'AIR');
|
||
CSX = SetMaterialProperty(CSX,'AIR','Epsilon',1,'Mue',1);
|
||
|
||
CSX = AddMaterial(CSX,'PTFE');
|
||
CSX = SetMaterialProperty(CSX,'PTFE','Epsilon',2.1);
|
||
|
||
%% ---------------- Geometry (X=b, Y=a, Z=length) ----------------
|
||
x0 = -b/2; x1 = b/2;
|
||
y0 = -a/2; y1 = a/2;
|
||
z0 = z_min; z1 = z_max;
|
||
|
||
% alumina fill
|
||
CSX = AddBox(CSX,'ALUMINA',1,[x0 y0 z0],[x1 y1 z1]);
|
||
|
||
% PEC walls + backshort
|
||
CSX = AddBox(CSX,'PEC',10,[x1 y0 z0],[x1+twall y1 z1]); % top
|
||
CSX = AddBox(CSX,'PEC',10,[x0-twall y0 z0],[x0 y1 z1]); % bottom
|
||
CSX = AddBox(CSX,'PEC',10,[x0 y1 z0],[x1 y1+twall z1]); % +Y
|
||
CSX = AddBox(CSX,'PEC',10,[x0 y0-twall z0],[x1 y0 z1]); % -Y
|
||
CSX = AddBox(CSX,'PEC',12,[x0 y0 z_short],[x1 y1 z_short+twall]); % backshort
|
||
|
||
% Kaiser taper -> slot offsets (alternate signs for dz=lambda_g/2)
|
||
w = kaiser(Nslots, kaiser_beta)'; w = w/max(w);
|
||
x_off_mag = x_off_base * sqrt(w);
|
||
x_off_mag = min(max(x_off_mag, x_off_minFrac*a), x_off_maxFrac*a);
|
||
signs = (-1).^(0:Nslots-1);
|
||
y_centers = signs .* x_off_mag;
|
||
|
||
% Drill each slot through top wall (AIR with higher priority)
|
||
for n = 1:Nslots
|
||
zc = z_first + (n-1)*dz;
|
||
y_c = y_centers(n);
|
||
y1s = y_c - slot_W/2;
|
||
y2s = y_c + slot_W/2;
|
||
z1s = zc - slot_L/2;
|
||
z2s = zc + slot_L/2;
|
||
CSX = AddBox(CSX,'AIR',20,[x1 y1s z1s],[x1+twall y2s z2s]);
|
||
end
|
||
|
||
%% ---------------- Mesh (STRUCT; BEFORE AddCoaxialPort) ----------------
|
||
if speedMode
|
||
maxCell = min([slot_W, r_in, lambda0/15]);
|
||
else
|
||
maxCell = min([slot_W, r_in, lambda0/20]);
|
||
end
|
||
|
||
% Y lines (include slot edges)
|
||
dy_lines = unique(sort([y0-2, y0, y1, y1+2, y_centers-slot_W/2, y_centers+slot_W/2]));
|
||
dy = SmoothMeshLines(dy_lines, maxCell, ratio);
|
||
|
||
% Z lines (include slot edges + feed/short)
|
||
slot_z_edges = zeros(1, 2*Nslots);
|
||
for n = 1:Nslots
|
||
zc = z_first + (n-1)*dz;
|
||
slot_z_edges(2*n-1) = zc - slot_L/2;
|
||
slot_z_edges(2*n) = zc + slot_L/2;
|
||
end
|
||
dz_lines = unique(sort([z0-5, z0, slot_z_edges, z_feed, z_short, z1, z1+5]));
|
||
dz = SmoothMeshLines(dz_lines, maxCell, ratio);
|
||
|
||
% X lines (thickness)
|
||
dx = SmoothMeshLines([x0-2, x0, x1, x1+2], maxCell, ratio);
|
||
|
||
mesh = struct('x', dx, 'y', dy, 'z', dz);
|
||
CSX = DefineRectGrid(CSX, unit, mesh);
|
||
|
||
%% ---------------- Coaxial probe feed (AFTER mesh) ----------------
|
||
% Coax along +Y from -Y wall at z = z_feed
|
||
start = [0, y0 - twall - 0.5, z_feed]; % outside PEC by 0.5 mm
|
||
stop = [0, y0 + 0.05, z_feed]; % 50 µm inside alumina
|
||
|
||
% 11-arg API: AddCoaxialPort(CSX,prio,portnr,'PEC','PTFE',start,stop,dir,r_i,r_o,r_os,...)
|
||
[CSX, port] = AddCoaxialPort( ...
|
||
CSX, 100, 1, 'PEC', 'PTFE', start, stop, 'y', r_in, r_out, r_os, ...
|
||
'ExciteAmp', 1.0);
|
||
|
||
%% ---------------- NF2FF box (robust across versions) ----------------
|
||
doNF2FF = doNF2FF || forcePatterns; % keep your toggles
|
||
|
||
if doNF2FF
|
||
nf_box_name = 'nf2ff';
|
||
nf_start = [x0 - (twall+1), y0 - (twall+1), z0 - 5];
|
||
nf_stop = [x1 + (twall+1), y1 + (twall+1), z1 + 5];
|
||
|
||
% Try modern creator that RETURNS a struct:
|
||
nf2ff = [];
|
||
try
|
||
nf2ff = CreateNF2FFBox(CSX, nf_box_name, nf_start, nf_stop);
|
||
catch
|
||
% Fallback: older helper that RETURNS nothing; synthesize the struct
|
||
try
|
||
AddNF2FFBox(CSX, nf_box_name, nf_start, nf_stop);
|
||
nf2ff = struct('name', nf_box_name);
|
||
catch
|
||
% Some very old versions add via AddBox/SetNF2FFBox—last resort:
|
||
error('NF2FF box creation failed. Your openEMS build lacks Create/AddNF2FFBox.');
|
||
end
|
||
end
|
||
end
|
||
|
||
|
||
|
||
%% ---------------- Write / View / Run ----------------
|
||
if speedMode, Sim_Path = 'sim_wg_alumina_fast';
|
||
else, Sim_Path = 'sim_wg_alumina_final'; end
|
||
Sim_CSX = 'wg_alumina.xml';
|
||
[~,~] = rmdir(Sim_Path, 's'); mkdir(Sim_Path);
|
||
WriteOpenEMS(fullfile(Sim_Path, Sim_CSX), FDTD, CSX);
|
||
|
||
if exist('AppCSXCAD','file'), AppCSXCAD(fullfile(Sim_Path, Sim_CSX)); end
|
||
RunOpenEMS(Sim_Path, Sim_CSX, '--engine=multithreaded');
|
||
|
||
%% ---------------- S11 (robust block) ----------------
|
||
f1 = f0*(1 - fbw/2);
|
||
f2 = f0*(1 + fbw/2);
|
||
Npts = 401;
|
||
fvec = linspace(f1, f2, Npts);
|
||
|
||
try
|
||
port = calcPort(port, Sim_Path, fvec, 'RefImpedance', 50);
|
||
catch
|
||
port = calcPort(port, Sim_Path, fvec);
|
||
end
|
||
|
||
if isfield(port, 'uf') && isfield(port.uf,'ref') && isfield(port.uf,'inc')
|
||
S11 = port.uf.ref ./ port.uf.inc;
|
||
elseif isfield(port, 'U') && isfield(port.U,'ref') && isfield(port.U,'inc')
|
||
S11 = port.U.ref ./ port.U.inc;
|
||
else
|
||
error('Unknown port struct fields. Inspect "port".');
|
||
end
|
||
|
||
figure; plot(fvec/1e9, 20*log10(abs(S11)), 'LineWidth',1.5);
|
||
grid on; xlabel('Frequency (GHz)'); ylabel('|S_{11}| (dB)');
|
||
if speedMode, title(sprintf('|S_{11}| alumina %d-slot (FAST)', Nslots));
|
||
else, title(sprintf('|S_{11}| alumina %d-slot (FINAL)', Nslots)); end
|
||
|
||
% S11 at f0 for mismatch efficiency
|
||
[~, idx0] = min(abs(fvec - f0));
|
||
Gamma0 = S11(idx0);
|
||
eta_m = 1 - abs(Gamma0)^2;
|
||
|
||
%% ---------------- Far-field, Directivity, Gain ----------------
|
||
if doNF2FF
|
||
theta = linspace(0,180,361); % deg
|
||
phi = linspace(0,360,721); % deg
|
||
|
||
% Ensure we pass a struct with .name; if a string slipped through, wrap it
|
||
if ischar(nf2ff)
|
||
nf2ff = struct('name', nf2ff);
|
||
end
|
||
|
||
try
|
||
ff3d = CalcNF2FF(nf2ff, Sim_Path, f0, theta, phi, 'OutDir','nf2ff');
|
||
catch
|
||
% Some builds don’t support 'OutDir' option
|
||
ff3d = CalcNF2FF(nf2ff, Sim_Path, f0, theta, phi);
|
||
end
|
||
|
||
|
||
% Radiation intensity U = |E|^2/(2*eta0) at r=1 m
|
||
eta0 = 376.730313668;
|
||
U = (E.^2) ./ (2*eta0); % W/sr
|
||
|
||
% Integrate over sphere to get Prad
|
||
th = deg2rad(theta(:)); % Nx1
|
||
ph = deg2rad(phi(:)).'; % 1xM
|
||
dth = [diff(th); th(end)-th(end-1)];
|
||
dph = [diff(ph), ph(end)-ph(end-1)];
|
||
[DTH, DPH] = ndgrid(dth, dph);
|
||
[SINE, ~] = ndgrid(sin(th), ones(size(ph)));
|
||
dOmega = SINE .* DTH .* DPH;
|
||
Prad = sum(U .* dOmega, 'all');
|
||
|
||
% Directivity
|
||
D = 4*pi * U / Prad;
|
||
[Dmax, imax] = max(D(:));
|
||
[i_th, i_ph] = ind2sub(size(D), imax);
|
||
th_max = theta(i_th); ph_max = phi(i_ph);
|
||
|
||
% Accepted power Pacc from port struct (best-effort across versions)
|
||
Pacc = NaN;
|
||
if isfield(port,'P_inc') && isfield(port,'P_ref')
|
||
Pinc0 = port.P_inc(idx0); Pref0 = port.P_ref(idx0);
|
||
Pacc = real(Pinc0 - Pref0);
|
||
elseif isfield(port,'uf') && isfield(port,'if') ...
|
||
&& isfield(port.uf,'inc') && isfield(port.if,'inc') ...
|
||
&& isfield(port.uf,'ref') && isfield(port.if,'ref')
|
||
Uinc = port.uf.inc(idx0); Iinc = port.if.inc(idx0);
|
||
Uref = port.uf.ref(idx0); Iref = port.if.ref(idx0);
|
||
Pacc = real(0.5*(Uinc*conj(Iinc) - Uref*conj(Iref)));
|
||
elseif isfield(port,'U') && isfield(port,'I') ...
|
||
&& isfield(port.U,'inc') && isfield(port.I,'inc') ...
|
||
&& isfield(port.U,'ref') && isfield(port.I,'ref')
|
||
Uinc = port.U.inc(idx0); Iinc = port.I.inc(idx0);
|
||
Uref = port.U.ref(idx0); Iref = port.I.ref(idx0);
|
||
Pacc = real(0.5*(Uinc*conj(Iinc) - Uref*conj(Iref)));
|
||
end
|
||
if ~isfinite(Pacc) || Pacc<=0
|
||
warning('Pacc not found in port struct; assuming Pacc = 1 W.');
|
||
Pacc = 1.0;
|
||
end
|
||
eta_rad = Prad / Pacc; eta_rad = max(min(eta_rad,1),0);
|
||
|
||
% Gains
|
||
Gmax = eta_rad * Dmax;
|
||
Greal_max = eta_m * Gmax;
|
||
Gmax_dBi = 10*log10(Gmax);
|
||
Greal_max_dBi = 10*log10(Greal_max);
|
||
|
||
fprintf('\n=== Gain summary @ %.3f GHz ===\n', f0/1e9);
|
||
fprintf('Prad = %.3f W, Pacc = %.3f W => eta_rad = %.1f %%\n', Prad, Pacc, 100*eta_rad);
|
||
fprintf('Mismatch efficiency eta_m = %.1f %% (|S11|=%.2f dB)\n', 100*eta_m, 20*log10(abs(Gamma0)));
|
||
fprintf('Dmax = %.2f (%.2f dBi) at (theta=%.1f°, phi=%.1f°)\n', Dmax, 10*log10(Dmax), th_max, ph_max);
|
||
fprintf('Gmax = %.2f dBi\n', Gmax_dBi);
|
||
fprintf('Greal_max = %.2f dBi (includes mismatch)\n\n', Greal_max_dBi);
|
||
|
||
% Quick E/H-plane plots (normalized)
|
||
[~,i_phi0] = min(abs(phi-0));
|
||
EnormE = E(:,i_phi0) / max(E(:,i_phi0));
|
||
figure; plot(theta, 20*log10(abs(EnormE)+1e-15),'LineWidth',1.5);
|
||
grid on; xlabel('\theta (deg)'); ylabel('Normalized |E| (dB)'); ylim([-60 0]);
|
||
title('E-plane cut (\phi \approx 0^\circ)');
|
||
|
||
[~,i_phiH] = min(abs(phi-90));
|
||
EnormH = E(:,i_phiH) / max(E(:,i_phiH));
|
||
figure; plot(theta, 20*log10(abs(EnormH)+1e-15),'LineWidth',1.5);
|
||
grid on; xlabel('\theta (deg)'); ylabel('Normalized |E| (dB)'); ylim([-60 0]);
|
||
title('H-plane cut (\phi \approx 90^\circ)');
|
||
end
|
||
|