From 9a6b8c4c7a8cafef178b132832b663c328fcf01a Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Tue, 16 Jan 2024 14:09:30 -0500 Subject: [PATCH 01/48] wip(skeletonizing): add support for fixing autapses in graphene volumes --- igneous/task_creation/skeleton.py | 20 +++++++++++++++++--- igneous/tasks/skeleton.py | 31 ++++++++++++++++++++++++++++--- 2 files changed, 45 insertions(+), 6 deletions(-) diff --git a/igneous/task_creation/skeleton.py b/igneous/task_creation/skeleton.py index bb7d4a6..1c665fd 100644 --- a/igneous/task_creation/skeleton.py +++ b/igneous/task_creation/skeleton.py @@ -71,7 +71,7 @@ def create_skeletonizing_tasks( parallel=1, fill_missing=False, sharded=False, frag_path=None, spatial_index=True, synapses=None, num_synapses=None, - dust_global=False, + dust_global=False, fix_autapses=False, cross_sectional_area=False, cross_sectional_area_smoothing_window=5, ): @@ -121,6 +121,17 @@ def create_skeletonizing_tasks( fix_borders: Allows trivial merging of single overlap tasks. You'll only want to set this to false if you're working on single or non-overlapping volumes. + fix_autapses: Only possible for graphene volumes. Uses PyChunkGraph (PCG) information + to fix autapses (when a neuron synapses onto itself). This requires splitting + contacts between the edges of two touching voxels. The algorithm for doing this + requires much more memory. + + This works by comparing the PYC L2 and root layers. L1 is watershed. L2 is the + connections only within an atomic chunk. The root layer provides the global + connectivity. Autapses can be distinguished at the L2 level, above that, they + may not be (and certainly not at the root level). We extract the voxel connectivity + graph from L2 and perform the overall trace at root connectivity. + dust_threshold: don't skeletonize labels smaller than this number of voxels as seen by a single task. dust_global: Use global voxel counts for the dust threshold instead of from @@ -159,6 +170,9 @@ def create_skeletonizing_tasks( shape = Vec(*shape) vol = CloudVolume(cloudpath, mip=mip, info=info) + if fix_autapses and vol.meta.path.format != "graphene": + raise ValueError("fix_autapses can only be performed on graphene volumes.") + if dust_threshold > 0 and dust_global: cf = CloudFiles(cloudpath) vxctfile = cf.join(vol.key, 'stats', 'voxel_counts.mb') @@ -247,8 +261,7 @@ def task(self, shape, offset): spatial_grid_shape=shape.clone(), # used for writing index filenames synapses=bbox_synapses, dust_global=dust_global, - cross_sectional_area=bool(cross_sectional_area), - cross_sectional_area_smoothing_window=int(cross_sectional_area_smoothing_window), + ) def synapses_for_bbox(self, shape, offset): @@ -292,6 +305,7 @@ def on_finish(self): 'spatial_index': bool(spatial_index), 'synapses': bool(synapses), 'dust_global': bool(dust_global), + 'fix_autapses': bool(fix_autapses), 'cross_sectional_area': bool(cross_sectional_area), 'cross_sectional_area_smoothing_window': int(cross_sectional_area_smoothing_window), }, diff --git a/igneous/tasks/skeleton.py b/igneous/tasks/skeleton.py index cfe80f7..fb7a2b5 100644 --- a/igneous/tasks/skeleton.py +++ b/igneous/tasks/skeleton.py @@ -1,4 +1,4 @@ -from typing import Optional,Sequence,Dict +from typing import Optional, Sequence, Dict, List from functools import reduce import itertools @@ -24,6 +24,7 @@ from cloudvolume.lib import Vec, Bbox, sip from cloudvolume.datasource.precomputed.sharding import synthesize_shard_files +import cc3d import fastremap import kimimaro @@ -89,6 +90,7 @@ def __init__( cross_sectional_area_shape_delta:int = 150, dry_run:bool = False, strip_integer_attributes:bool = True, + fix_autapses:bool = False, ): super().__init__( cloudpath, shape, offset, mip, @@ -101,7 +103,8 @@ def __init__( spatial_grid_shape, synapses, bool(dust_global), bool(cross_sectional_area), int(cross_sectional_area_smoothing_window), int(cross_sectional_area_shape_delta), - bool(dry_run), bool(strip_integer_attributes) + bool(dry_run), bool(strip_integer_attributes), + bool(fix_autapses), ) self.bounds = Bbox(offset, Vec(*shape) + Vec(*offset)) self.index_bounds = Bbox(offset, Vec(*spatial_grid_shape) + Vec(*offset)) @@ -136,6 +139,10 @@ def execute(self): dust_threshold = 0 all_labels = self.apply_global_dust_threshold(vol, all_labels) + voxel_graph = None + if self.fix_autapses: + voxel_graph = self.voxel_connectivity_graph(vol, bbox) + skeletons = kimimaro.skeletonize( all_labels, self.teasar_params, object_ids=self.object_ids, @@ -148,6 +155,7 @@ def execute(self): fill_holes=self.fill_holes, parallel=self.parallel, extra_targets_after=extra_targets_after.keys(), + voxel_graph=voxel_graph, ) del all_labels @@ -189,7 +197,24 @@ def execute(self): if self.spatial_index: self.upload_spatial_index(vol, path, index_bbox, skeletons) - + + def voxel_connectivity_graph(self, vol:CloudVolume, bbox:Bbox) -> np.ndarray: + layer_2 = vol.download(bbox, stop_layer=2, agglomerate=True) + + sgx, sgy, sgz = list(np.ceil(bbox.size3() / vol.meta.graph_chunk_size).astype(int)) + + vcg = np.zeros(layer_2.shape, dtype=np.uint32, order="F") + + for gx,gy,gz in xyzrange([sgx, sgy, sgz]): + bbx = Bbox((gx,gy,gz), (gx+1, gy+1, gz+1)) + bbx *= vol.meta.graph_chunk_size + + cutout = np.asfortranarray(layer_2[bbx.to_slices()]) + vcg_cutout = cc3d.voxel_connectivity_graph(cutout, connectivity=26) + vcg[bbx.to_slices()] = vcg_cutout + + return vcg + def compute_cross_sectional_area(self, vol, bbox, skeletons): # Why redownload a bigger image? In order to avoid clipping the # cross sectional areas on the edges. From 08a8151d005bf53880eb30ee0f5982bda12aa1f0 Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Tue, 16 Jan 2024 14:14:41 -0500 Subject: [PATCH 02/48] feat: add support for graphene timestamp --- igneous/task_creation/skeleton.py | 6 +++++- igneous/tasks/skeleton.py | 5 ++++- 2 files changed, 9 insertions(+), 2 deletions(-) diff --git a/igneous/task_creation/skeleton.py b/igneous/task_creation/skeleton.py index 1c665fd..4c0b8fe 100644 --- a/igneous/task_creation/skeleton.py +++ b/igneous/task_creation/skeleton.py @@ -74,6 +74,7 @@ def create_skeletonizing_tasks( dust_global=False, fix_autapses=False, cross_sectional_area=False, cross_sectional_area_smoothing_window=5, + timestamp=None, ): """ Assign tasks with one voxel overlap in a regular grid @@ -166,6 +167,7 @@ def create_skeletonizing_tasks( to the total computation.) cross_sectional_area_smoothing_window: Perform a rolling average of the normal vectors across these many vectors. + timestamp: for graphene volumes only, you can specify the timepoint to use """ shape = Vec(*shape) vol = CloudVolume(cloudpath, mip=mip, info=info) @@ -261,7 +263,8 @@ def task(self, shape, offset): spatial_grid_shape=shape.clone(), # used for writing index filenames synapses=bbox_synapses, dust_global=dust_global, - + fix_autapses=bool(fix_autapses), + timestamp=timestamp, ) def synapses_for_bbox(self, shape, offset): @@ -306,6 +309,7 @@ def on_finish(self): 'synapses': bool(synapses), 'dust_global': bool(dust_global), 'fix_autapses': bool(fix_autapses), + 'timestamp': timestamp, 'cross_sectional_area': bool(cross_sectional_area), 'cross_sectional_area_smoothing_window': int(cross_sectional_area_smoothing_window), }, diff --git a/igneous/tasks/skeleton.py b/igneous/tasks/skeleton.py index fb7a2b5..e1354d7 100644 --- a/igneous/tasks/skeleton.py +++ b/igneous/tasks/skeleton.py @@ -91,6 +91,7 @@ def __init__( dry_run:bool = False, strip_integer_attributes:bool = True, fix_autapses:bool = False, + timestamp:Optional[int] = None, ): super().__init__( cloudpath, shape, offset, mip, @@ -104,7 +105,7 @@ def __init__( bool(cross_sectional_area), int(cross_sectional_area_smoothing_window), int(cross_sectional_area_shape_delta), bool(dry_run), bool(strip_integer_attributes), - bool(fix_autapses), + bool(fix_autapses), timestamp, ) self.bounds = Bbox(offset, Vec(*shape) + Vec(*offset)) self.index_bounds = Bbox(offset, Vec(*spatial_grid_shape) + Vec(*offset)) @@ -115,6 +116,8 @@ def execute(self): info=self.info, cdn_cache=False, parallel=self.parallel, fill_missing=self.fill_missing, + agglomerate=True, + timestamp=timestamp, ) bbox = Bbox.clamp(self.bounds, vol.bounds) index_bbox = Bbox.clamp(self.index_bounds, vol.bounds) From 7f2e48762724374ddc1b34bb73c17fd3ae019be2 Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Wed, 17 Jan 2024 02:21:18 -0500 Subject: [PATCH 03/48] wip: ensure shape matches atomic chunks --- igneous/task_creation/skeleton.py | 11 +++++++++-- igneous/tasks/skeleton.py | 4 +++- 2 files changed, 12 insertions(+), 3 deletions(-) diff --git a/igneous/task_creation/skeleton.py b/igneous/task_creation/skeleton.py index 4c0b8fe..1f03325 100644 --- a/igneous/task_creation/skeleton.py +++ b/igneous/task_creation/skeleton.py @@ -172,8 +172,15 @@ def create_skeletonizing_tasks( shape = Vec(*shape) vol = CloudVolume(cloudpath, mip=mip, info=info) - if fix_autapses and vol.meta.path.format != "graphene": - raise ValueError("fix_autapses can only be performed on graphene volumes.") + if fix_autapses: + if vol.meta.path.format != "graphene": + raise ValueError("fix_autapses can only be performed on graphene volumes.") + + if not np.all(shape % vol.meta.graph_chunk_size == 0): + raise ValueError( + f"shape must be a multiple of the graph chunk size. Got: {shape}, " + f"{vol.meta.graph_chunk_size}" + ) if dust_threshold > 0 and dust_global: cf = CloudFiles(cloudpath) diff --git a/igneous/tasks/skeleton.py b/igneous/tasks/skeleton.py index e1354d7..2a110e5 100644 --- a/igneous/tasks/skeleton.py +++ b/igneous/tasks/skeleton.py @@ -204,13 +204,15 @@ def execute(self): def voxel_connectivity_graph(self, vol:CloudVolume, bbox:Bbox) -> np.ndarray: layer_2 = vol.download(bbox, stop_layer=2, agglomerate=True) - sgx, sgy, sgz = list(np.ceil(bbox.size3() / vol.meta.graph_chunk_size).astype(int)) + shape = bbox.size()[:3] + sgx, sgy, sgz = list(np.ceil(shape / vol.meta.graph_chunk_size).astype(int)) vcg = np.zeros(layer_2.shape, dtype=np.uint32, order="F") for gx,gy,gz in xyzrange([sgx, sgy, sgz]): bbx = Bbox((gx,gy,gz), (gx+1, gy+1, gz+1)) bbx *= vol.meta.graph_chunk_size + bbx = Bbox.clamp(bbx, (0,0,0), shape) cutout = np.asfortranarray(layer_2[bbx.to_slices()]) vcg_cutout = cc3d.voxel_connectivity_graph(cutout, connectivity=26) From 2b892c0b2e2a21c72490fa5b395b6d21dc61069c Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Mon, 13 May 2024 19:32:47 -0400 Subject: [PATCH 04/48] fix: add back in xyzrange --- igneous/tasks/skeleton.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/igneous/tasks/skeleton.py b/igneous/tasks/skeleton.py index 2a110e5..26f4c4a 100644 --- a/igneous/tasks/skeleton.py +++ b/igneous/tasks/skeleton.py @@ -21,7 +21,7 @@ import cloudvolume from cloudvolume import CloudVolume, Skeleton, paths -from cloudvolume.lib import Vec, Bbox, sip +from cloudvolume.lib import Vec, Bbox, sip, xyzrange from cloudvolume.datasource.precomputed.sharding import synthesize_shard_files import cc3d From dcfc93854e6a37e1cb35e7ee02dcf765c3147373 Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Tue, 14 May 2024 14:26:32 -0400 Subject: [PATCH 05/48] fix(SkeletonTask): generalize frag_path for windows Also interpret bare paths as "file://" --- igneous/tasks/skeleton.py | 12 +++++++----- 1 file changed, 7 insertions(+), 5 deletions(-) diff --git a/igneous/tasks/skeleton.py b/igneous/tasks/skeleton.py index 26f4c4a..febe1f5 100644 --- a/igneous/tasks/skeleton.py +++ b/igneous/tasks/skeleton.py @@ -107,6 +107,8 @@ def __init__( bool(dry_run), bool(strip_integer_attributes), bool(fix_autapses), timestamp, ) + if isinstance(self.frag_path, str): + self.frag_path = cloudfiles.paths.normalize(self.frag_path) self.bounds = Bbox(offset, Vec(*shape) + Vec(*offset)) self.index_bounds = Bbox(offset, Vec(*spatial_grid_shape) + Vec(*offset)) @@ -123,7 +125,10 @@ def execute(self): index_bbox = Bbox.clamp(self.index_bounds, vol.bounds) path = skeldir(self.cloudpath) - path = os.path.join(self.cloudpath, path) + if self.frag_path is None: + path = vol.meta.join(self.cloudpath, path) + else: + path = CloudFiles(self.frag_path).join(path) all_labels = vol[ bbox.to_slices() ] all_labels = all_labels[:,:,:,0] @@ -191,10 +196,7 @@ def execute(self): return skeletons if self.sharded: - if self.frag_path: - self.upload_batch(vol, os.path.join(self.frag_path, skeldir(self.cloudpath)), index_bbox, skeletons) - else: - self.upload_batch(vol, path, index_bbox, skeletons) + self.upload_batch(vol, path, index_bbox, skeletons) else: self.upload_individuals(vol, path, bbox, skeletons) From f3e138f64043704a3639e11aa08ca5874c762c3c Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Tue, 14 May 2024 14:27:06 -0400 Subject: [PATCH 06/48] fix(SkeletonTask): use agglomerate and timestamp correctly --- igneous/tasks/skeleton.py | 8 +++++--- 1 file changed, 5 insertions(+), 3 deletions(-) diff --git a/igneous/tasks/skeleton.py b/igneous/tasks/skeleton.py index febe1f5..a76acfc 100644 --- a/igneous/tasks/skeleton.py +++ b/igneous/tasks/skeleton.py @@ -118,8 +118,6 @@ def execute(self): info=self.info, cdn_cache=False, parallel=self.parallel, fill_missing=self.fill_missing, - agglomerate=True, - timestamp=timestamp, ) bbox = Bbox.clamp(self.bounds, vol.bounds) index_bbox = Bbox.clamp(self.index_bounds, vol.bounds) @@ -130,7 +128,11 @@ def execute(self): else: path = CloudFiles(self.frag_path).join(path) - all_labels = vol[ bbox.to_slices() ] + all_labels = vol.download( + bbox.to_slices(), + agglomerate=True, + timestamp=self.timestamp + ) all_labels = all_labels[:,:,:,0] if self.mask_ids: From 971eb146f7518aa14eaac2bb530d6985ae804741 Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Tue, 14 May 2024 14:38:37 -0400 Subject: [PATCH 07/48] fix: merge error removed cross_sectional_area params --- igneous/task_creation/skeleton.py | 2 ++ 1 file changed, 2 insertions(+) diff --git a/igneous/task_creation/skeleton.py b/igneous/task_creation/skeleton.py index 1f03325..b748d1f 100644 --- a/igneous/task_creation/skeleton.py +++ b/igneous/task_creation/skeleton.py @@ -272,6 +272,8 @@ def task(self, shape, offset): dust_global=dust_global, fix_autapses=bool(fix_autapses), timestamp=timestamp, + cross_sectional_area=bool(cross_sectional_area), + cross_sectional_area_smoothing_window=int(cross_sectional_area_smoothing_window), ) def synapses_for_bbox(self, shape, offset): From a3ab5bb3e4aba11a5a99e6e5c05355906ab1f566 Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Fri, 17 May 2024 23:41:30 -0400 Subject: [PATCH 08/48] fix(SkeletonTask): apply timestamp in correct location --- igneous/tasks/skeleton.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/igneous/tasks/skeleton.py b/igneous/tasks/skeleton.py index a76acfc..2cd61d7 100644 --- a/igneous/tasks/skeleton.py +++ b/igneous/tasks/skeleton.py @@ -206,7 +206,7 @@ def execute(self): self.upload_spatial_index(vol, path, index_bbox, skeletons) def voxel_connectivity_graph(self, vol:CloudVolume, bbox:Bbox) -> np.ndarray: - layer_2 = vol.download(bbox, stop_layer=2, agglomerate=True) + layer_2 = vol.download(bbox, stop_layer=2, agglomerate=True, timestamp=self.timestamp) shape = bbox.size()[:3] sgx, sgy, sgz = list(np.ceil(shape / vol.meta.graph_chunk_size).astype(int)) From 86a8a8f5dfb8721e7006824bd02c040f8fac47e0 Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Wed, 22 May 2024 02:08:25 -0400 Subject: [PATCH 09/48] fix: fixup the edges of the vcg with root vcg --- igneous/tasks/skeleton.py | 45 ++++++++++++++++++++++++++++++++++++--- 1 file changed, 42 insertions(+), 3 deletions(-) diff --git a/igneous/tasks/skeleton.py b/igneous/tasks/skeleton.py index 2cd61d7..b4da274 100644 --- a/igneous/tasks/skeleton.py +++ b/igneous/tasks/skeleton.py @@ -151,7 +151,7 @@ def execute(self): voxel_graph = None if self.fix_autapses: - voxel_graph = self.voxel_connectivity_graph(vol, bbox) + voxel_graph = self.voxel_connectivity_graph(vol, bbox, all_labels) skeletons = kimimaro.skeletonize( all_labels, self.teasar_params, @@ -205,8 +205,18 @@ def execute(self): if self.spatial_index: self.upload_spatial_index(vol, path, index_bbox, skeletons) - def voxel_connectivity_graph(self, vol:CloudVolume, bbox:Bbox) -> np.ndarray: - layer_2 = vol.download(bbox, stop_layer=2, agglomerate=True, timestamp=self.timestamp) + def voxel_connectivity_graph( + self, + vol:CloudVolume, + bbox:Bbox, + root_labels:np.ndarray, + ) -> np.ndarray: + layer_2 = vol.download( + bbox, + stop_layer=2, + agglomerate=True, + timestamp=self.timestamp + ) shape = bbox.size()[:3] sgx, sgy, sgz = list(np.ceil(shape / vol.meta.graph_chunk_size).astype(int)) @@ -221,6 +231,35 @@ def voxel_connectivity_graph(self, vol:CloudVolume, bbox:Bbox) -> np.ndarray: cutout = np.asfortranarray(layer_2[bbx.to_slices()]) vcg_cutout = cc3d.voxel_connectivity_graph(cutout, connectivity=26) vcg[bbx.to_slices()] = vcg_cutout + del vcg_cutout + + del layer_2 + + # the proper way to do this would be to get the lowest the L3..LN root + # as needed, but the lazy way to do this is to get the root labels + # which will retain a few errors, but overall the error rate should be + # over 100x less. We need to shade in the sides of the connectivity graph + # with edges that represent the connections between the adjacent boxes. + + root_vcg = cc3d.voxel_connectivity_graph(root_labels, connectivity=26) + + for gx,gy,gz in xyzrange([sgx, sgy, sgz]): + bbx = Bbox((gx,gy,gz), (gx+1, gy+1, gz+1)) + bbx *= vol.meta.graph_chunk_size + bbx = Bbox.clamp(bbx, (0,0,0), shape) + + slicearr = [] + for i in range(3): + bbx1 = bbx.clone() + bbx1.maxpt[i] = bbx1.minpt[i] + 1 + slicearr.append(bbx1) + + bbx1 = bbx.clone() + bbx1.minpt[i] = bbx1.maxpt[i] - 1 + slicearr.append(bbx1) + + for bbx1 in slicearr: + vcg[bbx1.to_slices()] = root_vcg[bbx1.to_slices()] return vcg From 7238b3aaef45f8d261863011162614f91b75c5c3 Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Thu, 23 May 2024 20:33:25 -0400 Subject: [PATCH 10/48] feat: remove requirements for a queue parameter --- igneous_cli/cli.py | 22 +++++++++++----------- 1 file changed, 11 insertions(+), 11 deletions(-) diff --git a/igneous_cli/cli.py b/igneous_cli/cli.py index 7e8803c..0d88b4f 100644 --- a/igneous_cli/cli.py +++ b/igneous_cli/cli.py @@ -891,7 +891,7 @@ def meshgroup(): @meshgroup.command("xfer") @click.argument("src", type=CloudPath()) @click.argument("dest", type=CloudPath()) -@click.option('--queue', required=True, help="AWS SQS queue or directory to be used for a task queue. e.g. sqs://my-queue or ./my-queue. See https://github.com/seung-lab/python-task-queue", type=str) +@click.option('--queue', help="AWS SQS queue or directory to be used for a task queue. e.g. sqs://my-queue or ./my-queue. See https://github.com/seung-lab/python-task-queue", type=str) @click.option("--sharded", is_flag=True, default=False, help="Generate shard fragments instead of outputing mesh fragments.", show_default=True) @click.option("--dir", "mesh_dir", type=str, default=None, help="Write meshes into this directory instead of the one indicated in the info file.") @click.option('--magnitude', default=2, help="Split up the work with 10^(magnitude) prefix based tasks.", show_default=True) @@ -928,7 +928,7 @@ def mesh_xfer( @meshgroup.command("forge") @click.argument("path", type=CloudPath()) -@click.option('--queue', required=True, help="AWS SQS queue or directory to be used for a task queue. e.g. sqs://my-queue or ./my-queue. See https://github.com/seung-lab/python-task-queue", type=str) +@click.option('--queue', help="AWS SQS queue or directory to be used for a task queue. e.g. sqs://my-queue or ./my-queue. See https://github.com/seung-lab/python-task-queue", type=str) @click.option('--mip', default=0, help="Perform meshing using this level of the image pyramid.", show_default=True) @click.option('--shape', type=Tuple3(), default=(448, 448, 448), help="Set the task shape in voxels.", show_default=True) @click.option('--simplify/--skip-simplify', is_flag=True, default=True, help="Enable mesh simplification.", show_default=True) @@ -979,7 +979,7 @@ def mesh_forge( @meshgroup.command("merge") @click.argument("path", type=CloudPath()) -@click.option('--queue', required=True, help="AWS SQS queue or directory to be used for a task queue. e.g. sqs://my-queue or ./my-queue. See https://github.com/seung-lab/python-task-queue", type=str) +@click.option('--queue', help="AWS SQS queue or directory to be used for a task queue. e.g. sqs://my-queue or ./my-queue. See https://github.com/seung-lab/python-task-queue", type=str) @click.option('--magnitude', default=2, help="Split up the work with 10^(magnitude) prefix based tasks. Default: 2 (100 tasks)") @click.option('--nlod', default=0, help="(multires) How many extra levels of detail to create.", show_default=True) @click.option('--vqb', default=16, help="(multires) Vertex quantization bits for stored model representation. 10 or 16 only.", show_default=True) @@ -1011,7 +1011,7 @@ def mesh_merge(ctx, path, queue, magnitude, nlod, vqb, dir, min_chunk_size): @meshgroup.command("merge-sharded") @click.argument("path", type=CloudPath()) -@click.option('--queue', required=True, help="AWS SQS queue or directory to be used for a task queue. e.g. sqs://my-queue or ./my-queue. See https://github.com/seung-lab/python-task-queue", type=str) +@click.option('--queue', help="AWS SQS queue or directory to be used for a task queue. e.g. sqs://my-queue or ./my-queue. See https://github.com/seung-lab/python-task-queue", type=str) @click.option('--nlod', default=1, help="Number of levels of detail to create.", type=int, show_default=True) @click.option('--vqb', default=16, help="Vertex quantization bits. Can be 10 or 16.", type=int, show_default=True) @click.option('--compress-level', default=7, help="Draco compression level.", type=int, show_default=True) @@ -1060,7 +1060,7 @@ def mesh_sharded_merge( @meshgroup.command("rm") @click.argument("path", type=CloudPath()) -@click.option('--queue', required=True, help="AWS SQS queue or directory to be used for a task queue. e.g. sqs://my-queue or ./my-queue. See https://github.com/seung-lab/python-task-queue", type=str) +@click.option('--queue', help="AWS SQS queue or directory to be used for a task queue. e.g. sqs://my-queue or ./my-queue. See https://github.com/seung-lab/python-task-queue", type=str) @click.option('--magnitude', default=2, help="Split up the work with 10^(magnitude) prefix based tasks. Default: 2 (100 tasks)") @click.option('--dir', 'mesh_dir', default=None, help="Target this directory instead of the one indicated in the info file.") @click.pass_context @@ -1109,7 +1109,7 @@ def spatialindexgroup(): @spatialindexgroup.command("create") @click.argument("path", type=CloudPath()) -@click.option('--queue', required=True, help="AWS SQS queue or directory to be used for a task queue. e.g. sqs://my-queue or ./my-queue. See https://github.com/seung-lab/python-task-queue", type=str) +@click.option('--queue', help="AWS SQS queue or directory to be used for a task queue. e.g. sqs://my-queue or ./my-queue. See https://github.com/seung-lab/python-task-queue", type=str) @click.option('--shape', default="448,448,448", type=Tuple3(), help="Shape in voxels of each indexing task.", show_default=True) @click.option('--mip', default=0, help="Perform indexing using this level of the image pyramid.", show_default=True) @click.option('--fill-missing', is_flag=True, default=False, help="Interpret missing image files as background instead of failing.", show_default=True) @@ -1278,7 +1278,7 @@ def skeleton_merge( @skeletongroup.command("merge-sharded") @click.argument("path", type=CloudPath()) -@click.option('--queue', required=True, help="AWS SQS queue or directory to be used for a task queue. e.g. sqs://my-queue or ./my-queue. See https://github.com/seung-lab/python-task-queue", type=str) +@click.option('--queue', help="AWS SQS queue or directory to be used for a task queue. e.g. sqs://my-queue or ./my-queue. See https://github.com/seung-lab/python-task-queue", type=str) @click.option('--min-cable-length', default=1000, help="Skip objects smaller than this physical path length.", type=float, show_default=True) @click.option('--max-cable-length', default=None, help="Skip objects larger than this physical path length. Default: no limit", type=float) @click.option('--tick-threshold', default=0, help="Remove small \"ticks\", or branches from the main skeleton one at a time from smallest to largest. Branches larger than this are preserved. Default: no elimination", type=float) @@ -1329,7 +1329,7 @@ def skeleton_sharded_merge( @imagegroup.command("rm") @click.argument("path", type=CloudPath()) -@click.option('--queue', required=True, help="AWS SQS queue or directory to be used for a task queue. e.g. sqs://my-queue or ./my-queue. See https://github.com/seung-lab/python-task-queue", type=str) +@click.option('--queue', help="AWS SQS queue or directory to be used for a task queue. e.g. sqs://my-queue or ./my-queue. See https://github.com/seung-lab/python-task-queue", type=str) @click.option('--mip', default=0, help="Which mip level to start deleting from. Default: 0") @click.option('--num-mips', default=5, help="The number of mip levels to delete at once. Default: 5") @click.option('--shape', default=None, help="The size of each deletion task as a comma delimited list. Must be a multiple of the chunk size.", type=Tuple3()) @@ -1353,7 +1353,7 @@ def delete_images( @skeletongroup.command("rm") @click.argument("path", type=CloudPath()) -@click.option('--queue', required=True, help="AWS SQS queue or directory to be used for a task queue. e.g. sqs://my-queue or ./my-queue. See https://github.com/seung-lab/python-task-queue", type=str) +@click.option('--queue', help="AWS SQS queue or directory to be used for a task queue. e.g. sqs://my-queue or ./my-queue. See https://github.com/seung-lab/python-task-queue", type=str) @click.option('--magnitude', default=2, help="Split up the work with 10^(magnitude) prefix based tasks. Default: 2 (100 tasks)") @click.option('--dir', 'skel_dir', default=None, help="Target this directory instead of the one indicated in the info file.") @click.pass_context @@ -1371,7 +1371,7 @@ def skeleton_rm(ctx, path, queue, magnitude, skel_dir): @skeletongroup.command("xfer") @click.argument("src", type=CloudPath()) @click.argument("dest", type=CloudPath()) -@click.option('--queue', required=True, help="AWS SQS queue or directory to be used for a task queue. e.g. sqs://my-queue or ./my-queue. See https://github.com/seung-lab/python-task-queue", type=str) +@click.option('--queue', help="AWS SQS queue or directory to be used for a task queue. e.g. sqs://my-queue or ./my-queue. See https://github.com/seung-lab/python-task-queue", type=str) @click.option("--sharded", is_flag=True, default=False, help="Generate shard fragments instead of outputing mesh fragments.", show_default=True) @click.option("--dir", "skel_dir", type=str, default=None, help="Write skeletons into this directory instead of the one indicated in the info file.") @click.option('--magnitude', default=2, help="Split up the work with 10^(magnitude) prefix based tasks.", show_default=True) @@ -1435,7 +1435,7 @@ def spatialindexgroupskel(): @spatialindexgroupskel.command("create") @click.argument("path", type=CloudPath()) -@click.option('--queue', required=True, help="AWS SQS queue or directory to be used for a task queue. e.g. sqs://my-queue or ./my-queue. See https://github.com/seung-lab/python-task-queue", type=str) +@click.option('--queue', help="AWS SQS queue or directory to be used for a task queue. e.g. sqs://my-queue or ./my-queue. See https://github.com/seung-lab/python-task-queue", type=str) @click.option('--shape', default="512,512,512", type=Tuple3(), help="Shape in voxels of each indexing task.", show_default=True) @click.option('--mip', default=0, help="Perform indexing using this level of the image pyramid.", show_default=True) @click.option('--fill-missing', is_flag=True, default=False, help="Interpret missing image files as background instead of failing.", show_default=True) From 3392897af7e74164e210e6f03902a11591880f4d Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Thu, 23 May 2024 22:33:48 -0400 Subject: [PATCH 11/48] redesign: use vol.info to get skeleton path --- igneous/tasks/skeleton.py | 11 +---------- 1 file changed, 1 insertion(+), 10 deletions(-) diff --git a/igneous/tasks/skeleton.py b/igneous/tasks/skeleton.py index b4da274..2634d40 100644 --- a/igneous/tasks/skeleton.py +++ b/igneous/tasks/skeleton.py @@ -40,15 +40,6 @@ def filename_to_segid(filename): segid, = matches.groups() return int(segid) -def skeldir(cloudpath): - cf = CloudFiles(cloudpath) - info = cf.get_json('info') - - skel_dir = 'skeletons/' - if 'skeletons' in info: - skel_dir = info['skeletons'] - return skel_dir - def strip_integer_attributes(skeletons): for skel in skeletons: skel.extra_attributes = [ @@ -122,7 +113,7 @@ def execute(self): bbox = Bbox.clamp(self.bounds, vol.bounds) index_bbox = Bbox.clamp(self.index_bounds, vol.bounds) - path = skeldir(self.cloudpath) + path = vol.info.get("skeletons", "skeletons") if self.frag_path is None: path = vol.meta.join(self.cloudpath, path) else: From 6ff5f7e3a41d576f166972296e4cddca0a09d81a Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Thu, 23 May 2024 22:34:28 -0400 Subject: [PATCH 12/48] fix: several errors --- igneous/tasks/skeleton.py | 10 ++++++---- 1 file changed, 6 insertions(+), 4 deletions(-) diff --git a/igneous/tasks/skeleton.py b/igneous/tasks/skeleton.py index 2634d40..38db81e 100644 --- a/igneous/tasks/skeleton.py +++ b/igneous/tasks/skeleton.py @@ -212,14 +212,16 @@ def voxel_connectivity_graph( shape = bbox.size()[:3] sgx, sgy, sgz = list(np.ceil(shape / vol.meta.graph_chunk_size).astype(int)) - vcg = np.zeros(layer_2.shape, dtype=np.uint32, order="F") + vcg = np.zeros(layer_2.shape[:3], dtype=np.uint32, order="F") + + clamp_box = Bbox([0,0,0], shape) for gx,gy,gz in xyzrange([sgx, sgy, sgz]): bbx = Bbox((gx,gy,gz), (gx+1, gy+1, gz+1)) bbx *= vol.meta.graph_chunk_size - bbx = Bbox.clamp(bbx, (0,0,0), shape) + bbx = Bbox.clamp(bbx, clamp_box) - cutout = np.asfortranarray(layer_2[bbx.to_slices()]) + cutout = np.asfortranarray(layer_2[bbx.to_slices()][:,:,:,0]) vcg_cutout = cc3d.voxel_connectivity_graph(cutout, connectivity=26) vcg[bbx.to_slices()] = vcg_cutout del vcg_cutout @@ -237,7 +239,7 @@ def voxel_connectivity_graph( for gx,gy,gz in xyzrange([sgx, sgy, sgz]): bbx = Bbox((gx,gy,gz), (gx+1, gy+1, gz+1)) bbx *= vol.meta.graph_chunk_size - bbx = Bbox.clamp(bbx, (0,0,0), shape) + bbx = Bbox.clamp(bbx, clamp_box) slicearr = [] for i in range(3): From 9838867630ef7c27a5f0f6750d4f9894c8e5cb07 Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Thu, 23 May 2024 22:34:47 -0400 Subject: [PATCH 13/48] feat: add fix_autapses to cli --- igneous_cli/cli.py | 7 +++++-- 1 file changed, 5 insertions(+), 2 deletions(-) diff --git a/igneous_cli/cli.py b/igneous_cli/cli.py index 0d88b4f..91a65da 100644 --- a/igneous_cli/cli.py +++ b/igneous_cli/cli.py @@ -1176,6 +1176,7 @@ def skeletongroup(): @click.option('--fix-branching', is_flag=True, default=True, help="Trades speed for quality of branching at forks.", show_default=True) @click.option('--fix-borders', is_flag=True, default=True, help="Allows trivial merging of single voxel overlap tasks. Only switch off for datasets that fit in a single task.", show_default=True) @click.option('--fix-avocados', is_flag=True, default=False, help="Fixes somata where nuclei and cytoplasm have separate segmentations.", show_default=True) +@click.option('--fix-autapses', is_flag=True, default=False, help="(graphene only) Fixes autapses by using the PyChunkGraph.", show_default=True) @click.option('--fill-holes', is_flag=True, default=False, help="Preprocess each cutout to eliminate background holes and holes caused by entirely contained inclusions. Warning: May remove labels that are considered inclusions.", show_default=True) @click.option('--dust-threshold', default=1000, help="Skip skeletonizing objects smaller than this number of voxels within a cutout.", type=int, show_default=True) @click.option('--dust-global/--dust-local', is_flag=True, default=False, help="Use global voxel counts for the dust threshold (when >0). To use this feature you must first compute the global voxel counts using the 'igneous image voxels' command.", show_default=True) @@ -1190,14 +1191,15 @@ def skeletongroup(): @click.option('--sharded', is_flag=True, default=False, help="Generate shard fragments instead of outputing skeleton fragments.", show_default=True) @click.option('--labels', type=TupleN(), default=None, help="Skeletonize only this comma separated list of labels.", show_default=True) @click.option('--cross-section', type=int, default=0, help="Compute the cross sectional area for each skeleton vertex. May add substantial computation time. Integer value is the normal vector rolling average smoothing window over vertices. 0 means off.", show_default=True) +@click.option('--output', '-o', type=CloudPath(), default=None, help="Output the results to a different place.", show_default=True) @click.pass_context def skeleton_forge( ctx, path, queue, mip, shape, fill_missing, dust_threshold, dust_global, spatial_index, - fix_branching, fix_borders, fix_avocados, + fix_branching, fix_borders, fix_avocados, fix_autapses, fill_holes, scale, const, soma_detect, soma_accept, soma_scale, soma_const, max_paths, sharded, labels, - cross_section, + cross_section, output, ): """ (1) Synthesize skeletons from segmentation cutouts. @@ -1242,6 +1244,7 @@ def skeleton_forge( dust_global=dust_global, object_ids=labels, cross_sectional_area=(cross_section > 0), cross_sectional_area_smoothing_window=int(cross_section), + frag_path=output, fix_autapses=fix_autapses, ) enqueue_tasks(ctx, queue, tasks) From cb3e5431e383c87011497123c9f50eae8d230311 Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Thu, 23 May 2024 23:31:29 -0400 Subject: [PATCH 14/48] fix: compensate for mip levels --- igneous/tasks/skeleton.py | 9 ++++++--- 1 file changed, 6 insertions(+), 3 deletions(-) diff --git a/igneous/tasks/skeleton.py b/igneous/tasks/skeleton.py index 38db81e..5077e6c 100644 --- a/igneous/tasks/skeleton.py +++ b/igneous/tasks/skeleton.py @@ -209,8 +209,11 @@ def voxel_connectivity_graph( timestamp=self.timestamp ) + graph_chunk_size = np.array(vol.meta.graph_chunk_size) / vol.meta.downsample_ratio(vol.mip) + graph_chunk_size = graph_chunk_size.astype(int) + shape = bbox.size()[:3] - sgx, sgy, sgz = list(np.ceil(shape / vol.meta.graph_chunk_size).astype(int)) + sgx, sgy, sgz = list(np.ceil(shape / graph_chunk_size).astype(int)) vcg = np.zeros(layer_2.shape[:3], dtype=np.uint32, order="F") @@ -218,7 +221,7 @@ def voxel_connectivity_graph( for gx,gy,gz in xyzrange([sgx, sgy, sgz]): bbx = Bbox((gx,gy,gz), (gx+1, gy+1, gz+1)) - bbx *= vol.meta.graph_chunk_size + bbx *= graph_chunk_size bbx = Bbox.clamp(bbx, clamp_box) cutout = np.asfortranarray(layer_2[bbx.to_slices()][:,:,:,0]) @@ -238,7 +241,7 @@ def voxel_connectivity_graph( for gx,gy,gz in xyzrange([sgx, sgy, sgz]): bbx = Bbox((gx,gy,gz), (gx+1, gy+1, gz+1)) - bbx *= vol.meta.graph_chunk_size + bbx *= graph_chunk_size bbx = Bbox.clamp(bbx, clamp_box) slicearr = [] From 3611bc8d0783b7a4297563b802bf84cf0694dfe1 Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Thu, 23 May 2024 23:35:24 -0400 Subject: [PATCH 15/48] fix: wrong path joining --- igneous/tasks/skeleton.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/igneous/tasks/skeleton.py b/igneous/tasks/skeleton.py index 5077e6c..9290766 100644 --- a/igneous/tasks/skeleton.py +++ b/igneous/tasks/skeleton.py @@ -117,7 +117,7 @@ def execute(self): if self.frag_path is None: path = vol.meta.join(self.cloudpath, path) else: - path = CloudFiles(self.frag_path).join(path) + path = CloudFiles(self.frag_path).join(self.frag_path, path) all_labels = vol.download( bbox.to_slices(), From 3bc97555242f1251df522e86bd2274a864611d44 Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Fri, 24 May 2024 02:27:16 -0400 Subject: [PATCH 16/48] fix: ensure mesh bounds are rectified after expansion --- igneous/task_creation/skeleton.py | 1 + 1 file changed, 1 insertion(+) diff --git a/igneous/task_creation/skeleton.py b/igneous/task_creation/skeleton.py index b748d1f..5b042b6 100644 --- a/igneous/task_creation/skeleton.py +++ b/igneous/task_creation/skeleton.py @@ -58,6 +58,7 @@ def bounds_from_mesh( bbxes.append(bounds) bounds = Bbox.expand(*bbxes) + bounds = bounds.expand_to_chunk_size(shape, offset=vol.voxel_offset) return Bbox.clamp(bounds, vol.bounds) def create_skeletonizing_tasks( From e6b917453432efa01bedbbe788b37bc1c00b29ba Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Sun, 9 Jun 2024 14:28:46 -0400 Subject: [PATCH 17/48] fix: wrong variable name in shard downsample --- igneous/tasks/image/image.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/igneous/tasks/image/image.py b/igneous/tasks/image/image.py index ae5964a..2e72cd6 100755 --- a/igneous/tasks/image/image.py +++ b/igneous/tasks/image/image.py @@ -673,7 +673,7 @@ def ImageShardDownsampleTask( output_img = np.zeros(shard_shape, dtype=src_vol.dtype, order="F") nz = int(math.ceil(bbox.dz / (chunk_size.z * factor[2]))) - dsfn = downsample_method_to_fn(method, sparse, vol) + dsfn = downsample_method_to_fn(method, sparse, src_vol) zbox = bbox.clone() zbox.maxpt.z = zbox.minpt.z + (chunk_size.z * factor[2]) From 329db4bbc9966906909c0caef7caaba6518f126e Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Wed, 24 Jul 2024 00:58:58 -0400 Subject: [PATCH 18/48] fix(shards/image): fix situation where chunk.z == dataset.z --- igneous/task_creation/image.py | 30 ++++++++++++++++++++---------- 1 file changed, 20 insertions(+), 10 deletions(-) diff --git a/igneous/task_creation/image.py b/igneous/task_creation/image.py index 9e06473..94ba786 100644 --- a/igneous/task_creation/image.py +++ b/igneous/task_creation/image.py @@ -406,14 +406,6 @@ def create_sharded_image_info( # maximum amount of information in the morton codes grid_size = np.ceil(Vec(*dataset_size) / Vec(*chunk_size)).astype(np.int64) max_bits = sum([ math.ceil(math.log2(size)) for size in grid_size ]) - if max_bits > 64: - raise ValueError( - f"{max_bits}, more than a 64-bit integer, " - "would be required to describe the chunk positions " - "in this dataset. Try increasing the chunk size or " - "increasing dataset bounds." - f"Dataset Size: {dataset_size} Chunk Size: {chunk_size}" - ) chunks_per_shard = math.ceil(uncompressed_shard_bytesize / (chunk_voxels * byte_width)) chunks_per_shard = 2 ** int(math.log2(chunks_per_shard)) @@ -423,7 +415,7 @@ def create_sharded_image_info( # approximate, would need to account for rounding effects to be exact # rounding is corrected for via max_bits - pre - mini below. - num_shards = num_chunks / chunks_per_shard + num_shards = num_chunks / chunks_per_shard def update_bits(): shard_bits = int(math.ceil(math.log2(num_shards))) @@ -465,7 +457,25 @@ def update_bits(): # in the morton codes, so if there's any slack from rounding, the # remainder goes into shard bits. preshift_bits = preshift_bits - minishard_bits - shard_bits = max_bits - preshift_bits - minishard_bits + if dataset_size[2] == chunk_size[2]: + additional_bits = (preshift_bits // 3) + i = 0 + while i < additional_bits: + max_bits += 1 + preshift_bits += 1 + if preshift_bits % 3 != 0: + i += 1 + + shard_bits = max(max_bits - preshift_bits - minishard_bits, 0) + + if max_bits > 64: + raise ValueError( + f"{max_bits}, more than a 64-bit integer, " + "would be required to describe the chunk positions " + "in this dataset. Try increasing the chunk size or " + "increasing dataset bounds." + f"Dataset Size: {dataset_size} Chunk Size: {chunk_size}" + ) if preshift_bits < 0: raise ValueError(f"Preshift bits cannot be negative. ({shard_bits}, {minishard_bits}, {preshift_bits}), total info: {max_bits} bits") From 6208048848e1252eb9671c8705d7bfe7bae83ae5 Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Tue, 16 Jan 2024 14:09:30 -0500 Subject: [PATCH 19/48] wip(skeletonizing): add support for fixing autapses in graphene volumes --- igneous/task_creation/skeleton.py | 20 +++++++++++++++++--- igneous/tasks/skeleton.py | 31 ++++++++++++++++++++++++++++--- 2 files changed, 45 insertions(+), 6 deletions(-) diff --git a/igneous/task_creation/skeleton.py b/igneous/task_creation/skeleton.py index bb7d4a6..1c665fd 100644 --- a/igneous/task_creation/skeleton.py +++ b/igneous/task_creation/skeleton.py @@ -71,7 +71,7 @@ def create_skeletonizing_tasks( parallel=1, fill_missing=False, sharded=False, frag_path=None, spatial_index=True, synapses=None, num_synapses=None, - dust_global=False, + dust_global=False, fix_autapses=False, cross_sectional_area=False, cross_sectional_area_smoothing_window=5, ): @@ -121,6 +121,17 @@ def create_skeletonizing_tasks( fix_borders: Allows trivial merging of single overlap tasks. You'll only want to set this to false if you're working on single or non-overlapping volumes. + fix_autapses: Only possible for graphene volumes. Uses PyChunkGraph (PCG) information + to fix autapses (when a neuron synapses onto itself). This requires splitting + contacts between the edges of two touching voxels. The algorithm for doing this + requires much more memory. + + This works by comparing the PYC L2 and root layers. L1 is watershed. L2 is the + connections only within an atomic chunk. The root layer provides the global + connectivity. Autapses can be distinguished at the L2 level, above that, they + may not be (and certainly not at the root level). We extract the voxel connectivity + graph from L2 and perform the overall trace at root connectivity. + dust_threshold: don't skeletonize labels smaller than this number of voxels as seen by a single task. dust_global: Use global voxel counts for the dust threshold instead of from @@ -159,6 +170,9 @@ def create_skeletonizing_tasks( shape = Vec(*shape) vol = CloudVolume(cloudpath, mip=mip, info=info) + if fix_autapses and vol.meta.path.format != "graphene": + raise ValueError("fix_autapses can only be performed on graphene volumes.") + if dust_threshold > 0 and dust_global: cf = CloudFiles(cloudpath) vxctfile = cf.join(vol.key, 'stats', 'voxel_counts.mb') @@ -247,8 +261,7 @@ def task(self, shape, offset): spatial_grid_shape=shape.clone(), # used for writing index filenames synapses=bbox_synapses, dust_global=dust_global, - cross_sectional_area=bool(cross_sectional_area), - cross_sectional_area_smoothing_window=int(cross_sectional_area_smoothing_window), + ) def synapses_for_bbox(self, shape, offset): @@ -292,6 +305,7 @@ def on_finish(self): 'spatial_index': bool(spatial_index), 'synapses': bool(synapses), 'dust_global': bool(dust_global), + 'fix_autapses': bool(fix_autapses), 'cross_sectional_area': bool(cross_sectional_area), 'cross_sectional_area_smoothing_window': int(cross_sectional_area_smoothing_window), }, diff --git a/igneous/tasks/skeleton.py b/igneous/tasks/skeleton.py index cfe80f7..fb7a2b5 100644 --- a/igneous/tasks/skeleton.py +++ b/igneous/tasks/skeleton.py @@ -1,4 +1,4 @@ -from typing import Optional,Sequence,Dict +from typing import Optional, Sequence, Dict, List from functools import reduce import itertools @@ -24,6 +24,7 @@ from cloudvolume.lib import Vec, Bbox, sip from cloudvolume.datasource.precomputed.sharding import synthesize_shard_files +import cc3d import fastremap import kimimaro @@ -89,6 +90,7 @@ def __init__( cross_sectional_area_shape_delta:int = 150, dry_run:bool = False, strip_integer_attributes:bool = True, + fix_autapses:bool = False, ): super().__init__( cloudpath, shape, offset, mip, @@ -101,7 +103,8 @@ def __init__( spatial_grid_shape, synapses, bool(dust_global), bool(cross_sectional_area), int(cross_sectional_area_smoothing_window), int(cross_sectional_area_shape_delta), - bool(dry_run), bool(strip_integer_attributes) + bool(dry_run), bool(strip_integer_attributes), + bool(fix_autapses), ) self.bounds = Bbox(offset, Vec(*shape) + Vec(*offset)) self.index_bounds = Bbox(offset, Vec(*spatial_grid_shape) + Vec(*offset)) @@ -136,6 +139,10 @@ def execute(self): dust_threshold = 0 all_labels = self.apply_global_dust_threshold(vol, all_labels) + voxel_graph = None + if self.fix_autapses: + voxel_graph = self.voxel_connectivity_graph(vol, bbox) + skeletons = kimimaro.skeletonize( all_labels, self.teasar_params, object_ids=self.object_ids, @@ -148,6 +155,7 @@ def execute(self): fill_holes=self.fill_holes, parallel=self.parallel, extra_targets_after=extra_targets_after.keys(), + voxel_graph=voxel_graph, ) del all_labels @@ -189,7 +197,24 @@ def execute(self): if self.spatial_index: self.upload_spatial_index(vol, path, index_bbox, skeletons) - + + def voxel_connectivity_graph(self, vol:CloudVolume, bbox:Bbox) -> np.ndarray: + layer_2 = vol.download(bbox, stop_layer=2, agglomerate=True) + + sgx, sgy, sgz = list(np.ceil(bbox.size3() / vol.meta.graph_chunk_size).astype(int)) + + vcg = np.zeros(layer_2.shape, dtype=np.uint32, order="F") + + for gx,gy,gz in xyzrange([sgx, sgy, sgz]): + bbx = Bbox((gx,gy,gz), (gx+1, gy+1, gz+1)) + bbx *= vol.meta.graph_chunk_size + + cutout = np.asfortranarray(layer_2[bbx.to_slices()]) + vcg_cutout = cc3d.voxel_connectivity_graph(cutout, connectivity=26) + vcg[bbx.to_slices()] = vcg_cutout + + return vcg + def compute_cross_sectional_area(self, vol, bbox, skeletons): # Why redownload a bigger image? In order to avoid clipping the # cross sectional areas on the edges. From 894b2fae5a7c1e29458f6a1df627cb1fc52be8a3 Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Tue, 16 Jan 2024 14:14:41 -0500 Subject: [PATCH 20/48] feat: add support for graphene timestamp --- igneous/task_creation/skeleton.py | 6 +++++- igneous/tasks/skeleton.py | 5 ++++- 2 files changed, 9 insertions(+), 2 deletions(-) diff --git a/igneous/task_creation/skeleton.py b/igneous/task_creation/skeleton.py index 1c665fd..4c0b8fe 100644 --- a/igneous/task_creation/skeleton.py +++ b/igneous/task_creation/skeleton.py @@ -74,6 +74,7 @@ def create_skeletonizing_tasks( dust_global=False, fix_autapses=False, cross_sectional_area=False, cross_sectional_area_smoothing_window=5, + timestamp=None, ): """ Assign tasks with one voxel overlap in a regular grid @@ -166,6 +167,7 @@ def create_skeletonizing_tasks( to the total computation.) cross_sectional_area_smoothing_window: Perform a rolling average of the normal vectors across these many vectors. + timestamp: for graphene volumes only, you can specify the timepoint to use """ shape = Vec(*shape) vol = CloudVolume(cloudpath, mip=mip, info=info) @@ -261,7 +263,8 @@ def task(self, shape, offset): spatial_grid_shape=shape.clone(), # used for writing index filenames synapses=bbox_synapses, dust_global=dust_global, - + fix_autapses=bool(fix_autapses), + timestamp=timestamp, ) def synapses_for_bbox(self, shape, offset): @@ -306,6 +309,7 @@ def on_finish(self): 'synapses': bool(synapses), 'dust_global': bool(dust_global), 'fix_autapses': bool(fix_autapses), + 'timestamp': timestamp, 'cross_sectional_area': bool(cross_sectional_area), 'cross_sectional_area_smoothing_window': int(cross_sectional_area_smoothing_window), }, diff --git a/igneous/tasks/skeleton.py b/igneous/tasks/skeleton.py index fb7a2b5..e1354d7 100644 --- a/igneous/tasks/skeleton.py +++ b/igneous/tasks/skeleton.py @@ -91,6 +91,7 @@ def __init__( dry_run:bool = False, strip_integer_attributes:bool = True, fix_autapses:bool = False, + timestamp:Optional[int] = None, ): super().__init__( cloudpath, shape, offset, mip, @@ -104,7 +105,7 @@ def __init__( bool(cross_sectional_area), int(cross_sectional_area_smoothing_window), int(cross_sectional_area_shape_delta), bool(dry_run), bool(strip_integer_attributes), - bool(fix_autapses), + bool(fix_autapses), timestamp, ) self.bounds = Bbox(offset, Vec(*shape) + Vec(*offset)) self.index_bounds = Bbox(offset, Vec(*spatial_grid_shape) + Vec(*offset)) @@ -115,6 +116,8 @@ def execute(self): info=self.info, cdn_cache=False, parallel=self.parallel, fill_missing=self.fill_missing, + agglomerate=True, + timestamp=timestamp, ) bbox = Bbox.clamp(self.bounds, vol.bounds) index_bbox = Bbox.clamp(self.index_bounds, vol.bounds) From 317612879c1d6cb7b0d1d28f238cce35607c2103 Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Wed, 17 Jan 2024 02:21:18 -0500 Subject: [PATCH 21/48] wip: ensure shape matches atomic chunks --- igneous/task_creation/skeleton.py | 11 +++++++++-- igneous/tasks/skeleton.py | 4 +++- 2 files changed, 12 insertions(+), 3 deletions(-) diff --git a/igneous/task_creation/skeleton.py b/igneous/task_creation/skeleton.py index 4c0b8fe..1f03325 100644 --- a/igneous/task_creation/skeleton.py +++ b/igneous/task_creation/skeleton.py @@ -172,8 +172,15 @@ def create_skeletonizing_tasks( shape = Vec(*shape) vol = CloudVolume(cloudpath, mip=mip, info=info) - if fix_autapses and vol.meta.path.format != "graphene": - raise ValueError("fix_autapses can only be performed on graphene volumes.") + if fix_autapses: + if vol.meta.path.format != "graphene": + raise ValueError("fix_autapses can only be performed on graphene volumes.") + + if not np.all(shape % vol.meta.graph_chunk_size == 0): + raise ValueError( + f"shape must be a multiple of the graph chunk size. Got: {shape}, " + f"{vol.meta.graph_chunk_size}" + ) if dust_threshold > 0 and dust_global: cf = CloudFiles(cloudpath) diff --git a/igneous/tasks/skeleton.py b/igneous/tasks/skeleton.py index e1354d7..2a110e5 100644 --- a/igneous/tasks/skeleton.py +++ b/igneous/tasks/skeleton.py @@ -204,13 +204,15 @@ def execute(self): def voxel_connectivity_graph(self, vol:CloudVolume, bbox:Bbox) -> np.ndarray: layer_2 = vol.download(bbox, stop_layer=2, agglomerate=True) - sgx, sgy, sgz = list(np.ceil(bbox.size3() / vol.meta.graph_chunk_size).astype(int)) + shape = bbox.size()[:3] + sgx, sgy, sgz = list(np.ceil(shape / vol.meta.graph_chunk_size).astype(int)) vcg = np.zeros(layer_2.shape, dtype=np.uint32, order="F") for gx,gy,gz in xyzrange([sgx, sgy, sgz]): bbx = Bbox((gx,gy,gz), (gx+1, gy+1, gz+1)) bbx *= vol.meta.graph_chunk_size + bbx = Bbox.clamp(bbx, (0,0,0), shape) cutout = np.asfortranarray(layer_2[bbx.to_slices()]) vcg_cutout = cc3d.voxel_connectivity_graph(cutout, connectivity=26) From 91bed59b5820f7ee8f8641a7d73c269bd83267e8 Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Mon, 13 May 2024 19:32:47 -0400 Subject: [PATCH 22/48] fix: add back in xyzrange --- igneous/tasks/skeleton.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/igneous/tasks/skeleton.py b/igneous/tasks/skeleton.py index 2a110e5..26f4c4a 100644 --- a/igneous/tasks/skeleton.py +++ b/igneous/tasks/skeleton.py @@ -21,7 +21,7 @@ import cloudvolume from cloudvolume import CloudVolume, Skeleton, paths -from cloudvolume.lib import Vec, Bbox, sip +from cloudvolume.lib import Vec, Bbox, sip, xyzrange from cloudvolume.datasource.precomputed.sharding import synthesize_shard_files import cc3d From bc215cea6441dea476a33b62a156e078388c75f6 Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Tue, 14 May 2024 14:26:32 -0400 Subject: [PATCH 23/48] fix(SkeletonTask): generalize frag_path for windows Also interpret bare paths as "file://" --- igneous/tasks/skeleton.py | 12 +++++++----- 1 file changed, 7 insertions(+), 5 deletions(-) diff --git a/igneous/tasks/skeleton.py b/igneous/tasks/skeleton.py index 26f4c4a..febe1f5 100644 --- a/igneous/tasks/skeleton.py +++ b/igneous/tasks/skeleton.py @@ -107,6 +107,8 @@ def __init__( bool(dry_run), bool(strip_integer_attributes), bool(fix_autapses), timestamp, ) + if isinstance(self.frag_path, str): + self.frag_path = cloudfiles.paths.normalize(self.frag_path) self.bounds = Bbox(offset, Vec(*shape) + Vec(*offset)) self.index_bounds = Bbox(offset, Vec(*spatial_grid_shape) + Vec(*offset)) @@ -123,7 +125,10 @@ def execute(self): index_bbox = Bbox.clamp(self.index_bounds, vol.bounds) path = skeldir(self.cloudpath) - path = os.path.join(self.cloudpath, path) + if self.frag_path is None: + path = vol.meta.join(self.cloudpath, path) + else: + path = CloudFiles(self.frag_path).join(path) all_labels = vol[ bbox.to_slices() ] all_labels = all_labels[:,:,:,0] @@ -191,10 +196,7 @@ def execute(self): return skeletons if self.sharded: - if self.frag_path: - self.upload_batch(vol, os.path.join(self.frag_path, skeldir(self.cloudpath)), index_bbox, skeletons) - else: - self.upload_batch(vol, path, index_bbox, skeletons) + self.upload_batch(vol, path, index_bbox, skeletons) else: self.upload_individuals(vol, path, bbox, skeletons) From f4a513fd94fbba0630bf344f02370633bc0a06e4 Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Tue, 14 May 2024 14:27:06 -0400 Subject: [PATCH 24/48] fix(SkeletonTask): use agglomerate and timestamp correctly --- igneous/tasks/skeleton.py | 8 +++++--- 1 file changed, 5 insertions(+), 3 deletions(-) diff --git a/igneous/tasks/skeleton.py b/igneous/tasks/skeleton.py index febe1f5..a76acfc 100644 --- a/igneous/tasks/skeleton.py +++ b/igneous/tasks/skeleton.py @@ -118,8 +118,6 @@ def execute(self): info=self.info, cdn_cache=False, parallel=self.parallel, fill_missing=self.fill_missing, - agglomerate=True, - timestamp=timestamp, ) bbox = Bbox.clamp(self.bounds, vol.bounds) index_bbox = Bbox.clamp(self.index_bounds, vol.bounds) @@ -130,7 +128,11 @@ def execute(self): else: path = CloudFiles(self.frag_path).join(path) - all_labels = vol[ bbox.to_slices() ] + all_labels = vol.download( + bbox.to_slices(), + agglomerate=True, + timestamp=self.timestamp + ) all_labels = all_labels[:,:,:,0] if self.mask_ids: From 2cfed740d63da067d3eefc2328b7a9aad8a5ec5b Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Tue, 14 May 2024 14:38:37 -0400 Subject: [PATCH 25/48] fix: merge error removed cross_sectional_area params --- igneous/task_creation/skeleton.py | 2 ++ 1 file changed, 2 insertions(+) diff --git a/igneous/task_creation/skeleton.py b/igneous/task_creation/skeleton.py index 1f03325..b748d1f 100644 --- a/igneous/task_creation/skeleton.py +++ b/igneous/task_creation/skeleton.py @@ -272,6 +272,8 @@ def task(self, shape, offset): dust_global=dust_global, fix_autapses=bool(fix_autapses), timestamp=timestamp, + cross_sectional_area=bool(cross_sectional_area), + cross_sectional_area_smoothing_window=int(cross_sectional_area_smoothing_window), ) def synapses_for_bbox(self, shape, offset): From 548fe0b5d3840e02a2bea5ae4e9eff5b3fd98180 Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Fri, 17 May 2024 23:41:30 -0400 Subject: [PATCH 26/48] fix(SkeletonTask): apply timestamp in correct location --- igneous/tasks/skeleton.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/igneous/tasks/skeleton.py b/igneous/tasks/skeleton.py index a76acfc..2cd61d7 100644 --- a/igneous/tasks/skeleton.py +++ b/igneous/tasks/skeleton.py @@ -206,7 +206,7 @@ def execute(self): self.upload_spatial_index(vol, path, index_bbox, skeletons) def voxel_connectivity_graph(self, vol:CloudVolume, bbox:Bbox) -> np.ndarray: - layer_2 = vol.download(bbox, stop_layer=2, agglomerate=True) + layer_2 = vol.download(bbox, stop_layer=2, agglomerate=True, timestamp=self.timestamp) shape = bbox.size()[:3] sgx, sgy, sgz = list(np.ceil(shape / vol.meta.graph_chunk_size).astype(int)) From 2fb01ac2b6a14668b4299b1407c4f48e80e3819c Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Wed, 22 May 2024 02:08:25 -0400 Subject: [PATCH 27/48] fix: fixup the edges of the vcg with root vcg --- igneous/tasks/skeleton.py | 45 ++++++++++++++++++++++++++++++++++++--- 1 file changed, 42 insertions(+), 3 deletions(-) diff --git a/igneous/tasks/skeleton.py b/igneous/tasks/skeleton.py index 2cd61d7..b4da274 100644 --- a/igneous/tasks/skeleton.py +++ b/igneous/tasks/skeleton.py @@ -151,7 +151,7 @@ def execute(self): voxel_graph = None if self.fix_autapses: - voxel_graph = self.voxel_connectivity_graph(vol, bbox) + voxel_graph = self.voxel_connectivity_graph(vol, bbox, all_labels) skeletons = kimimaro.skeletonize( all_labels, self.teasar_params, @@ -205,8 +205,18 @@ def execute(self): if self.spatial_index: self.upload_spatial_index(vol, path, index_bbox, skeletons) - def voxel_connectivity_graph(self, vol:CloudVolume, bbox:Bbox) -> np.ndarray: - layer_2 = vol.download(bbox, stop_layer=2, agglomerate=True, timestamp=self.timestamp) + def voxel_connectivity_graph( + self, + vol:CloudVolume, + bbox:Bbox, + root_labels:np.ndarray, + ) -> np.ndarray: + layer_2 = vol.download( + bbox, + stop_layer=2, + agglomerate=True, + timestamp=self.timestamp + ) shape = bbox.size()[:3] sgx, sgy, sgz = list(np.ceil(shape / vol.meta.graph_chunk_size).astype(int)) @@ -221,6 +231,35 @@ def voxel_connectivity_graph(self, vol:CloudVolume, bbox:Bbox) -> np.ndarray: cutout = np.asfortranarray(layer_2[bbx.to_slices()]) vcg_cutout = cc3d.voxel_connectivity_graph(cutout, connectivity=26) vcg[bbx.to_slices()] = vcg_cutout + del vcg_cutout + + del layer_2 + + # the proper way to do this would be to get the lowest the L3..LN root + # as needed, but the lazy way to do this is to get the root labels + # which will retain a few errors, but overall the error rate should be + # over 100x less. We need to shade in the sides of the connectivity graph + # with edges that represent the connections between the adjacent boxes. + + root_vcg = cc3d.voxel_connectivity_graph(root_labels, connectivity=26) + + for gx,gy,gz in xyzrange([sgx, sgy, sgz]): + bbx = Bbox((gx,gy,gz), (gx+1, gy+1, gz+1)) + bbx *= vol.meta.graph_chunk_size + bbx = Bbox.clamp(bbx, (0,0,0), shape) + + slicearr = [] + for i in range(3): + bbx1 = bbx.clone() + bbx1.maxpt[i] = bbx1.minpt[i] + 1 + slicearr.append(bbx1) + + bbx1 = bbx.clone() + bbx1.minpt[i] = bbx1.maxpt[i] - 1 + slicearr.append(bbx1) + + for bbx1 in slicearr: + vcg[bbx1.to_slices()] = root_vcg[bbx1.to_slices()] return vcg From 521c9d29e06fe39e2c7354e4a8ae376fc41c015d Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Thu, 23 May 2024 20:33:25 -0400 Subject: [PATCH 28/48] feat: remove requirements for a queue parameter --- igneous_cli/cli.py | 22 +++++++++++----------- 1 file changed, 11 insertions(+), 11 deletions(-) diff --git a/igneous_cli/cli.py b/igneous_cli/cli.py index 4800575..a7862e8 100644 --- a/igneous_cli/cli.py +++ b/igneous_cli/cli.py @@ -896,7 +896,7 @@ def meshgroup(): @meshgroup.command("xfer") @click.argument("src", type=CloudPath()) @click.argument("dest", type=CloudPath()) -@click.option('--queue', required=True, help="AWS SQS queue or directory to be used for a task queue. e.g. sqs://my-queue or ./my-queue. See https://github.com/seung-lab/python-task-queue", type=str) +@click.option('--queue', help="AWS SQS queue or directory to be used for a task queue. e.g. sqs://my-queue or ./my-queue. See https://github.com/seung-lab/python-task-queue", type=str) @click.option("--sharded", is_flag=True, default=False, help="Generate shard fragments instead of outputing mesh fragments.", show_default=True) @click.option("--dir", "mesh_dir", type=str, default=None, help="Write meshes into this directory instead of the one indicated in the info file.") @click.option('--magnitude', default=2, help="Split up the work with 10^(magnitude) prefix based tasks.", show_default=True) @@ -933,7 +933,7 @@ def mesh_xfer( @meshgroup.command("forge") @click.argument("path", type=CloudPath()) -@click.option('--queue', required=True, help="AWS SQS queue or directory to be used for a task queue. e.g. sqs://my-queue or ./my-queue. See https://github.com/seung-lab/python-task-queue", type=str) +@click.option('--queue', help="AWS SQS queue or directory to be used for a task queue. e.g. sqs://my-queue or ./my-queue. See https://github.com/seung-lab/python-task-queue", type=str) @click.option('--mip', default=0, help="Perform meshing using this level of the image pyramid.", show_default=True) @click.option('--shape', type=Tuple3(), default=(448, 448, 448), help="Set the task shape in voxels.", show_default=True) @click.option('--simplify/--skip-simplify', is_flag=True, default=True, help="Enable mesh simplification.", show_default=True) @@ -984,7 +984,7 @@ def mesh_forge( @meshgroup.command("merge") @click.argument("path", type=CloudPath()) -@click.option('--queue', required=True, help="AWS SQS queue or directory to be used for a task queue. e.g. sqs://my-queue or ./my-queue. See https://github.com/seung-lab/python-task-queue", type=str) +@click.option('--queue', help="AWS SQS queue or directory to be used for a task queue. e.g. sqs://my-queue or ./my-queue. See https://github.com/seung-lab/python-task-queue", type=str) @click.option('--magnitude', default=2, help="Split up the work with 10^(magnitude) prefix based tasks. Default: 2 (100 tasks)") @click.option('--nlod', default=0, help="(multires) How many extra levels of detail to create.", show_default=True) @click.option('--vqb', default=16, help="(multires) Vertex quantization bits for stored model representation. 10 or 16 only.", show_default=True) @@ -1016,7 +1016,7 @@ def mesh_merge(ctx, path, queue, magnitude, nlod, vqb, dir, min_chunk_size): @meshgroup.command("merge-sharded") @click.argument("path", type=CloudPath()) -@click.option('--queue', required=True, help="AWS SQS queue or directory to be used for a task queue. e.g. sqs://my-queue or ./my-queue. See https://github.com/seung-lab/python-task-queue", type=str) +@click.option('--queue', help="AWS SQS queue or directory to be used for a task queue. e.g. sqs://my-queue or ./my-queue. See https://github.com/seung-lab/python-task-queue", type=str) @click.option('--nlod', default=1, help="Number of levels of detail to create.", type=int, show_default=True) @click.option('--vqb', default=16, help="Vertex quantization bits. Can be 10 or 16.", type=int, show_default=True) @click.option('--compress-level', default=7, help="Draco compression level.", type=int, show_default=True) @@ -1065,7 +1065,7 @@ def mesh_sharded_merge( @meshgroup.command("rm") @click.argument("path", type=CloudPath()) -@click.option('--queue', required=True, help="AWS SQS queue or directory to be used for a task queue. e.g. sqs://my-queue or ./my-queue. See https://github.com/seung-lab/python-task-queue", type=str) +@click.option('--queue', help="AWS SQS queue or directory to be used for a task queue. e.g. sqs://my-queue or ./my-queue. See https://github.com/seung-lab/python-task-queue", type=str) @click.option('--magnitude', default=2, help="Split up the work with 10^(magnitude) prefix based tasks. Default: 2 (100 tasks)") @click.option('--dir', 'mesh_dir', default=None, help="Target this directory instead of the one indicated in the info file.") @click.pass_context @@ -1114,7 +1114,7 @@ def spatialindexgroup(): @spatialindexgroup.command("create") @click.argument("path", type=CloudPath()) -@click.option('--queue', required=True, help="AWS SQS queue or directory to be used for a task queue. e.g. sqs://my-queue or ./my-queue. See https://github.com/seung-lab/python-task-queue", type=str) +@click.option('--queue', help="AWS SQS queue or directory to be used for a task queue. e.g. sqs://my-queue or ./my-queue. See https://github.com/seung-lab/python-task-queue", type=str) @click.option('--shape', default="448,448,448", type=Tuple3(), help="Shape in voxels of each indexing task.", show_default=True) @click.option('--mip', default=0, help="Perform indexing using this level of the image pyramid.", show_default=True) @click.option('--fill-missing', is_flag=True, default=False, help="Interpret missing image files as background instead of failing.", show_default=True) @@ -1283,7 +1283,7 @@ def skeleton_merge( @skeletongroup.command("merge-sharded") @click.argument("path", type=CloudPath()) -@click.option('--queue', required=True, help="AWS SQS queue or directory to be used for a task queue. e.g. sqs://my-queue or ./my-queue. See https://github.com/seung-lab/python-task-queue", type=str) +@click.option('--queue', help="AWS SQS queue or directory to be used for a task queue. e.g. sqs://my-queue or ./my-queue. See https://github.com/seung-lab/python-task-queue", type=str) @click.option('--min-cable-length', default=1000, help="Skip objects smaller than this physical path length.", type=float, show_default=True) @click.option('--max-cable-length', default=None, help="Skip objects larger than this physical path length. Default: no limit", type=float) @click.option('--tick-threshold', default=0, help="Remove small \"ticks\", or branches from the main skeleton one at a time from smallest to largest. Branches larger than this are preserved. Default: no elimination", type=float) @@ -1334,7 +1334,7 @@ def skeleton_sharded_merge( @imagegroup.command("rm") @click.argument("path", type=CloudPath()) -@click.option('--queue', required=True, help="AWS SQS queue or directory to be used for a task queue. e.g. sqs://my-queue or ./my-queue. See https://github.com/seung-lab/python-task-queue", type=str) +@click.option('--queue', help="AWS SQS queue or directory to be used for a task queue. e.g. sqs://my-queue or ./my-queue. See https://github.com/seung-lab/python-task-queue", type=str) @click.option('--mip', default=0, help="Which mip level to start deleting from. Default: 0") @click.option('--num-mips', default=5, help="The number of mip levels to delete at once. Default: 5") @click.option('--shape', default=None, help="The size of each deletion task as a comma delimited list. Must be a multiple of the chunk size.", type=Tuple3()) @@ -1358,7 +1358,7 @@ def delete_images( @skeletongroup.command("rm") @click.argument("path", type=CloudPath()) -@click.option('--queue', required=True, help="AWS SQS queue or directory to be used for a task queue. e.g. sqs://my-queue or ./my-queue. See https://github.com/seung-lab/python-task-queue", type=str) +@click.option('--queue', help="AWS SQS queue or directory to be used for a task queue. e.g. sqs://my-queue or ./my-queue. See https://github.com/seung-lab/python-task-queue", type=str) @click.option('--magnitude', default=2, help="Split up the work with 10^(magnitude) prefix based tasks. Default: 2 (100 tasks)") @click.option('--dir', 'skel_dir', default=None, help="Target this directory instead of the one indicated in the info file.") @click.pass_context @@ -1376,7 +1376,7 @@ def skeleton_rm(ctx, path, queue, magnitude, skel_dir): @skeletongroup.command("xfer") @click.argument("src", type=CloudPath()) @click.argument("dest", type=CloudPath()) -@click.option('--queue', required=True, help="AWS SQS queue or directory to be used for a task queue. e.g. sqs://my-queue or ./my-queue. See https://github.com/seung-lab/python-task-queue", type=str) +@click.option('--queue', help="AWS SQS queue or directory to be used for a task queue. e.g. sqs://my-queue or ./my-queue. See https://github.com/seung-lab/python-task-queue", type=str) @click.option("--sharded", is_flag=True, default=False, help="Generate shard fragments instead of outputing mesh fragments.", show_default=True) @click.option("--dir", "skel_dir", type=str, default=None, help="Write skeletons into this directory instead of the one indicated in the info file.") @click.option('--magnitude', default=2, help="Split up the work with 10^(magnitude) prefix based tasks.", show_default=True) @@ -1440,7 +1440,7 @@ def spatialindexgroupskel(): @spatialindexgroupskel.command("create") @click.argument("path", type=CloudPath()) -@click.option('--queue', required=True, help="AWS SQS queue or directory to be used for a task queue. e.g. sqs://my-queue or ./my-queue. See https://github.com/seung-lab/python-task-queue", type=str) +@click.option('--queue', help="AWS SQS queue or directory to be used for a task queue. e.g. sqs://my-queue or ./my-queue. See https://github.com/seung-lab/python-task-queue", type=str) @click.option('--shape', default="512,512,512", type=Tuple3(), help="Shape in voxels of each indexing task.", show_default=True) @click.option('--mip', default=0, help="Perform indexing using this level of the image pyramid.", show_default=True) @click.option('--fill-missing', is_flag=True, default=False, help="Interpret missing image files as background instead of failing.", show_default=True) From f50059bd3e2a4948a1979f8e458ee001b7aced09 Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Thu, 23 May 2024 22:33:48 -0400 Subject: [PATCH 29/48] redesign: use vol.info to get skeleton path --- igneous/tasks/skeleton.py | 11 +---------- 1 file changed, 1 insertion(+), 10 deletions(-) diff --git a/igneous/tasks/skeleton.py b/igneous/tasks/skeleton.py index b4da274..2634d40 100644 --- a/igneous/tasks/skeleton.py +++ b/igneous/tasks/skeleton.py @@ -40,15 +40,6 @@ def filename_to_segid(filename): segid, = matches.groups() return int(segid) -def skeldir(cloudpath): - cf = CloudFiles(cloudpath) - info = cf.get_json('info') - - skel_dir = 'skeletons/' - if 'skeletons' in info: - skel_dir = info['skeletons'] - return skel_dir - def strip_integer_attributes(skeletons): for skel in skeletons: skel.extra_attributes = [ @@ -122,7 +113,7 @@ def execute(self): bbox = Bbox.clamp(self.bounds, vol.bounds) index_bbox = Bbox.clamp(self.index_bounds, vol.bounds) - path = skeldir(self.cloudpath) + path = vol.info.get("skeletons", "skeletons") if self.frag_path is None: path = vol.meta.join(self.cloudpath, path) else: From 6f5e0944076aae095f93ba78f8afe71a6d70c86f Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Thu, 23 May 2024 22:34:28 -0400 Subject: [PATCH 30/48] fix: several errors --- igneous/tasks/skeleton.py | 10 ++++++---- 1 file changed, 6 insertions(+), 4 deletions(-) diff --git a/igneous/tasks/skeleton.py b/igneous/tasks/skeleton.py index 2634d40..38db81e 100644 --- a/igneous/tasks/skeleton.py +++ b/igneous/tasks/skeleton.py @@ -212,14 +212,16 @@ def voxel_connectivity_graph( shape = bbox.size()[:3] sgx, sgy, sgz = list(np.ceil(shape / vol.meta.graph_chunk_size).astype(int)) - vcg = np.zeros(layer_2.shape, dtype=np.uint32, order="F") + vcg = np.zeros(layer_2.shape[:3], dtype=np.uint32, order="F") + + clamp_box = Bbox([0,0,0], shape) for gx,gy,gz in xyzrange([sgx, sgy, sgz]): bbx = Bbox((gx,gy,gz), (gx+1, gy+1, gz+1)) bbx *= vol.meta.graph_chunk_size - bbx = Bbox.clamp(bbx, (0,0,0), shape) + bbx = Bbox.clamp(bbx, clamp_box) - cutout = np.asfortranarray(layer_2[bbx.to_slices()]) + cutout = np.asfortranarray(layer_2[bbx.to_slices()][:,:,:,0]) vcg_cutout = cc3d.voxel_connectivity_graph(cutout, connectivity=26) vcg[bbx.to_slices()] = vcg_cutout del vcg_cutout @@ -237,7 +239,7 @@ def voxel_connectivity_graph( for gx,gy,gz in xyzrange([sgx, sgy, sgz]): bbx = Bbox((gx,gy,gz), (gx+1, gy+1, gz+1)) bbx *= vol.meta.graph_chunk_size - bbx = Bbox.clamp(bbx, (0,0,0), shape) + bbx = Bbox.clamp(bbx, clamp_box) slicearr = [] for i in range(3): From 80ed2c7033a9c1ec28ef6b2b12e0e86444413e39 Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Thu, 23 May 2024 22:34:47 -0400 Subject: [PATCH 31/48] feat: add fix_autapses to cli --- igneous_cli/cli.py | 7 +++++-- 1 file changed, 5 insertions(+), 2 deletions(-) diff --git a/igneous_cli/cli.py b/igneous_cli/cli.py index a7862e8..29224ba 100644 --- a/igneous_cli/cli.py +++ b/igneous_cli/cli.py @@ -1181,6 +1181,7 @@ def skeletongroup(): @click.option('--fix-branching', is_flag=True, default=True, help="Trades speed for quality of branching at forks.", show_default=True) @click.option('--fix-borders', is_flag=True, default=True, help="Allows trivial merging of single voxel overlap tasks. Only switch off for datasets that fit in a single task.", show_default=True) @click.option('--fix-avocados', is_flag=True, default=False, help="Fixes somata where nuclei and cytoplasm have separate segmentations.", show_default=True) +@click.option('--fix-autapses', is_flag=True, default=False, help="(graphene only) Fixes autapses by using the PyChunkGraph.", show_default=True) @click.option('--fill-holes', is_flag=True, default=False, help="Preprocess each cutout to eliminate background holes and holes caused by entirely contained inclusions. Warning: May remove labels that are considered inclusions.", show_default=True) @click.option('--dust-threshold', default=1000, help="Skip skeletonizing objects smaller than this number of voxels within a cutout.", type=int, show_default=True) @click.option('--dust-global/--dust-local', is_flag=True, default=False, help="Use global voxel counts for the dust threshold (when >0). To use this feature you must first compute the global voxel counts using the 'igneous image voxels' command.", show_default=True) @@ -1195,14 +1196,15 @@ def skeletongroup(): @click.option('--sharded', is_flag=True, default=False, help="Generate shard fragments instead of outputing skeleton fragments.", show_default=True) @click.option('--labels', type=TupleN(), default=None, help="Skeletonize only this comma separated list of labels.", show_default=True) @click.option('--cross-section', type=int, default=0, help="Compute the cross sectional area for each skeleton vertex. May add substantial computation time. Integer value is the normal vector rolling average smoothing window over vertices. 0 means off.", show_default=True) +@click.option('--output', '-o', type=CloudPath(), default=None, help="Output the results to a different place.", show_default=True) @click.pass_context def skeleton_forge( ctx, path, queue, mip, shape, fill_missing, dust_threshold, dust_global, spatial_index, - fix_branching, fix_borders, fix_avocados, + fix_branching, fix_borders, fix_avocados, fix_autapses, fill_holes, scale, const, soma_detect, soma_accept, soma_scale, soma_const, max_paths, sharded, labels, - cross_section, + cross_section, output, ): """ (1) Synthesize skeletons from segmentation cutouts. @@ -1247,6 +1249,7 @@ def skeleton_forge( dust_global=dust_global, object_ids=labels, cross_sectional_area=(cross_section > 0), cross_sectional_area_smoothing_window=int(cross_section), + frag_path=output, fix_autapses=fix_autapses, ) enqueue_tasks(ctx, queue, tasks) From 37e9dc68d33302f6d44b7c64916a74de4dac561f Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Thu, 23 May 2024 23:31:29 -0400 Subject: [PATCH 32/48] fix: compensate for mip levels --- igneous/tasks/skeleton.py | 9 ++++++--- 1 file changed, 6 insertions(+), 3 deletions(-) diff --git a/igneous/tasks/skeleton.py b/igneous/tasks/skeleton.py index 38db81e..5077e6c 100644 --- a/igneous/tasks/skeleton.py +++ b/igneous/tasks/skeleton.py @@ -209,8 +209,11 @@ def voxel_connectivity_graph( timestamp=self.timestamp ) + graph_chunk_size = np.array(vol.meta.graph_chunk_size) / vol.meta.downsample_ratio(vol.mip) + graph_chunk_size = graph_chunk_size.astype(int) + shape = bbox.size()[:3] - sgx, sgy, sgz = list(np.ceil(shape / vol.meta.graph_chunk_size).astype(int)) + sgx, sgy, sgz = list(np.ceil(shape / graph_chunk_size).astype(int)) vcg = np.zeros(layer_2.shape[:3], dtype=np.uint32, order="F") @@ -218,7 +221,7 @@ def voxel_connectivity_graph( for gx,gy,gz in xyzrange([sgx, sgy, sgz]): bbx = Bbox((gx,gy,gz), (gx+1, gy+1, gz+1)) - bbx *= vol.meta.graph_chunk_size + bbx *= graph_chunk_size bbx = Bbox.clamp(bbx, clamp_box) cutout = np.asfortranarray(layer_2[bbx.to_slices()][:,:,:,0]) @@ -238,7 +241,7 @@ def voxel_connectivity_graph( for gx,gy,gz in xyzrange([sgx, sgy, sgz]): bbx = Bbox((gx,gy,gz), (gx+1, gy+1, gz+1)) - bbx *= vol.meta.graph_chunk_size + bbx *= graph_chunk_size bbx = Bbox.clamp(bbx, clamp_box) slicearr = [] From fe9267fe4ee7ca78cb508f8c027583db148f33bd Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Thu, 23 May 2024 23:35:24 -0400 Subject: [PATCH 33/48] fix: wrong path joining --- igneous/tasks/skeleton.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/igneous/tasks/skeleton.py b/igneous/tasks/skeleton.py index 5077e6c..9290766 100644 --- a/igneous/tasks/skeleton.py +++ b/igneous/tasks/skeleton.py @@ -117,7 +117,7 @@ def execute(self): if self.frag_path is None: path = vol.meta.join(self.cloudpath, path) else: - path = CloudFiles(self.frag_path).join(path) + path = CloudFiles(self.frag_path).join(self.frag_path, path) all_labels = vol.download( bbox.to_slices(), From e7809da3fd5243876b85d46179dce551eb3dddf3 Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Fri, 24 May 2024 02:27:16 -0400 Subject: [PATCH 34/48] fix: ensure mesh bounds are rectified after expansion --- igneous/task_creation/skeleton.py | 1 + 1 file changed, 1 insertion(+) diff --git a/igneous/task_creation/skeleton.py b/igneous/task_creation/skeleton.py index b748d1f..5b042b6 100644 --- a/igneous/task_creation/skeleton.py +++ b/igneous/task_creation/skeleton.py @@ -58,6 +58,7 @@ def bounds_from_mesh( bbxes.append(bounds) bounds = Bbox.expand(*bbxes) + bounds = bounds.expand_to_chunk_size(shape, offset=vol.voxel_offset) return Bbox.clamp(bounds, vol.bounds) def create_skeletonizing_tasks( From 67e18cc6d6039a9d516ea5ed5d9ce4675c35f3da Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Mon, 29 Jul 2024 03:41:13 -0400 Subject: [PATCH 35/48] fix: don't import view from CloudVolume --- igneous/tasks/mesh/mesh.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/igneous/tasks/mesh/mesh.py b/igneous/tasks/mesh/mesh.py index 9f5878a..7502e14 100644 --- a/igneous/tasks/mesh/mesh.py +++ b/igneous/tasks/mesh/mesh.py @@ -12,7 +12,7 @@ from cloudfiles import CloudFiles, CloudFile import cloudfiles.paths -from cloudvolume import CloudVolume, view +from cloudvolume import CloudVolume from cloudvolume.lib import Vec, Bbox, jsonify import mapbuffer from mapbuffer import MapBuffer, IntMap From d8da603318d19c9c7240536df9f09b0313385a69 Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Mon, 29 Jul 2024 03:42:30 -0400 Subject: [PATCH 36/48] fix: don't import view from cloudvolume --- igneous/tasks/mesh/multires.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/igneous/tasks/mesh/multires.py b/igneous/tasks/mesh/multires.py index bc56136..b44423e 100644 --- a/igneous/tasks/mesh/multires.py +++ b/igneous/tasks/mesh/multires.py @@ -16,7 +16,7 @@ from cloudfiles import CloudFiles, CloudFile -from cloudvolume import CloudVolume, Mesh, view, paths +from cloudvolume import CloudVolume, Mesh, paths from cloudvolume.lib import Vec, Bbox, jsonify, sip, toiter, first from cloudvolume.datasource.precomputed.mesh.multilod \ import MultiLevelPrecomputedMeshManifest, to_stored_model_space From 9f9c0b096a8c3e7f109cee67195ec4ed8d88cfbe Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Mon, 29 Jul 2024 20:01:49 -0400 Subject: [PATCH 37/48] fix: intelligent frag path resolution --- igneous/tasks/skeleton.py | 9 ++++++++- 1 file changed, 8 insertions(+), 1 deletion(-) diff --git a/igneous/tasks/skeleton.py b/igneous/tasks/skeleton.py index 9290766..4faa5b5 100644 --- a/igneous/tasks/skeleton.py +++ b/igneous/tasks/skeleton.py @@ -117,7 +117,14 @@ def execute(self): if self.frag_path is None: path = vol.meta.join(self.cloudpath, path) else: - path = CloudFiles(self.frag_path).join(self.frag_path, path) + # if the path is to a volume root, follow the info instructions, + # otherwise place the files exactly where frag path says to + test_path = CloudFiles(self.frag_path).join(self.frag_path, "info") + test_info = CloudFile(test_path).get_json() + if test_info is None or 'scales' in test_info: + path = CloudFiles(self.frag_path).join(self.frag_path, path) + else: + path = self.frag_path all_labels = vol.download( bbox.to_slices(), From e02882071ca1bac49a58a46bbe4a22960ab56437 Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Tue, 30 Jul 2024 15:34:31 -0400 Subject: [PATCH 38/48] fix: handle binarization issue in repairing --- igneous/tasks/skeleton.py | 7 ++++++- 1 file changed, 6 insertions(+), 1 deletion(-) diff --git a/igneous/tasks/skeleton.py b/igneous/tasks/skeleton.py index 4faa5b5..1634c48 100644 --- a/igneous/tasks/skeleton.py +++ b/igneous/tasks/skeleton.py @@ -342,6 +342,11 @@ def reprocess_skel(pts, skel): diff = bbox.minpt - skel_bbx.minpt skel.vertices += diff * vol.resolution + # we binarized the label for memory's sake, + # so need to harmonize that with the skeleton ID + segid = skel.id + skel.id = 1 + kimimaro.cross_sectional_area( binary_image, skel, anisotropy=vol.resolution, @@ -351,7 +356,7 @@ def reprocess_skel(pts, skel): fill_holes=self.fill_holes, repair_contacts=True, ) - + skel.id = segid skel.vertices -= diff * vol.resolution for skel in repair_skels: From 52b5bb8d44b0c0b3aaccbb63be857f4c01170712 Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Tue, 30 Jul 2024 21:42:39 -0400 Subject: [PATCH 39/48] fix: smart path --- igneous/tasks/skeleton.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/igneous/tasks/skeleton.py b/igneous/tasks/skeleton.py index 1634c48..3d6d51c 100644 --- a/igneous/tasks/skeleton.py +++ b/igneous/tasks/skeleton.py @@ -121,7 +121,7 @@ def execute(self): # otherwise place the files exactly where frag path says to test_path = CloudFiles(self.frag_path).join(self.frag_path, "info") test_info = CloudFile(test_path).get_json() - if test_info is None or 'scales' in test_info: + if test_info is not None and 'scales' in test_info: path = CloudFiles(self.frag_path).join(self.frag_path, path) else: path = self.frag_path From 12c3117a75326b3f5ab1adf0ff4bbd0a7664bc0b Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Wed, 31 Jul 2024 01:35:19 -0400 Subject: [PATCH 40/48] perf: skip downloading the "huge" bbox if no skeletons --- igneous/tasks/skeleton.py | 3 +++ 1 file changed, 3 insertions(+) diff --git a/igneous/tasks/skeleton.py b/igneous/tasks/skeleton.py index 3d6d51c..59f1564 100644 --- a/igneous/tasks/skeleton.py +++ b/igneous/tasks/skeleton.py @@ -267,6 +267,9 @@ def voxel_connectivity_graph( return vcg def compute_cross_sectional_area(self, vol, bbox, skeletons): + if len(skeletons) == 0: + return skeletons + # Why redownload a bigger image? In order to avoid clipping the # cross sectional areas on the edges. delta = int(self.cross_sectional_area_shape_delta) From be2b4622a42ad2bb88c02957f1c42a1a026667a5 Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Thu, 1 Aug 2024 20:41:16 -0400 Subject: [PATCH 41/48] fix: ensure info file uploaded to frag path --- igneous/task_creation/skeleton.py | 11 ++++++++++- 1 file changed, 10 insertions(+), 1 deletion(-) diff --git a/igneous/task_creation/skeleton.py b/igneous/task_creation/skeleton.py index 5b042b6..5d9a016 100644 --- a/igneous/task_creation/skeleton.py +++ b/igneous/task_creation/skeleton.py @@ -15,7 +15,7 @@ from cloudvolume import CloudVolume from cloudvolume.lib import Vec, Bbox, max2, min2, xyzrange, find_closest_divisor, yellow, jsonify from cloudvolume.datasource.precomputed.sharding import ShardingSpecification -from cloudfiles import CloudFiles +from cloudfiles import CloudFiles, CloudFile from igneous.tasks import ( SkeletonTask, UnshardedSkeletonMergeTask, @@ -225,6 +225,15 @@ def create_skeletonizing_tasks( vol.skeleton.meta.commit_info() + if frag_path: + frag_info_path = CloudFiles(frag_path).join(frag_path, "info") + frag_info = CloudFile(frag_info_path).get_json() + if not frag_info: + CloudFile(frag_info_path).put_json(vol.skeleton.meta.info) + elif frag_info["scales"]: + frag_info_path = CloudFiles(frag_path).join(frag_path, vol.info["skeletons"], "info") + CloudFile(frag_info_path).put_json(vol.skeleton.meta.info) + will_postprocess = bool(np.any(vol.bounds.size3() > shape)) bounds = vol.bounds.clone() From 586bdf1198ae6f939f569c3c6d30772706a31ef2 Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Fri, 2 Aug 2024 01:53:20 -0400 Subject: [PATCH 42/48] fix: check for key --- igneous/task_creation/skeleton.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/igneous/task_creation/skeleton.py b/igneous/task_creation/skeleton.py index 5d9a016..1265466 100644 --- a/igneous/task_creation/skeleton.py +++ b/igneous/task_creation/skeleton.py @@ -230,7 +230,7 @@ def create_skeletonizing_tasks( frag_info = CloudFile(frag_info_path).get_json() if not frag_info: CloudFile(frag_info_path).put_json(vol.skeleton.meta.info) - elif frag_info["scales"]: + elif 'scales' in frag_info: frag_info_path = CloudFiles(frag_path).join(frag_path, vol.info["skeletons"], "info") CloudFile(frag_info_path).put_json(vol.skeleton.meta.info) From be59802ea8a72ffe775ce9dc1d520624afcdf29e Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Tue, 6 Aug 2024 17:34:50 -0400 Subject: [PATCH 43/48] feat: support agglomeration for graphene for cross sections --- igneous/tasks/skeleton.py | 10 +++++++--- 1 file changed, 7 insertions(+), 3 deletions(-) diff --git a/igneous/tasks/skeleton.py b/igneous/tasks/skeleton.py index 59f1564..9adbc19 100644 --- a/igneous/tasks/skeleton.py +++ b/igneous/tasks/skeleton.py @@ -282,9 +282,13 @@ def compute_cross_sectional_area(self, vol, bbox, skeletons): huge_bbox.grow(int(np.max(bbox.size()) / 2) + 1) huge_bbox = Bbox.clamp(huge_bbox, vol.bounds) - mem_vol = vol.image.memory_cutout( - huge_bbox, mip=vol.mip, - encoding="crackle", compress=False + mem_vol = vol.memory_cutout( + huge_bbox, + mip=vol.mip, + encoding="crackle", + compress=False, + agglomerate=True, + timestamp=self.timestamp, ) all_labels = mem_vol[big_bbox][...,0] From 32ed272e0ccc90dae451e1a498f68ae20ed92afd Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Tue, 6 Aug 2024 18:07:28 -0400 Subject: [PATCH 44/48] feat: support static root_ids for a graphene volume --- igneous/task_creation/skeleton.py | 6 ++++++ igneous/tasks/skeleton.py | 23 ++++++++++++++++++++++- igneous_cli/cli.py | 5 ++++- 3 files changed, 32 insertions(+), 2 deletions(-) diff --git a/igneous/task_creation/skeleton.py b/igneous/task_creation/skeleton.py index 1265466..e99ab06 100644 --- a/igneous/task_creation/skeleton.py +++ b/igneous/task_creation/skeleton.py @@ -76,6 +76,7 @@ def create_skeletonizing_tasks( cross_sectional_area=False, cross_sectional_area_smoothing_window=5, timestamp=None, + root_ids_cloudpath=None, ): """ Assign tasks with one voxel overlap in a regular grid @@ -169,6 +170,9 @@ def create_skeletonizing_tasks( cross_sectional_area_smoothing_window: Perform a rolling average of the normal vectors across these many vectors. timestamp: for graphene volumes only, you can specify the timepoint to use + root_ids_cloudpath: for graphene volumes, if you have a materialized archive + if your desired timepoint, you can use this path for fetching root ID + segmentation as it is far more efficient. """ shape = Vec(*shape) vol = CloudVolume(cloudpath, mip=mip, info=info) @@ -284,6 +288,7 @@ def task(self, shape, offset): timestamp=timestamp, cross_sectional_area=bool(cross_sectional_area), cross_sectional_area_smoothing_window=int(cross_sectional_area_smoothing_window), + root_ids_cloudpath=root_ids_cloudpath, ) def synapses_for_bbox(self, shape, offset): @@ -331,6 +336,7 @@ def on_finish(self): 'timestamp': timestamp, 'cross_sectional_area': bool(cross_sectional_area), 'cross_sectional_area_smoothing_window': int(cross_sectional_area_smoothing_window), + 'root_ids_cloudpath': root_ids_cloudpath, }, 'by': operator_contact(), 'date': strftime('%Y-%m-%d %H:%M %Z'), diff --git a/igneous/tasks/skeleton.py b/igneous/tasks/skeleton.py index 9adbc19..dd67ea4 100644 --- a/igneous/tasks/skeleton.py +++ b/igneous/tasks/skeleton.py @@ -83,6 +83,7 @@ def __init__( strip_integer_attributes:bool = True, fix_autapses:bool = False, timestamp:Optional[int] = None, + root_ids_cloudpath:Optional[str] = None, ): super().__init__( cloudpath, shape, offset, mip, @@ -97,6 +98,7 @@ def __init__( int(cross_sectional_area_shape_delta), bool(dry_run), bool(strip_integer_attributes), bool(fix_autapses), timestamp, + root_ids_cloudpath, ) if isinstance(self.frag_path, str): self.frag_path = cloudfiles.paths.normalize(self.frag_path) @@ -104,8 +106,15 @@ def __init__( self.index_bounds = Bbox(offset, Vec(*spatial_grid_shape) + Vec(*offset)) def execute(self): + # For graphene volumes, if we've materialized the root IDs + # into a static archive, let's use that because it's way more + # efficient for fetching root IDs. + cloudpath = self.cloudpath + if self.root_ids_cloudpath: + cloudpath = self.root_ids_cloudpath + vol = CloudVolume( - self.cloudpath, mip=self.mip, + cloudpath, mip=self.mip, info=self.info, cdn_cache=False, parallel=self.parallel, fill_missing=self.fill_missing, @@ -209,6 +218,18 @@ def voxel_connectivity_graph( bbox:Bbox, root_labels:np.ndarray, ) -> np.ndarray: + + if vol.meta.path.format != "graphene": + vol = CloudVolume( + self.cloudpath, mip=self.mip, + info=self.info, cdn_cache=False, + parallel=self.parallel, + fill_missing=self.fill_missing, + ) + + if vol.meta.path.format != "graphene": + raise ValueError("Can't extract a voxel connectivity graph from non-graphene volumes.") + layer_2 = vol.download( bbox, stop_layer=2, diff --git a/igneous_cli/cli.py b/igneous_cli/cli.py index 29224ba..77d1540 100644 --- a/igneous_cli/cli.py +++ b/igneous_cli/cli.py @@ -1197,6 +1197,8 @@ def skeletongroup(): @click.option('--labels', type=TupleN(), default=None, help="Skeletonize only this comma separated list of labels.", show_default=True) @click.option('--cross-section', type=int, default=0, help="Compute the cross sectional area for each skeleton vertex. May add substantial computation time. Integer value is the normal vector rolling average smoothing window over vertices. 0 means off.", show_default=True) @click.option('--output', '-o', type=CloudPath(), default=None, help="Output the results to a different place.", show_default=True) +@click.option('--timestamp', type=int, default=None, help="(graphene) Use the proofreading state at this UNIX timestamp.", show_default=True) +@click.option('--root-ids', type=CloudPath(), default=None, help="(graphene) If you have a materialization of graphene root ids for this timepoint, it's more efficient to use it than making requests to the graphene server.", show_default=True) @click.pass_context def skeleton_forge( ctx, path, queue, mip, shape, @@ -1204,7 +1206,7 @@ def skeleton_forge( fix_branching, fix_borders, fix_avocados, fix_autapses, fill_holes, scale, const, soma_detect, soma_accept, soma_scale, soma_const, max_paths, sharded, labels, - cross_section, output, + cross_section, output, timestamp, root_ids, ): """ (1) Synthesize skeletons from segmentation cutouts. @@ -1250,6 +1252,7 @@ def skeleton_forge( cross_sectional_area=(cross_section > 0), cross_sectional_area_smoothing_window=int(cross_section), frag_path=output, fix_autapses=fix_autapses, + timestamp=timestamp, root_ids_cloudpath=root_ids, ) enqueue_tasks(ctx, queue, tasks) From 7e915c0905342b396be5bdc33da90f0b0f8fb494 Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Thu, 8 Aug 2024 19:01:40 -0400 Subject: [PATCH 45/48] fix: regression in LuminanceLevelsTask due to CloudVolume upgrade --- igneous/tasks/image/image.py | 6 +++++- 1 file changed, 5 insertions(+), 1 deletion(-) diff --git a/igneous/tasks/image/image.py b/igneous/tasks/image/image.py index 2e72cd6..31443e5 100755 --- a/igneous/tasks/image/image.py +++ b/igneous/tasks/image/image.py @@ -327,6 +327,9 @@ def execute(self): cts = np.bincount(img2d) levels[0:len(cts)] += cts.astype(np.uint64) + if len(bboxes) == 0: + return + covered_area = sum([bbx.volume() for bbx in bboxes]) bboxes = [(bbox.volume(), bbox.size3()) for bbox in bboxes] @@ -376,7 +379,8 @@ def select_bounding_boxes(self, dataset_bounds): patch_start += self.offset bbox = Bbox(patch_start, patch_start + sample_shape.size3()) bbox = Bbox.clamp(bbox, dataset_bounds) - bboxes.append(bbox) + if not bbox.subvoxel(): + bboxes.append(bbox) return bboxes @queueable From 42bcde15d9c380447969aa79723048893b1b9498 Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Fri, 9 Aug 2024 21:13:22 -0400 Subject: [PATCH 46/48] refactor: simplify vcg code --- igneous/tasks/skeleton.py | 20 ++++---------------- 1 file changed, 4 insertions(+), 16 deletions(-) diff --git a/igneous/tasks/skeleton.py b/igneous/tasks/skeleton.py index dd67ea4..30fae46 100644 --- a/igneous/tasks/skeleton.py +++ b/igneous/tasks/skeleton.py @@ -234,8 +234,8 @@ def voxel_connectivity_graph( bbox, stop_layer=2, agglomerate=True, - timestamp=self.timestamp - ) + timestamp=self.timestamp, + )[...,0] graph_chunk_size = np.array(vol.meta.graph_chunk_size) / vol.meta.downsample_ratio(vol.mip) graph_chunk_size = graph_chunk_size.astype(int) @@ -243,20 +243,7 @@ def voxel_connectivity_graph( shape = bbox.size()[:3] sgx, sgy, sgz = list(np.ceil(shape / graph_chunk_size).astype(int)) - vcg = np.zeros(layer_2.shape[:3], dtype=np.uint32, order="F") - - clamp_box = Bbox([0,0,0], shape) - - for gx,gy,gz in xyzrange([sgx, sgy, sgz]): - bbx = Bbox((gx,gy,gz), (gx+1, gy+1, gz+1)) - bbx *= graph_chunk_size - bbx = Bbox.clamp(bbx, clamp_box) - - cutout = np.asfortranarray(layer_2[bbx.to_slices()][:,:,:,0]) - vcg_cutout = cc3d.voxel_connectivity_graph(cutout, connectivity=26) - vcg[bbx.to_slices()] = vcg_cutout - del vcg_cutout - + vcg = cc3d.voxel_connectivity_graph(layer_2, connectivity=26) del layer_2 # the proper way to do this would be to get the lowest the L3..LN root @@ -266,6 +253,7 @@ def voxel_connectivity_graph( # with edges that represent the connections between the adjacent boxes. root_vcg = cc3d.voxel_connectivity_graph(root_labels, connectivity=26) + clamp_box = Bbox([0,0,0], shape) for gx,gy,gz in xyzrange([sgx, sgy, sgz]): bbx = Bbox((gx,gy,gz), (gx+1, gy+1, gz+1)) From c85eee9ac46381c4b27e491ecee1fb8497df1a55 Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Mon, 12 Aug 2024 22:27:56 -0400 Subject: [PATCH 47/48] feat: add encoding level support for jpegxl --- igneous/task_creation/common.py | 2 ++ 1 file changed, 2 insertions(+) diff --git a/igneous/task_creation/common.py b/igneous/task_creation/common.py index 5600fe5..17ed399 100644 --- a/igneous/task_creation/common.py +++ b/igneous/task_creation/common.py @@ -226,6 +226,8 @@ def set_encoding(cv, mip, encoding, encoding_level): if encoding == "jpeg": scale["jpeg_quality"] = encoding_level + elif encoding == "jpegxl": + scale["jpegxl_quality"] = encoding_level elif encoding == "png": scale["png_level"] = encoding_level elif encoding == "fpzip": From 27ce91bb57ca498729521cf00b8283d38535f9da Mon Sep 17 00:00:00 2001 From: William Silversmith Date: Wed, 14 Aug 2024 01:29:16 -0400 Subject: [PATCH 48/48] fix: don't print error if max_mips is non-zero --- igneous/tasks/image/image.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/igneous/tasks/image/image.py b/igneous/tasks/image/image.py index 31443e5..b093148 100755 --- a/igneous/tasks/image/image.py +++ b/igneous/tasks/image/image.py @@ -70,7 +70,7 @@ def downsample_and_upload( if max_mips is not None: factors = factors[:max_mips] - if len(factors) == 0: + if len(factors) == 0 and max_mips: print("No factors generated. Image Shape: {}, Downsample Shape: {}, Volume Shape: {}, Bounds: {}".format( image.shape, ds_shape, vol.volume_size, bounds) )