Hello guys,
I'm seeking to calculate the polarimetric RCS for different 3D shapes. To get
started, I aimed to calculate the none-polarimetric RCS of a simple sphere. The
following code that I have implemented should be able to do so. However, the
result is not what I expected. After struggling for days, I've still not been
able to solve the issue. So I would appreciate it if you have a look at the
code and the corresponding result, so you may be able to help.
=============================================
from __future__ import division
import meep as mp
import math
import numpy
import matplotlib.pyplot as plt
import scipy.constants
# ------------------------------------------------------- #
# calculating the (non)-polarimetric RCS of a sphere
# ------------------------------------------------------- #
# Setting up the resolution
resolution = 5
# Define unit (mm)
unit = 1e-3
# ------------------------------------------------------- #
# Defining the parameters. All frequencies in Hz and lengths in mm
# ------------------------------------------------------- #
f_center = 75e9
f_bw = f_center * 2
f_steps = 151 # number of frequency samples in the considered bandwidth
# Define environment size
box_length = 10e-3
box_size = numpy.array([box_length, box_length, box_length])
# Define N2F Measurement box
near2far_box_length = 8e-3
near2far_box_size = numpy.array([near2far_box_length, near2far_box_length,
near2far_box_length])
# Define Perfectly Matched Layers (PML)
pml_thickness = 2e-3
# Sphere radius
radius = 1e-3
# Farfield position
farfield_pos = mp.Vector3(0, 0, -1)
# Environment with PML considered
total_box_size = box_size + 2 * pml_thickness
# ------------------------------------------------------- #
# converting all parameters to MEEP unit [f_meep = unit/λ] & [l_meep = l/unit]]
# ------------------------------------------------------- #
f_c_meep_unit = (unit / scipy.constants.c) * f_center
f_bw_meep_unit = (unit / scipy.constants.c) * f_bw
box_size_meep_unit = box_size/unit
near2far_box_size_meep_unit = near2far_box_size/unit
pml_thickness_meep_unit = pml_thickness/unit
total_box_size_meep_unit = total_box_size/unit
radius_meep_unit = radius/unit
farfield_pos_meep_unit = farfield_pos/unit
# ------------------------------------------------------- #
# Setting up the source's size and position
wave_source_center = mp.Vector3(z=(box_size[2] + near2far_box_size[2]) / -4.0 /
unit)
wave_source_size = mp.Vector3(x=0.9 * box_size_meep_unit[0], y=0.9 *
box_size_meep_unit[1])
# Build PML
pml_layers = [mp.PML(pml_thickness_meep_unit)]
# Specify source field component
src_cmpt = mp.Ex
# Define source wave in time domain
gaussian_source = mp.GaussianSource(f_c_meep_unit, fwidth=f_bw_meep_unit)
plane_wave_source = mp.Source(src=gaussian_source,
component=mp.Ex,
center=wave_source_center,
size=wave_source_size,
amp_func=lambda pos: 1.0) # plane wave
sources = [plane_wave_source]
geometry = [mp.Sphere(radius=radius_meep_unit,
center=mp.Vector3(),
material=mp.metal)]
# Initialising the simulation without the shape for calculating the scattered
field later on (calibration)
sim =
mp.Simulation(cell_size=(mp.Vector3(x=total_box_size_meep_unit[0],y=total_box_size_meep_unit[1],z=total_box_size_meep_unit[2])),
resolution=resolution,
sources=sources,
boundary_layers=pml_layers)
# Build a N2F box
positive_boundary_position = numpy.zeros(3) + near2far_box_size_meep_unit / 2
negative_boundary_position = numpy.zeros(3) - near2far_box_size_meep_unit / 2
size = near2far_box_size_meep_unit
nearfield_box = sim.add_near2far(f_c_meep_unit, f_bw_meep_unit, f_steps,
mp.Near2FarRegion(mp.Vector3(x=positive_boundary_position[0]),
size=mp.Vector3(y=size[1], z=size[2])),
mp.Near2FarRegion(mp.Vector3(x=negative_boundary_position[0]),
size=mp.Vector3(y=size[1], z=size[2]), weight=-1),
mp.Near2FarRegion(mp.Vector3(y=positive_boundary_position[1]),
size=mp.Vector3(x=size[0], z=size[2])),
mp.Near2FarRegion(mp.Vector3(y=negative_boundary_position[1]),
size=mp.Vector3(x=size[0], z=size[2]), weight=-1),
mp.Near2FarRegion(mp.Vector3(z=positive_boundary_position[2]),
size=mp.Vector3(x=size[0], y=size[1])),
mp.Near2FarRegion(mp.Vector3(z=negative_boundary_position[2]),
size=mp.Vector3(x=size[0], y=size[1]), weight=-1))
# Position where field data should be measured and evaluated as a stop condition
decayed_field_eval_position = mp.Vector3(z=((box_size[2] +
near2far_box_size[2]) / 4.0) / unit)
# 2D animating the calibrating simulation
v = mp.Volume(center=mp.Vector3(),
size=mp.Vector3(x=total_box_size_meep_unit[0], z=total_box_size_meep_unit[2]))
anim = mp.Animate2D(sim,
fields=mp.Ex,
normalize=True,
output_plane=v,
realtime=True)
# Running the simulation without the sphere
sim.run(mp.at_every(1, anim),
until_after_sources=mp.stop_when_fields_decayed(10, src_cmpt,
decayed_field_eval_position, 1e-9))
# Determine fields around the N2F box in frequency domain which are induced by
the source only
no_target_reference = sim.get_near2far_data(nearfield_box)
# From now let's consider the shape
sim.reset_meep()
sim = mp.Simulation(cell_size=(mp.Vector3(x=total_box_size_meep_unit[0],
y=total_box_size_meep_unit[1], z=total_box_size_meep_unit[2])),
resolution=resolution,
sources=sources,
boundary_layers=pml_layers,
geometry=geometry)
# Build a N2F box
positive_boundary_position = numpy.zeros(3) + near2far_box_size_meep_unit / 2
negative_boundary_position = numpy.zeros(3) - near2far_box_size_meep_unit / 2
size = near2far_box_size_meep_unit
nearfield_box = sim.add_near2far(f_c_meep_unit, f_bw_meep_unit, f_steps,
mp.Near2FarRegion(mp.Vector3(x=positive_boundary_position[0]),
size=mp.Vector3(y=size[1], z=size[2])),
mp.Near2FarRegion(mp.Vector3(x=negative_boundary_position[0]),
size=mp.Vector3(y=size[1], z=size[2]), weight=-1),
mp.Near2FarRegion(mp.Vector3(y=positive_boundary_position[1]),
size=mp.Vector3(x=size[0], z=size[2])),
mp.Near2FarRegion(mp.Vector3(y=negative_boundary_position[1]),
size=mp.Vector3(x=size[0], z=size[2]), weight=-1),
mp.Near2FarRegion(mp.Vector3(z=positive_boundary_position[2]),
size=mp.Vector3(x=size[0], y=size[1])),
mp.Near2FarRegion(mp.Vector3(z=negative_boundary_position[2]),
size=mp.Vector3(x=size[0], y=size[1]), weight=-1))
# Calculating the scattered field from the sphere
sim.load_minus_near2far_data(nearfield_box, no_target_reference)
# Frequency range of observed waves in f_steps samples (actually equal to that
for the source)
freqs_meep_unit = mp.get_near2far_freqs(nearfield_box)
# 2D animating the scattering simulation
anim = mp.Animate2D(sim,
fields=mp.Ex,
normalize=True,
output_plane=v,
realtime=True)
# Running the simulation with the sphere
sim.run(mp.at_every(1, anim),
until_after_sources=mp.stop_when_fields_decayed(10, src_cmpt,
decayed_field_eval_position, 1e-9))
# Calculate the farfield
target_farfield = numpy.array(sim.get_farfield(nearfield_box,
farfield_pos_meep_unit))
# Changing format from [Ex(f0),Ey(f0),Ez(f0),Hx(f0),Hy(f0),Hz(f0),Ex(f1)....]
to a 6 column matrix
target_farfield = target_farfield.reshape((-1, 6))
# ------------------------------------------------------- #
# Calculating and plotting the RCS
# ------------------------------------------------------- #
source_amplitude = [gaussian_source.fourier_transform(f) for f in
freqs_meep_unit]
# converting back the units
frequencies = numpy.array(freqs_meep_unit) * scipy.constants.c / unit
plt.figure()
plt.plot(2 * math.pi * radius * frequencies / scipy.constants.c,
4 * numpy.abs(target_farfield[:, 0])**2 /
(numpy.abs(source_amplitude)**2 * radius**2))
plt.grid(True)
plt.xlabel('Circumference / Wave length')
plt.ylabel('RCS / (pi * r**2)')
plt.show()
=============================================
Please note that my goal is to expand this code to a polarimetric simulation
later on. Therefore, I intentionally didn't use flux for RCS calculation. For
building an appropriate N2F structure, the sphere should be irradiated with a
plane wave. Interesting enough, the PMLs seem to lead to none-perfect plane
waves, as it is demonstrated by the animations. So maybe this could be the
problem.
Kind regards, Khashayar
_______________________________________________
meep-discuss mailing list
meep-discuss@ab-initio.mit.edu
http://ab-initio.mit.edu/cgi-bin/mailman/listinfo/meep-discuss
