Update thalamus (#11)

* Overhaul the thalamus into steps, and add docs

Load thalamus meshes from CLI input, clean CLI

* Improve Blender instructions

* Change thalamus region list regex

This updates the regular expression used for thalamus placement-hints to
be in a different format that has been tested successfully, excludes
habenular and peripeduncular subregions, and to be valid for the
hierarchy/annotation used at its appropriate step in the Atlas pipeline.
For information on which regions were chosen and this list was created,
see the internal BBP Confluence page located at "Circuits > Mouse
Thalamus > Atlas-based Whole-thalamus subregion selection". This regex
has been built from the region list of the desired and present thalamus
regions as of the "final" version of the hierarchy and annotation built
by the Atlas pipeline, which is the output of the rule
`split_barrel_ccfv3_l23split`.

This change is meant to go in tandem with
BlueBrain/atlas-direction-vectors#27 .

* Update layer names to be Atlas-Pipeline-compatible

* Fix formatting errors

* Attempt to update tests for new thal workflow

* Replace part of test anno with region in metadata

* fix test + format

* format

* Fix final linting issues

This does a lot of small things for passing the linting.

For mypy, I had to add additional ignores for the Trimesh returned
types, since the ignore on the module as a whole wasn't preventing mypy
from expecting `load_mesh` to return a Geometry object, which is a
grandparent of Trimesh objects. I don't know if Trimesh changed their
API, I couldn't figure it out from the docs, and I don't know why mypy
was raising this now. In all the cases I could test or see, a proper
"Trimesh" object was returned instead of the more generic Geometry. I
don't think we need to worry about this.

For the pylint disable W0613 (unused-argument), I needed some
polymorphism for the thalamus case, but I wasn't sure how to handle that
alongside the linters' type-checking. I think this is the simplest
solution and is harmless.

Everything else is minor.

* Make Alexis changes to CLI

* Apply MG code review changes

---------

Co-authored-by: Austin E. Soplata <[email protected]>
Co-authored-by: arnaudon <[email protected]>

Author: 3 people

atlas_placement_hints/app/metadata/thalamus_metadata.json

@@ -7,9 +7,10 @@
 
     },
     "layers": {
-        "names": ["RT", "THnotRT"],
-        "queries": ["RT", "@^\\bAD|AM|AMd|AMv|ATN|AV|CL|CM|DORpm|DORsm|EPI|Eth|GENd|GENv|IAD|IAM|IGL|ILM|IMD|IntG|LAT|LD|LGd|LGd-co|LGd-ip|LGd-sh|LGv|LGvl|LGvm|LH|LP|MD|MDc|MDl|MDm|MED|MG|MGd|MGm|MGv|MH|MTN|PCN|PF|PIL|PIN|PO|POL|PP|PR|PT|PVT|PoT|RE|REth|RH|SGN|SMT|SPA|SPF|SPFm|SPFp|SubG|TH|VAL|VENT|VM|VP|VPL|VPLpc|VPM|VPMpc|Xi\\b$"],
+        "names": ["layer_1", "layer_2"],
+        "queries": ["RT",
+"@^(?:AD|AMd|AMv|AV|CL|CM|Eth|IAD|IAM|IGL|IMD|IntG|LD|LGd-co|LGd-ip|LGd-sh|LGv_O|LP|MD_O|MGd|MGm|MGv|PCN|PF|PIL|PO|POL|PR|PT|PVT|PoT|RE|RH|SGN|SMT|SPA|SPFm|SPFp|SubG|TH_O|VAL|VM|VPL|VPLpc|VPM|VPMpc|Xi)$"],
         "attribute": "acronym",
         "with_descendants": false
     }
-}
+}

atlas_placement_hints/app/placement_hints.py

@@ -74,6 +74,7 @@ def _placement_hints(  # pylint: disable=too-many-locals
     max_thicknesses: Optional[List[float]] = None,
     flip_direction_vectors: bool = False,
     has_hemispheres: bool = False,
+    thalamus_meshes_dir: str = "",
 ) -> None:
     """
     Compute the placement hints for a laminar region of the mouse brain.
@@ -93,6 +94,10 @@ def _placement_hints(  # pylint: disable=too-many-locals
         has_hemispheres: (optional) If True, split the volume into halves along the z-axis and
             handle each of theses 'hemispheres' separately. Otherwise the whole volume is handled.
             Defaults to True.
+        thalamus_meshes_dir: (optional) Path of the directory to load thalamus meshes
+            from. Currently only used for thalamus. Required if you are producing thalamus
+            placement-hints. Defaults to None.
+
     """
     direction_vectors = voxcell.VoxelData.load_nrrd(direction_vectors_path)
     assert_meta_properties([direction_vectors, atlas.region])
@@ -106,6 +111,7 @@ def _placement_hints(  # pylint: disable=too-many-locals
         max_thicknesses,
         flip_direction_vectors=flip_direction_vectors,
         has_hemispheres=has_hemispheres,
+        thalamus_meshes_dir=thalamus_meshes_dir,
     )
     if not Path(output_dir).exists():
         os.makedirs(output_dir)
@@ -305,41 +311,84 @@ def isocortex(
     required=True,
     help="path of the directory to write. It will be created if it doesn't exist.",
 )
+@click.option(
+    "--thalamus-meshes-dir",
+    required=True,
+    help="""Path of the directory to use for either saving or loading thalamus
+    meshes. It will be created if it doesn't exist.""",
+)
+@click.option(
+    "--load-cut-thalamus-meshes",
+    required=False,
+    help="""(Optional) Flag to load your custom thalamus meshes, and then use
+    them to calculate placement-hints.""",
+    default=False,
+    is_flag=True,
+)
 @log_args(L)
