From 0d65a721d8b6bc581622ebf41ce67c8bfc9c4169 Mon Sep 17 00:00:00 2001
From: Chris Maffeo <cmaffeo2@illinois.edu>
Date: Sat, 8 Feb 2020 11:15:00 -0600
Subject: [PATCH] Added routines for getting bounding box

---
 mrdna/model/nbPot.py  |  4 +++-
 mrdna/segmentmodel.py | 36 ++++++++++++++++++++++++++++--------
 2 files changed, 31 insertions(+), 9 deletions(-)

diff --git a/mrdna/model/nbPot.py b/mrdna/model/nbPot.py
index c26c3b0..71a3c2a 100644
--- a/mrdna/model/nbPot.py
+++ b/mrdna/model/nbPot.py
@@ -80,7 +80,9 @@ class AbstractNbDnaScheme(NonbondedScheme):
         if self.debye_length != self.debye_length0:
             ## units "e**2/(4 * pi * 80 epsilon0 AA)" kcal_mol
             A = 4.1507964
-            du = (A/x) * (np.exp(-x/self.debye_length) - np.exp(-x/self.debye_length0))
+            with np.errstate(divide='ignore',invalid='ignore'):
+                du = (A/x) * (np.exp(-x/self.debye_length) - np.exp(-x/self.debye_length0))
+            
             du[x < 10] = du[x>=10][0]
             u = u + du
 
diff --git a/mrdna/segmentmodel.py b/mrdna/segmentmodel.py
index 8b85c3f..22f205d 100644
--- a/mrdna/segmentmodel.py
+++ b/mrdna/segmentmodel.py
@@ -453,6 +453,16 @@ class Segment(ConnectableElement, Group):
         tck, u = self.position_spline_params
         return np.mean(self.contour_to_position(u), axis=0)
 
+    def get_bounding_box( self, num_points=3 ):
+        positions = np.zeros( (num_points, 3) )
+        i = 0
+        for c in np.linspace(0,1,num_points):
+            positions[i] = (self.contour_to_position(c))
+            i += 1
+        min_ = np.array([np.min(positions[:,i]) for i in range(3)])
+        max_ = np.array([np.max(positions[:,i]) for i in range(3)])
+        return min_,max_
+
     def _get_location_positions(self):
         return [self.contour_to_nt_pos(l.address) for l in self.locations]
 
@@ -3241,17 +3251,27 @@ proc calcforces {} {
 }
 """)
 
-    def dimensions_from_structure( self, padding_factor=1.5, isotropic=False ):
-        positions = []
+    def get_bounding_box( self, num_points=3 ):
+        positions = np.zeros( (len(self.segments)*num_points, 3) )
+        i = 0
         for s in self.segments:
-            positions.append(s.contour_to_position(0))
-            positions.append(s.contour_to_position(0.5))
-            positions.append(s.contour_to_position(1))
-        positions = np.array(positions)
-        dx,dy,dz = [(np.max(positions[:,i])-np.min(positions[:,i])+30)*padding_factor for i in range(3)]
+            for c in np.linspace(0,1,num_points):
+                positions[i] = (s.contour_to_position(c))
+                i += 1
+        min_ = np.array([np.min(positions[:,i]) for i in range(3)])
+        max_ = np.array([np.max(positions[:,i]) for i in range(3)])
+        return min_,max_
+
+    def get_bounding_box_center( self, num_points=3 ):
+        min_,max_ = self.get_bounding_box(num_points)
+        return 0.5*(max_+min_)
+
+    def dimensions_from_structure( self, padding_factor=1.5, isotropic=False ):
+        min_,max_ = self.get_bounding_box()
+        dx,dy,dz = (max_-min_+30)*padding_factor
         if isotropic:
             dx = dy = dz = max((dx,dy,dz))
-        return [dx,dy,dz]
+        return np.array([dx,dy,dz])
 
     def add_grid_potential(self, grid_file, scale=1, per_nucleotide=True):
         grid_file = Path(grid_file)
-- 
GitLab