1 """@namespace IMP.pmi1.restraints.stereochemistry
2 Restraints for keeping correct stereochemistry.
5 from __future__
import print_function
15 from operator
import itemgetter
16 from math
import pi,log,sqrt
22 """ This class creates a restraint between consecutive TempResidue objects OR an entire
23 PMI MOlecule object. """
28 disorderedlength=
False,
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
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
52 raise Exception(
"ConnectivityRestraint: only pass stuff from one Molecule, please")
56 self.m = list(hiers)[0].get_model()
72 SortedSegments.append((start, end, startres))
73 SortedSegments = sorted(SortedSegments, key=itemgetter(2))
76 self.particle_pairs=[]
77 for x
in range(len(SortedSegments) - 1):
79 last = SortedSegments[x][1]
80 first = SortedSegments[x + 1][0]
82 apply_restraint =
True
96 apply_restraint =
False
105 residuegap = firstresn - lastresn - 1
106 if disorderedlength
and (nreslast / 2 + nresfirst / 2 + residuegap) > 20.0:
109 optdist = sqrt(5 / 3) * 1.93 * \
110 (nreslast / 2 + nresfirst / 2 + residuegap) ** 0.6
118 optdist = (0.0 + (float(residuegap) + 1.0) * 3.6) * scale
125 pt0 = last.get_particle()
126 pt1 = first.get_particle()
127 self.particle_pairs.append((pt0,pt1))
130 print(
"Adding sequence connectivity restraint between", pt0.get_name(),
" and ", pt1.get_name(),
'of distance', optdist)
131 self.rs.add_restraint(r)
133 def set_label(self, label):
136 def get_weight(self):
139 def add_to_model(self):
142 def get_restraint(self):
145 def set_weight(self, weight):
147 self.rs.set_weight(weight)
149 def get_output(self):
151 score = self.weight * self.rs.unprotected_evaluate(
None)
152 output[
"_TotalScore"] = str(score)
153 output[
"ConnectivityRestraint_" + self.label] = str(score)
157 """ Returns number of connectivity restraints """
158 return len(self.rs.get_restraints())
161 """ Returns the list of connected particles pairs """
162 return self.particle_pairs
165 return self.weight * self.rs.unprotected_evaluate(
None)
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)
176 included_objects=
None,
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
204 if other_objects
is not None:
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]
219 other_ps = [h.get_particle()
for h
in other_hierarchies]
221 self.mdl = representation.model
224 resolution=resolution,
225 hierarchies=hierarchies)
229 resolution=resolution,
230 hierarchies=other_hierarchies)
232 raise Exception(
"Must pass Representation or included_objects")
255 self.rs.add_restraint(evr)
257 def add_excluded_particle_pairs(self, excluded_particle_pairs):
259 inverted=[(p1,p0)
for p0,p1
in excluded_particle_pairs]
264 self.cpc.add_pair_filter(icpf)
266 def set_label(self, label):
269 def add_to_model(self):
272 def get_restraint(self):
275 def set_weight(self, weight):
277 self.rs.set_weight(weight)
279 def get_output(self):
281 score = self.weight * self.rs.unprotected_evaluate(
None)
282 output[
"_TotalScore"] = str(score)
283 output[
"ExcludedVolumeSphere_" + self.label] = str(score)
287 return self.weight * self.rs.unprotected_evaluate(
None)
290 """Enforce ideal Helix dihedrals and bonds for a selection at resolution 0"""
297 @param hierarchy the root node
298 @param selection_tuple (start, stop, molname, copynum=0)
301 self.mdl = hierarchy.get_model()
305 start = selection_tuple[0]
306 stop = selection_tuple[1]
307 mol = selection_tuple[2]
309 if len(selection_tuple)>3:
310 copy_index = selection_tuple[3]
313 residue_indexes=range(start,stop+1))
314 ps = sel.get_selected_particles(with_representation=
False)
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):
325 def get_weight(self):
328 def add_to_model(self):
331 def get_restraint(self):
334 def set_weight(self, weight):
336 self.rs.set_weight(weight)
338 def get_number_of_bonds(self):
339 return self.r.get_number_of_bonds()
341 def get_number_of_dihedrals(self):
342 return self.r.get_number_of_dihedrals()
345 return self.weight * self.rs.unprotected_evaluate(
None)
347 def get_output(self):
349 score = self.evaluate()
350 output[
"_TotalScore"] = str(score)
351 output[
"HelixRestraint_" + self.label] = str(score)
355 """ Add bond restraint between pair of consecutive
356 residues/beads to enforce the stereochemistry.
360 selection_tuple=
None,
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
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(
382 elif objects
is not None:
384 self.m = particles[0].get_model()
397 (distance - jitter, distance + jitter), strength)
402 raise ValueError(
"wrong length of pair")
405 raise TypeError(
"not a residue")
408 print(
"ResidueBondRestraint: adding a restraint between %s %s" % (pair[0].get_name(), pair[1].get_name()))
409 self.rs.add_restraint(
414 def set_label(self, label):
416 self.rs.set_name(label)
417 for r
in self.rs.get_restraints():
420 def add_to_model(self):
423 def get_restraint(self):
426 def set_weight(self, weight):
428 self.rs.set_weight(weight)
430 def get_excluded_pairs(self):
431 return self.pairslist
433 def get_output(self):
435 score = self.weight * self.rs.unprotected_evaluate(
None)
436 output[
"_TotalScore"] = str(score)
437 output[
"ResidueBondRestraint_" + self.label] = str(score)
442 """Add angular restraint between triplets of consecutive
443 residues/beads to enforce the stereochemistry.
447 selection_tuple=
None,
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(
460 elif objects
is not None:
462 self.m = particles[0].get_model()
470 (pi * anglemin / 180.0,
471 pi * anglemax / 180.0),
477 raise ValueError(
"wrong length of triplet")
480 raise TypeError(
"not a residue")
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(
492 def set_label(self, label):
494 self.rs.set_name(label)
495 for r
in self.rs.get_restraints():
498 def add_to_model(self):
501 def get_restraint(self):
504 def set_weight(self, weight):
506 self.rs.set_weight(weight)
508 def get_excluded_pairs(self):
509 return self.pairslist
511 def get_output(self):
513 score = self.weight * self.rs.unprotected_evaluate(
None)
514 output[
"_TotalScore"] = str(score)
515 output[
"ResidueAngleRestraint_" + self.label] = str(score)
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.
529 selection_tuple=
None,
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(
541 elif objects
is not None:
543 self.m = particles[0].get_model()
550 if stringsequence
is None:
551 stringsequence =
"T" * (len(particles) - 3)
556 raise ValueError(
"wrong length of quadruplet")
559 raise TypeError(
"not a residue")
562 dihedraltype = stringsequence[n]
563 if dihedraltype ==
"C":
567 (pi * anglemin / 180.0,
568 pi * anglemax / 180.0),
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),
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(
587 self.pairslist.append(
589 self.pairslist.append(
592 def set_label(self, label):
594 self.rs.set_name(label)
595 for r
in self.rs.get_restraints():
598 def add_to_model(self):
601 def get_restraint(self):
604 def set_weight(self, weight):
606 self.rs.set_weight(weight)
608 def get_excluded_pairs(self):
609 return self.pairslist
611 def get_output(self):
613 score = self.weight * self.rs.unprotected_evaluate(
None)
614 output[
"_TotalScore"] = str(score)
615 output[
"ResidueDihedralRestraint_" + self.label] = str(score)
619 """Experimental, requires isd_emxl for now"""
629 raise ValueError(
"IMP.isd_emxl is needed")
634 self.particles = IMP.pmi1.tools.select_by_tuple(
638 self.m = representation.prot.get_model()
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
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
657 (bondrslist, anglrslist, diherslist,
658 pairslist) = self.get_CA_force_field()
659 self.pairslist = pairslist
662 self.anglrs.add_restraints(anglrslist)
663 self.dihers.add_restraints(diherslist)
664 self.bondrs.add_restraints(bondrslist)
666 def set_label(self, label):
669 def add_to_model(self):
674 def get_CA_force_field(self):
680 for res
in range(0, len(self.particles) - 1):
682 ps = self.particles[res:res + 2]
685 br = self.get_distance_restraint(ps[0], ps[1], 3.78, 416.0)
686 br.set_name(
'Bond_restraint')
687 bondrslist.append(br)
689 for res
in range(0, len(self.particles) - 4):
694 ps = self.particles[res:res + 5]
697 score_dih] = self.read_potential_dihedral(
698 self.ssstring[res:res + 4],
704 dr = IMP.isd_emxl.CADihedralRestraint(
713 dr.set_name(
'Dihedral restraint')
714 diherslist.append(dr)
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)
722 dr = IMP.isd_emxl.CAAngleRestraint(
728 dr.set_name(
'Angle restraint')
729 anglrslist.append(dr)
730 return (bondrslist, anglrslist, diherslist, pairslist)
732 def read_potential_dihedral(self, string, mix=False):
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)
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]))
756 for i
in range(0, 36):
757 for j
in range(0, 36):
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 * \
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) * \
773 for i
in range(len(counts)):
774 score_dih[i] = -self.kt_caff * self.log(counts[i])
775 return [phi0, phi1, score_dih]
777 def read_potential_angle(self, string, mix=False):
781 for i
in range(0, 180):
782 psi.append(i / 180.0 * pi)
783 score_ang.append(0.0)
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):
791 score_ang[ii] = -self.kt_caff * self.log(float(riga[2]))
796 for i
in range(0, 180):
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):
804 counts[ii] += self.nativeness * float(riga[2])
805 if (len(riga) == 3
and riga[0] ==
"---"):
807 counts[ii] += (1.0 - self.nativeness) * float(riga[2])
809 for i
in range(0, 180):
810 score_ang[i] = -self.kt_caff * self.log(counts[i])
811 return [psi, score_ang]
813 def get_excluded_pairs(self):
814 return self.pairslist
816 def get_restraint(self):
818 tmprs.add_restraint(self.anglrs)
819 tmprs.add_restraint(self.dihers)
820 tmprs.add_restraint(self.bondrs)
823 def get_distance_restraint(self, p0, p1, d0, kappa):
829 def get_output(self):
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)
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)
844 """Add harmonic restraints between all pairs
846 def __init__(self,representation=None,
847 selection_tuples=
None,
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
864 if representation
is None and hierarchy
is not None:
865 self.m = hierarchy.get_model()
866 for st
in selection_tuples:
871 sel =
IMP.atom.Selection(hierarchy,molecule=st[2],residue_indexes=range(st[0],st[1]+1),
872 copy_index=copy_index)
874 sel =
IMP.atom.Selection(hierarchy,molecule=st[2],residue_indexes=range(st[0],st[1]+1),
875 copy_index=copy_index,
877 particles+=sel.get_selected_particles()
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):
886 particles.append(p.get_particle())
888 raise Exception(
"must pass representation or hierarchy")
896 for r
in self.rs.get_restraints():
897 a1,a2 = r.get_inputs()
900 print(
'ElasticNetwork: created',self.rs.get_number_of_restraints(),
'restraints')
902 def set_label(self, label):
904 self.rs.set_name(label)
905 for r
in self.rs.get_restraints():
908 def add_to_model(self):
911 def get_restraint(self):
914 def set_weight(self, weight):
916 self.rs.set_weight(weight)
918 def get_excluded_pairs(self):
919 return self.pairslist
921 def get_output(self):
923 score = self.weight * self.rs.unprotected_evaluate(
None)
924 output[
"_TotalScore"] = str(score)
925 output[
"ElasticNetworkRestraint_" + self.label] = str(score)
930 """ Enable CHARMM force field """
936 enable_nonbonded=
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
952 kB = (1.381 * 6.02214) / 4184.0
953 if representation
is not None:
954 root = representation.prot
956 self.mdl = root.get_model()
958 self.nonbonded_rs =
IMP.RestraintSet(self.mdl, 1.0 / (kB * ff_temp),
'NONBONDED')
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(
973 entire_residues=
True,
974 exclude_backbone=
False)
978 self.ps = limit_to_ps
982 print(
'init bonds score',r.unprotected_evaluate(
None))
983 self.bonds_rs.add_restraint(r)
985 ff.add_well_depths(root)
987 atoms = IMP.atom.get_by_type(root,IMP.atom.ATOM_TYPE)
990 if (zone_ps
is not None)
and zone_nonbonded:
991 print(
'stereochemistry: zone_nonbonded is True')
993 backbone_types=[
'C',
'N',
'CB',
'O']
995 for n
in backbone_types])
996 backbone_atoms = sel.get_selected_particles()
998 sel_ps=IMP.pmi1.topology.get_particles_within_zone(
1002 entire_residues=
True,
1003 exclude_backbone=
True)
1013 self.nbl.add_pair_filter(r.get_full_pair_filter())
1016 self.nonbonded_rs.add_restraint(pr)
1017 print(
'CHARMM is set up')
1019 def set_label(self, label):
1021 self.rs.set_name(label)
1022 for r
in self.rs.get_restraints():
1025 def add_to_model(self):
1029 def get_restraint(self):
1032 def get_close_pair_container(self):
1035 def set_weight(self, weight):
1036 self.weight = weight
1037 self.rs.set_weight(weight)
1039 def get_output(self):
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)
1051 """Add bonds and improper dihedral restraints for the CBs
1054 self, rnums, representation, selection_tuple, strength=10.0, kappa=1.0,
1055 jitter_angle=0.0, jitter_improper=0.0):
1059 ca-cb is a constraint, no restraint needed
1064 self.m = representation.prot.get_model()
1074 ca, cb = self.get_ca_cb(
1075 IMP.pmi1.tools.select_by_tuple(representation,
1076 (rnum, rnum,
'chainA'), resolution=0))
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))
1119 self.rset_angles.add_restraint(ar13u)
1120 self.rset_angles.add_restraint(ar13l)
1126 self.rset_angles.add_restraint(ar23u)
1127 self.rset_angles.add_restraint(ar23l)
1128 if not nter
and not cter:
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)
1157 def get_ca_cb(self, atoms):
1162 ca = a.get_particle()
1164 cb = a.get_particle()
1167 def set_label(self, label):
1169 self.rs.set_name(label)
1170 for r
in self.rs.get_restraints():
1173 def add_to_model(self):
1176 def get_restraint(self):
1179 def set_weight(self, weight):
1180 self.weight = weight
1181 self.rs.set_weight(weight)
1183 def get_excluded_pairs(self):
1184 return self.pairslist
1186 def get_output(self):
1188 score = self.weight * self.rs.unprotected_evaluate(
None)
1189 output[
"_TotalScore"] = str(score)
1190 output[
"PseudoAtomicRestraint_" + self.label] = str(score)
1194 """Create harmonic restraints between the reference and (transformed) clones.
1195 \note Wraps IMP::core::TransformedDistancePairScore with an IMP::core::Harmonic
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
1215 self.mdl = refs[0].get_model()
1219 if len(clones_list)!=len(transforms):
1220 raise Exception(
'Error: There should be as many clones as transforms')
1223 for tmp_clones,trans
in zip(clones_list,transforms):
1225 if len(clones)!=len(refs):
1226 raise Exception(
"Error: len(references)!=len(clones)")
1228 for p0,p1
in zip(refs,clones):
1232 p1.get_particle_index()])
1233 self.rs.add_restraint(r)
1234 print(
'created symmetry network with',self.rs.get_number_of_restraints(),
'restraints')
1236 def set_label(self, label):
1238 self.rs.set_name(label)
1239 for r
in self.rs.get_restraints():
1242 def add_to_model(self):
1245 def get_restraint(self):
1248 def set_weight(self, weight):
1249 self.weight = weight
1250 self.rs.set_weight(weight)
1252 def get_excluded_pairs(self):
1253 return self.pairslist
1255 def get_output(self):
1257 score = self.weight * self.rs.unprotected_evaluate(
None)
1258 output[
"SymmetryRestraint_" + self.label] = str(score)
1259 output[
"_TotalScore"] = str(score)
1264 """ This class creates a restraint between the termini two polypeptides, to simulate the sequence connectivity. """
1269 disorderedlength=
False,
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
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
1293 nter_lastres=ssn[-1][1]
1294 cter_firstres=ssc[0][0]
1295 self.m = nter_lastres.get_model()
1299 optdist = (3.6) * scale
1306 pt0 = nter_lastres.get_particle()
1307 pt1 = cter_firstres.get_particle()
1310 print(
"Adding fusion connectivity restraint between", pt0.get_name(),
" and ", pt1.get_name(),
'of distance', optdist)
1311 self.rs.add_restraint(r)
1313 def set_label(self, label):
1316 def get_weight(self):
1319 def add_to_model(self):
1322 def get_restraint(self):
1325 def set_weight(self, weight):
1326 self.weight = weight
1327 self.rs.set_weight(weight)
1329 def get_output(self):
1331 score = self.weight * self.rs.unprotected_evaluate(
None)
1332 output[
"_TotalScore"] = str(score)
1333 output[
"FusionRestraint_" + self.label] = str(score)
1337 return self.weight * self.rs.unprotected_evaluate(
None)
1344 """Restrain the dihedral between planes defined by three particles.
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.
1350 def __init__(self, particle_triplets, angle=0., k=1., label=None,
1353 @param particle_triplets List of lists of 3 particles. Each triplet
1354 defines a plane. Dihedrals of adjacent planes
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
1363 model = particle_triplets[0][0].get_model()
1364 super(PlaneDihedralRestraint, self).
__init__(model, label=label,
1367 angle = pi * angle / 180.
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(
1375 self.rs.add_restraint(
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)
Enforce CHARMM stereochemistry on the given Hierarchy.
def get_num_restraints
Returns number of connectivity restraints.
A member of a rigid body, it has internal (local) coordinates.
static bool get_is_setup(const IMP::ParticleAdaptor &p)
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.
def __init__
need to add: ca-ca bond ca-cb is a constraint, no restraint needed ca-ca-ca cb-ca-ca-cb ...
Add harmonic restraints between all pairs.
Dihedral restraint between four particles.
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)
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.
Distance restraint between two particles.
Enable CHARMM force field.
Object used to hold a set of restraints.
Store a list of ParticleIndexPairs.
A well with harmonic barriers.
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.
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).
Experimental, requires isd_emxl for now.
Representation of the system.
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.
A decorator for a residue.
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...
def __init__
Setup the CHARMM restraint on a selection.
The general base class for IMP exceptions.
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.
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.
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)
Restraint a set of residues to use ideal helix dihedrals and bonds.