+# pylint: disable=too-many-arguments
 def thalamus(
-    verbose, annotation_path, hierarchy_path, metadata_path, direction_vectors_path, output_dir
+    verbose,
+    annotation_path,
+    hierarchy_path,
+    metadata_path,
+    direction_vectors_path,
+    output_dir,
+    thalamus_meshes_dir,
+    load_cut_thalamus_meshes,
 ):
     """Generate and save the placement hints of the mouse thalamus.
 
-    Placement hints are saved under the names sepecified in `app/metadata/thalamus_metadata.json`.
-    Default to:
+    First, call this without passing '--load-cut-thalamus-meshes' to
+    create your region meshes, but not your placement-hints. Then, hand-cut
+    your meshes according to the documentation in
+    'atlas_placement_hints/layered_atlas.py::ThalamusAtlas'. Finally, call
+    this again while passing '--load-cut-thalamus-meshes'.
+
+    Placement hints are saved under the names specified in
+    `app/metadata/thalamus_metadata.json`. These default to:
 
     \b
     - `[PH]y.nrrd`
-    - `[PH]Rt.nrrd`, `[PH]VPL.nrrd`
+    - `[PH]layer_1.nrrd` (This is the "top-most" layer, equivalent to "RT".
+      This previously used the filename `[PH]RT.nrrd`)
+    - `[PH]layer_2.nrrd` (This is the "deepest" layer, equivalent to the
+      thalamus except the habenular, peripeduncular, and reticular regions.
+      This previously used the filename `[PH]THnotRT.nrrd`)
 
-    A report together with an nrrd volume on problematic distance computations are generated
-    in `output_dir` under the names:
+    A report together with an nrrd volume on problematic distance computations
+    are generated in `output_dir` under the names:
 
     \b
     - `distance_report.json`
-    - `<Thalamus>_problematic_voxel_mask.nrrd` (mask of the voxels for which the computed
-    placement hints cannot be trusted).  <Thalamus> is the region name specified in
-    thalamus_metadata.json. Defaults to "Thalamus".
-
-    The annotation file can contain the thalamus or a superset.
-    For the algorithm to work properly, some space should separate the boundary
-    of the thalamus from the boundary of its enclosing array.
+    - `<Thalamus>_problematic_voxel_mask.nrrd` (mask of the voxels for which
+      the computed placement hints cannot be trusted).  <Thalamus> is the
+      region name specified in thalamus_metadata.json. Defaults to "Thalamus".
+
+    The annotation file can contain the thalamus or a superset. For the
+    algorithm to work properly, some space should separate the boundary of the
+    thalamus from the boundary of its enclosing array.
+
+    For instructions on all steps necessary to generate the thalamus' placement
+    hints, see
+    'atlas-placement-hints/atlas_placement_hints/layered_atlas.py::ThalamusAtlas'
+    and its methods for details.
     """
     set_verbose(L, verbose)
 
     atlas = _create_layered_atlas(annotation_path, hierarchy_path, metadata_path)
-    _placement_hints(
-        atlas,
-        direction_vectors_path,
-        output_dir,
-        has_hemispheres=True,
-    )
+
+    if load_cut_thalamus_meshes:
+        _placement_hints(
+            atlas,
+            direction_vectors_path,
+            output_dir,
+            thalamus_meshes_dir=thalamus_meshes_dir,
+            has_hemispheres=True,
+        )
+    else:
+        Path(thalamus_meshes_dir).mkdir(parents=True, exist_ok=True)
+        atlas.create_uncut_thalamus_meshes(thalamus_meshes_dir)
 
 
 @app.command()

atlas_placement_hints/compute_placement_hints.py

@@ -27,6 +27,7 @@ def compute_placement_hints(
     max_thicknesses: Optional[List[float]] = None,
     flip_direction_vectors: bool = False,
     has_hemispheres: bool = True,
+    thalamus_meshes_dir: str = "",
 ) -> Tuple[DistanceInfo, Dict]:
     """
     Compute the placement hints for a laminar region of the mouse brain.
@@ -43,6 +44,9 @@ def compute_placement_hints(
         has_hemispheres: (optional) If True, split the volume into halves along the z-axis and
             handle each of theses 'hemispheres' separately. Otherwise the whole volume is
             handled. Defaults to True.
+        thalamus_meshes_dir: (optional) Path of the directory to load thalamus meshes
+            from. Currently only used for thalamus. Required if you are producing thalamus
+            placement-hints. Defaults to None.
 
     Returns:
         Tuple with the following items.
@@ -70,6 +74,7 @@ def compute_placement_hints(
         direction_vectors,
         flip_direction_vectors=flip_direction_vectors,
         has_hemispheres=has_hemispheres,
+        thalamus_meshes_dir=thalamus_meshes_dir,
     )
 
     distances_to_meshes = distances_info["distances_to_layer_boundaries"]

atlas_placement_hints/distances/create_watertight_mesh.py

@@ -210,7 +210,7 @@ def create_watertight_trimesh(
         optimized_mesh = trimesh.load_mesh(output_filepath_opt)
 
     if optimization_info:
-        log_mesh_optimization_info(optimized_mesh, unoptimized_mesh)
+        log_mesh_optimization_info(optimized_mesh, unoptimized_mesh)  # type: ignore
 
-    optimized_mesh.fix_normals()
-    return optimized_mesh
+    optimized_mesh.fix_normals()  # type: ignore
+    return optimized_mesh  # type: ignore