199 lines
6.7 KiB
Python
199 lines
6.7 KiB
Python
|
# -*- coding: utf-8 -*-
|
||
|
|
"""
|
||
|
|
Bent Patch Antenna Tutorial
|
||
|
|
|
||
|
|
Tested with
|
||
|
|
- python 3.4
|
||
|
|
- openEMS v0.0.33+
|
||
|
|
|
||
|
|
(C) 2016 Thorsten Liebig <thorsten.liebig@gmx.de>
|
||
|
|
|
||
|
|
"""
|
||
|
|
|
||
|
|
### Import Libraries
|
||
|
|
import os, tempfile
|
||
|
|
from pylab import *
|
||
|
|
from mpl_toolkits.mplot3d import Axes3D
|
||
|
|
|
||
|
|
from CSXCAD import CSXCAD
|
||
|
|
|
||
|
|
from openEMS.openEMS import openEMS
|
||
|
|
from openEMS.physical_constants import *
|
||
|
|
|
||
|
|
|
||
|
|
### Setup the simulation
|
||
|
|
Sim_Path = os.path.join(tempfile.gettempdir(), 'Bent_Patch')
|
||
|
|
|
||
|
|
post_proc_only = False
|
||
|
|
|
||
|
|
unit = 1e-3 # all length in mm
|
||
|
|
|
||
|
|
f0 = 2.4e9 # center frequency, frequency of interest!
|
||
|
|
lambda0 = round(C0/f0/unit) # wavelength in mm
|
||
|
|
fc = 0.5e9 # 20 dB corner frequency
|
||
|
|
|
||
|
|
# patch width in alpha-direction
|
||
|
|
patch_width = 32 # resonant length in alpha-direction
|
||
|
|
patch_radius = 50 # radius
|
||
|
|
patch_length = 40 # patch length in z-direction
|
||
|
|
|
||
|
|
#substrate setup
|
||
|
|
substrate_epsR = 3.38
|
||
|
|
substrate_kappa = 1e-3 * 2*pi*2.45e9 * EPS0*substrate_epsR
|
||
|
|
substrate_width = 80
|
||
|
|
substrate_length = 90
|
||
|
|
substrate_thickness = 1.524
|
||
|
|
substrate_cells = 4
|
||
|
|
|
||
|
|
#setup feeding
|
||
|
|
feed_pos = -5.5 #feeding position in x-direction
|
||
|
|
feed_width = 2 #feeding port width
|
||
|
|
feed_R = 50 #feed resistance
|
||
|
|
|
||
|
|
# size of the simulation box
|
||
|
|
SimBox_rad = 2*100
|
||
|
|
SimBox_height = 1.5*200
|
||
|
|
|
||
|
|
### Setup FDTD parameter & excitation function
|
||
|
|
FDTD = openEMS(CoordSystem=1) # init a cylindrical FDTD
|
||
|
|
f0 = 2e9 # center frequency
|
||
|
|
fc = 1e9 # 20 dB corner frequency
|
||
|
|
FDTD.SetGaussExcite(f0, fc)
|
||
|
|
FDTD.SetBoundaryCond(['MUR', 'MUR', 'MUR', 'MUR', 'MUR', 'MUR']) # boundary conditions
|
||
|
|
|
||
|
|
### Setup the Geometry & Mesh
|
||
|
|
# init a cylindrical mesh
|
||
|
|
CSX = CSXCAD.ContinuousStructure(CoordSystem=1)
|
||
|
|
FDTD.SetCSX(CSX)
|
||
|
|
mesh = CSX.GetGrid()
|
||
|
|
mesh.SetDeltaUnit(unit)
|
||
|
|
|
||
|
|
### Setup the geometry using cylindrical coordinates
|
||
|
|
# calculate some width as an angle in radiant
|
||
|
|
patch_ang_width = patch_width/(patch_radius+substrate_thickness)
|
||
|
|
substr_ang_width = substrate_width/patch_radius
|
||
|
|
feed_angle = feed_pos/patch_radius
|
||
|
|
|
||
|
|
# create patch
|
||
|
|
patch = CSX.AddMetal('patch') # create a perfect electric conductor (PEC)
|
||
|
|
start = [patch_radius+substrate_thickness, -patch_ang_width/2, -patch_length/2 ]
|
||
|
|
stop = [patch_radius+substrate_thickness, patch_ang_width/2, patch_length/2 ]
|
||
|
|
CSX.AddBox(patch, priority=10, start=start, stop=stop, edges2grid='all') # add a box-primitive to the metal property 'patch'
|
||
|
|
|
||
|
|
# create substrate
|
||
|
|
substrate = CSX.AddMaterial('substrate', epsilon=substrate_epsR, kappa=substrate_kappa )
|
||
|
|
start = [patch_radius , -substr_ang_width/2, -substrate_length/2]
|
||
|
|
stop = [patch_radius+substrate_thickness, substr_ang_width/2, substrate_length/2]
|
||
|
|
CSX.AddBox(substrate, start=start, stop=stop, edges2grid='all')
|
||
|
|
|
||
|
|
# save current density oon the patch
|
||
|
|
jt_patch = CSX.AddDump('Jt_patch', dump_type=3, file_type=1)
|
||
|
|
start = [patch_radius+substrate_thickness, -substr_ang_width/2, -substrate_length/2]
|
||
|
|
stop = [patch_radius+substrate_thickness, +substr_ang_width/2, substrate_length/2]
|
||
|
|
CSX.AddBox(jt_patch, start=start, stop=stop)
|
||
|
|
|
||
|
|
# create ground
|
||
|
|
gnd = CSX.AddMetal('gnd') # create a perfect electric conductor (PEC)
|
||
|
|
start = [patch_radius, -substr_ang_width/2, -substrate_length/2]
|
||
|
|
stop = [patch_radius, +substr_ang_width/2, +substrate_length/2]
|
||
|
|
CSX.AddBox(gnd, priority=10, start=start, stop=stop, edges2grid='all')
|
||
|
|
|
||
|
|
# apply the excitation & resist as a current source
|
||
|
|
start = [patch_radius , feed_angle, 0]
|
||
|
|
stop = [patch_radius+substrate_thickness, feed_angle, 0]
|
||
|
|
port = FDTD.AddLumpedPort(1 ,feed_R, start, stop, 'r', 1.0, priority=50, edges2grid='all')
|
||
|
|
|
||
|
|
### Finalize the Mesh
|
||
|
|
# add the simulation domain size
|
||
|
|
mesh.AddLine('r', patch_radius+np.array([-20, SimBox_rad]))
|
||
|
|
mesh.AddLine('a', [-0.75*pi, 0.75*pi])
|
||
|
|
mesh.AddLine('z', [-SimBox_height/2, SimBox_height/2])
|
||
|
|
|
||
|
|
# add some lines for the substrate
|
||
|
|
mesh.AddLine('r', patch_radius+np.linspace(0,substrate_thickness,substrate_cells))
|
||
|
|
|
||
|
|
# generate a smooth mesh with max. cell size: lambda_min / 20
|
||
|
|
max_res = C0 / (f0+fc) / unit / 20
|
||
|
|
max_ang = max_res/(SimBox_rad+patch_radius) # max res in radiant
|
||
|
|
mesh.SmoothMeshLines(0, max_res, 1.4)
|
||
|
|
mesh.SmoothMeshLines(1, max_ang, 1.4)
|
||
|
|
mesh.SmoothMeshLines(2, max_res, 1.4)
|
||
|
|
|
||
|
|
## Add the nf2ff recording box
|
||
|
|
nf2ff = FDTD.CreateNF2FFBox()
|
||
|
|
|
||
|
|
### Run the simulation
|
||
|
|
if 0: # debugging only
|
||
|
|
CSX_file = os.path.join(Sim_Path, 'bent_patch.xml')
|
||
|
|
if not os.path.exists(Sim_Path):
|
||
|
|
os.mkdir(Sim_Path)
|
||
|
|
CSX.Write2XML(CSX_file)
|
||
|
|
os.system(r'AppCSXCAD "{}"'.format(CSX_file))
|
||
|
|
|
||
|
|
|
||
|
|
if not post_proc_only:
|
||
|
|
FDTD.Run(Sim_Path, verbose=3, cleanup=True)
|
||
|
|
|
||
|
|
### Postprocessing & plotting
|
||
|
|
f = np.linspace(max(1e9,f0-fc),f0+fc,401)
|
||
|
|
port.CalcPort(Sim_Path, f)
|
||
|
|
Zin = port.uf_tot / port.if_tot
|
||
|
|
s11 = port.uf_ref/port.uf_inc
|
||
|
|
s11_dB = 20.0*np.log10(np.abs(s11))
|
||
|
|
|
||
|
|
figure()
|
||
|
|
plot(f/1e9, s11_dB)
|
||
|
|
grid()
|
||
|
|
ylabel('s11 (dB)')
|
||
|
|
xlabel('frequency (GHz)')
|
||
|
|
|
||
|
|
P_in = 0.5*np.real(port.uf_tot * np.conj(port.if_tot)) # antenna feed power
|
||
|
|
|
||
|
|
# plot feed point impedance
|
||
|
|
figure()
|
||
|
|
plot( f/1e6, real(Zin), 'k-', linewidth=2, label=r'$\Re(Z_{in})$' )
|
||
|
|
grid()
|
||
|
|
plot( f/1e6, imag(Zin), 'r--', linewidth=2, label=r'$\Im(Z_{in})$' )
|
||
|
|
title( 'feed point impedance' )
|
||
|
|
xlabel( 'frequency (MHz)' )
|
||
|
|
ylabel( 'impedance ($\Omega$)' )
|
||
|
|
legend( )
|
||
|
|
|
||
|
|
|
||
|
|
idx = np.where((s11_dB<-10) & (s11_dB==np.min(s11_dB)))[0]
|
||
|
|
if not len(idx)==1:
|
||
|
|
print('No resonance frequency found for far-field calulation')
|
||
|
|
else:
|
||
|
|
f_res = f[idx[0]]
|
||
|
|
theta = np.arange(-180.0, 180.0, 2.0)
|
||
|
|
print("Calculate NF2FF")
|
||
|
|
nf2ff_res_phi0 = nf2ff.CalcNF2FF(Sim_Path, f_res, theta, 0, center=np.array([patch_radius+substrate_thickness, 0, 0])*unit, read_cached=True, outfile='nf2ff_xz.h5')
|
||
|
|
|
||
|
|
figure(figsize=(15, 7))
|
||
|
|
ax = subplot(121, polar=True)
|
||
|
|
E_norm = 20.0*np.log10(nf2ff_res_phi0.E_norm/np.max(nf2ff_res_phi0.E_norm)) + nf2ff_res_phi0.Dmax
|
||
|
|
ax.plot(np.deg2rad(theta), 10**(np.squeeze(E_norm)/20), linewidth=2, label='xz-plane')
|
||
|
|
ax.grid(True)
|
||
|
|
ax.set_xlabel('theta (deg)')
|
||
|
|
ax.set_theta_zero_location('N')
|
||
|
|
ax.set_theta_direction(-1)
|
||
|
|
ax.legend(loc=3)
|
||
|
|
|
||
|
|
phi = theta
|
||
|
|
nf2ff_res_theta90 = nf2ff.CalcNF2FF(Sim_Path, f_res, 90, phi, center=np.array([patch_radius+substrate_thickness, 0, 0])*unit, read_cached=True, outfile='nf2ff_xy.h5')
|
||
|
|
|
||
|
|
ax = subplot(122, polar=True)
|
||
|
|
E_norm = 20.0*np.log10(nf2ff_res_theta90.E_norm/np.max(nf2ff_res_theta90.E_norm)) + nf2ff_res_theta90.Dmax
|
||
|
|
ax.plot(np.deg2rad(phi), 10**(np.squeeze(E_norm)/20), linewidth=2, label='xy-plane')
|
||
|
|
ax.grid(True)
|
||
|
|
ax.set_xlabel('phi (deg)')
|
||
|
|
suptitle('Bent Patch Anteanna Pattern\nFrequency: {} GHz'.format(f_res/1e9), fontsize=14)
|
||
|
|
ax.legend(loc=3)
|
||
|
|
|
||
|
|
print( 'radiated power: Prad = {:.2e} Watt'.format(nf2ff_res_theta90.Prad[0]))
|
||
|
|
print( 'directivity: Dmax = {:.1f} ({:.1f} dBi)'.format(nf2ff_res_theta90.Dmax[0], 10*np.log10(nf2ff_res_theta90.Dmax[0])))
|
||
|
|
print( 'efficiency: nu_rad = {:.1f} %'.format(100*nf2ff_res_theta90.Prad[0]/real(P_in[idx[0]])))
|
||
|
|
|
||
|
|
show()
|
||
|
|
|