Changeset 1325

Timestamp:

08/14/24 02:48:39 (5 weeks ago)

Author:

Maciej Komosinski

Message:

Minor optimizations

Location:

framspy/dissimilarity

Files:

• alignmodel.py (modified) (2 diffs)
• density_distribution.py (modified) (13 diffs)

framspy/dissimilarity/alignmodel.py

@@ -2 +2 @@
 import numpy as np
 
+
 def wcentre(matrix, weights):
-    sw = weights.sum()
-    swx = (matrix*weights).sum(axis=1)
-    swx /= sw
-    return (matrix.transpose()-swx).transpose()*np.sqrt(weights)
-    
+        sw = weights.sum()
+        swx = (matrix * weights).sum(axis=1)
+        swx /= sw
+        return (matrix.transpose() - swx).transpose() * np.sqrt(weights)
+
+
 def weightedMDS(distances, weights):
-    n = len(weights)
-    distances = distances**2
-    for i in range(2):
-        distances = wcentre(distances, weights)
-        distances = distances.T
-    distances *= -0.5
-    _, eigenvalues, vh = np.linalg.svd(distances)
-    W = (vh/np.sqrt(weights)).T
-    S = np.zeros((n,n))
-    np.fill_diagonal(S, eigenvalues)
-    S = S**0.5
-    dcoords = W.dot(S)
-    coords = np.zeros((n, 3))
-    coords[:,0]=dcoords[:,0]
-    for i in range(1,3):
-        if n>i:
-            coords[:,i]=dcoords[:,i]
-    return coords
-    
+        n = len(weights)
+        distances = distances ** 2
+        for i in range(2):
+                distances = wcentre(distances, weights)
+                distances = distances.T
+        distances *= -0.5
+        _, eigenvalues, vh = np.linalg.svd(distances)
+        W = (vh / np.sqrt(weights)).T
+        S = np.zeros((n, n))
+        np.fill_diagonal(S, eigenvalues)
+        S = S ** 0.5
+        dcoords = W.dot(S)
+        coords = np.zeros((n, 3))
+        coords[:, 0] = dcoords[:, 0]
+        for i in range(1, 3):
+                if n > i:
+                        coords[:, i] = dcoords[:, i]
+        return coords
+
 
 def align(model, fixedZaxis=False):
-        numparts=model.numparts._value()
+        numparts = model.numparts._value()
         distmatrix = np.zeros((numparts, numparts), dtype=float)
         for p1 in range(numparts):
-                for p2 in range(numparts): #TODO optimize, only calculate a triangle
-                        P1=model.getPart(p1)
-                        P2=model.getPart(p2)
+                for p2 in range(p1 + 1, numparts):  # only calculate a triangle since Euclidean distance is symmetrical
+                        P1 = model.getPart(p1)
+                        P2 = model.getPart(p2)
                         if fixedZaxis:
-                                #fixed vertical axis, so pretend all points are on the xy plane
+                                # fixed vertical axis, so pretend all points are on the xy plane
                                 z_dist = 0
                         else:
-                                z_dist = (P1.z._value()-P2.z._value())**2
-                        distmatrix[p1,p2]=math.sqrt((P1.x._value()-P2.x._value())**2+(P1.y._value()-P2.y._value())**2+z_dist)
-        
+                                z_dist = (P1.z._value() - P2.z._value()) ** 2
+                        distmatrix[p1, p2] = distmatrix[p2, p1] = math.sqrt((P1.x._value() - P2.x._value()) ** 2 + (P1.y._value() - P2.y._value()) ** 2 + z_dist)
+
         if model.numjoints._value() > 0:
-                weightvector=np.zeros((numparts), dtype=int)
+                weightvector = np.zeros((numparts), dtype=int)
         else:
-                weightvector=np.ones((numparts), dtype=int)
-        
+                weightvector = np.ones((numparts), dtype=int)
+
         for j in range(model.numjoints._value()):
-                J=model.getJoint(j)
-                weightvector[J.p1._value()]+=1
-                weightvector[J.p2._value()]+=1
-        weightvector=weightvector.astype(float) # convert to float once, since later it would be promoted to float so many times anyway...
+                J = model.getJoint(j)
+                weightvector[J.p1._value()] += 1
+                weightvector[J.p2._value()] += 1
+        weightvector = weightvector.astype(float)  # convert to float once, since later it would be promoted to float so many times anyway...
         coords = weightedMDS(distmatrix, weightvector)
 
@@ -66 +68 @@
                                 P.z = coords[p, 2]
 
-
         if fixedZaxis:
                 if np.shape(coords)[1] > 2:
-                #restore original z coordinate
+                        # restore original z coordinate
                         for p in range(numparts):
-                                P=model.getPart(p)
-                                coords[p,2]=P.z._value()
-
+                                P = model.getPart(p)
+                                coords[p, 2] = P.z._value()

