@@ -862,25 +862,6 @@ def compute(
862862 )
863863 return units
864864
865- unit_column_candidates = [
866- "UNITNAME" ,
867- "unitname" ,
868- "UNIT_NAME" ,
869- "UnitName" ,
870- "unit" ,
871- "Unit" ,
872- "UNIT" ,
873- "name" ,
874- "Name" ,
875- ]
876- unit_column = next ((col for col in unit_column_candidates if col in geology_data .columns ), None )
877- if unit_column is None :
878- logger .warning (
879- "AlongSection: Unable to identify a unit-name column in geology data; expected one of %s." ,
880- unit_column_candidates ,
881- )
882- return units
883-
884865 sections = self .sections .copy ()
885866 geology = geology_data .copy ()
886867
@@ -898,6 +879,44 @@ def compute(
898879 logger .error ("Failed to create spatial index for geology data: %s" , e , exc_info = True )
899880 geology_sindex = None
900881
882+ #TODO: check if sections and geology have same crs and reproject if needed. then check if sections overlap geology
883+ sections_bounds = sections .total_bounds
884+ geology_bounds = geology .total_bounds
885+ if (
886+ sections_bounds [2 ] < geology_bounds [0 ]
887+ or sections_bounds [0 ] > geology_bounds [2 ]
888+ or sections_bounds [3 ] < geology_bounds [1 ]
889+ or sections_bounds [1 ] > geology_bounds [3 ]
890+ ):
891+ logger .warning ("AlongSection: Sections do not overlap geology extent; skipping thickness calculation." )
892+ return units
893+
894+ if geology_sindex is None :
895+ overlaps = False
896+ for geom in sections .geometry :
897+ if geom is None or geom .is_empty :
898+ continue
899+ if geology .intersects (geom ).any ():
900+ overlaps = True
901+ break
902+ if not overlaps :
903+ logger .warning ("AlongSection: Sections do not overlap geology; skipping thickness calculation." )
904+ return units
905+ else :
906+ overlaps = False
907+ for geom in sections .geometry :
908+ if geom is None or geom .is_empty :
909+ continue
910+ candidate_idx = list (geology_sindex .intersection (geom .bounds ))
911+ if not candidate_idx :
912+ continue
913+ if geology .iloc [candidate_idx ].intersects (geom ).any ():
914+ overlaps = True
915+ break
916+ if not overlaps :
917+ logger .warning ("AlongSection: Sections do not overlap geology; skipping thickness calculation." )
918+ return units
919+
901920 dip_column = next ((col for col in ("DIP" , "dip" , "Dip" ) if col in structure_data .columns ), None )
902921 orientation_tree = None
903922 orientation_dips = None
@@ -1006,7 +1025,7 @@ def compute(
10061025 )
10071026 default_dip_warning_emitted = True
10081027
1009- thickness = segment ["length" ] * abs ( math .sin (math .radians (dip_value ) ))
1028+ thickness = segment ["length" ] * math .sin (math .radians (dip_value ))
10101029 if thickness <= 0 :
10111030 continue
10121031
0 commit comments