object filenames WIP

02_CODE/cfg/hfe-nodaratis.yaml

@@ -9,7 +9,7 @@ hydra:
   mode: MULTIRUN
   sweeper:
     params:
-      simulations.grayscale_filenames: C0000000
+      simulations.grayscale_filenames: C0003093
 
 
 mesher:

02_CODE/cfg/paths-nodaratis.yaml

@@ -10,11 +10,11 @@ paths:
   odb_python_script: 02_CODE/src/hfe_abq/readODB_acc.py
 
 filenames:
-  filename_postfix_cort_mask: _V3_CORT_MASK.AIM
-  filename_postfix_trab_mask: _V3_TRAB_MASK.AIM
+  filename_postfix_cort_mask: _V3_CORT_MASK.mha
+  filename_postfix_trab_mask: _V3_TRAB_MASK.mha
   filename_postfix_mask: _MASK.AIM                                  # if mask separate == False, give this filename
-  filename_postfix_bmd: _V3_FINAL.AIM
-  filename_postfix_seg: _V3_SEG.AIM
+  filename_postfix_bmd: _UNCOMP.mha
+  filename_postfix_seg: _SEG.mha
   filename_postfix_common: _common_region_MASK.mhd                  # if registration is true -> file of common region
   filename_postfix_transform: _transformation.tfm                   # if registration is true -> transformation file
 

02_CODE/cfg/simulations-nodaratis.yaml

@@ -1,6 +1,6 @@
 simulations:
   grayscale_filenames:
-    C0002217
+    C0003093
 
   folder_id:
-    C0002217: XCT2_181
+    C0003093: 442_R_75_F

02_CODE/src/hfe_abq/aim2fe.py

@@ -132,7 +132,7 @@ def aim2fe_psl(cfg, sample):
         lock = threading.Lock()
         for item in image_list:
             t = threading.Thread(
-                target=imutils.read_aim, args=(item, filenames, bone, lock)
+                target=imutils.read_image, args=(item, filenames, bone, lock)
             )
             threads.append(t)
             sleep(0.1)  # to avoid overloading the print statements
@@ -142,7 +142,7 @@ def aim2fe_psl(cfg, sample):
     else:
         image_list = ["BMD", "SEG"]
         for _, item in enumerate(image_list):
-            bone = imutils.read_aim(item, filenames, bone, lock=None)
+            bone = imutils.read_image(item, filenames, bone, lock=None)
         bone = imutils.read_aim_mask_combined("MASK", filenames, bone)
 
     # image_list = ["BMD", "SEG", "CORTMASK", "TRABMASK"]

02_CODE/src/hfe_accurate/material_mapping.py

@@ -18,7 +18,7 @@
 from hfe_abq.write_abaqus import AbaqusWriter
 from scipy.spatial import KDTree  # type: ignore
 
-from ..hfe_utils.io_utils import FileConfig
+from hfe_utils.io_utils import FileConfig
 
 LOGGING_NAME = "HFE-ACCURATE"
 logger = logging.getLogger(LOGGING_NAME)

02_CODE/src/hfe_utils/imutils.py

@@ -11,7 +11,7 @@
 from mpl_toolkits.axes_grid1 import make_axes_locatable  # type: ignore
 from vtk.util.numpy_support import numpy_to_vtk, vtk_to_numpy  # type: ignore
 
-from io_utils import FileConfig
+from hfe_utils.io_utils import FileConfig
 
 # flake8: noqa: E501
 
@@ -279,16 +279,47 @@ def AIMreader(fileINname, spacing):
     imvtk = reader.GetOutput()
     return imvtk, spacing, scaling, slope, intercept, header
 
-def mhd_reader(filename):
-    imvtk = sitk.ReadImage(filename)
-    spacing = imvtk.GetSpacing()
+def mha_reader(filename):
+    imsitk = sitk.ReadImage(filename)
+    print(type(imsitk))
+    spacing = imsitk.GetSpacing()
     scaling = None
     slope = None
     intercept = None
     header = None
-    return imvtk, spacing, scaling, slope, intercept, header
-
-
+    return imsitk, spacing, scaling, slope, intercept, header
+
+def mha_header_reader(filename):
+    reader = sitk.ImageFileReader()
+    reader.SetFileName(filename)
+    reader.LoadPrivateTagsOn()
+    reader.ReadImageInformation()
+
+    # Extract metadata and replace _LINEBREAK_ with actual line breaks for SimpleITK compatibility
+    metadata = {}
+    for k in reader.GetMetaDataKeys():
+        v = reader.GetMetaData(k).replace('_LINEBREAK_', '\n')
+        metadata[k] = v
+        print("({0}) = = \"{1}\"".format(k, v))
+
+    # Extract Density: slope from the metadata
+    density_slope = None
+    density_intercept = None
+    mu_scaling = None
+    for _, value in metadata.items():
+        if 'Density: slope' in value or 'Density: intercept' in value or 'Mu_Scaling' in value:
+            lines = value.split('\n')
+            for line in lines:
+                if 'Density: slope' in line:
+                    density_slope = float(line.split()[-1])
+                elif 'Density: intercept' in line:
+                    density_intercept = float(line.split()[-1])
+                elif 'Mu_Scaling' in line:
+                    mu_scaling = int(line.split()[-1])
+
+    spacing = reader.GetSpacing()
+
+    return spacing, mu_scaling, density_slope, density_intercept
 
 def read_img_param(filenames: FileConfig):
     """
@@ -306,10 +337,10 @@ def read_img_param(filenames: FileConfig):
     Raises:
         Exception: If an error occurs while reading the AIM image.
     """
-    print("\n ... read AIM files")
+    print("\n ... read images")
 
     raw_filename = filenames.raw_name
-
+    print(f'raw_filename: {raw_filename}')
     if raw_filename.suffix == 'AIM*':
         try:
             _, spacing, scaling, slope, intercept, _ = AIMreader(
@@ -320,15 +351,19 @@ def read_img_param(filenames: FileConfig):
             raise
     if raw_filename.suffix in ('.mha', '.mhd'):
         try:
-            _, spacing, scaling, slope, intercept, _ = mhd_reader(raw_filename)
+            spacing, scaling, slope, intercept = mha_header_reader(raw_filename)
         except Exception as e:
             logger.exception(f"An error occurred while using mhd_reader: {e}")
             raise
 
     return spacing, scaling, slope, intercept
 
 
-def read_aim(name, filenames, bone, lock):
+def general_image_reader(filename):
+    return sitk.ReadImage(filename)
+
+
+def read_image(name, filenames, bone, lock):
     """
     Reads an AIM image and stores it in the bone dictionary.
 
@@ -351,11 +386,14 @@ def read_aim(name, filenames, bone, lock):
     print("\n ... read file: " + name)
 
     spacing = bone["Spacing"]
-    IMG_vtk = AIMreader(filenames[name + "name"], spacing)[0]
-    IMG_array = vtk2numpy(IMG_vtk)
+    if name.suffix == '*AIM*':
+        IMG_vtk = AIMreader(filenames[name + "name"], spacing)[0]
+        IMG_array = vtk2numpy(IMG_vtk)
+        IMG_sitk = sitk.GetImageFromArray(IMG_array)
+    else:
+        IMG_sitk = general_image_reader(filenames[name])
 
     # pad image to avoid having non-zero values at the boundary
-    IMG_sitk = sitk.GetImageFromArray(IMG_array)
     IMG_pad = pad_image(IMG_sitk, iso_pad_size=10)
     IMG_pad = sitk.Flip(IMG_pad, [True, False, False])