General Usage
The basic building blocks are introduced in the following simple example; more details are provided afterwards:
Introductory Example: Plane Wave
using SphericalScattering, StaticArrays
# define excitation: plane wave travelling in positive z-direction with x-polarization
ex = planeWave(frequency=10e6) # Hz
# define scatterer: PEC sphere
sp = PECSphere(radius = 1.0)
# define an observation point
point_cart = [SVector(2.0, 2.0, 3.2)]
# compute scattered fields
E = scatteredfield(sp, ex, ElectricField(point_cart))
H = scatteredfield(sp, ex, MagneticField(point_cart))
FF = scatteredfield(sp, ex, FarField(point_cart))
Defining Observation Points
In order to define points the StaticArrays package has to be used.
using StaticArrays
# defining a single point
point_cart = [SVector(2.0, 2.0, 3.2)]
# defining multiple points (along a line)
point_cart = [SVector(5.0, 5.0, z) for z in -2:0.2:2]
A Cartesian basis is used for all coordinates and field components!
Defining an Excitation
For all available excitations a simple constructor with keyword arguments and default values is available. For more details see the APIs of the
Defining a Scatterer
For all available scatteres a simple constructor with keyword arguments and default values is available. For more details see the APIs of the
- PEC sphere
- PMC sphere
- Dielectric sphere
- Multilayer dielectric sphere
- Multilayer dielectric sphere with PEC core
- Dielectric sphere with thin impedance layer
Computing Fields
The incident, scattered, and total fields can be computed.
Incident Fields
For each excitation the far-field, the electric, and the magnetic near-field without a scatterer can be determined:
E = field(ex, ElectricField(point_cart))
H = field(ex, MagneticField(point_cart))
FF = field(ex, FarField(point_cart))
For the uniform field excitation, only the electric field as well as the scalar potential can be calculated:
Φ = field(ex, ScalarPotential(point_cart))
Scattered Fields
For each excitation the scattered fields from a given sphere can be determined as:
E = scatteredfield(sp, ex, ElectricField(point_cart))
H = scatteredfield(sp, ex, MagneticField(point_cart))
FF = scatteredfield(sp, ex, FarField(point_cart))
For the uniform static field excitation, only the electric field as well as the scalar potential can be calculated:
Φ = scatteredfield(sp, ex, ScalarPotential(point_cart))
For the dielectric sphere with thin impedance layer two additional quantities are available:
Φ = scatteredfield(sp, ex, ScalarPotentialJump(point_cart))
E = scatteredfield(sp, ex, DisplacementField(point_cart))
Total Fields
For each excitation the total fields in the presence of a given sphere can be determined as:
E = field(sp, ex, ElectricField(point_cart))
H = field(sp, ex, MagneticField(point_cart))
FF = field(sp, ex, FarField(point_cart))
For the uniform field excitation, only the electric field as well as the scalar potential can be calculated:
Φ = field(sp, ex, ScalarPotential(point_cart))
Radar Cross Section
To compute the bistatic radar cross section (RCS), the function
σ = rcs(sp, ex, points_cart)
is provided. For the monostatic RCS, the function
σ = rcs(sp, ex)
is provided.
The RCS is (so far) only defined for a plane wave excitation.
Conversion Between Bases
Methods are provided to convert between Cartesian and spherical coordinates:
point_cart = SphericalScattering.sph2cart.(point_sph)
point_sph = SphericalScattering.cart2sph.(point_cart)
Converting fields:
F_cart = SphericalScattering.convertSpherical2Cartesian.(F_sph, point_sph)
F_sph = SphericalScattering.convertCartesian2Spherical.(F_cart, point_sph)
Plotting Fields
To visualize far-fields and near-fields the functions
SphericalScattering.sphericalGridPoints
— FunctionpatternPoints(;r=1.0, resolution=5)
Convenience function returning points (in Cartesian and spherical coordinates) on a spherical grid at a distance 'r' and a resolution in degrees.
SphericalScattering.phiCutPoints
— FunctionphiCutPoints(ϕ; r=1.0, resolution=1)
Convenience function returning points (in Cartesian coordinates) on a phi cut at a distance 'r' and a resolution in degrees.
SphericalScattering.thetaCutPoints
— FunctionthetaCutPoints(ϑ; r=1.0, resolution=1)
Convenience function returning points (in Cartesian coordinates) on a theta cut at a distance 'r' and a resolution in degrees.
as well as the functions
plotff(F, points_sph; scale="log", normalize=true, type="abs")
plotffcut(F, points; scale="log", normalize=true, format="polar")
are provided (after loading the PlotlyJS package). For more details see the visualization of fields examples.
Issues have been reported with PlotlyJS if an installation via Conda is employed. See, e.g., this thread for issues with Julia and Conda.