custEM.meshgen package
Submodules
custEM.meshgen.bathy_tools module
@author: Rochlitz.R
- class custEM.meshgen.bathy_tools.Bathy(bathy_file, x=None, y=None, centering=False, easting_shift=None, northing_shift=None, min_dist=50.0, min_points=10, coast_refinement=None, do_nothing=False)[source]
Bases:
object
Coastline identification class for digital elevation models. Please note that the basic methods are working properly for a couple of specific models, but a more general application of this class needs some debugging and further extensions.
- add_point_to_coast_line(idx, direction, clockwise, tol=None)[source]
Add a coordinate to a coastline path if the distance to the previous position on the coastline is larger than min_dist.
Required arguments
- idx, type list of [int, int]
x and y index of the DEM grid
- direction, type str
direction string, either ‘w’, ‘a’, ‘s’, or ‘d’
- clockwise, type bool
go in clock- or counterclockwise direction
Keyword arguments
- tol = None, type float
custom search tolerance if required, not recommended to be changed
- add_point_towards_east(idx)[source]
Add a coordinate to a coastline path in eastern direction.
Required arguments
- idx, type list of [int, int]
x and y index of the DEM grid
- add_point_towards_north(idx)[source]
Add a coordinate to a coastline path in northern direction.
Required arguments
- idx, type list of [int, int]
x and y index of the DEM grid
- add_point_towards_south(idx)[source]
Add a coordinate to a coastline path in southern direction.
Required arguments
- idx, type list of [int, int]
x and y index of the DEM grid
- add_point_towards_west(idx)[source]
Add a coordinate to a coastline path in western direction.
Required arguments
- idx, type list of [int, int]
x and y index of the DEM grid
- build_closing_frame(idx, clockwise)[source]
Continue coastline path in clockwise or counterclockwise direction if the horizontal mesh boundary is reached.
Required arguments
- idx, type list of [int, int]
x and y index of the DEM grid
- clockwise, type bool
go in clock- or counterclockwise direction
- check_distance(x, y)[source]
Calculate distance between a coordinate [x, y] and the previous position on the coastline.
Required arguments
- x, y, typefloat
x and y coordinates for the horizontal distance evaluation
- check_neighbours(kk, jj)[source]
Check z-values of the eight neighbors around a position in the DEM. Utility function that needs an update.
Required arguments
- kk, jj, type int
x/y indices of the DEM grid
- find_direction(idx, from_boundary=False)[source]
Identifiy clock- or counterclockwise search path direction based on regular grid DEM data automatically. Does not work in all cases so far.
Required arguments
- idx, type list of [int, int]
x and y index of the DEM grid to start the direction search from
Keyword arguments
- from_boundary = False, type bool
internally used flag to check if a coastline path was following the horizontal mesh boundaries before entering the grid again
- find_last_dir(pre_idx, idx)[source]
Evaluate previous direction during the coastline search.
Required arguments
- pre_idx, type list of [int, int]
previous x and y index of the DEM grid
- idx, type list of [int, int]
x and y index of the DEM grid
- find_start_point(jjj=0)[source]
Evaluate start indices for a new coastline search in the DEM.
Keyword arguments
- jjj, type int
counter to continute until all indices in the DEM-grid were considered
- go_clockwise(side, idx, from_corner=False)[source]
Continue coastline path in clockwise direction if the horizontal mesh boundary is reached.
Required arguments
- side, type str
string specifying one of the four sides
- idx, type list of [int, int]
x and y index of the DEM grid
- from_corner = False, type bool
internally used flag to continue a path after a corner of the horizontal mesh extent was reached
- go_counter_clockwise(side, idx, from_corner=False)[source]
Continue coastline path in counterclockwise direction if the horizontal mesh boundary is reached.
Required arguments
- side, type str
string specifying one of the four sides
- idx, type list of [int, int]
x and y index of the DEM grid
- from_corner = False, type bool
internally used flag to continue a path after a corner of the horizontal mesh extent was reached
- loadASC(topo_file, centering, easting_shift, northing_shift)[source]
Load asc file (DEM matrix with location header). This function is redudant to an identical function in the DEM class and will be removed in a major update of the bathy-tools.
Required arguments
see Bathy class description
- loadTXT(topo_file, centering, easting_shift, northing_shift)[source]
Load txt file with column-based list of DEM coordinates. This function is redudant to an identical function in the DEM class and will be removed in a major update of the bathy-tools.
Required arguments
see Bathy class description
custEM.meshgen.dem_interpolator module
@authors: Rochlitz.R & Günther.T
- class custEM.meshgen.dem_interpolator.DEM(demfile, x=None, y=None, centering=False, easting_shift=None, northing_shift=None, revert_z=False)[source]
Bases:
object
Interpolation class for digital elevation models (DEM). Automatically creates interpolation objects during the class initialization, based on DEM information stored as list of coordinates in txt format or as matrix in asc format.
- loadASC(ascfile, centering, easting_shift, northing_shift)[source]
Load ASC file (DEM matrix with location header).
Required arguments
see DEM class description
- loadTXT(demfile, centering, easting_shift, northing_shift)[source]
Load txt file with column-based list of DEM coordinates.
Required arguments
see DEM class description
- show(cmap='terrain', cbar=True, ax=None, **kwargs)[source]
Show digital elevation model. Needs an update.
Keyword arguments
- cmap = “terrain”, type str ()
matplotlib colormap definiton
- cbar = True, type bool
add colorbar to the plot or not
- ax = None, type matplotlib figure axes object
add the plot to a given axes object or create a new one
- **kwargs, type keyword arguments
add additional keyword arguments for the plot style (e.g., lw)
- translate(centering=False, easting_shift=None, northing_shift=None, back=False)[source]
Function to translate coordinates in horizontal direction.
Keyword arguments
see DEM class description
- back = False, type bool
set True to translate back in the opposite direction; useful for rotating coordinates around a specific origin
custEM.meshgen.invmesh_tools module
@author: Rochlitz.R
- class custEM.meshgen.invmesh_tools.PrismWorld(name, x_extent, x_reduction=100.0, y_depth=300.0, z_depth=300.0, prism_quality=32.0, prism_area=10000.0, n_prisms=50.0, orthogonal_txs=None, txs=[], surface_rxs=[], **blankworld_kwargs)[source]
Bases:
object
Specific class for automated generation of pseudo-2D / 3D inversion meshes discretized by prisms.
- add_prisms(back_nodes, front_nodes)[source]
Add nodes and triangle faces to include horizontal prisms as conformal PLC.
- add_triangles(back_nodes, front_nodes)[source]
Add front and back triangle faces of the prisms. All the if statements are required to account for potential nodes on the x-directed boundary edges of the rectangular prism outcrops at the surface.
- create_2d_mesh(name, prism_quality, prism_area)[source]
Create and store 2D mesh in x-z slice with triangle, which will be the basis for elongating the triangles to prisms in y-direction later on.
- create_surface_rectangles()[source]
Build rectangular frames of the prisms at the surface to be included in the surface mesh later on for avoiding intersections. Check also if transmitter paths are intersecting these rectangles, and if yes, find the intersections and add the corresponding points to the frames.
custEM.meshgen.mesh_convert module
@author: modified by Rochlitz.R
- custEM.meshgen.mesh_convert.mesh2xml2h5(mesh_name, m_dir, rot)[source]
Converts TetGen meshes in Medit (.mesh) format to xml and h5 files for the usage in FEniCS. This function is called automatically during the intialization of the MOD class and converts a mesh if it does not already exist as xml and h5 file.
Required arguments
- mesh_name, type str
mesh name
- m_dir, type str
path to mesh directory
Functionality
Convert between .mesh and .xml, parser implemented as a state machine:
0 = read ‘Dimension’ 1 = read dimension 2 = read ‘Vertices’ 3 = read number of vertices 4 = read next vertex 5 = read ‘Triangles’ or ‘Tetrahedra’ 6 = read number of cells 7 = read next cell 8 = read next domains 9 = done
- custEM.meshgen.mesh_convert.write_cell_interval(ofile, cell, n0, n1)[source]
Write interval (1D) information to domain file.
- custEM.meshgen.mesh_convert.write_cell_tetrahedron(ofile, cell, n0, n1, n2, n3)[source]
Write tetrahedron 3D information to file.
- custEM.meshgen.mesh_convert.write_cell_triangle(ofile, cell, n0, n1, n2)[source]
Write triangle 2D information to domain file.
- custEM.meshgen.mesh_convert.write_domain_marker(ofile, cell, marker, dfile=False)[source]
Write domain markers to domain file.
Write footer cells to file.
Write footer to domain file.
Write footer edges to file.
Write footer graph to file.
Write footer of mesh file.
Write footer vertices to file.
- custEM.meshgen.mesh_convert.write_graph_edge(ofile, v1, v2, weight=1)[source]
Write edge graph to file.
- custEM.meshgen.mesh_convert.write_graph_vertex(ofile, vertex, num_edges, weight=1)[source]
Write vertex graph to file.
- custEM.meshgen.mesh_convert.write_header_cells(ofile, num_cells)[source]
Write header cells to file.
- custEM.meshgen.mesh_convert.write_header_domains(ofile, num_domains, dfile=False)[source]
Write header to domain file.
- custEM.meshgen.mesh_convert.write_header_edges(ofile, num_edges)[source]
Write header edges to file.
- custEM.meshgen.mesh_convert.write_header_graph(ofile, graph_type)[source]
Write header graph to file.
- custEM.meshgen.mesh_convert.write_header_mesh(ofile, cell_type, dim)[source]
Write header of mesh file.
- custEM.meshgen.mesh_convert.write_header_vertices(ofile, num_vertices)[source]
Write header vertices to file.
custEM.meshgen.meshgen_tools module
@author: Rochlitz.R
- class custEM.meshgen.meshgen_tools.BlankWorld(**raw_kwargs)[source]
Bases:
object
Toolbox for mesh generation using pyGIMLi functionalities to generate a world Omega, add structures such as interfaces, anomalies, transmitters etc. to this world and finally call TetGen to built a mesh based on the exported polyfile.
Important notice
This class was devloped over four years and provides useful functionalities for the mesh generation, but there are many parts which need an overhaul in the future. Maybe an updated version of these tools will be moved to pyGIMLi in future and updated to work with TetGen 1.6
- add_air_tri(paths)[source]
Add refinement trianles to Omega.
Required arguments:
- paths, type list
paths is a list of pointsets (arrays with shape (n, 3)) which will be connected from point to point.
- add_area_marker(poly, pos, face_size=0.0, marker=777)[source]
Add region marker to the surface or interfaces for 2D meshing, solely used internally to define area-constraints for cooresponding meshes.
Required arguments:
- poly, type pyGIMLi Polygon
polygon to which the marker is defined.
- pos, type list or array of length (3,)
x, y, z coordinates where the marker should be placed, note that this point must be located inside the correct domain!
Keyword arguments:
- face_size = 0., type int
maximum facet-area contstraint, default (0.) means no constraint.
- marker = 777, type int
region marker to define sub-areas within the poly.
- add_brick(start, stop, cell_size=0.0, marker=None)[source]
Add rectangular brick anomaly to Omega. Note that this is a special case of the add_plate() method, which can be used alternatively.
Required arguments:
- start, type list
list of back lower left corner [x_min,, y_min, z_min].
- stop, type list
list of front upper right corner [x_max,, y_max, z_max].
Keyword arguments:
- cell_size = 0., type float
maximum cell-volume contstraint, default (0.) means no constraint.
- add_cylinder(r, h, origin=[0.0, 0.0, -200.0], n_segs=16, dip=0.0, dip_azimuth=0.0, cell_size=0.0, marker=None, marker_position=None)[source]
Add cylinder-anomaly to Omega.
Required arguments
- r, type float
radius of the cylinder
- h, type float
height of the cylinder
Keyword arguments
- origin = [0., 0., -200.], type list/array of three floats
center (of mass) coordinate of the cylinder; it is rotated in case of non-zero dip or dip azimuth values around this point
- n_segs, type int
number of segments to approximate the circumference
- dip = 0., type float,
dip angle of the plate in degree
- dip_azimuth = 0., type float
dip azimuth angle in degree, zero points towards positive x-direction
- cell_size = 0., type float
maximum cell-volume contstraint, default (0.) means no constraint.
- marker = None, type int
specify custom marker for plate if required
- add_hollow_cylinder(r1, r2, h=100.0, origin=[0.0, 0.0, 0.0], rotation=None, cell_size=0.0, n_segs=16, preserve=[], close=False, mark=False, close1st=False)[source]
Add hollow cylinder-anomaly to Omega. In progress!
Required arguments:
- r1, type float
inner radius of cylinder.
- r2, type float
outer radius of cylinder.
Keyword arguments:
- cell_size = 0., type float
maximum cell-volume contstraint, default (0.) means no constraint.
- add_interface(frame=None, topo=None, id_offset=None, eval_zl=False)[source]
Add more or less horizontal 2D interface (layer interfaces) with or without topography to Omega.
- add_intersecting_anomaly(surface_intersection=False, surface_path=[], intersection_paths=[], intersecting_layers=None, cell_size=None, bottom=None, top=None, marker=None, marker_position=None)[source]
- add_inv_domains(depth, poly=None, rxs=None, cell_size=0.0, x_frame=1000.0, y_frame=1000.0, z_frame=1000.0, d=100.0, marker_positions=None, overwrite_markers=None, split_depths=None)[source]
- add_marker(label, pos, marker=None, cell_size=0.0)[source]
Add region marker to domains. This function is called internally for each layer or anomaly added. By default, domains are numbered consequtively. It is usually not recommended not to maunaly set domain marker, even though the possibility is given by overwriting the keyword argument marker=None when adding, e.g., anomalies.
Required arguments:
- label, type str
label for defining the domain information print output after calling TetGen, which is enabled by default.
- pos, type list or array of length (3,)
x, y, z coordinates where the marker should be placed, note that this point must be located inside the correct domain!
Keyword arguments:
- marker = None, type int
set manually defined region marker for this domain.
- cell_size = 0., type float
maximum cell-volume contstraint, default (0.) means no constraint.
- add_paths(paths, closed_paths=None, marker=None)[source]
Add (observation) lines to Omega.
Required arguments:
- paths, type list
paths is a list of pointsets (arrays with shape (n, 3)) which will be connected from point to point.
Keyword arguments:
- closed_paths = None, type list
a list of same length as paths is required with True or False values defining which path should be closed and which not.
- add_plate(dx=200.0, dy=200.0, dz=50.0, origin=[0.0, 0.0, -200.0], dip=0.0, dip_azimuth=0.0, cell_size=0.0, marker=None)[source]
Add rectangular plate anomaly to Omega. Rectangular bodys can be added by using the method add_brick().
Keyword arguments
- dx, dy = 200., type float
dimensions in x/y-directions (before rotation)
- dz = 50., type float
thickness of the plate
- origin = [0., 0., -200.], type list/array of three floats
center (of mass) coordinate of the plate; plate is rotated in case of non-zero dip or dip azimuth values around this point
- dip = 0., type float,
dip angle of the plate in degree
- dip_azimuth = 0., type float
dip azimuth angle in degree, zero points towards positive x-direction
- cell_size = 0., type float
maximum cell-volume contstraint, default (0.) means no constraint.
- marker = None, type int
specify custom marker for plate if required
- add_points(points, marker=None)[source]
Add points to Omega.
Required arguments
- points, type array with shape (n, 3)
arrqay of points to be added
Keyword arguments
- marker = None, type int
specify custom marker for points if required
- add_prism(poly, top, bottom, cell_size=0.0, marker=None, marker_pos=None)[source]
Add rectangular brick anomaly to Omega. Note that this is a special case of the add_plate() method, which can be used alternatively.
Required arguments:
- start, type list
list of back lower left corner [x_min,, y_min, z_min].
- stop, type list
list of front upper right corner [x_max,, y_max, z_max].
Keyword arguments:
- cell_size = 0., type float
maximum cell-volume contstraint, default (0.) means no constraint.
- add_quadrangle_face(n1, n2, n3, n4)[source]
Add quadrangle facet to Omega defined by nodes n1 - n4.
Required arguments:
- n1, n2, n3 n4, type int
four node id’s to build a quandrangle.
- add_rx(points, marker=None)[source]
Add recevier points for automated interpolation to Omega.
Required arguments
- points, type list or array with shape (n, 3)
list of receiver points to add
- add_surface_anomaly(insert_paths=[], depths=None, cell_sizes=None, marker=None, dips=None, dip_azimuths=None, split_depths=None, marker_positions=None, split_coords=None)[source]
- add_topo(node_positions, topo=None, eval_zl=False)[source]
Apply topography, either DEM information or synthetic descriptions as z = f(x, y) to the z-values of an interface mesh.
Required arguments:
- node_positions, type array of shape (n, 3)
positions of all nodes of the interface mesh
Keyword arguments:
- topo = None, type str or function
if None, self.topo will be used for the earth’s surface, for further interfaces, either additional DEM information can be used or a synthetic definition via a function f(x, y)
- add_tx(tx, grounded=True, marker=None)[source]
Add transmitter paths to Omega.
Required arguments:
- tx, type list
list of pointsets (arrays with shape (n, 3)) which will be connected from point to point.
Keyword arguments:
- grounded = False, type bool
close first and last point of path to build a grounded loop or not.
- build_fullspace_mesh(cell_size=None, boundary_marker=99)[source]
Create a fullspace mesh.
Keyword arguments:
- cell_size = 0., type float
create a fullspace mesh with cells of maximum volume ([m^3]) of max_cell_size
- build_halfspace_mesh(boundary_marker=99, **mesh_build_kwargs)[source]
Create halfspace mesh.
Keyword arguments:
- q = 34., type float
quality (see triangle documentation) of the surface mesh
- outer_area_cell_size = None, type float
if not specified when initializing the BlankWorld instance, set or owerwrite the maximum area of cells on the surface for the “outer area” here
- inner_area_cell_size = None, type float
if not specified when initializing the BlankWorld instance, set or owerwrite the maximum area of cells on the surface for the “inner area” here
- build_layered_earth_mesh(n_layers, layer_depths, insert_struct=None, closed_struct=None, tx_struct=None, struct_interface=0, boundary_marker=99, **le_kwargs)[source]
Create layered earth mesh.
Required arguments:
- n_layers, type int
number of subsurface layers
- layer_depths, type list of floats
list with len == n_layers - 1, specifying the depths of 1D subsurface layer interfaces. the z-axis is pointing upwards, so negative values correspond to increasing depths. This argument will be overwritten, if subsurface topography is applied
- build_surface(insert_line_tx=[], insert_loop_tx=[], insert_points=[], insert_lines=[], insert_loops=[], insert_paths=[], insert_ground_tri=[], closed_path=False, subsurface_interface=False, line_marker=None, loop_marker=None, path_marker=None, line_tx_marker=None, loop_tx_marker=None, point_marker=None, add_to_existing=False)[source]
Create pyGIMLi surface polygon for 2D meshing with “triangle” later on. Note! Structures at the surface such as CSEM transmitters or observation lines / points have to incorporated here.
Keyword arguments
- insert_line_tx = [], type list
list of grounded transmitter paths, numbered successively
- insert_loop_tx = [], type list
list of ungrounded transmitter paths; note, loop_tx will be always added and numbered after all line tx
- insert_points = [], type list
list of meshgen_utils “pointSet” definitions
- insert_lines = [], type list
list of meshgen_utils “line” definitions
- insert_loops = [], type list
list of meshgen_utils “loop” definitions
- insert_paths = [], type list
list of paths descriped by coordinates
- closed_path = False, type bool or list of len == len(insert_paths)
close the segsment between the first and last point of a given “PATH” to either generate a grounded (line) or ungrounded (loop) CSEM source or in general, open or closed polygons
- subsurface_interface = False, type bool
for internal usage only, specifies if the 2D polygon will be on the surface or a subsurface interface
- call_tetgen(tet_param='default', ignore_tx_z=False, export_before=True, export_format=' -k', name=None, copy=True, print_infos=True, suppress='NEF')[source]
Wrapper-function to call TetGen version 1.5 and copy the resulting mesh into the default “m_dir”/_mesh directory.
Keyword arguments:
- tet_param = ‘default’, type str
specify your parameters which should be passed to the TetGen call, e.g., ‘**’-pq1.4aYAO2/3’** without including flags for the export. Alternative keyword arguments are ‘raw’ (‘**’-p’) or ‘faces’ (’’-pYA’**). It is also possible to directly pass the TetGen- options (aside from export), e.g., tet_param=’-pq1.2aMAO2/3’.
- export_before = True, type bool
automatically export the pyGIMLi mesh Omega as “.poly”-file in the working directory self.w_dir before calling TetGen.
- export_vtk = True, type bool
set False for not exporting .vtk-file to be viewed in Paraview.
- name = None, type str
specify custom export name for the mesh, defualt is self.name.
- copy = True, type bool
automatically copy the mesh in .mesh (Medit)-format to the save- directory self.s_dir for automated conversion on the first computation run.
- print_infos = True, type bool
set False if no domain information should be printed at the end.
- create_box()[source]
Create bounding box from given limits and initialize inner_area and outer_area properties, if a separation into a fine-discretized inner (observation) area at the surface and outer area (boundary-area) was specified during intitializing the BlankWorld class.
- extend_anomaly_outcrops()[source]
Function for extending anomaly-bodies into a layer. Completely automated for internal usage. For more details, it is referred to the add_surface_anomaly() method, where users can define the anomalies!
- extend_intersecting_anomaly()[source]
Function for extending anomaly-bodies which interesect layers.
- extend_world(x_fact=10.0, y_fact=10.0, z_fact=10.0, marker=[0, 1], cell_sizes=None, boundary_marker=99, initial=True)[source]
alias for there_is_always_a_bigger_world() method
- get_topo_vals(node_positions, topo=None, z=0.0)[source]
Apply topography, either DEM information or synthetic descriptions as z = f(x, y) to the z-values of an interface mesh.
Required arguments:
- node_positions, type array of shape (n, 3)
positions of all nodes of the interface mesh
Keyword arguments:
- topo = None, type str or function
if None, self.topo will be used for the earth’s surface, for further interfaces, either additional DEM information can be used or a synthetic definition via a function f(x, y)
- init_dem(surface_mesh_coords, topo=None)[source]
Initialize DEM instance for interpolating real-world topography from digital-elevation-model-files and check if the extent of the dem-file coordinates is large enough for the specified mesh dimensions.
Required arguments:
- surface_mesh_coords, type array of shape (n, 3)
array containing all coordinates of the surface mesh for comparing the spatial extent in x- and y-direction.
- remove_duplicate_poly_faces()[source]
Check which self.Omega.n_poly faces were added multiple times if subsurface structures touch each other and remove the duplicates.
- there_is_always_a_bigger_world(x_fact=10.0, y_fact=10.0, z_fact=10.0, marker=[0, 1], cell_sizes=None, boundary_marker=99, initial=True)[source]
Appends a bigger world (halfspace) to the existing world Omega. This bigger world is synonymous to a “tetrahedron boundary”.
Keyword arguments:
- x_fact, y_fact, z_fact = 10., type float
factor in x-, y- and z-direction multiplied to the original domain size specified by the **x/y/zlim* class attributes.
- marker = [0, 1], list of two int
specify your own marker values for the bounding halfspace mesh.
custEM.meshgen.meshgen_utils module
@author: Rochlitz.R
- custEM.meshgen.meshgen_utils.add_edges(surface, edges, nn, size, marker, closed=False)[source]
Add edge to polygon.
- custEM.meshgen.meshgen_utils.assign_topography(nodes, t_dir=None, topo=None, z=0.0, centering=True, easting_shift=None, northing_shift=None, revert_z=False, rotation=None)[source]
Assign topography on nodes.
Required arguments
- nodes, type list/array with shape (:, 3)
list or array of 3D coordinates
Keyword arguments
- t_dir = None, type str
directory of topography file
- topo = None, type function or str
topography function or file
- z = 0., type float
add constant offset to topography in z-direction
- centering = True, type bool
Set True if extent of topography file should be centered horizontally around origin (0, 0)
- easting_shift = None, type float
add specified shift to x-coordinates of topography file
- northing_shift = None, type float
add specified shift to y-coordinates of topography file
- rotation = None, type float
rotate complete topography file about specified angle(DEG) horizontally
- custEM.meshgen.meshgen_utils.closest_node(coord, coords)[source]
Find the closest coordinate in a list of 2D/3D corrdinates.
Required arguments
- node, list/array with length 3
reference node
- nodes, list/array with shape (:, 2) or (:, 3)
list or array of nodes to be checked
- custEM.meshgen.meshgen_utils.create_circle(r, n, h)[source]
Create a list of n (number of) 3D coordinates on the circumference of a circle with radius r and height h.
Required arguments
- r, type float
the radius of the circle
- n, type int
number of nodes on the cirumference of the circle
- h, type float
z-coordinates of the circle
- custEM.meshgen.meshgen_utils.create_ellipse(a, b, n, h)[source]
Create a list of n (number of) 3D coordinates on the circumference of an ellipsum with ‘radii’ a and b, and height h.
Required arguments
- a, b, type float
the two radius of the ellipse
- n, type int
number of nodes on the cirumference of the ellipse
- h, type float
z-coordinates of the ellipse
- custEM.meshgen.meshgen_utils.export_tetgen_edge_file(preserve_edges, edge_marker, filename)[source]
Write a given set of edges which should be preserved into a Ascii file in Tetgen .edge format.
Required arguments
- preserve_edges, type list
edges which should be exported to the .edge file
- edge_marker, type list
list of markers assigned to the edges
- filename, type str
name of output file
- custEM.meshgen.meshgen_utils.export_tetgen_node_file(preserve_nodes, filename)[source]
Write a given set of nodes which should be preserved into an Ascii file in Tetgen .edge format.
Required arguments
- preserve_nodes, type list
edges which should be exported to the .node file
- filename, type str
name of output file
- custEM.meshgen.meshgen_utils.export_tetgen_poly_file(poly, filename, float_format='.3f', poly_marker=None, **kwargs)[source]
Write a given piecewise linear complex (mesh/poly) into an Ascii file in Tetgen .poly format.
Required arguments
- poly, type pyGIMLi mesh
piecewise linear complex which should be exported to the .poly file
- filename, type str
name of output file
Keyword arguments
- float_format = ‘.12e’, type string
format that will be used to write float values in the Ascii file, default is the exponential float form with a precision of 12 digits
- custEM.meshgen.meshgen_utils.find_closest(l1, l2)[source]
Find the closest entry in a list of 1D corrdinates in comparison to all entries of another list.
Required arguments
- l1, type list
reference list of coordinates
- l2, type list
another list of coordinate, usually with a different length
- custEM.meshgen.meshgen_utils.find_edge_node_ids(frame, mesh, id_offset)[source]
Internal utility function for meshgen tools to find the edge node ids of a horizontal 2D mesh for the extension of the mesh in z-direction.
Required arguments
- frame, type array with shape (:, 3)
references frame defined by coordinates
- mesh, type pyGIMLi mesh
mesh, whose coordinates should be checked
- id_offset, type int
offset of ids between 2D surface meshes and 3D mesh (should be removed if the meshgen_tools get an overhaul)
- custEM.meshgen.meshgen_utils.find_nodes_on_segment(p0, p1, mesh, topo_f)[source]
Find nodes on a segment.
Required arguments
mesh
- custEM.meshgen.meshgen_utils.find_on_frame(frame_coords, coords)[source]
Check if coordinates are located on a frame.
Required arguments
- frame_coords, type array with shape (:, 3)
array of frame coordinates
- coords, type array with shape (:, 3)
array of coordinates to be checked
- custEM.meshgen.meshgen_utils.get_closest(coord, coords, n=None, max_dist=1.0)[source]
Find the ids the closest coordinates in a list of 2D/3D corrdinates, either using a distance threshold or by specifying the number of closests coordinates.
Required arguments
- coord, list/array with length 3
reference node
- coords, list/array with shape (:, 2) or (:, 3)
list or array of nodes to be checked
Keyword arguments
- n = None, type int
number of closests coordinates, overwrites max_dist
- max_dist = 1., type float
return all coordinate ids which are closer than max_dists
- custEM.meshgen.meshgen_utils.get_intersect(a1, a2, b1, b2)[source]
Return intersection coordinate of two segments a and b.
- custEM.meshgen.meshgen_utils.get_unique_poly_ids(polys, max_poly_length)[source]
Get unique ids of nodes of multiple polygons.
Required arguments
- polys, type list
list of polygons
- max_poly_length
maximum number of nodes of within all polygons
- custEM.meshgen.meshgen_utils.inside_poly(x, y, path)[source]
Return True if a coordinate (x, y) is inside a polygon defined by a list of verticies [(x1, y1), (x2, x2), … , (xN, yN)].
Reference: http://www.ariel.com.au/a/python-point-int-poly.html
Required arguments
- x, type float
x coordinate of point to be checked
- y, type float
y coordinate of point to be checked
- path, type list/array with shape (:, 3)
list or array of 3D coordinates
- custEM.meshgen.meshgen_utils.is_between(a, b, c, a_tol=0.01, print_points=False)[source]
Search if a Nedelec-dof is located on the edge between two points in 3D, needed for crooked sources defined as path. Search if point c is in between a and b.
Required arguments:
- a, b, c, type list/array of length 3
three points to be evaluated
Keyword arguments:
- a_tol = 1e-2, type float
absolute tolerance to account for numerical inaccuracies
- print_points = False, type bool
set True to enable print of coordinates
- custEM.meshgen.meshgen_utils.is_between_1D(a, b, c, a_tol=0.01, print_points=False)[source]
Search if a x-coordinate is located in an interval in in horizontal direction. Search if point c is in between a and b.
Required arguments
- a, b, c, type list/array of length 3
three points to be evaluated
Keyword arguments
- a_tol = 1e-2, type float
absolute tolerance to account for numerical inaccuracies
- print_points = False, type bool
set True to enable print of coordinates
- custEM.meshgen.meshgen_utils.is_between_2D(a, b, c, a_tol=0.01, print_points=False)[source]
Search if a Nedelec-dof is located on the edge between two points in 2D in horizontal direction. Search if point c is in between a and b.
Required arguments
- a, b, c, type list/array of length 3
three points to be evaluated
Keyword arguments
- a_tol = 1e-2, type float
absolute tolerance to account for numerical inaccuracies
- print_points = False, type bool
set True to enable print of coordinates
- custEM.meshgen.meshgen_utils.is_between_2D_exact(a, b, c, a_tol=0.01, print_points=False)[source]
Search if a Nedelec-dof is located on the edge between two points in 2D in horizontal direction. Search if point c is in between a and b.
Required arguments
- a, b, c, type list/array of length 3
three points to be evaluated
Keyword arguments
- a_tol = 1e-2, type float
absolute tolerance to account for numerical inaccuracies
- print_points = False, type bool
set True to enable print of coordinates
- custEM.meshgen.meshgen_utils.is_between_2Dvert(a, b, c, a_tol=0.01, print_points=False)[source]
Search if a Nedelec-dof is located on the edge between two points in 2D in vertical direction. Search if point c is in between a and b.
Required arguments
- a, b, c, type list/array of length 3
three points to be evaluated
Keyword arguments
- a_tol = 1e-2, type float
absolute tolerance to account for numerical inaccuracies
- print_points = False, type bool
set True to enable print of coordinates
- custEM.meshgen.meshgen_utils.line(start=[0.0, 0.0, 0.0], stop=[0.0, 0.0, 0.0], n_segs=100, z=0.0, topo=None, t_dir=None, **topo_kwargs)[source]
Return 3D coordinates of points on a line in arbitrary direction.
Keyword arguments
- start, stop = [0., 0., 0.], type list or array
start and stop points
- n_segs = 100, type int
number of segments on the line
- topo_kwargs
listed below
- z = 0., type float
offset of the line in z-direction
- topo = None, type str or python function
DEM-topography file or synthetic topography (python) function.
- t_dir = None, type str
directory of topography (DEM) files
- centering = False, type bool
set True, if the given DEM should be centered relative to the computational domain coordinate frame
- easting_shift, northing_shift = None, type float
if centering is False, a custom offset of a given DEM relaive to the computational domain in x- and y- direction can be set
- revert_z = False, type bool
revert z-values of the DEM if the original height information was provided in a donward-oriented coorindate system
- rotation (RAD) = None, type float
rotation of the DEM around the center coordinate (set via the centering flag or easting/northing_shift) relative to the computational domain, given in rad, not degree
- custEM.meshgen.meshgen_utils.line_x(start=-100.0, stop=100.0, n_segs=100, y=0.0, z=0.0, topo=None, t_dir=None, **topo_kwargs)[source]
Return 3D coordinates of points on a line in x-direction.
Keyword arguments
- start, stop = +-1e2, type float
start and stop values in x-direction
- n_segs = 100, type int
number of segments on the line
- y, z = 0., type float
offset of the line in y- and z-direction
- topo_kwargs
see description of line method or meshgen_tools/dem_interpolator
- custEM.meshgen.meshgen_utils.line_y(start=-100.0, stop=100.0, n_segs=100, x=0.0, z=0.0, topo=None, t_dir=None, **topo_kwargs)[source]
Return 3D coordinates of points on a line in y-direction.
Keyword arguments
- start, stop = +-1e2, type float
start and stop values in x-direction
- n_segs = 100, type int
number of segments on the line
- x, z = 0., type float
offset of the line in x- and z-direction
- topo_kwargs
see description of line method or meshgen_tools/dem_interpolator
- custEM.meshgen.meshgen_utils.line_z(start=-100.0, stop=100.0, n_segs=100, x=0.0, y=0.0)[source]
Return 3D coordinates of points on a line in z-direction.
Keyword arguments
- start, stop = +-1e2, type float
start and stop values in x-direction
- n_segs = 100, type int
number of segments on the line
- x, y = 0., type float
offset of the line in x- and y-direction
- custEM.meshgen.meshgen_utils.loop_c(origin=[0.0, 0.0, 0.0], r=50.0, n_segs=100, z=0.0, topo=None, t_dir=None, **topo_kwargs)[source]
Return 3D coordinates of points on a circular loop.
Keyword arguments
- origin = [0., 0., 0.], type list of floats
x-, y-, z- coordinates of the center of the loop
- r = 5e1, type float
radius of the loop
- n_segs = 100, type int
number of segments on the loop circumference
- topo_kwargs
see description of line method or meshgen_tools/dem_interpolator
- custEM.meshgen.meshgen_utils.loop_e(origin=[0.0, 0.0, 0.0], a=50.0, b=50.0, n_segs=100, z=0.0, topo=None, t_dir=None, **topo_kwargs)[source]
Return 3D coordinates of points on a ellipsoidal loop.
Keyword arguments
- origin = [0., 0., 0.], type list of floats
x-, y-, z- coordinates of the center of the loop
- a, b = 5e1, type float
radius of the loop
- n_segs = 100, type int
number of segments on the loop circumference
- topo_kwargs
see description of line method or meshgen_tools/dem_interpolator
- custEM.meshgen.meshgen_utils.loop_r(start=[-100.0, -100.0], stop=[100.0, 100.0], n_segs=4, nx=None, ny=None, z=0.0, topo=None, t_dir=None, loop_rotation=None, **topo_kwargs)[source]
Return 3D coordinates of points on a rectangular loop.
Keyword arguments
- start, stop = +-[1e2, 1e2], type list: [float, float]
opposite corner coordinates of rectangular loop
- n_segs = 4, type int
number of overall segments of the loop, this keyword argument is ignored if nx and ny are specified
- nx, ny = 100, type int
number of segments on the x- and y-directed edges of the loop
- loop_rotation = None,
rotate loop horizontally around angle in (RAD)
- topo_kwargs
see description of line method or meshgen_tools/dem_interpolator
- custEM.meshgen.meshgen_utils.npoly_to_triangles(world, ids)[source]
Add polygones with more than three edges in form of triangles to Omega.
- custEM.meshgen.meshgen_utils.order_list(ref, points, tmp)[source]
Sort values in list in increasing order. Deprecated and not used anymore.
- custEM.meshgen.meshgen_utils.pointset(pointset)[source]
Return a given set of 3D points, either list or array, as numpy array.
Required arguments
- pointset, type list
list of coordinates to be converted to an array of shape (:, 3)
- custEM.meshgen.meshgen_utils.poly_area(path)[source]
Caluclate enclosed area of polygon.
Required arguments
- path, type array with shape (:, 2) or (:, 3)
array of coordinates describing the polygon
- custEM.meshgen.meshgen_utils.poly_peri(path)[source]
Calculate perimeter of polygon.
Required arguments
- path, type array with shape (:, 2) or (:, 3)
array of coordinates describing the polygon
- custEM.meshgen.meshgen_utils.refine_adaptive(coords, txs, r1=20.0, r2=5.0, r3=1.0, d3=200.0, d2=500.0, min_tx_dist=100.0, rot=30.0, n_segs=3)[source]
Decription to add
- custEM.meshgen.meshgen_utils.refine_path(path, n_segs=None, length=None)[source]
Refine each section of a 3D path into n_segs equal segments.
Required arguments
- path, type list/array with shape (:, 3)
list or array of 3D coordinates
- n_segs, type int
number of nodes to be inserted within each segment of the path
- length, type float
maximum length of any segment of the path
- custEM.meshgen.meshgen_utils.refine_rx(coords, r=5.0, rot=None, n_segs=3)[source]
Refine rx position by enclosing them with hexagons of radius r. The smaller r, the better the refinement. Note that each hexagon will be subdivided into 6 equilateral triangles, which will support optimum tetrahedra quality during the meshing process.
Required arguments
- coords, type list/array with shape (:, 3)
list or array of 3D coordinates
Keyword arguments
- r = 5., type float
radius of surrounding equilateral (triangle if n_segs=3) shape
- rot = None, type float
horizontal rotation angle(DEG) between subsequent shapes
- n_segs = 3, type int
number of corners of surrounding equilateral polygon shape
- custEM.meshgen.meshgen_utils.rotate(array, alpha, axis='z', tensor=False)[source]
Rotate array about angle alpha(RAD) around given axis Default rotation axis is z, alternatively x or y.
Required arguments
- array, type array of shape (:, 3)
array of coordinates to be rotated
- alpha, type float
angle(RAD) specifying the rotation
Keyword arguments
- axis = ‘z’, type str
rotation axis
- custEM.meshgen.meshgen_utils.rotate_around_point(array, shift, angles, axes, pre_rotate=False)[source]
Rotate array around arbitrary point defined by shift, not around centroid.
Required arguments
- array, type array of shape (:, 3)
array of coordinates to be rotated
- shift, type list or array of length 3
x-, y-, and z- values defining the shift
- angles, type list of floats
list of angles(RAD) specifying the rotation
- axes, type list of str
list of rotation axes
Keyword arguments
- pre_rotate = False, type bool
conduct an initial roation using the second entry of axes/angles
- custEM.meshgen.meshgen_utils.rx_grid(xlinsp, ylinsp, z=0.0)[source]
Build Rx grid from two orthogonal vectors.
Required arguments
- xvec/yvec, type list/array of len 3
start, stop and number of nodes for linspace method, building the vectors in x and y direction.
- custEM.meshgen.meshgen_utils.translate(array, shift, back=False)[source]
Translate array about given shift
Required arguments
- array, type array with shape (:, 3)
array of coordinates to ne translated
- shift, type list or array of length 3
x-, y-, and z- values defining the shift
Keyword arguments
- back = False, type bool
set “True” is shift in opposite direction is required
- custEM.meshgen.meshgen_utils.write_synth_topo_to_asc(t_dir, file_name, topo, spacing=25.0, x_min=-10000.0, y_min=-10000.0)[source]
Write an array of topography data to an .asc file.
Required arguments
- t_dir, type str
name of output directory
- file_name, type str
name of exported topography file
- topo, type array with shape (:, 3)
2D array containing only z_coordinates with topography information
Keyword arguments:
- spacing = 25., type float
horizontal grid spacing in x- and y-direction
- x_min = -1e4, type float
x-coordinate of “lower left” corner
- y_min = -1e4, type float
y-coordiante of “lower_left” corner
- custEM.meshgen.meshgen_utils.write_synth_topo_to_xyz(t_dir, file_name, topo)[source]
Write an array of topography data to an .xyz file.
Required arguments
- t_dir, type str
name of output directory
- file_name, type str
name of exported topography file
- topo, type array with shape (:, 3)
array containing coordinates with topography information
custEM.meshgen.vtk_convert module
@author: modified by Rochlitz.R
- custEM.meshgen.vtk_convert.vtk2xml2h5(mesh_name, m_dir, rot, overwrite_markers)[source]
Converts TetGen meshes in VTK (.vtk) format to xml and h5 files for the usage in FEniCS. This function is called automatically during the intialization of the MOD class and converts a mesh if it does not already exist as xml and h5 file.
Required arguments
- mesh_name, type str
mesh name
- m_dir, type str
path to mesh directory
Functionality
Convert between .vtk and .xml, parser implemented as a state machine:
0 = read number of vertices 1 = read vertices 2 = read number of cells 3 = read cells 4 = read number of markers 5 = read cell markers 6 = done
- custEM.meshgen.vtk_convert.write_cell_tetrahedron(ofile, cell, n0, n1, n2, n3)[source]
Write tetrahedron 3D information to file.
- custEM.meshgen.vtk_convert.write_domain_marker(ofile, cell, marker, dfile=False)[source]
Write domain markers to domain file.
Write footer cells to file.
Write footer to domain file.
Write footer of mesh file.
Write footer vertices to file.
- custEM.meshgen.vtk_convert.write_header_cells(ofile, num_cells)[source]
Write header cells to file.
- custEM.meshgen.vtk_convert.write_header_domains(ofile, num_domains, dfile=False)[source]
Write header to domain file.
- custEM.meshgen.vtk_convert.write_header_mesh(ofile, cell_type, dim)[source]
Write header of mesh file.
- custEM.meshgen.vtk_convert.write_header_vertices(ofile, num_vertices)[source]
Write header vertices to file.
Module contents
meshgen
Submodules:
meshgen_tools for tetrahedral mesh generation
meshgen_utils providing utility functions for mesh generation
invmesh_tools providing specific classes for building inversion meshes
mesh_convert for automated conversion of .mesh files to .xml and .h5
dem_interpolator for interpolating real digital elevation data
bathy_tools for identifying coastline paths