Merge pull request #192 from francois-drielsma/develop

Multiple bug fixes for 2x2 prod.

analysis/classes/Interaction.py

@@ -50,6 +50,7 @@ class Interaction:
 
     # Attributes that specify coordinates
     _COORD_ATTRS = ['points', 'vertex']
+    _TRUTH_COORD_ATTRS = []
 
     def __init__(self,
                  interaction_id: int = -1,
@@ -335,7 +336,8 @@ def convert_to_cm(self, meta):
         '''
         assert self._units == 'px'
         for attr in self._COORD_ATTRS:
-            setattr(self, attr, pixel_to_cm(getattr(self, attr), meta))
+            center = not attr in self._TRUTH_COORD_ATTRS
+            setattr(self, attr, pixel_to_cm(getattr(self, attr), meta, center=center))
 
         self._units = 'cm'
 

analysis/classes/Particle.py

@@ -85,6 +85,7 @@ class Particle:
 
     # Attributes that specify coordinates
     _COORD_ATTRS = ['points', 'start_point', 'end_point']
+    _TRUTH_COORD_ATTRS = []
 
     def __init__(self,
                  group_id: int = -1,
@@ -465,7 +466,8 @@ def convert_to_cm(self, meta):
         '''
         assert self._units == 'px'
         for attr in self._COORD_ATTRS:
-            setattr(self, attr, pixel_to_cm(getattr(self, attr), meta))
+            center = not attr in self._TRUTH_COORD_ATTRS
+            setattr(self, attr, pixel_to_cm(getattr(self, attr), meta, center=center))
         self._units = 'cm'
 
     @property

analysis/classes/ParticleFragment.py

@@ -44,7 +44,9 @@ class ParticleFragment:
         compose a particle.
     '''
 
+    # Attributes that specify coordinates
     _COORD_ATTRS = ['points', 'start_point', 'end_point']
+    _TRUTH_COORD_ATTRS = []
 
     def __init__(self,
                  fragment_id: int = -1,
@@ -238,7 +240,8 @@ def convert_to_cm(self, meta):
         '''
         assert self._units == 'px'
         for attr in self._COORD_ATTRS:
-            setattr(self, attr, pixel_to_cm(getattr(self, attr), meta))
+            center = not attr in self._TRUTH_COORD_ATTRS
+            setattr(self, attr, pixel_to_cm(getattr(self, attr), meta, center=center))
         self._units = 'cm'
 
     @property

analysis/classes/TruthInteraction.py

@@ -8,7 +8,6 @@
 from . import Interaction, TruthParticle
 from .Interaction import _process_interaction_attributes
 
-from mlreco.utils import pixel_to_cm
 from mlreco.utils.globals import PID_LABELS
 from mlreco.utils.decorators import inherit_docstring
 
@@ -45,6 +44,7 @@ class TruthInteraction(Interaction):
     # Attributes that specify coordinates
     _COORD_ATTRS = Interaction._COORD_ATTRS +\
 	['truth_points', 'sed_points', 'truth_vertex']
+    _TRUTH_COORD_ATTRS = ['truth_vertex']
 
     # Define placeholder values (-np.inf for float, -sys.maxsize for int)
     _SCALAR_KEYS = {'bjorken_x': -np.inf, 

analysis/classes/TruthParticle.py

@@ -3,7 +3,6 @@
 
 from . import Particle
 
-from mlreco.utils import pixel_to_cm
 from mlreco.utils.globals import PDG_TO_PID, TRACK_SHP, SHAPE_LABELS, \
         PID_LABELS, PID_MASSES
 from mlreco.utils.decorators import inherit_docstring
@@ -45,6 +44,8 @@ class TruthParticle(Particle):
     _COORD_ATTRS = Particle._COORD_ATTRS +\
 	['truth_points', 'sed_points', 'position', 'end_position',\
 	'parent_position', 'ancestor_position', 'first_step', 'last_step']
+    _TRUTH_COORD_ATTRS = ['position', 'end_position', 'parent_position',
+                          'ancestor_position', 'first_step', 'last_step']
 
     def __init__(self,
                  *args,

analysis/classes/TruthParticleFragment.py

@@ -1,7 +1,6 @@
 import numpy as np
 from typing import Counter, List, Union
 
-from mlreco.utils import pixel_to_cm
 from mlreco.utils.globals import PDG_TO_PID, TRACK_SHP, SHAPE_LABELS, PID_LABELS
 from mlreco.utils.decorators import inherit_docstring
 
@@ -27,6 +26,8 @@ class TruthParticleFragment(ParticleFragment):
     _COORD_ATTRS = ParticleFragment._COORD_ATTRS + \
 	['truth_points', 'sed_points', 'position', 'end_position', \
 	'parent_position', 'ancestor_position', 'first_step', 'last_step']
+    _TRUTH_COORD_ATTRS = ['position', 'end_position', 'parent_position',
+                          'ancestor_position', 'first_step', 'last_step']
 
     def __init__(self, 
                  *args,

analysis/classes/builders.py

@@ -411,7 +411,6 @@ def _build_truth(self,
             simE_deposits   = None
 
         # point_labels   = data['point_labels'][entry]
-        # unit_convert = lambda x: pixel_to_cm_1d(x, meta) if self.convert_to_cm == True else x
 
         # For debugging
         voxel_counts = 0
@@ -926,10 +925,3 @@ def match_points_to_particles(ppn_points : np.ndarray,
         dist = cdist(ppn_coords, particle.points)
         matches = ppn_points_type[dist.min(axis=1) < ppn_distance_threshold]
         particle.ppn_candidates = matches.reshape(-1, 7)
-
-def pixel_to_cm_1d(vec, meta):
-    out = np.zeros_like(vec)
-    out[0] = meta[0] + meta[6] * vec[0]
-    out[1] = meta[1] + meta[7] * vec[1]
-    out[2] = meta[2] + meta[8] * vec[2]
-    return out

analysis/classes/matching.py

@@ -290,7 +290,8 @@ def match_particles_all(particles_x : Union[List[Particle], List[TruthParticle]]
                 key = (px.id, matched.id)
                 matches[key] = (px, matched)
 
-    out_counts = np.array(out_counts)
+    out_counts = np.empty(len(out_counts), dtype=object)
+    out_counts[:] = out_counts
 
     return matches, out_counts
 

analysis/manager.py

@@ -149,9 +149,9 @@ def initialize_base(self,
         # Load the full chain configuration, if it is provided
         self.chain_config = chain_config
         if chain_config is not None:
-            #cfg = yaml.safe_load(open(chain_config, 'r').read())
-            process_config(chain_config, verbose=False)
-            self.chain_config = chain_config
+            cfg = yaml.safe_load(open(chain_config, 'r').read())
+            process_config(cfg, verbose=False)
+            self.chain_config = cfg
 
         # Initialize data product builders
         self.builders = {}
@@ -366,6 +366,7 @@ def convert_pixels_to_cm(self, data, result):
             'input_data', 'segment_label',
             'particles_label', 'cluster_label', 'kinematics_label', 'sed'
         ])
+
         result_has_voxels = set([
             'input_rescaled',
             'cluster_label_adapted',
@@ -389,10 +390,10 @@ def convert_pixels_to_cm(self, data, result):
 
         for key, val in data.items():
             if key in data_has_voxels:
-                data[key] = [self._pixel_to_cm(arr, meta) for arr in val]
+                data[key] = [self._pixel_to_cm(arr, meta, center=True) for arr in val]
         for key, val in result.items():
             if key in result_has_voxels:
-                result[key] = [self._pixel_to_cm(arr, meta) for arr in val]
+                result[key] = [self._pixel_to_cm(arr, meta, center=True) for arr in val]
             if key in data_products:
                 for plist in val:
                     for p in plist:
@@ -673,7 +674,7 @@ def _set_iteration(self, dataset):
         assert self.max_iteration <= len(dataset)
 
     @staticmethod
-    def _pixel_to_cm(arr, meta):
+    def _pixel_to_cm(arr, meta, center=False):
         '''
         Converts tensor pixel coordinates to detector coordinates
 
@@ -683,6 +684,9 @@ def _pixel_to_cm(arr, meta):
             Tensor of which to convert the coordinate columns
         meta : np.ndarray
             Metadata information to operate the translation
+        center : bool, default False
+            Whether to place the coordinates at the center of the pixel or not.
+            Provides a unbiased estimate for true pixel coordinates
         '''
-        arr[:, COORD_COLS] = pixel_to_cm(arr[:, COORD_COLS], meta)
+        arr[:, COORD_COLS] = pixel_to_cm(arr[:, COORD_COLS], meta, center=center)
         return arr

mlreco/iotools/parsers/cluster.py

@@ -70,6 +70,7 @@ def parse_cluster3d(cluster_event,
                     type_include_secondary = False,
                     primary_include_mpr = True,
                     break_clusters = False,
+                    break_clusters_eps = 1.1,
                     min_size = -1):
     """
     a function to retrieve a 3D clusters tensor
@@ -92,6 +93,7 @@ def parse_cluster3d(cluster_event,
               type_include_secondary: false
               primary_include_mpr: true
               break_clusters: false
+              break_clusters_eps: 1.1
 
     Configuration
     -------------
@@ -107,6 +109,7 @@ def parse_cluster3d(cluster_event,
     type_include_secondary: bool
     primary_include_mpr: bool
     break_clusters: bool
+    break_clusters_eps: float
 
     Returns
     -------
@@ -142,7 +145,7 @@ def parse_cluster3d(cluster_event,
     if add_particle_info:
         assert particle_event is not None,\
                 'Must provide particle tree if particle information is included'
-        num_particles = particle_event.size()
+        num_particles = particle_event.as_vector().size()
         assert num_particles == num_clusters or num_particles == num_clusters-1,\
                 'The number of particles must be aligned with the number of clusters'
 
@@ -191,7 +194,7 @@ def parse_cluster3d(cluster_event,
 
             # If requested, break cluster into pieces that do not touch each other
             if break_clusters:
-                dbscan = DBSCAN(eps=1.1, min_samples=1, metric='chebyshev')
+                dbscan = DBSCAN(eps=break_clusters_eps, min_samples=1, metric='chebyshev')
                 frag_labels = np.unique(dbscan.fit(voxels).labels_, return_inverse=True)[-1]
                 features[1] = id_offset + frag_labels
                 id_offset += max(frag_labels) + 1

mlreco/iotools/readers.py

@@ -374,7 +374,7 @@ def load_key(self, in_file, event, data_blob, result_blob, key, nested):
             # If the reference points at a group, unpack
             el_refs = group[key]['index'][region_ref].flatten()
             if len(group[key]['index'].shape) == 1:
-                ret = np.empty(len(el_refs), dtype=np.object)
+                ret = np.empty(len(el_refs), dtype=object)
                 ret[:] = [group[key]['elements'][r] for r in el_refs]
                 if len(group[key]['elements'].shape) > 1:
                     for i in range(len(el_refs)):

mlreco/iotools/writers.py

@@ -101,6 +101,7 @@ def __init__(self,
                  skip_input_keys: list = [],
                  result_keys: list = None,
                  skip_result_keys: list = [],
+                 dummy_keys: list = [],
                  append_file: bool = False,
                  merge_groups: bool = False):
         '''
@@ -129,6 +130,7 @@ def __init__(self,
         self.skip_input_keys  = skip_input_keys
         self.result_keys      = result_keys
         self.skip_result_keys = skip_result_keys
+        self.dummy_keys       = dummy_keys
         self.append_file      = append_file
         self.merge_groups     = merge_groups
         self.ready            = False
@@ -174,6 +176,13 @@ def create(self, data_blob, result_blob=None, cfg=None):
         for key in self.result_keys:
             self.register_key(result_blob, key, 'result')
 
+        # If requested, add dummy datasets for some requested keys
+        for key in self.dummy_keys:
+            assert key not in self.key_dict, (
+                    "Dummy key exists in the requested keys already, abort.")
+            dummy_blob = {key: [np.empty(0, dtype=np.float32)]}
+            self.register_key(dummy_blob, key, 'result')
+
         # Initialize the output HDF5 file
         with h5py.File(self.file_name, 'w') as out_file:
             # Initialize the info dataset that stores top-level description of what is stored
@@ -250,7 +259,7 @@ def register_key(self, blob, key, category):
                     # List containing a single list of scalars per batch ID
                     self.key_dict[key]['dtype'] = type(ref_obj)
 
-                elif not isinstance(blob[key][ref_id], list) and not blob[key][ref_id].dtype == np.object:
+                elif not isinstance(blob[key][ref_id], list) and not blob[key][ref_id].dtype == object:
                     # List containing a single ndarray of scalars per batch ID
                     self.key_dict[key]['dtype'] = blob[key][ref_id].dtype
                     self.key_dict[key]['width'] = blob[key][ref_id].shape[1] if len(blob[key][ref_id].shape) == 2 else 0
@@ -433,8 +442,14 @@ def append(self, data_blob=None, result_blob=None, cfg=None):
             self.create(data_blob, result_blob, cfg)
             self.ready = True
 
-        # Append file
+        # Create a dummy blob to fill dummy keys with
         self.batch_size = len(data_blob['index'])
+        if self.dummy_keys:
+            dummy_blob = {}
+            for key in self.dummy_keys:
+                dummy_blob[key] = [np.empty(0, dtype=np.float32) for b in range(self.batch_size)]
+
+        # Append file
         with h5py.File(self.file_name, 'a') as out_file:
             # Loop over batch IDs
             for batch_id in range(self.batch_size):
@@ -448,6 +463,8 @@ def append(self, data_blob=None, result_blob=None, cfg=None):
                     self.append_key(out_file, event, data_blob, key, batch_id)
                 for key in self.result_keys:
                     self.append_key(out_file, event, result_blob, key, batch_id)
+                for key in self.dummy_keys:
+                    self.append_key(out_file, event, dummy_blob, key, batch_id)
 
                 # Append event
                 event_id  = len(out_file['events'])

mlreco/utils/globals.py

@@ -60,9 +60,10 @@
 }
 
 # Invalid larcv.Particle labels
-INVAL_ID   = 9223372036854775807 # larcv.kINVALID_INSTANCEID
-INVAL_TID  = 4294967295          # larcv.kINVALID_UINT
-INVAL_PDG  = 0                   # Invalid particle PDG code
+OLD_INVAL_ID = 65535               # larcv.kINVALID_INSTANCEID, pre-update
+INVAL_ID     = 9223372036854775807 # larcv.kINVALID_INSTANCEID
+INVAL_TID    = 4294967295          # larcv.kINVALID_UINT
+INVAL_PDG    = 0                   # Invalid particle PDG code
 
 # Particle ID of each recognized particle species
 PHOT_PID = 0