Merge pull request #84 from Hendrik-code/new_features_robert

New features robert

Author: robert-graf

.gitignore

@@ -163,4 +163,6 @@ poetry.lock
 tutorials/tutorial_data_processing/*
 tutorials/*PixelPandemonium/*
 tutorials/dataset-PixelPandemonium/*
-*.html
+*.html
+_*.py
+dicom_select

TPTBox/core/bids_files.py

@@ -1394,7 +1394,7 @@ def __getitem__(self, item: str) -> list[BIDS_FILE]:
         try:
             return self.data_dict[item]
         except KeyError as e:
-            raise KeyError(f"BIDS_Family does not contain key {item}, only {self.keys()}") from e
+            raise KeyError(f"BIDS_Family does not contain key '{item}', only {list(self.keys())}; {self.family_id=}") from e
 
     def __setitem__(self, key, value):
         self.data_dict[key] = value

TPTBox/core/dicom/dicom2nii_utils.py

@@ -201,9 +201,11 @@ def save_json(json_ob, file, check_exist=False):
     """
 
     def convert(obj):
-        if isinstance(obj, np.int64):
+        if isinstance(obj, np.integer):
             return int(obj)
-        raise TypeError
+        if isinstance(obj, np.floating):
+            return float(obj)
+        raise TypeError(type(obj))
 
     if check_exist and test_name_conflict(json_ob, file):
         raise FileExistsError(file)

TPTBox/core/nii_poi_abstract.py

@@ -56,7 +56,7 @@ def to_gird(self) -> Grid:
     @property
     def shape_int(self):
         assert self.shape is not None, "need shape information"
-        return tuple(np.rint(list(self.shape)).astype(int))
+        return tuple(np.rint(list(self.shape)).astype(int).tolist())
 
     @property
     def spacing(self):
@@ -68,11 +68,11 @@ def spacing(self, value: ZOOMS):
 
     def __str__(self) -> str:
         try:
-            origin = tuple(np.around(self.origin, decimals=2))
+            origin = tuple(np.around(self.origin, decimals=2).tolist())
         except Exception:
             origin = self.origin
         try:
-            zoom = tuple(np.around(self.zoom, decimals=2))
+            zoom = tuple(np.around(self.zoom, decimals=2).tolist())
         except Exception:
             zoom = self.zoom
 

TPTBox/core/nii_wrapper.py

@@ -437,6 +437,16 @@ def rotation(self)->np.ndarray:
         rotation = rotation_zoom / zoom
         return rotation
 
+    @property
+    def direction(self) -> list:
+        direction_matrix = self.rotation
+        direction_flat = direction_matrix.flatten().tolist()
+        return direction_flat
+    @property
+    def direction_itk(self) -> list:
+        a = np.array(self.direction)
+        a[:len(a)//3*2]*=-1
+        return a.tolist()
 
     @orientation.setter
     def orientation(self, value: AX_CODES):
@@ -488,7 +498,7 @@ def get_array(self) -> np.ndarray:
         return self._arr.copy()
     def numpy(self, *_args):
         return self.get_array()
-    def set_array(self,arr:np.ndarray|Self, inplace=False,verbose:logging=False)-> Self:  # noqa: ARG002
+    def set_array(self,arr:np.ndarray|Self, inplace=False,verbose:logging=False,seg=None)-> Self:  # noqa: ARG002
         """Creates a NII where the array is replaces with the input array.
 
         Note: This function works "Out-of-place" by default, like all other methods.
@@ -500,6 +510,7 @@ def set_array(self,arr:np.ndarray|Self, inplace=False,verbose:logging=False)-> S
         Returns:
             self
         """
+
         if hasattr(arr,"get_array"):
             arr = arr.get_array() # type: ignore
         if arr.dtype == bool:
@@ -520,10 +531,12 @@ def set_array(self,arr:np.ndarray|Self, inplace=False,verbose:logging=False)-> S
         #if all(a is None for a in self.header.get_slope_inter()):
         #    nii.header.set_slope_inter(1,self.get_c_val()) # type: ignore
         if inplace:
+            if seg is not None:
+                self.seg = seg
             self.nii = nii
             return self
         else:
-            return self.copy(nii) # type: ignore
+            return self.copy(nii,seg=seg) # type: ignore
 
     def set_array_(self,arr:np.ndarray,verbose:logging=True):
         return self.set_array(arr,inplace=True,verbose=verbose)
@@ -1255,13 +1268,13 @@ def filter_connected_components(self, labels: int |list[int]|None=None,min_volum
         """
         assert self.seg, "This only works on segmentations"
         arr = np_filter_connected_components(self.get_seg_array(), largest_k_components=max_count_component,label_ref=labels,connectivity=connectivity,return_original_labels=keep_label,min_volume=min_volume,max_volume=max_volume,removed_to_label=removed_to_label,)
-        #if keep_label and labels is not None:
-        #    if isinstance(labels,int):
-        #        labels = [labels]
-        #    old_labels = [i for i in self.unique() if i not in labels]
-        #    if len(old_labels) != 0:
-        #        s = self.extract_label(old_labels,keep_label=True)
-        #        nii[s != 0] = s[s!=0]
+        if keep_label and labels is not None:
+            if isinstance(labels,int):
+                labels = [labels]
+            old_labels = [i for i in self.unique() if i not in labels]
+            if len(old_labels) != 0:
+                s = self.extract_label(old_labels,keep_label=True).get_array()
+                arr[s != 0] = s[s!=0]
         #print("filter",nii.unique())
         #assert max_count_component is None or nii.max() <= max_count_component, nii.unique()
         return self.set_array(arr, inplace=inplace)
@@ -1648,10 +1661,10 @@ def map_labels(self, label_map:LABEL_MAP , verbose:logging=True, inplace=False):
     def map_labels_(self, label_map: LABEL_MAP, verbose:logging=True):
         return self.map_labels(label_map,verbose=verbose,inplace=True)
 
-    def copy(self, nib:Nifti1Image|_unpacked_nii|None = None):
+    def copy(self, nib:Nifti1Image|_unpacked_nii|None = None,seg=None):
         if nib is None:
             nib = (self.get_array().copy(), self.affine.copy(), self.header.copy())
-        return NII(nib,seg=self.seg,c_val = self.c_val,info = self.info)
+        return NII(nib,seg=self.seg if seg is None else seg,c_val = self.c_val,info = self.info)
 
     def clone(self):
         return self.copy()
@@ -1745,9 +1758,10 @@ def __array__(self,dtype=None):
                 return self._arr
             else:
                 return self._arr.astype(dtype, copy=False)
-    def __array_wrap__(self, array):
+    def __array_wrap__(self, array,context=None, return_scalar=False):
+        assert not return_scalar,context
         if array.shape != self.shape:
-            raise SyntaxError(f"Function call induce a shape change of nii image. Before {self.shape} after {array.shape}.")
+            raise SyntaxError(f"Function call induce a shape change of nii image. Before {self.shape} after {array.shape}. {context}")
         return self.set_array(array)
     def __getitem__(self, key)-> Any:
         if isinstance(key,Sequence):

TPTBox/core/np_utils.py

@@ -74,6 +74,22 @@ def np_extract_label(
 
 
 def cc3dstatistics(arr: UINTARRAY, use_crop: bool = True) -> dict:
+    """
+    Computes connected component statistics for a labeled array using connected components 3D (cc3d).
+
+    Args:
+        arr (UINTARRAY): A 3D array of unsigned integers or booleans where each connected component
+                         is labeled with a unique integer. Typically output from a labeling function.
+        use_crop (bool): If True, the function attempts to crop the input array around non-zero regions
+                         to improve performance and focus statistics on the area of interest. Defaults to True.
+
+    Returns:
+        dict: A dictionary containing statistics of the connected components, such as their sizes,
+              bounding boxes, and possibly centroids, depending on implementation of `_cc3dstats`.
+
+    Raises:
+        AssertionError: If the input array is not of an unsigned integer or boolean dtype.
+    """
     assert np.issubdtype(arr.dtype, np.unsignedinteger) or np.issubdtype(arr.dtype, np.bool_), (
         f"cc3dstatistics expects uint type, got {arr.dtype}"
     )
@@ -164,7 +180,7 @@ def np_unique_withoutzero(arr: UINTARRAY) -> list[int]:
     return [i for i in np_unique(arr) if i != 0]
 
 
-def np_center_of_mass(arr: UINTARRAY) -> dict[int, COORDINATE]:
+def np_center_of_mass(arr: UINTARRAY) -> dict[int, np.ndarray]:
     """Calculates center of mass, mapping label in array to a coordinate (float) (exluding zero)
 
     Args:
@@ -662,7 +678,7 @@ def np_connected_components(
         verbose: If true, will print out if the array does not have any CC
 
     Returns:
-        arr_cc: UINTARRAY, N: int
+        arr_cc: UINTARRAY, N: number of cc
     """
     assert np.min(arr) == 0, f"min value of mask not zero, got {np.min(arr)}"
     assert np.max(arr) >= 0, f"wrong normalization, max value is not >= 0, got {np_unique(arr)}"
@@ -757,7 +773,7 @@ def np_filter_connected_components(
 
     arr2 = arr.copy()
     labels: Sequence[int] = _to_labels(arr, label_ref)
-    arr2[np.isin(arr, labels, invert=True)] = 0  # type:ignore
+    arr2[np.isin(arr2, labels, invert=True)] = 0  # type:ignore
 
     labels_out, n = _connected_components(arr2, connectivity=connectivity, return_N=True)
     if largest_k_components is None:
@@ -831,7 +847,7 @@ def np_translate_to_center_of_array(image: np.ndarray) -> np.ndarray:
 
 
 def np_translate_arr(arr: np.ndarray, translation_vector: tuple[int, int] | tuple[int, int, int]) -> np.ndarray:
-    """Translates nonzero values of an input array according to a 2D or 3D translation vector. Values that would be shifted beyond the boundary are removed!
+    """Translates values of an input array according to a 2D or 3D translation vector. Values that would be shifted beyond the boundary are removed!
 
     Args:
         arr: input array

TPTBox/core/poi_fun/poi_global.py

@@ -1,5 +1,6 @@
 from __future__ import annotations
 
+import json
 from copy import deepcopy
 from pathlib import Path
 
@@ -10,6 +11,7 @@
 from TPTBox.core.poi_fun.poi_abstract import Abstract_POI, POI_Descriptor
 from TPTBox.core.poi_fun.save_load import FORMAT_GLOBAL, load_poi, save_poi
 from TPTBox.core.vert_constants import Abstract_lvl, logging
+from TPTBox.logger.log_file import log
 
 ###### GLOBAL POI #####
 
@@ -138,7 +140,7 @@ def to_other(self, msk: Has_Grid, verbose=False) -> poi.POI:
                 v = (-v[0], -v[1], v[2])  # noqa: PLW2901
             v_out = msk.global_to_local(v)
             if verbose:
-                print(v, "-->", v_out)
+                log.print(v, "-->", v_out)
             out[k1, k2] = tuple(v_out)
 
         return poi.POI(centroids=out, **msk._extract_affine(), info=self.info, format=self.format)
@@ -155,8 +157,8 @@ def copy(self, centroids: POI_Descriptor | None = None) -> Self:
     @classmethod
     def load(cls, poi: poi.POI_Reference, itk_coords: bool | None = None) -> Self:
         poi_obj = load_poi(poi)
-        if not isinstance(poi_obj, POI_Global):
-            poi_obj = poi_obj.to_global(itk_coords if itk_coords is not None else False)
+        if not poi_obj.is_global:
+            poi_obj = poi_obj.to_global(itk_coords if itk_coords is not None else False)  # type: ignore
         if itk_coords is not None:
             assert itk_coords == poi_obj.itk_coords, "not implemented swichting to/from itk_coords to nii "
         return poi_obj  # type: ignore
@@ -173,3 +175,79 @@ def save(
         return save_poi(
             self, out_path, make_parents, additional_info, save_hint=save_hint, resample_reference=resample_reference, verbose=verbose
         )
+
+    def save_mrk(self, filepath: str | Path, color=None, split_by_region=True, split_by_subregion=False):
+        """
+        Save the POI data to a .mrk.json file in Slicer Markups format.
+        Automatically sets coordinate system based on itk_coords.
+        Includes level_one_info and level_two_info in the description.
+        Preserves metadata from `info` dictionary.
+        """
+        if color is None:
+            color = self.info.get("color", [1.0, 0.0, 0.0])
+        filepath = Path(filepath)
+        if not filepath.name.endswith(".mrk.json"):
+            filepath = filepath.parent / (filepath.stem + ".mrk.json")
+        coordinate_system = "LPS" if self.itk_coords else "RAS"
+
+        # Create list of control points
+        from TPTBox import NII
+        from TPTBox.mesh3D.mesh_colors import get_color_by_label
+
+        list_markups = {}
+        for region, subregion, coords in self.centroids.items():
+            try:
+                name = self.level_two_info(subregion).name
+            except Exception:
+                name = subregion
+            try:
+                name2 = self.level_one_info(region).name
+            except Exception:
+                name2 = region
+            key = "P"
+            color2 = color
+            if split_by_region:
+                key += str(region) + "_"
+                color2 = get_color_by_label(region).rgb.tolist()
+            if split_by_subregion:
+                key += str(subregion)
+                color2 = get_color_by_label(region).rgb.tolist()
+            if key not in list_markups:
+                list_markups[key] = {
+                    "type": "Fiducial",
+                    "coordinateSystem": coordinate_system,
+                    "locked": False,
+                    "labelFormat": "%N-%d",
+                    "controlPoints": [],
+                    "display": {
+                        "visibility": True,
+                        "opacity": 1.0,
+                        "color": color2.copy(),
+                        "propertiesLabelVisibility": False,
+                    },
+                    "description": "",  # self.info,
+                }
+
+            list_markups[key]["controlPoints"].append(
+                {
+                    "id": f"{region}-{subregion}",
+                    "label": f"{region}-{subregion}",
+                    "description": name,
+                    "associatedNodeID": name2,
+                    "position": list(coords),
+                    "orientation": [1.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 1.0],
+                    "selected": False,
+                    "locked": False,
+                    "visibility": True,
+                    "positionStatus": "defined",
+                }
+            )
+        mrk_data = {
+            "markups": list(list_markups.values()),
+            "schema": "https://raw.githubusercontent.com/slicer/slicer/master/Modules/Loadable/Markups/Resources/Schema/markups-schema-v1.0.3.json#",
+            # "coordinateSystem": coordinate_system,
+        }
+
+        with open(filepath, "w") as f:
+            json.dump(mrk_data, f, indent=2)
+        log.on_save(f"Saved .mrk.json to {filepath}")