Skip to content

fix: use first line from shapely and add all dip/thickness estimates #246

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
rabii-chaarani merged 2 commits into from
Jan 31, 2026
Merged

fix: use first line from shapely and add all dip/thickness estimates #246

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
6edef4b
fix: use first line from shapely and add all dip/thickness estimates
lachlangrose Jan 27, 2026
c190e51
fix: add call methods
lachlangrose Jan 28, 2026
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
8 changes: 6 additions & 2 deletions map2loop/contact_extractor.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -2,7 +2,7 @@
import pandas
import shapely
from .logging import getLogger
from typing import Optional
from typing import Any, Optional

logger = getLogger(__name__)

Expand All @@ -13,7 +13,11 @@ def __init__(self, geology: geopandas.GeoDataFrame, faults: Optional[geopandas.G
self.contacts = None
self.basal_contacts = None
self.all_basal_contacts = None

def __call__(self, **kwds: Any) -> Any:
if 'stratigraphic_column' in kwds:
return self.extract_basal_contacts(kwds['stratigraphic_column'], save_contacts=kwds.get('save_contacts', True))
else:
return self.extract_all_contacts(save_contacts=kwds.get('save_contacts', True))
def extract_all_contacts(self, save_contacts: bool = True) -> geopandas.GeoDataFrame:
logger.info("Extracting contacts")
geology = self.geology.copy()
Expand Down
3 changes: 2 additions & 1 deletion map2loop/sampler.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -136,7 +136,8 @@ def __init__(self, spacing: float = 50.0, dtm_data: Optional[gdal.Dataset] = Non
spacing = max(spacing, 1.0)
self.spacing = spacing


def __call__(self, **kwargs):
return self.sample(**kwargs)
@beartype.beartype
def sample(
self, spatial_data: geopandas.GeoDataFrame
Expand Down
11 changes: 6 additions & 5 deletions map2loop/thickness_calculator.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -380,10 +380,8 @@ def compute(
for _, row in basal_contact.iterrows():
# find the shortest line between the basal contact points and top contact points
short_line = shapely.shortest_line(row.geometry, top_contact_geometry)
_lines.append(short_line[0])

# check if the short line is
if self.max_line_length is not None and short_line.length > self.max_line_length:
if self.max_line_length is not None and short_line[0].length > self.max_line_length:
continue
if self.dtm_data is not None:
inv_geotransform = gdal.InvGeoTransform(self.dtm_data.GetGeoTransform())
Expand All @@ -408,12 +406,14 @@ def compute(
# find the indices of the points that are within 5% of the length of the shortest line
try:
# GEOS 3.10.0+
indices = shapely.dwithin(short_line, interp_points, line_length * 0.25)
indices = shapely.dwithin(short_line[0], interp_points, line_length * 0.25)
except UnsupportedGEOSVersionError:
Comment on lines 406 to 410
Copy link

Copilot AI Jan 27, 2026

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The comment says points are selected within 5% of the shortest-line length, but the threshold used is line_length * 0.25 (25%). Please update the comment or adjust the factor so they match.

Copilot uses AI. Check for mistakes.
indices= numpy.array([shapely.distance(short_line[0],point)<= (line_length * 0.25) for point in interp_points])
Copy link

Copilot AI Jan 27, 2026

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This fallback branch is not formatted according to the repo’s Black configuration (e.g., missing spaces around = and after commas). Please reformat this line to match the rest of the codebase.

Suggested change
indices= numpy.array([shapely.distance(short_line[0],point)<= (line_length * 0.25) for point in interp_points])
indices = numpy.array([shapely.distance(short_line[0], point) <= (line_length * 0.25) for point in interp_points])

Copilot uses AI. Check for mistakes.
# get the dip of the points that are within
_dip = numpy.deg2rad(dip[indices])
_dips.append(_dip)
if len(_dip) > 0:
_lines.extend([short_line[0]]*len(_dip))
_dips.extend(_dip)
Comment on lines +414 to +416
Copy link

Copilot AI Jan 27, 2026

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This change alters the shape/semantics of the debug output (self.lines) by emitting one row per dip sample (_lines.extend(...); _dips.extend(...)). There are existing tests for InterpolatedStructure.compute, but none assert the contents/lengths of thickness_calculator.lines; please add a regression test to cover the new behavior (e.g., that len(lines) == len(lines['dip']) and multiple dip samples produce multiple line rows).

Copilot uses AI. Check for mistakes.
# calculate the true thickness t = L * sin(dip)
thickness = line_length * numpy.sin(_dip)

Expand Down Expand Up @@ -475,6 +475,7 @@ def compute(
else:
self.location_tracking = geopandas.GeoDataFrame(columns=['p1_x', 'p1_y', 'p1_z', 'p2_x', 'p2_y', 'p2_z', 'thickness', 'unit', 'geometry'], crs=basal_contacts.crs)
# Create GeoDataFrame for lines

Copy link

Copilot AI Jan 27, 2026

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

There appears to be trailing whitespace on this blank line. Please remove the extra spaces (or delete the blank line) to avoid lint/format churn.

Suggested change

Copilot uses AI. Check for mistakes.
self.lines = geopandas.GeoDataFrame(geometry=_lines, crs=basal_contacts.crs)
self.lines['dip'] = _dips

Expand Down
Loading