1+ from typing import List , Dict , Optional
2+
3+ import numpy as np
4+
5+ from ...core .data .dual_contouring_mesh import DualContouringMesh
6+ from ...core .data .stacks_structure import StacksStructure
7+
8+
9+ def _apply_fault_relations_to_overlaps (
10+ all_meshes : List [DualContouringMesh ],
11+ voxel_overlaps : Dict [str , dict ],
12+ stacks_structure : StacksStructure
13+ ) -> None :
14+ """
15+ Apply fault relations to voxel overlaps by updating mesh vertices.
16+
17+ Args:
18+ all_meshes: List of dual contouring meshes
19+ voxel_overlaps: Dictionary containing overlap information between stacks
20+ stacks_structure: Structure containing fault relations and stack information
21+ """
22+ if stacks_structure .faults_relations is None :
23+ return
24+
25+ faults_relations = stacks_structure .faults_relations
26+ n_stacks = stacks_structure .n_stacks
27+ surfaces_per_stack = stacks_structure .number_of_surfaces_per_stack_vector
28+
29+ # Process fault relations
30+ for origin_stack , destination_stack in _get_fault_pairs (faults_relations , n_stacks ):
31+ surface_range = _get_surface_range (surfaces_per_stack , destination_stack )
32+
33+ for surface_n in surface_range :
34+ overlap_key = f"stack_{ origin_stack } { surface_n }
35+
36+ if overlap_key in voxel_overlaps :
37+ _apply_vertex_sharing (
38+ all_meshes ,
39+ origin_stack ,
40+ surface_n ,
41+ voxel_overlaps [overlap_key ]
42+ )
43+
44+
45+ def _get_fault_pairs (faults_relations : np .ndarray , n_stacks : int ):
46+ """Generate pairs of stacks that have fault relations."""
47+ for origin_stack in range (n_stacks ):
48+ for destination_stack in range (n_stacks ):
49+ if faults_relations [origin_stack , destination_stack ]:
50+ yield origin_stack , destination_stack
51+
52+
53+ def _get_surface_range (surfaces_per_stack : np .ndarray , stack_index : int ) -> range :
54+ """Get the range of surfaces for a given stack."""
55+ return range (
56+ surfaces_per_stack [stack_index ],
57+ surfaces_per_stack [stack_index + 1 ]
58+ )
59+
60+
61+ def _apply_vertex_sharing (
62+ all_meshes : List [DualContouringMesh ],
63+ origin_mesh_idx : int ,
64+ destination_mesh_idx : int ,
65+ overlap_data : dict
66+ ) -> None :
67+ """
68+ Apply vertex sharing between origin and destination meshes based on overlap data.
69+
70+ Args:
71+ all_meshes: List of dual contouring meshes
72+ origin_mesh_idx: Index of mesh that serves as the source of vertices
73+ destination_mesh_idx: Index of mesh that receives vertices from origin
74+ overlap_data: Dictionary containing indices and overlap information
75+ """
76+ if not _are_valid_mesh_indices (all_meshes , origin_mesh_idx , destination_mesh_idx ):
77+ return
78+
79+ origin_mesh = all_meshes [origin_mesh_idx ]
80+ destination_mesh = all_meshes [destination_mesh_idx ]
81+
82+ # Share vertices from origin to destination
83+ origin_indices = overlap_data ["indices_in_stack_i" ]
84+ destination_indices = overlap_data ["indices_in_stack_j" ]
85+
86+ destination_mesh .vertices [destination_indices ] = origin_mesh .vertices [origin_indices ]
87+
88+
89+ def _are_valid_mesh_indices (all_meshes : List [DualContouringMesh ], * indices : int ) -> bool :
90+ """Check if all provided mesh indices are valid."""
91+ return all (0 <= idx < len (all_meshes ) for idx in indices )
92+
93+
94+ def find_repeated_voxels_across_stacks (all_left_right_codes : List [np .ndarray ]) -> Dict [str , dict ]:
95+ """
96+ Find repeated voxels using NumPy operations for efficient processing of large arrays.
97+
98+ Args:
99+ all_left_right_codes: List of left_right_codes arrays, one per stack
100+
101+ Returns:
102+ Dictionary with detailed overlap analysis between stack pairs
103+ """
104+ if not all_left_right_codes :
105+ return {}
106+
107+ stack_codes = _generate_voxel_codes (all_left_right_codes )
108+ return _find_overlaps_between_stacks (stack_codes , all_left_right_codes )
109+
110+
111+ def _generate_voxel_codes (all_left_right_codes : List [np .ndarray ]) -> List [np .ndarray ]:
112+ """Generate voxel codes for each stack using packed bit directions."""
113+ from gempy_engine .modules .dual_contouring .fancy_triangulation import _StaticTriangulationData
114+
115+ pack_directions = _StaticTriangulationData .get_pack_directions_into_bits ()
116+ stack_codes = []
117+
118+ for left_right_codes in all_left_right_codes :
119+ if left_right_codes .size > 0 :
120+ voxel_codes = (left_right_codes * pack_directions ).sum (axis = 1 )
121+ stack_codes .append (voxel_codes )
122+ else :
123+ stack_codes .append (np .array ([]))
124+
125+ return stack_codes
126+
127+
128+ def _find_overlaps_between_stacks (
129+ stack_codes : List [np .ndarray ],
130+ all_left_right_codes : List [np .ndarray ]
131+ ) -> Dict [str , dict ]:
132+ """Find overlaps between all pairs of stacks."""
133+ overlaps = {}
134+
135+ for i in range (len (stack_codes )):
136+ for j in range (i + 1 , len (stack_codes )):
137+ overlap_data = _process_stack_pair (
138+ stack_codes [i ], stack_codes [j ],
139+ all_left_right_codes [i ], all_left_right_codes [j ],
140+ i , j
141+ )
142+
143+ if overlap_data :
144+ overlaps [f"stack_{ i } { j } ] = overlap_data
145+
146+ return overlaps
147+
148+
149+ def _process_stack_pair (
150+ codes_i : np .ndarray ,
151+ codes_j : np .ndarray ,
152+ left_right_i : np .ndarray ,
153+ left_right_j : np .ndarray ,
154+ stack_i : int ,
155+ stack_j : int
156+ ) -> Optional [dict ]:
157+ """Process a pair of stacks to find overlapping voxels."""
158+ if codes_i .size == 0 or codes_j .size == 0 :
159+ return None
160+
161+ common_codes = np .intersect1d (codes_i , codes_j )
162+
163+ if len (common_codes ) == 0 :
164+ return None
165+
166+ # Find indices of common voxels in each stack
167+ indices_i = np .isin (codes_i , common_codes )
168+ indices_j = np .isin (codes_j , common_codes )
169+
170+ return {
171+ 'common_voxel_codes' : common_codes ,
172+ 'count' : len (common_codes ),
173+ 'indices_in_stack_i' : np .where (indices_i )[0 ],
174+ 'indices_in_stack_j' : np .where (indices_j )[0 ],
175+ 'common_binary_codes_i' : left_right_i [indices_i ],
176+ 'common_binary_codes_j' : left_right_j [indices_j ]
177+ }
0 commit comments