bump version

to dev

.github/workflows/ci-docs.yml

@@ -0,0 +1,30 @@
+name: ci-docs
+on:
+  push:
+    branches:
+      - master
+      - main
+permissions:
+  contents: write
+jobs:
+  deploy:
+    runs-on: ubuntu-latest
+    steps:
+      - uses: actions/checkout@v4
+      - name: Configure Git Credentials
+        run: |
+          git config user.name github-actions[bot]
+          git config user.email 41898282+github-actions[bot]@users.noreply.github.com
+      - uses: actions/setup-python@v5
+        with:
+          python-version: 3.12.9
+      - run: echo "cache_id=$(date --utc '+%V')" >> $GITHUB_ENV
+      - uses: actions/cache@v4
+        with:
+          key: mkdocs-material-${{ env.cache_id }}
+          path: .cache
+          restore-keys: |
+            mkdocs-material-
+      - run: pip install -e .
+      - run: pip install mkdocs-material mkdocstrings-python
+      - run: mkdocs gh-deploy --force

.gitignore

@@ -1,3 +1,6 @@
+# Custom
+dev-tools/
+
 # Byte-compiled / optimized / DLL files
 __pycache__/
 *.py[codz]

README.md

@@ -1,2 +1,49 @@
 # pg-rad
-Primary Gamma RADiation landscape
+Primary Gamma RADiation landscape - Development
+
+## Clone
+```
+git clone https://github.com/pim-n/pg-rad
+cd pg-rad
+git checkout dev
+```
+
+or
+
+```
+git@github.com:pim-n/pg-rad.git
+cd pg-rad
+git checkout dev
+```
+
+## Dependencies / venv
+
+With Python verion `>=3.12.4` and `<3.13`, create a virtual environment and install pg-rad.
+
+```
+python3 -m venv .venv
+source .venv/bin/activate
+```
+
+With the virtual environment activated, run:
+
+```
+pip install -e .[dev]
+```
+
+## Tests
+
+Tests can be run with `pytest` from the root directory of the repository. With the virtual environment activated, run:
+
+```
+pytest
+```
+
+## Local viewing of documentation
+
+PG-RAD uses [Material for MkDocs](https://squidfunk.github.io/mkdocs-material/) for generating documentation. It can be locally viewed by (in the venv) running:
+```
+mkdocs serve
+```
+
+where you can add the `--livereload` flag to automatically update the documentation as you write to the Markdown files.

docs/API/landscape/create_landscape_from_path.md

@@ -0,0 +1,4 @@
+---
+title: pg_rad.landscape.create_landscape_from_path
+---
+::: pg_rad.landscape.create_landscape_from_path

docs/API/landscape/landscape.md

@@ -0,0 +1,4 @@
+---
+title: pg_rad.landscape.Landscape
+---
+::: pg_rad.landscape.Landscape

docs/API/objects/object.md

@@ -0,0 +1,5 @@
+---
+title: pg_rad.objects.Object
+---
+
+::: pg_rad.objects.Object

docs/API/path/path.md

@@ -0,0 +1,4 @@
+---
+title: pg_rad.path.Path
+---
+::: pg_rad.path.Path

docs/API/path/path_from_RT90.md

@@ -0,0 +1,5 @@
+---
+title: pg_rad.path.path_from_RT90
+---
+::: pg_rad.path.path_from_RT90
+

docs/API/path/simplify_path.md

@@ -0,0 +1,4 @@
+---
+title: pg_rad.path.simplify_path
+---
+::: pg_rad.path.simplify_path

docs/API/sources/point_source.md

@@ -0,0 +1,5 @@
+---
+title: pg_rad.sources.PointSource
+---
+
+::: pg_rad.sources.PointSource

docs/index.md

@@ -0,0 +1,45 @@
+# Welcome!
+
+Primary Gamma RADiation Landscapes (PG-RAD) is a Python package for research in source localization. It can simulate mobile gamma spectrometry data acquired from vehicle-borne detectors along a predefined path (e.g. a road).
+
+## Requirements
+
+PG-RAD requires Python `3.12`. The guides below assume a unix-like system.
+
+## Installation (CLI)
+
+<!--pipx seems like a possible option to install python package in a contained environment on unix-->
+
+Lorem ipsum
+
+## Installation (Python module)
+
+If you are interested in using PG-RAD in another Python project, create a virtual environment first:
+
+```
+python3 -m venv .venv
+```
+
+Then install PG-RAD in it:
+
+```
+source .venv/bin/activate
+(.venv) pip install git+https://github.com/pim-n/pg-rad
+```
+
+See how to get started with PG-RAD with your own Python code [here](pg-rad-in-python).
+
+## For developers
+```
+git clone https://github.com/pim-n/pg-rad
+cd pg-rad
+git checkout dev
+```
+
+or
+
+```
+git@github.com:pim-n/pg-rad.git
+cd pg-rad
+git checkout dev
+```

docs/javascripts/mathjax.js

@@ -0,0 +1,18 @@
+window.MathJax = {
+  tex: {
+    inlineMath: [['$', '$'], ["\\(", "\\)"]],
+    displayMath: [['$$', '$$'], ["\\[", "\\]"]],
+    processEscapes: true,
+    processEnvironments: true
+  },
+  options: {
+    processHtmlClass: "arithmatex"
+  }
+};
+
+document$.subscribe(() => { 
+  MathJax.startup.output.clearCache()
+  MathJax.typesetClear()
+  MathJax.texReset()
+  MathJax.typesetPromise()
+})

docs/pg-rad-in-cli.md

@@ -0,0 +1,4 @@
+---
+title: Using PG-RAD in CLI
+---
+Lorem ipsum.

mkdocs.yml

@@ -0,0 +1,49 @@
+site_name: PG-RAD Documentation
+theme:
+  name: material
+  palette:
+    # Dark Mode
+    - scheme: slate
+      toggle:
+        icon: material/weather-sunny
+        name: Dark mode
+      primary: blue
+      accent: deep purple
+
+    # Light Mode
+    - scheme: default
+      toggle:
+        icon: material/weather-night
+        name: Light mode
+      primary: light blue
+      accent: blue
+  features:
+    - content.code.copy
+
+markdown_extensions:
+  - pymdownx.highlight:
+      anchor_linenums: true
+      line_spans: __span
+      pygments_lang_class: true
+  - pymdownx.inlinehilite
+  - pymdownx.snippets
+  - pymdownx.superfences
+  - pymdownx.arithmatex:
+      generic: true
+
+extra_javascript:
+  - javascripts/mathjax.js
+  - https://cdn.jsdelivr.net/npm/mathjax@3/es5/tex-mml-chtml.js
+
+plugins:
+- mkdocs-jupyter:
+    execute: false
+- mkdocstrings:
+    enabled: !ENV [ENABLE_MKDOCSTRINGS, true]
+    default_handler: python
+    locale: en
+    handlers:
+      python:
+        options:
+          show_source: false
+          show_root_heading: false

pyproject.toml

@@ -0,0 +1,32 @@
+[build-system]
+requires = ["setuptools>=64.0", "wheel"]
+build-backend = "setuptools.build_meta"
+
+[tool.setuptools.packages.find]
+where = ["src"]
+
+[project]
+name = "pg-rad"
+version = "0.2.1"
+authors = [
+  { name="Pim Nelissen", email="pi0274ne-s@student.lu.se" },
+]
+description = "Primary Gamma RADiation Landscape"
+readme = "README.md"
+requires-python = ">=3.12.4,<3.13"
+dependencies = [
+    "matplotlib>=3.9.2",
+    "numpy>=2",
+    "pandas>=2.3.1",
+    "piecewise_regression==1.5.0",
+    "pyyaml>=6.0.2"
+]
+license = "MIT"
+license-files = ["LICEN[CS]E*"]
+
+[project.urls]
+Homepage = "https://github.com/pim-n/pg-rad"
+Issues = "https://github.com/pim-n/pg-rad/issues"
+
+[project.optional-dependencies]
+dev = ["pytest", "mkinit", "notebook", "mkdocs-material", "mkdocstrings-python", "mkdocs-jupyter"]

requirements.txt

@@ -0,0 +1,6 @@
+matplotlib>=3.9.2
+notebook>=7.2.1
+numpy>=2
+pandas>=2.3.1
+piecewise_regression==1.5.0
+pyyaml>=6.0.2

src/pg_rad/configs/logging.yml

@@ -0,0 +1,15 @@
+version: 1
+disable_existing_loggers: false
+formatters:
+  simple:
+    format: '%(asctime)s - %(levelname)s: %(message)s'
+handlers:
+  stdout:
+    class: logging.StreamHandler
+    formatter: simple
+    stream: ext://sys.stdout
+loggers:
+  root:
+    level: INFO
+    handlers:
+      - stdout

src/pg_rad/dataloader/__init__.py

@@ -0,0 +1,8 @@
+# do not expose internal logger when running mkinit
+__ignore__ = ["logger"]
+
+from pg_rad.dataloader import dataloader
+
+from pg_rad.dataloader.dataloader import (load_data,)
+
+__all__ = ['dataloader', 'load_data']

src/pg_rad/dataloader/dataloader.py

@@ -0,0 +1,28 @@
+import logging
+
+import pandas as pd
+
+from pg_rad.exceptions import DataLoadError, InvalidCSVError
+
+logger = logging.getLogger(__name__)
+
+def load_data(filename: str) -> pd.DataFrame:
+    logger.debug(f"Attempting to load file: {filename}")
+
+    try:
+        df = pd.read_csv(filename, delimiter=',')
+
+    except FileNotFoundError as e:
+        logger.error(f"File not found: {filename}")
+        raise DataLoadError(f"File does not exist: {filename}") from e
+
+    except pd.errors.ParserError as e:
+        logger.error(f"Invalid CSV format: {filename}")
+        raise InvalidCSVError(f"Invalid CSV file: {filename}") from e
+
+    except Exception as e:
+        logger.exception(f"Unexpected error while loading {filename}")
+        raise DataLoadError("Unexpected error while loading data") from e
+
+    logger.debug(f"File loaded: {filename}")
+    return df

src/pg_rad/exceptions/__init__.py

@@ -0,0 +1,10 @@
+# do not expose internal logger when running mkinit
+__ignore__ = ["logger"]
+
+from pg_rad.exceptions import exceptions
+
+from pg_rad.exceptions.exceptions import (ConvergenceError, DataLoadError,
+                                          InvalidCSVError,)
+
+__all__ = ['ConvergenceError', 'DataLoadError', 'InvalidCSVError',
+           'exceptions']

src/pg_rad/exceptions/exceptions.py

@@ -0,0 +1,8 @@
+class ConvergenceError(Exception):
+    """Raised when an algorithm fails to converge."""
+
+class DataLoadError(Exception):
+    """Base class for data loading errors."""
+
+class InvalidCSVError(DataLoadError):
+    """Raised when a file is not a valid CSV."""

src/pg_rad/isotopes/__init__.py

@@ -0,0 +1,8 @@
+# do not expose internal logger when running mkinit
+__ignore__ = ["logger"]
+
+from pg_rad.isotopes import isotope
+
+from pg_rad.isotopes.isotope import (Isotope,)
+
+__all__ = ['Isotope', 'isotope']

src/pg_rad/isotopes/isotope.py

@@ -0,0 +1,23 @@
+class Isotope:
+    """Represents the essential information of an isotope.
+
+    Args:
+        name (str): Full name (e.g. Caesium-137).
+        E (float): Energy of the primary gamma in keV.
+        b (float): Branching ratio for the gamma at energy E.
+    """  
+    def __init__(
+            self,
+            name: str,
+            E: float,
+            b: float
+    ):        
+        if E <= 0:
+            raise ValueError("primary_gamma must be a positive energy (keV).")
+
+        if not (0 <= b <= 1):
+            raise ValueError("branching_ratio_pg must be a ratio (0 <= b <= 1)")
+
+        self.name = name
+        self.E = E
+        self.b = b

src/pg_rad/landscape/__init__.py

@@ -0,0 +1,8 @@
+# do not expose internal logger when running mkinit
+__ignore__ = ["logger"]
+
+from pg_rad.landscape import landscape
+
+from pg_rad.landscape.landscape import (Landscape, create_landscape_from_path,)
+
+__all__ = ['Landscape', 'create_landscape_from_path', 'landscape']

src/pg_rad/landscape/landscape.py

@@ -0,0 +1,131 @@
+import logging
+
+from matplotlib import pyplot as plt
+from matplotlib.patches import Circle
+import numpy as np
+
+from pg_rad.path import Path
+from pg_rad.objects import PointSource
+
+logger = logging.getLogger(__name__)
+
+class Landscape:
+    """A generic Landscape that can contain a Path and sources.
+
+    Args:
+        air_density (float, optional): Air density in kg / m^3. Defaults to 1.243.
+        size (int | tuple[int, int, int], optional): Size of the world. Defaults to 500.
+        scale (str, optional): The scale of the size argument passed. Defaults to 'meters'.
+    """    
+    def __init__(
+            self,
+            air_density: float = 1.243,
+            size: int | tuple[int, int, int] = 500,
+            scale = 'meters',
+            ):    
+        if isinstance(size, int):
+            self.world = np.zeros((size, size, size))
+        elif isinstance(size, tuple) and len(size) == 3:
+            self.world = np.zeros(size)
+        else:
+            raise TypeError("size must be an integer or a tuple of 3 integers.")
+
+        self.air_density = air_density
+        self.scale = scale
+
+        self.path: Path = None
+        self.sources: list[PointSource] = []
+        logger.debug("Landscape initialized.")
+    
+    def plot(self, z = 0):
+        """Plot a slice of the world at a height `z`.
+
+        Args:
+            z (int, optional): Height of slice. Defaults to 0.
+
+        Returns:
+            fig, ax: Matplotlib figure objects.
+        """        
+        x_lim, y_lim, _ = self.world.shape
+
+        fig, ax = plt.subplots()
+        ax.set_xlim(right=x_lim)
+        ax.set_ylim(top=y_lim)
+        ax.set_xlabel(f"X [{self.scale}]")
+        ax.set_ylabel(f"Y [{self.scale}]")
+
+        if not self.path == None:
+            ax.plot(self.path.x_list, self.path.y_list, 'bo-')
+        
+        for s in self.sources:
+            if np.isclose(s.z, z):
+                dot = Circle(
+                        (s.x, s.y),
+                        radius=5,
+                        color=s.color,
+                        zorder=5
+                )
+
+                ax.text(
+                    s.x + 0.06,
+                    s.y + 0.06,
+                    s.name,
+                    color=s.color,
+                    fontsize=10,
+                    ha="left",
+                    va="bottom",
+                    zorder=6
+                )
+
+                ax.add_patch(dot)
+                
+        return fig, ax
+
+    def add_sources(self, *sources: PointSource):
+        """Add one or more point sources to the world.
+
+        Args:
+            *sources (pg_rad.sources.PointSource): One or more sources, passed as
+            Source1, Source2, ...
+        Raises:
+            ValueError: If the source is outside the boundaries of the landscape.
+        """        
+
+        max_x, max_y, max_z = self.world.shape[:3]
+
+        if any(
+            not (0 <= source.x < max_x and
+                 0 <= source.y < max_y and
+                 0 <= source.z < max_z)
+            for source in sources
+        ):
+            raise ValueError("One or more sources are outside the landscape!")
+
+        self.sources.extend(sources)
+
+    def set_path(self, path: Path):
+        """
+        Set the path in the landscape.
+        """
+        self.path = path
+    
+def create_landscape_from_path(path: Path, max_z = 500):
+    """Generate a landscape from a path, using its dimensions to determine
+    the size of the landscape.
+
+    Args:
+        path (Path): A Path object describing the trajectory.
+        max_z (int, optional): Height of the world. Defaults to 500 meters.
+
+    Returns:
+        landscape (pg_rad.landscape.Landscape): A landscape with dimensions based on the provided Path.
+    """    
+    max_x = np.ceil(max(path.x_list))
+    max_y = np.ceil(max(path.y_list))
+
+    landscape = Landscape(
+        size = (max_x, max_y, max_z)
+        )
+    
+    landscape.path = path
+    return landscape

src/pg_rad/logging/__init__.py

@@ -0,0 +1,5 @@
+from pg_rad.logging import logger
+
+from pg_rad.logging.logger import (setup_logger,)
+
+__all__ = ['logger', 'setup_logger']

src/pg_rad/logging/logger.py

@@ -0,0 +1,20 @@
+import logging
+import pathlib
+
+import yaml
+
+def setup_logger(log_level: str = "WARNING"):
+    levels = ["DEBUG", "INFO", "WARNING", "ERROR", "CRITICAL"]
+    
+    if not log_level in levels:
+        raise ValueError(f"Log level must be one of {levels}.")
+    
+    base_dir = pathlib.Path(__file__).resolve().parent
+    config_file = base_dir / "configs" / "logging.yml"
+
+    with open(config_file) as f:
+        config = yaml.safe_load(f)
+    
+    config["loggers"]["root"]["level"] = log_level
+
+    logging.config.dictConfig(config)

src/pg_rad/objects/__init__.py

@@ -0,0 +1,13 @@
+# do not expose internal logger when running mkinit
+__ignore__ = ["logger"]
+
+from pg_rad.objects import detectors
+from pg_rad.objects import objects
+from pg_rad.objects import sources
+
+from pg_rad.objects.detectors import (Detector,)
+from pg_rad.objects.objects import (Object,)
+from pg_rad.objects.sources import (PointSource,)
+
+__all__ = ['Detector', 'Object', 'PointSource', 'detectors', 'objects',
+           'sources']

src/pg_rad/objects/objects.py

@@ -0,0 +1,33 @@
+import math
+from typing import Self
+
+class Object:
+    def __init__(
+            self,
+            x: float,
+            y: float,
+            z: float,
+            name: str = "Unnamed object",
+            color: str = 'grey'):
+        """
+        A generic object.
+
+        Args:
+            x (float): X coordinate.
+            y (float): Y coordinate.
+            z (float): Z coordinate.
+            name (str, optional): Name for the object. Defaults to "Unnamed object".
+            color (str, optional): Matplotlib compatible color string. Defaults to "red".
+        """
+        
+        self.x = x
+        self.y = y
+        self.z = z
+        self.name = name
+        self.color = color
+    
+    def distance_to(self, other: Self) -> float:
+        return math.dist(
+            (self.x, self.y, self.z),
+            (other.x, other.y, other.z),
+        )

src/pg_rad/objects/sources.py

@@ -0,0 +1,49 @@
+import logging
+
+from .objects import Object
+from pg_rad.isotopes import Isotope
+
+logger = logging.getLogger(__name__)
+
+class PointSource(Object):
+    _id_counter = 1
+    def __init__(
+            self,
+            x: float,
+            y: float,
+            z: float,
+            activity: int,
+            isotope: Isotope,
+            name: str | None = None,
+            color: str = "red"):
+        """A point source.
+
+        Args:
+            x (float): X coordinate.
+            y (float): Y coordinate.
+            z (float): Z coordinate.
+            activity (int): Activity A in MBq.
+            isotope (Isotope): The isotope.
+            name (str | None, optional): Can give the source a unique name.
+            Defaults to None, making the name sequential.
+            (Source-1, Source-2, etc.).
+            color (str, optional): Matplotlib compatible color string. Defaults to "red".
+        """
+        
+        self.id = PointSource._id_counter
+        PointSource._id_counter += 1
+
+        # default name derived from ID if not provided
+        if name is None:
+            name = f"Source {self.id}"
+
+        super().__init__(x, y, z, name, color)
+
+        self.activity = activity
+        self.isotope = isotope
+        self.color = color
+        
+        logger.debug(f"Source created: {self.name}")
+
+    def __repr__(self):
+        return f"PointSource(name={self.name}, pos={(self.x, self.y, self.z)}, isotope={self.isotope.name}, A={self.activity} MBq)"

src/pg_rad/path/__init__.py

@@ -0,0 +1,9 @@
+# do not expose internal logger when running mkinit
+__ignore__ = ["logger"]
+
+from pg_rad.path import path
+
+from pg_rad.path.path import (Path, PathSegment, path_from_RT90,
+                              simplify_path,)
+
+__all__ = ['Path', 'PathSegment', 'path', 'path_from_RT90', 'simplify_path']

src/pg_rad/path/path.py

@@ -0,0 +1,190 @@
+from collections.abc import Sequence
+import logging
+import math
+
+from matplotlib import pyplot as plt
+import numpy as np
+import pandas as pd
+import piecewise_regression
+
+from pg_rad.exceptions import ConvergenceError
+
+logger = logging.getLogger(__name__)
+
+class PathSegment:
+    def __init__(self, a: tuple[float, float], b: tuple[float, float]):
+        """A straight Segment of a Path, from (x_a, y_a) to (x_b, y_b).
+
+        Args:
+            a (tuple[float, float]): The starting point (x_a, y_a).
+            b (tuple[float, float]): The final point (x_b, y_b).
+        """        
+        self.a = a
+        self.b = b
+
+    def get_length(self) -> float:
+        return math.dist(self.a, self.b)
+    
+    length = property(get_length)
+
+    def __str__(self) -> str:
+        return str(f"({self.a}, {self.b})")
+    
+    def __getitem__(self, index) -> float:
+        if index == 0:
+            return self.a
+        elif index == 1:
+            return self.b
+        else:
+            raise IndexError
+
+class Path:
+    def __init__(
+            self,
+            coord_list: Sequence[tuple[float, float]],
+            z: float = 0,
+            path_simplify = False
+                 ):
+        """Construct a path of sequences based on a list of coordinates.
+
+        Args:
+            coord_list (Sequence[tuple[float, float]]): List of x,y coordinates.
+            z (float, optional): Height of the path. Defaults to 0.
+            path_simplify (bool, optional): Whether to pg_rad.path.simplify_path(). Defaults to False.
+        """        
+        
+        if len(coord_list) < 2:
+            raise ValueError("Must provide at least two coordinates as a list of tuples, e.g. [(x1, y1), (x2, y2)]")
+
+        x, y = tuple(zip(*coord_list))
+    
+        if path_simplify:
+            try:
+                x, y = simplify_path(list(x), list(y))
+            except ConvergenceError:
+                logger.warning("Continuing without simplifying path.")
+
+        self.x_list = list(x)
+        self.y_list = list(y)
+
+        coord_list = list(zip(x, y))
+
+        self.segments = [PathSegment(i, ip1) for i, ip1 in zip(coord_list, coord_list[1:])]
+
+        self.z = z
+
+        logger.debug("Path created.")
+
+    def get_length(self) -> float:
+        return sum([s.length for s in self.segments])
+    
+    length = property(get_length)
+
+    def __getitem__(self, index) -> PathSegment:
+        return self.segments[index]
+
+    def __str__(self) -> str:
+        return str([str(s) for s in self.segments])
+    
+    def plot(self, **kwargs):
+        """
+        Plot the path using matplotlib.
+        """  
+        plt.plot(self.x_list, self.y_list, **kwargs)
+
+def simplify_path(
+        x: Sequence[float],
+        y: Sequence[float],
+        keep_endpoints_equal: bool = False,
+        n_breakpoints: int = 3
+        ):
+    """From full resolution x and y arrays, return a piecewise linearly approximated/simplified pair of x and y arrays.
+
+    This function uses the `piecewise_regression` package. From a full set of
+    coordinate pairs, the function fits linear sections, automatically finding
+    the number of breakpoints and their positions.
+
+    On why the default value of n_breakpoints is 3, from the `piecewise_regression`
+    docs:
+    "If you do not have (or do not want to use) initial guesses for the number
+    of breakpoints, you can set it to n_breakpoints=3, and the algorithm will
+    randomly generate start_values. With a 50% chance, the bootstrap restarting
+    algorithm will either use the best currently converged breakpoints or
+    randomly generate new start_values, escaping the local optima in two ways in
+    order to find better global optima."
+
+    Args:
+        x (Sequence[float]): Full list of x coordinates.
+        y (Sequence[float]): Full list of y coordinates.
+        keep_endpoints_equal (bool, optional): Whether or not to force start
+        and end to be exactly equal to the original. This will worsen the linear
+        approximation at the beginning and end of path. Defaults to False.
+        n_breakpoints (int, optional): Number of breakpoints. Defaults to 3.
+
+    Returns:
+        x (list[float]): Reduced list of x coordinates.
+        y (list[float]): Reduced list of y coordinates.
+
+    Raises:
+        ConvergenceError: If the fitting algorithm failed to simplify the path.
+
+    Reference:
+        Pilgrim, C., (2021). piecewise-regression (aka segmented regression) in Python. Journal of Open Source Software, 6(68), 3859, https://doi.org/10.21105/joss.03859.
+    """
+    
+    logger.debug(f"Attempting piecewise regression on path.")
+
+    pw_fit = piecewise_regression.Fit(x, y, n_breakpoints=n_breakpoints)
+    pw_res = pw_fit.get_results()
+
+    if pw_res == None:
+        logger.warning("Piecewise regression failed to converge.")
+        raise ConvergenceError("Piecewise regression failed to converge.") 
+    
+    est = pw_res['estimates']
+
+    # extract and sort breakpoints
+    breakpoints_x = sorted(
+        v['estimate'] for k, v in est.items() if k.startswith('breakpoint')
+    )
+ 
+    x_points = [x[0]] + breakpoints_x + [x[-1]]
+
+    y_points = pw_fit.predict(x_points)
+
+    if keep_endpoints_equal:
+        logger.debug("Forcing endpoint equality.")
+        y_points[0] = y[0]
+        y_points[-1] = y[-1]
+
+    logger.info(
+        f"Piecewise regression reduced path from {len(x)-1} to {len(x_points)-1} segments."
+    )
+
+    return x_points, y_points
+
+def path_from_RT90(
+        df: pd.DataFrame,
+        east_col: str = "East",
+        north_col: str = "North",
+        **kwargs
+        ) -> Path:
+    """Construct a path from East and North formatted coordinates (RT90) in a Pandas DataFrame.
+
+    Args:
+        df (pandas.DataFrame): DataFrame containing at least the two columns noted in the cols argument.
+        east_col (str): The column name for the East coordinates.
+        north_col (str): The column name for the North coordinates.
+
+    Returns:
+        Path: A Path object built from the aquisition coordinates in the DataFrame.
+    """    
+    
+    east_arr = np.array(df[east_col]) - min(df[east_col])
+    north_arr = np.array(df[north_col]) - min(df[north_col])
+
+    coord_pairs = list(zip(east_arr, north_arr))
+
+    path = Path(coord_pairs, **kwargs)
+    logger.debug("Loaded path from provided RT90 coordinates.")
+    return path

tests/data/coordinates.csv

@@ -0,0 +1,107 @@
+,East,North
+0,1324671.2,6187244.9
+1,1324671.8,6187239.9
+2,1324672.7,6187235.0
+3,1324673.5,6187230.1
+4,1324675.1,6187225.4
+5,1324677.4,6187220.9
+6,1324679.2,6187216.3
+7,1324681.5,6187211.8
+8,1324683.9,6187207.4
+9,1324686.5,6187203.2
+10,1324689.0,6187198.9
+11,1324692.3,6187195.1
+12,1324695.6,6187191.5
+13,1324698.1,6187187.1
+14,1324701.9,6187184.1
+15,1324704.2,6187179.6
+16,1324707.7,6187176.0
+17,1324710.2,6187171.7
+18,1324712.8,6187167.4
+19,1324715.1,6187163.0
+20,1324718.3,6187159.2
+21,1324721.3,6187155.3
+22,1324725.0,6187151.9
+23,1324728.0,6187147.9
+24,1324732.0,6187145.0
+25,1324736.5,6187142.8
+26,1324741.0,6187140.7
+27,1324745.9,6187140.1
+28,1324750.6,6187138.6
+29,1324755.2,6187136.8
+30,1324760.1,6187136.4
+31,1324765.1,6187136.1
+32,1324770.1,6187135.9
+33,1324774.8,6187134.2
+34,1324779.8,6187133.7
+35,1324784.8,6187133.5
+36,1324789.8,6187133.3
+37,1324794.5,6187132.1
+38,1324799.3,6187131.1
+39,1324804.0,6187129.8
+40,1324808.0,6187126.7
+41,1324811.7,6187123.3
+42,1324815.2,6187119.8
+43,1324819.1,6187116.6
+44,1324822.3,6187112.9
+45,1324825.5,6187109.0
+46,1324828.6,6187105.1
+47,1324832.3,6187101.8
+48,1324836.6,6187099.3
+49,1324840.3,6187095.9
+50,1324843.9,6187092.4
+51,1324847.2,6187088.7
+52,1324851.6,6187086.5
+53,1324856.3,6187084.6
+54,1324860.1,6187081.4
+55,1324864.8,6187079.7
+56,1324868.9,6187076.9
+57,1324872.9,6187073.9
+58,1324876.6,6187070.4
+59,1324880.4,6187067.3
+60,1324884.3,6187064.1
+61,1324887.6,6187060.4
+62,1324891.1,6187056.8
+63,1324894.8,6187053.5
+64,1324898.1,6187049.8
+65,1324901.7,6187046.3
+66,1324905.5,6187043.1
+67,1324909.3,6187039.9
+68,1324913.0,6187036.4
+69,1324916.4,6187032.8
+70,1324919.6,6187029.0
+71,1324923.3,6187025.6
+72,1324926.3,6187021.6
+73,1324929.4,6187017.7
+74,1324933.0,6187014.2
+75,1324936.6,6187010.8
+76,1324939.8,6187007.0
+77,1324942.7,6187002.9
+78,1324945.9,6186999.1
+79,1324948.4,6186994.8
+80,1324951.9,6186991.2
+81,1324954.9,6186987.2
+82,1324957.4,6186982.8
+83,1324960.4,6186978.9
+84,1324962.7,6186974.4
+85,1324965.8,6186970.5
+86,1324968.2,6186966.1
+87,1324971.1,6186962.0
+88,1324973.7,6186957.8
+89,1324976.4,6186953.6
+90,1324978.8,6186949.2
+91,1324981.8,6186945.2
+92,1324984.3,6186940.9
+93,1324987.0,6186936.8
+94,1324989.3,6186932.3
+95,1324992.1,6186928.1
+96,1324994.2,6186923.6
+97,1324996.7,6186919.3
+98,1324998.5,6186914.8
+99,1325001.4,6186910.7
+100,1325003.7,6186906.2
+101,1325006.8,6186902.3
+102,1325009.9,6186898.4
+103,1325012.9,6186894.4
+104,1325015.3,6186890.0
+105,1325018.9,6186886.5

tests/test_path_functionality.py

@@ -0,0 +1,47 @@
+import pathlib
+
+import numpy as np
+import pytest
+
+from pg_rad.dataloader import load_data
+from pg_rad.path import Path, path_from_RT90
+
+@pytest.fixture
+def test_df():
+    csv_path = pathlib.Path(__file__).parent / "data/coordinates.csv"
+    return load_data(csv_path)
+
+def test_piecewise_regression(test_df):
+    """_Verify whether the intermediate points deviate less than 0.1 SD._"""
+        
+    p_full = path_from_RT90(
+        test_df,
+        east_col="East",
+        north_col="North",
+        simplify_path=False
+    )
+ 
+    p_simpl = path_from_RT90(
+        test_df,
+        east_col="East",
+        north_col="North",
+        simplify_path=True
+    )
+
+    x_f = np.array(p_full.x_list)
+    y_f = np.array(p_full.y_list)
+
+    x_s = np.array(p_simpl.x_list)
+    y_s = np.array(p_simpl.y_list)
+
+    sd = np.std(y_f)
+
+    for xb, yb in zip(x_s[1:-1], y_s[1:-1]):
+        # find nearest original x index
+        idx = np.argmin(np.abs(x_f - xb))
+        deviation = abs(yb - y_f[idx])
+
+        assert deviation < 0.1 * sd, (
+            f"Breakpoint deviation too large: {deviation:.4f} "
+            f"(threshold {0.1 * sd:.4f}) at x={xb:.2f}"
+            )

tests/test_sources.py

@@ -0,0 +1,35 @@
+import numpy as np
+import pytest
+
+from pg_rad.sources import PointSource
+
+@pytest.fixture
+def test_sources():
+    pos_a = np.random.rand(3)
+    pos_b = np.random.rand(3)
+
+    a = PointSource(*tuple(pos_a), strength = None)
+    b = PointSource(*tuple(pos_b), strength = None)
+
+    return pos_a, pos_b, a, b
+
+def test_if_distances_equal(test_sources):
+    """Verify whether from PointSource A to PointSource B is the same as B to A."""
+
+    _, _, a, b = test_sources
+
+    assert a.distance_to(b) == b.distance_to(a)
+
+def test_distance_calculation(test_sources):
+    """Verify whether distance between two PointSources is calculated correctly."""
+
+    pos_a, pos_b, a, b = test_sources
+
+    dx = pos_b[0] - pos_a[0]
+    dy = pos_b[1] - pos_a[1]
+    dz = pos_b[2] - pos_a[2]
+
+    assert np.isclose(
+        a.distance_to(b),
+        np.sqrt(dx**2 + dy**2 + dz**2),
+        rtol = 1e-12)