1 """@namespace IMP.pmi1.mmcif
2 @brief Support for the mmCIF file format.
4 IMP has basic support for writing out files in mmCIF format, for
5 deposition in [PDB-Dev](https://pdb-dev.wwpdb.org/).
6 mmCIF files are currently generated by creating an
7 IMP.pmi1.mmcif.ProtocolOutput class, and attaching it to an
8 IMP.pmi1.representation.Representation object, after which any
9 generated models and metadata are collected and output as mmCIF.
12 from __future__
import print_function
41 import ihm.representation
43 import ihm.cross_linkers
45 def _assign_id(obj, seen_objs, obj_by_id):
46 """Assign a unique ID to obj, and track all ids in obj_by_id."""
47 if obj
not in seen_objs:
48 if not hasattr(obj,
'id'):
50 obj.id = len(obj_by_id)
51 seen_objs[obj] = obj.id
53 obj.id = seen_objs[obj]
56 def _get_by_residue(p):
57 """Determine whether the given particle represents a specific residue
58 or a more coarse-grained object."""
62 class _ComponentMapper(object):
63 """Map a Particle to a component name"""
64 def __init__(self, prot):
67 self.name =
'cif-output'
68 self.o.dictionary_pdbs[self.name] = self.prot
69 self.o._init_dictchain(self.name, self.prot, multichar_chain=
True)
71 def __getitem__(self, p):
72 protname, is_a_bead = self.o.get_prot_name_from_particle(self.name, p)
76 class _AsymMapper(object):
77 """Map a Particle to an asym_unit"""
78 def __init__(self, simo, prot):
80 self._cm = _ComponentMapper(prot)
82 def __getitem__(self, p):
83 protname = self._cm[p]
84 return self.simo.asym_units[protname]
86 def get_feature(self, ps):
87 """Get an ihm.restraint.Feature that covers the given particles"""
95 rngs.append(asym(rind, rind))
99 rngs.append(asym(rinds[0], rinds[-1]))
101 raise ValueError(
"Unsupported particle type %s" % str(p))
102 return ihm.restraint.ResidueFeature(rngs)
105 class _AllSoftware(object):
106 def __init__(self, system):
108 self.modeller_used = self.phyre2_used =
False
110 name=
"Integrative Modeling Platform (IMP)",
111 version=IMP.__version__,
112 classification=
"integrative model building",
113 description=
"integrative model building",
114 location=
'https://integrativemodeling.org')
116 name=
"IMP PMI module",
117 version=IMP.pmi1.__version__,
118 classification=
"integrative model building",
119 description=
"integrative model building",
120 location=
'https://integrativemodeling.org')
121 self.system.software.extend([imp, pmi])
124 if hasattr(imp,
'citation'):
125 if sys.version_info[0] > 2:
127 javi =
'Vel\u00e1zquez-Muriel J'
129 javi =
'Velazquez-Muriel J'
130 imp.citation = ihm.Citation(
132 title=
'Putting the pieces together: integrative modeling '
133 'platform software for structure determination of '
134 'macromolecular assemblies',
135 journal=
'PLoS Biol', volume=10, page_range=
'e1001244',
137 authors=[
'Russel D',
'Lasker K',
'Webb B', javi,
'Tjioe E',
138 'Schneidman-Duhovny D',
'Peterson B',
'Sali A'],
139 doi=
'10.1371/journal.pbio.1001244')
140 pmi.citation = ihm.Citation(
142 title=
'Modeling Biological Complexes Using Integrative '
143 'Modeling Platform.',
144 journal=
'Methods Mol Biol', volume=2022, page_range=(353, 377),
146 authors=[
'Saltzberg D',
'Greenberg CH',
'Viswanath S',
147 'Chemmama I',
'Webb B',
'Pellarin R',
'Echeverria I',
149 doi=
'10.1007/978-1-4939-9608-7_15')
151 def set_modeller_used(self, version, date):
152 if self.modeller_used:
154 self.modeller_used =
True
156 name=
'MODELLER', classification=
'comparative modeling',
157 description=
'Comparative modeling by satisfaction '
158 'of spatial restraints, build ' + date,
159 location=
'https://salilab.org/modeller/', version=version)
160 self.system.software.append(s)
161 if hasattr(s,
'citation'):
162 s.citation = ihm.Citation(
164 title=
'Comparative protein modelling by satisfaction of '
165 'spatial restraints.',
166 journal=
'J Mol Biol', volume=234, page_range=(779, 815),
167 year=1993, authors=[
'Sali A',
'Blundell TL'],
168 doi=
'10.1006/jmbi.1993.1626')
170 def set_phyre2_used(self):
173 self.phyre2_used =
True
175 name=
'Phyre2', classification=
'protein homology modeling',
176 description=
'Protein Homology/analogY Recognition Engine V 2.0',
177 version=
'2.0', location=
'http://www.sbg.bio.ic.ac.uk/~phyre2/')
178 if hasattr(s,
'citation'):
179 s.citation = ihm.Citation(
181 title=
'The Phyre2 web portal for protein modeling, '
182 'prediction and analysis.',
183 journal=
'Nat Protoc', volume=10, page_range=(845, 858),
184 authors=[
'Kelley LA',
'Mezulis S',
'Yates CM',
'Wass MN',
186 year=2015, doi=
'10.1038/nprot.2015.053')
187 self.system.software.append(s)
189 def _get_fragment_is_rigid(fragment):
190 """Determine whether a fragment is modeled rigidly"""
195 class _PDBFragment(ihm.representation.ResidueSegment):
196 """Record details about part of a PDB file used as input
198 def __init__(self, state, component, start, end, pdb_offset,
199 pdbname, chain, hier, asym_unit):
201 super(_PDBFragment, self).__init__(
202 asym_unit=asym_unit.pmi_range(start, end),
203 rigid=
None, primitive=
'sphere')
204 self.component, self.start, self.end, self.offset, self.pdbname \
205 = component, start, end, pdb_offset, pdbname
206 self.state, self.chain, self.hier = state, chain, hier
211 rigid = property(
lambda self: _get_fragment_is_rigid(self),
212 lambda self, val:
None)
214 def combine(self, other):
218 class _BeadsFragment(ihm.representation.FeatureSegment):
219 """Record details about beads used to represent part of a component."""
221 def __init__(self, state, component, start, end, count, hier, asym_unit):
222 super(_BeadsFragment, self).__init__(asym_unit=asym_unit(start, end),
223 rigid=
None, primitive=
'sphere', count=count)
224 self.state, self.component, self.hier \
225 = state, component, hier
227 rigid = property(
lambda self: _get_fragment_is_rigid(self),
228 lambda self, val:
None)
230 def combine(self, other):
232 if type(other) == type(self)
and \
233 other.asym_unit.seq_id_range[0] == self.asym_unit.seq_id_range[1] + 1:
234 self.asym_unit.seq_id_range = (self.asym_unit.seq_id_range[0],
235 other.asym_unit.seq_id_range[1])
236 self.count += other.count
240 class _AllModelRepresentations(object):
241 def __init__(self, simo):
245 self.fragments = OrderedDict()
246 self._all_representations = {}
248 def copy_component(self, state, name, original, asym_unit):
249 """Copy all representation for `original` in `state` to `name`"""
252 newf.asym_unit = asym_unit(*f.asym_unit.seq_id_range)
254 for rep
in self.fragments:
255 if original
in self.fragments[rep]:
256 if name
not in self.fragments[rep]:
257 self.fragments[rep][name] = OrderedDict()
258 self.fragments[rep][name][state] = [copy_frag(f)
259 for f
in self.fragments[rep][original][state]]
262 first_state = list(self.fragments[rep][name].keys())[0]
263 if state
is first_state:
264 representation = self._all_representations[rep]
265 representation.extend(self.fragments[rep][name][state])
267 def add_fragment(self, state, representation, fragment):
268 """Add a model fragment."""
269 comp = fragment.component
270 id_rep = id(representation)
271 self._all_representations[id_rep] = representation
272 if id_rep
not in self.fragments:
273 self.fragments[id_rep] = OrderedDict()
274 if comp
not in self.fragments[id_rep]:
275 self.fragments[id_rep][comp] = OrderedDict()
276 if state
not in self.fragments[id_rep][comp]:
277 self.fragments[id_rep][comp][state] = []
278 fragments = self.fragments[id_rep][comp][state]
279 if len(fragments) == 0
or not fragments[-1].combine(fragment):
280 fragments.append(fragment)
283 first_state = list(self.fragments[id_rep][comp].keys())[0]
284 if state
is first_state:
285 representation.append(fragment)
288 class _AllDatasets(object):
289 """Track all datasets generated by PMI and add them to the ihm.System"""
290 def __init__(self, system):
293 self._datasets_by_state = {}
294 self._restraints_by_state = {}
296 def get_all_group(self, state):
297 """Get a DatasetGroup encompassing all datasets so far in this state"""
301 g = ihm.dataset.DatasetGroup(self._datasets_by_state.get(state, [])
302 + [r.dataset
for r
in self._restraints_by_state.get(state, [])
306 def add(self, state, dataset):
307 """Add a new dataset."""
308 self._datasets.append(dataset)
309 if state
not in self._datasets_by_state:
310 self._datasets_by_state[state] = []
311 self._datasets_by_state[state].append(dataset)
314 self.system.orphan_datasets.append(dataset)
318 """Add the dataset for a restraint"""
319 if state
not in self._restraints_by_state:
320 self._restraints_by_state[state] = []
321 self._restraints_by_state[state].append(restraint)
324 class _CrossLinkRestraint(ihm.restraint.CrossLinkRestraint):
325 """Restrain to a set of cross links"""
328 _label_map = {
'wtDSS':
'DSS',
'scDSS':
'DSS',
'scEDC':
'EDC'}
329 _descriptor_map = {
'DSS': ihm.cross_linkers.dss,
330 'EDC': ihm.cross_linkers.edc}
332 def __init__(self, pmi_restraint):
333 self.pmi_restraint = pmi_restraint
334 label = self.pmi_restraint.label
336 label = self._label_map.get(label, label)
338 super(_CrossLinkRestraint, self).__init__(
339 dataset=self.pmi_restraint.dataset,
340 linker=self._get_chem_descriptor(label))
343 def _get_chem_descriptor(cls, label):
344 if label
not in cls._descriptor_map:
348 d = ihm.ChemDescriptor(label)
349 cls._descriptor_map[label] = d
350 return cls._descriptor_map[label]
352 def _set_psi_sigma(self, model):
355 if model.m != self.pmi_restraint.m:
357 for resolution
in self.pmi_restraint.sigma_dictionary:
358 statname =
'ISDCrossLinkMS_Sigma_%s_%s' % (resolution, self.label)
359 if model.stats
and statname
in model.stats:
360 sigma = float(model.stats[statname])
361 p = self.pmi_restraint.sigma_dictionary[resolution][0]
362 old_values.append((p, p.get_scale()))
364 for psiindex
in self.pmi_restraint.psi_dictionary:
365 statname =
'ISDCrossLinkMS_Psi_%s_%s' % (psiindex, self.label)
366 if model.stats
and statname
in model.stats:
367 psi = float(model.stats[statname])
368 p = self.pmi_restraint.psi_dictionary[psiindex][0]
369 old_values.append((p, p.get_scale()))
373 return list(reversed(old_values))
375 def add_fits_from_model_statfile(self, model):
377 old_values = self._set_psi_sigma(model)
381 for xl
in self.cross_links:
383 xl.fits[model] = ihm.restraint.CrossLinkFit(
384 psi=xl.psi, sigma1=xl.sigma1, sigma2=xl.sigma2)
386 for p, scale
in old_values:
390 def __set_dataset(self, val):
391 self.pmi_restraint.dataset = val
392 dataset = property(
lambda self: self.pmi_restraint.dataset,
396 def get_asym_mapper_for_state(simo, state, asym_map):
397 asym = asym_map.get(state,
None)
399 asym = _AsymMapper(simo, state.prot)
400 asym_map[state] = asym
403 class _PMICrossLink(object):
406 psi = property(
lambda self: self.psi_p.get_scale(),
407 lambda self, val:
None)
408 sigma1 = property(
lambda self: self.sigma1_p.get_scale(),
409 lambda self, val:
None)
410 sigma2 = property(
lambda self: self.sigma2_p.get_scale(),
411 lambda self, val:
None)
414 class _ResidueCrossLink(ihm.restraint.ResidueCrossLink, _PMICrossLink):
418 class _FeatureCrossLink(ihm.restraint.FeatureCrossLink, _PMICrossLink):
422 class _EM2DRestraint(ihm.restraint.EM2DRestraint):
423 def __init__(self, state, pmi_restraint, image_number, resolution,
424 pixel_size, image_resolution, projection_number,
426 self.pmi_restraint, self.image_number = pmi_restraint, image_number
427 super(_EM2DRestraint, self).__init__(
428 dataset=pmi_restraint.datasets[image_number],
429 assembly=state.modeled_assembly,
430 segment=
False, number_raw_micrographs=micrographs_number,
431 pixel_size_width=pixel_size, pixel_size_height=pixel_size,
432 image_resolution=image_resolution,
433 number_of_projections=projection_number)
436 def __get_dataset(self):
437 return self.pmi_restraint.datasets[self.image_number]
438 def __set_dataset(self, val):
439 self.pmi_restraint.datasets[self.image_number] = val
440 dataset = property(__get_dataset, __set_dataset)
442 def add_fits_from_model_statfile(self, model):
443 ccc = self._get_cross_correlation(model)
444 transform = self._get_transformation(model)
445 rot = transform.get_rotation()
446 rm = [[e
for e
in rot.get_rotation_matrix_row(i)]
for i
in range(3)]
447 self.fits[model] = ihm.restraint.EM2DRestraintFit(
448 cross_correlation_coefficient=ccc,
450 tr_vector=transform.get_translation())
452 def _get_transformation(self, model):
453 """Get the transformation that places the model on the image"""
454 stats = model.em2d_stats
or model.stats
455 prefix =
'ElectronMicroscopy2D_%s_Image%d' % (self.pmi_restraint.label,
456 self.image_number + 1)
457 r = [float(stats[prefix +
'_Rotation%d' % i])
for i
in range(4)]
458 t = [float(stats[prefix +
'_Translation%d' % i])
462 inv = model.transform.get_inverse()
466 def _get_cross_correlation(self, model):
467 """Get the cross correlation coefficient between the model projection
469 stats = model.em2d_stats
or model.stats
470 return float(stats[
'ElectronMicroscopy2D_%s_Image%d_CCC'
471 % (self.pmi_restraint.label,
472 self.image_number + 1)])
475 class _EM3DRestraint(ihm.restraint.EM3DRestraint):
477 def __init__(self, simo, state, pmi_restraint, target_ps, densities):
478 self.pmi_restraint = pmi_restraint
479 super(_EM3DRestraint, self).__init__(
480 dataset=pmi_restraint.dataset,
481 assembly=self._get_assembly(densities, simo, state),
482 fitting_method=
'Gaussian mixture models',
483 number_of_gaussians=len(target_ps))
486 def __set_dataset(self, val):
487 self.pmi_restraint.dataset = val
488 dataset = property(
lambda self: self.pmi_restraint.dataset,
491 def _get_assembly(self, densities, simo, state):
492 """Get the Assembly that this restraint acts on"""
493 cm = _ComponentMapper(state.prot)
496 components[cm[d]] =
None
497 a = simo._get_subassembly(components,
498 name=
"EM subassembly",
499 description=
"All components that fit the EM map")
502 def add_fits_from_model_statfile(self, model):
503 ccc = self._get_cross_correlation(model)
504 self.fits[model] = ihm.restraint.EM3DRestraintFit(
505 cross_correlation_coefficient=ccc)
507 def _get_cross_correlation(self, model):
508 """Get the cross correlation coefficient between the model
510 if model.stats
is not None:
511 return float(model.stats[
'GaussianEMRestraint_%s_CCC'
512 % self.pmi_restraint.label])
515 class _GeometricRestraint(ihm.restraint.GeometricRestraint):
517 def __init__(self, simo, state, pmi_restraint, geometric_object,
518 feature, distance, sigma):
519 self.pmi_restraint = pmi_restraint
520 super(_GeometricRestraint, self).__init__(
521 dataset=pmi_restraint.dataset,
522 geometric_object=geometric_object, feature=feature,
523 distance=distance, harmonic_force_constant = 1. / sigma,
527 def __set_dataset(self, val):
528 self.pmi_restraint.dataset = val
529 dataset = property(
lambda self: self.pmi_restraint.dataset,
533 class _ReplicaExchangeProtocolStep(ihm.protocol.Step):
534 def __init__(self, state, rex):
535 if rex.monte_carlo_sample_objects
is not None:
536 method =
'Replica exchange monte carlo'
538 method =
'Replica exchange molecular dynamics'
539 super(_ReplicaExchangeProtocolStep, self).__init__(
540 assembly=state.modeled_assembly,
542 method=method, name=
'Sampling',
543 num_models_begin=
None,
544 num_models_end=rex.vars[
"number_of_frames"],
545 multi_scale=
True, multi_state=
False, ordered=
False,
549 class _SimpleProtocolStep(ihm.protocol.Step):
550 def __init__(self, state, num_models_end, method):
551 super(_SimpleProtocolStep, self).__init__(
552 assembly=state.modeled_assembly,
554 method=method, name=
'Sampling',
555 num_models_begin=
None,
556 num_models_end=num_models_end,
557 multi_scale=
True, multi_state=
False, ordered=
False,
561 class _Chain(object):
562 """Represent a single chain in a Model"""
563 def __init__(self, pmi_chain_id, asym_unit):
564 self.pmi_chain_id, self.asym_unit = pmi_chain_id, asym_unit
568 def add(self, xyz, atom_type, residue_type, residue_index,
569 all_indexes, radius):
570 if atom_type
is None:
571 self.spheres.append((xyz, residue_type, residue_index,
572 all_indexes, radius))
574 self.atoms.append((xyz, atom_type, residue_type, residue_index,
575 all_indexes, radius))
576 orig_comp = property(
lambda self: self.comp)
578 class _TransformedChain(object):
579 """Represent a chain that is a transformed version of another"""
580 def __init__(self, orig_chain, asym_unit, transform):
581 self.orig_chain, self.asym_unit = orig_chain, asym_unit
582 self.transform = transform
584 def __get_spheres(self):
585 for (xyz, residue_type, residue_index, all_indexes,
586 radius)
in self.orig_chain.spheres:
587 yield (self.transform * xyz, residue_type, residue_index,
589 spheres = property(__get_spheres)
591 def __get_atoms(self):
592 for (xyz, atom_type, residue_type, residue_index, all_indexes,
593 radius)
in self.orig_chain.atoms:
594 yield (self.transform * xyz, atom_type, residue_type, residue_index,
596 atoms = property(__get_atoms)
598 entity = property(
lambda self: self.orig_chain.entity)
599 orig_comp = property(
lambda self: self.orig_chain.comp)
601 class _Excluder(object):
602 def __init__(self, component, simo):
603 self._seqranges = simo._exclude_coords.get(component, [])
605 def is_excluded(self, indexes):
606 """Return True iff the given sequence range is excluded."""
607 for seqrange
in self._seqranges:
608 if indexes[0] >= seqrange[0]
and indexes[-1] <= seqrange[1]:
612 class _Model(ihm.model.Model):
613 def __init__(self, prot, simo, protocol, assembly, representation):
614 super(_Model, self).__init__(assembly=assembly, protocol=protocol,
615 representation=representation)
616 self.simo = weakref.proxy(simo)
622 self.em2d_stats =
None
625 self._is_restrained =
True
628 self.m = prot.get_model()
629 o.dictionary_pdbs[name] = prot
630 o._init_dictchain(name, prot, multichar_chain=
True)
631 (particle_infos_for_pdb,
632 self.geometric_center) = o.get_particle_infos_for_pdb_writing(name)
634 self._make_spheres_atoms(particle_infos_for_pdb, o, name, simo)
637 def all_chains(self, simo):
638 """Yield all chains, including transformed ones"""
640 for c
in self.chains:
642 chain_for_comp[c.comp] = c
643 for tc
in simo._transformed_components:
644 orig_chain = chain_for_comp.get(tc.original,
None)
646 asym = simo.asym_units[tc.name]
647 c = _TransformedChain(orig_chain, asym, tc.transform)
651 def _make_spheres_atoms(self, particle_infos_for_pdb, o, name, simo):
652 entity_for_chain = {}
655 for protname, chain_id
in o.dictchain[name].items():
656 entity_for_chain[chain_id] = simo.entities[protname]
657 comp_for_chain[chain_id] = protname
661 correct_asym[chain_id] = simo.asym_units[protname]
668 for (xyz, atom_type, residue_type, chain_id, residue_index,
669 all_indexes, radius)
in particle_infos_for_pdb:
670 if chain
is None or chain.pmi_chain_id != chain_id:
671 chain = _Chain(chain_id, correct_asym[chain_id])
672 chain.entity = entity_for_chain[chain_id]
673 chain.comp = comp_for_chain[chain_id]
674 self.chains.append(chain)
675 excluder = _Excluder(chain.comp, simo)
676 if not excluder.is_excluded(all_indexes
if all_indexes
677 else [residue_index]):
678 chain.add(xyz, atom_type, residue_type, residue_index,
681 def parse_rmsf_file(self, fname, component):
682 self.rmsf[component] = rmsf = {}
683 with open(fname)
as fh:
685 resnum, blocknum, val = line.split()
686 rmsf[int(resnum)] = (int(blocknum), float(val))
688 def get_rmsf(self, component, indexes):
689 """Get the RMSF value for the given residue indexes."""
692 rmsf = self.rmsf[component]
693 blocknums = dict.fromkeys(rmsf[ind][0]
for ind
in indexes)
694 if len(blocknums) != 1:
695 raise ValueError(
"Residue indexes %s aren't all in the same block"
697 return rmsf[indexes[0]][1]
700 for chain
in self.all_chains(self.simo):
701 pmi_offset = chain.asym_unit.entity.pmi_offset
702 for atom
in chain.atoms:
703 (xyz, atom_type, residue_type, residue_index,
704 all_indexes, radius) = atom
705 pt = self.transform * xyz
706 yield ihm.model.Atom(asym_unit=chain.asym_unit,
707 seq_id=residue_index - pmi_offset,
708 atom_id=atom_type.get_string(),
710 x=pt[0], y=pt[1], z=pt[2])
712 def get_spheres(self):
713 for chain
in self.all_chains(self.simo):
714 pmi_offset = chain.asym_unit.entity.pmi_offset
715 for sphere
in chain.spheres:
716 (xyz, residue_type, residue_index,
717 all_indexes, radius) = sphere
718 if all_indexes
is None:
719 all_indexes = (residue_index,)
720 pt = self.transform * xyz
721 yield ihm.model.Sphere(asym_unit=chain.asym_unit,
722 seq_id_range=(all_indexes[0] - pmi_offset,
723 all_indexes[-1] - pmi_offset),
724 x=pt[0], y=pt[1], z=pt[2], radius=radius,
725 rmsf=self.get_rmsf(chain.orig_comp, all_indexes))
728 class _AllProtocols(object):
729 def __init__(self, simo):
730 self.simo = weakref.proxy(simo)
732 self.protocols = OrderedDict()
734 def add_protocol(self, state):
735 """Add a new Protocol"""
736 if state
not in self.protocols:
737 self.protocols[state] = []
738 p = ihm.protocol.Protocol()
739 self.simo.system.orphan_protocols.append(p)
740 self.protocols[state].append(p)
742 def add_step(self, step, state):
743 """Add a ProtocolStep to the last Protocol of the given State"""
744 if state
not in self.protocols:
745 self.add_protocol(state)
746 protocol = self.get_last_protocol(state)
747 if len(protocol.steps) == 0:
748 step.num_models_begin = 0
750 step.num_models_begin = protocol.steps[-1].num_models_end
751 protocol.steps.append(step)
752 step.id = len(protocol.steps)
754 step.dataset_group = self.simo.all_datasets.get_all_group(state)
756 def add_postproc(self, step, state):
757 """Add a postprocessing step to the last protocol"""
758 protocol = self.get_last_protocol(state)
759 if len(protocol.analyses) == 0:
760 protocol.analyses.append(ihm.analysis.Analysis())
761 protocol.analyses[-1].steps.append(step)
763 def get_last_protocol(self, state):
764 """Return the most recently-added _Protocol"""
765 return self.protocols[state][-1]
768 class _AllStartingModels(object):
769 def __init__(self, simo):
773 self.models = OrderedDict()
776 def add_pdb_fragment(self, fragment):
777 """Add a starting model PDB fragment."""
778 comp = fragment.component
779 state = fragment.state
780 if comp
not in self.models:
781 self.models[comp] = OrderedDict()
782 if state
not in self.models[comp]:
783 self.models[comp][state] = []
784 models = self.models[comp][state]
785 if len(models) == 0 \
786 or models[-1].fragments[0].pdbname != fragment.pdbname:
787 model = self._add_model(fragment)
791 models[-1].fragments.append(weakref.proxy(fragment))
793 pmi_offset = models[-1].asym_unit.entity.pmi_offset
794 sid_begin = min(fragment.start + fragment.offset - pmi_offset,
795 models[-1].asym_unit.seq_id_range[0])
796 sid_end = max(fragment.end + fragment.offset - pmi_offset,
797 models[-1].asym_unit.seq_id_range[1])
798 models[-1].asym_unit = fragment.asym_unit.asym(sid_begin, sid_end)
799 fragment.starting_model = models[-1]
801 def _add_model(self, f):
802 parser = ihm.metadata.PDBParser()
803 r = parser.parse_file(f.pdbname)
805 self.simo._add_dataset(r[
'dataset'])
807 templates = r[
'templates'].get(f.chain, [])
810 self.simo.system.locations.append(t.alignment_file)
812 self.simo._add_dataset(t.dataset)
813 pmi_offset = f.asym_unit.entity.pmi_offset
815 asym_unit=f.asym_unit.asym.pmi_range(f.start, f.end),
816 dataset=r[
'dataset'], asym_id=f.chain,
817 templates=templates, offset=f.offset + pmi_offset,
818 metadata=r[
'metadata'],
819 software=r[
'software'][0]
if r[
'software']
else None,
820 script_file=r[
'script'])
821 m.fragments = [weakref.proxy(f)]
825 class _StartingModel(ihm.startmodel.StartingModel):
826 def get_seq_dif(self):
830 pmi_offset = self.asym_unit.entity.pmi_offset
831 mh = IMP.mmcif.data._StartingModelAtomHandler(self.templates,
833 for f
in self.fragments:
835 residue_indexes=list(range(f.start - f.offset,
836 f.end - f.offset + 1)))
837 for a
in mh.get_ihm_atoms(sel.get_selected_particles(),
838 f.offset - pmi_offset):
840 self._seq_dif = mh._seq_dif
843 class _ReplicaExchangeAnalysisPostProcess(ihm.analysis.ClusterStep):
844 """Post processing using AnalysisReplicaExchange0 macro"""
846 def __init__(self, rex, num_models_begin):
849 for fname
in self.get_all_stat_files():
850 with open(fname)
as fh:
851 num_models_end += len(fh.readlines())
852 super(_ReplicaExchangeAnalysisPostProcess, self).__init__(
853 feature=
'RMSD', num_models_begin=num_models_begin,
854 num_models_end=num_models_end)
856 def get_stat_file(self, cluster_num):
857 return os.path.join(self.rex._outputdir,
"cluster.%d" % cluster_num,
860 def get_all_stat_files(self):
861 for i
in range(self.rex._number_of_clusters):
862 yield self.get_stat_file(i)
865 class _ReplicaExchangeAnalysisEnsemble(ihm.model.Ensemble):
866 """Ensemble generated using AnalysisReplicaExchange0 macro"""
868 num_models_deposited =
None
870 def __init__(self, pp, cluster_num, model_group, num_deposit):
871 with open(pp.get_stat_file(cluster_num))
as fh:
872 num_models = len(fh.readlines())
873 super(_ReplicaExchangeAnalysisEnsemble, self).__init__(
874 num_models=num_models,
875 model_group=model_group, post_process=pp,
876 clustering_feature=pp.feature,
877 name=model_group.name)
878 self.cluster_num = cluster_num
879 self.num_models_deposited = num_deposit
881 def get_rmsf_file(self, component):
882 return os.path.join(self.post_process.rex._outputdir,
883 'cluster.%d' % self.cluster_num,
884 'rmsf.%s.dat' % component)
886 def load_rmsf(self, model, component):
887 fname = self.get_rmsf_file(component)
888 if os.path.exists(fname):
889 model.parse_rmsf_file(fname, component)
891 def get_localization_density_file(self, fname):
892 return os.path.join(self.post_process.rex._outputdir,
893 'cluster.%d' % self.cluster_num,
896 def load_localization_density(self, state, fname, select_tuple, asym_units):
897 fullpath = self.get_localization_density_file(fname)
898 if os.path.exists(fullpath):
899 details =
"Localization density for %s %s" \
900 % (fname, self.model_group.name)
901 local_file = ihm.location.OutputFileLocation(fullpath,
903 for s
in select_tuple:
904 if isinstance(s, tuple)
and len(s) == 3:
905 asym = asym_units[s[2]].pmi_range(s[0], s[1])
908 den = ihm.model.LocalizationDensity(file=local_file,
910 self.densities.append(den)
912 def load_all_models(self, simo, state):
913 stat_fname = self.post_process.get_stat_file(self.cluster_num)
915 with open(stat_fname)
as fh:
916 stats = ast.literal_eval(fh.readline())
918 rmf_file = os.path.join(os.path.dirname(stat_fname),
919 "%d.rmf3" % model_num)
920 for c
in state.all_modeled_components:
922 state._pmi_object.set_coordinates_from_rmf(c, rmf_file, 0,
923 force_rigid_update=
True)
927 if model_num >= self.num_models_deposited:
931 def _get_precision(self):
932 precfile = os.path.join(self.post_process.rex._outputdir,
933 "precision.%d.%d.out" % (self.cluster_num,
935 if not os.path.exists(precfile):
938 r = re.compile(
'All .*/cluster.%d/ average centroid distance ([\d\.]+)'
940 with open(precfile)
as fh:
944 return float(m.group(1))
946 precision = property(
lambda self: self._get_precision(),
947 lambda self, val:
None)
949 class _SimpleEnsemble(ihm.model.Ensemble):
950 """Simple manually-created ensemble"""
952 num_models_deposited =
None
954 def __init__(self, pp, model_group, num_models, drmsd,
955 num_models_deposited, ensemble_file):
956 super(_SimpleEnsemble, self).__init__(
957 model_group=model_group, post_process=pp, num_models=num_models,
958 file=ensemble_file, precision=drmsd, name=model_group.name,
959 clustering_feature=
'dRMSD')
960 self.num_models_deposited = num_models_deposited
962 def load_localization_density(self, state, component, asym, local_file):
963 den = ihm.model.LocalizationDensity(file=local_file, asym_unit=asym)
964 self.densities.append(den)
967 class _EntityMapper(dict):
968 """Handle mapping from IMP components to CIF entities.
969 Multiple components may map to the same entity if they share sequence."""
970 def __init__(self, system):
971 super(_EntityMapper, self).__init__()
972 self._sequence_dict = {}
976 def add(self, component_name, sequence, offset):
977 def entity_seq(sequence):
980 return [
'UNK' if s ==
'X' else s
for s
in sequence]
983 if sequence
not in self._sequence_dict:
986 d = component_name.split(
"@")[0].split(
".")[0]
987 entity = Entity(entity_seq(sequence), description=d,
989 self.system.entities.append(entity)
990 self._sequence_dict[sequence] = entity
991 self[component_name] = self._sequence_dict[sequence]
994 class _TransformedComponent(object):
995 def __init__(self, name, original, transform):
996 self.name, self.original, self.transform = name, original, transform
999 class _SimpleRef(object):
1000 """Class with similar interface to weakref.ref, but keeps a strong ref"""
1001 def __init__(self, ref):
1007 class _State(ihm.model.State):
1008 """Representation of a single state in the system."""
1010 def __init__(self, pmi_object, po):
1015 self._pmi_object = weakref.proxy(pmi_object)
1016 if hasattr(pmi_object,
'state'):
1019 self._pmi_state = _SimpleRef(pmi_object.state)
1021 self._pmi_state = weakref.ref(pmi_object)
1023 old_name = self.name
1024 super(_State, self).__init__(experiment_type=
'Fraction of bulk')
1025 self.name = old_name
1029 self.modeled_assembly = ihm.Assembly(
1030 name=
"Modeled assembly",
1031 description=self.get_postfixed_name(
1032 "All components modeled by IMP"))
1033 po.system.orphan_assemblies.append(self.modeled_assembly)
1035 self.all_modeled_components = []
1038 return hash(self._pmi_state())
1039 def __eq__(self, other):
1040 return self._pmi_state() == other._pmi_state()
1042 def add_model_group(self, group):
1046 def get_prefixed_name(self, name):
1047 """Prefix the given name with the state name, if available."""
1049 return self.short_name +
' ' + name
1053 return name[0].upper() + name[1:]
if name
else ''
1055 def get_postfixed_name(self, name):
1056 """Postfix the given name with the state name, if available."""
1058 return "%s in state %s" % (name, self.short_name)
1062 short_name = property(
lambda self: self._pmi_state().short_name)
1063 long_name = property(
lambda self: self._pmi_state().long_name)
1064 def __get_name(self):
1065 return self._pmi_state().long_name
1066 def __set_name(self, val):
1067 self._pmi_state().long_name = val
1068 name = property(__get_name, __set_name)
1072 """A single entity in the system.
1073 This functions identically to the base ihm.Entity class, but it
1074 allows identifying residues by either the IHM numbering scheme
1075 (seq_id, which is always contiguous starting at 1) or by the PMI
1076 scheme (which is similar, but may not start at 1 if the sequence in
1077 the FASTA file has one or more N-terminal gaps"""
1078 def __init__(self, sequence, pmi_offset, *args, **keys):
1081 self.pmi_offset = pmi_offset
1082 super(Entity, self).__init__(sequence, *args, **keys)
1085 """Return a single IHM residue indexed using PMI numbering"""
1086 return self.residue(res_id - self.pmi_offset)
1089 """Return a range of IHM residues indexed using PMI numbering"""
1090 off = self.pmi_offset
1091 return self(res_id_begin - off, res_id_end - off)
1095 """A single asymmetric unit in the system.
1096 This functions identically to the base ihm.AsymUnit class, but it
1097 allows identifying residues by either the IHM numbering scheme
1098 (seq_id, which is always contiguous starting at 1) or by the PMI
1099 scheme (which is similar, but may not start at 1 if the sequence in
1100 the FASTA file has one or more N-terminal gaps"""
1102 def __init__(self, entity, *args, **keys):
1103 super(AsymUnit, self).__init__(
1104 entity, auth_seq_id_map=entity.pmi_offset, *args, **keys)
1107 """Return a single IHM residue indexed using PMI numbering"""
1108 return self.residue(res_id - self.entity.pmi_offset)
1111 """Return a range of IHM residues indexed using PMI numbering"""
1112 off = self.entity.pmi_offset
1113 return self(res_id_begin - off, res_id_end - off)
1117 """Class to encode a modeling protocol as mmCIF.
1119 IMP has basic support for writing out files in mmCIF format, for
1120 deposition in [PDB-Dev](https://pdb-dev.wwpdb.org/).
1121 After creating an instance of this class, attach it to an
1122 IMP.pmi1.representation.Representation object. After this, any
1123 generated models and metadata are automatically collected and
1126 def __init__(self, fh):
1128 self.system = ihm.System()
1129 self._state_group = ihm.model.StateGroup()
1130 self.system.state_groups.append(self._state_group)
1133 self._state_ensemble_offset = 0
1134 self._each_metadata = []
1135 self._file_datasets = []
1136 self._main_script = os.path.abspath(sys.argv[0])
1139 loc = ihm.location.WorkflowFileLocation(path=self._main_script,
1140 details=
"The main integrative modeling script")
1141 self.system.locations.append(loc)
1144 self.__asym_states = {}
1145 self._working_directory = os.getcwd()
1147 "Default representation")
1148 self.entities = _EntityMapper(self.system)
1150 self.asym_units = {}
1151 self._all_components = {}
1152 self.all_modeled_components = []
1153 self._transformed_components = []
1154 self.sequence_dict = {}
1157 self._xy_plane = ihm.geometry.XYPlane()
1158 self._xz_plane = ihm.geometry.XZPlane()
1159 self._z_axis = ihm.geometry.ZAxis()
1160 self._center_origin = ihm.geometry.Center(0,0,0)
1161 self._identity_transform = ihm.geometry.Transformation.identity()
1164 self._exclude_coords = {}
1166 self.all_representations = _AllModelRepresentations(self)
1167 self.all_protocols = _AllProtocols(self)
1168 self.all_datasets = _AllDatasets(self.system)
1169 self.all_starting_models = _AllStartingModels(self)
1171 self.all_software = _AllSoftware(self.system)
1174 """Create a new Representation and return it. This can be
1175 passed to add_model(), add_bead_element() or add_pdb_element()."""
1176 r = ihm.representation.Representation()
1177 self.system.orphan_representations.append(r)
1181 """Don't record coordinates for the given domain.
1182 Coordinates for the given domain (specified by a component name
1183 and a 2-element tuple giving the start and end residue numbers)
1184 will be excluded from the mmCIF file. This can be used to exclude
1185 parts of the structure that weren't well resolved in modeling.
1186 Any bead or residue that lies wholly within this range will be
1187 excluded. Multiple ranges for a given component can be excluded
1188 by calling this method multiple times."""
1189 if component
not in self._exclude_coords:
1190 self._exclude_coords[component] = []
1191 self._exclude_coords[component].append(seqrange)
1193 def _is_excluded(self, component, start, end):
1194 """Return True iff this chunk of sequence should be excluded"""
1195 for seqrange
in self._exclude_coords.get(component, ()):
1196 if start >= seqrange[0]
and end <= seqrange[1]:
1199 def _add_state(self, state):
1200 """Create a new state and return a pointer to it."""
1201 self._state_ensemble_offset = len(self.system.ensembles)
1202 s = _State(state, self)
1203 self._state_group.append(s)
1204 self._last_state = s
1207 def get_file_dataset(self, fname):
1208 for d
in self._file_datasets:
1209 fd = d.get(os.path.abspath(fname),
None)
1213 def _get_chain_for_component(self, name, output):
1214 """Get the chain ID for a component, if any."""
1216 if name
in self.asym_units:
1217 return self.asym_units[name]._id
1222 def _get_assembly_comps(self, assembly):
1223 """Get the names of the components in the given assembly"""
1227 comps[ca.details] =
None
1231 """Make a new component that's a transformed copy of another.
1232 All representation for the existing component is copied to the
1234 assembly_comps = self._get_assembly_comps(state.modeled_assembly)
1235 if name
in assembly_comps:
1236 raise ValueError(
"Component %s already exists" % name)
1237 elif original
not in assembly_comps:
1238 raise ValueError(
"Original component %s does not exist" % original)
1239 self.create_component(state, name,
True)
1240 self.add_component_sequence(state, name, self.sequence_dict[original])
1241 self._transformed_components.append(_TransformedComponent(
1242 name, original, transform))
1243 self.all_representations.copy_component(state, name, original,
1244 self.asym_units[name])
1246 def create_component(self, state, name, modeled):
1247 new_comp = name
not in self._all_components
1248 self._all_components[name] =
None
1250 state.all_modeled_components.append(name)
1252 self.asym_units[name] =
None
1253 self.all_modeled_components.append(name)
1255 def add_component_sequence(self, state, name, seq):
1256 def get_offset(seq):
1258 for i
in range(len(seq)):
1261 seq, offset = get_offset(seq)
1262 if name
in self.sequence_dict:
1263 if self.sequence_dict[name] != seq:
1264 raise ValueError(
"Sequence mismatch for component %s" % name)
1266 self.sequence_dict[name] = seq
1267 self.entities.add(name, seq, offset)
1268 if name
in self.asym_units:
1269 if self.asym_units[name]
is None:
1271 entity = self.entities[name]
1272 asym =
AsymUnit(entity, details=name)
1273 self.system.asym_units.append(asym)
1274 self.asym_units[name] = asym
1275 state.modeled_assembly.append(self.asym_units[name])
1277 def _add_restraint_model_fits(self):
1278 """Add fits to restraints for all known models"""
1279 for group, m
in self.system._all_models():
1280 if m._is_restrained:
1281 for r
in self.system.restraints:
1282 if hasattr(r,
'add_fits_from_model_statfile'):
1283 r.add_fits_from_model_statfile(m)
1285 def _add_em2d_raw_micrographs(self):
1286 """If the deprecated metadata.EMMicrographsDataset class was used,
1287 transfer its data to the EM2D restraint."""
1288 for r
in self.system.restraints:
1289 if isinstance(r, _EM2DRestraint):
1290 for d
in r.dataset.parents:
1291 if isinstance(d, IMP.pmi1.metadata.EMMicrographsDataset):
1292 r.number_raw_micrographs = d.number
1293 if isinstance(d.location,
1294 IMP.pmi1.metadata.FileLocation):
1295 d.location.content_type = \
1296 ihm.location.InputFileLocation.content_type
1300 ihm.dumper.write(self.fh, [self.system])
1303 """Do any final processing on the class hierarchy.
1304 After calling this method, the `system` member (an instance
1305 of `ihm.System`) completely reproduces the PMI modeling, and
1306 can be written out to an mmCIF file with `ihm.dumper.write`,
1307 and/or modified using the ihm API."""
1308 self._add_restraint_model_fits()
1310 self._add_em2d_raw_micrographs()
1313 self.system.software.extend(m
for m
in self._metadata
1314 if isinstance(m, ihm.Software))
1315 self.system.citations.extend(m
for m
in self._metadata
1316 if isinstance(m, ihm.Citation))
1317 self.system.locations.extend(m
for m
in self._metadata
1318 if isinstance(m, ihm.location.FileLocation))
1321 all_repos = [m
for m
in self._metadata
1322 if isinstance(m, ihm.location.Repository)]
1323 self.system.update_locations_in_repositories(all_repos)
1325 def add_pdb_element(self, state, name, start, end, offset, pdbname,
1326 chain, hier, representation=
None):
1327 if self._is_excluded(name, start, end):
1329 if representation
is None:
1330 representation = self.default_representation
1331 asym = self.asym_units[name]
1332 p = _PDBFragment(state, name, start, end, offset, pdbname, chain,
1334 self.all_representations.add_fragment(state, representation, p)
1335 self.all_starting_models.add_pdb_fragment(p)
1337 def add_bead_element(self, state, name, start, end, num, hier,
1338 representation=
None):
1339 if self._is_excluded(name, start, end):
1341 if representation
is None:
1342 representation = self.default_representation
1343 asym = self.asym_units[name]
1344 pmi_offset = asym.entity.pmi_offset
1345 b = _BeadsFragment(state, name, start - pmi_offset, end - pmi_offset,
1347 self.all_representations.add_fragment(state, representation, b)
1349 def get_cross_link_group(self, pmi_restraint):
1350 r = _CrossLinkRestraint(pmi_restraint)
1351 self.system.restraints.append(r)
1352 self._add_restraint_dataset(r)
1355 def add_experimental_cross_link(self, r1, c1, r2, c2, rsr):
1356 if c1
not in self._all_components
or c2
not in self._all_components:
1362 e1 = self.entities[c1]
1363 e2 = self.entities[c2]
1364 xl = ihm.restraint.ExperimentalCrossLink(residue1=e1.pmi_residue(r1),
1365 residue2=e2.pmi_residue(r2))
1366 rsr.experimental_cross_links.append([xl])
1369 def add_cross_link(self, state, ex_xl, p1, p2, length, sigma1_p, sigma2_p,
1372 asym = get_asym_mapper_for_state(self, state, self.__asym_states)
1373 d = ihm.restraint.UpperBoundDistanceRestraint(length)
1375 if _get_by_residue(p1)
and _get_by_residue(p2):
1376 cls = _ResidueCrossLink
1378 cls = _FeatureCrossLink
1379 xl = cls(ex_xl, asym1=asym[p1], asym2=asym[p2], distance=d,
1382 xl.psi_p, xl.sigma1_p, xl.sigma2_p = psi_p, sigma1_p, sigma2_p
1383 rsr.cross_links.append(xl)
1385 def add_replica_exchange(self, state, rex):
1390 self.all_protocols.add_step(_ReplicaExchangeProtocolStep(state, rex),
1393 def _add_simple_dynamics(self, num_models_end, method):
1395 state = self._last_state
1396 self.all_protocols.add_step(_SimpleProtocolStep(state, num_models_end,
1399 def _add_protocol(self):
1401 state = self._last_state
1402 self.all_protocols.add_protocol(state)
1404 def _add_dataset(self, dataset):
1405 return self.all_datasets.add(self._last_state, dataset)
1407 def _add_restraint_dataset(self, restraint):
1408 return self.all_datasets.add_restraint(self._last_state, restraint)
1410 def _add_simple_postprocessing(self, num_models_begin, num_models_end):
1412 state = self._last_state
1413 pp = ihm.analysis.ClusterStep(
'RMSD', num_models_begin, num_models_end)
1414 self.all_protocols.add_postproc(pp, state)
1417 def _add_no_postprocessing(self, num_models):
1419 state = self._last_state
1420 pp = ihm.analysis.EmptyStep()
1421 pp.num_models_begin = pp.num_models_end = num_models
1422 self.all_protocols.add_postproc(pp, state)
1425 def _add_simple_ensemble(self, pp, name, num_models, drmsd,
1426 num_models_deposited, localization_densities,
1428 """Add an ensemble generated by ad hoc methods (not using PMI).
1429 This is currently only used by the Nup84 system."""
1431 state = self._last_state
1432 group = ihm.model.ModelGroup(name=state.get_postfixed_name(name))
1433 state.add_model_group(group)
1436 if isinstance(ensemble_file, IMP.pmi1.metadata.FileLocation):
1437 ensemble_file.content_type = \
1438 ihm.location.OutputFileLocation.content_type
1439 self.system.locations.append(ensemble_file)
1440 e = _SimpleEnsemble(pp, group, num_models, drmsd, num_models_deposited,
1442 self.system.ensembles.append(e)
1443 for c
in state.all_modeled_components:
1444 den = localization_densities.get(c,
None)
1447 if isinstance(den, IMP.pmi1.metadata.FileLocation):
1448 den.content_type = \
1449 ihm.location.OutputFileLocation.content_type
1450 e.load_localization_density(state, c, self.asym_units[c], den)
1454 """Point a previously-created ensemble to an 'all-models' file.
1455 This could be a trajectory such as DCD, an RMF, or a multimodel
1458 if isinstance(location, IMP.pmi1.metadata.FileLocation):
1459 location.content_type = ihm.location.OutputFileLocation.content_type
1460 self.system.locations.append(location)
1462 ind = i + self._state_ensemble_offset
1463 self.system.ensembles[ind].file = location
1465 def add_replica_exchange_analysis(self, state, rex, density_custom_ranges):
1471 protocol = self.all_protocols.get_last_protocol(state)
1472 num_models = protocol.steps[-1].num_models_end
1473 pp = _ReplicaExchangeAnalysisPostProcess(rex, num_models)
1474 self.all_protocols.add_postproc(pp, state)
1475 for i
in range(rex._number_of_clusters):
1476 group = ihm.model.ModelGroup(name=state.get_prefixed_name(
1477 'cluster %d' % (i + 1)))
1478 state.add_model_group(group)
1480 e = _ReplicaExchangeAnalysisEnsemble(pp, i, group, 1)
1481 self.system.ensembles.append(e)
1483 for fname, stuple
in sorted(density_custom_ranges.items()):
1484 e.load_localization_density(state, fname, stuple,
1486 for stats
in e.load_all_models(self, state):
1487 m = self.add_model(group)
1490 m.name =
'Best scoring model'
1493 for c
in state.all_modeled_components:
1496 def _get_subassembly(self, comps, name, description):
1497 """Get an Assembly consisting of the given components.
1498 `compdict` is a dictionary of the components to add, where keys
1499 are the component names and values are the sequence ranges (or
1500 None to use all residues in the component)."""
1502 for comp, seqrng
in comps.items():
1503 a = self.asym_units[comp]
1504 asyms.append(a
if seqrng
is None else a(*seqrng))
1506 a = ihm.Assembly(asyms, name=name, description=description)
1509 def _add_foxs_restraint(self, model, comp, seqrange, dataset, rg, chi,
1511 """Add a basic FoXS fit. This is largely intended for use from the
1513 assembly = self._get_subassembly({comp:seqrange},
1514 name=
"SAXS subassembly",
1515 description=
"All components that fit SAXS data")
1516 r = ihm.restraint.SASRestraint(dataset, assembly, segment=
False,
1517 fitting_method=
'FoXS', fitting_atom_type=
'Heavy atoms',
1518 multi_state=
False, radius_of_gyration=rg, details=details)
1519 r.fits[model] = ihm.restraint.SASRestraintFit(chi_value=chi)
1520 self.system.restraints.append(r)
1521 self._add_restraint_dataset(r)
1523 def add_em2d_restraint(self, state, r, i, resolution, pixel_size,
1524 image_resolution, projection_number,
1525 micrographs_number):
1526 r = _EM2DRestraint(state, r, i, resolution, pixel_size,
1527 image_resolution, projection_number,
1529 self.system.restraints.append(r)
1530 self._add_restraint_dataset(r)
1532 def add_em3d_restraint(self, state, target_ps, densities, pmi_restraint):
1534 r = _EM3DRestraint(self, state, pmi_restraint, target_ps, densities)
1535 self.system.restraints.append(r)
1536 self._add_restraint_dataset(r)
1538 def add_zaxial_restraint(self, state, ps, lower_bound, upper_bound,
1539 sigma, pmi_restraint):
1540 self._add_geometric_restraint(state, ps, lower_bound, upper_bound,
1541 sigma, pmi_restraint, self._xy_plane)
1543 def add_yaxial_restraint(self, state, ps, lower_bound, upper_bound,
1544 sigma, pmi_restraint):
1545 self._add_geometric_restraint(state, ps, lower_bound, upper_bound,
1546 sigma, pmi_restraint, self._xz_plane)
1548 def add_xyradial_restraint(self, state, ps, lower_bound, upper_bound,
1549 sigma, pmi_restraint):
1550 self._add_geometric_restraint(state, ps, lower_bound, upper_bound,
1551 sigma, pmi_restraint, self._z_axis)
1553 def _add_geometric_restraint(self, state, ps, lower_bound, upper_bound,
1554 sigma, pmi_restraint, geom):
1555 asym = get_asym_mapper_for_state(self, state, self.__asym_states)
1556 r = _GeometricRestraint(self, state, pmi_restraint, geom,
1557 asym.get_feature(ps),
1558 ihm.restraint.LowerUpperBoundDistanceRestraint(
1559 lower_bound, upper_bound),
1561 self.system.restraints.append(r)
1562 self._add_restraint_dataset(r)
1564 def _get_membrane(self, tor_R, tor_r, tor_th):
1565 """Get an object representing a half-torus membrane"""
1566 if not hasattr(self,
'_seen_membranes'):
1567 self._seen_membranes = {}
1570 membrane_id = tuple(int(x * 100.)
for x
in (tor_R, tor_r, tor_th))
1571 if membrane_id
not in self._seen_membranes:
1572 m = ihm.geometry.HalfTorus(center=self._center_origin,
1573 transformation=self._identity_transform,
1574 major_radius=tor_R, minor_radius=tor_r, thickness=tor_th,
1575 inner=
True, name=
'Membrane')
1576 self._seen_membranes[membrane_id] = m
1577 return self._seen_membranes[membrane_id]
1579 def add_membrane_surface_location_restraint(
1580 self, state, ps, tor_R, tor_r, tor_th, sigma, pmi_restraint):
1581 self._add_membrane_restraint(state, ps, tor_R, tor_r, tor_th, sigma,
1582 pmi_restraint, ihm.restraint.UpperBoundDistanceRestraint(0.))
1584 def add_membrane_exclusion_restraint(
1585 self, state, ps, tor_R, tor_r, tor_th, sigma, pmi_restraint):
1586 self._add_membrane_restraint(state, ps, tor_R, tor_r, tor_th, sigma,
1587 pmi_restraint, ihm.restraint.LowerBoundDistanceRestraint(0.))
1589 def _add_membrane_restraint(self, state, ps, tor_R, tor_r, tor_th,
1590 sigma, pmi_restraint, rsr):
1591 asym = get_asym_mapper_for_state(self, state, self.__asym_states)
1592 r = _GeometricRestraint(self, state, pmi_restraint,
1593 self._get_membrane(tor_R, tor_r, tor_th),
1594 asym.get_feature(ps), rsr, sigma)
1595 self.system.restraints.append(r)
1596 self._add_restraint_dataset(r)
1598 def add_model(self, group, assembly=None, representation=None):
1599 state = self._last_state
1600 if representation
is None:
1601 representation = self.default_representation
1602 protocol = self.all_protocols.get_last_protocol(state)
1603 m = _Model(state.prot, self, protocol,
1604 assembly
if assembly
else state.modeled_assembly,
1609 def _update_locations(self, filelocs):
1610 """Update FileLocation to point to a parent repository, if any"""
1611 all_repos = [m
for m
in self._metadata
1612 if isinstance(m, ihm.location.Repository)]
1613 for fileloc
in filelocs:
1614 ihm.location.Repository._update_in_repos(fileloc, all_repos)
1616 _metadata = property(
lambda self:
1617 itertools.chain.from_iterable(self._each_metadata))
Select non water and non hydrogen atoms.
static bool get_is_setup(const IMP::ParticleAdaptor &p)
A decorator to associate a particle with a part of a protein/DNA/RNA.
static bool get_is_setup(const IMP::ParticleAdaptor &p)
def create_representation
Create a new Representation and return it.
def finalize
Do any final processing on the class hierarchy.
Class for easy writing of PDBs, RMFs, and stat files.
Class to encode a modeling protocol as mmCIF.
A single asymmetric unit in the system.
def exclude_coordinates
Don't record coordinates for the given domain.
Base class for capturing a modeling protocol.
def pmi_residue
Return a single IHM residue indexed using PMI numbering.
static bool get_is_setup(const IMP::ParticleAdaptor &p)
void read_pdb(TextInput input, int model, Hierarchy h)
Classes to represent data structures used in mmCIF.
void add_restraint(RMF::FileHandle fh, Restraint *hs)
static bool get_is_setup(Model *m, ParticleIndex pi)
Ints get_index(const ParticlesTemp &particles, const Subset &subset, const Subsets &excluded)
def set_ensemble_file
Point a previously-created ensemble to an 'all-models' file.
def pmi_range
Return a range of IHM residues indexed using PMI numbering.
Basic utilities for handling cryo-electron microscopy 3D density maps.
def pmi_residue
Return a single IHM residue indexed using PMI numbering.
def create_transformed_component
Make a new component that's a transformed copy of another.
Transformation3D get_identity_transformation_3d()
Return a transformation that does not do anything.
Representation of the system.
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...
A single entity in the system.
def pmi_range
Return a range of IHM residues indexed using PMI numbering.
Classes for writing output files and processing them.
Functionality for loading, creating, manipulating and scoring atomic structures.
Hierarchies get_leaves(const Selection &h)
Select hierarchy particles identified by the biological name.
Select all ATOM and HETATM records with the given chain ids.
Inferential scoring building on methods developed as part of the Inferential Structure Determination ...