Tutorial: Routing Pipes Through a CAD Assembly
===============================================

This tutorial walks through a complete routing workflow using a window frame assembly — a
real mechanical component with two nested volumes that pipes must navigate through. By the
end, you will have routed multiple pipes of different types through the assembly and
visualized the results.

**What you will learn:**

- Loading a STEP CAD file and preparing volumes
- Defining port connection points
- Creating pipe types with physical properties
- Building a CadMap (the voxelized routing space)
- Running the pathfinding step
- Accessing and visualizing routing results
- Running the optional geometry optimization step

**Source script**: ``scripts/simple_test_cases/window_test_case.py``

.. contents:: Steps
   :local:
   :depth: 1


Prerequisites
-------------

Make sure ``routing`` is installed:

.. code-block:: bash

   pip install routing

You also need a STEP file of your assembly. This tutorial uses
``data/step/Cas_test_windows_V5.step`` (included in the repository).


Step 1 — Load the CAD Model
----------------------------

Load the STEP file using volmdlr and convert it to a ``VolumeModel``:

.. code-block:: python

   import os
   from pathlib import Path
   from volmdlr.step import Step
   from volmdlr.model import VolumeModel

   step_file = "data/step/Cas_test_windows_V5.step"
   model = Step.from_file(step_file)
   volume_model = model.to_volume_model()

   # Optional: set transparency for visualization
   volume_model.primitives[0].alpha = 0.5  # border volume
   volume_model.primitives[1].alpha = 0.5  # outer volume

The window model has two volumes:

- ``primitives[0]`` — the border frame (inner structure)
- ``primitives[1]`` — the outer shell

Pipes will be routed through and around both of these.


Step 2 — Define Weighting Rules
---------------------------------

:class:`~routing.core.design_rules.WeightingRule` tells the router how to weight paths
based on their distance from a CAD volume. A higher weight means the router avoids that
region; lower weight means it prefers it.

.. code-block:: python

   from routing.core.design_rules import WeightingRule

   weighting_rule = WeightingRule(
       function="linear",    # weight increases linearly with distance
       min_length=0.0,       # apply from distance 0 m
       max_length=0.2,       # apply up to distance 0.2 m
       min_value=1,          # weight = 1 at distance 0 (near surface)
       max_value=10,         # weight = 10 at distance 0.2 m (far from surface)
   )

This rule makes paths prefer to stay close to the CAD surfaces rather than floating in
open space.

**Available functions**: ``"linear"``, ``"sqrt"`` (square root), ``"constant"``.


Step 3 — Wrap Volumes as RuledVolumes
---------------------------------------

:class:`~routing.core.core.RuledVolume` attaches design rules to a CAD volume. Each
volume independently contributes its rules to the overall routing cost:

.. code-block:: python

   from routing.core.core import RuledVolume

   outer_volume = RuledVolume.from_volume_model(
       VolumeModel([volume_model.primitives[1]]),
       design_rules=[weighting_rule],
       name="outer",
   )
   border_volume = RuledVolume.from_volume_model(
       VolumeModel([volume_model.primitives[0]]),
       design_rules=[weighting_rule],
       name="border",
   )

.. tip::

   Different volumes can have different rules. For example, you might use a stricter
   ``WeightingRule`` near moving parts and a gentler one near static structures.


Step 4 — Define Pipe Types
----------------------------

:class:`~routing.core.core.PipeType` defines the physical properties of a pipe.
The most important parameters are the outer radius and the bend radius ratio
(minimum bend radius expressed as a multiple of the pipe diameter):

.. code-block:: python

   from piping_3d import piping
   from routing.core.core import PipeType

   # Large pipe: 30 mm outer diameter
   big_section = piping.Section(radius_equivalent=0.030 / 2)

   pipe_type_large = PipeType(
       section=big_section,
       name="LargePipe",
       radius_of_curvature_ratio_min=1.5,   # bend radius ≥ 1.5 × diameter
       length_before_turn_min=0.0,           # no mandatory straight before turn
       color=(0, 150, 200),                  # display color (RGB 0–255)
       type_="TypeA",                        # group label (used for pooling)
   )

   # Small pipe: 10 mm outer diameter, needs 30 mm straight before each turn
   small_section = piping.Section(radius_equivalent=0.010 / 2)

   pipe_type_small = PipeType(
       section=small_section,
       name="SmallPipe",
       radius_of_curvature_ratio_min=1.5,
       length_before_turn_min=0.030,         # 30 mm straight before turns
       color=(200, 0, 0),
       type_="TypeB",
   )

**Key parameters:**

- ``radius_of_curvature_ratio_min``: The pipe cannot bend tighter than this multiple of
  its diameter. Common values: 1.5× for flexible hoses, 3–5× for rigid pipes.
- ``length_before_turn_min``: Minimum straight length at both ends of every turn. Set to
  0 for fully flexible tubing.
- ``type_``: A string label used to group pipes for pooling (bundling). ``type_`` stored in the same
  ``PoolingSpecification`` will follow the same ``PoolingMode``.
