你已经派生过 gprMax
镜像自地址
https://gitee.com/sunhf/gprMax.git
已同步 2025-08-07 15:10:13 +08:00
Move fractal objects into seperate files
这个提交包含在:
171
gprMax/fractals/fractal_surface.py
普通文件
171
gprMax/fractals/fractal_surface.py
普通文件
@@ -0,0 +1,171 @@
|
||||
import numpy as np
|
||||
import numpy.typing as npt
|
||||
from scipy import fftpack
|
||||
|
||||
from gprMax import config
|
||||
from gprMax.cython.fractals_generate import generate_fractal2D
|
||||
|
||||
np.seterr(divide="raise")
|
||||
|
||||
|
||||
class FractalSurface:
|
||||
"""Fractal surfaces."""
|
||||
|
||||
surfaceIDs = ["xminus", "xplus", "yminus", "yplus", "zminus", "zplus"]
|
||||
|
||||
def __init__(self, xs, xf, ys, yf, zs, zf, dimension, seed):
|
||||
"""
|
||||
Args:
|
||||
xs, xf, ys, yf, zs, zf: floats for the extent of the fractal surface
|
||||
(one pair of coordinates must be equal
|
||||
to correctly define a surface).
|
||||
dimension: float for the fractal dimension that controls the fractal
|
||||
distribution.
|
||||
seed: int for seed value for random number generator.
|
||||
"""
|
||||
|
||||
self.ID = None
|
||||
self.surfaceID = None
|
||||
self.start = np.array([xs, ys, zs], dtype=np.int32)
|
||||
self.stop = np.array([xf, yf, zf], dtype=np.int32)
|
||||
self.dtype = np.dtype(np.complex128)
|
||||
self.seed = seed
|
||||
self.dimension = (
|
||||
dimension # Fractal dimension from: http://dx.doi.org/10.1017/CBO9781139174695
|
||||
)
|
||||
self.weighting = np.array([1, 1], dtype=np.float64)
|
||||
self.fractalrange = (0, 0)
|
||||
self.filldepth = 0
|
||||
self.grass = []
|
||||
|
||||
@property
|
||||
def xs(self) -> int:
|
||||
return self.start[0]
|
||||
|
||||
@xs.setter
|
||||
def xs(self, value: int):
|
||||
self.start[0] = value
|
||||
|
||||
@property
|
||||
def ys(self) -> int:
|
||||
return self.start[1]
|
||||
|
||||
@ys.setter
|
||||
def ys(self, value: int):
|
||||
self.start[1] = value
|
||||
|
||||
@property
|
||||
def zs(self) -> int:
|
||||
return self.start[2]
|
||||
|
||||
@zs.setter
|
||||
def zs(self, value: int):
|
||||
self.start[2] = value
|
||||
|
||||
@property
|
||||
def xf(self) -> int:
|
||||
return self.stop[0]
|
||||
|
||||
@xf.setter
|
||||
def xf(self, value: int):
|
||||
self.stop[0] = value
|
||||
|
||||
@property
|
||||
def yf(self) -> int:
|
||||
return self.stop[1]
|
||||
|
||||
@yf.setter
|
||||
def yf(self, value: int):
|
||||
self.stop[1] = value
|
||||
|
||||
@property
|
||||
def zf(self) -> int:
|
||||
return self.stop[2]
|
||||
|
||||
@zf.setter
|
||||
def zf(self, value: int):
|
||||
self.stop[2] = value
|
||||
|
||||
@property
|
||||
def size(self) -> npt.NDArray[np.int32]:
|
||||
return self.stop - self.start
|
||||
|
||||
@property
|
||||
def nx(self) -> int:
|
||||
return self.xf - self.xs
|
||||
|
||||
@property
|
||||
def ny(self) -> int:
|
||||
return self.yf - self.ys
|
||||
|
||||
@property
|
||||
def nz(self) -> int:
|
||||
return self.zf - self.zs
|
||||
|
||||
def get_surface_dims(self):
|
||||
"""Gets the dimensions of the fractal surface based on surface plane."""
|
||||
|
||||
if self.xs == self.xf:
|
||||
surfacedims = (self.ny, self.nz)
|
||||
elif self.ys == self.yf:
|
||||
surfacedims = (self.nx, self.nz)
|
||||
elif self.zs == self.zf:
|
||||
surfacedims = (self.nx, self.ny)
|
||||
|
||||
return surfacedims
|
||||
|
||||
def generate_fractal_surface(self):
|
||||
"""Generate a 2D array with a fractal distribution."""
|
||||
|
||||
surfacedims = self.get_surface_dims()
|
||||
|
||||
self.fractalsurface = np.zeros(surfacedims, dtype=self.dtype)
|
||||
|
||||
# Positional vector at centre of array, scaled by weighting
|
||||
v1 = np.array(
|
||||
[
|
||||
self.weighting[0] * (surfacedims[0]) / 2,
|
||||
self.weighting[1] * (surfacedims[1]) / 2,
|
||||
]
|
||||
)
|
||||
|
||||
# 2D array of random numbers to be convolved with the fractal function
|
||||
rng = np.random.default_rng(seed=self.seed)
|
||||
A = rng.standard_normal(size=(surfacedims[0], surfacedims[1]))
|
||||
|
||||
# 2D FFT
|
||||
A = fftpack.fftn(A)
|
||||
|
||||
# Generate fractal
|
||||
generate_fractal2D(
|
||||
surfacedims[0],
|
||||
surfacedims[1],
|
||||
0,
|
||||
0,
|
||||
surfacedims[0],
|
||||
surfacedims[1],
|
||||
config.get_model_config().ompthreads,
|
||||
self.dimension,
|
||||
self.weighting,
|
||||
v1,
|
||||
A,
|
||||
self.fractalsurface,
|
||||
)
|
||||
|
||||
# Set DC component of FFT to zero
|
||||
self.fractalsurface[0, 0] = 0
|
||||
# Take the real part (numerical errors can give rise to an imaginary part)
|
||||
# of the IFFT, and convert type to floattype. N.B calculation of fractals
|
||||
# must always be carried out at double precision, i.e. float64, complex128
|
||||
self.fractalsurface = np.real(fftpack.ifftn(self.fractalsurface)).astype(
|
||||
config.sim_config.dtypes["float_or_double"], copy=False
|
||||
)
|
||||
# Scale the fractal volume according to requested range
|
||||
fractalmin = np.amin(self.fractalsurface)
|
||||
fractalmax = np.amax(self.fractalsurface)
|
||||
fractalrange = fractalmax - fractalmin
|
||||
self.fractalsurface = (
|
||||
self.fractalsurface * ((self.fractalrange[1] - self.fractalrange[0]) / fractalrange)
|
||||
+ self.fractalrange[0]
|
||||
- ((self.fractalrange[1] - self.fractalrange[0]) / fractalrange) * fractalmin
|
||||
)
|
@@ -17,7 +17,7 @@
|
||||
# along with gprMax. If not, see <http://www.gnu.org/licenses/>.
|
||||
|
||||
import logging
|
||||
from typing import List, Optional, Tuple, Union
|
||||
from typing import List, Optional, Tuple
|
||||
|
||||
import numpy as np
|
||||
import numpy.typing as npt
|
||||
@@ -26,10 +26,8 @@ from mpi4py_fft import PFFT, DistArray, newDistArray
|
||||
from mpi4py_fft.pencil import Subcomm
|
||||
from scipy import fftpack
|
||||
|
||||
import gprMax.config as config
|
||||
from gprMax.cython.fractals_generate import generate_fractal2D, generate_fractal3D
|
||||
from gprMax.grid.mpi_grid import MPIGrid
|
||||
from gprMax.output_controllers.grid_view import MPIGridView
|
||||
from gprMax import config
|
||||
from gprMax.cython.fractals_generate import generate_fractal3D
|
||||
from gprMax.utilities.mpi import Dim, Dir, get_relative_neighbour
|
||||
from gprMax.utilities.utilities import round_value
|
||||
|
||||
@@ -37,169 +35,6 @@ logger = logging.getLogger(__name__)
|
||||
np.seterr(divide="raise")
|
||||
|
||||
|
||||
class FractalSurface:
|
||||
"""Fractal surfaces."""
|
||||
|
||||
surfaceIDs = ["xminus", "xplus", "yminus", "yplus", "zminus", "zplus"]
|
||||
|
||||
def __init__(self, xs, xf, ys, yf, zs, zf, dimension, seed):
|
||||
"""
|
||||
Args:
|
||||
xs, xf, ys, yf, zs, zf: floats for the extent of the fractal surface
|
||||
(one pair of coordinates must be equal
|
||||
to correctly define a surface).
|
||||
dimension: float for the fractal dimension that controls the fractal
|
||||
distribution.
|
||||
seed: int for seed value for random number generator.
|
||||
"""
|
||||
|
||||
self.ID = None
|
||||
self.surfaceID = None
|
||||
self.start = np.array([xs, ys, zs], dtype=np.int32)
|
||||
self.stop = np.array([xf, yf, zf], dtype=np.int32)
|
||||
self.dtype = np.dtype(np.complex128)
|
||||
self.seed = seed
|
||||
self.dimension = (
|
||||
dimension # Fractal dimension from: http://dx.doi.org/10.1017/CBO9781139174695
|
||||
)
|
||||
self.weighting = np.array([1, 1], dtype=np.float64)
|
||||
self.fractalrange = (0, 0)
|
||||
self.filldepth = 0
|
||||
self.grass = []
|
||||
|
||||
@property
|
||||
def xs(self) -> int:
|
||||
return self.start[0]
|
||||
|
||||
@xs.setter
|
||||
def xs(self, value: int):
|
||||
self.start[0] = value
|
||||
|
||||
@property
|
||||
def ys(self) -> int:
|
||||
return self.start[1]
|
||||
|
||||
@ys.setter
|
||||
def ys(self, value: int):
|
||||
self.start[1] = value
|
||||
|
||||
@property
|
||||
def zs(self) -> int:
|
||||
return self.start[2]
|
||||
|
||||
@zs.setter
|
||||
def zs(self, value: int):
|
||||
self.start[2] = value
|
||||
|
||||
@property
|
||||
def xf(self) -> int:
|
||||
return self.stop[0]
|
||||
|
||||
@xf.setter
|
||||
def xf(self, value: int):
|
||||
self.stop[0] = value
|
||||
|
||||
@property
|
||||
def yf(self) -> int:
|
||||
return self.stop[1]
|
||||
|
||||
@yf.setter
|
||||
def yf(self, value: int):
|
||||
self.stop[1] = value
|
||||
|
||||
@property
|
||||
def zf(self) -> int:
|
||||
return self.stop[2]
|
||||
|
||||
@zf.setter
|
||||
def zf(self, value: int):
|
||||
self.stop[2] = value
|
||||
|
||||
@property
|
||||
def size(self) -> npt.NDArray[np.int32]:
|
||||
return self.stop - self.start
|
||||
|
||||
@property
|
||||
def nx(self) -> int:
|
||||
return self.xf - self.xs
|
||||
|
||||
@property
|
||||
def ny(self) -> int:
|
||||
return self.yf - self.ys
|
||||
|
||||
@property
|
||||
def nz(self) -> int:
|
||||
return self.zf - self.zs
|
||||
|
||||
def get_surface_dims(self):
|
||||
"""Gets the dimensions of the fractal surface based on surface plane."""
|
||||
|
||||
if self.xs == self.xf:
|
||||
surfacedims = (self.ny, self.nz)
|
||||
elif self.ys == self.yf:
|
||||
surfacedims = (self.nx, self.nz)
|
||||
elif self.zs == self.zf:
|
||||
surfacedims = (self.nx, self.ny)
|
||||
|
||||
return surfacedims
|
||||
|
||||
def generate_fractal_surface(self):
|
||||
"""Generate a 2D array with a fractal distribution."""
|
||||
|
||||
surfacedims = self.get_surface_dims()
|
||||
|
||||
self.fractalsurface = np.zeros(surfacedims, dtype=self.dtype)
|
||||
|
||||
# Positional vector at centre of array, scaled by weighting
|
||||
v1 = np.array(
|
||||
[
|
||||
self.weighting[0] * (surfacedims[0]) / 2,
|
||||
self.weighting[1] * (surfacedims[1]) / 2,
|
||||
]
|
||||
)
|
||||
|
||||
# 2D array of random numbers to be convolved with the fractal function
|
||||
rng = np.random.default_rng(seed=self.seed)
|
||||
A = rng.standard_normal(size=(surfacedims[0], surfacedims[1]))
|
||||
|
||||
# 2D FFT
|
||||
A = fftpack.fftn(A)
|
||||
|
||||
# Generate fractal
|
||||
generate_fractal2D(
|
||||
surfacedims[0],
|
||||
surfacedims[1],
|
||||
0,
|
||||
0,
|
||||
surfacedims[0],
|
||||
surfacedims[1],
|
||||
config.get_model_config().ompthreads,
|
||||
self.dimension,
|
||||
self.weighting,
|
||||
v1,
|
||||
A,
|
||||
self.fractalsurface,
|
||||
)
|
||||
|
||||
# Set DC component of FFT to zero
|
||||
self.fractalsurface[0, 0] = 0
|
||||
# Take the real part (numerical errors can give rise to an imaginary part)
|
||||
# of the IFFT, and convert type to floattype. N.B calculation of fractals
|
||||
# must always be carried out at double precision, i.e. float64, complex128
|
||||
self.fractalsurface = np.real(fftpack.ifftn(self.fractalsurface)).astype(
|
||||
config.sim_config.dtypes["float_or_double"], copy=False
|
||||
)
|
||||
# Scale the fractal volume according to requested range
|
||||
fractalmin = np.amin(self.fractalsurface)
|
||||
fractalmax = np.amax(self.fractalsurface)
|
||||
fractalrange = fractalmax - fractalmin
|
||||
self.fractalsurface = (
|
||||
self.fractalsurface * ((self.fractalrange[1] - self.fractalrange[0]) / fractalrange)
|
||||
+ self.fractalrange[0]
|
||||
- ((self.fractalrange[1] - self.fractalrange[0]) / fractalrange) * fractalmin
|
||||
)
|
||||
|
||||
|
||||
class FractalVolume:
|
||||
"""Fractal volumes."""
|
||||
|
||||
@@ -680,83 +515,3 @@ class MPIFractalVolume(FractalVolume):
|
||||
mpi_type = MPI.FLOAT.Create_subarray(shape.tolist(), subshape.tolist(), starts.tolist())
|
||||
mpi_type.Commit()
|
||||
return mpi_type
|
||||
|
||||
|
||||
class Grass:
|
||||
"""Geometry information for blades of grass."""
|
||||
|
||||
def __init__(self, numblades, seed):
|
||||
"""
|
||||
Args:
|
||||
numblades: int for the number of blades of grass.
|
||||
seed: int for seed value for random number generator.
|
||||
"""
|
||||
|
||||
self.numblades = numblades
|
||||
self.geometryparams = np.zeros(
|
||||
(self.numblades, 6), dtype=config.sim_config.dtypes["float_or_double"]
|
||||
)
|
||||
self.seed = seed
|
||||
self.set_geometry_parameters()
|
||||
|
||||
def set_geometry_parameters(self):
|
||||
"""Sets randomly defined parameters that will be used to calculate
|
||||
blade and root geometries.
|
||||
"""
|
||||
|
||||
self.R1 = np.random.default_rng(seed=self.seed)
|
||||
self.R2 = np.random.default_rng(seed=self.seed)
|
||||
self.R3 = np.random.default_rng(seed=self.seed)
|
||||
self.R4 = np.random.default_rng(seed=self.seed)
|
||||
self.R5 = np.random.default_rng(seed=self.seed)
|
||||
self.R6 = np.random.default_rng(seed=self.seed)
|
||||
|
||||
for i in range(self.numblades):
|
||||
self.geometryparams[i, 0] = 10 + 20 * self.R1.random()
|
||||
self.geometryparams[i, 1] = 10 + 20 * self.R2.random()
|
||||
self.geometryparams[i, 2] = self.R3.choice([-1, 1])
|
||||
self.geometryparams[i, 3] = self.R4.choice([-1, 1])
|
||||
|
||||
def calculate_blade_geometry(self, blade, height):
|
||||
"""Calculates the x and y coordinates for a given height of grass blade.
|
||||
|
||||
Args:
|
||||
blade: int for the numeric ID of grass blade.
|
||||
height: float for the height of grass blade.
|
||||
|
||||
Returns:
|
||||
x, y: floats for the x and y coordinates of grass blade.
|
||||
"""
|
||||
|
||||
x = (
|
||||
self.geometryparams[blade, 2]
|
||||
* (height / self.geometryparams[blade, 0])
|
||||
* (height / self.geometryparams[blade, 0])
|
||||
)
|
||||
y = (
|
||||
self.geometryparams[blade, 3]
|
||||
* (height / self.geometryparams[blade, 1])
|
||||
* (height / self.geometryparams[blade, 1])
|
||||
)
|
||||
x = round_value(x)
|
||||
y = round_value(y)
|
||||
|
||||
return x, y
|
||||
|
||||
def calculate_root_geometry(self, root, depth):
|
||||
"""Calculates the x and y coordinates for a given depth of grass root.
|
||||
|
||||
Args:
|
||||
root: int for the umeric ID of grass root.
|
||||
depth: float for the depth of grass root.
|
||||
|
||||
Returns:
|
||||
x, y: floats for the x and y coordinates of grass root.
|
||||
"""
|
||||
|
||||
self.geometryparams[root, 4] += -1 + 2 * self.R5.random()
|
||||
self.geometryparams[root, 5] += -1 + 2 * self.R6.random()
|
||||
x = round(self.geometryparams[root, 4])
|
||||
y = round(self.geometryparams[root, 5])
|
||||
|
||||
return x, y
|
86
gprMax/fractals/grass.py
普通文件
86
gprMax/fractals/grass.py
普通文件
@@ -0,0 +1,86 @@
|
||||
import numpy as np
|
||||
|
||||
from gprMax import config
|
||||
from gprMax.utilities.utilities import round_value
|
||||
|
||||
np.seterr(divide="raise")
|
||||
|
||||
|
||||
class Grass:
|
||||
"""Geometry information for blades of grass."""
|
||||
|
||||
def __init__(self, numblades, seed):
|
||||
"""
|
||||
Args:
|
||||
numblades: int for the number of blades of grass.
|
||||
seed: int for seed value for random number generator.
|
||||
"""
|
||||
|
||||
self.numblades = numblades
|
||||
self.geometryparams = np.zeros(
|
||||
(self.numblades, 6), dtype=config.sim_config.dtypes["float_or_double"]
|
||||
)
|
||||
self.seed = seed
|
||||
self.set_geometry_parameters()
|
||||
|
||||
def set_geometry_parameters(self):
|
||||
"""Sets randomly defined parameters that will be used to calculate
|
||||
blade and root geometries.
|
||||
"""
|
||||
|
||||
self.R1 = np.random.default_rng(seed=self.seed)
|
||||
self.R2 = np.random.default_rng(seed=self.seed)
|
||||
self.R3 = np.random.default_rng(seed=self.seed)
|
||||
self.R4 = np.random.default_rng(seed=self.seed)
|
||||
self.R5 = np.random.default_rng(seed=self.seed)
|
||||
self.R6 = np.random.default_rng(seed=self.seed)
|
||||
|
||||
for i in range(self.numblades):
|
||||
self.geometryparams[i, 0] = 10 + 20 * self.R1.random()
|
||||
self.geometryparams[i, 1] = 10 + 20 * self.R2.random()
|
||||
self.geometryparams[i, 2] = self.R3.choice([-1, 1])
|
||||
self.geometryparams[i, 3] = self.R4.choice([-1, 1])
|
||||
|
||||
def calculate_blade_geometry(self, blade, height):
|
||||
"""Calculates the x and y coordinates for a given height of grass blade.
|
||||
|
||||
Args:
|
||||
blade: int for the numeric ID of grass blade.
|
||||
height: float for the height of grass blade.
|
||||
|
||||
Returns:
|
||||
x, y: floats for the x and y coordinates of grass blade.
|
||||
"""
|
||||
|
||||
x = (
|
||||
self.geometryparams[blade, 2]
|
||||
* (height / self.geometryparams[blade, 0])
|
||||
* (height / self.geometryparams[blade, 0])
|
||||
)
|
||||
y = (
|
||||
self.geometryparams[blade, 3]
|
||||
* (height / self.geometryparams[blade, 1])
|
||||
* (height / self.geometryparams[blade, 1])
|
||||
)
|
||||
x = round_value(x)
|
||||
y = round_value(y)
|
||||
|
||||
return x, y
|
||||
|
||||
def calculate_root_geometry(self, root, depth):
|
||||
"""Calculates the x and y coordinates for a given depth of grass root.
|
||||
|
||||
Args:
|
||||
root: int for the umeric ID of grass root.
|
||||
depth: float for the depth of grass root.
|
||||
|
||||
Returns:
|
||||
x, y: floats for the x and y coordinates of grass root.
|
||||
"""
|
||||
|
||||
self.geometryparams[root, 4] += -1 + 2 * self.R5.random()
|
||||
self.geometryparams[root, 5] += -1 + 2 * self.R6.random()
|
||||
x = round(self.geometryparams[root, 4])
|
||||
y = round(self.geometryparams[root, 5])
|
||||
|
||||
return x, y
|
@@ -20,7 +20,8 @@ import logging
|
||||
|
||||
import numpy as np
|
||||
|
||||
from gprMax.fractals.fractals import FractalSurface, Grass
|
||||
from gprMax.fractals.fractal_surface import FractalSurface
|
||||
from gprMax.fractals.grass import Grass
|
||||
from gprMax.grid.fdtd_grid import FDTDGrid
|
||||
from gprMax.materials import create_grass
|
||||
from gprMax.user_objects.rotatable import RotatableMixin
|
||||
|
@@ -20,7 +20,7 @@ import logging
|
||||
|
||||
import numpy as np
|
||||
|
||||
from gprMax.fractals.fractals import FractalSurface
|
||||
from gprMax.fractals.fractal_surface import FractalSurface
|
||||
from gprMax.grid.fdtd_grid import FDTDGrid
|
||||
from gprMax.user_objects.rotatable import RotatableMixin
|
||||
from gprMax.user_objects.user_objects import GeometryUserObject
|
||||
|
@@ -22,7 +22,7 @@ import numpy as np
|
||||
|
||||
import gprMax.config as config
|
||||
from gprMax.cython.geometry_primitives import build_voxels_from_array, build_voxels_from_array_mask
|
||||
from gprMax.fractals.fractals import FractalVolume, MPIFractalVolume
|
||||
from gprMax.fractals.fractal_volume import FractalVolume, MPIFractalVolume
|
||||
from gprMax.grid.fdtd_grid import FDTDGrid
|
||||
from gprMax.grid.mpi_grid import MPIGrid
|
||||
from gprMax.materials import ListMaterial
|
||||
|
在新工单中引用
屏蔽一个用户