diff --git a/.pre-commit-config.yaml b/.pre-commit-config.yaml index e8f62503..9ea4dd7b 100644 --- a/.pre-commit-config.yaml +++ b/.pre-commit-config.yaml @@ -34,7 +34,11 @@ repos: hooks: - id: nbstripout exclude: ^docs/ - args: [ '--extra-keys=metadata.kernelspec metadata.language_info' ] + args: [ + '--extra-keys=metadata.kernelspec metadata.language_info', + '--keep-metadata-keys=cell.metadata.tags cell.metadata.jupyter', + '--keep-id', + ] - repo: https://github.com/pre-commit/pre-commit-hooks rev: v6.0.0 diff --git a/notebooks/user_guide.ipynb b/notebooks/user_guide.ipynb index 1841ae47..0ff88eb1 100644 --- a/notebooks/user_guide.ipynb +++ b/notebooks/user_guide.ipynb @@ -3499,25 +3499,12 @@ "metadata": {}, "outputs": [], "source": [ - "from pedpy import (\n", - " compute_grid_cell_polygon_intersection_area,\n", - " get_grid_cells,\n", - ")\n", - "\n", "grid_size = 0.4\n", - "grid_cells, _, _ = get_grid_cells(walkable_area=walkable_area, grid_size=grid_size)\n", "\n", - "min_frame_profiles = 250 # We use here just an excerpt of the\n", - "max_frame_profiles = 400 # trajectory data to reduce compute time\n", + "min_frame_profiles = 500 # We use here just an excerpt of the\n", + "max_frame_profiles = 1000 # trajectory data to reduce compute time\n", "\n", - "profile_data = profile_data[profile_data.frame.between(min_frame_profiles, max_frame_profiles)]\n", - "\n", - "# Compute the grid intersection area for the resorted profile data (they have the same sorting)\n", - "# for usage in multiple calls to not run the compute heavy operation multiple times\n", - "(\n", - " grid_cell_intersection_area,\n", - " resorted_profile_data,\n", - ") = compute_grid_cell_polygon_intersection_area(data=profile_data, grid_cells=grid_cells)" + "profile_data = profile_data[profile_data.frame.between(min_frame_profiles, max_frame_profiles)]" ] }, { @@ -3567,44 +3554,31 @@ }, "outputs": [], "source": [ - "from pedpy import plot_measurement_setup\n", + "from pedpy import PEDPY_ORANGE, get_grid_cells, plot_measurement_setup\n", "\n", - "plot_measurement_setup(\n", + "ax = plot_measurement_setup(\n", " walkable_area=walkable_area,\n", " measurement_areas=[profile_measurement_area],\n", " ma_line_width=2,\n", " ma_alpha=0.2,\n", - ").set_aspect(\"equal\")\n", - "plt.show()" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "from pedpy import (\n", - " compute_grid_cell_polygon_intersection_area,\n", - " get_grid_cells,\n", ")\n", "\n", + "grid_size = 0.4\n", "grid_cells_measurement_area, _, _ = get_grid_cells(\n", " axis_aligned_measurement_area=profile_measurement_area,\n", " grid_size=grid_size,\n", ")\n", + "grid_cells = [MeasurementArea(cell) for cell in grid_cells_measurement_area]\n", + "plot_measurement_setup(\n", + " measurement_areas=grid_cells,\n", + " ma_line_color=PEDPY_ORANGE,\n", + " ma_line_width=0.3,\n", + " ma_alpha=0,\n", + " axes=ax,\n", + ")\n", + "ax.set_aspect(\"equal\")\n", "\n", - "min_frame_profiles = 250 # We use here just an excerpt of the\n", - "max_frame_profiles = 400 # trajectory data to reduce compute time\n", - "\n", - "profile_data = profile_data[profile_data.frame.between(min_frame_profiles, max_frame_profiles)]\n", - "\n", - "# Compute the grid intersection area for the resorted profile data (they have the same sorting)\n", - "# for usage in multiple calls to not run the compute heavy operation multiple times\n", - "(\n", - " grid_cell_intersection_area_measurement_area,\n", - " resorted_profile_data_measurement_area,\n", - ") = compute_grid_cell_polygon_intersection_area(data=profile_data, grid_cells=grid_cells_measurement_area)" + "plt.show()" ] }, { @@ -3667,17 +3641,15 @@ "from pedpy import SpeedMethod, compute_speed_profile\n", "\n", "voronoi_speed_profile = compute_speed_profile(\n", - " data=resorted_profile_data,\n", + " data=profile_data,\n", " walkable_area=walkable_area,\n", - " grid_intersections_area=grid_cell_intersection_area,\n", " grid_size=grid_size,\n", " speed_method=SpeedMethod.VORONOI,\n", ")\n", "\n", "arithmetic_speed_profile = compute_speed_profile(\n", - " data=resorted_profile_data,\n", + " data=profile_data,\n", " walkable_area=walkable_area,\n", - " grid_intersections_area=grid_cell_intersection_area,\n", " grid_size=grid_size,\n", " speed_method=SpeedMethod.ARITHMETIC,\n", ")\n", @@ -3770,19 +3742,17 @@ "from pedpy import AxisAlignedMeasurementArea, SpeedMethod, compute_speed_profile\n", "\n", "voronoi_speed_profile = compute_speed_profile(\n", - " data=resorted_profile_data_measurement_area,\n", + " data=profile_data,\n", " walkable_area=walkable_area,\n", " axis_aligned_measurement_area=profile_measurement_area,\n", - " grid_intersections_area=grid_cell_intersection_area_measurement_area,\n", " grid_size=grid_size,\n", " speed_method=SpeedMethod.VORONOI,\n", ")\n", "\n", "arithmetic_speed_profile = compute_speed_profile(\n", - " data=resorted_profile_data_measurement_area,\n", + " data=profile_data,\n", " walkable_area=walkable_area,\n", " axis_aligned_measurement_area=profile_measurement_area,\n", - " grid_intersections_area=grid_cell_intersection_area_measurement_area,\n", " grid_size=grid_size,\n", " speed_method=SpeedMethod.ARITHMETIC,\n", ")\n", @@ -3916,16 +3886,13 @@ "source": [ "from pedpy import DensityMethod, compute_density_profile\n", "\n", - "# here it is important to use the resorted data, as it needs to be in the same ordering as \"grid_cell_intersection_area\"\n", "voronoi_density_profile = compute_density_profile(\n", - " data=resorted_profile_data,\n", + " data=profile_data,\n", " walkable_area=walkable_area,\n", - " grid_intersections_area=grid_cell_intersection_area,\n", " grid_size=grid_size,\n", " density_method=DensityMethod.VORONOI,\n", ")\n", "\n", - "# here the unsorted data can be used\n", "classic_density_profile = compute_density_profile(\n", " data=profile_data,\n", " walkable_area=walkable_area,\n", @@ -4005,17 +3972,14 @@ "source": [ "from pedpy import DensityMethod, compute_density_profile\n", "\n", - "# here it is important to use the resorted data, as it needs to be in the same ordering as \"grid_cell_intersection_area\"\n", "voronoi_density_profile = compute_density_profile(\n", - " data=resorted_profile_data_measurement_area,\n", + " data=profile_data,\n", " walkable_area=walkable_area,\n", " axis_aligned_measurement_area=profile_measurement_area,\n", - " grid_intersections_area=grid_cell_intersection_area_measurement_area,\n", " grid_size=grid_size,\n", " density_method=DensityMethod.VORONOI,\n", ")\n", "\n", - "# here the unsorted data can be used\n", "classic_density_profile = compute_density_profile(\n", " data=profile_data,\n", " walkable_area=walkable_area,\n", diff --git a/pedpy/methods/profile_calculator.py b/pedpy/methods/profile_calculator.py index c22150ab..e1ee2911 100644 --- a/pedpy/methods/profile_calculator.py +++ b/pedpy/methods/profile_calculator.py @@ -28,7 +28,7 @@ WalkableArea, ) from pedpy.data.trajectory_data import TrajectoryData -from pedpy.errors import PedPyRuntimeError, PedPyTypeError, PedPyValueError +from pedpy.errors import PedPyTypeError, PedPyValueError from pedpy.internal.utils import alias @@ -212,8 +212,8 @@ def compute_profiles( # noqa: D417 density_method: density method to compute the density profile (default: :attr:`DensityMethod.VORONOI`) gaussian_width: full width at half maximum for Gaussian - approximation of the density, only needed when using - :attr:`DensityMethod.GAUSSIAN`. + approximation. Used when :attr:`DensityMethod.GAUSSIAN` or + :attr:`SpeedMethod.GAUSSIAN` is selected. axis_aligned_measurement_area (AxisAlignedMeasurementArea): Measurement area for which the profiles are computed. individual_voronoi_speed_data: deprecated alias for @@ -224,40 +224,77 @@ def compute_profiles( # noqa: D417 """ ( grid_cells, - _, - _, + rows, + cols, ) = get_grid_cells( walkable_area=walkable_area, axis_aligned_measurement_area=axis_aligned_measurement_area, grid_size=grid_size, ) - ( - grid_intersections_area, - internal_data, - ) = _compute_grid_polygon_intersection( - data=data, - grid_cells=grid_cells, + needs_intersection = density_method == DensityMethod.VORONOI or speed_method in ( + SpeedMethod.VORONOI, + SpeedMethod.ARITHMETIC, ) - density_profiles = compute_density_profile( - data=internal_data, - grid_intersections_area=grid_intersections_area, - density_method=density_method, - walkable_area=walkable_area, - axis_aligned_measurement_area=axis_aligned_measurement_area, - grid_size=grid_size, - gaussian_width=gaussian_width, - ) + # Pre-compute all frame-invariant quantities once. + bounds = None + if density_method == DensityMethod.CLASSIC or speed_method == SpeedMethod.MEAN: + if axis_aligned_measurement_area is not None: + bounds = axis_aligned_measurement_area.bounds + elif walkable_area is not None: + bounds = walkable_area.bounds + + x_center = y_center = None + if density_method == DensityMethod.GAUSSIAN or speed_method == SpeedMethod.GAUSSIAN: + grid_center = np.vectorize(shapely.centroid)(grid_cells) + x_center = shapely.get_x(grid_center[:cols]) + y_center = shapely.get_y(grid_center[::cols]) + + if density_method == DensityMethod.GAUSSIAN and gaussian_width is None: + raise PedPyValueError("Computing a Gaussian density profile needs a parameter 'gaussian_width'.") + if speed_method == SpeedMethod.GAUSSIAN and gaussian_width is None: + raise PedPyValueError("Computing a Gaussian speed profile needs a parameter 'gaussian_width'.") + + density_profiles: list[npt.NDArray[np.float64]] = [] + speed_profiles: list[npt.NDArray[np.float64]] = [] + + for _, frame_data in data.groupby(FRAME_COL): + # Precompute intersection once; shared by density and speed when both + # need it, avoiding a redundant shapely intersection call per frame. + grid_intersections_area_frame = None + if needs_intersection: + grid_intersections_area_frame = _compute_frame_grid_intersection( + frame_data=frame_data, + grid_cells=grid_cells, + ) - speed_profiles = compute_speed_profile( - data=internal_data, - grid_intersections_area=grid_intersections_area, - speed_method=speed_method, - walkable_area=walkable_area, - axis_aligned_measurement_area=axis_aligned_measurement_area, - grid_size=grid_size, - ) + density = _compute_density_for_frame( + frame_data=frame_data, + density_method=density_method, + grid_intersections_area_frame=grid_intersections_area_frame, + grid_cells=grid_cells, + bounds=bounds, + x_center=x_center, + y_center=y_center, + gaussian_width=gaussian_width, + grid_size=grid_size, + ) + density_profiles.append(density.reshape(rows, cols)) + + speed = _compute_speed_for_frame( + frame_data=frame_data, + speed_method=speed_method, + grid_intersections_area_frame=grid_intersections_area_frame, + grid_cells=grid_cells, + bounds=bounds, + x_center=x_center, + y_center=y_center, + gaussian_width=gaussian_width, + fill_value=np.nan, + grid_size=grid_size, + ) + speed_profiles.append(speed.reshape(rows, cols)) return ( density_profiles, @@ -299,9 +336,13 @@ def compute_density_profile( profiles density_method: density method to compute the density profile - grid_intersections_area: intersection of grid cells with the Voronoi - polygons (result from - :func:`compute_grid_cell_polygon_intersection_area`) + grid_intersections_area: (Optional) intersection of grid cells with + the Voronoi polygons (result from + :func:`compute_grid_cell_polygon_intersection_area`). If not + provided when using :attr:`DensityMethod.VORONOI`, the + intersections are computed on-the-fly per frame, + which avoids allocating the full (num_grid_cells x num_rows) matrix + at the cost of recomputing the per-frame intersection. gaussian_width: full width at half maximum for Gaussian approximation of the density, only needed when using :attr:`DensityMethod.GAUSSIAN`. @@ -319,63 +360,55 @@ def compute_density_profile( grid_size=grid_size, ) - grid_center = np.vectorize(shapely.centroid)(grid_cells) - x_center = shapely.get_x(grid_center[:cols]) - y_center = shapely.get_y(grid_center[::cols]) - data_grouped_by_frame = data.groupby(FRAME_COL) - density_profiles = [] - for ( - frame, - frame_data, - ) in data_grouped_by_frame: - if density_method == DensityMethod.VORONOI: - if grid_intersections_area is None: - raise PedPyRuntimeError( - "Computing a Voronoi density profile needs the parameter `grid_intersections_area`." - ) - - grid_intersections_area_frame = grid_intersections_area[ - :, - data_grouped_by_frame.indices[frame], - ] - - density = _compute_voronoi_density_profile( + # Fast path for VORONOI + pre-computed intersections: bypass the generic + # _compute_density_for_frame dispatch and call the implementation directly, + # since all pre-conditions are already met and no other invariants need + # pre-computation. + if density_method == DensityMethod.VORONOI and grid_intersections_area is not None: + _gia = grid_intersections_area + _grid_area = grid_cells[0].area + return [ + _compute_voronoi_density_profile( frame_data=frame_data, - grid_intersections_area=grid_intersections_area_frame, - grid_area=grid_cells[0].area, - ) - elif density_method == DensityMethod.CLASSIC: - if walkable_area is not None: - bounds = walkable_area.bounds - if axis_aligned_measurement_area is not None: - bounds = axis_aligned_measurement_area.bounds - - density = _compute_classic_density_profile( - frame_data=frame_data, - bounds=bounds, - grid_size=grid_size, - ) - elif density_method == DensityMethod.GAUSSIAN: - if gaussian_width is None: - raise PedPyValueError("Computing a Gaussian density profile needs a parameter 'gaussian_width'.") + grid_intersections_area=_gia[:, data_grouped_by_frame.indices[frame]], + grid_area=_grid_area, + ).reshape(rows, cols) + for frame, frame_data in data_grouped_by_frame + ] + + # Pre-compute frame-invariant quantities. + bounds = None + if density_method == DensityMethod.CLASSIC: + if axis_aligned_measurement_area is not None: + bounds = axis_aligned_measurement_area.bounds + elif walkable_area is not None: + bounds = walkable_area.bounds + + x_center = y_center = None + if density_method == DensityMethod.GAUSSIAN: + grid_center = np.vectorize(shapely.centroid)(grid_cells) + x_center = shapely.get_x(grid_center[:cols]) + y_center = shapely.get_y(grid_center[::cols]) + + if density_method == DensityMethod.GAUSSIAN and gaussian_width is None: + raise PedPyValueError("Computing a Gaussian density profile needs a parameter 'gaussian_width'.") - density = _compute_gaussian_density_profile( - frame_data=frame_data, - center_x=x_center, - center_y=y_center, - width=gaussian_width, - ) - else: - raise PedPyValueError("density method not accepted.") - - density_profiles.append( - density.reshape( - rows, - cols, - ) + density_profiles = [] + for _, frame_data in data_grouped_by_frame: + density = _compute_density_for_frame( + frame_data=frame_data, + density_method=density_method, + grid_intersections_area_frame=None, + grid_cells=grid_cells, + bounds=bounds, + x_center=x_center, + y_center=y_center, + gaussian_width=gaussian_width, + grid_size=grid_size, ) + density_profiles.append(density.reshape(rows, cols)) return density_profiles @@ -492,6 +525,86 @@ def _compute_gaussian_density( return np.array(gauss_density.T) +def _compute_density_for_frame( + *, + frame_data: pd.DataFrame, + density_method: DensityMethod, + grid_intersections_area_frame: Optional[npt.NDArray[np.float64]], + grid_cells: npt.NDArray[shapely.Polygon], + bounds: Optional[tuple[float, float, float, float]], + x_center: Optional[npt.NDArray[np.float64]], + y_center: Optional[npt.NDArray[np.float64]], + gaussian_width: Optional[float], + grid_size: float, +) -> npt.NDArray[np.float64]: + """Compute the density profile for a single frame. + + This is the single authoritative dispatch for all density methods. + All public profile functions (:func:`compute_profiles`, + :func:`compute_density_profile`) delegate here so that adding a new + :class:`DensityMethod` only requires changing this function. + + Args: + frame_data: DataFrame for a single frame. + density_method: Density method to use. + grid_intersections_area_frame: Pre-computed intersection areas for + this frame (shape: ``num_grid_cells x num_pedestrians``). + Pass ``None`` when using :attr:`DensityMethod.VORONOI` to have + the intersections computed on-the-fly. + grid_cells: Grid cells covering the area. Required for + :attr:`DensityMethod.VORONOI` (on-the-fly intersection and grid + area). + bounds: Bounding box ``(min_x, min_y, max_x, max_y)``; required for + :attr:`DensityMethod.CLASSIC`. + x_center: Grid cell centre x-coordinates; required for + :attr:`DensityMethod.GAUSSIAN`. + y_center: Grid cell centre y-coordinates; required for + :attr:`DensityMethod.GAUSSIAN`. + gaussian_width: FWHM of the Gaussian kernel; required for + :attr:`DensityMethod.GAUSSIAN`. + grid_size: Side length of one grid cell; required for + :attr:`DensityMethod.CLASSIC`. + + Returns: + Flat NumPy array of density values (one element per grid cell, + not yet reshaped into rows x cols). + """ + if density_method == DensityMethod.VORONOI: + if grid_intersections_area_frame is None: + grid_intersections_area_frame = _compute_frame_grid_intersection( + frame_data=frame_data, + grid_cells=grid_cells, + ) + return _compute_voronoi_density_profile( + frame_data=frame_data, + grid_intersections_area=grid_intersections_area_frame, + grid_area=grid_cells[0].area, + ) + elif density_method == DensityMethod.CLASSIC: + if bounds is None: + raise PedPyValueError("bounds is required for DensityMethod.CLASSIC") + return _compute_classic_density_profile( + frame_data=frame_data, + bounds=bounds, + grid_size=grid_size, + ) + elif density_method == DensityMethod.GAUSSIAN: + if gaussian_width is None: + raise PedPyValueError("gaussian_width is required for DensityMethod.GAUSSIAN") + if x_center is None: + raise PedPyValueError("x_center is required for DensityMethod.GAUSSIAN") + if y_center is None: + raise PedPyValueError("y_center is required for DensityMethod.GAUSSIAN") + return _compute_gaussian_density_profile( + frame_data=frame_data, + center_x=x_center, + center_y=y_center, + width=gaussian_width, + ) + else: + raise PedPyValueError("density method not accepted.") + + def compute_speed_profile( *, data: pd.DataFrame, @@ -500,7 +613,7 @@ def compute_speed_profile( speed_method: SpeedMethod, grid_intersections_area: Optional[npt.NDArray[np.float64]] = None, fill_value: float = np.nan, - gaussian_width: float = 0.5, + gaussian_width: Optional[float] = 0.5, axis_aligned_measurement_area: Optional[AxisAlignedMeasurementArea] = None, # pylint: disable=too-many-arguments ) -> Sequence[npt.NDArray[np.float64]]: @@ -539,7 +652,10 @@ def compute_speed_profile( speed profile grid_intersections_area: (Optional) intersection areas of grid cells with Voronoi polygons (result from - :func:`compute_grid_cell_polygon_intersection_area`) + :func:`compute_grid_cell_polygon_intersection_area`). If not + provided when using :attr:`SpeedMethod.VORONOI` or + :attr:`SpeedMethod.ARITHMETIC`, the intersections are computed + on-the-fly per frame. fill_value: fill value for cells with no pedestrians inside when using :attr:`SpeedMethod.MEAN` (default = `np.nan`) gaussian_width: (Optional) The full width at half maximum (FWHM) for @@ -567,72 +683,51 @@ def compute_speed_profile( data_grouped_by_frame = data.groupby(FRAME_COL) - speed_profiles = [] + # Pre-compute frame-invariant quantities. + bounds = None + if speed_method == SpeedMethod.MEAN: + if axis_aligned_measurement_area is not None: + bounds = axis_aligned_measurement_area.bounds + elif walkable_area is not None: + bounds = walkable_area.bounds - for ( - frame, - frame_data, - ) in data_grouped_by_frame: - if speed_method == SpeedMethod.VORONOI: - if grid_intersections_area is None: - raise PedPyRuntimeError( - "Computing a Arithmetic speed profile needs the parameter `grid_intersections_area`." - ) - grid_intersections_area_frame = grid_intersections_area[ - :, - data_grouped_by_frame.indices[frame], - ] + x_center = y_center = None + if speed_method == SpeedMethod.GAUSSIAN: + grid_center = np.vectorize(shapely.centroid)(grid_cells) + x_center = shapely.get_x(grid_center[:cols]) + y_center = shapely.get_y(grid_center[::cols]) - speed = _compute_voronoi_speed_profile( - frame_data=frame_data, - grid_intersections_area=grid_intersections_area_frame, - grid_area=grid_cells[0].area, - ) - elif speed_method == SpeedMethod.ARITHMETIC: - if grid_intersections_area is None: - raise PedPyRuntimeError( - "Computing a Arithmetic speed profile needs the parameter `grid_intersections_area`." - ) + if speed_method == SpeedMethod.GAUSSIAN and gaussian_width is None: + raise PedPyValueError("Computing a Gaussian speed profile needs a parameter 'gaussian_width'.") + + speed_profiles = [] + for frame, frame_data in data_grouped_by_frame: + # Resolve the frame-level intersection from the global matrix when the + # caller pre-computed it; otherwise pass None and let the helper + # compute it on-the-fly. + grid_intersections_area_frame = None + if grid_intersections_area is not None and speed_method in ( + SpeedMethod.VORONOI, + SpeedMethod.ARITHMETIC, + ): grid_intersections_area_frame = grid_intersections_area[ :, data_grouped_by_frame.indices[frame], ] - speed = _compute_arithmetic_voronoi_speed_profile( - frame_data=frame_data, - grid_intersections_area=grid_intersections_area_frame, - ) - elif speed_method == SpeedMethod.MEAN: - if walkable_area is not None: - bounds = walkable_area.bounds - if axis_aligned_measurement_area is not None: - bounds = axis_aligned_measurement_area.bounds - - speed = _compute_mean_speed_profile( - frame_data=frame_data, - bounds=bounds, - grid_size=grid_size, - fill_value=fill_value, - ) - elif speed_method == SpeedMethod.GAUSSIAN: - grid_center = np.vectorize(shapely.centroid)(grid_cells) - center_x = shapely.get_x(grid_center[:cols]) - center_y = shapely.get_y(grid_center[::cols]) - speed = _compute_gaussian_speed_profile( - frame_data=frame_data, - center_x=center_x, - center_y=center_y, - fwhm=gaussian_width, - ) - else: - raise PedPyValueError("Speed method not accepted.") - - speed_profiles.append( - speed.reshape( - rows, - cols, - ) + speed = _compute_speed_for_frame( + frame_data=frame_data, + speed_method=speed_method, + grid_intersections_area_frame=grid_intersections_area_frame, + grid_cells=grid_cells, + bounds=bounds, + x_center=x_center, + y_center=y_center, + gaussian_width=gaussian_width, + fill_value=fill_value, + grid_size=grid_size, ) + speed_profiles.append(speed.reshape(rows, cols)) return speed_profiles @@ -876,6 +971,101 @@ def _compute_mean_speed_profile( return speed +def _compute_speed_for_frame( + *, + frame_data: pd.DataFrame, + speed_method: SpeedMethod, + grid_intersections_area_frame: Optional[npt.NDArray[np.float64]], + grid_cells: npt.NDArray[shapely.Polygon], + bounds: Optional[tuple[float, float, float, float]], + x_center: Optional[npt.NDArray[np.float64]], + y_center: Optional[npt.NDArray[np.float64]], + gaussian_width: Optional[float], + fill_value: float, + grid_size: float, +) -> npt.NDArray[np.float64]: + """Compute the speed profile for a single frame. + + This helper centralizes per-frame dispatch among the supported + :class:`SpeedMethod` values for callers that route speed + computation through it. Some public profile functions may still use + method-specific optimized paths directly. + + + Args: + frame_data: DataFrame for a single frame. + speed_method: Speed method to use. + grid_intersections_area_frame: Pre-computed intersection areas for + this frame (shape: ``num_grid_cells x num_pedestrians``). + Pass ``None`` when using :attr:`SpeedMethod.VORONOI` or + :attr:`SpeedMethod.ARITHMETIC` to have the intersections computed + on-the-fly. + grid_cells: Grid cells covering the area. Required for + :attr:`SpeedMethod.VORONOI` and :attr:`SpeedMethod.ARITHMETIC` + (on-the-fly intersection and grid area). + bounds: Bounding box ``(min_x, min_y, max_x, max_y)``; required for + :attr:`SpeedMethod.MEAN`. + x_center: Grid cell centre x-coordinates; required for + :attr:`SpeedMethod.GAUSSIAN`. + y_center: Grid cell centre y-coordinates; required for + :attr:`SpeedMethod.GAUSSIAN`. + gaussian_width: FWHM of the Gaussian kernel; required for + :attr:`SpeedMethod.GAUSSIAN`. + fill_value: Value for empty cells; used by :attr:`SpeedMethod.MEAN`. + grid_size: Side length of one grid cell; required for + :attr:`SpeedMethod.MEAN`. + + Returns: + Flat NumPy array of speed values (one element per grid cell, + not yet reshaped into rows x cols). + """ + if speed_method == SpeedMethod.VORONOI: + if grid_intersections_area_frame is None: + grid_intersections_area_frame = _compute_frame_grid_intersection( + frame_data=frame_data, + grid_cells=grid_cells, + ) + return _compute_voronoi_speed_profile( + frame_data=frame_data, + grid_intersections_area=grid_intersections_area_frame, + grid_area=grid_cells[0].area, + ) + elif speed_method == SpeedMethod.ARITHMETIC: + if grid_intersections_area_frame is None: + grid_intersections_area_frame = _compute_frame_grid_intersection( + frame_data=frame_data, + grid_cells=grid_cells, + ) + return _compute_arithmetic_voronoi_speed_profile( + frame_data=frame_data, + grid_intersections_area=grid_intersections_area_frame, + ) + elif speed_method == SpeedMethod.MEAN: + if bounds is None: + raise PedPyValueError("bounds is required for SpeedMethod.MEAN") + return _compute_mean_speed_profile( + frame_data=frame_data, + bounds=bounds, + grid_size=grid_size, + fill_value=fill_value, + ) + elif speed_method == SpeedMethod.GAUSSIAN: + if gaussian_width is None: + raise PedPyValueError("gaussian_width is required for SpeedMethod.GAUSSIAN") + if x_center is None: + raise PedPyValueError("x_center is required for SpeedMethod.GAUSSIAN") + if y_center is None: + raise PedPyValueError("y_center is required for SpeedMethod.GAUSSIAN") + return _compute_gaussian_speed_profile( + frame_data=frame_data, + center_x=x_center, + center_y=y_center, + fwhm=gaussian_width, + ) + else: + raise PedPyValueError("Speed method not accepted.") + + def compute_grid_cell_polygon_intersection_area( *, data: pd.DataFrame, @@ -929,27 +1119,52 @@ def compute_grid_cell_polygon_intersection_area( ) +def _compute_frame_grid_intersection( + *, + frame_data: pd.DataFrame, + grid_cells: npt.NDArray[shapely.Polygon], +) -> npt.NDArray[np.float64]: + """Compute grid cell-polygon intersection areas for a single frame. + + Args: + frame_data: DataFrame for a single frame containing a 'polygon' column. + grid_cells: Grid cells used for computing the profiles. + + Returns: + Intersection areas with shape (num_grid_cells, num_pedestrians_in_frame). + """ + return shapely.area( + shapely.intersection( + np.asarray(grid_cells)[:, np.newaxis], + np.asarray(frame_data.polygon)[np.newaxis, :], + ) + ) + + def _compute_grid_polygon_intersection( *, - data, - grid_cells, -): + data: pd.DataFrame, + grid_cells: npt.NDArray[shapely.Polygon], +) -> Tuple[npt.NDArray[np.float64], pd.DataFrame]: internal_data = data.copy(deep=True) internal_data = internal_data.sort_values(by=FRAME_COL) internal_data = internal_data.reset_index(drop=True) - grid_intersections_area = shapely.area( - shapely.intersection( - np.array(grid_cells)[ - :, - np.newaxis, - ], - np.array(internal_data.polygon)[ - np.newaxis, - :, - ], - ) - ) + # Process frame-by-frame to avoid materializing a + # (num_grid_cells x num_total_rows) array of Shapely geometries, + # which causes excessive memory usage for large datasets. + # Preallocate the output array and fill it frame-by-frame to avoid + # doubling peak memory usage from concatenation. + num_rows = len(internal_data) + grid_intersections_area = np.empty((len(grid_cells), num_rows), dtype=np.float64) + + col_offset = 0 + for _, frame_data in internal_data.groupby(FRAME_COL): + frame_result = _compute_frame_grid_intersection(frame_data=frame_data, grid_cells=grid_cells) + frame_size = frame_result.shape[1] + grid_intersections_area[:, col_offset : col_offset + frame_size] = frame_result + col_offset += frame_size + return ( grid_intersections_area, internal_data, diff --git a/tests/unit_tests/methods/test_profile_calculator.py b/tests/unit_tests/methods/test_profile_calculator.py index 314192aa..0664a6ac 100644 --- a/tests/unit_tests/methods/test_profile_calculator.py +++ b/tests/unit_tests/methods/test_profile_calculator.py @@ -1,14 +1,26 @@ +import numpy as np +import pandas as pd import pytest import shapely from shapely.geometry import box +import pedpy.methods.profile_calculator as _profile_calc_module +from pedpy.column_identifier import FRAME_COL, ID_COL from pedpy.data.geometry import ( AxisAlignedMeasurementArea, MeasurementArea, WalkableArea, ) from pedpy.errors import PedPyTypeError, PedPyValueError -from pedpy.methods.profile_calculator import get_grid_cells +from pedpy.methods.profile_calculator import ( + DensityMethod, + SpeedMethod, + compute_density_profile, + compute_grid_cell_polygon_intersection_area, + compute_profiles, + compute_speed_profile, + get_grid_cells, +) def test_get_grid_cells_walkable_area(): @@ -77,3 +89,476 @@ def test_get_grid_cells_axis_aligned_measurement_area_not_instance(): match="`axis_aligned_measurement_area` must be an instance of AxisAlignedMeasurementArea", ): get_grid_cells(axis_aligned_measurement_area="not_an_area", grid_size=1.0) + + +def _make_profile_test_data(): + """Create test data with Voronoi polygons across two frames.""" + walkable_area = WalkableArea(shapely.box(0, 0, 2, 2)) + grid_size = 1.0 + grid_cells, rows, cols = get_grid_cells(walkable_area=walkable_area, grid_size=grid_size) + + # Two frames, each with 2 pedestrians and their Voronoi polygons + data = pd.DataFrame( + { + ID_COL: [1, 2, 1, 2], + FRAME_COL: [0, 0, 1, 1], + "x": [0.5, 1.5, 0.5, 1.5], + "y": [0.5, 1.5, 1.5, 0.5], + "speed": [1.0, 2.0, 1.5, 2.5], + "polygon": [ + shapely.box(0, 0, 1, 1), + shapely.box(1, 1, 2, 2), + shapely.box(0, 1, 1, 2), + shapely.box(1, 0, 2, 1), + ], + } + ) + return data, walkable_area, grid_size, grid_cells, rows, cols + + +# ── on-the-fly vs precomputed parity ───────────────────────────────────────── +# The three Voronoi/Arithmetic methods all support passing a pre-computed +# intersection matrix. One parametrized test covers all three paths. + + +@pytest.mark.parametrize( + "method_kwargs", + [ + {"density_method": DensityMethod.VORONOI}, + {"speed_method": SpeedMethod.VORONOI}, + {"speed_method": SpeedMethod.ARITHMETIC}, + ], + ids=["density-voronoi", "speed-voronoi", "speed-arithmetic"], +) +def test_precomputed_intersection_matches_on_the_fly(method_kwargs): + """Profiles computed with a pre-computed intersection matrix equal those computed on-the-fly.""" + data, walkable_area, grid_size, grid_cells, _, _ = _make_profile_test_data() + intersection_area, resorted_data = compute_grid_cell_polygon_intersection_area(data=data, grid_cells=grid_cells) + + is_density = "density_method" in method_kwargs + compute_fn = compute_density_profile if is_density else compute_speed_profile + + profiles_precomputed = compute_fn( + data=resorted_data, + walkable_area=walkable_area, + grid_size=grid_size, + grid_intersections_area=intersection_area, + **method_kwargs, + ) + profiles_on_the_fly = compute_fn( + data=data, + walkable_area=walkable_area, + grid_size=grid_size, + **method_kwargs, + ) + + assert len(profiles_precomputed) == len(profiles_on_the_fly) + for pre, fly in zip(profiles_precomputed, profiles_on_the_fly, strict=True): + np.testing.assert_allclose(pre, fly) + + +def test_grid_cell_polygon_intersection_frame_chunked(): + """Verify _compute_grid_polygon_intersection produces correct results with frame-by-frame processing.""" + data, _, _, grid_cells, _, _ = _make_profile_test_data() + + intersection_area, _resorted_data = compute_grid_cell_polygon_intersection_area(data=data, grid_cells=grid_cells) + + # 4 grid cells x 4 total rows + assert intersection_area.shape == (4, 4) + + # Each polygon exactly covers one grid cell, so each column should have + # exactly one non-zero entry equal to 1.0 (area of the 1x1 grid cell) + for col_idx in range(4): + col = intersection_area[:, col_idx] + assert np.isclose(np.sum(col), 1.0) + assert np.count_nonzero(col > 0) == 1 + + +# ── helpers ─────────────────────────────────────────────────────────────────── + + +def _single_cell_area() -> tuple[WalkableArea, float]: + """1×1 m walkable area with grid_size=1 → single 1 m² grid cell.""" + return WalkableArea(shapely.box(0, 0, 1, 1)), 1.0 + + +def _single_ped_data(n_frames: int = 1, speed: float = 1.0) -> pd.DataFrame: + """N frames, each with one ped whose Voronoi polygon covers the full 1×1 cell.""" + return pd.DataFrame( + { + ID_COL: list(range(n_frames)), + FRAME_COL: list(range(n_frames)), + "x": [0.5] * n_frames, + "y": [0.5] * n_frames, + "speed": [speed] * n_frames, + "polygon": [shapely.box(0, 0, 1, 1)] * n_frames, + } + ) + + +# ── compute_profiles ───────────────────────────────────────────────────────── + + +def test_compute_profiles_voronoi_returns_one_profile_per_frame(): + wa, gs = _single_cell_area() + n_frames = 3 + density_profiles, speed_profiles = compute_profiles( + data=_single_ped_data(n_frames=n_frames), + walkable_area=wa, + grid_size=gs, + speed_method=SpeedMethod.VORONOI, + ) + assert len(density_profiles) == n_frames + assert len(speed_profiles) == n_frames + # Each profile must be reshaped to (rows, cols) — i.e. (1, 1) for a single cell. + assert density_profiles[0].shape == (1, 1) + assert speed_profiles[0].shape == (1, 1) + + +def test_compute_profiles_classic_density_mean_speed_walkable_area(): + """Exercise bounds computation path via walkable_area.""" + wa, gs = _single_cell_area() + density_profiles, speed_profiles = compute_profiles( + data=_single_ped_data(n_frames=2), + walkable_area=wa, + grid_size=gs, + density_method=DensityMethod.CLASSIC, + speed_method=SpeedMethod.MEAN, + ) + assert len(density_profiles) == 2 + assert len(speed_profiles) == 2 + + +def test_compute_profiles_classic_density_mean_speed_axis_aligned(): + """Exercise bounds computation path via axis_aligned_measurement_area.""" + aa = AxisAlignedMeasurementArea(0, 0, 1, 1) + density_profiles, speed_profiles = compute_profiles( + data=_single_ped_data(n_frames=2), + axis_aligned_measurement_area=aa, + grid_size=1.0, + density_method=DensityMethod.CLASSIC, + speed_method=SpeedMethod.MEAN, + ) + assert len(density_profiles) == 2 + assert len(speed_profiles) == 2 + + +def test_compute_profiles_gaussian_density_and_speed(): + """Exercise x_center/y_center path via both GAUSSIAN methods.""" + wa = WalkableArea(shapely.box(0, 0, 2, 2)) + data = pd.DataFrame( + { + ID_COL: [1], + FRAME_COL: [0], + "x": [1.0], + "y": [1.0], + "speed": [1.0], + "polygon": [shapely.box(0, 0, 2, 2)], + } + ) + density_profiles, speed_profiles = compute_profiles( + data=data, + walkable_area=wa, + grid_size=1.0, + density_method=DensityMethod.GAUSSIAN, + speed_method=SpeedMethod.GAUSSIAN, + gaussian_width=1.0, + ) + assert len(density_profiles) == 1 + assert len(speed_profiles) == 1 + + +def test_compute_profiles_gaussian_density_raises_without_gaussian_width(): + wa, gs = _single_cell_area() + with pytest.raises(PedPyValueError, match="gaussian_width"): + compute_profiles( + data=_single_ped_data(), + walkable_area=wa, + grid_size=gs, + density_method=DensityMethod.GAUSSIAN, + speed_method=SpeedMethod.VORONOI, + ) + + +def test_compute_profiles_gaussian_speed_raises_without_gaussian_width(): + wa, gs = _single_cell_area() + with pytest.raises(PedPyValueError, match="gaussian_width"): + compute_profiles( + data=_single_ped_data(), + walkable_area=wa, + grid_size=gs, + density_method=DensityMethod.VORONOI, + speed_method=SpeedMethod.GAUSSIAN, + ) + + +# ── compute_density_profile missing branches ────────────────────────────────── + + +def test_classic_density_with_axis_aligned_measurement_area(): + aa = AxisAlignedMeasurementArea(0, 0, 1, 1) + profiles = compute_density_profile( + data=_single_ped_data(), + axis_aligned_measurement_area=aa, + grid_size=1.0, + density_method=DensityMethod.CLASSIC, + ) + np.testing.assert_allclose(profiles[0], [[1.0]]) + + +def test_gaussian_density_raises_without_gaussian_width(): + wa, gs = _single_cell_area() + with pytest.raises(PedPyValueError, match="gaussian_width"): + compute_density_profile( + data=_single_ped_data(), + walkable_area=wa, + grid_size=gs, + density_method=DensityMethod.GAUSSIAN, + ) + + +# ── _compute_density_for_frame else-branch ──────────────────────────────────── + + +def test_density_for_frame_invalid_method_raises(): + """Defensive else-branch: passing a non-DensityMethod value raises PedPyValueError.""" + wa, gs = _single_cell_area() + grid_cells, _, _ = get_grid_cells(walkable_area=wa, grid_size=gs) + with pytest.raises(PedPyValueError, match="density method not accepted"): + _profile_calc_module._compute_density_for_frame( + frame_data=_single_ped_data(n_frames=1), + density_method=None, # not a DensityMethod member → hits else branch + grid_intersections_area_frame=None, + grid_cells=grid_cells, + bounds=None, + x_center=None, + y_center=None, + gaussian_width=None, + grid_size=gs, + ) + + +# ── compute_speed_profile missing branches ──────────────────────────────────── + + +def test_gaussian_speed_raises_without_gaussian_width(): + wa, gs = _single_cell_area() + with pytest.raises(PedPyValueError, match="gaussian_width"): + compute_speed_profile( + data=_single_ped_data(), + walkable_area=wa, + grid_size=gs, + speed_method=SpeedMethod.GAUSSIAN, + gaussian_width=None, + ) + + +def test_mean_speed_with_axis_aligned_measurement_area(): + aa = AxisAlignedMeasurementArea(0, 0, 1, 1) + profiles = compute_speed_profile( + data=_single_ped_data(speed=2.0), + axis_aligned_measurement_area=aa, + grid_size=1.0, + speed_method=SpeedMethod.MEAN, + ) + np.testing.assert_allclose(profiles[0], [[2.0]]) + + +# ── _compute_speed_for_frame else-branch ────────────────────────────────────── + + +def test_speed_for_frame_invalid_method_raises(): + """Defensive else-branch: passing a non-SpeedMethod value raises PedPyValueError.""" + wa, gs = _single_cell_area() + grid_cells, _, _ = get_grid_cells(walkable_area=wa, grid_size=gs) + with pytest.raises(PedPyValueError, match="Speed method not accepted"): + _profile_calc_module._compute_speed_for_frame( + frame_data=_single_ped_data(n_frames=1), + speed_method=None, # not a SpeedMethod member → hits else branch + grid_intersections_area_frame=None, + grid_cells=grid_cells, + bounds=None, + x_center=None, + y_center=None, + gaussian_width=None, + fill_value=np.nan, + grid_size=gs, + ) + + +# ── compute_grid_cell_polygon_intersection_area: empty data ─────────────────── + + +def test_grid_cell_polygon_intersection_empty_data_returns_zero_columns(): + """Empty input → result shape (n_cells, 0) without concatenation error.""" + wa, gs = _single_cell_area() + grid_cells, _, _ = get_grid_cells(walkable_area=wa, grid_size=gs) + empty_data = pd.DataFrame( + { + ID_COL: pd.Series([], dtype=int), + FRAME_COL: pd.Series([], dtype=int), + "polygon": pd.Series([], dtype=object), + } + ) + intersection, _ = compute_grid_cell_polygon_intersection_area(data=empty_data, grid_cells=grid_cells) + assert intersection.shape == (len(grid_cells), 0) + + +# ── analytical correctness tests ────────────────────────────────────────────── + + +def test_voronoi_density_single_ped_covering_full_cell(): + wa, gs = _single_cell_area() + profiles = compute_density_profile( + data=_single_ped_data(), + walkable_area=wa, + grid_size=gs, + density_method=DensityMethod.VORONOI, + ) + assert len(profiles) == 1 + np.testing.assert_allclose(profiles[0], [[1.0]]) + + +def test_classic_density_single_ped_in_only_cell(): + wa, gs = _single_cell_area() + profiles = compute_density_profile( + data=_single_ped_data(), + walkable_area=wa, + grid_size=gs, + density_method=DensityMethod.CLASSIC, + ) + np.testing.assert_allclose(profiles[0], [[1.0]]) + + +def test_voronoi_speed_single_ped_covering_full_cell(): + wa, gs = _single_cell_area() + profiles = compute_speed_profile( + data=_single_ped_data(speed=2.5), + walkable_area=wa, + grid_size=gs, + speed_method=SpeedMethod.VORONOI, + ) + np.testing.assert_allclose(profiles[0], [[2.5]]) + + +def test_arithmetic_speed_two_peds_in_same_cell(): + """2 peds each covering half the cell → arithmetic mean of their speeds.""" + wa, gs = _single_cell_area() + data = pd.DataFrame( + { + ID_COL: [1, 2], + FRAME_COL: [0, 0], + "x": [0.25, 0.75], + "y": [0.5, 0.5], + "speed": [2.0, 4.0], + "polygon": [shapely.box(0, 0, 0.5, 1), shapely.box(0.5, 0, 1, 1)], + } + ) + profiles = compute_speed_profile( + data=data, + walkable_area=wa, + grid_size=gs, + speed_method=SpeedMethod.ARITHMETIC, + ) + np.testing.assert_allclose(profiles[0], [[3.0]]) + + +def test_mean_speed_two_peds_in_same_cell(): + """2 peds in the same cell → histogram-weighted mean of their speeds.""" + wa, gs = _single_cell_area() + data = pd.DataFrame( + { + ID_COL: [1, 2], + FRAME_COL: [0, 0], + "x": [0.3, 0.7], + "y": [0.4, 0.6], + "speed": [2.0, 4.0], + "polygon": [shapely.box(0, 0, 0.5, 1), shapely.box(0.5, 0, 1, 1)], + } + ) + profiles = compute_speed_profile( + data=data, + walkable_area=wa, + grid_size=gs, + speed_method=SpeedMethod.MEAN, + ) + np.testing.assert_allclose(profiles[0], [[3.0]]) + + +def test_gaussian_density_ped_equidistant_from_all_cell_centers(): + """Ped placed equidistant from all 4 cell centers → all 4 cells identical (symmetry).""" + wa = WalkableArea(shapely.box(0, 0, 2, 2)) + # Cell centers: (0.5,0.5), (1.5,0.5), (0.5,1.5), (1.5,1.5). + # Ped at (1,1) has |Δx|=|Δy|=0.5 to every center → identical Gaussian weight. + data = pd.DataFrame( + { + ID_COL: [1], + FRAME_COL: [0], + "x": [1.0], + "y": [1.0], + "polygon": [shapely.box(0, 0, 2, 2)], + } + ) + profiles = compute_density_profile( + data=data, + walkable_area=wa, + grid_size=1.0, + density_method=DensityMethod.GAUSSIAN, + gaussian_width=1.0, + ) + assert profiles[0].shape == (2, 2) + np.testing.assert_allclose(profiles[0], np.full((2, 2), profiles[0][0, 0]), rtol=1e-10) + + +# ── dispatch / program-flow tests ───────────────────────────────────────────── +# Strategy: spy on the private low-level implementations and assert they are +# called exactly once per frame. This verifies the dispatch logic in +# _compute_density_for_frame / _compute_speed_for_frame without needing to +# replicate expected output values. + + +@pytest.mark.parametrize( + "density_method,spy_target", + [ + (DensityMethod.VORONOI, "_compute_voronoi_density_profile"), + (DensityMethod.CLASSIC, "_compute_classic_density_profile"), + (DensityMethod.GAUSSIAN, "_compute_gaussian_density_profile"), + ], +) +def test_density_dispatch_calls_correct_impl_per_frame(mocker, density_method, spy_target): + """compute_density_profile routes to the correct implementation once per frame.""" + spy = mocker.spy(_profile_calc_module, spy_target) + wa, gs = _single_cell_area() + n_frames = 3 + compute_density_profile( + data=_single_ped_data(n_frames=n_frames), + walkable_area=wa, + grid_size=gs, + density_method=density_method, + gaussian_width=1.0, # ignored unless GAUSSIAN + ) + assert spy.call_count == n_frames + + +@pytest.mark.parametrize( + "speed_method,spy_target", + [ + (SpeedMethod.VORONOI, "_compute_voronoi_speed_profile"), + (SpeedMethod.ARITHMETIC, "_compute_arithmetic_voronoi_speed_profile"), + (SpeedMethod.MEAN, "_compute_mean_speed_profile"), + (SpeedMethod.GAUSSIAN, "_compute_gaussian_speed_profile"), + ], +) +def test_speed_dispatch_calls_correct_impl_per_frame(mocker, speed_method, spy_target): + """compute_speed_profile routes to the correct implementation once per frame.""" + spy = mocker.spy(_profile_calc_module, spy_target) + wa, gs = _single_cell_area() + n_frames = 3 + compute_speed_profile( + data=_single_ped_data(n_frames=n_frames), + walkable_area=wa, + grid_size=gs, + speed_method=speed_method, + gaussian_width=1.0, # ignored unless GAUSSIAN + ) + assert spy.call_count == n_frames