SciPy

pynibs.mesh package

Submodules

pynibs.mesh.mesh_struct module

class pynibs.mesh.mesh_struct.Mesh(mesh_name, subject_id, subject_folder)

Bases: object

” Mesh class to initialize default attributes.

fill_defaults(approach)

Initializes attributes for a headreco mesh.

Parameters:

approach (str) – ‘headreco’ ‘mri2mesh’ ‘charm’

print()

Print self information.

write_to_hdf5(fn_hdf5, check_file_exist=False, verbose=False)

Write this mesh’ attributes to .hdf5 file.

Parameters:

  • fn_hdf5 (str) –

  • check_file_exist (bool) – Check if provided filenames exist, warn if not.

  • verbose (bool) – Print self information

class pynibs.mesh.mesh_struct.ROI(subject_id, roi_name, mesh_name)

Bases: object

” Region of interest class to initialize default attributes.

print()

Print self information.

write_to_hdf5(fn_hdf5, check_file_exist=False, verbose=False)

Write this mesh’ attributes to .hdf5 file.

Parameters:

  • fn_hdf5 (str) –

  • check_file_exist (bool) – Check if provided filenames exist, warn if not.

  • verbose (bool) – Print self information

class pynibs.mesh.mesh_struct.TetrahedraLinear(points, triangles, triangles_regions, tetrahedra, tetrahedra_regions)

Bases: object

Mesh, consisting of linear tetrahedra.

Parameters:

  • points (array of float [N_points x 3]) – Vertices of FE mesh

  • triangles (np.ndarray of int [N_tri x 3]) – Connectivity of points forming triangles

  • triangles_regions (np.ndarray of int [N_tri x 1]) – Region identifiers of triangles

  • tetrahedra (np.ndarray of int [N_tet x 4]) – Connectivity of points forming tetrahedra

  • tetrahedra_regions (np.ndarray of int [N_tet x 1]) – Region identifiers of tetrahedra

N_points

Number of vertices

Type:

int

N_tet

Number of tetrahedra

Type:

int

N_tri

Number of triangles

Type:

int

N_region

Number of regions

Type:

int

region

Region labels

Type:

np.ndarray of int

tetrahedra_volume

Volumes of tetrahedra

Type:

np.ndarray of float [N_tet x 1]

tetrahedra_center

Center of tetrahedra

Type:

np.ndarray of float [N_tet x 1]

triangles_center

Center of triangles

Type:

np.ndarray of float [N_tri x 1]

triangles_normal

Normal components of triangles pointing outwards

Type:

np.ndarray of float [N_tri x 3]

calc_E(grad_phi, omegaA)

Calculate electric field with gradient of electric potential and omega-scaled magnetic vector potential A.

\mathbf{E}=-\nabla\varphi-\omega\mathbf{A}

Parameters:

  • grad_phi (np.ndarray of float [N_tet x 3]) – Gradient of Scalar DOF in tetrahedra center

  • omegaA (np.ndarray of float [N_tet x 3]) – Magnetic vector potential in tetrahedra center (scaled with angular frequency omega)

Returns:

E – Electric field in tetrahedra center

Return type:

np.ndarray of float [N_tet x 3]

calc_E_normal_tangential_surface(E, fname)

Calculate normal and tangential component of electric field on given surfaces of mesh instance.

Parameters:

  • E (np.ndarray of float [N_tri x 3]) – Electric field data on surfaces

  • fname (str) – Filename of the .txt file containing the tetrahedra indices, which are adjacent to the surface triangles generated by the method “calc_surface_adjacent_tetrahedra_idx_list(self, fname)”

Returns:

  • En_pos (np.ndarray of float [N_tri x 3]) – Normal component of electric field of top side (outside) of surface

  • En_neg (np.ndarray of float [N_tri x 3]) – Normal component of electric field of bottom side (inside) of surface

  • n (np.ndarray of float [N_tri x 3]) – Normal vector

  • Et (np.ndarray of float [N_tri x 3]) – Tangential component of electric field lying in surface

  • t (np.ndarray of float [N_tri x 3]) – Tangential vector

calc_E_on_GM_WM_surface(E, roi)

Determines the normal and tangential component of the induced electric field on a GM-WM surface using nearest neighbour principle.

Parameters:

  • E (np.ndarray of float [N_tri x 3]) – Induced electric field given in the tetrahedra centre of the mesh instance

  • roi (pynibs.roi.RegionOfInterestSurface) – RegionOfInterestSurface object class instance

