diff --git a/mdpow/analysis/hbond.py b/mdpow/analysis/hbond.py new file mode 100644 index 00000000..d137a620 --- /dev/null +++ b/mdpow/analysis/hbond.py @@ -0,0 +1,86 @@ +# MDPOW: hbond.py +# 2022 Cade Duckworth + +#Review and update old docs + +import pandas as pd + +from MDAnalysis.analysis.hydrogenbonds.hbond_analysis import HydrogenBondAnalysis as HBA + +from mdpow.analysis.ensemble import Ensemble, EnsembleAtomGroup, EnsembleAnalysis + +import logging + +logger = logging.getLogger('mdpow.analysis.hbond') + +class HBondAnalysis(EnsembleAnalysis): + """Analyzes potential hydrogen bonds of solute from a single + :class:`~mdpow.analysis.ensemble.Ensemble` . + :keywords: + *solute* + single :class:`~mdpow.analysis.ensemble.EnsembleAtomGroup` + selected by resname provided for Universe. Selections must be from the same + :class:`~mdpow.analysis.ensemble.Ensemble` . + Data is returned in a :class:`pandas.DataFrame` with observations sorted by + solvent, interaction, lambda, frame, donor index, + hydrogen index, acceptor index, distance, angle. + .. ruberic:: Example + Typical Workflow:: + ens = Ensemble(dirname='Mol', solvents=('water', 'octanol'), + interactions=('Coulomb', 'VDW')) + solute = ens.select_atoms('resname UNK') + hb = HBondAnalysis(solute, acceptors_sel='resname UNK', d_a_cutoff=2.5) + hb.run(start=0, stop=1000, step=10) + + HydrogenBondAnalysis + numpy array results format: + , , + , + , + , + + HydrogenBondAnalysis + :keywords: + acceptors_sel=None + donors_sel=None + hydrogens_sel=None + d_h_cutoff=1.2 + d_a_cutoff=3.0 + d_h_a_angle_cutoff=150 + update_selections=True + + https://docs.mdanalysis.org/stable/documentation_pages/analysis/hydrogenbonds.html + + ***NEED TO REFINE THIS AND ADD MORE DETAIL*** + """ + + def __init__(self, solute: EnsembleAtomGroup, **kwargs): + super(HBondAnalysis, self).__init__(solute.ensemble) + self._solute = solute + self.hb_kwargs = kwargs + + def _prepare_ensemble(self): + self._col = ['solvent', 'interaction', + 'lambda', 'time', 'frame', + 'donor index', 'hydrogen index', + 'acceptor index', 'distance', 'angle'] + self.results = pd.DataFrame(columns=self._col) + self._res_dict = {key: [] for key in self._col} + + def _single_universe(self): + solute = self._solute[self._key] + hb = HBA(universe=self._solute.ensemble[self._key], **self.hb_kwargs) + hb.run(verbose=True, start=self.start, stop=self.stop, step=self.step) + results = hb.results.hbonds + + for h in results: + self._time = self._system.trajectory[int(h[0])].time + result = [self._key[0], self._key[1],self._key[2], + self._time, int(h[0]), int(h[1]), + int(h[2]), int(h[3]), h[4], h[5]] + for i in range(len(self._col)): + self._res_dict[self._col[i]].append(result[i]) + + def _conclude_ensemble(self): + for k in self._col: + self.results[k] = self._res_dict[k]