protocaller/ProtoCaller

minimise not take effect when protonated was set to True

Closed this issue · 2 comments

kexul commented

Here is the code to reproduce:

lig = Ligand('test1.sdf', protonated=True, minimise=True)
lig.parametrise()

Here is the sdf file I used, which is 2D structure.
test1.zip

The generated inpcrd file and mol2 file is 2D, though minimise is set to True (3D structure was expected).

In Ligand.py, the parametrise function only take effect on self.protonated_filename . It may be better if self.molecule is parametrised.

def parametrise(self, params=None, molecule_type="ligand", id=None, reparametrise=False):
        """
        Parametrises the ligand using ProtoCaller.Parametrise.

        Parameters
        ----------
        params : ProtoCaller.Parametrise.Params
            Force field parameters.
        molecule_type : str
            The type of the molecule. One of: "ligand" and "cofactor".
        id : str
            The name of the molecule. Default: equal to the ligand name.
        reparametrise : bool
            Whether to reparametrise an already parametrised ligand.
        """
        with self.workdir:
            if self._parametrised and not reparametrise:
                _logging.debug("Ligand %s is already parametrised." % self.name)
                return

            _logging.info("Parametrising ligand %s..." % self.name)
            if not self.protonated:
                _logging.warning("Cannot parametrise unprotonated ligand. Protonating first with default parameters...")
                self.protonate()

            if params is None:
                params = _parametrise.Params()

            # we convert the protonated file into a pdb so that antechamber can read it
            filename = _babel.babelTransform(self.protonated_filename, "pdb")
            if id is None: id = self.name

            charge = _rdmolops.GetFormalCharge(self.molecule)
            self.parametrised_files = _parametrise.parametriseFile(params=params, filename=filename,
                                                                   molecule_type=molecule_type, id=id, charge=charge)

The code for this class is a bit convoluted and probably needs a more thorough refactoring, although I don't have the time to do so. In any case, I just pushed a commit that should take care of this behaviour. Let me know if there are any issues with it.

kexul commented

Thanks, the new version worked like a charm.