Returns:

  • E_normal (np.ndarray of float [N_points x 3]) – Normal vector of electric field on GM-WM surface

  • E_tangential (np.ndarray of float [N_points x 3]) – Tangential vector of electric field on GM-WM surface

calc_E_on_GM_WM_surface3(phi, dAdt, roi, verbose=True, mode='components')

Determines the normal and tangential component of the induced electric field on a GM-WM surface by recalculating phi and dA/dt in an epsilon environment around the GM/WM surface (upper and lower GM-WM surface).

Parameters:

  • phi (np.ndarray of float) – (N_nodes, 1) Scalar electric potential given in the nodes of the mesh.

  • dAdt (np.ndarray of float) – (N_nodes, 3) Magnetic vector potential given in the nodes of the mesh.

  • roi (pynibs.mesh.mesh_struct.ORI) – RegionOfInterestSurface object class instance.

  • verbose (bool) – Print information to stdout.

  • mode (str) – Select mode of output: - “components” : return x, y, and z component of tangential and normal components - “magnitude” : return magnitude of tangential and normal component (normal with sign for direction)

Returns:

  • E_normal (np.ndarray of float) – (N_nodes, 3) Normal vector of electric field on GM-WM surface.

  • E_tangential (np.ndarray of float) – (N_nodes, 3) Tangential vector of electric field on GM-WM surface.

calc_E_on_GM_WM_surface_simnibs(phi, dAdt, roi, subject, verbose=False, mesh_idx=0)

Determines the normal and tangential component of the induced electric field on a GM-WM surface by recalculating phi and dA/dt in an epsilon environment around the GM/WM surface (upper and lower GM-WM surface) or by using the Simnibs interpolation function.

Parameters:

  • phi (np.ndarray of float) – (N_nodes, 1) Scalar electric potential given in the nodes of the mesh.

  • dAdt (np.ndarray of float) – (N_nodes, 3) Magnetic vector potential given in the nodes of the mesh.

  • roi (pynibs.mesh.mesh_struct.ROI) – RegionOfInterestSurface object class instance.

  • subject (pynibs.subject.Subject) – Subject object loaded from .hdf5 file.

  • verbose (bool) – Print information to stdout.

  • mesh_idx (int) – Mesh index.

Returns:

  • E_normal (np.ndarray of float) – (N_points, 3) Normal vector of electric field on GM-WM surface.

  • E_tangential (np.ndarray of float) – (N_points, 3) Tangential vector of electric field on GM-WM surface.

calc_E_on_GM_WM_surface_simnibs_KW(phi, dAdt, roi, subject, verbose=False, mesh_idx=0)

Determines the normal and tangential component of the induced electric field on a GM-WM surface by recalculating phi and dA/dt in an epsilon environment around the GM/WM surface (upper and lower GM-WM surface) or by using the Simnibs interpolation function.

Parameters:

  • phi (np.ndarray of float) – (N_nodes, 1) Scalar electric potential given in the nodes of the mesh.

  • dAdt (np.ndarray of float) – (N_nodes, 1) Magnetic vector potential given in the nodes of the mesh.

  • roi (pynibs.mesh.mesh_struct.ROI) – RegionOfInterestSurface object class instance.

  • subject (pynibs.subject.Subject) – Subject object loaded from .hdf5 file.

  • verbose (bool) – Print information to stdout.

  • mesh_idx (int) – Mesh index.

Returns:

  • E_normal (np.ndarray of float) – (N_points, 3) Normal vector of electric field on GM-WM surface.

  • E_tangential (np.ndarray of float) – (N_points, 3) Tangential vector of electric field on GM-WM surface.

calc_J(E, sigma)

Calculate current density J. The conductivity sigma is a list of np.arrays containing conductivities of regions (scalar and/or tensor).

\mathbf{J} = [\sigma]\mathbf{E}

Parameters:

  • E (np.ndarray of float [N_tet x 3]) – Electric field in tetrahedra center

  • sigma (list of np.ndarray of float [N_regions][3 x 3]) – Conductivities of regions (scalar and/or tensor).

Returns:

E – Electric field in tetrahedra center

Return type:

np.ndarray of float [N_tet x 3]

calc_QOI_in_points(qoi, points_out)

Calculate QOI_out in points_out using the mesh instance and the quantity of interest (QOI).

Parameters:

  • qoi (np.ndarray of float) – Quantity of interest in nodes of tetrahedra mesh instance

  • points_out (np.ndarray of float) – Point coordinates (x, y, z) where the qoi is going to be interpolated by linear basis functions

Returns:

qoi_out – Quantity of interest in points_out

