Add functionality to create centroids and centroid connector to model network

[lmz]

This PR adds end-to-end support for creating centroids and centroid connectors in a model network, including the core APIs, geometry/link generation logic, and tests/documentation updates needed to use and validate the workflow.

What this PR does

This PR introduces functionality to build model-zone access structures by generating:

• Centroid nodes (zone representation in the network)
• Centroid connector links between centroids and nearby network nodes
• The feature is designed to support common modeling workflows where zones/TAZs need to be integrated into the routable network with generated access/egress links.

Key changes

• Added APIs/utilities to create centroid nodes from input zone geometry/data.
• Added logic to create centroid connector links from each centroid to eligible nearby network nodes.
• Added validation and guardrails around connector creation inputs and outputs.
• Added/updated tests covering centroid and connector creation behavior.
• Updated docs/examples to demonstrate expected usage and outputs.

Why this change is needed

Many demand and assignment workflows require zone-level connectors to represent how trips enter/exit the physical network. This PR provides a native, repeatable way to generate those structures directly in network_wrangler, reducing manual preprocessing and ensuring consistent connector generation behavior.

Implementation notes

Centroid generation and connector generation are implemented as explicit network-building operations so they can be run as part of a reproducible wrangling pipeline.
Connector creation is based on proximity/eligibility rules to attach centroids to nearby nodes.
The new functionality is intended to be configurable for model-specific connector policies.

Testing

This PR includes tests for:

• Correct centroid creation from valid inputs
• Connector generation for expected nearest/eligible nodes
• Error/edge handling for invalid or missing inputs
• Basic integrity checks on produced nodes/links
• Backward compatibility

Existing workflows are unaffected unless the new centroid/connector functionality is explicitly invoked. No breaking changes are intended for users not using these new APIs.

Limitations / follow-ups

Connector selection policies may still need tuning for specific model conventions (e.g., max connectors per zone, facility/type filtering, directional behavior).
Future enhancements can expand configurability and performance for very large zone sets.

---

• To see the specific tasks where the Asana app for GitHub is being used, see below:
  • https://app.asana.com/0/0/1215523318902339

[lmz]

cc @i-am-sijia

[e-lo]

Some things to address :-)

[e-lo]

most of this belongs in: https://github.com/network-wrangler/network_wrangler/blob/main/docs/how_to.md

[e-lo]

1. Do we want ML_geometry to be mandatory?
2. This coercing should be its own function (e.g. _fill_missing_link_geometries_from_nodes())
3. All the apply() is very SLOW...especially when we are using it just for counting.

Suggest
(a) just doing the counting to make sure things work for a test on a tiny network (not here).
(b) i think you can do a select just on ones that don't work (if any) and then apply? If its null should be skipping entirely.

[e-lo]

1. this should be an option rather than default - in many (most?) cases this indicates an error which shouldn't be obfuscated.
2. should encapsulate in a function "merge nodes"
3. should probably error/warn if distances between are really big (i.e. almost definitely indicates an error)

[e-lo]

Add usage

[e-lo]

Should this be in a sub-module similar to /nodes /links?

[e-lo]

This helper can be removed...see below

[e-lo]

I was pretty sure there was a way to do this without the extra applies and this is what Claude came up with which I think is ok. Its basically just keeping a running tab of the min angle of already selected connectors for each zone and updating as they are selected. I only 68% understand this (but didn't burn a ton of time on it either) but it should be faster.

fit_col = f"{mode}_centroid_fit"
# rows are already sorted by [zone_id, fit_col, distance_from_centroid], so within a
# zone the first row is the best fit + closest.
for zone_num in mode_node_df[zone_id].unique():
    idx = mode_node_df.index[mode_node_df[zone_id] == zone_num]
    if len(idx) == 0:
        continue

    angles = mode_node_df.loc[idx, "centroid_angle"].to_numpy()
    fit = mode_node_df.loc[idx, fit_col].to_numpy()
    conn = np.zeros(len(idx), dtype=int)

    def circ_sep(a):  # wraparound angular distance of every angle from scalar a
        d = np.abs(angles - a)
        return np.minimum(d, 360 - d)

    conn[0] = 1                       # connector 1: best fit + closest
    min_sep = circ_sep(angles[0])     # running min separation to the selected set

    for k in range(2, num_centroid_connectors + 1):
        unsel = conn == 0
        if not unsel.any():
            break
        best_fit = fit[unsel].min()               # best remaining fit level
        eligible = unsel & (fit == best_fit)
        pick = np.argmax(np.where(eligible, min_sep, -np.inf))
        conn[pick] = k
        min_sep = np.minimum(min_sep, circ_sep(angles[pick]))   # incremental, O(n)

    mode_node_df.loc[idx, "connector_num"] = conn

[e-lo]

this should probably be a parameter and not hard coded...right? there are different limits in various softwares.

[e-lo]

if do this here - should do it to transit!?

[e-lo]

do vector operations applied to series...

def bearings_to_cardinal_directions(bearings: pd.Series, cardinal_only: bool = False) -> pd.Series:
    """Vectorized bearing (deg) -> cardinal/intercardinal label. NaN bearings map to None."""

    _CARDINAL_LABELS = np.array(
    [d.value for d in (CardinalDirection.N, CardinalDirection.E, CardinalDirection.S, CardinalDirection.W)]
) #also do intercardinal if want that....

    labels, sector_deg = (
        (_CARDINAL_LABELS, 90.0) if cardinal_only else (_INTERCARDINAL_LABELS, 45.0)
    )
    # Shift by half a sector so North straddles 0°, then floor-divide into sector indices.
    sector = ((bearings + sector_deg / 2) % 360) // sector_deg

    directions = pd.Series(pd.NA, index=bearings.index, dtype="object")
    valid = sector.notna()
    directions[valid] = labels[sector[valid].astype(int)]
    return directions

[lmz]

Thank you @e-lo ! I'll go through these ASAP, hopefully by end of week!