From ac8c85eda0e2f4ae63ff698724e971d53c63d7a4 Mon Sep 17 00:00:00 2001 From: Diego Prada-Gracia Date: Wed, 21 Feb 2024 23:00:46 -0600 Subject: [PATCH] In process --- molsysmt/_private/digestion/argument/item.py | 8 ++ .../group/terminal_capping/__init__.py | 2 +- .../{names.py => group_names.py} | 2 +- .../terminal_capping/is_terminal_capping.py | 4 +- molsysmt/form/file_crd/__init__.py | 6 + molsysmt/form/file_crd/to_molsysmt_MolSys.py | 16 +++ .../form/file_crd/to_molsysmt_Structures.py | 52 +++++++++ .../form/file_crd/to_molsysmt_Topology.py | 104 ++++++++++++++++++ .../form/molsysmt_MolSys/append_structures.py | 18 ++- .../form/molsysmt_MolSys/to_file_h5msm.py | 5 +- molsysmt/form/parmed_Structure/get.py | 3 + molsysmt/form/string_aminoacids3/is_form.py | 4 +- molsysmt/native/h5msm_file_handler.py | 2 +- molsysmt/native/structures.py | 31 ++++-- 14 files changed, 234 insertions(+), 23 deletions(-) rename molsysmt/element/group/terminal_capping/{names.py => group_names.py} (62%) create mode 100644 molsysmt/form/file_crd/to_molsysmt_MolSys.py create mode 100644 molsysmt/form/file_crd/to_molsysmt_Structures.py create mode 100644 molsysmt/form/file_crd/to_molsysmt_Topology.py diff --git a/molsysmt/_private/digestion/argument/item.py b/molsysmt/_private/digestion/argument/item.py index a0401ccbb..6b1bac59e 100644 --- a/molsysmt/_private/digestion/argument/item.py +++ b/molsysmt/_private/digestion/argument/item.py @@ -25,6 +25,7 @@ def digest_item(item, form=None, caller=None): >>> digest_item(item=sequence, form="biopython.Seq") """ + from molsysmt.basic import get_form try: in_form = get_form(item) @@ -32,6 +33,13 @@ def digest_item(item, form=None, caller=None): except: output = False + if caller.startswith('molsysmt.form') and caller.endswith('append_structures'): + + if item is None: + return None + elif output==True: + return item + if output: if form is not None: if in_form!=form: diff --git a/molsysmt/element/group/terminal_capping/__init__.py b/molsysmt/element/group/terminal_capping/__init__.py index ac835e101..99a961336 100644 --- a/molsysmt/element/group/terminal_capping/__init__.py +++ b/molsysmt/element/group/terminal_capping/__init__.py @@ -1,2 +1,2 @@ -from .names import n_terminal_capping_names, c_terminal_capping_names, names +from .group_names import n_terminal_capping_names, c_terminal_capping_names, group_names from .is_terminal_capping import is_terminal_capping diff --git a/molsysmt/element/group/terminal_capping/names.py b/molsysmt/element/group/terminal_capping/group_names.py similarity index 62% rename from molsysmt/element/group/terminal_capping/names.py rename to molsysmt/element/group/terminal_capping/group_names.py index f5d779fcc..012bed049 100644 --- a/molsysmt/element/group/terminal_capping/names.py +++ b/molsysmt/element/group/terminal_capping/group_names.py @@ -6,5 +6,5 @@ 'NME' ] -names = n_terminal_capping_names + c_terminal_capping_names +group_names = n_terminal_capping_names + c_terminal_capping_names diff --git a/molsysmt/element/group/terminal_capping/is_terminal_capping.py b/molsysmt/element/group/terminal_capping/is_terminal_capping.py index b35c81fdd..cdc183010 100644 --- a/molsysmt/element/group/terminal_capping/is_terminal_capping.py +++ b/molsysmt/element/group/terminal_capping/is_terminal_capping.py @@ -1,9 +1,9 @@ -from .names import names +from .group_names import group_names def is_terminal_capping(name): """ To be written soon... """ - return (name in names) + return (name in group_names) diff --git a/molsysmt/form/file_crd/__init__.py b/molsysmt/form/file_crd/__init__.py index 13409845a..f40825515 100644 --- a/molsysmt/form/file_crd/__init__.py +++ b/molsysmt/form/file_crd/__init__.py @@ -16,6 +16,9 @@ from .set import * from .iterators import StructuresIterator, TopologyIterator +from .to_molsysmt_MolSys import to_molsysmt_MolSys +from .to_molsysmt_Topology import to_molsysmt_Topology +from .to_molsysmt_Structures import to_molsysmt_Structures from .to_molsysmt_MolSysOld import to_molsysmt_MolSysOld from .to_molsysmt_TopologyOld import to_molsysmt_TopologyOld from .to_molsysmt_StructuresOld import to_molsysmt_StructuresOld @@ -23,6 +26,9 @@ _convert_to={ 'file:crd': extract, + 'molsysmt.MolSys': to_molsysmt_MolSys, + 'molsysmt.Topology': to_molsysmt_Topology, + 'molsysmt.Structures': to_molsysmt_Structures, 'molsysmt.MolSysOld': to_molsysmt_MolSysOld, 'molsysmt.TopologyOld': to_molsysmt_TopologyOld, 'molsysmt.StructuresOld': to_molsysmt_StructuresOld, diff --git a/molsysmt/form/file_crd/to_molsysmt_MolSys.py b/molsysmt/form/file_crd/to_molsysmt_MolSys.py new file mode 100644 index 000000000..6f86d7352 --- /dev/null +++ b/molsysmt/form/file_crd/to_molsysmt_MolSys.py @@ -0,0 +1,16 @@ +from molsysmt._private.digestion import digest + +@digest(form='file:crd') +def to_molsysmt_MolSys(item, atom_indices='all', structure_indices='all', get_missing_bonds=False, + skip_digestion=False): + + from molsysmt.native.molsys import MolSys + from .to_molsysmt_Topology import to_molsysmt_Topology + from .to_molsysmt_Structures import to_molsysmt_Structures + + tmp_item = MolSys() + tmp_item.topology = to_molsysmt_Topology(item, atom_indices=atom_indices, structure_indices=structure_indices) + tmp_item.structures = to_molsysmt_Structures(item, atom_indices=atom_indices, structure_indices=structure_indices) + + return tmp_item + diff --git a/molsysmt/form/file_crd/to_molsysmt_Structures.py b/molsysmt/form/file_crd/to_molsysmt_Structures.py new file mode 100644 index 000000000..60fc8a07a --- /dev/null +++ b/molsysmt/form/file_crd/to_molsysmt_Structures.py @@ -0,0 +1,52 @@ +from molsysmt._private.digestion import digest +from molsysmt import pyunitwizard as puw +import numpy as np + +@digest(form='file:crd') +def to_molsysmt_Structures(item, atom_indices='all', structure_indices='all', skip_digestion=False): + + # EXT: + # (i10,2x,a) natoms,'EXT' + # (2I10,2X,A8,2X,A8,3F20.10,2X,A8,2X,A8,F20.10) + # iatom,ires,resn,typr,x,y,z,segid,rid,wmain + # standard: + # (i5) natoms + # (2I5,1X,A4,1X,A4,3F10.5,1X,A4,1X,A4,F10.5) + # iatom,ires,resn,typr,x,y,z,segid,orig_resid,wmain + + from molsysmt.native.structures import Structures + + tmp_item = Structures() + + x = [] + y = [] + z = [] + + extended = False + + with open(item) as fff: + for line in fff: + if line.strip().startswith('*') or line.strip() == "": + continue + field = line.split() + if len(field)==1: + n_atoms = int(field[0]) + elif len(field)==2: + n_atoms = int(field[0]) + extended = True + else: + x.append(float(field[4])) + y.append(float(field[5])) + z.append(float(field[6])) + + if len(x)!=n_atoms: + raise ValueError + + coordinates = puw.quantity(np.expand_dims(np.column_stack([x,y,z]), axis=0), unit='angstroms', standardized=True) + + tmp_item.append(structure_id=[0], coordinates=coordinates) + + del(x,y,z,coordinates) + + return tmp_item + diff --git a/molsysmt/form/file_crd/to_molsysmt_Topology.py b/molsysmt/form/file_crd/to_molsysmt_Topology.py new file mode 100644 index 000000000..db729c6f4 --- /dev/null +++ b/molsysmt/form/file_crd/to_molsysmt_Topology.py @@ -0,0 +1,104 @@ +from molsysmt._private.digestion import digest +from molsysmt import pyunitwizard as puw +from molsysmt.element.atom.get_atom_type import _get_atom_type_from_atom_name +from molsysmt.element.group.get_group_type import _get_group_type_from_group_name +import numpy as np + +@digest(form='file:crd') +def to_molsysmt_Topology(item, atom_indices='all', structure_indices='all', skip_digestion=False): + + # EXT: + # (i10,2x,a) natoms,'EXT' + # (2I10,2X,A8,2X,A8,3F20.10,2X,A8,2X,A8,F20.10) + # iatom,ires,resn,typr,x,y,z,segid,rid,wmain + # standard: + # (i5) natoms + # (2I5,1X,A4,1X,A4,3F10.5,1X,A4,1X,A4,F10.5) + # iatom,ires,resn,typr,x,y,z,segid,orig_resid,wmain + + from molsysmt.native.topology import Topology + + tmp_item = Topology() + + atom_id = [] + atom_name = [] + atom_type = [] + group_index = [] + group_id = [] + group_name = [] + group_type = [] + chain_index = [] + chain_name = [] + + extended = False + + former_group_id = -1 + former_chain_name = '' + + aux_group_index = -1 + aux_chain_index = -1 + + with open(item) as fff: + for line in fff: + if line.strip().startswith('*') or line.strip() == "": + continue + field = line.split() + if len(field)==1: + n_atoms = int(field[0]) + elif len(field)==2: + n_atoms = int(field[0]) + extended = True + else: + atom_id.append(int(field[0])) + atom_name.append(field[3]) + + aux_group_id = int(field[1]) + aux_chain_name = field[7] + + if former_group_id!=aux_group_id: + former_group_id=aux_group_id + group_id.append(aux_group_id) + group_name.append(field[2]) + aux_group_index += 1 + + if former_chain_name!=aux_chain_name: + former_chain_name=aux_chain_name + chain_name.append(aux_chain_name) + aux_chain_index += 1 + + group_index.append(aux_group_index) + chain_index.append(aux_chain_index) + + if len(atom_id)!=n_atoms: + raise ValueError + + for ii in atom_name: + atom_type.append(_get_atom_type_from_atom_name(ii)) + + for ii in group_name: + group_type.append(_get_group_type_from_group_name(ii)) + + n_groups = len(group_name) + n_chains = len(chain_name) + + tmp_item.reset_atoms(n_atoms=n_atoms) + tmp_item.reset_groups(n_groups=n_groups) + tmp_item.reset_chains(n_chains=n_chains) + + tmp_item.atoms.atom_id = np.array(atom_id, dtype=int) + tmp_item.atoms.atom_name = np.array(atom_name, dtype=object) + tmp_item.atoms.atom_type = np.array(atom_type, dtype=object) + tmp_item.atoms.group_index = np.array(group_index, dtype=int) + tmp_item.atoms.chain_index = np.array(chain_index, dtype=int) + tmp_item.groups.group_id = np.array(group_id, dtype=int) + tmp_item.groups.group_name = np.array(group_name, dtype=object) + tmp_item.groups.group_type = np.array(group_type, dtype=object) + tmp_item.chains.chain_id = np.arange(n_chains, dtype=int) + tmp_item.chains.chain_name = np.array(chain_name, dtype=object) + + del(atom_id, atom_name, atom_type, + group_index, group_id, group_name, group_type, + chain_index, chain_name) + + return tmp_item + diff --git a/molsysmt/form/molsysmt_MolSys/append_structures.py b/molsysmt/form/molsysmt_MolSys/append_structures.py index 91d395987..b77018544 100644 --- a/molsysmt/form/molsysmt_MolSys/append_structures.py +++ b/molsysmt/form/molsysmt_MolSys/append_structures.py @@ -1,10 +1,22 @@ from molsysmt._private.exceptions import NotImplementedMethodError from molsysmt._private.digestion import digest +from molsysmt._private.variables import is_all -@digest(form='molsysmt.MolSys') -def append_structures(item, atom_indices='all', structure_indices='all', skip_digestion=False): - item.append_structures(item, atom_indices=atom_indices, structure_indices=structure_indices, skip_digestion=True) +@digest(form='molsysmt.MolSys', to_form='molsysmt.MolSys') +def append_structures(to_item, item=None, structure_id=None, time=None, coordinates=None, velocities=None, + box=None, temperature=None, potential_energy=None, kinetic_energy=None, + atom_indices='all', structure_indices='all', skip_digestion=False): + + if item is not None: + to_item.structures.append_structures(item, atom_indices=atom_indices, structure_indices=structure_indices, + skip_digestion=True) + else: + to_item.structures.append(structure_id=structure_id, time=time, coordinates=coordinates, + velocities=velocities, box=box, temperature=temperature, + potential_energy=potential_energy, kinetic_energy=kinetic_energy, + atom_indices=atom_indices, structure_indices=structure_indices, + skip_digestion=True) pass diff --git a/molsysmt/form/molsysmt_MolSys/to_file_h5msm.py b/molsysmt/form/molsysmt_MolSys/to_file_h5msm.py index ab0583e81..d922e5fc5 100644 --- a/molsysmt/form/molsysmt_MolSys/to_file_h5msm.py +++ b/molsysmt/form/molsysmt_MolSys/to_file_h5msm.py @@ -13,9 +13,8 @@ def to_file_h5msm(item, atom_indices='all', structure_indices='all', output_file compression_opts=compression_opts, int_precision=int_precision, float_precision=float_precision, closed=False, skip_digestion=True) - _add_topology_to_h5msm(item.topology, handler, atom_indices=atom_indices, skip_digestion=True) - _add_structures_to_h5msm(item.structures, handler, atom_indices=atom_indices, structure_indices=structure_indices, - skip_digestion=True) + _add_topology_to_h5msm(item.topology, handler, atom_indices=atom_indices) + _add_structures_to_h5msm(item.structures, handler, atom_indices=atom_indices, structure_indices=structure_indices) handler.close() diff --git a/molsysmt/form/parmed_Structure/get.py b/molsysmt/form/parmed_Structure/get.py index c539e656e..e5f53802d 100644 --- a/molsysmt/form/parmed_Structure/get.py +++ b/molsysmt/form/parmed_Structure/get.py @@ -66,6 +66,9 @@ def get_coordinates_from_atom(item, indices='all', structure_indices='all', skip output = item.coordinates + if len(output.shape)==2 and output.shape[1]==3: + output = output[np.newaxis,:,:] + if output is None: return output diff --git a/molsysmt/form/string_aminoacids3/is_form.py b/molsysmt/form/string_aminoacids3/is_form.py index fc2bff8d3..e7c326fb1 100644 --- a/molsysmt/form/string_aminoacids3/is_form.py +++ b/molsysmt/form/string_aminoacids3/is_form.py @@ -13,8 +13,8 @@ def is_form(item): elif not ' ' in item: - from molsysmt.element.group.amino_acid import names as aminoacid_names - from molsysmt.element.group.terminal_capping import names as terminal_capping_names + from molsysmt.element.group.amino_acid import group_names as aminoacid_names + from molsysmt.element.group.terminal_capping import group_names as terminal_capping_names tmp_item = item.upper() diff --git a/molsysmt/native/h5msm_file_handler.py b/molsysmt/native/h5msm_file_handler.py index 01c00fe45..bc6ec054d 100644 --- a/molsysmt/native/h5msm_file_handler.py +++ b/molsysmt/native/h5msm_file_handler.py @@ -8,7 +8,7 @@ class H5MSMFileHandler(): def __init__(self, filename, io_mode='r', creator='MolSysMT', compression="gzip", compression_opts=4, int_precision='single', float_precision='single', length_unit=None, time_unit=None, energy_unit=None, - temperature_unit=None, charge_unit=None, mass_unit=None, closed=False): + temperature_unit=None, charge_unit=None, mass_unit=None, closed=False, skip_digestion=False): self.file = None diff --git a/molsysmt/native/structures.py b/molsysmt/native/structures.py index 42afab4c0..6e57a58fd 100644 --- a/molsysmt/native/structures.py +++ b/molsysmt/native/structures.py @@ -132,7 +132,7 @@ def _append_coordinates(self, coordinates, atom_indices='all', structure_indices raise ValueError( "The coordinates to be appended in the system " "need to have the same number of atoms.") - self.coordinates = self._concatenate_arrays(self.coordinates, coordinates) + self.coordinates = _concatenate_arrays(self.coordinates, coordinates) if is_all(structure_indices): if is_all(atom_indices): self.coordinates = np.concatenate([self.coordinates, coordinates]) @@ -166,7 +166,7 @@ def _append_velocities(self, velocities, atom_indices='all', structure_indices=' raise ValueError( "The velocities to be appended in the system " "need to have the same number of atoms.") - self.velocities = self._concatenate_arrays(self.velocities, velocities) + self.velocities = _concatenate_arrays(self.velocities, velocities) if is_all(structure_indices): if is_all(atom_indices): self.velocities = np.concatenate([self.velocities, velocities]) @@ -376,31 +376,31 @@ def append(self, structure_id=None, time=None, coordinates=None, velocities=None else: if (self.structure_id is not None) and (structure_id is not None): - self._append_structure_id(self, structure_id, structure_indices=structure_indices, + self._append_structure_id(structure_id, structure_indices=structure_indices, skip_digestion=True) if (self.time is not None) and (time is not None): - self._append_time(self, time, structure_indices=structure_indices, skip_digestion=True) + self._append_time(time, structure_indices=structure_indices, skip_digestion=True) if (self.coordinates is not None) and (coordinates is not None): - self._append_coordinates(self, coordinates, atom_indices=atom_indices, + self._append_coordinates(coordinates, atom_indices=atom_indices, structure_indices=structure_indices, skip_digestion=True) if (self.velocities is not None) and (velocities is not None): - self._append_velocities(self, velocities, atom_indices=atom_indices, + self._append_velocities(velocities, atom_indices=atom_indices, structure_indices=structure_indices, skip_digestion=True) if (self.box is not None) and (box is not None): - self._append_box(self, box, structure_indices=structure_indices, skip_digestion=True) + self._append_box(box, structure_indices=structure_indices, skip_digestion=True) if (self.temperature is not None) and (temperature is not None): - self._append_temperature(self, temperature, structure_indices=structure_indices, skip_digestion=True) + self._append_temperature(temperature, structure_indices=structure_indices, skip_digestion=True) if (self.potential_energy is not None) and (potential_energy is not None): - self._append_potential_energy(self, potential_energy, structure_indices=structure_indices) + self._append_potential_energy(potential_energy, structure_indices=structure_indices) if (self.kinetic_energy is not None) and (kinetic_energy is not None): - self._append_kinetic_energy(self, kinetic_energy, structure_indices=structure_indices, + self._append_kinetic_energy(kinetic_energy, structure_indices=structure_indices, skip_digestion=True) if self.n_structures==0 and self.n_atoms==0: @@ -625,3 +625,14 @@ def copy(self): """ Returns a copy of the structures.""" return deepcopy(self) + +@staticmethod +def _concatenate_arrays(array_1, array_2): + """ Concatenates two arrays provided that they are not null.""" + if array_2 is not None: + if array_1 is None: + raise ValueError( + f"The trajectory has no array to append the new frame.") + else: + return np.concatenate([array_1, array_2]) +