1 """@namespace IMP.pmi.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):
152 score = self.weight * self.rs.unprotected_evaluate(
None)
153 output[
"_TotalScore"] = str(score)
154 output[
"ConnectivityRestraint_" + self.label] = str(score)
158 """ Returns number of connectivity restraints """
159 return len(self.rs.get_restraints())
162 """ Returns the list of connected particles pairs """
163 return self.particle_pairs
166 return self.weight * self.rs.unprotected_evaluate(
None)
169 """A class to create an excluded volume restraint for a set of particles at a given resolution.
170 Can be initialized as a bipartite restraint between two sets of particles.
171 # Potential additional function: Variable resolution for each PMI object. Perhaps passing selection_tuples
172 with (PMI_object, resolution)
177 included_objects=
None,
182 @param representation DEPRECATED - just pass objects
183 @param included_objects Can be one of the following inputs:
184 IMP Hierarchy, PMI System/State/Molecule/TempResidue, or a list/set of them
185 @param other_objects Initializes a bipartite restraint between included_objects and other_objects
186 Same format as included_objects
187 @param resolution The resolution particles at which to impose the restraint.
188 By default, the coarsest particles will be chosen.
189 If a number is chosen, for each particle, the closest
190 resolution will be used (see IMP.atom.Selection).
191 @param kappa Restraint strength
205 if other_objects
is not None:
214 if representation
is None:
215 if hierarchies
is None:
216 raise Exception(
"Must at least pass included objects")
217 self.mdl = hierarchies[0].get_model()
218 included_ps = [h.get_particle()
for h
in hierarchies]
220 other_ps = [h.get_particle()
for h
in other_hierarchies]
222 self.mdl = representation.m
225 resolution=resolution,
226 hierarchies=hierarchies)
230 resolution=resolution,
231 hierarchies=other_hierarchies)
233 raise Exception(
"Must pass Representation or included_objects")
256 self.rs.add_restraint(evr)
258 def add_excluded_particle_pairs(self, excluded_particle_pairs):
260 inverted=[(p1,p0)
for p0,p1
in excluded_particle_pairs]
265 self.cpc.add_pair_filter(icpf)
267 def set_label(self, label):
270 def add_to_model(self):
273 def get_restraint(self):
276 def set_weight(self, weight):
278 self.rs.set_weight(weight)
280 def get_output(self):
283 score = self.weight * self.rs.unprotected_evaluate(
None)
284 output[
"_TotalScore"] = str(score)
285 output[
"ExcludedVolumeSphere_" + self.label] = str(score)
289 return self.weight * self.rs.unprotected_evaluate(
None)
292 """Enforce ideal Helix dihedrals and bonds for a selection at resolution 0"""
299 @param hierarchy the root node
300 @param selection_tuple (start, stop, molname, copynum=0)
303 self.mdl = hierarchy.get_model()
307 start = selection_tuple[0]
308 stop = selection_tuple[1]
309 mol = selection_tuple[2]
311 if len(selection_tuple)>3:
312 copy_index = selection_tuple[3]
315 residue_indexes=range(start,stop+1))
316 ps = sel.get_selected_particles(with_representation=
False)
320 self.rs.add_restraint(self.r)
321 print(
'Created helix %s.%i.%i-%i with %i dihedrals and %i bonds'%(
322 mol,copy_index,start,stop,
323 self.get_number_of_bonds(),self.get_number_of_dihedrals()))
324 def set_label(self, label):
327 def get_weight(self):
330 def add_to_model(self):
333 def get_restraint(self):
336 def set_weight(self, weight):
338 self.rs.set_weight(weight)
340 def get_number_of_bonds(self):
341 return self.r.get_number_of_bonds()
343 def get_number_of_dihedrals(self):
344 return self.r.get_number_of_dihedrals()
347 return self.weight * self.rs.unprotected_evaluate(
None)
349 def get_output(self):
352 score = self.evaluate()
353 output[
"_TotalScore"] = str(score)
354 output[
"HelixRestraint_" + self.label] = str(score)
358 """ Add bond restraint between pair of consecutive
359 residues/beads to enforce the stereochemistry.
363 selection_tuple=
None,
369 @param representation (PMI1)
370 @param selection_tuple Requested selection (PMI1)
371 @param objects (PMI2)
372 @param distance Resting distance for restraint
373 @param strength Bond constant
374 @param jitter Defines the +- added to the optimal distance in the harmonic well restraint
375 used to increase the tolerance
378 if representation
is not None and selection_tuple
is not None:
379 self.m = representation.prot.get_model()
380 particles = IMP.pmi.tools.select_by_tuple(
385 elif objects
is not None:
387 self.m = particles[0].get_model()
400 (distance - jitter, distance + jitter), strength)
405 raise ValueError(
"wrong length of pair")
408 raise TypeError(
"not a residue")
411 print(
"ResidueBondRestraint: adding a restraint between %s %s" % (pair[0].get_name(), pair[1].get_name()))
412 self.rs.add_restraint(
417 def set_label(self, label):
419 self.rs.set_name(label)
420 for r
in self.rs.get_restraints():
423 def add_to_model(self):
426 def get_restraint(self):
429 def set_weight(self, weight):
431 self.rs.set_weight(weight)
433 def get_excluded_pairs(self):
434 return self.pairslist
436 def get_output(self):
439 score = self.weight * self.rs.unprotected_evaluate(
None)
440 output[
"_TotalScore"] = str(score)
441 output[
"ResidueBondRestraint_" + self.label] = str(score)
446 """Add angular restraint between triplets of consecutive
447 residues/beads to enforce the stereochemistry.
451 selection_tuple=
None,
457 if representation
is not None and selection_tuple
is not None:
458 self.m = representation.prot.get_model()
459 particles = IMP.pmi.tools.select_by_tuple(
464 elif objects
is not None:
466 self.m = particles[0].get_model()
474 (pi * anglemin / 180.0,
475 pi * anglemax / 180.0),
481 raise ValueError(
"wrong length of triplet")
484 raise TypeError(
"not a residue")
487 print(
"ResidueAngleRestraint: adding a restraint between %s %s %s" % (triplet[0].get_name(), triplet[1].get_name(), triplet[2].get_name()))
488 self.rs.add_restraint(
496 def set_label(self, label):
498 self.rs.set_name(label)
499 for r
in self.rs.get_restraints():
502 def add_to_model(self):
505 def get_restraint(self):
508 def set_weight(self, weight):
510 self.rs.set_weight(weight)
512 def get_excluded_pairs(self):
513 return self.pairslist
515 def get_output(self):
518 score = self.weight * self.rs.unprotected_evaluate(
None)
519 output[
"_TotalScore"] = str(score)
520 output[
"ResidueAngleRestraint_" + self.label] = str(score)
525 """Add dihedral restraints between quadruplet of consecutive
526 residues/beads to enforce the stereochemistry.
527 Give as input a string of "C" and "T", meaning cys (0+-40) or trans (180+-40)
528 dihedral. The length of the string must be \#residue-3.
529 Without the string, the dihedral will be assumed trans.
534 selection_tuple=
None,
539 if representation
is not None and selection_tuple
is not None:
540 self.m = representation.prot.get_model()
541 particles = IMP.pmi.tools.select_by_tuple(
546 elif objects
is not None:
548 self.m = particles[0].get_model()
555 if stringsequence
is None:
556 stringsequence =
"T" * (len(particles) - 3)
561 raise ValueError(
"wrong length of quadruplet")
564 raise TypeError(
"not a residue")
567 dihedraltype = stringsequence[n]
568 if dihedraltype ==
"C":
572 (pi * anglemin / 180.0,
573 pi * anglemax / 180.0),
575 print(
"ResidueDihedralRestraint: adding a CYS restraint between %s %s %s %s" % (quadruplet[0].get_name(), quadruplet[1].get_name(),
576 quadruplet[2].get_name(), quadruplet[3].get_name()))
577 if dihedraltype ==
"T":
578 anglemin = 180 - 70.0
579 anglemax = 180 + 70.0
581 (pi * anglemin / 180.0,
582 pi * anglemax / 180.0),
584 print(
"ResidueDihedralRestraint: adding a TRANS restraint between %s %s %s %s" % (quadruplet[0].get_name(), quadruplet[1].get_name(),
585 quadruplet[2].get_name(), quadruplet[3].get_name()))
586 self.rs.add_restraint(
592 self.pairslist.append(
594 self.pairslist.append(
597 def set_label(self, label):
599 self.rs.set_name(label)
600 for r
in self.rs.get_restraints():
603 def add_to_model(self):
606 def get_restraint(self):
609 def set_weight(self, weight):
611 self.rs.set_weight(weight)
613 def get_excluded_pairs(self):
614 return self.pairslist
616 def get_output(self):
619 score = self.weight * self.rs.unprotected_evaluate(
None)
620 output[
"_TotalScore"] = str(score)
621 output[
"ResidueDihedralRestraint_" + self.label] = str(score)
625 """Experimental, requires isd_emxl for now"""
635 raise ValueError(
"IMP.isd_emxl is needed")
640 self.particles = IMP.pmi.tools.select_by_tuple(
644 self.m = representation.prot.get_model()
648 self.anglfilename = IMP.isd_emxl.get_data_path(
"CAAngleRestraint.dat")
649 self.dihefilename = IMP.isd_emxl.get_data_path(
650 "CADihedralRestraint.dat")
651 self.nativeness = nativeness
652 self.kt_caff = kt_caff
658 if len(self.particles) != len(ssstring):
659 print(len(self.particles), len(ssstring))
660 print(
"SecondaryStructure: residue range and SS string incompatible")
661 self.ssstring = ssstring
663 (bondrslist, anglrslist, diherslist,
664 pairslist) = self.get_CA_force_field()
665 self.pairslist = pairslist
668 self.anglrs.add_restraints(anglrslist)
669 self.dihers.add_restraints(diherslist)
670 self.bondrs.add_restraints(bondrslist)
672 def set_label(self, label):
675 def add_to_model(self):
680 def get_CA_force_field(self):
686 for res
in range(0, len(self.particles) - 1):
688 ps = self.particles[res:res + 2]
691 br = self.get_distance_restraint(ps[0], ps[1], 3.78, 416.0)
692 br.set_name(
'Bond_restraint')
693 bondrslist.append(br)
695 for res
in range(0, len(self.particles) - 4):
700 ps = self.particles[res:res + 5]
703 score_dih] = self.read_potential_dihedral(
704 self.ssstring[res:res + 4],
710 dr = IMP.isd_emxl.CADihedralRestraint(
719 dr.set_name(
'Dihedral restraint')
720 diherslist.append(dr)
722 for res
in range(0, len(self.particles) - 2):
723 ps = self.particles[res:res + 3]
724 [psi, score_ang] = self.read_potential_angle(
725 self.ssstring[res:res + 2],
True)
728 dr = IMP.isd_emxl.CAAngleRestraint(
734 dr.set_name(
'Angle restraint')
735 anglrslist.append(dr)
736 return (bondrslist, anglrslist, diherslist, pairslist)
738 def read_potential_dihedral(self, string, mix=False):
743 for i
in range(0, 36):
744 phi0.append(i * 10.0 / 180.0 * pi)
745 phi1.append(i * 10.0 / 180.0 * pi)
746 for j
in range(0, 36):
747 score_dih.append(0.0)
750 f = open(self.dihefilename,
'r')
751 for line
in f.readlines():
752 riga = (line.strip()).split()
753 if (len(riga) == 4
and riga[0] == string):
754 ii = int(float(riga[1]) / 10.0)
755 jj = int(float(riga[2]) / 10.0)
756 score_dih[ii * len(phi0) + jj] = - \
757 self.kt_caff * self.log(float(riga[3]))
762 for i
in range(0, 36):
763 for j
in range(0, 36):
765 f = open(self.dihefilename,
'r')
766 for line
in f.readlines():
767 riga = (line.strip()).split()
768 if (len(riga) == 4
and riga[0] == string):
769 ii = int(float(riga[1]) / 10.0)
770 jj = int(float(riga[2]) / 10.0)
771 counts[ii * len(phi0) + jj] += self.nativeness * \
773 if (len(riga) == 4
and riga[0] ==
"-----"):
774 ii = int(float(riga[1]) / 10.0)
775 jj = int(float(riga[2]) / 10.0)
776 counts[ii * len(phi0) + jj] += (1.0 - self.nativeness) * \
779 for i
in range(len(counts)):
780 score_dih[i] = -self.kt_caff * self.log(counts[i])
781 return [phi0, phi1, score_dih]
783 def read_potential_angle(self, string, mix=False):
787 for i
in range(0, 180):
788 psi.append(i / 180.0 * pi)
789 score_ang.append(0.0)
792 f = open(self.anglfilename,
'r')
793 for line
in f.readlines():
794 riga = (line.strip()).split()
795 if (len(riga) == 3
and riga[0] == string):
797 score_ang[ii] = -self.kt_caff * self.log(float(riga[2]))
802 for i
in range(0, 180):
805 f = open(self.anglfilename,
'r')
806 for line
in f.readlines():
807 riga = (line.strip()).split()
808 if (len(riga) == 3
and riga[0] == string):
810 counts[ii] += self.nativeness * float(riga[2])
811 if (len(riga) == 3
and riga[0] ==
"---"):
813 counts[ii] += (1.0 - self.nativeness) * float(riga[2])
815 for i
in range(0, 180):
816 score_ang[i] = -self.kt_caff * self.log(counts[i])
817 return [psi, score_ang]
819 def get_excluded_pairs(self):
820 return self.pairslist
822 def get_restraint(self):
824 tmprs.add_restraint(self.anglrs)
825 tmprs.add_restraint(self.dihers)
826 tmprs.add_restraint(self.bondrs)
829 def get_distance_restraint(self, p0, p1, d0, kappa):
835 def get_output(self):
838 score_angle = self.anglrs.unprotected_evaluate(
None)
839 score_dihers = self.dihers.unprotected_evaluate(
None)
840 score_bondrs = self.bondrs.unprotected_evaluate(
None)
841 output[
"_TotalScore"] = str(score_angle + score_dihers + score_bondrs)
843 output[
"SecondaryStructure_Angles_" + self.label] = str(score_angle)
844 output[
"SecondaryStructure_Dihedrals_" +
845 self.label] = str(score_dihers)
846 output[
"SecondaryStructure_Bonds_" + self.label] = str(score_bondrs)
851 """Add harmonic restraints between all pairs
853 def __init__(self,representation=None,
854 selection_tuples=
None,
861 @param representation Representation object
862 @param selection_tuples Selecting regions for the restraint [[start,stop,molname,copy_index=0],...]
863 @param resolution Resolution for applying restraint
864 @param strength Bond strength
865 @param dist_cutoff Cutoff for making restraints
866 @param ca_only Selects only CAlphas. Only matters if resolution=0.
867 @param hierarchy Root hierarchy to select from, use this instead of representation
871 if representation
is None and hierarchy
is not None:
872 self.m = hierarchy.get_model()
873 for st
in selection_tuples:
878 sel =
IMP.atom.Selection(hierarchy,molecule=st[2],residue_indexes=range(st[0],st[1]+1),
879 copy_index=copy_index)
881 sel =
IMP.atom.Selection(hierarchy,molecule=st[2],residue_indexes=range(st[0],st[1]+1),
882 copy_index=copy_index,
884 particles+=sel.get_selected_particles()
886 self.m = representation.mdl
887 for st
in selection_tuples:
888 print(
'selecting with',st)
889 for p
in IMP.pmi.tools.select_by_tuple(representation,st,resolution=resolution):
893 particles.append(p.get_particle())
895 raise Exception(
"must pass representation or hierarchy")
903 for r
in self.rs.get_restraints():
904 a1,a2 = r.get_inputs()
907 print(
'ElasticNetwork: created',self.rs.get_number_of_restraints(),
'restraints')
909 def set_label(self, label):
911 self.rs.set_name(label)
912 for r
in self.rs.get_restraints():
915 def add_to_model(self):
918 def get_restraint(self):
921 def set_weight(self, weight):
923 self.rs.set_weight(weight)
925 def get_excluded_pairs(self):
926 return self.pairslist
928 def get_output(self):
931 score = self.weight * self.rs.unprotected_evaluate(
None)
932 output[
"_TotalScore"] = str(score)
933 output[
"ElasticNetworkRestraint_" + self.label] = str(score)
938 """ Enable CHARMM force field """
944 enable_nonbonded=
True,
946 zone_nonbonded=
False,
947 representation=
None):
948 """Setup the CHARMM restraint on a selection. Expecting atoms.
949 @param root The node at which to apply the restraint
950 @param ff_temp The temperature of the force field
951 @param zone_ps Create a zone around this set of particles
952 Automatically includes the entire residue (incl. backbone)
953 @param zone_size The size for looking for neighbor residues
954 @param enable_nonbonded Allow the repulsive restraint
955 @param enable_bonded Allow the bonded restraint
956 @param zone_nonbonded EXPERIMENTAL: exclude from nonbonded all sidechains that aren't in zone!
957 @param representation Legacy representation object
960 kB = (1.381 * 6.02214) / 4184.0
961 if representation
is not None:
962 root = representation.prot
964 self.mdl = root.get_model()
966 self.nonbonded_rs =
IMP.RestraintSet(self.mdl, 1.0 / (kB * ff_temp),
'NONBONDED')
973 topology = ff.create_topology(root)
974 topology.apply_default_patches()
975 topology.setup_hierarchy(root)
976 if zone_ps
is not None:
977 limit_to_ps=IMP.pmi.topology.get_particles_within_zone(
981 entire_residues=
True,
982 exclude_backbone=
False)
986 self.ps = limit_to_ps
990 print(
'init bonds score',r.unprotected_evaluate(
None))
991 self.bonds_rs.add_restraint(r)
993 ff.add_well_depths(root)
995 atoms = IMP.atom.get_by_type(root,IMP.atom.ATOM_TYPE)
998 if (zone_ps
is not None)
and zone_nonbonded:
999 print(
'stereochemistry: zone_nonbonded is True')
1001 backbone_types=[
'C',
'N',
'CB',
'O']
1003 for n
in backbone_types])
1004 backbone_atoms = sel.get_selected_particles()
1006 sel_ps=IMP.pmi.topology.get_particles_within_zone(
1010 entire_residues=
True,
1011 exclude_backbone=
True)
1021 self.nbl.add_pair_filter(r.get_full_pair_filter())
1024 self.nonbonded_rs.add_restraint(pr)
1025 print(
'CHARMM is set up')
1027 def set_label(self, label):
1029 self.rs.set_name(label)
1030 for r
in self.rs.get_restraints():
1033 def add_to_model(self):
1037 def get_restraint(self):
1040 def get_close_pair_container(self):
1043 def set_weight(self, weight):
1044 self.weight = weight
1045 self.rs.set_weight(weight)
1047 def get_output(self):
1050 bonds_score = self.weight * self.bonds_rs.unprotected_evaluate(
None)
1051 nonbonded_score = self.weight * self.nonbonded_rs.unprotected_evaluate(
None)
1052 score=bonds_score+nonbonded_score
1053 output[
"_TotalScore"] = str(score)
1054 output[
"CHARMM_BONDS"] = str(bonds_score)
1055 output[
"CHARMM_NONBONDED"] = str(nonbonded_score)
1060 """Add bonds and improper dihedral restraints for the CBs
1063 self, rnums, representation, selection_tuple, strength=10.0, kappa=1.0,
1064 jitter_angle=0.0, jitter_improper=0.0):
1068 ca-cb is a constraint, no restraint needed
1073 self.m = representation.prot.get_model()
1083 ca, cb = self.get_ca_cb(
1084 IMP.pmi.tools.select_by_tuple(representation,
1085 (rnum, rnum,
'chainA'), resolution=0))
1089 ca_prev, cb_prev = self.get_ca_cb(
1090 IMP.pmi.tools.select_by_tuple(representation,
1091 (rnum - 1, rnum - 1,
'chainA'), resolution=0))
1092 ca_next, cb_next = self.get_ca_cb(
1093 IMP.pmi.tools.select_by_tuple(representation,
1094 (rnum + 1, rnum + 1,
'chainA'), resolution=0))
1128 self.rset_angles.add_restraint(ar13u)
1129 self.rset_angles.add_restraint(ar13l)
1135 self.rset_angles.add_restraint(ar23u)
1136 self.rset_angles.add_restraint(ar23l)
1137 if not nter
and not cter:
1161 self.rset_angles.add_restraint(idru)
1162 self.rset_angles.add_restraint(idrl)
1163 self.rs.add_restraint(self.rset_bonds)
1164 self.rs.add_restraint(self.rset_angles)
1166 def get_ca_cb(self, atoms):
1171 ca = a.get_particle()
1173 cb = a.get_particle()
1176 def set_label(self, label):
1178 self.rs.set_name(label)
1179 for r
in self.rs.get_restraints():
1182 def add_to_model(self):
1185 def get_restraint(self):
1188 def set_weight(self, weight):
1189 self.weight = weight
1190 self.rs.set_weight(weight)
1192 def get_excluded_pairs(self):
1193 return self.pairslist
1195 def get_output(self):
1198 score = self.weight * self.rs.unprotected_evaluate(
None)
1199 output[
"_TotalScore"] = str(score)
1200 output[
"PseudoAtomicRestraint_" + self.label] = str(score)
1204 """Create harmonic restraints between the reference and (transformed) clones.
1205 \note Wraps IMP::core::TransformedDistancePairScore with an IMP::core::Harmonic
1214 @param references List of particles for symmetry reference
1215 @param clones_list List of lists of clone particles
1216 @param transforms Transforms moving each selection to the first selection
1217 @param label Label for output
1218 @param strength The elastic bond strength
1219 \note You will have to perform an IMP::atom::Selection to get the particles you need.
1220 Will check to make sure the numbers match.
1222 self.mdl = root.get_model()
1226 if len(clones_list)!=len(transforms):
1227 raise Exception(
'Error: There should be as many clones as transforms')
1230 for clones,trans
in zip(clones_list,transforms):
1231 if len(clones)!=len(references):
1232 raise Exception(
"Error: len(references)!=len(clones)")
1234 for p0,p1
in zip(references,clones):
1236 self.rs.add_restraint(r)
1238 print(
'created symmetry network with',self.rs.get_number_of_restraints(),
'restraints')
1240 def set_label(self, label):
1242 self.rs.set_name(label)
1243 for r
in self.rs.get_restraints():
1246 def add_to_model(self):
1249 def get_restraint(self):
1252 def set_weight(self, weight):
1253 self.weight = weight
1254 self.rs.set_weight(weight)
1256 def get_excluded_pairs(self):
1257 return self.pairslist
1259 def get_output(self):
1262 score = self.weight * self.rs.unprotected_evaluate(
None)
1263 output[
"SymmetryRestraint_" + self.label] = str(score)
1264 output[
"_TotalScore"] = str(score)
1269 """ This class creates a restraint between the termini two polypeptides, to simulate the sequence connectivity. """
1274 disorderedlength=
False,
1279 @param nterminal - single PMI2 Hierarchy/molecule at the nterminal
1280 @param cterminal - single PMI2 Hierarchy/molecule at the cterminal
1281 @param scale Scale the maximal distance between the beads by this factor when disorderedlength is False.
1282 The maximal distance is calculated as ((float(residuegap) + 1.0) * 3.6) * scale.
1283 @param disorderedlength - This flag uses either disordered length
1284 calculated for random coil peptides (True) or zero
1285 surface-to-surface distance between beads (False)
1286 as optimal distance for the sequence connectivity
1288 @param upperharmonic - This flag uses either harmonic (False)
1289 or upperharmonic (True) in the intra-pair
1290 connectivity restraint.
1291 @param resolution - The resolution to connect things at - only used if you pass PMI objects
1292 @param label - A string to identify this restraint in the output/stat file
1298 nter_lastres=ssn[-1][1]
1299 cter_firstres=ssc[0][0]
1300 self.m = nter_lastres.get_model()
1304 optdist = (3.6) * scale
1311 pt0 = nter_lastres.get_particle()
1312 pt1 = cter_firstres.get_particle()
1315 print(
"Adding fusion connectivity restraint between", pt0.get_name(),
" and ", pt1.get_name(),
'of distance', optdist)
1316 self.rs.add_restraint(r)
1318 def set_label(self, label):
1321 def get_weight(self):
1324 def add_to_model(self):
1327 def get_restraint(self):
1330 def set_weight(self, weight):
1331 self.weight = weight
1332 self.rs.set_weight(weight)
1334 def get_output(self):
1337 score = self.weight * self.rs.unprotected_evaluate(
None)
1338 output[
"_TotalScore"] = str(score)
1339 output[
"FusionRestraint_" + self.label] = str(score)
1343 return self.weight * self.rs.unprotected_evaluate(
None)
1350 """Restrain the dihedral between planes defined by three particles.
1352 This restraint is useful for restraining the twist of a string of
1353 more or less identical rigid bodies, so long as the curvature is mild.
1356 def __init__(self, particle_triplets, angle=0., k=1., label=None,
1359 @param particle_triplets List of lists of 3 particles. Each triplet
1360 defines a plane. Dihedrals of adjacent planes
1362 @param angle Angle of plane dihedral in degrees
1363 @param k Strength of restraint
1364 @param label Label for output
1365 @param weight Weight of restraint
1366 \note Particles defining planes should be rigid and more or less
1367 parallel for proper behavior
1369 m = particle_triplets[0][0].get_model()
1370 super(PlaneDihedralRestraint, self).
__init__(m, label=label,
1373 angle = pi * angle / 180.
1375 for i, t1
in enumerate(particle_triplets[:-1]):
1376 t2 = particle_triplets[i + 1]
1377 q1 = [t1[1], t1[0], t2[0], t2[1]]
1378 q2 = [t1[2], t1[0], t2[0], t2[2]]
A filter which returns true if a container containers the Pair.
Add dihedral restraints between quadruplet of consecutive residues/beads to enforce the stereochemist...
CHARMMParameters * get_heavy_atom_CHARMM_parameters()
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.
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)
Enable CHARMM force field.
A class to store an fixed array of same-typed values.
Enforce ideal Helix dihedrals and bonds for a selection at resolution 0.
This class creates a restraint between the termini two polypeptides, to simulate the sequence connect...
Add harmonic restraints between all pairs.
Restrain the dihedral between planes defined by three particles.
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-proximals 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.
Representation of the system.
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)
Add bond restraint between pair of consecutive residues/beads to enforce the stereochemistry.
The standard decorator for manipulating molecular structures.
RestraintSet * create_elastic_network(const Particles &ps, Float dist_cutoff, Float strength)
Create an elastic network restraint set.
Add bonds and improper dihedral restraints for the CBs.
Store a list of ParticleIndexes.
A decorator for a particle representing an atom.
def get_particle_pairs
Returns the list of connected particles pairs.
def __init__
Setup the CHARMM restraint on a selection.
Create harmonic restraints between the reference and (transformed) clones.
Add angular restraint between triplets of consecutive residues/beads to enforce the stereochemistry...
def get_num_restraints
Returns number of connectivity restraints.
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...
This class creates a restraint between consecutive TempResidue objects OR an entire PMI MOlecule obje...
A class to create an excluded volume restraint for a set of particles at a given resolution.
The general base class for IMP exceptions.
def __init__
need to add: ca-ca bond ca-cb is a constraint, no restraint needed ca-ca-ca cb-ca-ca-cb ...
Applies a PairScore to a Pair.
Experimental, requires isd_emxl for now.
Set up the representation of all proteins and nucleic acid macromolecules.
Functionality for loading, creating, manipulating and scoring atomic structures.
Select hierarchy particles identified by the biological name.
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.
Base class for PMI restraints, which wrap IMP.Restraint(s).
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.