diff --git a/driptorch/firing/strip_contour.py b/driptorch/firing/strip_contour.py index 1944e7a..a6ea15e 100644 --- a/driptorch/firing/strip_contour.py +++ b/driptorch/firing/strip_contour.py @@ -112,7 +112,7 @@ def _init_paths(self, paths: dict, **kwargs) -> dict: source_line, neighborhood_size=1, z_multiplier=elevation_influence) # Determine level set values for ignition path slicing - if heat_depth == depth: + if not heat_depth or heat_depth == depth: levels = range(depth, int(np.max(cost_distance.data)), depth) else: levels = [depth] diff --git a/driptorch/io.py b/driptorch/io.py index 6eb1fe3..6011534 100644 --- a/driptorch/io.py +++ b/driptorch/io.py @@ -15,6 +15,7 @@ from shapely.geometry.base import BaseGeometry from shapely.ops import transform from typing import Union +import copy class Projector: @@ -183,7 +184,7 @@ def read_geojson_polygon(geojson: dict) -> Polygon: def write_geojson(geometries: list[BaseGeometry], src_epsg: int, dst_epsg: int = 4326, properties={}, - style={}, elapsed_time=None) -> dict: + style={}, elapsed_time=None, max_line_segment_time=None, path_properties: list[dict]=None) -> dict: """Write a list of shapely geometries to GeoJSON Parameters @@ -200,7 +201,10 @@ def write_geojson(geometries: list[BaseGeometry], src_epsg: int, dst_epsg: int = Rendering style applied to all features. Defaults to {}. elapsed_time : float, optional Time elapsed during the firing operation. Defaults to None. - + max_line_segment_time : int, optional + Maximum time (in milliseconds) for a line segment. + path_properties : list[dict], optional + Properties for each segment of the path. Returns ------- dict @@ -218,17 +222,95 @@ def write_geojson(geometries: list[BaseGeometry], src_epsg: int, dst_epsg: int = for i, geometry in enumerate(geometries): props = {} + # add properties for name in property_names: props[name] = properties[name][i] + # add style + props.update(style) + + # add any path_properties for this segment of the path + if path_properties is not None: + props.update(path_properties[i]) + features.append( { 'type': 'Feature', - 'properties': props | style, + 'properties': props, 'geometry': mapping(projector.forward(geometry)) } ) + if max_line_segment_time is not None and 'times' in features[0]['properties']: + + new_features = [] + + for feature in features: + + # only split lines + if feature['geometry']['type'] != 'LineString': + new_features.append(copy.deepcopy(feature)) + continue + + times = feature['properties']['times'] + coords = feature['geometry']['coordinates'] + prev_idx = 0 + cur_idx = 0 + + while cur_idx < len(times): + + while times[cur_idx] - times[prev_idx] > max_line_segment_time: + + # create a new feature ending at previous time + max time + end_time = times[prev_idx] + max_line_segment_time + new_times = times[prev_idx:cur_idx] + [end_time] + + fraction = max_line_segment_time / (times[cur_idx] - times[prev_idx]) + + # calculate distance for accumulated line + distance = LineString(coords[prev_idx:cur_idx+1]).length + + max_distance = distance * fraction + + # find the segment where the distance along accumulated line + # is past the max distance + accum_idx = prev_idx + accum_distance = 0 + while accum_distance < max_distance: + next_distance = Point(coords[accum_idx]).distance(Point(coords[accum_idx+1])) + accum_distance += next_distance + accum_idx += 1 + + subfraction = (max_distance - (accum_distance - next_distance)) / next_distance + + new_x = coords[accum_idx-1][0] + (coords[accum_idx][0] - coords[accum_idx-1][0]) * subfraction + new_y = coords[accum_idx-1][1] + (coords[accum_idx][1] - coords[accum_idx-1][1]) * subfraction + end_coords = (new_x, new_y) + new_coords = coords[prev_idx:accum_idx] + (end_coords,) + + new_feature = copy.deepcopy(feature) + new_feature['properties']['times'] = new_times + new_feature['geometry']['coordinates'] = new_coords + new_features.append(new_feature) + + # modify current feature to start at previous time + max_line_segment_time + feature['properties']['times'] = [end_time] + times[accum_idx:] + feature['geometry']['coordinates'] = (end_coords,) + coords[accum_idx:] + times = feature['properties']['times'] + coords = feature['geometry']['coordinates'] + + prev_idx = 0 + cur_idx = 0 + + cur_idx += 1 + + + if len(times) > 0: + new_features.append(copy.deepcopy(feature)) + + + features = new_features + # Compile the features in a feature collection geojson = { 'type': 'FeatureCollection', diff --git a/driptorch/pattern.py b/driptorch/pattern.py index 7c4bb66..891f885 100644 --- a/driptorch/pattern.py +++ b/driptorch/pattern.py @@ -4,7 +4,6 @@ # Core imports from __future__ import annotations -import pdb import copy from time import time as unix_time import warnings @@ -56,10 +55,12 @@ class Pattern: Path geometry epsg : int EPSG code for path geometry + properties: list[dict], optional + Set of properties of the path """ def __init__(self, heat: list[int], igniter: list[int], leg: list[int], times: list[list[float]], - geometry: list[LineString], epsg: int, + geometry: list[LineString], epsg: int, properties: list[dict] = None ): self.heat = heat @@ -69,6 +70,10 @@ def __init__(self, heat: list[int], igniter: list[int], leg: list[int], times: l self.geometry = geometry self.epsg = epsg + self.properties = properties + if self.properties is None: + self.properties = [{} for _ in range(len(heat))] + # Compute the total elapsed time for the ignition crew times_ak = ak.Array(self.times) min_time = ak.min(times_ak) @@ -76,7 +81,7 @@ def __init__(self, heat: list[int], igniter: list[int], leg: list[int], times: l self.elapsed_time = max_time - min_time @classmethod - def from_dict(cls, paths_dict: dict, epsg: int) -> Pattern: + def from_dict(cls, paths_dict: dict, epsg: int = None) -> Pattern: """Alternative constructor for initializing a Pattern object with a dictionary of path parameters @@ -85,7 +90,7 @@ def from_dict(cls, paths_dict: dict, epsg: int) -> Pattern: paths_dict : dict Dictionary of path parameters epsg : Int - EPSG code of path geometries + EPSG code of path geometries. If not specified, epsg must by in paths_dict Returns ------- @@ -93,7 +98,18 @@ def from_dict(cls, paths_dict: dict, epsg: int) -> Pattern: A new instance of Pattern """ + if epsg is None: + if 'epsg' not in paths_dict: + raise ValueError('epsg not specified and not in paths_dict') + epsg = paths_dict['epsg'] + paths_dict["geometry"] = [shape(x) for x in paths_dict["geometry"]] + + if 'properties' in paths_dict: + properties = paths_dict["properties"] + else: + properties = [{} for _ in range(len(paths_dict["heat"]))] + return cls( paths_dict["heat"], paths_dict["igniter"], @@ -101,6 +117,7 @@ def from_dict(cls, paths_dict: dict, epsg: int) -> Pattern: paths_dict["times"], paths_dict["geometry"], epsg, + properties, ) def to_dict(self) -> dict: @@ -117,6 +134,8 @@ def to_dict(self) -> dict: "igniter": self.igniter, "leg": self.leg, "times": self.times, + "epsg": self.epsg, + "properties": self.properties, # convert to geoJSON for storage "geometry": [x.__geo_interface__ for x in self.geometry], } @@ -133,9 +152,14 @@ def empty_path_dict() -> dict: return {"heat": [], "igniter": [], "leg": [], "geometry": []} - def to_json(self) -> dict: + def to_json(self, max_line_segment_time=None) -> dict: """Write the Pattern to a GeoJSON dictionary + Parameters + ---------- + max_line_segment_time : int, optional + Maximum time (in milliseconds) for a line segment. + Returns ------- dict @@ -176,6 +200,8 @@ def to_json(self) -> dict: properties=props, style=style, elapsed_time=self.elapsed_time, + max_line_segment_time=max_line_segment_time, + path_properties=self.properties ) def translate(self, x_off: float, y_off: float) -> Pattern: @@ -275,15 +301,19 @@ def to_quicfire(self, burn_unit: BurnUnit, filename: str = None, time_offset=0, geometry, times, self.elapsed_time, resolution=resolution ) - def merge(self, pattern: Pattern, time_offset: float = 0, inplace: bool = False) -> Pattern: + def merge(self, pattern: Pattern, time_offset: float = 0, time_offset_end: bool = True, inplace: bool = False) -> Pattern: """Merge an input pattern with self Parameters ---------- pattern : Pattern Input pattern to be merged into self - time_offset : float - Time offset between patterns (seconds) + time_offset : float, optional + Time offset for other pattern. See time_offset_end for when other pattern executes in parallel or + sequentially this self pattern. (seconds) + time_offset_end : bool, optional + If true, patterns are executed sequentially and time_offset is the time between patterns. + Otherwise, patterns are executed in parallel and time_offset is the delayed start time for the other pattern. inplace : bool, optional Overwrites Pattern object if true, otherwise return new Pattern object. Defaults to True. @@ -304,10 +334,16 @@ def merge(self, pattern: Pattern, time_offset: float = 0, inplace: bool = False) # Get an empty path dictionary merged_dict = self.empty_path_dict() - # Delay times in the second pattern by the elapsed time of the first pattern and offset - delayed_times = ( - ak.max(self.times) + time_offset + ak.Array(pattern.times) - ).to_list() + if time_offset_end: + # Delay times in the second pattern by the elapsed time of the first pattern and offset + delayed_times = ( + ak.max(self.times) + time_offset + ak.Array(pattern.times) + ).to_list() + else: + # Delay times in the second pattern by the offset + delayed_times = ( + time_offset + ak.Array(pattern.times) + ).to_list() # Build merged pattern attributes merged_dict["heat"] = self.heat + pattern.heat @@ -315,14 +351,20 @@ def merge(self, pattern: Pattern, time_offset: float = 0, inplace: bool = False) merged_dict["leg"] = self.leg + pattern.leg merged_dict["times"] = self.times + delayed_times merged_dict["geometry"] = self.geometry + pattern.geometry + merged_dict["properties"] = self.properties + pattern.properties # Instantiate a new Pattern object merged_pattern = self.from_dict(merged_dict, self.epsg) # Compute summed elapse time - merged_pattern.elapsed_time = ( - self.elapsed_time + time_offset + pattern.elapsed_time - ) + if time_offset_end: + merged_pattern.elapsed_time = ( + self.elapsed_time + time_offset + pattern.elapsed_time + ) + else: + merged_pattern.elapsed_time = ( + max(self.elapsed_time, time_offset + pattern.elapsed_time) + ) # If inplace then overwrite self with the merged pattern if inplace: diff --git a/driptorch/unit.py b/driptorch/unit.py index 425dc32..20d765c 100644 --- a/driptorch/unit.py +++ b/driptorch/unit.py @@ -15,7 +15,7 @@ # External imports import numpy as np -from shapely.geometry import Polygon, LineString, Point +from shapely.geometry import Polygon, LineString, Point, MultiPolygon from shapely import affinity @@ -183,6 +183,17 @@ def buffer_downwind(self, width: float) -> BurnUnit: # buffer to cut out the downfiring blackline area buffered_polygon = self.polygon.difference(fore_line_buffer) + if isinstance(buffered_polygon, MultiPolygon): + print('WARNING: downwind buffered geometry is multipolygon; using only largest polygon.') + multi = list(buffered_polygon) + max_area = -999999 + for poly in multi: + if poly.area > max_area: + max_area = poly.area + max_poly = poly + + buffered_polygon = max_poly + return BurnUnit(buffered_polygon, self.firing_direction, utm_epsg=self.utm_epsg) def difference(self, burn_unit: BurnUnit) -> BurnUnit: diff --git a/setup.py b/setup.py index 53e054a..b8b14ca 100644 --- a/setup.py +++ b/setup.py @@ -35,7 +35,7 @@ 'pandas==1.4.2', 'pyproj==3.3.1', 'folium==0.12.1.post1', - 'gcsfs==2022.1.0', + 'gcsfs==2024.2.0', 'zarr==2.13.3', 'scipy==1.9.3', 'scikit-image==0.19.3'