Automatically detect properties which can be computed

SOAP/compute_halo_properties.py

@@ -164,6 +164,9 @@ def compute_halo_properties():
     cellgrid.snapshot_datasets.setup_defined_constants(
         parameter_file.get_defined_constants()
     )
+    # Tell the parameter file which datasets are in the input files, so that
+    # properties which cannot be computed can be skipped or reported
+    parameter_file.set_available_datasets(cellgrid.snapshot_datasets.datasets_in_file)
     parameter_file.record_property_timings = args.record_property_timings
 
     # Try to load parameters for RecentlyHeatedGasFilter. If a property that uses the
@@ -473,6 +476,7 @@ def compute_halo_properties():
         if args.record_property_timings:
             print("Storing processing time for each property")
         parameter_file.print_unregistered_properties(halo_prop_list, dmo=args.dmo)
+        parameter_file.print_skipped_properties(halo_prop_list, dmo=args.dmo)
         parameter_file.print_invalid_properties(halo_prop_list)
         parameter_file.print_variation_warnings()
         if not parameter_file.renclose_enabled():
@@ -481,6 +485,12 @@ def compute_halo_properties():
             )
         category_filter.print_filters()
 
+        # Properties enabled in the parameter file must be computed, so abort
+        # if the input files do not contain the datasets they require
+        parameter_file.print_uncomputable_properties()
+        if len(parameter_file.uncomputable_properties):
+            comm_world.Abort(1)
+
     # Ensure output dir exists
     if comm_world_rank == 0:
         try:

SOAP/core/parameter_file.py

@@ -15,6 +15,28 @@
 
 from SOAP import property_table
 
+# Lazily built map from the output name of a property to its entry in the
+# property table. The property table itself is keyed on an internal name.
+_PROPERTY_BY_NAME = None
+
+
+def _property_by_name(name: str):
+    """
+    Look up a property in the property table by its output name.
+
+    Returns None if there is no property with this name, which happens for
+    properties that are not defined by the property table (e.g. the dummy
+    properties used when generating the documentation).
+    """
+    global _PROPERTY_BY_NAME
+    if _PROPERTY_BY_NAME is None:
+        _PROPERTY_BY_NAME = {
+            prop.name: prop
+            for prop in property_table.PropertyTable.full_property_list.values()
+        }
+    return _PROPERTY_BY_NAME.get(name)
+
+
 # Known parameter file structure, used by check_schema to flag typos. A value
 # of None means the keys directly under that section are user-named or free-form
 # and are not checked; a set lists the only keys allowed directly under that