Return type:

np.ndarray of float

calc_QOI_in_points_tet_idx(qoi, points_out, tet_idx)

Calculate QOI_out in points_out sitting in tet_idx using the mesh instance and the quantity of interest (QOI).

Parameters:

  • qoi (np.ndarray of float) – Quantity of interest in nodes of tetrahedra mesh instance

  • points_out (np.ndarray of float) – Point coordinates (x, y, z) where the qoi is going to be interpolated by linear basis functions

  • tet_idx (np.ndarray of int) – Element indices where the points_out are sitting

Returns:

qoi_out – Quantity of interest in points_out

Return type:

np.ndarray of float

calc_gradient(phi)

Calculate gradient of scalar DOF in tetrahedra center.

Parameters:

phi (np.ndarray of float [N_nodes]) – Scalar DOF the gradient is calculated for

Returns:

grad_phi – Gradient of Scalar DOF in tetrahedra center

Return type:

np.ndarray of float [N_tet x 3]

calc_surface_adjacent_tetrahedra_idx_list(fname)

Determine the indices of the tetrahedra touching the surfaces and save the indices into a .txt file specified with fname.

Parameters:

fname (str) – Filename of output .txt file

Returns:

<File> – Element indices of the tetrahedra touching the surfaces (outer-most elements)

Return type:

.txt file

data_elements2nodes(data)

Transforms an data in tetrahedra into the nodes after Zienkiewicz et al. (1992) [1]. Can only transform volume data, i.e. needs the data in the surrounding tetrahedra to average it to the nodes. Will not work well for discontinuous fields (like E, if several tissues are used).

Parameters:

data (np.ndarray [N_elements x N_data]) – Data in tetrahedra

Returns:

data_nodes – Data in nodes

Return type:

np.ndarray [N_nodes x N_data]

Notes

[1]

Zienkiewicz, Olgierd Cecil, and Jian Zhong Zhu. “The superconvergent patch recovery and a posteriori error estimates. Part 1: The recovery technique.” International Journal for Numerical Methods in Engineering 33.7 (1992): 1331-1364.

data_nodes2elements(data)

Interpolate data given in the nodes to the tetrahedra center.

Parameters:

data (np.ndarray [N_nodes x N_data]) – Data in nodes

Returns:

data_elements – Data in elements

Return type:

np.ndarray [N_elements x N_data]

get_faces(tetrahedra_indexes=None)

Creates a list of nodes in each face and a list of faces in each tetrahedra.

Parameters:

tetrahedra_indexes (np.ndarray) – Indices of the tetrehedra where the faces are to be determined (default: all tetrahedra)

Returns:

  • faces (np.ndarray) – List of nodes in faces, in arbitrary order

  • th_faces (np.ndarray) – List of faces in each tetrahedra, starts at 0, order=((0, 2, 1), (0, 1, 3), (0, 3, 2), (1, 2, 3))

  • face_adjacency_list (np.ndarray) – List of tetrahedron adjacent to each face, filled with -1 if a face is in a single tetrahedron. Not in the normal element ordering, but only in the order the tetrahedra are presented

get_outside_faces(tetrahedra_indexes=None)

Creates a list of nodes in each face that are in the outer volume.

Parameters:

tetrahedra_indices (np.ndarray) – Indices of the tetrehedra where the outer volume is to be determined (default: all tetrahedra)

Returns:

faces – List of nodes in faces in arbitrary order

Return type:

np.ndarray

pynibs.mesh.transformations module

pynibs.mesh.transformations.cell_data_to_point_data(tris, data_tris, nodes, method='nearest')

A wrapper for scipy.interpolate.griddata to interpolate cell data to node data.

Parameters:

  • tris (np.ndarray) – (n_tri, 3) element number list.

  • data_tris (np.ndarray) – (n_tri x 3) data in tris.

  • nodes (np.ndarray) – (n_nodes, 3) nodes coordinates.

  • method (str, default: 'nearest') – Which method to use for interpolation. Default uses NearestNDInterpolator.

Returns:

data_nodes – Data in nodes

Return type:

np.ndarray

pynibs.mesh.transformations.data_elements2nodes(data, con, precise=False)

Transforms data in elements (triangles or tetrahedra) to nodes. Data can be list of multiple data arrays.

Parameters:

  • data (np.ndarray of float or list of np.ndarray) – (N_elements, N_data) Data given in the elements (multiple datasets who fit to con may be passed in a list).

  • con (np.ndarray of int) – triangles: (N_elements. 3). tetrahedra: (N_elements, 4). Connectivity index list forming the elements.

  • precise (bool, default: False) – Compute data transformation precisely but slow. Better for near-0 values..