- ``min_slope``: Optional minimum slope (metres/metre) for gravity-drained pipes.


Step 5 — Define Ports
-----------------------

A :class:`~routing.core.core.Port` is a connection point in 3D space. Each port has a
position, an orientation (the direction the pipe must exit or enter), and an optional
minimum straight length before the first turn:

.. code-block:: python

   from routing.core.core import Port

   # A port at (50 mm, -400 mm, 100 mm) pointing in +Y direction
   port_in = Port(
       coordinates=(0.05, -0.4, 0.1),   # (x, y, z) in metres
       direction=(0, 1, 0),              # pipe exits in +Y
       length=0.1,                       # 100 mm mandatory straight at the port
       name="inlet",
   )

   port_out = Port(
       coordinates=(0.05, 0.8, 0.1),
       direction=(0, -1, 0),             # pipe enters from +Y (points inward)
       length=0.1,
       name="outlet",
   )

.. note::

   The ``direction`` does not need to be a unit vector — it is normalized automatically.
   The convention is: the direction points **away from** the port along the pipe.


Step 6 — Create Specifications
--------------------------------

A :class:`~routing.core.core.Specification` pairs a start port, end port, and pipe type.
Create one per pipe to route:

.. code-block:: python

   from routing.core.core import Specification

   spec_large = Specification(
       start=port_in,
       end=port_out,
       pipe_type=pipe_type_large,
       name="route_large",
   )

For many pipes, it is convenient to define ports as ``PortCouple`` objects and build
specifications in a loop:

.. code-block:: python

   from routing.core.core import PortCouple

   # All port pairs for the assembly
   port_pairs = {
       "PipeA": PortCouple(
           Port(coordinates=(0.05, -0.4, 0.1), direction=(0, 1, 0), length=0.1),
           Port(coordinates=(0.05,  0.8, 0.1), direction=(0, -1, 0), length=0.1),
       ),
       "PipeB": PortCouple(
           Port(coordinates=(0.37, -0.4, 0.3), direction=(0, 1, 0), length=0.08),
           Port(coordinates=(0.45,  0.8, 0.18), direction=(0, -1, 0), length=0.04),
       ),
   }

   pipe_types = {"PipeA": pipe_type_large, "PipeB": pipe_type_small}

   specifications = [
       Specification(port_pairs[name].start, port_pairs[name].end, pipe_types[name], name)
       for name in port_pairs
   ]


Step 7 — Build the CadMap
--------------------------

The :class:`~routing.core.cadmap.CadMap` is the central object: it voxelizes the routing
space, computes distance fields from the CAD surfaces, and holds all design rules.

.. code-block:: python

   from routing.core.cadmap import CadMap
   from routing.core.design_rules import TurnPenalizer

   voxel_size = 0.021   # 21 mm voxels — about 0.7× the large pipe diameter

   # Global rules applied to all pipe types
   global_rules = [
       TurnPenalizer(cost=3.0),  # each 90° turn costs 3× a straight voxel
   ]

   # Automatically compute a tight bounding box around all volumes
   bounding_box = CadMap.build_bounding_volume(
       ruled_volumes=[outer_volume, border_volume],
       border_rules=[],
       voxel_size=voxel_size,
       voxel_margin=1,           # add 1 voxel of margin around the bounding box
   )

   cadmap = CadMap.from_project_inputs(
       design_rules=global_rules,
       ruled_volumes=[outer_volume, border_volume],
       bounding_box=bounding_box,
       voxel_size=voxel_size,
       exact_distances=True,     # use exact CAD mesh distances (slower but accurate)
       build_cad=False,          # do not rebuild CAD during voxelization
   )

**Choosing voxel_size:**

The voxel size controls the trade-off between accuracy and speed:

- **Too large**: paths may miss narrow passages or violate bend radius constraints
- **Too small**: voxelization becomes slow and memory-intensive
- A good starting point: ``voxel_size = 0.5 × (largest pipe radius)``

**exact_distances**: When ``True``, distances are computed from the actual CAD mesh
triangles (accurate but slower). When ``False``, distances use the voxel grid
approximation (fast but slightly less accurate for weighting rules).


Step 8 — Find Routing Paths
----------------------------

Create a :class:`~routing.core.route_planner.RoutePlanner` with your CadMap and a
pathfinding algorithm, then call ``generate()``:

.. code-block:: python

   from routing.core.finders import SmartFinder
   from routing.core.route_planner import RoutePlanner

   # A* with Manhattan distance heuristic, no diagonal movement
   finder = SmartFinder("AStar", "manhattan", "never")

   route_planner = RoutePlanner(cadmap, finder)

   routing_result = route_planner.generate(
       specifications=specifications,
       n_iterations=1,     # 0 = single best path (default); N = best + N alternatives
       steady_mode=False,  # False = stop at first failure; True = continue on failures
   )

**Available finders:**

- ``"AStar"`` — standard A* (fast, grid-aligned paths)
- ``"ThetaStar"`` — Theta* (smoother paths, slightly slower)
- ``"Dijkstra"`` — Dijkstra's algorithm (slower, always finds shortest path)