@@ -104,6 +126,24 @@ def __init__(
         # on a single rank via print_variation_warnings()
         self.variation_warnings = []
 
+        # Datasets present in the input files, as {particle type: set of dataset
+        # names}. While this is None no property is considered uncomputable,
+        # which is the case when building a parameter file from a dictionary.
+        self.available_datasets = None
+
+        # Properties which are not in the parameter file, and which are not
+        # calculated because the input files lack the datasets they need.
+        # Only used when calculate_missing_properties is True.
+        self.skipped_properties = set()
+
+        # Properties which are enabled in the parameter file, but which cannot
+        # be calculated because the input files lack the datasets they need.
+        self.uncomputable_properties = {}
+
+        # Names of the properties used by the filters defined in the parameter
+        # file, generated on demand by _filter_property_names()
+        self.filter_property_names = None
+
     def get_parameters(self) -> Dict:
         """
         Get a copy of the parameter dictionary.
@@ -148,6 +188,63 @@ def _validate_filter_name(self, filter_name, context: str) -> None:
                 f'the "filters" section of the parameter file'
             )
 
+    def set_available_datasets(self, datasets_in_file: Dict) -> None:
+        """
+        Record which datasets are present in the input files
+
+        Parameters:
+         - datasets_in_file: Dict
+           Dictionary of the datasets present in the snapshot and extra-input
+           files, as {particle type: set of dataset names}.
+        """
+        self.available_datasets = datasets_in_file
+
+    def missing_datasets(self, property_name: str) -> List[str]:
+        """
+        Get the datasets which are required to compute the given property, but
+        which are not present in the input files.
+
+        Particle types which are absent from the input files entirely are not
+        considered. This matches the check done in SWIFTCellGrid.check_datasets_exist.
+
+        Returns the (aliased) names of the missing datasets. An empty list is
+        returned if the property can be computed.
+        """
+        if self.available_datasets is None:
+            return []
+        prop = _property_by_name(property_name)
+        if prop is None:
+            return []
+        missing = []
+        for dataset in prop.particle_properties:
+            ptype, name = self.get_particle_property(dataset)
+            # Skip particle types which are not in the input files at all
+            if ptype not in self.available_datasets:
+                continue
+            dataset_name = f"{ptype}/{name}"
+            if (name not in self.available_datasets[ptype]) and (
+                dataset_name not in missing
+            ):
+                missing.append(dataset_name)
+        return missing
+
+    def _filter_property_names(self) -> set:
+        """
+        Get the names of the properties used by the filters defined in the
+        parameter file.
+
+        These properties are never skipped, since the filters that use them
+        would otherwise be undefined.
+        """
+        if self.filter_property_names is None:
+            self.filter_property_names = set()
+            for filter_info in self.get_filters().values():
+                for prop in filter_info.get("properties", []):
+                    # Filters name properties by their full path in the output,
+                    # e.g. BoundSubhalo/NumberOfGasParticles
+                    self.filter_property_names.add(prop.split("/")[-1])
+        return self.filter_property_names
+
     def get_property_filters(self, base_halo_type: str, full_list: List[str]) -> Dict:
         """
         Get a dictionary with the filter that should be applied to each
@@ -176,18 +273,18 @@ def get_property_filters(self, base_halo_type: str, full_list: List[str]) -> Dic
 
         if not base_halo_type in self.parameters:
             self.parameters[base_halo_type] = {}
-        # Handle the case where no properties are listed for the halo type
+        # Handle the case where no properties are listed for the halo type. Each
+        # property is then treated as if it were missing from the parameter file.
         if not "properties" in self.parameters[base_halo_type]:
             self.parameters[base_halo_type]["properties"] = {}
-            for property in full_list:
-                self.parameters[base_halo_type]["properties"][
-                    property
-                ] = self.calculate_missing_properties()
+        listed = self.parameters[base_halo_type]["properties"]
         filters = {}
         for property in full_list:
+            # Datasets this property needs which are not in the input files
+            missing = self.missing_datasets(property)
             # Check if property is listed in the parameter file for this base_halo_type
-            if property in self.parameters[base_halo_type]["properties"]:
-                filter_name = self.parameters[base_halo_type]["properties"][property]
+            if property in listed:
+                filter_name = listed[property]
                 # filter_name will a dict if we want different behaviour
                 # for snapshots/snipshots
                 if isinstance(filter_name, dict):
@@ -200,15 +297,27 @@ def get_property_filters(self, base_halo_type: str, full_list: List[str]) -> Dic
                 if filter_name == True:
                     filter_name = "basic"
                 filters[property] = filter_name
+                # An uncomputable property enabled in the parameter file was asked for
+                # explicitly, so we abort rather than quietly omitting it
+                if filter_name and missing:
+                    self.uncomputable_properties[property] = missing
             # Property is not listed in the parameter file for this base_halo_type
+            elif not self.calculate_missing_properties():
+                filters[property] = False
+            elif missing and property not in self._filter_property_names():
+                # The property was not asked for explicitly and cannot be
+                # computed, so it is skipped. Properties used by a filter are
+                # never skipped, since the filter would then be undefined.
+                filters[property] = False
+                listed[property] = False
+                self.skipped_properties.add(property)
             else:
-                if self.calculate_missing_properties():
-                    filters[property] = "basic"
-                    self.parameters[base_halo_type]["properties"][property] = "basic"
-                    if self.unregistered_parameters is not None:
-                        self.unregistered_parameters.add((base_halo_type, property))
-                else:
-                    filters[property] = False
+                filters[property] = "basic"
+                listed[property] = "basic"
+                if self.unregistered_parameters is not None:
+                    self.unregistered_parameters.add((base_halo_type, property))
+                if missing:
+                    self.uncomputable_properties[property] = missing
             if isinstance(filters[property], str):
                 self._validate_filter_name(
                     filters[property], f"{base_halo_type}/{property}"
@@ -241,12 +350,9 @@ def print_unregistered_properties(
         # In a DMO run, drop properties that will be skipped because they are
         # not DMO properties, so the printed list matches the output
         if dmo and halo_prop_list is not None:
-            dmo_flag = {}
-            for halo_type in halo_prop_list:
-                for prop in halo_type.property_list.values():
-                    dmo_flag[(halo_type.base_halo_type, prop.name)] = prop.dmo_property
+            non_dmo_names = self._non_dmo_property_names(halo_prop_list)
             unregistered = {
-                entry for entry in unregistered if dmo_flag.get(entry, True)
+                entry for entry in unregistered if entry[1] not in non_dmo_names
             }
 
         if len(unregistered):
@@ -256,6 +362,72 @@ def print_unregistered_properties(
             for base_halo_type, property in sorted(unregistered):
                 print(f"  {base_halo_type.ljust(30)}{property}")
 
+    def _non_dmo_property_names(self, halo_prop_list) -> set:
+        """
+        Get the names of the properties which are not calculated in a DMO run.
+
+        Properties which are not found in halo_prop_list are not included, since
+        we cannot tell whether they would be calculated.
+
+        Parameters:
+         - halo_prop_list: List
+           List of the halo property calculations that are enabled.
+        """
+        names = set()
+        for halo_type in halo_prop_list:
+            for prop in halo_type.property_list.values():
+                if not prop.dmo_property:
+                    names.add(prop.name)
+        return names
+
+    def print_skipped_properties(self, halo_prop_list=None, dmo: bool = False) -> None:
+        """
+        Print a list of the properties which are not in the parameter file, and
+        which are not calculated because the input files lack the datasets they
+        need.
+
+        Each property is listed once, rather than once per halo type, since
+        whether a dataset is present does not depend on the halo type.
+
+        In a DMO run the property calculators skip any property that is not
+        flagged as a DMO property.
+        """
+
+        # A property that is explicitly enabled in the parameter file
+        # is reported as uncomputable instead
+        skipped = set(self.skipped_properties) - set(self.uncomputable_properties)
+
+        # In a DMO run, drop properties that would be skipped anyway because
+        # they are not DMO properties
+        if dmo and halo_prop_list is not None:
+            skipped -= self._non_dmo_property_names(halo_prop_list)
+
+        if len(skipped):
+            print(
+                "Not computing the following properties as required datasets are missing:"
+            )
+            for property in sorted(skipped):
+                print(f"  {property}")
+
+    def print_uncomputable_properties(self) -> None:
+        """
+        Print a list of the properties which are enabled in the parameter file,
+        but which cannot be calculated because the input files lack the datasets
+        they need. Only the missing datasets are listed for each property, not
+        all of the datasets it requires.
+        """
+        if not len(self.uncomputable_properties):
+            return
+        print(
+            f"Cannot compute {len(self.uncomputable_properties)} properties enabled in "
+            f"the parameter file (only the missing datasets are listed for each one):",
+            flush=True,
+        )
+        for property in sorted(self.uncomputable_properties):
+            print(f"  {property}", flush=True)
+            for dataset in self.uncomputable_properties[property]:
+                print(f"    {dataset}", flush=True)
+
     def print_invalid_properties(self, halo_prop_list) -> None:
         """
         Print a list of any properties in the parameter file that are not present in

SOAP/core/swift_cells.py

@@ -451,31 +451,50 @@ def verify_extra_input(self, comm):
     def check_datasets_exist(self, required_datasets, halo_prop_list):
         # Check we have all the fields needed for each property
         # Doing it at this point rather than in masked cells since we want
-        # to output a list of properties that require the missing fields
-        for ptype in set(self.ptypes).intersection(set(required_datasets.keys())):
-            for name in required_datasets[ptype]:
+        # to output a list of properties that require the missing fields.
+        # Properties which cannot be computed are normally skipped or reported
+        # when the parameter file is resolved, so reaching this point indicates
+        # a problem with an alias or with a declared dependency.
+        missing_datasets = []
+        for ptype in sorted(
+            set(self.ptypes).intersection(set(required_datasets.keys()))
+        ):
+            for name in sorted(required_datasets[ptype]):
                 # Note that the field names in required_datasets have already had
                 # any aliases applied, so we can check the raw files themselves
                 in_extra = (self.extra_filenames is not None) and (
                     name in self.extra_metadata_combined[ptype]
                 )
                 in_snap = name in self.snap_metadata[ptype]
                 if not (in_extra or in_snap):
-                    dataset = f"{ptype}/{name}"
-                    print(f"The following properties require {dataset}:")
-                    full_property_list = property_table.PropertyTable.full_property_list
-                    for k, v in full_property_list.items():
-                        # Skip property if it doesn't require this dataset
-                        if dataset not in v.particle_properties:
-                            continue
-                        # Only print if the property is being calculated for some halo type
-                        for halo_prop in halo_prop_list:
-                            if halo_prop.property_filters.get(v.name, False):
-                                print(f"  {v.name}")
-                                break
-                    raise KeyError(
-                        f"Can't find required dataset {dataset} in input file(s)!"
-                    )
+                    missing_datasets.append(f"{ptype}/{name}")
+
+        if not missing_datasets:
+            return
+
+        # Report the missing datasets
+        full_property_list = property_table.PropertyTable.full_property_list
+        for dataset in missing_datasets:
+            print(f"The following properties require {dataset}:")
+            n_properties = 0
+            for k, v in full_property_list.items():
+                # Skip property if it doesn't require this dataset
+                if dataset not in v.particle_properties:
+                    continue
+                # Only print if the property is being calculated for some halo type
+                for halo_prop in halo_prop_list:
+                    if halo_prop.property_filters.get(v.name, False):
+                        print(f"  {v.name}")
+                        n_properties += 1
+                        break
+            if n_properties == 0:
+                # No enabled property lists this dataset, so it is one that SOAP
+                # always reads. Disabling properties will not help here.
+                print(f"  (none, {dataset} is required for every calculation)")
+        raise KeyError(
+            "Can't find required dataset(s) "
+            f"{', '.join(missing_datasets)} in input file(s)!"
+        )
 
     def prepare_read(self, ptype, mask):
         """