Returns:

out – (N_nodes, N_data) Data in nodes.

Return type:

np.ndarray of float or list of np.ndarray

pynibs.mesh.transformations.data_nodes2elements(data, con)

Transforms data in nodes to elements (triangles or tetrahedra).

Parameters:

  • data (np.ndarray of float) – (N_nodes, N_data) Data given in the nodes.

  • con (np.ndarray of int) – triangles: (N_elements, 3). tetrahedra: (N_elements, 4). Connectivity index list forming the elements.

Returns:

out – (N_elements, N_data) Data given in the element scenters.

Return type:

np.ndarray of float

pynibs.mesh.transformations.map_data_to_surface(datasets, points_datasets, con_datasets, fname_fsl_gm, fname_fsl_wm, fname_midlayer=None, delta=0.5, input_data_in_center=True, return_data_in_center=True, data_substitute=-1)

Maps data from ROI of fsl surface (wm, gm, or midlayer) to given Freesurfer brain surface (wm, gm, inflated).

Parameters:

  • datasets (np.ndarray of float [N_points x N_data] or list of np.ndarray) – Data in nodes or center of triangles in ROI (specify this in “data_in_center”)

  • points_datasets (np.ndarray of float [N_points x 3] or list of np.ndarray) – Point coordinates (x,y,z) of ROI where data in datasets list is given, the points have to be a subset of the GM/WM surface (has to be provided for each dataset)

  • con_datasets (np.ndarray of int [N_tri x 3] or list of np.ndarray) – Connectivity matrix of dataset points (has to be provided for each dataset)

  • fname_fsl_gm (str or list of str or list of None) – Filename of pial surface fsl file(s) (one or two hemispheres) e.g. in mri2msh: …/fs_ID/surf/lh.pial

  • fname_fsl_wm (str or list of str or list of None) – Filename of wm surface fsl file(s) (one or two hemispheres) e.g. in mri2msh: …/fs_ID/surf/lh.white

  • fname_midlayer (str or list of str) – Filename of midlayer surface fsl file(s) (one or two hemispheres) e.g. in headreco: …/fs_ID/surf/lh.central

  • delta (float) – Distance parameter where gm-wm surface was generated 0…1 (default: 0.5) 0 -> WM surface 1 -> GM surface

  • input_data_in_center (bool) – Flag if data in datasets in given in triangle centers or in points (Default: True)

  • return_data_in_center (bool) – Flag if data should be returned in nodes or in elements (Default: True)

  • data_substitute (float) – Data substitute with this number for all points in the inflated brain, which do not belong to the given data set

Returns:

data_mapped – Mapped data to target brain surface. In points or elements

Return type:

np.ndarray of float [N_points_inf x N_data]

pynibs.mesh.transformations.midlayer_2_surf(midlayer_data, coords_target, coords_midlayer, midlayer_con=None, midlayer_data_in_nodes=False, max_dist=5, outside_roi_val=0, precise_map=True)

Convert midlayer data to whole-brain surface data, e.g. grey matter. Output is returned as data in nodes.

Parameters:

  • midlayer_data (np.ndarray of float) – (n_elm_midlayer,) or (n_nodes_midlayer,), the data in the midlayer.

  • coords_target (np.ndarray of float) – (n_nodes_target, 3) Coordinates of the nodes of the target surface.

  • coords_midlayer (np.ndarray of float) – (n_nodes_midlayer, 3) Coordinates of the nodes of the midlayer surface.

  • midlayer_con (np.ndarray of int, optional) – (n_elm_midlayer, 3) Connectivity of the midlayer elements. Provide if data_in_points == True.

  • midlayer_data_in_nodes (bool, default=False) – If midlayer data is provided in nodes, set to True and provide midlayer_con.

  • max_dist (float, default=5) – Maximum distance between target and midlayer nodes to pull data from midlayer_data for.

  • outside_roi_val (float, default=0) – Areas outside of max_dist are filled with outside_roi_val.

  • precise_map (bool, default=True) – If elements to nodes mapping is done, perform this precise and slow or not.

Returns:

data_target – (n_nodes_target, 1) The data in nodes of the target surface.

Return type:

np.ndarray

pynibs.mesh.transformations.project_on_scalp(coords, mesh, scalp_tag=1005)

Find the node in the scalp closest to each coordinate

