IMP logo
IMP Reference Guide  2.8.0
The Integrative Modeling Platform
tools.py
1 #!/usr/bin/env python
2 
3 """@namespace IMP.pmi.tools
4  Miscellaneous utilities.
5 """
6 
7 from __future__ import print_function
8 import IMP
9 import IMP.algebra
10 import IMP.isd
11 import IMP.pmi
12 import IMP.pmi.topology
13 import collections
14 import itertools
15 from math import log,pi,sqrt,exp
16 import sys,os
17 import random
18 import ast
19 import time
20 import RMF
21 import IMP.rmf
22 from collections import defaultdict
23 try:
24  from collections import OrderedDict
25 except ImportError:
26  from IMP.pmi._compat_collections import OrderedDict
27 
28 def _get_restraint_set_key():
29  if not hasattr(_get_restraint_set_key, 'pmi_rs_key'):
30  _get_restraint_set_key.pmi_rs_key = IMP.ModelKey("PMI restraints")
31  return _get_restraint_set_key.pmi_rs_key
32 
33 def _add_restraint_set(model, mk):
34  rs = IMP.RestraintSet(model, "All PMI restraints")
35  model.add_data(mk, rs)
36  return rs
37 
38 def add_restraint_to_model(model, restraint):
39  """Add a PMI restraint to the model.
40  Since Model.add_restraint() no longer exists (in modern IMP restraints
41  should be added to a ScoringFunction instead) store them instead in
42  a RestraintSet, and keep a reference to it in the Model."""
43  mk = _get_restraint_set_key()
44  if model.get_has_data(mk):
45  rs = IMP.RestraintSet.get_from(model.get_data(mk))
46  else:
47  rs = _add_restraint_set(model, mk)
48  rs.add_restraint(restraint)
49 
50 def get_restraint_set(model):
51  """Get a RestraintSet containing all PMI restraints added to the model"""
52  mk = _get_restraint_set_key()
53  if not model.get_has_data(mk):
54  print("WARNING: no restraints added to model yet")
55  _add_restraint_set(model, mk)
56  return IMP.RestraintSet.get_from(model.get_data(mk))
57 
58 class Stopwatch(object):
59  """Collect timing information.
60  Add an instance of this class to outputobjects to get timing information
61  in a stat file."""
62 
63  def __init__(self, isdelta=True):
64  """Constructor.
65  @param isdelta if True (the default) then report the time since the
66  last use of this class; if False, report cumulative time."""
67  self.starttime = time.clock()
68  self.label = "None"
69  self.isdelta = isdelta
70 
71  def set_label(self, labelstr):
72  self.label = labelstr
73 
74  def get_output(self):
75  output = {}
76  if self.isdelta:
77  newtime = time.clock()
78  output["Stopwatch_" + self.label + "_delta_seconds"] \
79  = str(newtime - self.starttime)
80  self.starttime = newtime
81  else:
82  output["Stopwatch_" + self.label + "_elapsed_seconds"] \
83  = str(time.clock() - self.starttime)
84  return output
85 
86 
87 class SetupNuisance(object):
88 
89  def __init__(self, m, initialvalue, minvalue, maxvalue, isoptimized=True):
90 
91  nuisance = IMP.isd.Scale.setup_particle(IMP.Particle(m), initialvalue)
92  if minvalue:
93  nuisance.set_lower(minvalue)
94  if maxvalue:
95  nuisance.set_upper(maxvalue)
96 
97  # m.add_score_state(IMP.core.SingletonConstraint(IMP.isd.NuisanceRangeModifier(),None,nuisance))
98  nuisance.set_is_optimized(nuisance.get_nuisance_key(), isoptimized)
99  self.nuisance = nuisance
100 
101  def get_particle(self):
102  return self.nuisance
103 
104 
105 class SetupWeight(object):
106 
107  def __init__(self, m, isoptimized=True):
108  pw = IMP.Particle(m)
109  self.weight = IMP.isd.Weight.setup_particle(pw)
110  self.weight.set_weights_are_optimized(True)
111 
112  def get_particle(self):
113  return self.weight
114 
115 
116 class SetupSurface(object):
117 
118  def __init__(self, m, center, normal, isoptimized=True):
119  p = IMP.Particle(m)
120  self.surface = IMP.core.Surface.setup_particle(p, center, normal)
121  self.surface.set_coordinates_are_optimized(isoptimized)
122  self.surface.set_normal_is_optimized(isoptimized)
123 
124  def get_particle(self):
125  return self.surface
126 
127 
129  "If you use this class please let the PMI developers know.")
130 class ParticleToSampleFilter(object):
131  def __init__(self, sampled_objects):
132  self.sampled_objects=sampled_objects
133  self.filters=[]
134 
135  def add_filter(self,filter_string):
136  self.filters.append(filter_string)
137 
138  def get_particles_to_sample(self):
139  particles_to_sample={}
140  for so in self.sampled_objects:
141  ps_dict=so.get_particles_to_sample()
142  for key in ps_dict:
143  for f in self.filters:
144  if f in key:
145  if key not in particles_to_sample:
146  particles_to_sample[key]=ps_dict[key]
147  else:
148  particles_to_sample[key]+=ps_dict[key]
149  return particles_to_sample
150 
151 class ParticleToSampleList(object):
152 
153  def __init__(self, label="None"):
154 
155  self.dictionary_particle_type = {}
156  self.dictionary_particle_transformation = {}
157  self.dictionary_particle_name = {}
158  self.label = label
159 
160  def add_particle(
161  self,
162  particle,
163  particle_type,
164  particle_transformation,
165  name):
166  if not particle_type in ["Rigid_Bodies", "Floppy_Bodies", "Nuisances", "X_coord", "Weights", "Surfaces"]:
167  raise TypeError("not the right particle type")
168  else:
169  self.dictionary_particle_type[particle] = particle_type
170  if particle_type == "Rigid_Bodies":
171  if type(particle_transformation) == tuple and len(particle_transformation) == 2 and type(particle_transformation[0]) == float and type(particle_transformation[1]) == float:
172  self.dictionary_particle_transformation[
173  particle] = particle_transformation
174  self.dictionary_particle_name[particle] = name
175  else:
176  raise TypeError("ParticleToSampleList: not the right transformation format for Rigid_Bodies, should be a tuple of floats")
177  elif particle_type == "Surfaces":
178  if type(particle_transformation) == tuple and len(particle_transformation) == 3 and all(isinstance(x, float) for x in particle_transformation):
179  self.dictionary_particle_transformation[
180  particle] = particle_transformation
181  self.dictionary_particle_name[particle] = name
182  else:
183  raise TypeError("ParticleToSampleList: not the right transformation format for Surfaces, should be a tuple of floats")
184  else:
185  if type(particle_transformation) == float:
186  self.dictionary_particle_transformation[
187  particle] = particle_transformation
188  self.dictionary_particle_name[particle] = name
189  else:
190  raise TypeError("ParticleToSampleList: not the right transformation format, should be a float")
191 
192  def get_particles_to_sample(self):
193  ps = {}
194  for particle in self.dictionary_particle_type:
195  key = self.dictionary_particle_type[
196  particle] + "ParticleToSampleList_" + self.dictionary_particle_name[particle] + "_" + self.label
197  value = (
198  [particle],
199  self.dictionary_particle_transformation[particle])
200  ps[key] = value
201  return ps
202 
203 
204 class Variance(object):
205 
206  def __init__(self, model, tau, niter, prot, th_profile, write_data=False):
207 
208  self.model = model
209  self.write_data = write_data
210  self.tau = tau
211  self.niter = niter
212  #! select particles from the model
213  particles = IMP.atom.get_by_type(prot, IMP.atom.ATOM_TYPE)
214  self.particles = particles
215  # store reference coordinates and theoretical profile
216  self.refpos = [IMP.core.XYZ(p).get_coordinates() for p in particles]
217  self.model_profile = th_profile
218 
219  def perturb_particles(self, perturb=True):
220  for i, p in enumerate(self.particles):
221  newpos = array(self.refpos[i])
222  if perturb:
223  newpos += random.normal(0, self.tau, 3)
224  newpos = IMP.algebra.Vector3D(newpos)
225  IMP.core.XYZ(p).set_coordinates(newpos)
226 
227  def get_profile(self):
228  model_profile = self.model_profile
229  p = model_profile.calculate_profile(self.particles, IMP.saxs.CA_ATOMS)
230  return array([model_profile.get_intensity(i) for i in
231  range(model_profile.size())])
232 
233  def init_variances(self):
234  # create placeholders
235  N = self.model_profile.size()
236  a = self.profiles[0][:]
237  self.m = matrix(a).T # Nx1
238  self.V = self.m * self.m.T
239  self.normm = linalg.norm(self.m)
240  self.normV = linalg.norm(self.V)
241 
242  def update_variances(self):
243  a = matrix(self.profiles[-1]) # 1xN
244  n = float(len(self.profiles))
245  self.m = a.T / n + (n - 1) / n * self.m
246  self.V = a.T * a + self.V
247  self.oldnormm = self.normm
248  self.oldnormV = self.normV
249  self.normm = linalg.norm(self.m)
250  self.normV = linalg.norm(self.V)
251  self.diffm = (self.oldnormm - self.normm) / self.oldnormm
252  self.diffV = (self.oldnormV - self.normV) / self.oldnormV
253 
254  def get_direct_stats(self, a):
255  nq = len(a[0])
256  nprof = len(a)
257  m = [0] * nq
258  for prof in a:
259  for q, I in enumerate(prof):
260  m[q] += I
261  m = array(m) / nprof
262  V = matrix(a)
263  V = V.T * V
264  Sigma = (matrix(a - m))
265  Sigma = Sigma.T * Sigma / (nprof - 1)
266  mi = matrix(diag(1. / m))
267  Sigmarel = mi.T * Sigma * mi
268  return m, V, Sigma, Sigmarel
269 
270  def store_data(self):
271  if not os.path.isdir('data'):
272  os.mkdir('data')
273  profiles = matrix(self.profiles)
274  self.directm, self.directV, self.Sigma, self.Sigmarel = \
275  self.get_direct_stats(array(profiles))
276  directV = self.directV
277  # print "V comparison",(linalg.norm(directV-self.V)/self.normV)
278  save('data/profiles', profiles)
279  # absolute profile differences
280  fl = open('data/profiles.dat', 'w')
281  for i, l in enumerate(array(profiles).T):
282  self.model_profile.get_q(i)
283  fl.write('%s ' % i)
284  for k in l:
285  fl.write('%s ' % (k - self.directm[i]))
286  fl.write('\n')
287  # relative profile differences
288  fl = open('data/profiles_rel.dat', 'w')
289  for i, l in enumerate(array(profiles).T):
290  self.model_profile.get_q(i)
291  fl.write('%s ' % i)
292  for k in l:
293  fl.write('%s ' % ((k - self.directm[i]) / self.directm[i]))
294  fl.write('\n')
295  save('data/m', self.directm)
296  save('data/V', self.directV)
297  Sigma = self.Sigma
298  save('data/Sigma', Sigma)
299  # Sigma matrix
300  fl = open('data/Sigma.dat', 'w')
301  model_profile = self.model_profile
302  for i in range(model_profile.size()):
303  qi = model_profile.get_q(i)
304  for j in range(model_profile.size()):
305  qj = model_profile.get_q(j)
306  vij = self.Sigma[i, j]
307  fl.write('%s %s %s\n' % (qi, qj, vij))
308  fl.write('\n')
309  # Sigma eigenvalues
310  fl = open('data/eigenvals', 'w')
311  for i in linalg.eigvalsh(Sigma):
312  fl.write('%s\n' % i)
313  Sigmarel = self.Sigmarel
314  save('data/Sigmarel', Sigmarel)
315  # Sigmarel matrix
316  fl = open('data/Sigmarel.dat', 'w')
317  model_profile = self.model_profile
318  for i in range(model_profile.size()):
319  qi = model_profile.get_q(i)
320  for j in range(model_profile.size()):
321  qj = model_profile.get_q(j)
322  vij = self.Sigmarel[i, j]
323  fl.write('%s %s %s\n' % (qi, qj, vij))
324  fl.write('\n')
325  # Sigma eigenvalues
326  fl = open('data/eigenvals_rel', 'w')
327  for i in linalg.eigvalsh(Sigmarel):
328  fl.write('%s\n' % i)
329  # mean profile
330  fl = open('data/mean.dat', 'w')
331  for i in range(len(self.directm)):
332  qi = self.model_profile.get_q(i)
333  fl.write('%s ' % qi)
334  fl.write('%s ' % self.directm[i])
335  fl.write('%s ' % sqrt(self.Sigma[i, i]))
336  fl.write('\n')
337 
338  def try_chol(self, jitter):
339  Sigma = self.Sigma
340  try:
341  linalg.cholesky(Sigma + matrix(eye(len(Sigma))) * jitter)
342  except linalg.LinAlgError:
343  print("Decomposition failed with jitter =", jitter)
344  return
345  print("Successful decomposition with jitter =", jitter)
346 
347  def run(self):
348  self.profiles = [self.get_profile()]
349  # self.init_variances()
350  for n in range(self.niter):
351  self.perturb_particles()
352  self.profiles.append(self.get_profile())
353  # self.update_variances()
354  #profiles = matrix(self.profiles)
355  # print n,self.diffm,self.diffV
356  # print
357  #
358  if self.write_data:
359  self.store_data()
360  # self.try_chol(0.)
361  # for i in logspace(-7,0,num=8):
362  # self.try_chol(i)
363 
364  def get_cov(self, relative=True):
365  if not relative:
366  return self.Sigma
367  else:
368  return self.Sigmarel
369 
370  #-------------------------------
371 
372 def get_random_cross_link_dataset(representation,
373  resolution=1.0,
374  number_of_cross_links=10,
375  ambiguity_probability=0.1,
376  confidence_score_range=[0,100],
377  avoid_same_particles=False):
378  '''Return a random cross-link dataset as a string.
379  Every line is a residue pair, together with UniqueIdentifier
380  and XL score.'''
381 
382  residue_pairs=get_random_residue_pairs(representation, resolution, number_of_cross_links, avoid_same_particles=avoid_same_particles)
383 
384  unique_identifier=0
385  cmin=float(min(confidence_score_range))
386  cmax=float(max(confidence_score_range))
387 
388  dataset="#\n"
389 
390  for (name1, r1, name2, r2) in residue_pairs:
391  if random.random() > ambiguity_probability:
392  unique_identifier+=1
393  score=random.random()*(cmax-cmin)+cmin
394  dataset+=str(name1)+" "+str(name2)+" "+str(r1)+" "+str(r2)+" "+str(score)+" "+str(unique_identifier)+"\n"
395 
396  return dataset
397 
398 
399  #-------------------------------
400 
401 def get_cross_link_data(directory, filename, dist, omega, sigma,
402  don=None, doff=None, prior=0, type_of_profile="gofr"):
403 
404  (distmin, distmax, ndist) = dist
405  (omegamin, omegamax, nomega) = omega
406  (sigmamin, sigmamax, nsigma) = sigma
407 
408  filen = IMP.isd.get_data_path("CrossLinkPMFs.dict")
409  with open(filen) as xlpot:
410  dictionary = ast.literal_eval(xlpot.readline())
411 
412  xpot = dictionary[directory][filename]["distance"]
413  pot = dictionary[directory][filename][type_of_profile]
414 
415  dist_grid = get_grid(distmin, distmax, ndist, False)
416  omega_grid = get_log_grid(omegamin, omegamax, nomega)
417  sigma_grid = get_log_grid(sigmamin, sigmamax, nsigma)
418 
419  if not don is None and not doff is None:
420  xlmsdata = IMP.isd.CrossLinkData(
421  dist_grid,
422  omega_grid,
423  sigma_grid,
424  xpot,
425  pot,
426  don,
427  doff,
428  prior)
429  else:
430  xlmsdata = IMP.isd.CrossLinkData(
431  dist_grid,
432  omega_grid,
433  sigma_grid,
434  xpot,
435  pot)
436  return xlmsdata
437 
438  #-------------------------------
439 
440 
441 def get_cross_link_data_from_length(length, xxx_todo_changeme3, xxx_todo_changeme4, xxx_todo_changeme5):
442  (distmin, distmax, ndist) = xxx_todo_changeme3
443  (omegamin, omegamax, nomega) = xxx_todo_changeme4
444  (sigmamin, sigmamax, nsigma) = xxx_todo_changeme5
445 
446  dist_grid = get_grid(distmin, distmax, ndist, False)
447  omega_grid = get_log_grid(omegamin, omegamax, nomega)
448  sigma_grid = get_log_grid(sigmamin, sigmamax, nsigma)
449 
450  xlmsdata = IMP.isd.CrossLinkData(dist_grid, omega_grid, sigma_grid, length)
451  return xlmsdata
452 
453 
454 def get_grid(gmin, gmax, ngrid, boundaries):
455  grid = []
456  dx = (gmax - gmin) / float(ngrid)
457  for i in range(0, ngrid + 1):
458  if(not boundaries and i == 0):
459  continue
460  if(not boundaries and i == ngrid):
461  continue
462  grid.append(gmin + float(i) * dx)
463  return grid
464 
465  #-------------------------------
466 
467 
468 def get_log_grid(gmin, gmax, ngrid):
469  grid = []
470  for i in range(0, ngrid + 1):
471  grid.append(gmin * exp(float(i) / ngrid * log(gmax / gmin)))
472  return grid
473 
474  #-------------------------------
475 
476 
478  '''
479  example '"{ID_Score}" > 28 AND "{Sample}" ==
480  "%10_1%" OR ":Sample}" == "%10_2%" OR ":Sample}"
481  == "%10_3%" OR ":Sample}" == "%8_1%" OR ":Sample}" == "%8_2%"'
482  '''
483 
484  import pyparsing as pp
485 
486  operator = pp.Regex(">=|<=|!=|>|<|==|in").setName("operator")
487  value = pp.QuotedString(
488  '"') | pp.Regex(
489  r"[+-]?\d+(:?\.\d*)?(:?[eE][+-]?\d+)?")
490  identifier = pp.Word(pp.alphas, pp.alphanums + "_")
491  comparison_term = identifier | value
492  condition = pp.Group(comparison_term + operator + comparison_term)
493 
494  expr = pp.operatorPrecedence(condition, [
495  ("OR", 2, pp.opAssoc.LEFT, ),
496  ("AND", 2, pp.opAssoc.LEFT, ),
497  ])
498 
499  parsedstring = str(expr.parseString(inputstring)) \
500  .replace("[", "(") \
501  .replace("]", ")") \
502  .replace(",", " ") \
503  .replace("'", " ") \
504  .replace("%", "'") \
505  .replace("{", "float(entry['") \
506  .replace("}", "'])") \
507  .replace(":", "str(entry['") \
508  .replace("}", "'])") \
509  .replace("AND", "and") \
510  .replace("OR", "or")
511  return parsedstring
512 
513 
514 def open_file_or_inline_text(filename):
515  try:
516  fl = open(filename, "r")
517  except IOError:
518  fl = filename.split("\n")
519  return fl
520 
521 
522 def get_drmsd(prot0, prot1):
523  drmsd = 0.
524  npairs = 0.
525  for i in range(0, len(prot0) - 1):
526  for j in range(i + 1, len(prot0)):
527  dist0 = IMP.core.get_distance(prot0[i], prot0[j])
528  dist1 = IMP.core.get_distance(prot1[i], prot1[j])
529  drmsd += (dist0 - dist1) ** 2
530  npairs += 1.
531  return math.sqrt(drmsd / npairs)
532 
533  #-------------------------------
534 
535 
536 def get_ids_from_fasta_file(fastafile):
537  ids = []
538  with open(fastafile) as ff:
539  for l in ff:
540  if l[0] == ">":
541  ids.append(l[1:-1])
542  return ids
543 
544 
545 def get_closest_residue_position(hier, resindex, terminus="N"):
546  '''
547  this function works with plain hierarchies, as read from the pdb,
548  no multi-scale hierarchies
549  '''
550  p = []
551  niter = 0
552  while len(p) == 0:
553  niter += 1
554  sel = IMP.atom.Selection(hier, residue_index=resindex,
555  atom_type=IMP.atom.AT_CA)
556 
557  if terminus == "N":
558  resindex += 1
559  if terminus == "C":
560  resindex -= 1
561 
562  if niter >= 10000:
563  print("get_closest_residue_position: exiting while loop without result")
564  break
565  p = sel.get_selected_particles()
566 
567  if len(p) == 1:
568  return IMP.core.XYZ(p[0]).get_coordinates()
569  elif len(p) == 0:
570  print("get_closest_residue_position: got NO residues for hierarchy %s and residue %i" % (hier, resindex))
571  raise Exception("get_closest_residue_position: got NO residues for hierarchy %s and residue %i" % (
572  hier, resindex))
573  else:
574  raise ValueError("got multiple residues for hierarchy %s and residue %i; the list of particles is %s" % (hier, resindex, str([pp.get_name() for pp in p])))
575 
576 def get_terminal_residue(representation, hier, terminus="C", resolution=1):
577  '''
578  Get the particle of the terminal residue at the GIVEN resolution
579  (NOTE: not the closest resolution!).
580  To get the terminal residue at the closest resolution use:
581  particles=IMP.pmi.tools.select_by_tuple(representation,molecule_name)
582  particles[0] and particles[-1] will be the first and last particles
583  corresponding to the two termini.
584  It is needed for instance to determine the last residue of a pdb.
585  @param hier hierarchy containing the terminal residue
586  @param terminus either 'N' or 'C'
587  @param resolution resolution to use.
588  '''
589  termresidue = None
590  termparticle = None
591 
592  ps=select(representation,
593  resolution=resolution,
594  hierarchies=[hier])
595 
596  for p in ps:
597  if IMP.pmi.Resolution(p).get_resolution() == resolution:
599  if terminus == "C":
600  if termresidue is None:
601  termresidue = max(residues)
602  termparticle = p
603  elif max(residues) >= termresidue:
604  termresidue = max(residues)
605  termparticle = p
606  elif terminus == "N":
607  if termresidue is None:
608  termresidue = min(residues)
609  termparticle = p
610  elif min(residues) <= termresidue:
611  termresidue = min(residues)
612  termparticle = p
613  else:
614  raise ValueError("terminus argument should be either N or C")
615  return termparticle
616 
617 def get_terminal_residue_position(representation, hier, terminus="C",
618  resolution=1):
619  """Get XYZ coordinates of the terminal residue at the GIVEN resolution"""
620  p = get_terminal_residue(representation, hier, terminus, resolution)
621  return IMP.core.XYZ(p).get_coordinates()
622 
623 def get_residue_gaps_in_hierarchy(hierarchy, start, end):
624  '''
625  Return the residue index gaps and contiguous segments in the hierarchy.
626 
627  @param hierarchy hierarchy to examine
628  @param start first residue index
629  @param end last residue index
630 
631  @return A list of lists of the form
632  [[1,100,"cont"],[101,120,"gap"],[121,200,"cont"]]
633  '''
634  gaps = []
635  for n, rindex in enumerate(range(start, end + 1)):
636  sel = IMP.atom.Selection(hierarchy, residue_index=rindex,
637  atom_type=IMP.atom.AT_CA)
638 
639  if len(sel.get_selected_particles()) == 0:
640  if n == 0:
641  # set the initial condition
642  rindexgap = start
643  rindexcont = start - 1
644  if rindexgap == rindex - 1:
645  # residue is contiguous with the previously discovered gap
646  gaps[-1][1] += 1
647  else:
648  # residue is not contiguous with the previously discovered gap
649  # hence create a new gap tuple
650  gaps.append([rindex, rindex, "gap"])
651  # update the index of the last residue gap
652  rindexgap = rindex
653  else:
654  if n == 0:
655  # set the initial condition
656  rindexgap = start - 1
657  rindexcont = start
658  if rindexcont == rindex - 1:
659  # residue is contiguous with the previously discovered
660  # continuous part
661  gaps[-1][1] += 1
662  else:
663  # residue is not contiguous with the previously discovered continuous part
664  # hence create a new cont tuple
665  gaps.append([rindex, rindex, "cont"])
666  # update the index of the last residue gap
667  rindexcont = rindex
668  return gaps
669 
670 
671 class map(object):
672 
673  def __init__(self):
674  self.map = {}
675 
676  def set_map_element(self, xvalue, yvalue):
677  self.map[xvalue] = yvalue
678 
679  def get_map_element(self, invalue):
680  if type(invalue) == float:
681  n = 0
682  mindist = 1
683  for x in self.map:
684  dist = (invalue - x) * (invalue - x)
685 
686  if n == 0:
687  mindist = dist
688  minx = x
689  if dist < mindist:
690  mindist = dist
691  minx = x
692  n += 1
693  return self.map[minx]
694  elif type(invalue) == str:
695  return self.map[invalue]
696  else:
697  raise TypeError("wrong type for map")
698 
699 
700 def select(representation,
701  resolution=None,
702  hierarchies=None,
703  selection_arguments=None,
704  name=None,
705  name_is_ambiguous=False,
706  first_residue=None,
707  last_residue=None,
708  residue=None,
709  representation_type=None):
710  '''
711  this function uses representation=SimplifiedModel
712  it returns the corresponding selected particles
713  representation_type="Beads", "Res:X", "Densities", "Representation", "Molecule"
714  '''
715 
716  if resolution is None:
717  allparticles = IMP.atom.get_leaves(representation.prot)
718  resolution_particles = None
719  hierarchies_particles = None
720  names_particles = None
721  residue_range_particles = None
722  residue_particles = None
723  representation_type_particles = None
724 
725  if not resolution is None:
726  resolution_particles = []
727  hs = representation.get_hierarchies_at_given_resolution(resolution)
728  for h in hs:
729  resolution_particles += IMP.atom.get_leaves(h)
730 
731  if not hierarchies is None:
732  hierarchies_particles = []
733  for h in hierarchies:
734  hierarchies_particles += IMP.atom.get_leaves(h)
735 
736  if not name is None:
737  names_particles = []
738  if name_is_ambiguous:
739  for namekey in representation.hier_dict:
740  if name in namekey:
741  names_particles += IMP.atom.get_leaves(
742  representation.hier_dict[namekey])
743  elif name in representation.hier_dict:
744  names_particles += IMP.atom.get_leaves(representation.hier_dict[name])
745  else:
746  print("select: component %s is not there" % name)
747 
748  if not first_residue is None and not last_residue is None:
749  sel = IMP.atom.Selection(representation.prot,
750  residue_indexes=range(first_residue, last_residue + 1))
751  residue_range_particles = [IMP.atom.Hierarchy(p)
752  for p in sel.get_selected_particles()]
753 
754  if not residue is None:
755  sel = IMP.atom.Selection(representation.prot, residue_index=residue)
756  residue_particles = [IMP.atom.Hierarchy(p)
757  for p in sel.get_selected_particles()]
758 
759  if not representation_type is None:
760  representation_type_particles = []
761  if representation_type == "Molecule":
762  for name in representation.hier_representation:
763  for repr_type in representation.hier_representation[name]:
764  if repr_type == "Beads" or "Res:" in repr_type:
765  h = representation.hier_representation[name][repr_type]
766  representation_type_particles += IMP.atom.get_leaves(h)
767 
768  elif representation_type == "PDB":
769  for name in representation.hier_representation:
770  for repr_type in representation.hier_representation[name]:
771  if repr_type == "Res:" in repr_type:
772  h = representation.hier_representation[name][repr_type]
773  representation_type_particles += IMP.atom.get_leaves(h)
774 
775  else:
776  for name in representation.hier_representation:
777  h = representation.hier_representation[
778  name][
779  representation_type]
780  representation_type_particles += IMP.atom.get_leaves(h)
781 
782  selections = [hierarchies_particles, names_particles,
783  residue_range_particles, residue_particles, representation_type_particles]
784 
785  if resolution is None:
786  selected_particles = set(allparticles)
787  else:
788  selected_particles = set(resolution_particles)
789 
790  for s in selections:
791  if not s is None:
792  selected_particles = (set(s) & selected_particles)
793 
794  return list(selected_particles)
795 
796 
797 def select_by_tuple(
798  representation,
799  tupleselection,
800  resolution=None,
801  name_is_ambiguous=False):
802  if isinstance(tupleselection, tuple) and len(tupleselection) == 3:
803  particles = IMP.pmi.tools.select(representation, resolution=resolution,
804  name=tupleselection[2],
805  first_residue=tupleselection[0],
806  last_residue=tupleselection[1],
807  name_is_ambiguous=name_is_ambiguous)
808  elif isinstance(tupleselection, str):
809  particles = IMP.pmi.tools.select(representation, resolution=resolution,
810  name=tupleselection,
811  name_is_ambiguous=name_is_ambiguous)
812  else:
813  raise ValueError('you passed something bad to select_by_tuple()')
814  # now order the result by residue number
815  particles = IMP.pmi.tools.sort_by_residues(particles)
816 
817  return particles
818 
819 def select_by_tuple_2(hier,tuple_selection,resolution):
820  """New tuple format: molname OR (start,stop,molname,copynum,statenum)
821  Copy and state are optional. Can also use 'None' for them which will get all.
822  You can also pass -1 for stop which will go to the end.
823  Returns the particles
824  """
825  kwds = {} # going to accumulate keywords
826  kwds['resolution'] = resolution
827  if type(tuple_selection) is str:
828  kwds['molecule'] = tuple_selection
829  elif type(tuple_selection) is tuple:
830  rbegin = tuple_selection[0]
831  rend = tuple_selection[1]
832  kwds['molecule'] = tuple_selection[2]
833  try:
834  copynum = tuple_selection[3]
835  if copynum is not None:
836  kwds['copy_index'] = copynum
837  except:
838  pass
839  try:
840  statenum = tuple_selection[4]
841  if statenum is not None:
842  kwds['state_index'] = statenum
843  except:
844  pass
845  if rend==-1:
846  if rbegin>1:
847  s = IMP.atom.Selection(hier,**kwds)
848  s -= IMP.atom.Selection(hier,
849  residue_indexes=range(1,rbegin),
850  **kwds)
851  return s.get_selected_particles()
852  else:
853  kwds['residue_indexes'] = range(rbegin,rend+1)
854  s = IMP.atom.Selection(hier,**kwds)
855  return s.get_selected_particles()
856 
857 
858 
859 def get_db_from_csv(csvfilename):
860  import csv
861  outputlist = []
862  with open(csvfilename) as fh:
863  csvr = csv.DictReader(fh)
864  for l in csvr:
865  outputlist.append(l)
866  return outputlist
867 
868 
869 class HierarchyDatabase(object):
870  """Store the representations for a system."""
871 
872  def __init__(self):
873  self.db = {}
874  # this dictionary map a particle to its root hierarchy
875  self.root_hierarchy_dict = {}
876  self.preroot_fragment_hierarchy_dict = {}
877  self.particle_to_name = {}
878  self.model = None
879 
880  def add_name(self, name):
881  if name not in self.db:
882  self.db[name] = {}
883 
884  def add_residue_number(self, name, resn):
885  resn = int(resn)
886  self.add_name(name)
887  if resn not in self.db[name]:
888  self.db[name][resn] = {}
889 
890  def add_resolution(self, name, resn, resolution):
891  resn = int(resn)
892  resolution = float(resolution)
893  self.add_name(name)
894  self.add_residue_number(name, resn)
895  if resolution not in self.db[name][resn]:
896  self.db[name][resn][resolution] = []
897 
898  def add_particles(self, name, resn, resolution, particles):
899  resn = int(resn)
900  resolution = float(resolution)
901  self.add_name(name)
902  self.add_residue_number(name, resn)
903  self.add_resolution(name, resn, resolution)
904  self.db[name][resn][resolution] += particles
905  for p in particles:
906  (rh, prf) = self.get_root_hierarchy(p)
907  self.root_hierarchy_dict[p] = rh
908  self.preroot_fragment_hierarchy_dict[p] = prf
909  self.particle_to_name[p] = name
910  if self.model is None:
911  self.model = particles[0].get_model()
912 
913  def get_model(self):
914  return self.model
915 
916  def get_names(self):
917  names = list(self.db.keys())
918  names.sort()
919  return names
920 
921  def get_particles(self, name, resn, resolution):
922  resn = int(resn)
923  resolution = float(resolution)
924  return self.db[name][resn][resolution]
925 
926  def get_particles_at_closest_resolution(self, name, resn, resolution):
927  resn = int(resn)
928  resolution = float(resolution)
929  closestres = min(self.get_residue_resolutions(name, resn),
930  key=lambda x: abs(float(x) - float(resolution)))
931  return self.get_particles(name, resn, closestres)
932 
933  def get_residue_resolutions(self, name, resn):
934  resn = int(resn)
935  resolutions = list(self.db[name][resn].keys())
936  resolutions.sort()
937  return resolutions
938 
939  def get_molecule_resolutions(self, name):
940  resolutions = set()
941  for resn in self.db[name]:
942  resolutions.update(list(self.db[name][resn].keys()))
943  resolutions.sort()
944  return resolutions
945 
946  def get_residue_numbers(self, name):
947  residue_numbers = list(self.db[name].keys())
948  residue_numbers.sort()
949  return residue_numbers
950 
951  def get_particles_by_resolution(self, name, resolution):
952  resolution = float(resolution)
953  particles = []
954  for resn in self.get_residue_numbers(name):
955  result = self.get_particles_at_closest_resolution(
956  name,
957  resn,
958  resolution)
959  pstemp = [p for p in result if p not in particles]
960  particles += pstemp
961  return particles
962 
963  def get_all_particles_by_resolution(self, resolution):
964  resolution = float(resolution)
965  particles = []
966  for name in self.get_names():
967  particles += self.get_particles_by_resolution(name, resolution)
968  return particles
969 
970  def get_root_hierarchy(self, particle):
971  prerootfragment = particle
972  while IMP.atom.Atom.get_is_setup(particle) or \
973  IMP.atom.Residue.get_is_setup(particle) or \
975  if IMP.atom.Atom.get_is_setup(particle):
976  p = IMP.atom.Atom(particle).get_parent()
977  elif IMP.atom.Residue.get_is_setup(particle):
978  p = IMP.atom.Residue(particle).get_parent()
979  elif IMP.atom.Fragment.get_is_setup(particle):
980  p = IMP.atom.Fragment(particle).get_parent()
981  prerootfragment = particle
982  particle = p
983  return (
984  (IMP.atom.Hierarchy(particle), IMP.atom.Hierarchy(prerootfragment))
985  )
986 
987  def get_all_root_hierarchies_by_resolution(self, resolution):
988  hierarchies = []
989  resolution = float(resolution)
990  particles = self.get_all_particles_by_resolution(resolution)
991  for p in particles:
992  rh = self.root_hierarchy_dict[p]
993  if rh not in hierarchies:
994  hierarchies.append(IMP.atom.Hierarchy(rh))
995  return hierarchies
996 
997  def get_preroot_fragments_by_resolution(self, name, resolution):
998  fragments = []
999  resolution = float(resolution)
1000  particles = self.get_particles_by_resolution(name, resolution)
1001  for p in particles:
1002  fr = self.preroot_fragment_hierarchy_dict[p]
1003  if fr not in fragments:
1004  fragments.append(fr)
1005  return fragments
1006 
1007  def show(self, name):
1008  print(name)
1009  for resn in self.get_residue_numbers(name):
1010  print(resn)
1011  for resolution in self.get_residue_resolutions(name, resn):
1012  print("----", resolution)
1013  for p in self.get_particles(name, resn, resolution):
1014  print("--------", p.get_name())
1015 
1016 
1017 def get_prot_name_from_particle(p, list_of_names):
1018  '''Return the component name provided a particle and a list of names'''
1019  root = p
1020  protname = root.get_name()
1021  is_a_bead = False
1022  while not protname in list_of_names:
1023  root0 = root.get_parent()
1024  if root0 == IMP.atom.Hierarchy():
1025  return (None, None)
1026  protname = root0.get_name()
1027 
1028  # check if that is a bead
1029  # this piece of code might be dangerous if
1030  # the hierarchy was called Bead :)
1031  if "Beads" in protname:
1032  is_a_bead = True
1033  root = root0
1034  return (protname, is_a_bead)
1035 
1036 
1038  '''
1039  Retrieve the residue indexes for the given particle.
1040 
1041  The particle must be an instance of Fragment,Residue or Atom
1042  or else returns an empty list
1043  '''
1044  resind = []
1046  resind = IMP.atom.Fragment(hier).get_residue_indexes()
1047  elif IMP.atom.Residue.get_is_setup(hier):
1048  resind = [IMP.atom.Residue(hier).get_index()]
1049  elif IMP.atom.Atom.get_is_setup(hier):
1050  a = IMP.atom.Atom(hier)
1051  resind = [IMP.atom.Residue(a.get_parent()).get_index()]
1052  else:
1053  resind = []
1054  return resind
1055 
1056 
1057 def sort_by_residues(particles):
1058  particles_residues = [(p, IMP.pmi.tools.get_residue_indexes(p))
1059  for p in particles]
1060  sorted_particles_residues = sorted(
1061  particles_residues,
1062  key=lambda tup: tup[1])
1063  particles = [p[0] for p in sorted_particles_residues]
1064  return particles
1065 
1066 
1067 def get_residue_to_particle_map(particles):
1068  # this function returns a dictionary that map particles to residue indexes
1069  particles = sort_by_residues(particles)
1070  particles_residues = [(p, IMP.pmi.tools.get_residue_indexes(p))
1071  for p in particles]
1072  return dict(zip(particles_residues, particles))
1073 
1074 #
1075 # Parallel Computation
1076 #
1077 
1078 
1080  """Synchronize data over a parallel run"""
1081  from mpi4py import MPI
1082  comm = MPI.COMM_WORLD
1083  rank = comm.Get_rank()
1084  number_of_processes = comm.size
1085  comm.Barrier()
1086  if rank != 0:
1087  comm.send(data, dest=0, tag=11)
1088 
1089  elif rank == 0:
1090  for i in range(1, number_of_processes):
1091  data_tmp = comm.recv(source=i, tag=11)
1092  if type(data) == list:
1093  data += data_tmp
1094  elif type(data) == dict:
1095  data.update(data_tmp)
1096  else:
1097  raise TypeError("data not supported, use list or dictionaries")
1098 
1099  for i in range(1, number_of_processes):
1100  comm.send(data, dest=i, tag=11)
1101 
1102  if rank != 0:
1103  data = comm.recv(source=0, tag=11)
1104  return data
1105 
1107  """Synchronize data over a parallel run"""
1108  from mpi4py import MPI
1109  comm = MPI.COMM_WORLD
1110  rank = comm.Get_rank()
1111  number_of_processes = comm.size
1112  comm.Barrier()
1113  if rank != 0:
1114  comm.send(data, dest=0, tag=11)
1115  elif rank == 0:
1116  for i in range(1, number_of_processes):
1117  data_tmp = comm.recv(source=i, tag=11)
1118  for k in data:
1119  data[k]+=data_tmp[k]
1120 
1121  for i in range(1, number_of_processes):
1122  comm.send(data, dest=i, tag=11)
1123 
1124  if rank != 0:
1125  data = comm.recv(source=0, tag=11)
1126  return data
1127 
1128 
1129 #
1130 ### Lists and iterators
1131 #
1132 
1133 def sublist_iterator(l, lmin=None, lmax=None):
1134  '''
1135  Yield all sublists of length >= lmin and <= lmax
1136  '''
1137  if lmin is None:
1138  lmin = 0
1139  if lmax is None:
1140  lmax = len(l)
1141  n = len(l) + 1
1142  for i in range(n):
1143  for j in range(i + 1, n):
1144  if len(l[i:j]) <= lmax and len(l[i:j]) >= lmin:
1145  yield l[i:j]
1146 
1147 
1148 def flatten_list(l):
1149  return [item for sublist in l for item in sublist]
1150 
1151 
1152 def list_chunks_iterator(list, length):
1153  """ Yield successive length-sized chunks from a list.
1154  """
1155  for i in range(0, len(list), length):
1156  yield list[i:i + length]
1157 
1158 
1159 def chunk_list_into_segments(seq, num):
1160  seq = list(seq)
1161  avg = len(seq) / float(num)
1162  out = []
1163  last = 0.0
1164 
1165  while last < len(seq):
1166  out.append(seq[int(last):int(last + avg)])
1167  last += avg
1168 
1169  return out
1170 
1171 
1172 class Segments(object):
1173 
1174  ''' This class stores integers
1175  in ordered compact lists eg:
1176  [[1,2,3],[6,7,8]]
1177  the methods help splitting and merging the internal lists
1178  Example:
1179  s=Segments([1,2,3]) is [[1,2,3]]
1180  s.add(4) is [[1,2,3,4]] (add right)
1181  s.add(3) is [[1,2,3,4]] (item already existing)
1182  s.add(7) is [[1,2,3,4],[7]] (new list)
1183  s.add([8,9]) is [[1,2,3,4],[7,8,9]] (add item right)
1184  s.add([5,6]) is [[1,2,3,4,5,6,7,8,9]] (merge)
1185  s.remove(3) is [[1,2],[4,5,6,7,8,9]] (split)
1186  etc.
1187  '''
1188 
1189  def __init__(self,index):
1190  '''index can be a integer or a list of integers '''
1191  if type(index) is int:
1192  self.segs=[[index]]
1193  elif type(index) is list:
1194  self.segs=[[index[0]]]
1195  for i in index[1:]:
1196  self.add(i)
1197 
1198  def add(self,index):
1199  '''index can be a integer or a list of integers '''
1200  if type(index) is int:
1201  mergeleft=None
1202  mergeright=None
1203  for n,s in enumerate(self.segs):
1204  if index in s:
1205  return 0
1206  else:
1207  if s[0]-index==1:
1208  mergeleft=n
1209  if index-s[-1]==1:
1210  mergeright=n
1211  if mergeright is None and mergeleft is None:
1212  self.segs.append([index])
1213  if not mergeright is None and mergeleft is None:
1214  self.segs[mergeright].append(index)
1215  if not mergeleft is None and mergeright is None:
1216  self.segs[mergeleft]=[index]+self.segs[mergeleft]
1217  if not mergeleft is None and not mergeright is None:
1218  self.segs[mergeright]=self.segs[mergeright]+[index]+self.segs[mergeleft]
1219  del self.segs[mergeleft]
1220 
1221  for n in range(len(self.segs)):
1222  self.segs[n].sort()
1223 
1224  self.segs.sort(key=lambda tup: tup[0])
1225 
1226  elif type(index) is list:
1227  for i in index:
1228  self.add(i)
1229 
1230  def remove(self,index):
1231  '''index can be a integer'''
1232  for n,s in enumerate(self.segs):
1233  if index in s:
1234  if s[0]==index:
1235  self.segs[n]=s[1:]
1236  elif s[-1]==index:
1237  self.segs[n]=s[:-1]
1238  else:
1239  i=self.segs[n].index(index)
1240  self.segs[n]=s[:i]
1241  self.segs.append(s[i+1:])
1242  for n in range(len(self.segs)):
1243  self.segs[n].sort()
1244  if len(self.segs[n])==0:
1245  del self.segs[n]
1246  self.segs.sort(key=lambda tup: tup[0])
1247 
1248  def get_flatten(self):
1249  ''' Returns a flatten list '''
1250  return [item for sublist in self.segs for item in sublist]
1251 
1252 
1253 
1254 #
1255 # COORDINATE MANIPULATION
1256 #
1257 
1258 
1259 def translate_hierarchy(hierarchy, translation_vector):
1260  '''
1261  Apply a translation to a hierarchy along the input vector.
1262  '''
1263  rbs = set()
1264  xyzs = set()
1265  if type(translation_vector) == list:
1266  transformation = IMP.algebra.Transformation3D(
1267  IMP.algebra.Vector3D(translation_vector))
1268  else:
1269  transformation = IMP.algebra.Transformation3D(translation_vector)
1270  for p in IMP.atom.get_leaves(hierarchy):
1272  rbs.add(IMP.core.RigidBody(p))
1274  rb = IMP.core.RigidMember(p).get_rigid_body()
1275  rbs.add(rb)
1276  else:
1277  xyzs.add(p)
1278  for xyz in xyzs:
1279  IMP.core.transform(IMP.core.XYZ(xyz), transformation)
1280  for rb in rbs:
1281  IMP.core.transform(rb, transformation)
1282 
1283 
1284 def translate_hierarchies(hierarchies, translation_vector):
1285  for h in hierarchies:
1286  IMP.pmi.tools.translate_hierarchy(h, translation_vector)
1287 
1288 
1289 def translate_hierarchies_to_reference_frame(hierarchies):
1290  xc = 0
1291  yc = 0
1292  zc = 0
1293  nc = 0
1294  for h in hierarchies:
1295  for p in IMP.atom.get_leaves(h):
1296  coor = IMP.core.XYZ(p).get_coordinates()
1297  nc += 1
1298  xc += coor[0]
1299  yc += coor[1]
1300  zc += coor[2]
1301  xc = xc / nc
1302  yc = yc / nc
1303  zc = zc / nc
1304  IMP.pmi.tools.translate_hierarchies(hierarchies, (-xc, -yc, -zc))
1305 
1306 
1307 #
1308 # Tools to simulate data
1309 #
1310 
1311 def normal_density_function(expected_value, sigma, x):
1312  return (
1313  1 / math.sqrt(2 * math.pi) / sigma *
1314  math.exp(-(x - expected_value) ** 2 / 2 / sigma / sigma)
1315  )
1316 
1317 
1318 def log_normal_density_function(expected_value, sigma, x):
1319  return (
1320  1 / math.sqrt(2 * math.pi) / sigma / x *
1321  math.exp(-(math.log(x / expected_value) ** 2 / 2 / sigma / sigma))
1322  )
1323 
1324 
1325 def get_random_residue_pairs(representation, resolution,
1326  number,
1327  max_distance=None,
1328  avoid_same_particles=False,
1329  names=None):
1330 
1331  particles = []
1332  if names is None:
1333  names=list(representation.hier_dict.keys())
1334 
1335  for name in names:
1336  prot = representation.hier_dict[name]
1337  particles += select(representation,name=name,resolution=resolution)
1338  random_residue_pairs = []
1339  while len(random_residue_pairs)<=number:
1340  p1 = random.choice(particles)
1341  p2 = random.choice(particles)
1342  if max_distance is not None and \
1343  core.get_distance(core.XYZ(p1), core.XYZ(p2)) > max_distance:
1344  continue
1345  r1 = random.choice(IMP.pmi.tools.get_residue_indexes(p1))
1346  r2 = random.choice(IMP.pmi.tools.get_residue_indexes(p2))
1347  if r1==r2 and avoid_same_particles: continue
1348  name1 = representation.get_prot_name_from_particle(p1)
1349  name2 = representation.get_prot_name_from_particle(p2)
1350  random_residue_pairs.append((name1, r1, name2, r2))
1351 
1352  return random_residue_pairs
1353 
1354 
1355 def get_random_data_point(
1356  expected_value,
1357  ntrials,
1358  sensitivity,
1359  sigma,
1360  outlierprob,
1361  begin_end_nbins_tuple,
1362  log=False,
1363  loggrid=False):
1364  import random
1365 
1366  begin = begin_end_nbins_tuple[0]
1367  end = begin_end_nbins_tuple[1]
1368  nbins = begin_end_nbins_tuple[2]
1369 
1370  if not loggrid:
1371  fmod_grid = get_grid(begin, end, nbins, True)
1372  else:
1373  fmod_grid = get_log_grid(begin, end, nbins)
1374 
1375  norm = 0
1376  cumul = []
1377  cumul.append(0)
1378 
1379  a = []
1380  for i in range(0, ntrials):
1381  a.append([random.random(), True])
1382 
1383  if sigma != 0.0:
1384  for j in range(1, len(fmod_grid)):
1385  fj = fmod_grid[j]
1386  fjm1 = fmod_grid[j - 1]
1387  df = fj - fjm1
1388 
1389  if not log:
1390  pj = normal_density_function(expected_value, sigma, fj)
1391  pjm1 = normal_density_function(expected_value, sigma, fjm1)
1392  else:
1393  pj = log_normal_density_function(expected_value, sigma, fj)
1394  pjm1 = log_normal_density_function(expected_value, sigma, fjm1)
1395 
1396  norm += (pj + pjm1) / 2.0 * df
1397  cumul.append(norm)
1398  # print fj, pj
1399 
1400  random_points = []
1401 
1402  for i in range(len(cumul)):
1403  # print i,a, cumul[i], norm
1404  for aa in a:
1405  if (aa[0] <= cumul[i] / norm and aa[1]):
1406  random_points.append(
1407  int(fmod_grid[i] / sensitivity) * sensitivity)
1408  aa[1] = False
1409 
1410  else:
1411  random_points = [expected_value] * ntrials
1412 
1413  for i in range(len(random_points)):
1414  if random.random() < outlierprob:
1415  a = random.uniform(begin, end)
1416  random_points[i] = int(a / sensitivity) * sensitivity
1417  print(random_points)
1418  '''
1419  for i in range(ntrials):
1420  if random.random() > OUTLIERPROB_:
1421  r=truncnorm.rvs(0.0,1.0,expected_value,BETA_)
1422  if r>1.0: print r,expected_value,BETA_
1423  else:
1424  r=random.random()
1425  random_points.append(int(r/sensitivity)*sensitivity)
1426  '''
1427 
1428  rmean = 0.
1429  rmean2 = 0.
1430  for r in random_points:
1431  rmean += r
1432  rmean2 += r * r
1433 
1434  rmean /= float(ntrials)
1435  rmean2 /= float(ntrials)
1436  stddev = math.sqrt(max(rmean2 - rmean * rmean, 0.))
1437  return rmean, stddev
1438 
1439 def print_multicolumn(list_of_strings, ncolumns=2, truncate=40):
1440 
1441  l = list_of_strings
1442 
1443  cols = ncolumns
1444  # add empty entries after l
1445  for i in range(len(l) % cols):
1446  l.append(" ")
1447 
1448  split = [l[i:i + len(l) / cols] for i in range(0, len(l), len(l) / cols)]
1449  for row in zip(*split):
1450  print("".join(str.ljust(i, truncate) for i in row))
1451 
1452 class ColorChange(object):
1453  '''Change color code to hexadecimal to rgb'''
1454  def __init__(self):
1455  self._NUMERALS = '0123456789abcdefABCDEF'
1456  self._HEXDEC = dict((v, int(v, 16)) for v in (x+y for x in self._NUMERALS for y in self._NUMERALS))
1457  self.LOWERCASE, self.UPPERCASE = 'x', 'X'
1458 
1459  def rgb(self,triplet):
1460  return float(self._HEXDEC[triplet[0:2]]), float(self._HEXDEC[triplet[2:4]]), float(self._HEXDEC[triplet[4:6]])
1461 
1462  def triplet(self,rgb, lettercase=None):
1463  if lettercase is None: lettercase=self.LOWERCASE
1464  return format(rgb[0]<<16 | rgb[1]<<8 | rgb[2], '06'+lettercase)
1465 
1466 # -------------- Collections --------------- #
1467 
1468 class OrderedSet(collections.MutableSet):
1469 
1470  def __init__(self, iterable=None):
1471  self.end = end = []
1472  end += [None, end, end] # sentinel node for doubly linked list
1473  self.map = {} # key --> [key, prev, next]
1474  if iterable is not None:
1475  self |= iterable
1476 
1477  def __len__(self):
1478  return len(self.map)
1479 
1480  def __contains__(self, key):
1481  return key in self.map
1482 
1483  def add(self, key):
1484  if key not in self.map:
1485  end = self.end
1486  curr = end[1]
1487  curr[2] = end[1] = self.map[key] = [key, curr, end]
1488 
1489  def discard(self, key):
1490  if key in self.map:
1491  key, prev, next = self.map.pop(key)
1492  prev[2] = next
1493  next[1] = prev
1494 
1495  def __iter__(self):
1496  end = self.end
1497  curr = end[2]
1498  while curr is not end:
1499  yield curr[0]
1500  curr = curr[2]
1501 
1502  def __reversed__(self):
1503  end = self.end
1504  curr = end[1]
1505  while curr is not end:
1506  yield curr[0]
1507  curr = curr[1]
1508 
1509  def pop(self, last=True):
1510  if not self:
1511  raise KeyError('set is empty')
1512  if last:
1513  key = self.end[1][0]
1514  else:
1515  key = self.end[2][0]
1516  self.discard(key)
1517  return key
1518 
1519  def __repr__(self):
1520  if not self:
1521  return '%s()' % (self.__class__.__name__,)
1522  return '%s(%r)' % (self.__class__.__name__, list(self))
1523 
1524  def __eq__(self, other):
1525  if isinstance(other, OrderedSet):
1526  return len(self) == len(other) and list(self) == list(other)
1527  return set(self) == set(other)
1528 
1529 
1530 class OrderedDefaultDict(OrderedDict):
1531  """Store objects in order they were added, but with default type.
1532  Source: http://stackoverflow.com/a/4127426/2608793
1533  """
1534  def __init__(self, *args, **kwargs):
1535  if not args:
1536  self.default_factory = None
1537  else:
1538  if not (args[0] is None or callable(args[0])):
1539  raise TypeError('first argument must be callable or None')
1540  self.default_factory = args[0]
1541  args = args[1:]
1542  super(OrderedDefaultDict, self).__init__(*args, **kwargs)
1543 
1544  def __missing__ (self, key):
1545  if self.default_factory is None:
1546  raise KeyError(key)
1547  self[key] = default = self.default_factory()
1548  return default
1549 
1550  def __reduce__(self): # optional, for pickle support
1551  args = (self.default_factory,) if self.default_factory else ()
1552  return self.__class__, args, None, None, self.iteritems()
1553 
1554 # -------------- PMI2 Tools --------------- #
1555 
1556 def set_coordinates_from_rmf(hier,rmf_fn,frame_num=0):
1557  """Extract frame from RMF file and fill coordinates. Must be identical topology.
1558  @param hier The (System) hierarchy to fill (e.g. after you've built it)
1559  @param rmf_fn The file to extract from
1560  @param frame_num The frame number to extract
1561  """
1562  rh = RMF.open_rmf_file_read_only(rmf_fn)
1563  IMP.rmf.link_hierarchies(rh,[hier])
1564  IMP.rmf.load_frame(rh, RMF.FrameID(frame_num))
1565  del rh
1566 
1567 def input_adaptor(stuff,
1568  pmi_resolution=0,
1569  flatten=False,
1570  selection_tuple=None,
1571  warn_about_slices=True):
1572  """Adapt things for PMI (degrees of freedom, restraints, ...)
1573  Returns list of list of hierarchies, separated into Molecules if possible.
1574  The input can be a list, or a list of lists (iterable of ^1 or iterable of ^2)
1575  (iterable of ^2) Hierarchy -> returns input as list of list of hierarchies,
1576  only one entry, not grouped by molecules.
1577  (iterable of ^2) PMI::System/State/Molecule/TempResidue ->
1578  returns residue hierarchies, grouped in molecules, at requested resolution
1579  @param stuff Can be one of the following inputs:
1580  IMP Hierarchy, PMI System/State/Molecule/TempResidue, or a list/set (of list/set) of them.
1581  Must be uniform input, however. No mixing object types.
1582  @param pmi_resolution For selecting, only does it if you pass PMI objects. Set it to "all"
1583  if you want all resolutions!
1584  @param flatten Set to True if you just want all hierarchies in one list.
1585  @param warn_about_slices Print a warning if you are requesting only part of a bead.
1586  Sometimes you just don't care!
1587  \note since this relies on IMP::atom::Selection, this will not return any objects if they weren't built!
1588  But there should be no problem if you request unbuilt residues, they should be ignored.
1589  """
1590 
1591  if stuff is None:
1592  return stuff
1593 
1594  if hasattr(stuff,'__iter__'):
1595  if len(stuff)==0:
1596  return stuff
1597  thelist=list(stuff)
1598 
1599  # iter of iter of should be ok
1600  if all(hasattr(el,'__iter__') for el in thelist):
1601  thelist = [i for sublist in thelist for i in sublist]
1602  elif any(hasattr(el,'__iter__') for el in thelist):
1603  raise Exception('input_adaptor: input_object must be a list or a list of lists')
1604 
1605  stuff = thelist
1606  else:
1607  stuff = [stuff]
1608 
1609  # check that it is a hierarchy homogenously:
1610  try:
1611  is_hierarchy=all(IMP.atom.Hierarchy.get_is_setup(s) for s in stuff)
1612  except NotImplementedError:
1613  is_hierarchy=False
1614  # get the other types homogenously
1615  is_system=all(isinstance(s, IMP.pmi.topology.System) for s in stuff)
1616  is_state=all(isinstance(s, IMP.pmi.topology.State) for s in stuff)
1617  is_molecule=all(isinstance(s, IMP.pmi.topology.Molecule) for s in stuff)
1618  is_temp_residue=all(isinstance(s, IMP.pmi.topology.TempResidue) for s in stuff)
1619 
1620  # now that things are ok, do selection if requested
1621  hier_list = []
1622  pmi_input = False
1623  if is_system or is_state or is_molecule or is_temp_residue:
1624  # if PMI, perform selection using gathered indexes
1625  pmi_input = True
1626  indexes_per_mol = OrderedDefaultDict(list) #key is Molecule object, value are residues
1627  if is_system:
1628  for system in stuff:
1629  for state in system.get_states():
1630  mdict = state.get_molecules()
1631  for molname in mdict:
1632  for copy in mdict[molname]:
1633  indexes_per_mol[copy] += [r.get_index() for r in copy.get_residues()]
1634  elif is_state:
1635  for state in stuff:
1636  mdict = state.get_molecules()
1637  for molname in mdict:
1638  for copy in mdict[molname]:
1639  indexes_per_mol[copy] += [r.get_index() for r in copy.get_residues()]
1640  elif is_molecule:
1641  for molecule in stuff:
1642  indexes_per_mol[molecule] += [r.get_index() for r in molecule.get_residues()]
1643  elif is_temp_residue:
1644  for tempres in stuff:
1645  indexes_per_mol[tempres.get_molecule()].append(tempres.get_index())
1646  for mol in indexes_per_mol:
1647  if pmi_resolution=='all':
1648  # because you select from the molecule,
1649  # this will start the search from the base resolution
1650  ps = select_at_all_resolutions(mol.get_hierarchy(),
1651  residue_indexes=indexes_per_mol[mol])
1652  else:
1653  sel = IMP.atom.Selection(mol.get_hierarchy(),
1654  resolution=pmi_resolution,
1655  residue_indexes=indexes_per_mol[mol])
1656  ps = sel.get_selected_particles()
1657 
1658  # check that you don't have any incomplete fragments!
1659  if warn_about_slices:
1660  rset = set(indexes_per_mol[mol])
1661  for p in ps:
1663  fset = set(IMP.atom.Fragment(p).get_residue_indexes())
1664  if not fset <= rset:
1665  print('BAD')
1666  minset = min(fset)
1667  maxset = max(fset)
1668  found = fset&rset
1669  minf = min(found)
1670  maxf = max(found)
1671  resbreak = maxf if minf==minset else minset-1
1672  print('WARNING: You are trying to select only part of the bead %s:%i-%i.\n'
1673  'The residues you requested are %i-%i. You can fix this by:\n'
1674  '1) requesting the whole bead/none of it or\n'
1675  '2) break the bead up by passing bead_extra_breaks=[\'%i\'] in '
1676  'molecule.add_representation()'
1677  %(mol.get_name(),minset,maxset,minf,maxf,resbreak))
1678  hier_list.append([IMP.atom.Hierarchy(p) for p in ps])
1679  elif is_hierarchy:
1680  #check
1681  ps=[]
1682  if pmi_resolution=='all':
1683  for h in stuff:
1685  else:
1686  for h in stuff:
1687  ps+=IMP.atom.Selection(h,resolution=pmi_resolution).get_selected_particles()
1688  hier_list=[IMP.atom.Hierarchy(p) for p in ps]
1689  if not flatten:
1690  hier_list = [hier_list]
1691  else:
1692  raise Exception('input_adaptor: you passed something of wrong type or a list with mixed types')
1693 
1694  if flatten and pmi_input:
1695  return [h for sublist in hier_list for h in sublist]
1696  else:
1697  return hier_list
1698 
1699 
1701  """Returns sequence-sorted segments array, each containing the first particle
1702  the last particle and the first residue index."""
1703 
1704  from operator import itemgetter
1705  hiers=IMP.pmi.tools.input_adaptor(mol)
1706  if len(hiers)>1:
1707  raise Exception("IMP.pmi.tools.get_sorted_segments: only pass stuff from one Molecule, please")
1708  hiers = hiers[0]
1709  SortedSegments = []
1710  for h in hiers:
1711  try:
1712  start = IMP.atom.Hierarchy(h).get_children()[0]
1713  except:
1714  start = IMP.atom.Hierarchy(h)
1715 
1716  try:
1717  end = IMP.atom.Hierarchy(h).get_children()[-1]
1718  except:
1719  end = IMP.atom.Hierarchy(h)
1720 
1721  startres = IMP.pmi.tools.get_residue_indexes(start)[0]
1722  endres = IMP.pmi.tools.get_residue_indexes(end)[-1]
1723  SortedSegments.append((start, end, startres))
1724  SortedSegments = sorted(SortedSegments, key=itemgetter(2))
1725  return SortedSegments
1726 
1727 def display_bonds(mol):
1728  """Decorate the sequence-consecutive particles from a PMI2 molecule with a bond,
1729  so that they appear connected in the rmf file"""
1730  SortedSegments=get_sorted_segments(mol)
1731  for x in range(len(SortedSegments) - 1):
1732 
1733  last = SortedSegments[x][1]
1734  first = SortedSegments[x + 1][0]
1735 
1736  p1 = last.get_particle()
1737  p2 = first.get_particle()
1738  if not IMP.atom.Bonded.get_is_setup(p1):
1740  if not IMP.atom.Bonded.get_is_setup(p2):
1742 
1745  IMP.atom.Bonded(p1),
1746  IMP.atom.Bonded(p2),1)
1747 
1748 
1749 class ThreeToOneConverter(defaultdict):
1750  """This class converts three to one letter codes, and return X for any unknown codes"""
1751  def __init__(self,is_nucleic=False):
1752 
1753  if not is_nucleic:
1754  threetoone = {'ALA': 'A', 'ARG': 'R', 'ASN': 'N', 'ASP': 'D',
1755  'CYS': 'C', 'GLU': 'E', 'GLN': 'Q', 'GLY': 'G',
1756  'HIS': 'H', 'ILE': 'I', 'LEU': 'L', 'LYS': 'K',
1757  'MET': 'M', 'PHE': 'F', 'PRO': 'P', 'SER': 'S',
1758  'THR': 'T', 'TRP': 'W', 'TYR': 'Y', 'VAL': 'V', 'UNK': 'X'}
1759  else:
1760  threetoone = {'ADE': 'A', 'URA': 'U', 'CYT': 'C', 'GUA': 'G',
1761  'THY': 'T', 'UNK': 'X'}
1762 
1763  defaultdict.__init__(self,lambda: "X", threetoone)
1764 
1765 
1766 
1767 
1768 def get_residue_type_from_one_letter_code(code,is_nucleic=False):
1769  threetoone=ThreeToOneConverter(is_nucleic)
1770  one_to_three={}
1771  for k in threetoone:
1772  one_to_three[threetoone[k]] = k
1773  return IMP.atom.ResidueType(one_to_three[code])
1774 
1775 
1776 def get_all_leaves(list_of_hs):
1777  """ Just get the leaves from a list of hierarchies """
1778  lvs = list(itertools.chain.from_iterable(IMP.atom.get_leaves(item) for item in list_of_hs))
1779  return lvs
1780 
1781 
1783  hiers=None,
1784  **kwargs):
1785  """Perform selection using the usual keywords but return ALL resolutions (BEADS and GAUSSIANS).
1786  Returns in flat list!
1787  """
1788 
1789  if hiers is None:
1790  hiers = []
1791  if hier is not None:
1792  hiers.append(hier)
1793  if len(hiers)==0:
1794  print("WARNING: You passed nothing to select_at_all_resolutions()")
1795  return []
1796  ret = OrderedSet()
1797  for hsel in hiers:
1798  try:
1799  htest = IMP.atom.Hierarchy.get_is_setup(hsel)
1800  except:
1801  raise Exception('select_at_all_resolutions: you have to pass an IMP Hierarchy')
1802  if not htest:
1803  raise Exception('select_at_all_resolutions: you have to pass an IMP Hierarchy')
1804  if 'resolution' in kwargs or 'representation_type' in kwargs:
1805  raise Exception("don't pass resolution or representation_type to this function")
1806  selB = IMP.atom.Selection(hsel,resolution=IMP.atom.ALL_RESOLUTIONS,
1807  representation_type=IMP.atom.BALLS,**kwargs)
1808  selD = IMP.atom.Selection(hsel,resolution=IMP.atom.ALL_RESOLUTIONS,
1809  representation_type=IMP.atom.DENSITIES,**kwargs)
1810  ret |= OrderedSet(selB.get_selected_particles())
1811  ret |= OrderedSet(selD.get_selected_particles())
1812  return list(ret)
1813 
1814 
1816  target_ps,
1817  sel_zone,
1818  entire_residues,
1819  exclude_backbone):
1820  """Utility to retrieve particles from a hierarchy within a
1821  zone around a set of ps.
1822  @param hier The hierarchy in which to look for neighbors
1823  @param target_ps The particles for zoning
1824  @param sel_zone The maximum distance
1825  @param entire_residues If True, will grab entire residues
1826  @param exclude_backbone If True, will only return sidechain particles
1827  """
1828 
1829  test_sel = IMP.atom.Selection(hier)
1830  backbone_types=['C','N','CB','O']
1831  if exclude_backbone:
1832  test_sel -= IMP.atom.Selection(hier,atom_types=[IMP.atom.AtomType(n)
1833  for n in backbone_types])
1834  test_ps = test_sel.get_selected_particles()
1835  nn = IMP.algebra.NearestNeighbor3D([IMP.core.XYZ(p).get_coordinates()
1836  for p in test_ps])
1837  zone = set()
1838  for target in target_ps:
1839  zone|=set(nn.get_in_ball(IMP.core.XYZ(target).get_coordinates(),sel_zone))
1840  zone_ps = [test_ps[z] for z in zone]
1841  if entire_residues:
1842  final_ps = set()
1843  for z in zone_ps:
1844  final_ps|=set(IMP.atom.Hierarchy(z).get_parent().get_children())
1845  zone_ps = [h.get_particle() for h in final_ps]
1846  return zone_ps
1847 
1848 
1850  """Returns unique objects in original order"""
1851  rbs = set()
1852  beads = []
1853  rbs_ordered = []
1854  if not hasattr(hiers,'__iter__'):
1855  hiers = [hiers]
1856  for p in get_all_leaves(hiers):
1858  rb = IMP.core.RigidMember(p).get_rigid_body()
1859  if rb not in rbs:
1860  rbs.add(rb)
1861  rbs_ordered.append(rb)
1863  rb = IMP.core.NonRigidMember(p).get_rigid_body()
1864  if rb not in rbs:
1865  rbs.add(rb)
1866  rbs_ordered.append(rb)
1867  beads.append(p)
1868  else:
1869  beads.append(p)
1870  return rbs_ordered,beads
1871 
1872 def get_molecules(input_objects):
1873  "This function returns the parent molecule hierarchies of given objects"
1874  stuff=input_adaptor(input_objects, pmi_resolution='all',flatten=True)
1875  molecules=set()
1876  for h in stuff:
1877  is_root=False
1878  is_molecule=False
1879  while not (is_molecule or is_root):
1880  root=IMP.atom.get_root(h)
1881  if root == h:
1882  is_root=True
1883  is_molecule=IMP.atom.Molecule.get_is_setup(h)
1884  if is_molecule:
1885  molecules.add(IMP.atom.Molecule(h))
1886  h=h.get_parent()
1887  return list(molecules)
1888 
1889 def get_densities(input_objects):
1890  """Given a list of PMI objects, returns all density hierarchies within
1891  these objects. The output of this function can be inputted into
1892  things such as EM restraints. This function is intended to gather density particles
1893  appended to molecules (and not other hierarchies which might have been appended to the root node directly).
1894  """
1895  # Note that Densities can only be selected at the Root or Molecule level and not at the Leaves level.
1896  # we'll first get all molecule hierarchies corresponding to the leaves.
1897  molecules=get_molecules(input_objects)
1898  densities=[]
1899  for i in molecules:
1900  densities+=IMP.atom.Selection(i,representation_type=IMP.atom.DENSITIES).get_selected_particles()
1901  return densities
1902 
1903 def shuffle_configuration(objects,
1904  max_translation=300., max_rotation=2.0 * pi,
1905  avoidcollision_rb=True, avoidcollision_fb=False,
1906  cutoff=10.0, niterations=100,
1907  bounding_box=None,
1908  excluded_rigid_bodies=[],
1909  hierarchies_excluded_from_collision=[],
1910  hierarchies_included_in_collision=[],
1911  verbose=False,
1912  return_debug=False):
1913  """Shuffle particles. Used to restart the optimization.
1914  The configuration of the system is initialized by placing each
1915  rigid body and each bead randomly in a box with a side of
1916  max_translation angstroms, and far enough from each other to
1917  prevent any steric clashes. The rigid bodies are also randomly rotated.
1918  @param objects Can be one of the following inputs:
1919  IMP Hierarchy, PMI System/State/Molecule/TempResidue, or a list/set of them
1920  @param max_translation Max translation (rbs and flexible beads)
1921  @param max_rotation Max rotation (rbs only)
1922  @param avoidcollision_rb check if the particle/rigid body was
1923  placed close to another particle; uses the optional
1924  arguments cutoff and niterations
1925  @param avoidcollision_fb Advanced. Generally you want this False because it's hard to shuffle beads.
1926  @param cutoff Distance less than this is a collision
1927  @param niterations How many times to try avoiding collision
1928  @param bounding_box Only shuffle particles within this box. Defined by ((x1,y1,z1),(x2,y2,z2)).
1929  @param excluded_rigid_bodies Don't shuffle these rigid body objects
1930  @param hierarchies_excluded_from_collision Don't count collision with these bodies
1931  @param hierarchies_included_in_collision Hierarchies that are not shuffled, but should be included in collision calculation (for fixed regions)
1932  @param verbose Give more output
1933  \note Best to only call this function after you've set up degrees of freedom
1934  For debugging purposes, returns: <shuffled indexes>, <collision avoided indexes>
1935  """
1936 
1937  ### checking input
1938  hierarchies = IMP.pmi.tools.input_adaptor(objects,
1939  pmi_resolution='all',
1940  flatten=True)
1941  rigid_bodies,flexible_beads = get_rbs_and_beads(hierarchies)
1942  if len(rigid_bodies)>0:
1943  mdl = rigid_bodies[0].get_model()
1944  elif len(flexible_beads)>0:
1945  mdl = flexible_beads[0].get_model()
1946  else:
1947  raise Exception("Could not find any particles in the hierarchy")
1948  if len(rigid_bodies) == 0:
1949  print("shuffle_configuration: rigid bodies were not intialized")
1950 
1951  ### gather all particles
1953  gcpf.set_distance(cutoff)
1954 
1955  # Add particles from excluded hierarchies to excluded list
1956  collision_excluded_hierarchies = IMP.pmi.tools.input_adaptor(hierarchies_excluded_from_collision,
1957  pmi_resolution='all',
1958  flatten=True)
1959 
1960  collision_included_hierarchies = IMP.pmi.tools.input_adaptor(hierarchies_included_in_collision,
1961  pmi_resolution='all',
1962  flatten=True)
1963 
1964  collision_excluded_idxs = set([l.get_particle().get_index() for h in collision_excluded_hierarchies \
1965  for l in IMP.core.get_leaves(h)])
1966 
1967  collision_included_idxs = set([l.get_particle().get_index() for h in collision_included_hierarchies \
1968  for l in IMP.core.get_leaves(h)])
1969 
1970  # Excluded collision with Gaussians
1971  all_idxs = [] #expand to representations?
1972  for p in IMP.pmi.tools.get_all_leaves(hierarchies):
1974  all_idxs.append(p.get_particle_index())
1976  collision_excluded_idxs.add(p.get_particle_index())
1977 
1978  if bounding_box is not None:
1979  ((x1, y1, z1), (x2, y2, z2)) = bounding_box
1980  ub = IMP.algebra.Vector3D(x1, y1, z1)
1981  lb = IMP.algebra.Vector3D(x2, y2, z2)
1982  bb = IMP.algebra.BoundingBox3D(ub, lb)
1983 
1984  all_idxs = set(all_idxs) | collision_included_idxs
1985  all_idxs = all_idxs - collision_excluded_idxs
1986  debug = []
1987  print('shuffling', len(rigid_bodies), 'rigid bodies')
1988  for rb in rigid_bodies:
1989  if rb not in excluded_rigid_bodies:
1990  # gather particles to avoid with this transform
1991  if avoidcollision_rb:
1992  rb_idxs = set(rb.get_member_particle_indexes()) - \
1993  collision_excluded_idxs
1994  other_idxs = all_idxs - rb_idxs
1995  if not other_idxs:
1996  continue
1997 
1998  # iterate, trying to avoid collisions
1999  niter = 0
2000  while niter < niterations:
2001  rbxyz = (rb.get_x(), rb.get_y(), rb.get_z())
2002 
2003  # local transform
2004  if bounding_box:
2005  translation = IMP.algebra.get_random_vector_in(bb)
2006  # First move to origin
2008  -IMP.core.XYZ(rb).get_coordinates())
2009  IMP.core.transform(rb, transformation_orig)
2010  old_coord = IMP.core.XYZ(rb).get_coordinates()
2012  transformation = IMP.algebra.Transformation3D(rotation,
2013  translation)
2014 
2015  else:
2017  rbxyz,
2018  max_translation,
2019  max_rotation)
2020 
2021  debug.append([rb, other_idxs if avoidcollision_rb else set()])
2022  IMP.core.transform(rb, transformation)
2023 
2024  # check collisions
2025  if avoidcollision_rb:
2026  mdl.update()
2027  npairs = len(gcpf.get_close_pairs(mdl,
2028  list(other_idxs),
2029  list(rb_idxs)))
2030  #print("NPAIRS:", npairs)
2031  if npairs==0:
2032  break
2033  else:
2034  niter += 1
2035  if verbose:
2036  print("shuffle_configuration: rigid body placed close to other %d particles, trying again..." % npairs)
2037  print("shuffle_configuration: rigid body name: " + rb.get_name())
2038  if niter == niterations:
2039  raise ValueError("tried the maximum number of iterations to avoid collisions, increase the distance cutoff")
2040  else:
2041  break
2042 
2043  print('shuffling', len(flexible_beads), 'flexible beads')
2044  for fb in flexible_beads:
2045  # gather particles to avoid
2046  if avoidcollision_fb:
2047  fb_idxs = set(IMP.get_indexes([fb]))
2048  other_idxs = all_idxs - fb_idxs
2049  if not other_idxs:
2050  continue
2051 
2052  # iterate, trying to avoid collisions
2053  niter = 0
2054  while niter < niterations:
2055  if bounding_box:
2056  translation = IMP.algebra.get_random_vector_in(bb)
2057  transformation = IMP.algebra.Transformation3D(translation)
2058  else:
2059  fbxyz = IMP.core.XYZ(fb).get_coordinates()
2060  transformation = IMP.algebra.get_random_local_transformation(fbxyz,
2061  max_translation,
2062  max_rotation)
2063 
2064  # For gaussians, treat this fb as an rb
2066  xyz=[fb.get_value(IMP.FloatKey(4)), fb.get_value(IMP.FloatKey(5)), fb.get_value(IMP.FloatKey(6))]
2067  xyz_transformed=transformation.get_transformed(xyz)
2068  if bounding_box:
2069  # Translate to origin
2070  fb.set_value(IMP.FloatKey(4),xyz_transformed[0]-xyz[0])
2071  fb.set_value(IMP.FloatKey(5),xyz_transformed[1]-xyz[1])
2072  fb.set_value(IMP.FloatKey(6),xyz_transformed[2]-xyz[2])
2073  debug.append([xyz,other_idxs if avoidcollision_fb else set()])
2074  xyz2=[fb.get_value(IMP.FloatKey(4)), fb.get_value(IMP.FloatKey(5)), fb.get_value(IMP.FloatKey(6))]
2075  else:
2076  fb.set_value(IMP.FloatKey(4),xyz_transformed[0])
2077  fb.set_value(IMP.FloatKey(5),xyz_transformed[1])
2078  fb.set_value(IMP.FloatKey(6),xyz_transformed[2])
2079  debug.append([xyz,other_idxs if avoidcollision_fb else set()])
2080  else:
2081  d =IMP.core.XYZ(fb)
2082  if bounding_box:
2083  # Translate to origin first
2084  IMP.core.transform(d, -d.get_coordinates())
2085  d =IMP.core.XYZ(fb)
2086  debug.append([d,other_idxs if avoidcollision_fb else set()])
2087  IMP.core.transform(d, transformation)
2088 
2089  if avoidcollision_fb:
2090  mdl.update()
2091  npairs = len(gcpf.get_close_pairs(mdl,
2092  list(other_idxs),
2093  list(fb_idxs)))
2094 
2095  if npairs==0:
2096  break
2097  else:
2098  niter += 1
2099  print("shuffle_configuration: floppy body placed close to other %d particles, trying again..." % npairs)
2100  if niter == niterations:
2101  raise ValueError("tried the maximum number of iterations to avoid collisions, increase the distance cutoff")
2102  else:
2103  break
2104  if return_debug:
2105  return debug
2106 
2107 def color2rgb(colorname):
2108  """Given a chimera color name, return RGB"""
2109  d = {'aquamarine': (0.4980392156862745, 1.0, 0.8313725490196079),
2110  'black': (0.0, 0.0, 0.0),
2111  'blue': (0.0, 0.0, 1.0),
2112  'brown': (0.6470588235294118, 0.16470588235294117, 0.16470588235294117),
2113  'chartreuse': (0.4980392156862745, 1.0, 0.0),
2114  'coral': (1.0, 0.4980392156862745, 0.3137254901960784),
2115  'cornflower blue': (0.39215686274509803, 0.5843137254901961, 0.9294117647058824),
2116  'cyan': (0.0, 1.0, 1.0),
2117  'dark cyan': (0.0, 0.5450980392156862, 0.5450980392156862),
2118  'dark gray': (0.6627450980392157, 0.6627450980392157, 0.6627450980392157),
2119  'dark green': (0.0, 0.39215686274509803, 0.0),
2120  'dark khaki': (0.7411764705882353, 0.7176470588235294, 0.4196078431372549),
2121  'dark magenta': (0.5450980392156862, 0.0, 0.5450980392156862),
2122  'dark olive green': (0.3333333333333333, 0.4196078431372549, 0.1843137254901961),
2123  'dark red': (0.5450980392156862, 0.0, 0.0),
2124  'dark slate blue': (0.2823529411764706, 0.23921568627450981, 0.5450980392156862),
2125  'dark slate gray': (0.1843137254901961, 0.30980392156862746, 0.30980392156862746),
2126  'deep pink': (1.0, 0.0784313725490196, 0.5764705882352941),
2127  'deep sky blue': (0.0, 0.7490196078431373, 1.0),
2128  'dim gray': (0.4117647058823529, 0.4117647058823529, 0.4117647058823529),
2129  'dodger blue': (0.11764705882352941, 0.5647058823529412, 1.0),
2130  'firebrick': (0.6980392156862745, 0.13333333333333333, 0.13333333333333333),
2131  'forest green': (0.13333333333333333, 0.5450980392156862, 0.13333333333333333),
2132  'gold': (1.0, 0.8431372549019608, 0.0),
2133  'goldenrod': (0.8549019607843137, 0.6470588235294118, 0.12549019607843137),
2134  'gray': (0.7450980392156863, 0.7450980392156863, 0.7450980392156863),
2135  'green': (0.0, 1.0, 0.0),
2136  'hot pink': (1.0, 0.4117647058823529, 0.7058823529411765),
2137  'khaki': (0.9411764705882353, 0.9019607843137255, 0.5490196078431373),
2138  'light blue': (0.6784313725490196, 0.8470588235294118, 0.9019607843137255),
2139  'light gray': (0.8274509803921568, 0.8274509803921568, 0.8274509803921568),
2140  'light green': (0.5647058823529412, 0.9333333333333333, 0.5647058823529412),
2141  'light sea green': (0.12549019607843137, 0.6980392156862745, 0.6666666666666666),
2142  'lime green': (0.19607843137254902, 0.803921568627451, 0.19607843137254902),
2143  'magenta': (1.0, 0.0, 1.0),
2144  'medium blue': (0.19607843137254902, 0.19607843137254902, 0.803921568627451),
2145  'medium purple': (0.5764705882352941, 0.4392156862745098, 0.8588235294117647),
2146  'navy blue': (0.0, 0.0, 0.5019607843137255),
2147  'olive drab': (0.4196078431372549, 0.5568627450980392, 0.13725490196078433),
2148  'orange red': (1.0, 0.27058823529411763, 0.0),
2149  'orange': (1.0, 0.4980392156862745, 0.0),
2150  'orchid': (0.8549019607843137, 0.4392156862745098, 0.8392156862745098),
2151  'pink': (1.0, 0.7529411764705882, 0.796078431372549),
2152  'plum': (0.8666666666666667, 0.6274509803921569, 0.8666666666666667),
2153  'purple': (0.6274509803921569, 0.12549019607843137, 0.9411764705882353),
2154  'red': (1.0, 0.0, 0.0),
2155  'rosy brown': (0.7372549019607844, 0.5607843137254902, 0.5607843137254902),
2156  'salmon': (0.9803921568627451, 0.5019607843137255, 0.4470588235294118),
2157  'sandy brown': (0.9568627450980393, 0.6431372549019608, 0.3764705882352941),
2158  'sea green': (0.1803921568627451, 0.5450980392156862, 0.3411764705882353),
2159  'sienna': (0.6274509803921569, 0.3215686274509804, 0.17647058823529413),
2160  'sky blue': (0.5294117647058824, 0.807843137254902, 0.9215686274509803),
2161  'slate gray': (0.4392156862745098, 0.5019607843137255, 0.5647058823529412),
2162  'spring green': (0.0, 1.0, 0.4980392156862745),
2163  'steel blue': (0.27450980392156865, 0.5098039215686274, 0.7058823529411765),
2164  'tan': (0.8235294117647058, 0.7058823529411765, 0.5490196078431373),
2165  'turquoise': (0.25098039215686274, 0.8784313725490196, 0.8156862745098039),
2166  'violet red': (0.8156862745098039, 0.12549019607843137, 0.5647058823529412),
2167  'white': (1.0, 1.0, 1.0),
2168  'yellow': (1.0, 1.0, 0.0)}
2169  return d[colorname]
2170 
2171 class Colors(object):
2172  def __init__(self):
2173  self.colors={
2174  "reds":[("maroon","#800000",(128,0,0)),("dark red","#8B0000",(139,0,0)),
2175  ("brown","#A52A2A",(165,42,42)),("firebrick","#B22222",(178,34,34)),
2176  ("crimson","#DC143C",(220,20,60)),("red","#FF0000",(255,0,0)),
2177  ("tomato","#FF6347",(255,99,71)),("coral","#FF7F50",(255,127,80)),
2178  ("indian red","#CD5C5C",(205,92,92)),("light coral","#F08080",(240,128,128)),
2179  ("dark salmon","#E9967A",(233,150,122)),("salmon","#FA8072",(250,128,114)),
2180  ("light salmon","#FFA07A",(255,160,122)),("orange red","#FF4500",(255,69,0)),
2181  ("dark orange","#FF8C00",(255,140,0))],
2182  "yellows":[("orange","#FFA500",(255,165,0)),("gold","#FFD700",(255,215,0)),
2183  ("dark golden rod","#B8860B",(184,134,11)),("golden rod","#DAA520",(218,165,32)),
2184  ("pale golden rod","#EEE8AA",(238,232,170)),("dark khaki","#BDB76B",(189,183,107)),
2185  ("khaki","#F0E68C",(240,230,140)),("olive","#808000",(128,128,0)),
2186  ("yellow","#FFFF00",(255,255,0)),("antique white","#FAEBD7",(250,235,215)),
2187  ("beige","#F5F5DC",(245,245,220)),("bisque","#FFE4C4",(255,228,196)),
2188  ("blanched almond","#FFEBCD",(255,235,205)),("wheat","#F5DEB3",(245,222,179)),
2189  ("corn silk","#FFF8DC",(255,248,220)),("lemon chiffon","#FFFACD",(255,250,205)),
2190  ("light golden rod yellow","#FAFAD2",(250,250,210)),("light yellow","#FFFFE0",(255,255,224))],
2191  "greens":[("yellow green","#9ACD32",(154,205,50)),("dark olive green","#556B2F",(85,107,47)),
2192  ("olive drab","#6B8E23",(107,142,35)),("lawn green","#7CFC00",(124,252,0)),
2193  ("chart reuse","#7FFF00",(127,255,0)),("green yellow","#ADFF2F",(173,255,47)),
2194  ("dark green","#006400",(0,100,0)),("green","#008000",(0,128,0)),
2195  ("forest green","#228B22",(34,139,34)),("lime","#00FF00",(0,255,0)),
2196  ("lime green","#32CD32",(50,205,50)),("light green","#90EE90",(144,238,144)),
2197  ("pale green","#98FB98",(152,251,152)),("dark sea green","#8FBC8F",(143,188,143)),
2198  ("medium spring green","#00FA9A",(0,250,154)),("spring green","#00FF7F",(0,255,127)),
2199  ("sea green","#2E8B57",(46,139,87)),("medium aqua marine","#66CDAA",(102,205,170)),
2200  ("medium sea green","#3CB371",(60,179,113)),("light sea green","#20B2AA",(32,178,170)),
2201  ("dark slate gray","#2F4F4F",(47,79,79)),("teal","#008080",(0,128,128)),
2202  ("dark cyan","#008B8B",(0,139,139))],
2203  "blues":[("dark turquoise","#00CED1",(0,206,209)),
2204  ("turquoise","#40E0D0",(64,224,208)),("medium turquoise","#48D1CC",(72,209,204)),
2205  ("pale turquoise","#AFEEEE",(175,238,238)),("aqua marine","#7FFFD4",(127,255,212)),
2206  ("powder blue","#B0E0E6",(176,224,230)),("cadet blue","#5F9EA0",(95,158,160)),
2207  ("steel blue","#4682B4",(70,130,180)),("corn flower blue","#6495ED",(100,149,237)),
2208  ("deep sky blue","#00BFFF",(0,191,255)),("dodger blue","#1E90FF",(30,144,255)),
2209  ("light blue","#ADD8E6",(173,216,230)),("sky blue","#87CEEB",(135,206,235)),
2210  ("light sky blue","#87CEFA",(135,206,250)),("midnight blue","#191970",(25,25,112)),
2211  ("navy","#000080",(0,0,128)),("dark blue","#00008B",(0,0,139)),
2212  ("medium blue","#0000CD",(0,0,205)),("blue","#0000FF",(0,0,255)),("royal blue","#4169E1",(65,105,225)),
2213  ("aqua","#00FFFF",(0,255,255)),("cyan","#00FFFF",(0,255,255)),("light cyan","#E0FFFF",(224,255,255))],
2214  "violets":[("blue violet","#8A2BE2",(138,43,226)),("indigo","#4B0082",(75,0,130)),
2215  ("dark slate blue","#483D8B",(72,61,139)),("slate blue","#6A5ACD",(106,90,205)),
2216  ("medium slate blue","#7B68EE",(123,104,238)),("medium purple","#9370DB",(147,112,219)),
2217  ("dark magenta","#8B008B",(139,0,139)),("dark violet","#9400D3",(148,0,211)),
2218  ("dark orchid","#9932CC",(153,50,204)),("medium orchid","#BA55D3",(186,85,211)),
2219  ("purple","#800080",(128,0,128)),("thistle","#D8BFD8",(216,191,216)),
2220  ("plum","#DDA0DD",(221,160,221)),("violet","#EE82EE",(238,130,238)),
2221  ("magenta / fuchsia","#FF00FF",(255,0,255)),("orchid","#DA70D6",(218,112,214)),
2222  ("medium violet red","#C71585",(199,21,133)),("pale violet red","#DB7093",(219,112,147)),
2223  ("deep pink","#FF1493",(255,20,147)),("hot pink","#FF69B4",(255,105,180)),
2224  ("light pink","#FFB6C1",(255,182,193)),("pink","#FFC0CB",(255,192,203))],
2225  "browns":[("saddle brown","#8B4513",(139,69,19)),("sienna","#A0522D",(160,82,45)),
2226  ("chocolate","#D2691E",(210,105,30)),("peru","#CD853F",(205,133,63)),
2227  ("sandy brown","#F4A460",(244,164,96)),("burly wood","#DEB887",(222,184,135)),
2228  ("tan","#D2B48C",(210,180,140)),("rosy brown","#BC8F8F",(188,143,143)),
2229  ("moccasin","#FFE4B5",(255,228,181)),("navajo white","#FFDEAD",(255,222,173)),
2230  ("peach puff","#FFDAB9",(255,218,185)),("misty rose","#FFE4E1",(255,228,225)),
2231  ("lavender blush","#FFF0F5",(255,240,245)),("linen","#FAF0E6",(250,240,230)),
2232  ("old lace","#FDF5E6",(253,245,230)),("papaya whip","#FFEFD5",(255,239,213)),
2233  ("sea shell","#FFF5EE",(255,245,238))],
2234  "greys":[("black","#000000",(0,0,0)),("dim gray / dim grey","#696969",(105,105,105)),
2235  ("gray / grey","#808080",(128,128,128)),("dark gray / dark grey","#A9A9A9",(169,169,169)),
2236  ("silver","#C0C0C0",(192,192,192)),("light gray / light grey","#D3D3D3",(211,211,211)),
2237  ("gainsboro","#DCDCDC",(220,220,220)),("white smoke","#F5F5F5",(245,245,245)),
2238  ("white","#FFFFFF",(255,255,255))]}
2239 
2240  def assign_color_group(self,color_group,representation,component_names):
2241  for n,p in enumerate(component_names):
2242  s=IMP.atom.Selection(representation.prot,molecule=p)
2243  psel=s.get_selected_particles()
2244  ctuple=self.colors[color_group][n]
2245  print("Assigning "+p+" to color "+ctuple[0])
2246  c=ctuple[2]
2247  color=IMP.display.Color(float(c[0])/255,float(c[1])/255,float(c[2])/255)
2248  for part in psel:
2250  IMP.display.Colored(part).set_color(color)
2251  else:
2253 
2254  def get_list_distant_colors(self):
2255  cnames = ['#F0F8FF', '#FAEBD7', '#00FFFF', '#7FFFD4', '#F0FFFF', '#F5F5DC',
2256  '#FFE4C4', '#000000', '#FFEBCD', '#0000FF', '#8A2BE2', '#A52A2A', '#DEB887',
2257  '#5F9EA0', '#7FFF00', '#D2691E', '#FF7F50', '#6495ED', '#FFF8DC', '#DC143C',
2258  '#00FFFF', '#00008B', '#008B8B', '#B8860B', '#A9A9A9', '#006400', '#BDB76B',
2259  '#8B008B', '#556B2F', '#FF8C00', '#9932CC', '#8B0000', '#E9967A', '#8FBC8F',
2260  '#483D8B', '#2F4F4F', '#00CED1', '#9400D3', '#FF1493', '#00BFFF', '#696969',
2261  '#1E90FF', '#B22222', '#FFFAF0', '#228B22', '#FF00FF', '#DCDCDC', '#F8F8FF',
2262  '#FFD700', '#DAA520', '#808080', '#008000', '#ADFF2F', '#F0FFF0', '#FF69B4',
2263  '#CD5C5C', '#4B0082', '#FFFFF0', '#F0E68C', '#E6E6FA', '#FFF0F5', '#7CFC00',
2264  '#FFFACD', '#ADD8E6', '#F08080', '#E0FFFF', '#FAFAD2', '#90EE90', '#D3D3D3',
2265  '#FFB6C1', '#FFA07A', '#20B2AA', '#87CEFA', '#778899', '#B0C4DE', '#FFFFE0',
2266  '#00FF00', '#32CD32', '#FAF0E6', '#FF00FF', '#800000', '#66CDAA', '#0000CD',
2267  '#BA55D3', '#9370DB', '#3CB371', '#7B68EE', '#00FA9A', '#48D1CC', '#C71585',
2268  '#191970', '#F5FFFA', '#FFE4E1', '#FFE4B5', '#FFDEAD', '#000080', '#FDF5E6',
2269  '#808000', '#6B8E23', '#FFA500', '#FF4500', '#DA70D6', '#EEE8AA', '#98FB98',
2270  '#AFEEEE', '#DB7093', '#FFEFD5', '#FFDAB9', '#CD853F', '#FFC0CB', '#DDA0DD',
2271  '#B0E0E6', '#800080', '#FF0000', '#BC8F8F', '#4169E1', '#8B4513', '#FA8072',
2272  '#FAA460', '#2E8B57', '#FFF5EE', '#A0522D', '#C0C0C0', '#87CEEB', '#6A5ACD',
2273  '#708090', '#FFFAFA', '#00FF7F', '#4682B4', '#D2B48C', '#008080', '#D8BFD8',
2274  '#FF6347', '#40E0D0', '#EE82EE', '#F5DEB3', '#FFFFFF', '#F5F5F5', '#FFFF00',
2275  '#9ACD32']
2276  return cnames
static bool get_is_setup(const IMP::ParticleAdaptor &p)
Definition: Molecule.h:35
def list_chunks_iterator
Yield successive length-sized chunks from a list.
Definition: tools.py:1152
Simple 3D transformation class.
static bool get_is_setup(const IMP::ParticleAdaptor &p)
Definition: Residue.h:155
def select_at_all_resolutions
Perform selection using the usual keywords but return ALL resolutions (BEADS and GAUSSIANS).
Definition: tools.py:1782
Represent an RGB color.
Definition: Color.h:24
A decorator to associate a particle with a part of a protein/DNA/RNA.
Definition: Fragment.h:20
def display_bonds
Decorate the sequence-consecutive particles from a PMI2 molecule with a bond, so that they appear con...
Definition: tools.py:1727
def get_restraint_set
Get a RestraintSet containing all PMI restraints added to the model.
Definition: tools.py:50
Store the representations for a system.
Definition: tools.py:869
def remove
index can be a integer
Definition: tools.py:1230
def get_random_cross_link_dataset
Return a random cross-link dataset as a string.
Definition: tools.py:372
def shuffle_configuration
Shuffle particles.
Definition: tools.py:1917
static bool get_is_setup(const IMP::ParticleAdaptor &p)
Definition: atom/Atom.h:241
void add_particles(RMF::FileHandle fh, const ParticlesTemp &hs)
def get_particles_within_zone
Utility to retrieve particles from a hierarchy within a zone around a set of ps.
Definition: tools.py:1815
Set of python classes to create a multi-state, multi-resolution IMP hierarchy.
def __init__
Constructor.
Definition: tools.py:63
static Weight setup_particle(Model *m, ParticleIndex pi)
Definition: Weight.h:31
A decorator for a particle which has bonds.
double get_drmsd(const Vector3DsOrXYZs0 &m0, const Vector3DsOrXYZs1 &m1)
Calculate distance the root mean square deviation between two sets of 3D points.
Definition: atom/distance.h:49
Rotation3D get_random_rotation_3d(const Rotation3D &center, double distance)
Pick a rotation at random near the provided one.
def scatter_and_gather_dict_append
Synchronize data over a parallel run.
Definition: tools.py:1106
def get_molecules
This function returns the parent molecule hierarchies of given objects.
Definition: tools.py:1872
Change color code to hexadecimal to rgb.
Definition: tools.py:1452
IMP::Vector< Color > Colors
Definition: Color.h:188
ParticlesTemp get_particles(Model *m, const ParticleIndexes &ps)
This class stores integers in ordered compact lists eg: [[1,2,3],[6,7,8]] the methods help splitting ...
Definition: tools.py:1172
The type of an atom.
static Surface setup_particle(Model *m, ParticleIndex pi)
Definition: Surface.h:45
void add_particle(RMF::FileHandle fh, Particle *hs)
log
Definition: log.py:1
static bool get_is_setup(const IMP::ParticleAdaptor &p)
Definition: rigid_bodies.h:619
static bool get_is_setup(const IMP::ParticleAdaptor &p)
Definition: XYZ.h:49
GenericHierarchies get_leaves(Hierarchy mhd)
Get all the leaves of the bit of hierarchy.
def get_prot_name_from_particle
Return the component name provided a particle and a list of names.
Definition: tools.py:1017
This class initializes the root node of the global IMP.atom.Hierarchy.
double get_distance(XYZR a, XYZR b)
Compute the sphere distance between a and b.
Definition: XYZR.h:89
def get_residue_gaps_in_hierarchy
Return the residue index gaps and contiguous segments in the hierarchy.
Definition: tools.py:629
def add
index can be a integer or a list of integers
Definition: tools.py:1198
Vector3D get_random_vector_in(const Cylinder3D &c)
Generate a random vector in a cylinder with uniform density.
def get_all_leaves
Just get the leaves from a list of hierarchies.
Definition: tools.py:1776
Add resolution to a particle.
Definition: Resolution.h:24
Bond create_bond(Bonded a, Bonded b, Bond o)
Connect the two wrapped particles by a custom bond.
Object used to hold a set of restraints.
Definition: RestraintSet.h:36
algebra::GridD< 3, algebra::DenseGridStorageD< 3, float >, float > get_grid(DensityMap *in_map)
Return a dense grid containing the voxels of the passed density map.
Stores a named protein chain.
def color2rgb
Given a chimera color name, return RGB.
Definition: tools.py:2107
def input_adaptor
Adapt things for PMI (degrees of freedom, restraints, ...) Returns list of list of hierarchies...
Definition: tools.py:1579
def get_terminal_residue
Get the particle of the terminal residue at the GIVEN resolution (NOTE: not the closest resolution!)...
Definition: tools.py:576
static bool get_is_setup(Model *m, ParticleIndex pi)
Definition: Fragment.h:46
def translate_hierarchy
Apply a translation to a hierarchy along the input vector.
Definition: tools.py:1259
ParticleIndexPairs get_indexes(const ParticlePairsTemp &ps)
static bool get_is_setup(Model *m, ParticleIndex pi)
Definition: Colored.h:48
The standard decorator for manipulating molecular structures.
Ints get_index(const ParticlesTemp &particles, const Subset &subset, const Subsets &excluded)
def select_by_tuple_2
New tuple format: molname OR (start,stop,molname,copynum,statenum) Copy and state are optional...
Definition: tools.py:819
A decorator for a particle representing an atom.
Definition: atom/Atom.h:234
def scatter_and_gather
Synchronize data over a parallel run.
Definition: tools.py:1079
void transform(XYZ a, const algebra::Transformation3D &tr)
Apply a transformation to the particle.
def add_restraint_to_model
Add a PMI restraint to the model.
Definition: tools.py:38
This class converts three to one letter codes, and return X for any unknown codes.
Definition: tools.py:1749
static Bonded setup_particle(Model *m, ParticleIndex pi)
The type for a residue.
def get_terminal_residue_position
Get XYZ coordinates of the terminal residue at the GIVEN resolution.
Definition: tools.py:617
def __init__
index can be a integer or a list of integers
Definition: tools.py:1189
void load_frame(RMF::FileConstHandle file, RMF::FrameID frame)
Load the given RMF frame into the state of the linked objects.
A decorator for a particle with x,y,z coordinates.
Definition: XYZ.h:30
def get_flatten
Returns a flatten list.
Definition: tools.py:1248
static Scale setup_particle(Model *m, ParticleIndex pi)
Definition: Scale.h:30
A base class for Keys.
Definition: Key.h:44
Collect timing information.
Definition: tools.py:58
static Colored setup_particle(Model *m, ParticleIndex pi, Color color)
Definition: Colored.h:62
A decorator for a particle that is part of a rigid body but not rigid.
Definition: rigid_bodies.h:638
def get_closest_residue_position
this function works with plain hierarchies, as read from the pdb, no multi-scale hierarchies ...
Definition: tools.py:545
std::ostream & show(Hierarchy h, std::ostream &out=std::cout)
Print the hierarchy using a given decorator to display each node.
Find all nearby pairs by testing all pairs.
static bool get_is_setup(const IMP::ParticleAdaptor &p)
Definition: rigid_bodies.h:147
def deprecated_object
Python decorator to mark a class as deprecated.
Definition: __init__.py:9583
A decorator for a residue.
Definition: Residue.h:134
General purpose algebraic and geometric methods that are expected to be used by a wide variety of IMP...
static bool get_is_setup(const IMP::ParticleAdaptor &p)
static bool get_is_setup(Model *m, ParticleIndex p)
Check if the particle has the needed attributes for a cast to succeed.
The general base class for IMP exceptions.
Definition: exception.h:49
VectorD< 3 > Vector3D
Definition: VectorD.h:395
Rotation3D get_identity_rotation_3d()
Return a rotation that does not do anything.
Definition: Rotation3D.h:244
def get_densities
Given a list of PMI objects, returns all density hierarchies within these objects.
Definition: tools.py:1889
void link_hierarchies(RMF::FileConstHandle fh, const atom::Hierarchies &hs)
Class to handle individual particles of a Model object.
Definition: Particle.h:41
Bond get_bond(Bonded a, Bonded b)
Get the bond between two particles.
Stores a list of Molecules all with the same State index.
std::string get_data_path(std::string file_name)
Return the full path to one of this module's data files.
def set_coordinates_from_rmf
Extract frame from RMF file and fill coordinates.
Definition: tools.py:1556
def select
this function uses representation=SimplifiedModel it returns the corresponding selected particles rep...
Definition: tools.py:700
Python classes to represent, score, sample and analyze models.
double get_resolution(Model *m, ParticleIndex pi)
Estimate the resolution of the hierarchy as used by Representation.
A decorator for a rigid body.
Definition: rigid_bodies.h:75
def cross_link_db_filter_parser
example '"{ID_Score}" > 28 AND "{Sample}" == "%10_1%" OR ":Sample}" == "%10_2%" OR ":Sample...
Definition: tools.py:477
static bool get_is_setup(const IMP::ParticleAdaptor &p)
Definition: core/Gaussian.h:47
def get_sorted_segments
Returns sequence-sorted segments array, each containing the first particle the last particle and the ...
Definition: tools.py:1700
def get_rbs_and_beads
Returns unique objects in original order.
Definition: tools.py:1849
Hierarchies get_leaves(const Selection &h)
A decorator for a molecule.
Definition: Molecule.h:24
Select hierarchy particles identified by the biological name.
Definition: Selection.h:66
Support for the RMF file format for storing hierarchical molecular data and markup.
def get_residue_indexes
Retrieve the residue indexes for the given particle.
Definition: tools.py:1037
static bool get_is_setup(const IMP::ParticleAdaptor &p)
Definition: rigid_bodies.h:640
Transformation3D get_random_local_transformation(Vector3D origin, double max_translation=5., double max_angle_in_rad=0.26)
Get a local transformation.
Temporarily stores residue information, even without structure available.
Store objects in order they were added, but with default type.
Definition: tools.py:1530
Inferential scoring building on methods developed as part of the Inferential Structure Determination ...
def sublist_iterator
Yield all sublists of length >= lmin and <= lmax.
Definition: tools.py:1133
A particle with a color.
Definition: Colored.h:23