Directly report properties which throw errors

SOAP/compute_halo_properties.py

@@ -12,6 +12,7 @@
 
 import os
 import os.path
+import sys
 import time
 import traceback
 import numpy as np
@@ -564,6 +565,8 @@ def compute_halo_properties():
         )
     except Exception as e:
         traceback.print_exc()
+        sys.stdout.flush()
+        sys.stderr.flush()
         comm_world.Abort(1)
 
     # Can stop the halo request thread now that all chunk tasks have executed

SOAP/core/halo_tasks.py

@@ -141,7 +141,9 @@ def process_single_halo(
                     # Calculation caused an unexpected error.
                     # Output the halo ID so we can debug this.
                     print(
-                        f"Object with HaloCatalogueIndex={input_halo['index']} encountered an error"
+                        f"Object with HaloCatalogueIndex={input_halo['index']} encountered an error "
+                        f"while calculating {halo_prop.name}",
+                        flush=True,
                     )
                     raise
                 else:

SOAP/particle_selection/SO_properties.py

@@ -3674,7 +3674,14 @@ def calculate(
                         unit = unit * unyt.Unit("a", registry=registry) ** a_exponent
                     if do_calculation[filter_name]:
                         t0_calc = time.time()
-                        val = getattr(part_props, name)
+                        try:
+                            val = getattr(part_props, name)
+                        except Exception as e:
+                            e.add_note(
+                                f"Error calculating {prop.name} for subhalo "
+                                f"{input_halo['index']}"
+                            )
+                            raise
                         if val is not None:
                             assert (
                                 SO[name].shape == val.shape

SOAP/particle_selection/aperture_properties.py

@@ -4176,7 +4176,14 @@ def calculate(
                     unit = unit * unyt.Unit("a", registry=registry) ** a_exponent
                 if do_calculation[filter_name]:
                     t0_calc = time.time()
-                    val = getattr(part_props, name)
+                    try:
+                        val = getattr(part_props, name)
+                    except Exception as e:
+                        e.add_note(
+                            f"Error calculating {prop.name} for subhalo "
+                            f"{input_halo['index']}"
+                        )
+                        raise
                     if val is not None:
                         assert (
                             aperture_sphere[name].shape == val.shape

SOAP/particle_selection/projected_aperture_properties.py

@@ -1950,7 +1950,14 @@ def calculate(
                         unit = unit * unyt.Unit("a", registry=registry) ** a_exponent
                     if do_calculation[filter_name]:
                         t0_calc = time.time()
-                        val = getattr(proj_part_props, name)
+                        try:
+                            val = getattr(proj_part_props, name)
+                        except Exception as e:
+                            e.add_note(
+                                f"Error calculating {prop.name} ({projname}) "
+                                f"for subhalo {input_halo['index']}"
+                            )
+                            raise
                         if val is not None:
                             assert (
                                 projected_aperture[projname][name].shape == val.shape

SOAP/particle_selection/subhalo_properties.py

@@ -2678,7 +2678,14 @@ def calculate(self, input_halo, search_radius, data, halo_result):
             )
             if do_calculation[filter_name]:
                 t0_calc = time.time()
-                val = getattr(part_props, name)
+                try:
+                    val = getattr(part_props, name)
+                except Exception as e:
+                    e.add_note(
+                        f"Error calculating {prop.name} for subhalo "
+                        f"{input_halo['index']}"
+                    )
+                    raise
                 if val is not None:
                     assert (
                         subhalo[name].shape == val.shape

SOAP/property_calculation/inertia_tensors.py

@@ -335,6 +335,9 @@ def get_weighted_projected_inertia_tensor(
         tensor = tensor.sum(axis=0)
         eig_val, eig_vec = np.linalg.eigh(tensor.value)
 
+        # Handle overflows into negative values of very small eigenvalues
+        eig_val = np.abs(eig_val)
+
         # Handle cases where there is only one particle after iterating.
         if q == 0:
             tensor.fill(0)