**Heuristics** (second argument): ``"manhattan"``, ``"euclidean"``, ``"octile"``

**Diagonal movement** (third argument): ``"never"``, ``"always"``, ``"only_when_no_obstacle"``


Step 9 — Access Results
------------------------

The ``routing_result`` object contains all found paths:

.. code-block:: python

   # Visualize the generated (grid-aligned) paths
   routing_result.plot_data_generated().plot()

   # Collect pipe bundles (packs) and display them in 3D
   packs = routing_result.collect_packs(only_optimized=False)
   routing_result.specific_packs_cad_view(packs)

   # Find bifurcation points (where bundled pipes split)
   bifurcation_pts = routing_result.cadmap.get_bifurcation_points(
       packs=packs,
       tolerance=voxel_size * 2,
       with_ports=False,
   )


Step 10 — Optimize Path Geometry (optional)
--------------------------------------------

The generated paths follow the voxel grid and can look "staircase-like". The optimization
step refines them into smooth, tangent-continuous geometry:

.. code-block:: python

   from routing.core.optimizer import Costs, Settings

   costs = Costs(
       length=1000,         # penalize total pipe length
       sideways=100,        # penalize lateral deviations
       short_line=5000,     # penalize very short pipe segments
       gravity=5000,        # penalize slope violations (for gravity-drained pipes)
       bend=250,            # penalize tight bends
       constraints=5000,    # penalize constraint violations
       interferences=1e8,   # strongly penalize CAD interferences
   )

   settings = Settings(
       voxel_size=0.01,               # fine voxel size for optimization domain
       max_shift=0.1,                 # max displacement of any control point (m)
       stabilization_threshold=5.0,  # stop when cost improvement < this threshold
       n_iterations=3,               # number of refinement passes
       max_iter=1000,                # max optimizer iterations per pipe
   )

   routing_result = route_planner.optimize(
       costs=costs,
       settings=settings,
       picked_path="best",    # optimize the best-found path (vs "last")
       max_refinements=0,     # 0 disables re-routing of problem segments (default is 5)
   )

   # Display optimized paths
   routing_result.plot_data_optimized().plot()

   packs_opt = routing_result.collect_packs(only_optimized=True)
   routing_result.specific_packs_cad_view(packs_opt)


Full Script
-----------

Here is the complete script:

.. code-block:: python

   import os
   from pathlib import Path
   from volmdlr.step import Step
   from volmdlr.model import VolumeModel
   from piping_3d import piping

   from routing.core.cadmap import CadMap
   from routing.core.core import Port, PortCouple, PipeType, RuledVolume, Specification
   from routing.core.design_rules import TurnPenalizer, WeightingRule
   from routing.core.finders import SmartFinder
   from routing.core.optimizer import Costs, Settings
   from routing.core.route_planner import RoutePlanner
   from routing.core.verbose import VerboseLevel, verbose

   verbose.set_level(VerboseLevel.INFO)

   # Load CAD
   model = Step.from_file("data/step/Cas_test_windows_V5.step")
   volume_model = model.to_volume_model()

   # Design rules
   weighting_rule = WeightingRule("linear", 0.0, 0.2, 1, 10)
   outer_volume = RuledVolume.from_volume_model(
       VolumeModel([volume_model.primitives[1]]), [weighting_rule], name="outer"
   )
   border_volume = RuledVolume.from_volume_model(
       VolumeModel([volume_model.primitives[0]]), [weighting_rule], name="border"
   )

   # Pipe types
   section = piping.Section(radius_equivalent=0.015)
   pipe_type = PipeType(section=section, name="PipeA", radius_of_curvature_ratio_min=1.5,
                        color=(0, 150, 200), type_="TypeA")

   # Ports and specifications
   spec = Specification(
       start=Port((0.05, -0.4, 0.1), (0, 1, 0), length=0.1),
       end=Port((0.05,  0.8, 0.1), (0, -1, 0), length=0.1),
       pipe_type=pipe_type,
       name="route_A",
   )

   # CadMap
   voxel_size = 0.021
   bounding_box = CadMap.build_bounding_volume([outer_volume, border_volume], [], voxel_size)
   cadmap = CadMap.from_project_inputs(
       design_rules=[TurnPenalizer(cost=3.0)],
       ruled_volumes=[outer_volume, border_volume],
       bounding_box=bounding_box,
       voxel_size=voxel_size,
       exact_distances=True,
   )

   # Route
   route_planner = RoutePlanner(cadmap, SmartFinder("AStar", "manhattan", "never"))
   result = route_planner.generate(specifications=[spec], n_iterations=1)  # best + 1 alternative

   # Display
   result.plot_data_generated().plot()


Next Steps
----------

- :doc:`design_rules` — add engineering constraints (gravity slope, clamping, straight lengths)
- :doc:`shaped_zones` — define zones that routing must avoid, prefer, or pass through
- :doc:`../user_guide/pooling` — bundle compatible pipes into organized packs
- :doc:`../user_guide/optimization` — tune the geometry optimizer