Parameters:

  • coords (nx3 np.ndarray) – Vectors to be transformed

  • mesh (pynibs.TetrahedraLinear or simnibs.msh.mesh_io.Msh) – Mesh structure in simnibs or pynibs format

  • scalp_tag (int, default: 1005) – Tag in the mesh where the scalp is to be set.

Returns:

points_closest – (n, 3) coordinates projected scalp (closest skin points)

Return type:

np.ndarry

pynibs.mesh.transformations.project_on_scalp_hdf5(coords, mesh, scalp_tag=1005)

Find the node in the scalp closest to each coordinate.

Parameters:

  • coords (np.ndarray) – (n, 3) Vectors to be transformed.

  • mesh (str or pynibs.TetrahedraLinear) – Filename of mesh in .hdf5 format or Mesh structure.

  • scalp_tag (int, default: 1005) – Tag in the mesh where the scalp is to be set.

Returns:

points_closest – (n, 3) coordinates projected scalp (closest skin points).

Return type:

np.ndarray

pynibs.mesh.transformations.refine_surface(fn_surf, fn_surf_refined, center, radius, repair=True, remesh=True, verbose=True)

Refines surface (.stl) in spherical ROI an saves as .stl file.

Parameters:

  • fn_surf (str) – Input filename (.stl).

  • fn_surf_refined (str) – Output filename (.stl).

  • center (np.ndarray of float) –

    1. Center of spherical ROI (x,y,z).

  • radius (float) – Radius of ROI.

  • repair (bool, default: True) – Repair surface mesh to ensure that it is watertight and forms a volume.

  • remesh (bool, default: False) – Perform remeshing with meshfix (also removes possibly overlapping facets and intersections).

  • verbose (bool, default: True) – Print output messages.

Returns:

<file>

Return type:

.stl file

pynibs.mesh.utils module

pynibs.mesh.utils.calc_distances(coords, mesh_fn, tissues=None)

Calculates the distances between coords and tissue types.

Parameters:

  • coords (list of list or list or np.ndarray) – Coordinates (X, Y, Z) to compute depths for.

  • mesh_fn (str) – pynibs.Mesh hdf5 filename.

  • tissues (list of int, optional) – Which tissue types to compute depths for. If none, distances to all tissue types are computed.

Returns:

distances – colunms: coorrd, tissue_type, distance

Return type:

pd.Dataframe()

pynibs.mesh.utils.calc_gradient_surface(phi, points, triangles)

Calculate gradient of potential phi on surface (i.e. tangential component) given in vertices of a triangular mesh forming a 2D surface.

Parameters:

  • phi (np.ndarray of float [N_points x 1]) – Potential in nodes

  • points (np.ndarray of float [N_points x 3]) – Coordinates of nodes (x,y,z)

  • triangles (np.ndarray of int32 [N_tri x 3]) – Connectivity of triangular mesh

Returns:

grad_phi – Gradient of potential phi on surface

Return type:

np.ndarray of float [N_tri x 3]

pynibs.mesh.utils.calc_tet_volume(points, abs=True)

Calculate tetrahedra volumes.

Parameters:

  • points (np.ndarray) – shape: (n_tets,4,3)

  • abs (bool, default: true) – Return magnitude

Returns:

volume – shape: (n_tets)

Return type:

np.ndarray

pynibs.mesh.utils.calc_tetrahedra_volume_cross(P1, P2, P3, P4)

Calculates volume of tetrahedra specified by the 4 points P1…P4 multiple tetrahedra can be defined by P1…P4 as 2-D np.ndarrays using the cross and vector dot product.

P1=\begin{bmatrix} x_{{tet}_1} & y_{{tet}_1} & z_{{tet}_1} \\ x_{{tet}_2} & y_{{tet}_2} & z_{{tet}_2} \\ ... & ... & ... \\ x_{{tet}_N} & y_{{tet}_N} & z_{{tet}_N} \\ \end{bmatrix}

Parameters:

  • P1 (np.ndarray of float [N_tet x 3]) – Coordinates of first point of tetrahedra

  • P2 (np.ndarray of float [N_tet x 3]) – Coordinates of second point of tetrahedra

  • P3 (np.ndarray of float [N_tet x 3]) – Coordinates of third point of tetrahedra

  • P4 (np.ndarray of float [N_tet x 3]) – Coordinates of fourth point of tetrahedra

Returns:

tetrahedra_volume – Volumes of tetrahedra

Return type:

np.ndarray of float [N_tet x 1]

pynibs.mesh.utils.calc_tetrahedra_volume_det(P1, P2, P3, P4)

Calculate volume of tetrahedron specified by 4 points P1…P4 multiple tetrahedra can be defined by P1…P4 as 2-D np.arrays using the determinant.