framspy/dissimilarity/density_distribution.py

@@ -41 +41 @@
 
 
-    def calculateNeighberhood(self,array,mean_coords):
+    def calculateNeighborhood(self,array,mean_coords):
         """ Calculates number of elements for given sample and set ups the center of this sample
         to the center of mass (calculated by mean of every coordinate)
@@ -55 +55 @@
         weight = len(array)
         if weight > 0:
-            point = [np.mean(array[:,0]),np.mean(array[:,1]),np.mean(array[:,2])]
+            point = np.mean(array, axis=0)  # equivalent to [np.mean(array[:,0]),np.mean(array[:,1]),np.mean(array[:,2])]
             return weight, point
         else:
@@ -62 +62 @@
 
     def calculateDistPoints(self,point1, point2):
-        """ Returns euclidean distance between two points
+        """ Returns Euclidean distance between two points
         Args (distribution):
             point1 ([float,float,float]) - coordinates of first point
@@ -71 +71 @@
 
         Returns:
-            [float]: euclidean distance
+            [float]: Euclidean distance
         """
         if self.frequency:
-            return abs(point1-point2)
+            return abs(point1-point2) # TODO vector instead of scalar returned?
         else:
-            return np.sqrt(np.sum(np.square(point1-point2)))
+            return np.linalg.norm(point1-point2, ord=2)
 
 
@@ -107 +107 @@
         """
         lens = len(s1)
-        indices = []
-        for i in range(lens):
-            if s1[i]==0 and s2[i]==0:
-                    indices.append(i)
+        indices = [i for i in range(lens) if s1[i]==0 and s2[i]==0]
 
         return np.delete(s1, indices), np.delete(s2, indices)
@@ -126 +123 @@
         """
         lens = len(s1[0])
-        indices = []
-        for i in range(lens):
-            if s1[1][i]==0 and s2[1][i]==0:
-                indices.append(i)
+        indices = [i for i in range(lens) if s1[1][i]==0 and s2[1][i]==0]
 
         s1 = [np.delete(s1[0], indices, axis=0), np.delete(s1[1], indices, axis=0)]
@@ -168 +162 @@
                         feature_array.append(len(array[rows]))
                     else:
-                        weight, point = self.calculateNeighberhood(array[rows],[edges_x[x]+step_x_half,edges_y[y]+step_y_half,edges_z[z]+step_z_half])
+                        weight, point = self.calculateNeighborhood(array[rows],[edges_x[x]+step_x_half,edges_y[y]+step_y_half,edges_z[z]+step_z_half])
                         feature_array.append(point)
                         weight_array.append(weight)
@@ -180 +174 @@
 
     def getSignaturesForPair(self,array1,array2):
-        """Generates signatures for given pair of models represented by array of voxels.
+        """Generates signatures for a given pair of models represented by array of voxels.
         We calculate space for given models by taking the extremas for each axis and dividing the space by the resolution.
         This divided space generate us samples which contains points. Each sample will have new coordinates which are mean of all points from it and weight which equals to the number of points.
@@ -248 +242 @@
 
     def calculateDissimforVoxels(self, voxels1, voxels2):
-        """Calculates EMD for pair of voxels representing models.
+        """Calculates EMD for a pair of voxels representing models.
         Args:
             voxels1 np.array([np.array(,dtype=float)]: list of voxels representing model1.
@@ -254 +248 @@
 
         Returns:
-            float: dissim for pair of list of voxels
+            float: dissim for a pair of list of voxels
         """
         numvox1 = len(voxels1)
@@ -278 +272 @@
         if self.metric == 'l1':
             if self.frequency:
-                out = np.linalg.norm((s1-s2), ord=1)
+                out = np.linalg.norm(s1-s2, ord=1)
             else:
-                out = np.linalg.norm((s1[1]-s2[1]), ord=1)
+                out = np.linalg.norm(s1[1]-s2[1], ord=1)
 
         elif self.metric == 'l2':
             if self.frequency:
-                out = np.linalg.norm((s1-s2))
+                out = np.linalg.norm(s1-s2)
             else:
-                out = np.linalg.norm((s1[1]-s2[1]))
+                out = np.linalg.norm(s1[1]-s2[1])
 
         elif self.metric == 'emd':
@@ -320 +314 @@
 
         Returns:
-            float: dissim for pair of strings representing models.
+            float: dissim for a pair of strings representing models.
         """     
 
@@ -349 +343 @@
         listOfVoxels = [self.getVoxels(g) for g in listOfGeno]
         for i in range(numOfGeno):
-            for j in range(numOfGeno):
+            for j in range(numOfGeno): # could only calculate a triangle if the definition of similarity and its calculation guarantees symmetry
                 dissimMatrix[i,j] = self.calculateDissimforVoxels(listOfVoxels[i], listOfVoxels[j])
         return dissimMatrix