IMP logo
IMP Reference Guide  2.20.0
The Integrative Modeling Platform
/restraints/stereochemistry.py
1 """@namespace IMP.pmi1.restraints.stereochemistry
2 Restraints for keeping correct stereochemistry.
3 """
4 
5 from __future__ import print_function
6 import IMP
7 import IMP.core
8 import IMP.algebra
9 import IMP.atom
10 import IMP.container
11 import IMP.isd
12 import itertools
13 import IMP.pmi1.tools
15 from operator import itemgetter
16 from math import pi,log,sqrt
17 import sys
18 
19 
20 
21 class ConnectivityRestraint(object):
22  """ This class creates a restraint between consecutive TempResidue objects OR an entire
23  PMI MOlecule object. """
24 
25  def __init__(self,
26  objects,
27  scale=1.0,
28  disorderedlength=False,
29  upperharmonic=True,
30  resolution=1,
31  label="None"):
32  """
33  @param objects - a list of hierarchies, PMI TempResidues OR a single Molecule
34  @param scale Scale the maximal distance between the beads by this factor when disorderedlength is False.
35  The maximal distance is calculated as ((float(residuegap) + 1.0) * 3.6) * scale.
36  @param disorderedlength - This flag uses either disordered length
37  calculated for random coil peptides (True) or zero
38  surface-to-surface distance between beads (False)
39  as optimal distance for the sequence connectivity
40  restraint.
41  @param upperharmonic - This flag uses either harmonic (False)
42  or upperharmonic (True) in the intra-pair
43  connectivity restraint.
44  @param resolution - The resolution to connect things at - only used if you pass PMI objects
45  @param label - A string to identify this restraint in the output/stat file
46  """
47  self.label = label
48  self.weight = 1.0
49 
50  hiers = IMP.pmi1.tools.input_adaptor(objects,resolution)
51  if len(hiers)>1:
52  raise Exception("ConnectivityRestraint: only pass stuff from one Molecule, please")
53  hiers = hiers[0]
54 
55  self.kappa = 10 # spring constant used for the harmonic restraints
56  self.m = list(hiers)[0].get_model()
57  SortedSegments = []
58  self.rs = IMP.RestraintSet(self.m, "connectivity_restraint")
59  for h in hiers:
60  try:
61  start = IMP.atom.Hierarchy(h).get_children()[0]
62  except:
63  start = IMP.atom.Hierarchy(h)
64 
65  try:
66  end = IMP.atom.Hierarchy(h).get_children()[-1]
67  except:
68  end = IMP.atom.Hierarchy(h)
69 
70  startres = IMP.pmi1.tools.get_residue_indexes(start)[0]
71  endres = IMP.pmi1.tools.get_residue_indexes(end)[-1]
72  SortedSegments.append((start, end, startres))
73  SortedSegments = sorted(SortedSegments, key=itemgetter(2))
74 
75  # connect the particles
76  self.particle_pairs=[]
77  for x in range(len(SortedSegments) - 1):
78 
79  last = SortedSegments[x][1]
80  first = SortedSegments[x + 1][0]
81 
82  apply_restraint = True
83 
84  # Apply connectivity runless ALL of the following are true:
85  # - first and last both have RigidBodyMember decorators
86  # - first and last are both RigidMembers
87  # - first and last are part of the same RigidBody object
88 
89  # Check for both in a rigid body
92 
93  #Check if the rigid body objects for each particle are the same object.
94  # if so, skip connectivity restraint
95  if IMP.core.RigidBodyMember(first).get_rigid_body() == IMP.core.RigidBodyMember(last).get_rigid_body():
96  apply_restraint = False
97 
98  if apply_restraint:
99 
100  nreslast = len(IMP.pmi1.tools.get_residue_indexes(last))
101  lastresn = IMP.pmi1.tools.get_residue_indexes(last)[-1]
102  nresfirst = len(IMP.pmi1.tools.get_residue_indexes(first))
103  firstresn = IMP.pmi1.tools.get_residue_indexes(first)[0]
104 
105  residuegap = firstresn - lastresn - 1
106  if disorderedlength and (nreslast / 2 + nresfirst / 2 + residuegap) > 20.0:
107  # calculate the distance between the sphere centers using Kohn
108  # PNAS 2004
109  optdist = sqrt(5 / 3) * 1.93 * \
110  (nreslast / 2 + nresfirst / 2 + residuegap) ** 0.6
111  # optdist2=sqrt(5/3)*1.93*((nreslast)**0.6+(nresfirst)**0.6)/2
112  if upperharmonic:
113  hu = IMP.core.HarmonicUpperBound(optdist, self.kappa)
114  else:
115  hu = IMP.core.Harmonic(optdist, self.kappa)
117  else: # default
118  optdist = (0.0 + (float(residuegap) + 1.0) * 3.6) * scale
119  if upperharmonic: # default
120  hu = IMP.core.HarmonicUpperBound(optdist, self.kappa)
121  else:
122  hu = IMP.core.Harmonic(optdist, self.kappa)
124 
125  pt0 = last.get_particle()
126  pt1 = first.get_particle()
127  self.particle_pairs.append((pt0,pt1))
128  r = IMP.core.PairRestraint(self.m, dps, (pt0.get_index(), pt1.get_index()))
129 
130  print("Adding sequence connectivity restraint between", pt0.get_name(), " and ", pt1.get_name(), 'of distance', optdist)
131  self.rs.add_restraint(r)
132 
133  def set_label(self, label):
134  self.label = label
135 
136  def get_weight(self):
137  return self.weight
138 
139  def add_to_model(self):
141 
142  def get_restraint(self):
143  return self.rs
144 
145  def set_weight(self, weight):
146  self.weight = weight
147  self.rs.set_weight(weight)
148 
149  def get_output(self):
150  output = {}
151  score = self.weight * self.rs.unprotected_evaluate(None)
152  output["_TotalScore"] = str(score)
153  output["ConnectivityRestraint_" + self.label] = str(score)
154  return output
155 
157  """ Returns number of connectivity restraints """
158  return len(self.rs.get_restraints())
159 
161  """ Returns the list of connected particles pairs """
162  return self.particle_pairs
163 
164  def evaluate(self):
165  return self.weight * self.rs.unprotected_evaluate(None)
166 
167 class ExcludedVolumeSphere(object):
168  """A class to create an excluded volume restraint for a set of particles at a given resolution.
169  Can be initialized as a bipartite restraint between two sets of particles.
170  # Potential additional function: Variable resolution for each PMI object. Perhaps passing selection_tuples
171  with (PMI_object, resolution)
172  """
173 
174  def __init__(self,
175  representation=None,
176  included_objects=None,
177  other_objects=None,
178  resolution=1000,
179  kappa=1.0):
180  """Constructor.
181  @param representation DEPRECATED - just pass objects
182  @param included_objects Can be one of the following inputs:
183  IMP Hierarchy, PMI System/State/Molecule/TempResidue, or a list/set of them
184  @param other_objects Initializes a bipartite restraint between included_objects and other_objects
185  Same format as included_objects
186  @param resolution The resolution particles at which to impose the restraint.
187  By default, the coarsest particles will be chosen.
188  If a number is chosen, for each particle, the closest
189  resolution will be used (see IMP.atom.Selection).
190  @param kappa Restraint strength
191  """
192 
193  self.weight = 1.0
194  self.kappa = kappa
195  self.label = "None"
196  self.cpc = None
197  bipartite = False
198 
199  # gather IMP hierarchies from input objects
200  hierarchies = IMP.pmi1.tools.input_adaptor(included_objects,
201  resolution,
202  flatten=True)
203  included_ps = []
204  if other_objects is not None:
205  bipartite = True
206  other_hierarchies = IMP.pmi1.tools.input_adaptor(other_objects,
207  resolution,
208  flatten=True)
209  other_ps = []
210 
211 
212  # perform selection
213  if representation is None:
214  if hierarchies is None:
215  raise Exception("Must at least pass included objects")
216  self.mdl = hierarchies[0].get_model()
217  included_ps = [h.get_particle() for h in hierarchies]
218  if bipartite:
219  other_ps = [h.get_particle() for h in other_hierarchies]
220  elif isinstance(representation, IMP.pmi1.representation.Representation):
221  self.mdl = representation.model
222  included_ps = IMP.pmi1.tools.select(
223  representation,
224  resolution=resolution,
225  hierarchies=hierarchies)
226  if bipartite:
227  other_ps = IMP.pmi1.tools.select(
228  representation,
229  resolution=resolution,
230  hierarchies=other_hierarchies)
231  else:
232  raise Exception("Must pass Representation or included_objects")
233 
234  # setup score
235  self.rs = IMP.RestraintSet(self.mdl, 'excluded_volume')
236  ssps = IMP.core.SoftSpherePairScore(self.kappa)
238  lsa.add(IMP.get_indexes(included_ps))
239 
240  # setup close pair container
241  if not bipartite:
243  self.cpc = IMP.container.ClosePairContainer(lsa, 0.0, rbcpf, 10.0)
244  evr = IMP.container.PairsRestraint(ssps, self.cpc)
245  else:
246  other_lsa = IMP.container.ListSingletonContainer(self.mdl)
247  other_lsa.add(IMP.get_indexes(other_ps))
249  lsa,
250  other_lsa,
251  0.0,
252  10.0)
253  evr = IMP.container.PairsRestraint(ssps, self.cpc)
254 
255  self.rs.add_restraint(evr)
256 
257  def add_excluded_particle_pairs(self, excluded_particle_pairs):
258  # add pairs to be filtered when calculating the score
259  inverted=[(p1,p0)for p0,p1 in excluded_particle_pairs]
260  lpc = IMP.container.ListPairContainer(self.mdl)
261  lpc.add(IMP.get_indexes(excluded_particle_pairs))
262  lpc.add(IMP.get_indexes(inverted))
264  self.cpc.add_pair_filter(icpf)
265 
266  def set_label(self, label):
267  self.label = label
268 
269  def add_to_model(self):
270  IMP.pmi1.tools.add_restraint_to_model(self.mdl, self.rs)
271 
272  def get_restraint(self):
273  return self.rs
274 
275  def set_weight(self, weight):
276  self.weight = weight
277  self.rs.set_weight(weight)
278 
279  def get_output(self):
280  output = {}
281  score = self.weight * self.rs.unprotected_evaluate(None)
282  output["_TotalScore"] = str(score)
283  output["ExcludedVolumeSphere_" + self.label] = str(score)
284  return output
285 
286  def evaluate(self):
287  return self.weight * self.rs.unprotected_evaluate(None)
288 
289 class HelixRestraint(object):
290  """Enforce ideal Helix dihedrals and bonds for a selection at resolution 0"""
291  def __init__(self,
292  hierarchy,
293  selection_tuple,
294  weight=1.0,
295  label=''):
296  """Constructor
297  @param hierarchy the root node
298  @param selection_tuple (start, stop, molname, copynum=0)
299  @param weight
300  """
301  self.mdl = hierarchy.get_model()
302  self.rs = IMP.RestraintSet(self.mdl)
303  self.weight = weight
304  self.label = label
305  start = selection_tuple[0]
306  stop = selection_tuple[1]
307  mol = selection_tuple[2]
308  copy_index = 0
309  if len(selection_tuple)>3:
310  copy_index = selection_tuple[3]
311 
312  sel = IMP.atom.Selection(hierarchy,molecule=mol,copy_index=copy_index,
313  residue_indexes=range(start,stop+1))
314  ps = sel.get_selected_particles(with_representation=False)
315  res = [IMP.atom.Residue(p) for p in ps]
316  self.rs = IMP.RestraintSet(self.mdl,self.weight)
317  self.r = IMP.atom.HelixRestraint(res)
318  self.rs.add_restraint(self.r)
319  print('Created helix %s.%i.%i-%i with %i dihedrals and %i bonds'%(
320  mol,copy_index,start,stop,
321  self.get_number_of_bonds(),self.get_number_of_dihedrals()))
322  def set_label(self, label):
323  self.label = label
324 
325  def get_weight(self):
326  return self.weight
327 
328  def add_to_model(self):
329  IMP.pmi1.tools.add_restraint_to_model(self.mdl, self.rs)
330 
331  def get_restraint(self):
332  return self.rs
333 
334  def set_weight(self, weight):
335  self.weight = weight
336  self.rs.set_weight(weight)
337 
338  def get_number_of_bonds(self):
339  return self.r.get_number_of_bonds()
340 
341  def get_number_of_dihedrals(self):
342  return self.r.get_number_of_dihedrals()
343 
344  def evaluate(self):
345  return self.weight * self.rs.unprotected_evaluate(None)
346 
347  def get_output(self):
348  output = {}
349  score = self.evaluate()
350  output["_TotalScore"] = str(score)
351  output["HelixRestraint_" + self.label] = str(score)
352  return output
353 
354 class ResidueBondRestraint(object):
355  """ Add bond restraint between pair of consecutive
356  residues/beads to enforce the stereochemistry.
357  """
358  def __init__(self,
359  representation=None,
360  selection_tuple=None,
361  objects=None,
362  distance=3.78,
363  strength=10.0,
364  jitter=None):
365  """Constructor
366  @param representation (PMI1)
367  @param selection_tuple Requested selection (PMI1)
368  @param objects (PMI2)
369  @param distance Resting distance for restraint
370  @param strength Bond constant
371  @param jitter Defines the +- added to the optimal distance in the harmonic well restraint
372  used to increase the tolerance
373  """
374 
375  if representation is not None and selection_tuple is not None:
376  self.m = representation.prot.get_model()
377  particles = IMP.pmi1.tools.select_by_tuple(
378  representation,
379  selection_tuple,
380  resolution=1)
381 
382  elif objects is not None:
383  particles = IMP.pmi1.tools.input_adaptor(objects,1,flatten=True)
384  self.m = particles[0].get_model()
385 
386  self.rs = IMP.RestraintSet(self.m, "Bonds")
387  self.weight = 1
388  self.label = "None"
389  self.pairslist = []
390 
391 
392 
393  if not jitter:
394  ts = IMP.core.Harmonic(distance, strength)
395  else:
397  (distance - jitter, distance + jitter), strength)
398 
399  for ps in IMP.pmi1.tools.sublist_iterator(particles, 2, 2):
400  pair = []
401  if len(ps) != 2:
402  raise ValueError("wrong length of pair")
403  for p in ps:
405  raise TypeError("not a residue")
406  else:
407  pair.append(p)
408  print("ResidueBondRestraint: adding a restraint between %s %s" % (pair[0].get_name(), pair[1].get_name()))
409  self.rs.add_restraint(
410  IMP.core.DistanceRestraint(self.m, ts, pair[0], pair[1]))
411  self.pairslist.append(IMP.ParticlePair(pair[0], pair[1]))
412  self.pairslist.append(IMP.ParticlePair(pair[1], pair[0]))
413 
414  def set_label(self, label):
415  self.label = label
416  self.rs.set_name(label)
417  for r in self.rs.get_restraints():
418  r.set_name(label)
419 
420  def add_to_model(self):
422 
423  def get_restraint(self):
424  return self.rs
425 
426  def set_weight(self, weight):
427  self.weight = weight
428  self.rs.set_weight(weight)
429 
430  def get_excluded_pairs(self):
431  return self.pairslist
432 
433  def get_output(self):
434  output = {}
435  score = self.weight * self.rs.unprotected_evaluate(None)
436  output["_TotalScore"] = str(score)
437  output["ResidueBondRestraint_" + self.label] = str(score)
438  return output
439 
440 
441 class ResidueAngleRestraint(object):
442  """Add angular restraint between triplets of consecutive
443  residues/beads to enforce the stereochemistry.
444  """
445  def __init__(self,
446  representation=None,
447  selection_tuple=None,
448  objects=None,
449  anglemin=100.0,
450  anglemax=140.0,
451  strength=10.0):
452 
453  if representation is not None and selection_tuple is not None:
454  self.m = representation.prot.get_model()
455  particles = IMP.pmi1.tools.select_by_tuple(
456  representation,
457  selection_tuple,
458  resolution=1)
459 
460  elif objects is not None:
461  particles = IMP.pmi1.tools.input_adaptor(objects,1,flatten=True)
462  self.m = particles[0].get_model()
463 
464  self.rs = IMP.RestraintSet(self.m, "Angles")
465  self.weight = 1
466  self.label = "None"
467  self.pairslist = []
468 
470  (pi * anglemin / 180.0,
471  pi * anglemax / 180.0),
472  strength)
473 
474  for ps in IMP.pmi1.tools.sublist_iterator(particles, 3, 3):
475  triplet = []
476  if len(ps) != 3:
477  raise ValueError("wrong length of triplet")
478  for p in ps:
480  raise TypeError("not a residue")
481  else:
482  triplet.append(p)
483  print("ResidueAngleRestraint: adding a restraint between %s %s %s" % (triplet[0].get_name(), triplet[1].get_name(), triplet[2].get_name()))
484  self.rs.add_restraint(
485  IMP.core.AngleRestraint(triplet[0].get_model(), ts,
486  triplet[0],
487  triplet[1],
488  triplet[2]))
489  self.pairslist.append(IMP.ParticlePair(triplet[0], triplet[2]))
490  self.pairslist.append(IMP.ParticlePair(triplet[2], triplet[0]))
491 
492  def set_label(self, label):
493  self.label = label
494  self.rs.set_name(label)
495  for r in self.rs.get_restraints():
496  r.set_name(label)
497 
498  def add_to_model(self):
500 
501  def get_restraint(self):
502  return self.rs
503 
504  def set_weight(self, weight):
505  self.weight = weight
506  self.rs.set_weight(weight)
507 
508  def get_excluded_pairs(self):
509  return self.pairslist
510 
511  def get_output(self):
512  output = {}
513  score = self.weight * self.rs.unprotected_evaluate(None)
514  output["_TotalScore"] = str(score)
515  output["ResidueAngleRestraint_" + self.label] = str(score)
516  return output
517 
518 
520  """Add dihedral restraints between quadruplet of consecutive
521  residues/beads to enforce the stereochemistry.
522  Give as input a string of "C" and "T", meaning cys (0+-40) or trans (180+-40)
523  dihedral. The length of the string must be \#residue-3.
524  Without the string, the dihedral will be assumed trans.
525  """
526  def __init__(
527  self,
528  representation=None,
529  selection_tuple=None,
530  objects=None,
531  stringsequence=None,
532  strength=10.0):
533 
534  if representation is not None and selection_tuple is not None:
535  self.m = representation.prot.get_model()
536  particles = IMP.pmi1.tools.select_by_tuple(
537  representation,
538  selection_tuple,
539  resolution=1)
540 
541  elif objects is not None:
542  particles = IMP.pmi1.tools.input_adaptor(objects,1,flatten=True)
543  self.m = particles[0].get_model()
544 
545  self.rs = IMP.RestraintSet(self.m, "Angles")
546  self.weight = 1
547  self.label = "None"
548  self.pairslist = []
549 
550  if stringsequence is None:
551  stringsequence = "T" * (len(particles) - 3)
552 
553  for n, ps in enumerate(IMP.pmi1.tools.sublist_iterator(particles, 4, 4)):
554  quadruplet = []
555  if len(ps) != 4:
556  raise ValueError("wrong length of quadruplet")
557  for p in ps:
559  raise TypeError("not a residue")
560  else:
561  quadruplet.append(p)
562  dihedraltype = stringsequence[n]
563  if dihedraltype == "C":
564  anglemin = -20.0
565  anglemax = 20.0
567  (pi * anglemin / 180.0,
568  pi * anglemax / 180.0),
569  strength)
570  print("ResidueDihedralRestraint: adding a CYS restraint between %s %s %s %s" % (quadruplet[0].get_name(), quadruplet[1].get_name(),
571  quadruplet[2].get_name(), quadruplet[3].get_name()))
572  if dihedraltype == "T":
573  anglemin = 180 - 70.0
574  anglemax = 180 + 70.0
576  (pi * anglemin / 180.0,
577  pi * anglemax / 180.0),
578  strength)
579  print("ResidueDihedralRestraint: adding a TRANS restraint between %s %s %s %s" % (quadruplet[0].get_name(), quadruplet[1].get_name(),
580  quadruplet[2].get_name(), quadruplet[3].get_name()))
581  self.rs.add_restraint(
582  IMP.core.DihedralRestraint(self.m,ts,
583  quadruplet[0],
584  quadruplet[1],
585  quadruplet[2],
586  quadruplet[3]))
587  self.pairslist.append(
588  IMP.ParticlePair(quadruplet[0], quadruplet[3]))
589  self.pairslist.append(
590  IMP.ParticlePair(quadruplet[3], quadruplet[0]))
591 
592  def set_label(self, label):
593  self.label = label
594  self.rs.set_name(label)
595  for r in self.rs.get_restraints():
596  r.set_name(label)
597 
598  def add_to_model(self):
600 
601  def get_restraint(self):
602  return self.rs
603 
604  def set_weight(self, weight):
605  self.weight = weight
606  self.rs.set_weight(weight)
607 
608  def get_excluded_pairs(self):
609  return self.pairslist
610 
611  def get_output(self):
612  output = {}
613  score = self.weight * self.rs.unprotected_evaluate(None)
614  output["_TotalScore"] = str(score)
615  output["ResidueDihedralRestraint_" + self.label] = str(score)
616  return output
617 
618 class SecondaryStructure(object):
619  """Experimental, requires isd_emxl for now"""
620  def __init__(self,
621  representation,
622  selection_tuple,
623  ssstring,
624  mixture=False,
625  nativeness=1.0,
626  kt_caff=0.1):
627 
628  if no_isd_emxl:
629  raise ValueError("IMP.isd_emxl is needed")
630 
631  # check that the secondary structure string
632  # is compatible with the ssstring
633 
634  self.particles = IMP.pmi1.tools.select_by_tuple(
635  representation,
636  selection_tuple,
637  resolution=1)
638  self.m = representation.prot.get_model()
639  self.dihe_dict = {}
640  self.ang_dict = {}
641  self.do_mix = {}
642  self.anglfilename = IMP.isd_emxl.get_data_path("CAAngleRestraint.dat")
643  self.dihefilename = IMP.isd_emxl.get_data_path(
644  "CADihedralRestraint.dat")
645  self.nativeness = nativeness
646  self.kt_caff = kt_caff
647  self.anglrs = IMP.RestraintSet(self.m, "Angles")
648  self.dihers = IMP.RestraintSet(self.m, "Dihedrals")
649  self.bondrs = IMP.RestraintSet(self.m, "Bonds")
650  self.label = "None"
651 
652  if len(self.particles) != len(ssstring):
653  print(len(self.particles), len(ssstring))
654  print("SecondaryStructure: residue range and SS string incompatible")
655  self.ssstring = ssstring
656 
657  (bondrslist, anglrslist, diherslist,
658  pairslist) = self.get_CA_force_field()
659  self.pairslist = pairslist
660 
661  # print anglrslist, diherslist, bondrslist, self.particles
662  self.anglrs.add_restraints(anglrslist)
663  self.dihers.add_restraints(diherslist)
664  self.bondrs.add_restraints(bondrslist)
665 
666  def set_label(self, label):
667  self.label = label
668 
669  def add_to_model(self):
670  IMP.pmi1.tools.add_restraint_to_model(self.m, self.anglrs)
671  IMP.pmi1.tools.add_restraint_to_model(self.m, self.dihers)
672  IMP.pmi1.tools.add_restraint_to_model(self.m, self.bondrs)
673 
674  def get_CA_force_field(self):
675  bondrslist = []
676  anglrslist = []
677  diherslist = []
678  pairslist = []
679  # add bonds
680  for res in range(0, len(self.particles) - 1):
681 
682  ps = self.particles[res:res + 2]
683  pairslist.append(IMP.ParticlePair(ps[0], ps[1]))
684  pairslist.append(IMP.ParticlePair(ps[1], ps[0]))
685  br = self.get_distance_restraint(ps[0], ps[1], 3.78, 416.0)
686  br.set_name('Bond_restraint')
687  bondrslist.append(br)
688  # add dihedrals
689  for res in range(0, len(self.particles) - 4):
690 
691  # if res not in dihe_dict: continue
692  # get the appropriate parameters
693  # get the particles
694  ps = self.particles[res:res + 5]
695  [phi0,
696  phi1,
697  score_dih] = self.read_potential_dihedral(
698  self.ssstring[res:res + 4],
699  True)
700  pairslist.append(IMP.ParticlePair(ps[0], ps[3]))
701  pairslist.append(IMP.ParticlePair(ps[3], ps[0]))
702  pairslist.append(IMP.ParticlePair(ps[1], ps[4]))
703  pairslist.append(IMP.ParticlePair(ps[4], ps[1]))
704  dr = IMP.isd_emxl.CADihedralRestraint(
705  ps[0],
706  ps[1],
707  ps[2],
708  ps[3],
709  ps[4],
710  phi0,
711  phi1,
712  score_dih)
713  dr.set_name('Dihedral restraint')
714  diherslist.append(dr)
715  # add angles
716  for res in range(0, len(self.particles) - 2):
717  ps = self.particles[res:res + 3]
718  [psi, score_ang] = self.read_potential_angle(
719  self.ssstring[res:res + 2], True)
720  pairslist.append(IMP.ParticlePair(ps[0], ps[2]))
721  pairslist.append(IMP.ParticlePair(ps[2], ps[0]))
722  dr = IMP.isd_emxl.CAAngleRestraint(
723  ps[0],
724  ps[1],
725  ps[2],
726  psi,
727  score_ang)
728  dr.set_name('Angle restraint')
729  anglrslist.append(dr)
730  return (bondrslist, anglrslist, diherslist, pairslist)
731 
732  def read_potential_dihedral(self, string, mix=False):
733  # read potentials for dihedral
734  score_dih = []
735  phi0 = []
736  phi1 = []
737  for i in range(0, 36):
738  phi0.append(i * 10.0 / 180.0 * pi)
739  phi1.append(i * 10.0 / 180.0 * pi)
740  for j in range(0, 36):
741  score_dih.append(0.0)
742  # open file
743  if not mix:
744  f = open(self.dihefilename, 'r')
745  for line in f.readlines():
746  riga = (line.strip()).split()
747  if (len(riga) == 4 and riga[0] == string):
748  ii = int(float(riga[1]) / 10.0)
749  jj = int(float(riga[2]) / 10.0)
750  score_dih[ii * len(phi0) + jj] = - \
751  self.kt_caff * self.log(float(riga[3]))
752  f.close()
753  if mix:
754  # mix random coil and native secondary structure
755  counts = []
756  for i in range(0, 36):
757  for j in range(0, 36):
758  counts.append(1.0)
759  f = open(self.dihefilename, 'r')
760  for line in f.readlines():
761  riga = (line.strip()).split()
762  if (len(riga) == 4 and riga[0] == string):
763  ii = int(float(riga[1]) / 10.0)
764  jj = int(float(riga[2]) / 10.0)
765  counts[ii * len(phi0) + jj] += self.nativeness * \
766  float(riga[3])
767  if (len(riga) == 4 and riga[0] == "-----"):
768  ii = int(float(riga[1]) / 10.0)
769  jj = int(float(riga[2]) / 10.0)
770  counts[ii * len(phi0) + jj] += (1.0 - self.nativeness) * \
771  float(riga[3])
772  f.close()
773  for i in range(len(counts)):
774  score_dih[i] = -self.kt_caff * self.log(counts[i])
775  return [phi0, phi1, score_dih]
776 
777  def read_potential_angle(self, string, mix=False):
778  # read potentials for angles
779  score_ang = []
780  psi = []
781  for i in range(0, 180):
782  psi.append(i / 180.0 * pi)
783  score_ang.append(0.0)
784  # read file
785  if not mix:
786  f = open(self.anglfilename, 'r')
787  for line in f.readlines():
788  riga = (line.strip()).split()
789  if (len(riga) == 3 and riga[0] == string):
790  ii = int(riga[1])
791  score_ang[ii] = -self.kt_caff * self.log(float(riga[2]))
792  f.close()
793  if mix:
794  # mix random coil and native secondary structure
795  counts = []
796  for i in range(0, 180):
797  counts.append(1.0)
798 
799  f = open(self.anglfilename, 'r')
800  for line in f.readlines():
801  riga = (line.strip()).split()
802  if (len(riga) == 3 and riga[0] == string):
803  ii = int(riga[1])
804  counts[ii] += self.nativeness * float(riga[2])
805  if (len(riga) == 3 and riga[0] == "---"):
806  ii = int(riga[1])
807  counts[ii] += (1.0 - self.nativeness) * float(riga[2])
808  f.close()
809  for i in range(0, 180):
810  score_ang[i] = -self.kt_caff * self.log(counts[i])
811  return [psi, score_ang]
812 
813  def get_excluded_pairs(self):
814  return self.pairslist
815 
816  def get_restraint(self):
817  tmprs = IMP.RestraintSet(self.m, 'tmp')
818  tmprs.add_restraint(self.anglrs)
819  tmprs.add_restraint(self.dihers)
820  tmprs.add_restraint(self.bondrs)
821  return tmprs
822 
823  def get_distance_restraint(self, p0, p1, d0, kappa):
824  h = IMP.core.Harmonic(d0, kappa)
826  pr = IMP.core.PairRestraint(self.m, dps, IMP.ParticlePair(p0, p1))
827  return pr
828 
829  def get_output(self):
830  output = {}
831  score_angle = self.anglrs.unprotected_evaluate(None)
832  score_dihers = self.dihers.unprotected_evaluate(None)
833  score_bondrs = self.bondrs.unprotected_evaluate(None)
834  output["_TotalScore"] = str(score_angle + score_dihers + score_bondrs)
835 
836  output["SecondaryStructure_Angles_" + self.label] = str(score_angle)
837  output["SecondaryStructure_Dihedrals_" +
838  self.label] = str(score_dihers)
839  output["SecondaryStructure_Bonds_" + self.label] = str(score_bondrs)
840  return output
841 
842 
844  """Add harmonic restraints between all pairs
845  """
846  def __init__(self,representation=None,
847  selection_tuples=None,
848  resolution=1,
849  strength=0.01,
850  dist_cutoff=10.0,
851  ca_only=True,
852  hierarchy=None):
853  """Constructor
854  @param representation Representation object
855  @param selection_tuples Selecting regions for the restraint [[start,stop,molname,copy_index=0],...]
856  @param resolution Resolution for applying restraint
857  @param strength Bond strength
858  @param dist_cutoff Cutoff for making restraints
859  @param ca_only Selects only CAlphas. Only matters if resolution=0.
860  @param hierarchy Root hierarchy to select from, use this instead of representation
861  """
862 
863  particles = []
864  if representation is None and hierarchy is not None:
865  self.m = hierarchy.get_model()
866  for st in selection_tuples:
867  copy_index=0
868  if len(st)>3:
869  copy_index=st[3]
870  if not ca_only:
871  sel = IMP.atom.Selection(hierarchy,molecule=st[2],residue_indexes=range(st[0],st[1]+1),
872  copy_index=copy_index)
873  else:
874  sel = IMP.atom.Selection(hierarchy,molecule=st[2],residue_indexes=range(st[0],st[1]+1),
875  copy_index=copy_index,
876  atom_type=IMP.atom.AtomType("CA"))
877  particles+=sel.get_selected_particles()
878  elif representation is not None and type(representation)==IMP.pmi1.representation.Representation:
879  self.m = representation.model
880  for st in selection_tuples:
881  print('selecting with',st)
882  for p in IMP.pmi1.tools.select_by_tuple(representation,st,resolution=resolution):
883  if (resolution==0 and ca_only and IMP.atom.Atom(p).get_atom_type()!=IMP.atom.AtomType("CA")):
884  continue
885  else:
886  particles.append(p.get_particle())
887  else:
888  raise Exception("must pass representation or hierarchy")
889 
890  self.weight = 1
891  self.label = "None"
892  self.pairslist = []
893 
894  # create score
895  self.rs = IMP.pmi1.create_elastic_network(particles,dist_cutoff,strength)
896  for r in self.rs.get_restraints():
897  a1,a2 = r.get_inputs()
898  self.pairslist.append(IMP.ParticlePair(a1,a2))
899  self.pairslist.append(IMP.ParticlePair(a2,a1))
900  print('ElasticNetwork: created',self.rs.get_number_of_restraints(),'restraints')
901 
902  def set_label(self, label):
903  self.label = label
904  self.rs.set_name(label)
905  for r in self.rs.get_restraints():
906  r.set_name(label)
907 
908  def add_to_model(self):
910 
911  def get_restraint(self):
912  return self.rs
913 
914  def set_weight(self, weight):
915  self.weight = weight
916  self.rs.set_weight(weight)
917 
918  def get_excluded_pairs(self):
919  return self.pairslist
920 
921  def get_output(self):
922  output = {}
923  score = self.weight * self.rs.unprotected_evaluate(None)
924  output["_TotalScore"] = str(score)
925  output["ElasticNetworkRestraint_" + self.label] = str(score)
926  return output
927 
928 
930  """ Enable CHARMM force field """
931  def __init__(self,
932  root=None,
933  ff_temp=300.0,
934  zone_ps=None,
935  zone_size=10.0,
936  enable_nonbonded=True,
937  enable_bonded=True,
938  zone_nonbonded=False,
939  representation=None):
940  """Setup the CHARMM restraint on a selection. Expecting atoms.
941  @param root The node at which to apply the restraint
942  @param ff_temp The temperature of the force field
943  @param zone_ps Create a zone around this set of particles
944  Automatically includes the entire residue (incl. backbone)
945  @param zone_size The size for looking for neighbor residues
946  @param enable_nonbonded Allow the repulsive restraint
947  @param enable_bonded Allow the bonded restraint
948  @param zone_nonbonded EXPERIMENTAL: exclude from nonbonded all sidechains that aren't in zone!
949  @param representation Legacy representation object
950  """
951 
952  kB = (1.381 * 6.02214) / 4184.0
953  if representation is not None:
954  root = representation.prot
955 
956  self.mdl = root.get_model()
957  self.bonds_rs = IMP.RestraintSet(self.mdl, 1.0 / (kB * ff_temp), 'BONDED')
958  self.nonbonded_rs = IMP.RestraintSet(self.mdl, 1.0 / (kB * ff_temp), 'NONBONDED')
959  self.weight = 1.0
960  self.label = ""
961 
962  ### setup topology and bonds etc
963  if enable_bonded:
965  topology = ff.create_topology(root)
966  topology.apply_default_patches()
967  topology.setup_hierarchy(root)
968  if zone_ps is not None:
969  limit_to_ps=IMP.pmi1.topology.get_particles_within_zone(
970  root,
971  zone_ps,
972  zone_size,
973  entire_residues=True,
974  exclude_backbone=False)
976  topology,
977  limit_to_ps)
978  self.ps = limit_to_ps
979  else:
980  r = IMP.atom.CHARMMStereochemistryRestraint(root, topology)
981  self.ps = IMP.core.get_leaves(root)
982  print('init bonds score',r.unprotected_evaluate(None))
983  self.bonds_rs.add_restraint(r)
984  ff.add_radii(root)
985  ff.add_well_depths(root)
986 
987  atoms = IMP.atom.get_by_type(root,IMP.atom.ATOM_TYPE)
988  ### non-bonded forces
989  if enable_nonbonded:
990  if (zone_ps is not None) and zone_nonbonded:
991  print('stereochemistry: zone_nonbonded is True')
992  # atoms list should only include backbone plus zone_ps!
993  backbone_types=['C','N','CB','O']
994  sel = IMP.atom.Selection(root,atom_types=[IMP.atom.AtomType(n)
995  for n in backbone_types])
996  backbone_atoms = sel.get_selected_particles()
997  #x = list(set(backbone_atoms)|set(limit_to_ps))
998  sel_ps=IMP.pmi1.topology.get_particles_within_zone(
999  root,
1000  zone_ps,
1001  zone_size,
1002  entire_residues=True,
1003  exclude_backbone=True)
1004 
1006  IMP.container.ListSingletonContainer(backbone_atoms),
1008  4.0)
1009  else:
1010  cont = IMP.container.ListSingletonContainer(self.mdl,atoms)
1011  self.nbl = IMP.container.ClosePairContainer(cont, 4.0)
1012  if enable_bonded:
1013  self.nbl.add_pair_filter(r.get_full_pair_filter())
1014  pairscore = IMP.isd.RepulsiveDistancePairScore(0,1)
1015  pr=IMP.container.PairsRestraint(pairscore, self.nbl)
1016  self.nonbonded_rs.add_restraint(pr)
1017  print('CHARMM is set up')
1018 
1019  def set_label(self, label):
1020  self.label = label
1021  self.rs.set_name(label)
1022  for r in self.rs.get_restraints():
1023  r.set_name(label)
1024 
1025  def add_to_model(self):
1026  IMP.pmi1.tools.add_restraint_to_model(self.mdl, self.bonds_rs)
1027  IMP.pmi1.tools.add_restraint_to_model(self.mdl, self.nonbonded_rs)
1028 
1029  def get_restraint(self):
1030  return self.rs
1031 
1032  def get_close_pair_container(self):
1033  return self.nbl
1034 
1035  def set_weight(self, weight):
1036  self.weight = weight
1037  self.rs.set_weight(weight)
1038 
1039  def get_output(self):
1040  output = {}
1041  bonds_score = self.weight * self.bonds_rs.unprotected_evaluate(None)
1042  nonbonded_score = self.weight * self.nonbonded_rs.unprotected_evaluate(None)
1043  score=bonds_score+nonbonded_score
1044  output["_TotalScore"] = str(score)
1045  output["CHARMM_BONDS"] = str(bonds_score)
1046  output["CHARMM_NONBONDED"] = str(nonbonded_score)
1047  return output
1048 
1049 
1051  """Add bonds and improper dihedral restraints for the CBs
1052  """
1053  def __init__(
1054  self, rnums, representation, selection_tuple, strength=10.0, kappa=1.0,
1055  jitter_angle=0.0, jitter_improper=0.0):
1056  '''
1057  need to add:
1058  ca-ca bond
1059  ca-cb is a constraint, no restraint needed
1060  ca-ca-ca
1061  cb-ca-ca-cb
1062  '''
1063 
1064  self.m = representation.prot.get_model()
1065  self.rs = IMP.RestraintSet(self.m, "PseudoAtomic")
1066  self.rset_angles = IMP.RestraintSet(self.m, "PseudoAtomic_Angles")
1067  self.rset_bonds = IMP.RestraintSet(self.m, "PseudoAtomic_Bonds")
1068  self.weight = 1
1069  self.label = "None"
1070  self.pairslist = []
1071 
1072  # residues=IMP.pmi1.tools.select_by_tuple(representation,selection_tuple,resolution=1)
1073  for rnum in rnums:
1074  ca, cb = self.get_ca_cb(
1075  IMP.pmi1.tools.select_by_tuple(representation,
1076  (rnum, rnum, 'chainA'), resolution=0))
1077  if not cb is None:
1078  nter = False
1079  cter = False
1080  ca_prev, cb_prev = self.get_ca_cb(
1081  IMP.pmi1.tools.select_by_tuple(representation,
1082  (rnum - 1, rnum - 1, 'chainA'), resolution=0))
1083  ca_next, cb_next = self.get_ca_cb(
1084  IMP.pmi1.tools.select_by_tuple(representation,
1085  (rnum + 1, rnum + 1, 'chainA'), resolution=0))
1086  if ca_prev is None:
1087  nter = True
1088  else:
1089  if ca_next is None:
1090  cter = True
1091  else:
1092  if (nter and cter):
1093  continue
1094 
1095  # adding a bond restraint between CA and CB
1096  # h=IMP.core.Harmonic(6.0,kappa)
1097  # dps=IMP.core.DistancePairScore(h)
1098  # pr=IMP.core.PairRestraint(dps,IMP.ParticlePair(ca,cb))
1099  # self.pairslist.append((ca,cb))
1100  # self.rset_bonds.add_restraint(pr)
1101 
1102  # creating improper dihedral restraint
1103  # hus=IMP.core.Harmonic(2.09,kappa)
1105  2.09 +
1106  jitter_angle /
1107  0.5,
1108  kappa)
1110  2.09 -
1111  jitter_angle /
1112  0.5,
1113  kappa)
1114  if not nter:
1115  # ar13=IMP.core.AngleRestraint(hus,ca_prev,ca,cb)
1116  # self.rset_angles.add_restraint(ar13)
1117  ar13u = IMP.core.AngleRestraint(hupp, ca_prev, ca, cb)
1118  ar13l = IMP.core.AngleRestraint(hlow, ca_prev, ca, cb)
1119  self.rset_angles.add_restraint(ar13u)
1120  self.rset_angles.add_restraint(ar13l)
1121  if not cter:
1122  # ar23=IMP.core.AngleRestraint(hus,ca_next,ca,cb)
1123  # self.rset_angles.add_restraint(ar23)
1124  ar23u = IMP.core.AngleRestraint(hupp, ca_next, ca, cb)
1125  ar23l = IMP.core.AngleRestraint(hlow, ca_next, ca, cb)
1126  self.rset_angles.add_restraint(ar23u)
1127  self.rset_angles.add_restraint(ar23l)
1128  if not nter and not cter:
1129  # hus2=IMP.core.Harmonic(0,kappa)
1130  # idr=IMP.core.DihedralRestraint(hus2,ca,ca_prev,ca_next,cb)
1131  # self.rset_angles.add_restraint(idr)
1132 
1133  hus2upp = IMP.core.HarmonicUpperBound(
1134  jitter_improper,
1135  kappa)
1136  hus2low = IMP.core.HarmonicLowerBound(
1137  -
1138  jitter_improper,
1139  kappa)
1141  hus2upp,
1142  ca,
1143  ca_prev,
1144  ca_next,
1145  cb)
1147  hus2low,
1148  ca,
1149  ca_prev,
1150  ca_next,
1151  cb)
1152  self.rset_angles.add_restraint(idru)
1153  self.rset_angles.add_restraint(idrl)
1154  self.rs.add_restraint(self.rset_bonds)
1155  self.rs.add_restraint(self.rset_angles)
1156 
1157  def get_ca_cb(self, atoms):
1158  ca = None
1159  cb = None
1160  for a in atoms:
1161  if IMP.atom.Atom(a).get_atom_type() == IMP.atom.AtomType("CA"):
1162  ca = a.get_particle()
1163  elif IMP.atom.Atom(a).get_atom_type() == IMP.atom.AtomType("CB"):
1164  cb = a.get_particle()
1165  return ca, cb
1166 
1167  def set_label(self, label):
1168  self.label = label
1169  self.rs.set_name(label)
1170  for r in self.rs.get_restraints():
1171  r.set_name(label)
1172 
1173  def add_to_model(self):
1174  IMP.pmi1.tools.add_restraint_to_model(self.m, self.rs)
1175 
1176  def get_restraint(self):
1177  return self.rs
1178 
1179  def set_weight(self, weight):
1180  self.weight = weight
1181  self.rs.set_weight(weight)
1182 
1183  def get_excluded_pairs(self):
1184  return self.pairslist
1185 
1186  def get_output(self):
1187  output = {}
1188  score = self.weight * self.rs.unprotected_evaluate(None)
1189  output["_TotalScore"] = str(score)
1190  output["PseudoAtomicRestraint_" + self.label] = str(score)
1191  return output
1192 
1193 class SymmetryRestraint(object):
1194  """Create harmonic restraints between the reference and (transformed) clones.
1195  \note Wraps IMP::core::TransformedDistancePairScore with an IMP::core::Harmonic
1196  """
1197  def __init__(self,
1198  references,
1199  clones_list,
1200  transforms,
1201  label='',
1202  strength=10.0,
1203  ca_only=False):
1204  """Constructor
1205  @param references Can be one of the following inputs:
1206  IMP Hierarchy, PMI System/State/Molecule/TempResidue, or a list/set of them
1207  @param clones_list List of lists of the above
1208  @param transforms Transforms moving each selection to the first selection
1209  @param label Label for output
1210  @param strength The elastic bond strength
1211  @param ca_only Optionally select so only CAlpha particles are used
1212  """
1213 
1214  refs = IMP.pmi1.tools.input_adaptor(references,flatten=True)
1215  self.mdl = refs[0].get_model()
1216  self.rs = IMP.RestraintSet(self.mdl, "Symmetry")
1217  self.weight = 1
1218  self.label = label
1219  if len(clones_list)!=len(transforms):
1220  raise Exception('Error: There should be as many clones as transforms')
1221 
1222  harmonic = IMP.core.Harmonic(0.,strength)
1223  for tmp_clones,trans in zip(clones_list,transforms):
1224  clones = IMP.pmi1.tools.input_adaptor(tmp_clones,flatten=True)
1225  if len(clones)!=len(refs):
1226  raise Exception("Error: len(references)!=len(clones)")
1227  pair_score = IMP.core.TransformedDistancePairScore(harmonic,trans)
1228  for p0,p1 in zip(refs,clones):
1229  if not ca_only or (IMP.atom.Atom(p0).get_atom_type()==IMP.atom.AtomType("CA")
1230  and IMP.atom.Atom(p1).get_atom_type()==IMP.atom.AtomType("CA")):
1231  r = IMP.core.PairRestraint(self.mdl, pair_score, [p0.get_particle_index(),
1232  p1.get_particle_index()])
1233  self.rs.add_restraint(r)
1234  print('created symmetry network with',self.rs.get_number_of_restraints(),'restraints')
1235 
1236  def set_label(self, label):
1237  self.label = label
1238  self.rs.set_name(label)
1239  for r in self.rs.get_restraints():
1240  r.set_name(label)
1241 
1242  def add_to_model(self):
1243  IMP.pmi1.tools.add_restraint_to_model(self.mdl, self.rs)
1244 
1245  def get_restraint(self):
1246  return self.rs
1247 
1248  def set_weight(self, weight):
1249  self.weight = weight
1250  self.rs.set_weight(weight)
1251 
1252  def get_excluded_pairs(self):
1253  return self.pairslist
1254 
1255  def get_output(self):
1256  output = {}
1257  score = self.weight * self.rs.unprotected_evaluate(None)
1258  output["SymmetryRestraint_" + self.label] = str(score)
1259  output["_TotalScore"] = str(score)
1260  return output
1261 
1262 
1263 class FusionRestraint(object):
1264  """ This class creates a restraint between the termini two polypeptides, to simulate the sequence connectivity. """
1265  def __init__(self,
1266  nterminal,
1267  cterminal,
1268  scale=1.0,
1269  disorderedlength=False,
1270  upperharmonic=True,
1271  resolution=1,
1272  label="None"):
1273  """
1274  @param nterminal - single PMI2 Hierarchy/molecule at the nterminal
1275  @param cterminal - single PMI2 Hierarchy/molecule at the cterminal
1276  @param scale Scale the maximal distance between the beads by this factor when disorderedlength is False.
1277  The maximal distance is calculated as ((float(residuegap) + 1.0) * 3.6) * scale.
1278  @param disorderedlength - This flag uses either disordered length
1279  calculated for random coil peptides (True) or zero
1280  surface-to-surface distance between beads (False)
1281  as optimal distance for the sequence connectivity
1282  restraint.
1283  @param upperharmonic - This flag uses either harmonic (False)
1284  or upperharmonic (True) in the intra-pair
1285  connectivity restraint.
1286  @param resolution - The resolution to connect things at - only used if you pass PMI objects
1287  @param label - A string to identify this restraint in the output/stat file
1288  """
1289  self.label = label
1290  self.weight = 1.0
1291  ssn=IMP.pmi1.tools.get_sorted_segments(nterminal)
1292  ssc=IMP.pmi1.tools.get_sorted_segments(cterminal)
1293  nter_lastres=ssn[-1][1]
1294  cter_firstres=ssc[0][0]
1295  self.m = nter_lastres.get_model()
1296 
1297  self.kappa = 10 # spring constant used for the harmonic restraints
1298 
1299  optdist = (3.6) * scale
1300  if upperharmonic: # default
1301  hu = IMP.core.HarmonicUpperBound(optdist, self.kappa)
1302  else:
1303  hu = IMP.core.Harmonic(optdist, self.kappa)
1305 
1306  pt0 = nter_lastres.get_particle()
1307  pt1 = cter_firstres.get_particle()
1308  r = IMP.core.PairRestraint(self.m, dps, (pt0.get_index(), pt1.get_index()))
1309  self.rs = IMP.RestraintSet(self.m, "fusion_restraint")
1310  print("Adding fusion connectivity restraint between", pt0.get_name(), " and ", pt1.get_name(), 'of distance', optdist)
1311  self.rs.add_restraint(r)
1312 
1313  def set_label(self, label):
1314  self.label = label
1315 
1316  def get_weight(self):
1317  return self.weight
1318 
1319  def add_to_model(self):
1320  IMP.pmi1.tools.add_restraint_to_model(self.m, self.rs)
1321 
1322  def get_restraint(self):
1323  return self.rs
1324 
1325  def set_weight(self, weight):
1326  self.weight = weight
1327  self.rs.set_weight(weight)
1328 
1329  def get_output(self):
1330  output = {}
1331  score = self.weight * self.rs.unprotected_evaluate(None)
1332  output["_TotalScore"] = str(score)
1333  output["FusionRestraint_" + self.label] = str(score)
1334  return output
1335 
1336  def evaluate(self):
1337  return self.weight * self.rs.unprotected_evaluate(None)
1338 
1339 
1340 
1341 
1343 
1344  """Restrain the dihedral between planes defined by three particles.
1345 
1346  This restraint is useful for restraining the twist of a string of
1347  more or less identical rigid bodies, so long as the curvature is mild.
1348  """
1349 
1350  def __init__(self, particle_triplets, angle=0., k=1., label=None,
1351  weight=1.):
1352  """Constructor
1353  @param particle_triplets List of lists of 3 particles. Each triplet
1354  defines a plane. Dihedrals of adjacent planes
1355  in list are scored.
1356  @param angle Angle of plane dihedral in degrees
1357  @param k Strength of restraint
1358  @param label Label for output
1359  @param weight Weight of restraint
1360  \note Particles defining planes should be rigid and more or less
1361  parallel for proper behavior
1362  """
1363  model = particle_triplets[0][0].get_model()
1364  super(PlaneDihedralRestraint, self).__init__(model, label=label,
1365  weight=weight)
1366 
1367  angle = pi * angle / 180.
1368  ds = IMP.core.Cosine(.5 * k, 1., -angle)
1369  for i, t1 in enumerate(particle_triplets[:-1]):
1370  t2 = particle_triplets[i + 1]
1371  q1 = [t1[1], t1[0], t2[0], t2[1]]
1372  q2 = [t1[2], t1[0], t2[0], t2[2]]
1373  self.rs.add_restraint(
1374  IMP.core.DihedralRestraint(self.model, ds, *q1))
1375  self.rs.add_restraint(
1376  IMP.core.DihedralRestraint(self.model, ds, *q2))
def sublist_iterator
Yield all sublists of length >= lmin and <= lmax.
Definition: /tools.py:1168
A filter which returns true if a container contains the Pair.
CHARMMParameters * get_heavy_atom_CHARMM_parameters()
A class to create an excluded volume restraint for a set of particles at a given resolution.
Lower bound harmonic function (non-zero when feature < mean)
static bool get_is_setup(const IMP::ParticleAdaptor &p)
Definition: Residue.h:158
Enforce CHARMM stereochemistry on the given Hierarchy.
def input_adaptor
Adapt things for PMI (degrees of freedom, restraints, ...) Returns list of list of hierarchies...
Definition: /tools.py:1621
def get_num_restraints
Returns number of connectivity restraints.
A member of a rigid body, it has internal (local) coordinates.
Definition: rigid_bodies.h:634
static bool get_is_setup(const IMP::ParticleAdaptor &p)
Definition: rigid_bodies.h:635
Various classes to hold sets of particles.
Upper bound harmonic function (non-zero when feature > mean)
Enforce ideal Helix dihedrals and bonds for a selection at resolution 0.
A class to store a fixed array of same-typed values.
Definition: Array.h:35
def __init__
need to add: ca-ca bond ca-cb is a constraint, no restraint needed ca-ca-ca cb-ca-ca-cb ...
Cosine function.
Definition: Cosine.h:25
Apply a function to the distance between two particles after transforming the second.
Dihedral restraint between four particles.
The type of an atom.
A score on the distance between the surfaces of two spheres.
Return all close unordered pairs of particles taken from the SingletonContainer.
static bool get_is_setup(const IMP::ParticleAdaptor &p)
Definition: rigid_bodies.h:749
Return all spatially-proximal pairs of particles (a,b) from the two SingletonContainers A and B...
GenericHierarchies get_leaves(Hierarchy mhd)
Get all the leaves of the bit of hierarchy.
def get_residue_indexes
Retrieve the residue indexes for the given particle.
Definition: /tools.py:1064
Distance restraint between two particles.
Object used to hold a set of restraints.
Definition: RestraintSet.h:41
Store a list of ParticleIndexPairs.
Miscellaneous utilities.
Definition: /tools.py:1
A well with harmonic barriers.
Definition: HarmonicWell.h:25
Angle restraint between three particles.
ParticleIndexPairs get_indexes(const ParticlePairsTemp &ps)
Get the indexes from a list of particle pairs.
Set up the representation of all proteins and nucleic acid macromolecules.
def get_particle_pairs
Returns the list of connected particles pairs.
The standard decorator for manipulating molecular structures.
Store a list of ParticleIndexes.
Add angular restraint between triplets of consecutive residues/beads to enforce the stereochemistry...
A decorator for a particle representing an atom.
Definition: atom/Atom.h:238
Restrain the dihedral between planes defined by three particles.
Add bonds and improper dihedral restraints for the CBs.
Base class for PMI restraints, which wrap IMP.Restraint(s).
Representation of the system.
def get_sorted_segments
Returns sequence-sorted segments array, each containing the first particle the last particle and the ...
Definition: /tools.py:1742
Score a pair of particles based on the distance between their centers.
RestraintSet * create_elastic_network(const Particles &ps, Float dist_cutoff, Float strength)
Create an elastic network restraint set.
def select
this function uses representation=SimplifiedModel it returns the corresponding selected particles rep...
Definition: /tools.py:727
A decorator for a residue.
Definition: Residue.h:137
Basic functionality that is expected to be used by a wide variety of IMP users.
General purpose algebraic and geometric methods that are expected to be used by a wide variety of IMP...
The general base class for IMP exceptions.
Definition: exception.h:48
Add bond restraint between pair of consecutive residues/beads to enforce the stereochemistry.
Create harmonic restraints between the reference and (transformed) clones.
Applies a PairScore to a Pair.
Definition: PairRestraint.h:31
def add_restraint_to_model
Add a PMI restraint to the model.
Definition: /tools.py:53
Functionality for loading, creating, manipulating and scoring atomic structures.
This class creates a restraint between consecutive TempResidue objects OR an entire PMI MOlecule obje...
This class creates a restraint between the termini two polypeptides, to simulate the sequence connect...
Select hierarchy particles identified by the biological name.
Definition: Selection.h:70
Add dihedral restraints between quadruplet of consecutive residues/beads to enforce the stereochemist...
Applies a PairScore to each Pair in a list.
A repulsive potential on the distance between two atoms.
Perform more efficient close pair finding when rigid bodies are involved.
Inferential scoring building on methods developed as part of the Inferential Structure Determination ...
Harmonic function (symmetric about the mean)
Definition: core/Harmonic.h:27
Restraint a set of residues to use ideal helix dihedrals and bonds.