P1=\begin{bmatrix} x_{{tet}_1} & y_{{tet}_1} & z_{{tet}_1} \\ x_{{tet}_2} & y_{{tet}_2} & z_{{tet}_2} \\ ... & ... & ... \\ x_{{tet}_N} & y_{{tet}_N} & z_{{tet}_N} \\ \end{bmatrix}

Parameters:

  • P1 (np.ndarray of float [N_tet x 3]) – Coordinates of first point of tetrahedra

  • P2 (np.ndarray of float [N_tet x 3]) – Coordinates of second point of tetrahedra

  • P3 (np.ndarray of float [N_tet x 3]) – Coordinates of third point of tetrahedra

  • P4 (np.ndarray of float [N_tet x 3]) – Coordinates of fourth point of tetrahedra

Returns:

tetrahedra_volume – Volumes of tetrahedra

Return type:

np.ndarray of float [N_tet x 1]

pynibs.mesh.utils.calc_tri_surface(points)

Calculate triangle surface areas.

Parameters:

points (np.ndarray) – (n_triangles,3,3)

Returns:

triangle_area

Return type:

np.ndarray

pynibs.mesh.utils.check_islands_for_single_elm(source_elm, connectivity=None, adjacency=None, island_crit=1)

This identifies islands in a mesh for a given element. An island is a set of elements, that is only connect via a single node to another set of elements. These islands usually crash the FEM solver and should be removed.

  1. Find all elements connect to source_elm via one node (1-node-neighbor)

  2. Start with source_elm and visit all 2-node-neighbors (‘shared-edge)

  3. Continue recursively with all 2-node-neighbors and visit their 2-node-neighbors

  4. See if any 1-node-neighbors have not been visited with this strategy. If so, an island has been found

Parameters:

  • source_elm (int) – The source element to check

  • connectivity (np.ndarray, optional) – Connectivity (‘node_number_list’) starting with 0. Can be triangles or tetrahedra (n_elms, 3) or (n_elms_4).

  • adjacency (np.ndparray, optional) – Adjenceny matrix (n_elm, n_elm). Weights are supposed to be number of shared nodes. Computed from neighbors if not provided.

  • island_crit (int, default: 'any') – How many nodes to define islands? ‘any’ -> Elements connected via a single node or single edge are defined as an island. ‘node’ -> Elements connected via a single _node_ are defined as an island. ‘edge’ -> Elements connected via a single _edge_ are defined as an island.

Returns:

  • n_visited (int)

  • n_not_visited (int)

  • neighbors_visited (dict, which neighbors have been visited and which have not)

pynibs.mesh.utils.cortical_depth(mesh_fn, geo_fn=None, write_xdmf=True, skin_surface_id=1005, verbose=False)

Compute skin-cortex-distance (SCD) for surface and volume data in mesh_fn.

Visualized cortical depth.

Cortical depth computed against skin surface.

Parameters:

  • mesh_fn (str) – TetrahedraLinear mesh file.

  • geo_fn (str, optional) – TetrahedraLinear mesh file with geometric data. If provided, geometric information is read from here.

  • write_xdmf (bool, default: True) – Write .xdmf or not.

  • skin_surface_id (int, default: 1005) – Which tissue type nr to compute distance against.

  • verbose (bool, default: False) – Print some verbosity information.

Returns:

  • <file> (.hdf5) – mesh_fn or geo_fn with SCD information in /data/tris/Cortex_dist and /data/tets/Cortex_dist.

  • <file> (.xdmf) – Only if write_xdmf == True.

pynibs.mesh.utils.determine_e_midlayer(fn_e_results, fn_mesh_hdf5, subject, mesh_idx, roi_idx, n_cpu=4, midlayer_fun='simnibs', phi_scaling=1.0, verbose=False)

Parallel version to determine the midlayer e-fields from a list of .hdf5 results files

Parameters:

  • fn_e_results (list of str) – List of results filenames (.hdf5 format)

  • fn_mesh_hdf5 (str) – Filename of corresponding mesh file

  • subject (pynibs.Subject) – Subject object

  • mesh_idx (int) – Mesh index

  • roi_idx (int) – ROI index

  • n_cpu (int, default: 4) – Number of parallel computations

  • midlayer_fun (str, default: "simnibs") – Method to determine the midlayer e-fields (“pynibs” or “simnibs”)

  • phi_scaling (float, default: 1.0) – Scaling factor of scalar potential to change between “m” and “mm”

Returns:

Adds midlayer e-field results to ROI

Return type:

<File> .hdf5 file

pynibs.mesh.utils.determine_e_midlayer_workhorse(fn_e_results, subject, mesh_idx, midlayer_fun, fn_mesh_hdf5, roi_idx, phi_scaling=1.0, verbose=False)
phi_scaling: float

simnibs < 3.0 : 1000. simnibs >= 3.0 : 1. (Default)

pynibs.mesh.utils.find_element_idx_by_points(nodes, con, points)

Finds the tetrahedral element index of an arbitrary point in the FEM mesh.

Parameters:

  • nodes (np.ndarray [N_nodes x 3]) – Coordinates (x, y, z) of the nodes

  • con (np.ndarray [N_tet x 4]) – Connectivity matrix

  • points (np.ndarray [N_points x 3]) – Points for which the element indices are found.

Returns:

ele_idx – Element indices of tetrahedra where corresponding ‘points’ are lying in

Return type:

np.ndarray [N_points]

pynibs.mesh.utils.find_island_elms(connectivity=None, adjacency=None, verbose=False, island_crit='edge', decision='cumulative')

Searches for islands in a mesh and returns element indices of the smallest island. Island is defines as a set of elements, which share a single node and/or single edge with the rest of the mesh.

Parameters:

  • connectivity (np.ndarray, optional) – Connectivity (‘node_number_list’) starting with 0. Can be triangles or tetrahedra (n_elms, 3) or (n_elms_4).

  • adjacency (np.ndparray, optional) – Adjenceny matrix (n_elm, n_elm). Weights are supposed to be number of shared nodes. Computed from neighbors if not provided.

  • island_crit (int, default: 'edge') – How many nodes to define islands? ‘node’ -> Elements connected via a single _node_ are defined as an island. ‘edge’ -> Elements connected via a single _edge_ are defined as an island.

  • decision (str, default: cumulative) – ‘cumulative’ -> Return all element indices that are not visited any times ‘smallest’ -> Return smallest island.

  • verbose (bool, optional) – Print some verbosity information. Default: False

Returns:

island

Return type:

list of island-elms

pynibs.mesh.utils.find_islands(connectivity=None, adjacency=None, island_crit='any', verbose=False, largest=False)

This identifies islands in a mesh. An island is a set of elements, that is only connect via a single node to another set of elements. These islands usually crash the FEM solver and should be removed.

For each element:

  1. Find all elements connect to source_elm via one node (1-node-neighbor)

  2. Start with source_elm and visit all 2-node-neighbors (‘shared-edge)

  3. Continue recursively with all 2-node-neighbors and visit their 2-node-neighbors

  4. See if any 1-node-neighbors have not been visited with this strategy. If so, an island has been found

Island detection

Islands are groups of elements that are only connected via a single node/edge to another group.

Parameters:

  • connectivity (np.ndarray, optional) – Connectivity (‘node_number_list’) starting with 0. Can be triangles or tetrahedra (n_elms, 3) or (n_elms_4).

  • adjacency (np.ndparray, optional) – Adjenceny matrix (n_elm, n_elm). Weights are supposed to be number of shared nodes. Computed from neighbors if not provided.

  • island_crit (int or str, default: 'any') – How many nodes to define islands? ‘any’ -> Elements connected via a single node or single edge are defined as an island. ‘node’ -> Elements connected via a single _node_ are defined as an island. ‘edge’ -> Elements connected via a single _edge_ are defined as an island.

  • largest (book, default: False) – Only return largest island, speeds up computation quite a bit if only one large, and many small islands exist.

  • verbose (bool, optional) – Print some verbosity information. Default: False

Returns:

  • elms_with_island (list) – Elements with neighboring islands

  • counter_visited (np.ndarray) – shape = (n_elms). How often as each element been visited.

  • counter_not_visited (np.ndarray) – shape = (n_elms). How often as each element not been visited.

pynibs.mesh.utils.find_nearest(array, value)

Given an “array”, and given a “value” , returns an index j such that “value” is between array[j] and array[j+1]. “array” must be monotonic increasing. j=-1 or j=len(array) is returned to indicate that “value” is out of range below and above respectively.

Parameters:

  • array (np.ndarray of float) – Monotonic increasing array.

  • value (float) – Target value the nearest neighbor index in array is computed for.

Returns:

idx – Index j such that “value” is between array[j] and array[j+1].

Return type:

int

pynibs.mesh.utils.get_indices_discontinuous_data(data, con, neighbor=False, deviation_factor=2, min_val=None, not_fitted_elms=None, crit='median', neigh_style='point')

