Cleanup

snudda/place/bend_morphologies.py

@@ -23,82 +23,17 @@ def check_if_inside(self, morphology: NeuronMorphologyExtended):
 
         return inside_flag
 
-    def bend_morphology_OLD(self, morphology: NeuronMorphologyExtended, k=50e-6):
-
-        # TODO: Parent point idx is included if parent section is of same type as section!
-        #       So we should not rotate the first point!
-
-        # k -- decay constant
-        n_random = 10
-        candidate_pos = np.zeros((n_random, 3))
-
-        parent_rotation_matrices = dict()  # indexed by parent_point_idx in SectionMetaData
-
-        # Iterate over all parts of neuron
-        for section in morphology.section_iterator():
-
-            # We need to track the rotation of each point, in particular save rotations
-            # for each branch point, so that all children can start with that rotation
-
-            if section.parent_section_idx in parent_rotation_matrices:
-                rotation_matrix = parent_rotation_matrices[section.parent_section_idx]
-            else:
-                rotation_matrix = np.eye(3)
-
-            parent_coords = section.morphology_data.position[section.parent_section_idx, :]
-
-            # Loop over all points in section
-            for point_idx in section.point_idx:
-                coords = section.morphology_data.geometry[point_idx, :]
-                dist = self.region_mesh.distance_to_border(points=coords)
-
-                P = 1 / (1 + np.exp(-dist/k))
-
-                if self.rng.uniform(1) < P:
-                    # We need to randomize new rotation matrix
-                    # https://docs.scipy.org/doc/scipy/reference/generated/scipy.spatial.transform.Rotation.html
-                    angles = self.rng.uniform(size=(n_random, 3), low=-0.2, high=0.2)  # Angles in radians
-                    rots = Rotation.from_euler(seq="XYZ", angles=angles)
-
-                    for idx, rot in enumerate(rots):
-                        candidate_pos[idx, :] = parent_coords[idx, :] + rots.apply(vectors=coords-parent_coords)
-
-                    candidate_dist = self.region_mesh.distance_to_border(points=candidate_pos)
-
-                    P_candidate = np.divide(1, 1 + np.exp(-candidate_dist/k))
-                    picked_idx = self.rng.choice(n_random, p=P_candidate)
-
-                    new_coords = candidate_pos[picked_idx, :]
-
-                    picked_rotation = rots[picked_idx]
-                    rotation_point = parent_coords
-
-                    # TODO: We also need to apply this rotation to ALL the other points in this section,
-                    #       and to all child sections of this branch.
-
-            else:
-                # Keep the old coords
-                new_coords = coords
-
-
-
-                # Old point is next points parent
-                parent_coords = new_coords
-
-            # Gradient is calculated based on distance to mesh
-            # Calculate delta_gradient = gradient_self - gradient_parent
-            # Amount of angle to bend f(delta_gradient)
-            # Randomize the rotation matrix
-
-            # After section is done, store the last rotation matrix, so its children can get their parent rotation
-            # store in dictionary?
-
-            parent_rotation_matrices[point_idx] = rotation_matrix
-
     def bend_morphology(self, morphology: NeuronMorphologyExtended, k=30e-6, n_random=5):
 
         # k -- how early will the neuron start bending when it approaches the border
 
+        # Check distance to border for all points, negative distance means inside
+        all_original_dist = self.region_mesh.distance_to_border(morphology.geometry[:, :3])
+        if (all_original_dist < 0).all():
+            # Morphology entirely inside mesh, nothing to do
+            return None, False
+
+        # We randomize n_random points, and store candidates here, and pick "best" one
         candidate_pos = np.zeros((n_random, 3))
 
         parent_direction = dict()
@@ -121,12 +56,14 @@ def bend_morphology(self, morphology: NeuronMorphologyExtended, k=30e-6, n_rando
                 else:
                     parent_point = np.zeros((3, ))
 
-                parent_dist = self.region_mesh.distance_to_border(points=parent_point.reshape((1,3)))[0]
+                parent_dist = self.region_mesh.distance_to_border(points=parent_point.reshape((1, 3)))[0]
                 parent_moved = False
 
             rot_rep = old_rotation_representation[section.section_id, section.section_type]
             new_rot_rep = []
 
+            section_dist = all_original_dist[section.point_idx]
+
             # Loop over all points in section
             for idx, (rotation, length) in enumerate(rot_rep):
 
@@ -142,7 +79,12 @@ def bend_morphology(self, morphology: NeuronMorphologyExtended, k=30e-6, n_rando
                 # TODO: This check can be cached (assuming no parent segment have been rotated, track that)
 
                 # Check if point is too close to edge
-                dist = self.region_mesh.distance_to_border(points=putative_point)[0]
+                if parent_moved:
+                    dist = self.region_mesh.distance_to_border(points=putative_point)[0]
+                else:
+                    # Parent has not moved, so use stored original distance
+                    dist = section_dist[idx]
+
                 P_move = 1 / (1 + np.exp(-dist/k))
 
                 # Cache the random numbers for segments in the section...