IMP  2.0.0
The Integrative Modeling Platform
sampling.py
1 
2 import IMP
3 import IMP.algebra as alg
4 import IMP.em2d.solutions_io as solutions_io
5 import IMP.em2d.imp_general.io as io
6 import IMP.em2d.Database as Database
7 import itertools
8 import logging
9 import os
10 
11 log = logging.getLogger("sampling")
12 
13 #########################################
14 """
15 
16  Functions required for sampling in a 2D/3D grid
17 
18 """
19 #########################################
20 
21 
22 def create_sampling_grid_3d(diameter, n_axis_points):
23  """
24  Creates a grid of positions (Vector3Ds), centered at 0.
25  @param diameter The shape of the grid is a sphere with this diameter.
26  @param n_axis_points Number of points used alogn an axis for the grid.
27  The axis X Y and Z will contain n_axis_points, equispaced.
28  The other regions of space will contain only the points allowed by
29  the size of the spere.
30  """
31  radius = diameter/2.0
32  step = diameter/n_axis_points
33  points = [-radius + step*n for n in range(n_axis_points +1 )]
34  points = filter(lambda x: abs(x) > 1, points) # filter points closer than 1A
35  log.info("Points along the axis %s", points)
36  positions = []
37  for x, y, z in itertools.product(points, points, points):
38  d = (x**2. + y**2. + z**2.)**.5
39  if(d < (1.01 *radius)): # allow 1% margin. Avoids approximation problems
40  positions.append( alg.Vector3D(x,y,z))
41  return positions
42 
43 
44 def create_sampling_grid_2d(diameter, n_axis_points):
45  """
46  Creates a grid of positions (Vector3Ds), centered at 0.
47  The shape of the grid is a circle with diameter given by the parameter.
48  n_axis_points is the number of points along an axis: The axis X Y
49  will contain n_axis_points, equispaced. The other regions of space will
50  contain only the points allowed by the size of the circle.
51  """
52  if(diameter < 0 ):
53  raise ValueError("create_sampling_grid_2d: Negative diameter.")
54  if(n_axis_points < 1 ):
55  raise ValueError("create_sampling_grid_2d: Less than one axis point.")
56  radius = diameter/2.0
57  step = diameter/n_axis_points
58  points = [-radius + step*n for n in range(n_axis_points +1 )]
59  log.info("Points along the axis %s", points)
60  positions = []
61  for x, y in itertools.product(points, points):
62  d = (x**2. + y**2.)**.5
63  if(d < (1.01 *radius)): # allow 1% margin. Avoids approximation problems
64  positions.append( alg.Vector3D(x,y, 0.0))
65  return positions
66 
67 
68 def get_orientations_nearby(rotation, n, f):
69  """
70  Rotations nearby a given one. They are got intepolating with
71  the rotations of the uniform coverage. The parameter f
72  (0 <= f <= 1) must be close to 0 to get orientations
73  that are close to the given orientation
74  - Values from 0.1 (tight cluster) to 0.4 (loose)
75  seem to be ok
76  n - number of rotations requested
77  """
78  log.debug("Computing nearby rotations around %s", rotation)
79  unif = alg.get_uniform_cover_rotations_3d(n)
80  id_rot = alg.get_identity_rotation_3d()
81  oris = []
82  for rot in unif:
83  r = alg.get_interpolated(rot, id_rot, f)
84  r = alg.compose(r, rotation )
85  oris.append(r)
86  return oris
87 
88 
89 
90 class SamplingSchema:
91 
92  def __init__(self, n_components, fixed, anchored):
93  """
94  Build a set of orientations and positions for sampling with DOMINO
95  @param n_components Number of components in the assembly
96  @param fixed List of True/False values. If fixed == True
97  for a component its position is not changed.
98  @param anchor A list of True/False values. If anchor = True,
99  the component is set at (0,0,0) without rotating it.
100  """
101  if(n_components < 1 ):
102  raise ValueError("SamplingSchema: Requesting less than 1 components.")
103 
104  self.anchored = anchored
105  self.fixed = fixed
106  self.n_components = n_components
107 
108  def get_sampling_transformations(self, i):
109  return self.transformations[i]
110 
111  def read_from_database(self, fn_database, fields=["reference_frames"],
112  max_number=False, orderby=False ):
113  """
114  Read orientations and positions from a database file.
115  self.anchored and self.fixed overwrite
116  the positions and orientations read from the database
117  """
118  if not os.path.exists(fn_database):
119  raise IOError("read_from_database: Database file not found. " \
120  "Are you perhaps trying to run the DOMINO optimization without " \
121  "having the database yet?")
122 
123 
124  db = solutions_io.ResultsDB()
125  db.connect(fn_database)
126  data = db.get_solutions(fields, max_number, orderby)
127  db.close()
128  self.transformations = [ [] for T in range(self.n_components)]
129  # Each record contains a reference frame for each of the
130  # components. But here the states considered make sense as columns.
131  # Each rigidbody is considered to have all the states in a column.
132  for d in data:
133  texts = d[0].split("/")
134  for i, t in zip(range(self.n_components), texts):
135  T = io.TextToTransformation3D(t).get_transformation()
136  self.transformations[i].append(T)
137 
138 
139 
140  # Set the anchored components
141  for i in range(self.n_components):
142  if self.anchored[i]:
143  origin = alg.Transformation3D(alg.get_identity_rotation_3d(),
144  alg.Vector3D(0.0, 0.0, 0.))
145  self.transformations[i] = [origin]
146 
147  # If fixed, only the first state is kept
148  for i in range(self.n_components):
149  if self.fixed[i]:
150  if len(self.transformations[i]) == 0:
151  raise ValueError("There are positions to keep fixed")
152  self.transformations[i] = [self.transformations[i][0]]