Get element indices (and the best neighbor index), where the data is discontinuous

Parameters:

  • data (np.ndarray of float [n_data]) – Data array to analyze given in the element center

  • con (np.ndarray of float [n_data, 3 or 4]) – Connectivity matrix

  • neighbor (bool, default: False) – Return also the element index of the “best” neighbor (w.r.t. median of data)

  • deviation_factor (float) – Allows data deviation from 1/deviation_factor < data[i]/median < deviation_factor

  • min_val (float, optional) – If given, only return elements which have a neighbor with data higher than min_val.

  • not_fitted_elms (np.ndarray) – If given, these elements are not used as neighbors

  • crit (str, default: median) – Criterium for best neighbor. Either median or max value

  • neigh_style (str, default: 'point') – Should neighbors share point or ‘edge’

Returns:

  • idx_disc (list of int [n_disc]) – Index list containing the indices of the discontinuous elements

  • idx_neighbor (list of int [n_disc]) – Index list containing the indices of the “best” neighbors of the discontinuous elements

pynibs.mesh.utils.get_sphere(mesh=None, mesh_fn=None, target=None, radius=None, roi_idx=None, roi=None, elmtype='tris', domain=None)

Return element idx of elements within a certain distance to provided target. Element indices are 0-based (tris and tets start at 0, ‘pynibs’ style) Elements might be ‘tris’ (default) or ‘tets’

If roi object / idx and mesh fn is provided, the roi is expected to have midlayer information and the roi geometry is used.

Parameters:

  • mesh (pynibs.mesh.mesh_struct.TetrahedraLinear, optional) –

  • mesh_fn (str, optional) – Filename to SimNIBS .msh or pyNIBS .hdf5 mesh file.

  • target (np.ndarray of float or list of float) – (3,) X, Y, Z coordinates of target.

  • radius (float) – Sphere radius im mm.

  • roi_idx (str or int, optional) – ROI name.

  • elmtype (str, default: 'tris') – Return triangles or tetrahedra in sphere around target. One of (‘tris’, ‘tets’).

Returns:

elms_in_sphere – (n_elements): Indices of elements found in ROI

Return type:

np.ndarray

pynibs.mesh.utils.in_hull(points, hull)

Test if points in points are in hull. points should be a [N x K] coordinates of N points in K dimensions. hull is either a scipy.spatial.Delaunay object or the [M x K] array of the coordinates of M points in Kdimensions for which Delaunay triangulation will be computed.

Parameters:

  • points (np.ndarray) – (N_points x 3) Set of floating point data to test whether they are lying inside the hull or not.

  • hull (scipy.spatial.Delaunay or np.ndarray) – (M x K) Surface data.

Returns:

inside – TRUE: point inside the hull FALSE: point outside the hull

Return type:

np.ndarray of bool

pynibs.mesh.utils.sample_sphere(n_points, r)

Creates n_points evenly spread in a sphere of radius r.

Parameters:

  • n_points (int) – Number of points to be spread, must be odd.

  • r (float) – Radius of sphere.

Returns:

points – (N x 3), Evenly spread points in a unit sphere.

Return type:

np.ndarray of float

pynibs.mesh.utils.tets_in_sphere(mesh, target, radius, roi, domain=None)

Worker function for get_sphere()

Returns element idx of elements within a certain distance to provided target. If roi object / idx and mesh fn is provided, the roi is expected to have midlayer information and the roi geometry is used.

If radius is None or 0, the nearest element is returned.

Parameters:

  • mesh (pynibs.TetrahedraLinear, optional) –

  • target (np.ndarray of float, optional) – (3,) X, Y, Z coordinates of target

  • radius (float, optional) – Sphere radius im mm

  • roi (pynibs.mesh.ROI, optional) – Region of interest

Returns:

tets_in_sphere – (n_tets): Indices of elements found in ROI

Return type:

np.ndarray

pynibs.mesh.utils.tris_in_sphere(mesh, target, radius, roi)

Worker function for get_sphere().

Returns triangle idx of elements within a certain distance to provided target. If roi object / idx and mesh fn is provided, the roi is expected to have midlayer information and the roi geometry is used.

If radius is None or 0, the nearest element is returned.

Parameters:
Returns:

tris_in_sphere – (n_triangles): Indices of elements found in sphere

Return type:

np.ndarray

Table of Contents

Previous topic

pynibs.expio package

This Page

Edit page
  • © Copyright 2023, Konstantin Weise, Ole Numssen.
  • Created using Sphinx 6.2.1.