gochem: github.com/rmera/gochem Index | Files | Directories

package chem

import "github.com/rmera/gochem"

Package chem is the main package of the goChem library. It provides atom and molecule structures, facilities for reading and writing some files used in computational chemistry and some functions for geometric manipulations and shape indicators.

	**goChem Capabilities**

    Reads/writes PDB and XYZ files.

    Reads XTC and DCD trajectory files, both sequentially and concurrently.

    Superimposes molecules (especially adequate for non-proteins since doesn't
	use sequence information). The user specify what atoms to use for the
	superimposing transformation calculation. Then all the atoms will be
	superimposed accordingly. Thus, non-identical molecules can be superimposed.

    Calculates RMSD between sets of coordinates.

    Allows to select atoms and coordinates by using a go slice of indexes.

    Allows to replace selected coordinates for a new set.

    Calculates moment tensor and elipsoid of inertia--related properties

    The Molecule object implements the sort.Interface interface, so atoms
	can easily be sorted by b-factors.

    Axis manipulation.
        Align a vector with the Z axis.
        Rotate around the Z axis until the xy projection of a vector becomes
		the Y axis.
        Rotate a sub-group of atoms in a molecule using any 2 coordinates as
		the rotation axis.

    The latter is implemented using Clifford algebra and, as a legacy version, Euler
	angles and rotation matrices (math for the Cliffor algebra implementation by Doc.
	Dr. Janne Pesonen). The Clifford algebra implementation is concurrent. In general,
	Clifford algebra is mathematically better behaving than Euler angles, which are
	not defined for certain rotations.

    Calculates and draws Ramachandran plots (uses the Plotinum library). for an
	aminoacidic chain or a subset of it.

    Generates input for, run and recover results from QM calculations with Turbomole,
	Orca and MOPAC (which must be obtained independently from their respective
	distributors). Interfacing gochem to other QM codes is fairly simple.

    goChem data can be JSON encoded and transfered in such a way that PyMOL
	(http://www.pymol.org plugins can send and received data from/to goChem programs)
	thus enabling the build of PyMOL plugins using goChem.

goChem implements its own matrix type for coordinates, VecMatrix, based in github.com/gonum/matrix.

Currently, each row of a VecMatrix represents one point in space. As this could change if github.com/gonum/matrix changes to col major (which seem unlikely at this point, but still), I recomend prefering the Vec* methods over the Row* methods when manipulating a VecMatrix. Vec* methods will change from row to column following the upstream Gonum library.

Index

Package Files

chem.go clifford.go conversion.go doc.go files.go geometric.go handy.go hydrogens.go interfaces.go matrixhelp.go ramacalc.go

Constants

const (
    ErrNilData          = PanicMsg("goChem: Nil data given ")
    ErrInconsistentData = PanicMsg("goChem: Inconsistent data length ")
    ErrNilMatrix        = PanicMsg("goChem: Attempted to access nil v3.Matrix or gonum/mat64.Dense")
    ErrNilAtoms         = PanicMsg("goChem: Topology has a nil []*Atom slice")
    ErrNilAtom          = PanicMsg("goChem: Attempted to copy from or to a nil Atom")
    ErrAtomOutOfRange   = PanicMsg("goChem: Requested/Attempted setting Atom out of range")
    ErrNilFrame         = PanicMsg("goChem: Attempted to acces nil frame")
    ErrNotXx3Matrix     = PanicMsg("goChem: A v3.Matrix should have 3 columns")
    ErrCliffordRotation = PanicMsg("goChem-Clifford: Target and Result matrices must have the same dimensions. They cannot reference the same matrix") //the only panic that the Clifford functions throw.
)
const (
    Deg2Rad = 0.0174533
    Rad2Deg = 1 / 0.0174533
    H2Kcal  = 627.509 //HArtree 2 Kcal/mol
    Kcal2H  = 1 / 627.509
    KJ2Kcal = 1 / 4.184
    Kcal2KJ = 4.184
    A2Bohr  = 1.889725989
    Bohr2A  = 1 / 1.889725989
    EV2Kcal = 23.061
    Kcal2EV = 1 / 23.061
)

Conversions

const (
    CHDist = 1.098 //C(sp3)--H distance in A
)

Others

func Angle Uses

func Angle(v1, v2 *v3.Matrix) float64

Angle takes 2 vectors and calculate the angle in radians between them It does not check for correctness or return errors!

func AntiProjection Uses

func AntiProjection(test, ref *v3.Matrix) *v3.Matrix

AntiProjection returns a vector in the direction of ref with the magnitude of a vector A would have if |test| was the magnitude of its projection in the direction of test.

func BestPlane Uses

func BestPlane(coords *v3.Matrix, mol ...Masser) (*v3.Matrix, error)

BestPlane returns a row vector that is normal to the plane that best contains the molecule if passed a nil Masser, it will simply set all masses to 1. If more than one Masser is passed Only the first will be considered

func BestPlaneP Uses

func BestPlaneP(evecs *v3.Matrix) (*v3.Matrix, error)

BestPlaneP takes sorted evecs, according to the eval,s and returns a row vector that is normal to the Plane that best contains the molecule. Notice that the function can't possibly check that the vectors are sorted. The P at the end of the name is for Performance. If That is not an issue it is safer to use the BestPlane function that wraps this one.

func CenterOfMass Uses

func CenterOfMass(geometry *v3.Matrix, massS ...*mat.Dense) (*v3.Matrix, error)

CenterOfMass returns the center of mass the atoms represented by the coordinates in geometry and the masses in mass, and an error. If no mass is given, it calculates the geometric center

func Corrupted Uses

func Corrupted(X Traj, R Atomer) error

Corrupted is a convenience function to check that a reference and a trajectory have the same number of atoms

func CutAlphaRef Uses

func CutAlphaRef(r Atomer, chain []string, list []int) []int

CutAlphaRef will return a list with the atoms in the residues indicated by list, in the chains given. The carbonyl carbon and amide nitrogen for each residue will be transformer into hydrogens. The MolID of the other backbone atoms will be set to -1 so they are no longer considered.

func CutBackRef Uses

func CutBackRef(r Atomer, chains []string, list [][]int) ([]int, error)

CutBackRef takes a list of lists of residues and deletes from r all atoms not in the list or not belonging to the chain chain. It caps the N and C terminal of each list with -COH for the N terminal and NH2 for C terminal. the residues on each sublist should be contiguous to each other. for instance, {6,7,8} is a valid sublist, {6,8,9} is not. This is NOT currently checked by the function!. It returns the list of kept atoms

func CutBetaRef Uses

func CutBetaRef(r Atomer, chain []string, list []int) []int

CutLateralRef will return a list with the atom indexes of the lateral chains of the residues in list for each of these residues it will change the alpha carbon to oxygen and change the residue number of the rest of the backbone to -1.

func Dihedral Uses

func Dihedral(a, b, c, d *v3.Matrix) float64

Dihedral calculates the dihedral between the points a, b, c, d, where the first plane is defined by abc and the second by bcd.

func EasyShape Uses

func EasyShape(coords *v3.Matrix, epsilon float64, mol ...Masser) (float64, float64, error)

Easy shape takes a matrix of coordinates, a value for epsilon (a number close to zero, the close, the more strict the orthogonality requriements are) and an (optative) masser and returns two shape indicators based on the elipsoid of inertia (or it massless equivalent) a linear and circular distortion indicators, and an error or nil (in that order). if you give a negative number as epsilon, the default (quite strict) will be used.

func EulerRotateAbout Uses

func EulerRotateAbout(coordsorig, ax1, ax2 *v3.Matrix, angle float64) (*v3.Matrix, error)

EulerRotateAbout uses Euler angles to rotate the coordinates in coordsorig around by angle radians around the axis given by the vector axis. It returns the rotated coordsorig, since the original is not affected. It seems more clunky than the RotateAbout, which uses Clifford algebra. I leave it for benchmark, mostly, and might remove it later. There is no test for this function!

func FixGromacsPDB Uses

func FixGromacsPDB(mol Atomer)

FixGromacsPDB fixes the problem that Gromacs PDBs have when there are more than 10000 residues Gromacs simply restarts the numbering. Since solvents (where this is likely to happen) don't have chain ID in Gromacs, it's impossible to distinguish between the water 1 and the water 10001. FixGromacsPDB Adds a chain ID to the newly restrated residue that is the letter/symbol coming after the last seen chain ID in the constant allchains defined in this file. The current implementation does nothing if a chain ID is already defined, even if it is wrong (if 9999 and the following 0 residue have the same chain).

func FixNumbering Uses

func FixNumbering(r Atomer)

FixNumbering will put the internal numbering+1 in the atoms and residue fields, so they match the current residues/atoms in the molecule

func MakeWater Uses

func MakeWater(a1, a2 *v3.Matrix, distance, angle float64, oxygen bool) *v3.Matrix

MakeWater Creates a water molecule at distance Angstroms from a2, in a direction that is angle radians from the axis defined by a1 and a2. Notice that the exact position of the water is not well defined when angle is not zero. One can always use the RotateAbout function to move the molecule to the desired location. If oxygen is true, the oxygen will be pointing to a2. Otherwise, one of the hydrogens will.

func MassCenter Uses

func MassCenter(in, oref *v3.Matrix, massS ...*mat.Dense) (*v3.Matrix, *v3.Matrix, error)

MassCenter centers in in the center of mass of oref. Mass must be A column vector. Returns the centered matrix and the displacement matrix.

func MemRMSD Uses

func MemRMSD(test, templa, tmp *v3.Matrix, indexes ...[]int) (float64, error)

MemRMSD calculates the RMSD between test and template, considering only the atoms present in the slices of int slices indexes. The first indexes slices will be assumed to contain test indexes and the second, template indexes. If you give only one (it must be the first one), it will be assumed to correspond to whatever molecule that has more atoms than the elements in the slice. Giving a nil or 0-lenght first slice and a non-nil second slice will cause MemRMSD to not consider neither of them. The same number of atoms has to be considered for the calculation in both systems. It does not superimpose the objects. To save memory, it asks for the temporary matrix it needs to be supplied: tmp must be Nx3 where N is the number of elements in testlst and templalst

func MolIDNameChain2Index Uses

func MolIDNameChain2Index(mol Atomer, molID int, name, chain string) (int, error)

MolIDNameChain2Index takes a molID (residue number), atom name, chain index and a molecule Atomer. it returns the index associated with the atom in question in the Ref. The function returns also an error (if failure of warning) or nil (if succses and no warnings). Note that this function is not efficient to call several times to retrieve many atoms.

func Molecules2Atoms Uses

func Molecules2Atoms(mol Atomer, residues []int, chains []string) []int

Molecules2Atoms gets a selection list from a list of residues. It select all the atoms that form part of the residues in the list. It doesnt return errors, if a residue is out of range, no atom will be returned for it. Atoms are also required to be part of one of the chains specified in chains.

func MomentTensor Uses

func MomentTensor(A *v3.Matrix, massslice ...[]float64) (*v3.Matrix, error)

MomentTensor returns the moment tensor for a matrix A of coordinates and a column vector mass with the respective massess.

func MultiPDBWrite Uses

func MultiPDBWrite(out io.Writer, Coords []*v3.Matrix, mol Atomer, Bfactors [][]float64) error

MultiPDBWrite writes a multiPDB for the molecule mol and the various coordinate sets in CandB, to the io.Writer given. CandB is a list of lists of *matrix.DenseMatrix. If it has 2 elements or more, the second will be used as Bfactors. If it has one element, all b-factors will be zero. Returns an error if fails, or nil if succeeds.

func OnesMass Uses

func OnesMass(lenght int) *v3.Matrix

OnesMass returns a column matrix with lenght rosw. This matrix can be used as a dummy mass matrix for geometric calculations.

func PDBFileWrite Uses

func PDBFileWrite(pdbname string, coords *v3.Matrix, mol Atomer, Bfactors []float64) error

PDBFileWrite writes a PDB for the molecule mol and the coordinates Coords.

func PDBStringWrite Uses

func PDBStringWrite(coords *v3.Matrix, mol Atomer, bfact []float64) (string, error)

PDBStringWrite writes a string in PDB format for a given reference, coordinate set and bfactor set, which must match each other returns the written string and error or nil.

func PDBWrite Uses

func PDBWrite(out io.Writer, coords *v3.Matrix, mol Atomer, bfact []float64) error

PDBWrite writes a PDB formatted sequence of bytes to an io.Writer for a given reference, coordinate set and bfactor set, which must match each other returns error or nil.

func Projection Uses

func Projection(test, ref *v3.Matrix) *v3.Matrix

Projection returns the projection of test in ref.

func RMSD Uses

func RMSD(test, templa *v3.Matrix, indexes ...[]int) (float64, error)

RMSD calculates the RMSD between test and template, considering only the atoms present in the slices of int slices indexes. The first indexes slices will be assumed to contain test indexes and the second, template indexes. If you give only one, it will be assumed to correspondo to whatever molecule that has more atoms than the elements in the slice. The same number of atoms has to be considered for superposition in both systems. The objects are not superimposed before the calculation.

func RamaCalc Uses

func RamaCalc(M *v3.Matrix, dihedrals []RamaSet) ([][]float64, error)

RamaCalc Obtains the values for the phi and psi dihedrals indicated in []Ramaset, for the structure M. The angles are in *degrees*. It returns a slice of 2-element slices, one for the phi the next for the psi dihedral, a and an error or nil.

func RhoShapeIndexes Uses

func RhoShapeIndexes(rhos []float64) (float64, float64, error)

RhoShapeIndexes Get shape indices based on the axes of the elipsoid of inertia. linear and circular distortion, in that order, and error or nil. Based on the work of Taylor et al., .(1983), J Mol Graph, 1, 30 This function has NOT been tested thoroughly in the sense of the appropiateness of the indexes definitions.

func Rhos Uses

func Rhos(momentTensor *v3.Matrix, epsilon ...float64) ([]float64, error)

Rhos returns the semiaxis of the elipoid of inertia given the the moment of inertia tensor.

func Rotate Uses

func Rotate(Target, Res, axis *v3.Matrix, angle float64) *v3.Matrix

Rotate takes the matrix Target and uses Clifford algebra to _concurrently_ rotate each of its rows by angle radians around axis. Axis must be a 3D row vector. Target must be an N,3 matrix. The result is put in Res, which is also returned.

func RotateAbout Uses

func RotateAbout(coordsorig, ax1, ax2 *v3.Matrix, angle float64) (*v3.Matrix, error)

RotateAbout about rotates the coordinates in coordsorig around by angle radians around the axis given by the vector axis. It returns the rotated coordsorig, since the original is not affected. Uses Clifford algebra.

func RotateP Uses

func RotateP(Target, Res, axis, Rv, Rvrev, tmp1, tmp2, tmp3, tmp4, itmp1, itmp2, itmp3 *v3.Matrix, angle float64, gorut int)

Rotate takes the matrix Target and uses Clifford algebra to _concurrently_ rotate each of its rows by angle radians around axis. Axis must be a 3D row vector. Target must be an N,3 matrix.

func RotateSer Uses

func RotateSer(Target, ToRot, axis *v3.Matrix, angle float64) *v3.Matrix

RotateSer takes the matrix Target and uses Clifford algebra to rotate each of its rows by angle radians around axis. Axis must be a 3D row vector. Target must be an N,3 matrix. The Ser in the name is from "serial". ToRot will be overwritten and returned

func RotateSerP Uses

func RotateSerP(Target, ToRot, axis, tmpv, Rv, Rvrev, Rotatedv, itmp1, itmp2, itmp3, itmp4, itmp5, itmp6 *v3.Matrix, angle float64)

RotateSerP takes the matrix Target and uses Clifford algebra to rotate each of its rows by angle radians around axis. Axis must be a 3D row vector. Target must be an N,3 matrix. The Ser in the name is from "serial". ToRot will be overwritten and returned. RotateSerP only allocates some floats but not any v3.Matrix. Instead, it takes the needed intermediates as arguments, hence the "P" for "performance" If performance is not an issue, use RotateSer instead, it will perform the allocations for you and call this function. Notice that if you use this function directly you may have to zero at least some of the intermediates before reusing them.

func RotatorAroundZ Uses

func RotatorAroundZ(gamma float64) (*v3.Matrix, error)

RotatorAroundZ returns an operator that will rotate a set of coordinates by gamma radians around the z axis.

func RotatorAroundZToNewY Uses

func RotatorAroundZToNewY(newy *v3.Matrix) (*v3.Matrix, error)

RotatorToNewY takes a set of coordinates (mol) and a vector (y). It returns a rotation matrix that, when applied to mol, will rotate it around the Z axis in such a way that the projection of newy in the XY plane will be aligned with the Y axis.

func RotatorToNewZ Uses

func RotatorToNewZ(newz *v3.Matrix) *v3.Matrix

RotatorToNewZ takes a matrix a row vector (newz). It returns a linear operator such that, when applied to a matrix mol ( with the operator on the right side) it will rotate mol such that the z axis is aligned with newz.

func RotatorTranslatorToSuper Uses

func RotatorTranslatorToSuper(test, templa *v3.Matrix) (*v3.Matrix, *v3.Matrix, *v3.Matrix, *v3.Matrix, error)

RotatorTranslatorToSuper superimposes the set of cartesian coordinates given as the rows of the matrix test on the gnOnes of the rows of the matrix templa. Returns the transformed matrix, the rotation matrix, 2 translation row vectors For the superposition plus an error. In order to perform the superposition, without using the transformed the first translation vector has to be added first to the moving matrix, then the rotation must be performed and finally the second translation has to be added. This is a low level function, although one can use it directly since it returns the transformed matrix. The math for this function is by Prof. Veronica Jimenez-Curihual, University of Concepcion, Chile.

func ScaleBond Uses

func ScaleBond(C, H *v3.Matrix, bond float64)

ScaleBond takes a C-H bond and moves the H (in place) so the distance between them is the one given (bond). CAUTION: I have only tested it for the case where the original distance>bond, although I expect it to also work in the other case.

func ScaleBonds Uses

func ScaleBonds(coords *v3.Matrix, mol Atomer, n1, n2 string, finallenght float64)

ScaleBonds scales all bonds between atoms in the same residue with names n1, n2 to a final lenght finallengt, by moving the atoms n2. the operation is executed in place.

func SelCone Uses

func SelCone(B, selection *v3.Matrix, angle, distance, thickness, initial float64, whatcone int) []int

SelCone, Given a set of cartesian points in sellist, obtains a vector "plane" normal to the best plane passing through the points. It selects atoms from the set A that are inside a cone in the direction of "plane" that starts from the geometric center of the cartesian points, and has an angle of angle (radians), up to a distance distance. The cone is approximated by a set of radius-increasing cilinders with height thickness. If one starts from one given point, 2 cgnOnes, one in each direction, are possible. If whatcone is 0, both cgnOnes are considered. if whatcone<0, only the cone opposite to the plane vector direction. If whatcone>0, only the cone in the plane vector direction. the 'initial' argument allows the construction of a truncate cone with a radius of initial.

func Super Uses

func Super(test, templa *v3.Matrix, indexes ...[]int) (*v3.Matrix, error)

Super determines the best rotation and translations to superimpose the coords in test considering only the atoms present in the slices of int slices indexes. The first indexes slices will be assumed to contain test indexes and the second, template indexes. If you give only one (it must be the first one), it will be assumed to correspondo to whatever molecule that has more atoms than the elements in the slice. Giving a nil or 0-lenght first slice and a non-nil second slice will cause MemRMSD to not consider neither of them. The same number of atoms has to be considered for superposition in both systems. It applies those rotation and translations to the whole molecule test, in palce. testlst and templalst must have the same number of elements.

func TagAtomsByName Uses

func TagAtomsByName(r Atomer, name string, list []int) int

TagAtomsByName will tag all atoms with a given name in a given list of atoms. return the number of tagged atoms

func XYZFileWrite Uses

func XYZFileWrite(xyzname string, Coords *v3.Matrix, mol Atomer) error

XYZWrite writes the mol Ref and the Coord coordinates in an XYZ file with name xyzname which will be created fot that. If the file exist it will be overwritten.

func XYZStringWrite Uses

func XYZStringWrite(Coords *v3.Matrix, mol Atomer) (string, error)

XYZStringWrite writes the mol Ref and the Coord coordinates in an XYZ-formatted string.

func XYZWrite Uses

func XYZWrite(out io.Writer, Coords *v3.Matrix, mol Atomer) error

XYZStringWrite writes the mol Ref and the Coord coordinates in an XYZ-formatted string.

type Atom Uses

type Atom struct {
    Name      string  //PDB name of the atom
    ID        int     //The PDB index of the atom
    Tag       int     //Just added this for something that someone might want to keep that is not a float.
    Molname   string  //PDB name of the residue or molecule (3-letter code for residues)
    Molname1  byte    //the one letter name for residues and nucleotids
    Char16    byte    //Whatever is in the column 16 (counting from 0) in a PDB file, anything.
    MolID     int     //PDB index of the corresponding residue or molecule
    Chain     string  //One-character PDB name for a chain.
    Mass      float64 //hopefully all these float64 are not too much memory
    Occupancy float64 //a PDB crystallographic field, often used to store values of interest.
    Vdw       float64 //radius
    Charge    float64 //Partial charge on an atom
    Symbol    string
    Het       bool // is the atom an hetatm in the pdb file? (if applicable)
}

Atom contains the atoms read except for the coordinates, which will be in a matrix and the b-factors, which are in a separate slice of float64.

func (*Atom) Copy Uses

func (N *Atom) Copy(A *Atom)

Copy returns a copy of the Atom object. puts the copy into the

type AtomMultiCharger Uses

type AtomMultiCharger interface {
    Atomer

    //Charge gets the total charge of the topology
    Charge() int

    //Multi returns the multiplicity of the topology
    Multi() int
}

AtomChargerMultier is atomer but also gives a charge and multiplicity

type Atomer Uses

type Atomer interface {

    //Atom returns the Atom corresponding to the index i
    //of the Atom slice in the Topology. Should panic if
    //out of range.
    Atom(i int) *Atom

    Len() int
}

Atomer is the basic interface for a topology.

type CError Uses

type CError struct {
    // contains filtered or unexported fields
}

CError (Concrete Error) is the concrete error type for the chem package, that implements chem.Error

func (CError) Decorate Uses

func (err CError) Decorate(dec string) []string

Decorate will add the dec string to the decoration slice of strings of the error, and return the resulting slice.

func (CError) Error Uses

func (err CError) Error() string

type ConcTraj Uses

type ConcTraj interface {

    //Is the trajectory ready to be read?
    Readable() bool

    /*NextConc takes a slice of bools and reads as many frames as elements the list has
    form the trajectory. The frames are discarted if the corresponding elemetn of the slice
    is false. The function returns a slice of channels through each of each of which
    a *matrix.DenseMatrix will be transmited*/
    NextConc(frames []*v3.Matrix) ([]chan *v3.Matrix, error)

    //Returns the number of atoms per frame
    Len() int
}

ConcTraje is an interface for a trajectory that can be read concurrently.

type Error Uses

type Error interface {
    Error() string
    Decorate(string) []string //This is the new thing for errors. It allows you to add information when you pass it up. Each call also returns the "decoration" slice of strins resulting from the current call. If passed an empty string, it should just return the current value, not add the empty string to the slice.

}

Error is the interface for errors that all packages in this library implement. The Decorate method allows to add and retrieve info from the error, without changing it's type or wrapping it around something else.

type LastFrameError Uses

type LastFrameError interface {
    TrajError
    NormalLastFrameTermination() //does nothing, just to separate this interface from other TrajError's

}

LastFrameError has a useless function to distinguish the harmless errors (i.e. last frame) so they can be filtered in a typeswith that looks for this interface.

type Masser Uses

type Masser interface {

    //Returns a column vector with the massess of all atoms
    Masses() ([]float64, error)
}

Masser can return a slice with the masses of each atom in the reference.

type Molecule Uses

type Molecule struct {
    *Topology
    Coords   []*v3.Matrix
    Bfactors [][]float64
    // contains filtered or unexported fields
}

Molecule contains all the info for a molecule in many states. The info that is expected to change between states, Coordinates and b-factors are stored separately from other atomic info.

func NewMolecule Uses

func NewMolecule(coords []*v3.Matrix, ats Atomer, bfactors [][]float64) (*Molecule, error)

NewMolecule makes a molecule with ats atoms, coords coordinates, bfactors b-factors charge charge and unpaired unpaired electrons, and returns it. It doesnt check for consitency across slices or correct charge or unpaired electrons.

func PDBFileRead Uses

func PDBFileRead(pdbname string, read_additional bool) (*Molecule, error)

PDBFileRead reads the atomic entries for a PDB file, returns a bunch of without coordinates, and the coordinates in a separate array of arrays. If there is one frame in the PDB the coordinates array will be of lenght 1. It also returns an error which is not really well set up right now.

func PDBRead Uses

func PDBRead(pdb io.Reader, read_additional bool) (*Molecule, error)

PDBRRead reads a pdb file from an io.Reader. Returns a bunch of without coordinates, and the coordinates in a separate array of arrays. If there is one frame in the PDB the coordinates array will be of lenght 1. It also returns an error which is not really well set up right now.

func Reduce Uses

func Reduce(mol Atomer, coords *v3.Matrix, build int, report *os.File, executable string) (*Molecule, error)

*************** NOTE: This function breaks in Go=>1.5 for unknown reasons. There seem to be an issue with StdinPipe, or the program Reduce been never actually excecuted **************** Reduce uses the Reduce program (Word, et. al. (1999) J. Mol. Biol. 285, 1735-1747. For more information see http://kinemage.biochem.duke.edu) To protonate a protein and flip residues. It writes the report from Reduce to a file called Reduce.err in the current dir. The Reduce executable can be given, in "executable", otherwise it will be assumed that it is a file called "reduce" and it is in the PATH.

func XYZFileRead Uses

func XYZFileRead(xyzname string) (*Molecule, error)

XYZFileRead Reads an xyz or multixyz file (as produced by Turbomole). Returns a Molecule and error or nil.

func XYZRead Uses

func XYZRead(xyzp io.Reader) (*Molecule, error)

Reads an xyz or multixyz formatted bufio.Reader (as produced by Turbomole). Returns a Molecule and error or nil.

func (*Molecule) AddFrame Uses

func (M *Molecule) AddFrame(newframe *v3.Matrix)

AddFrame akes a matrix of coordinates and appends them at the end of the Coords. It checks that the number of coordinates matches the number of atoms.

func (*Molecule) AddManyFrames Uses

func (M *Molecule) AddManyFrames(newframes []*v3.Matrix)

AddManyFrames adds the array of matrices newfames to the molecule. It checks that the number of coordinates matches the number of atoms.

func (*Molecule) Coord Uses

func (M *Molecule) Coord(atom, frame int) *v3.Matrix

Coord returns the coords for the atom atom in the frame frame. panics if frame or coords are out of range.

func (*Molecule) Copy Uses

func (M *Molecule) Copy(A *Molecule)

Copy puts in the receiver a copy of the molecule A including coordinates

func (*Molecule) Corrupted Uses

func (M *Molecule) Corrupted() error

Corrupted checks whether the molecule is corrupted, i.e. the coordinates don't match the number of atoms. It also checks That the coordinate matrices have 3 columns.

func (*Molecule) Current Uses

func (M *Molecule) Current() int

Current returns the number of the next read frame

func (*Molecule) Del Uses

func (M *Molecule) Del(i int) error

Deletes atom i and its coordinates from the molecule.

func (*Molecule) DelCoord Uses

func (M *Molecule) DelCoord(i int) error

Deletes the coodinate i from every frame of the molecule.

func (*Molecule) InitRead Uses

func (M *Molecule) InitRead() error

Initializes molecule to be read as a traj (not tested!)

func (*Molecule) Less Uses

func (M *Molecule) Less(i, j int) bool

Less: Should the atom i be sorted before atom j?

func (*Molecule) NFrames Uses

func (M *Molecule) NFrames() int

NFrames returns the number of frames in the molecule

func (*Molecule) Next Uses

func (M *Molecule) Next(V *v3.Matrix) error

Returns the next frame and an error

func (*Molecule) NextConc Uses

func (M *Molecule) NextConc(frames []bool) ([]chan *v3.Matrix, error)

NextConc takes a slice of bools and reads as many frames as elements the list has form the trajectory. The frames are discarted if the corresponding elemetn of the slice * is false. The function returns a slice of channels through each of each of which * a *matrix.DenseMatrix will be transmited

func (*Molecule) Readable Uses

func (M *Molecule) Readable() bool

Checks that the molecule exists and has some existent Coordinates, in which case returns true. It returns false otherwise. The coordinates could still be empty

func (*Molecule) SetCurrent Uses

func (M *Molecule) SetCurrent(i int)

SetCurrent sets the value of the frame nex to be read to i.

func (*Molecule) Swap Uses

func (M *Molecule) Swap(i, j int)

Swap function, as demanded by sort.Interface. It swaps atoms, coordinates (all frames) and bfactors of the molecule.

type PanicMsg Uses

type PanicMsg string

PanicMsg is the type used for all the panics raised in the chem package so they can be easily recovered if needed. goChem only raises panics for programming errors: Attempts to to out of a matrice's bounds, dimension mismatches in matrices, etc.

func (PanicMsg) Error Uses

func (v PanicMsg) Error() string

Error returns a string with an error message

type RamaSet Uses

type RamaSet struct {
    Cprev   int
    N       int
    Ca      int
    C       int
    Npost   int
    MolID   int
    Molname string
}

func RamaList Uses

func RamaList(M Atomer, chains string, resran []int) ([]RamaSet, error)

RamaList takes a molecule and returns a slice of RamaSet, which contains the indexes for each dihedral to be included in a Ramachandran plot. It gets the dihedral indices for all residues in the range resran. if resran has 2 elements defining the boundaries. Otherwise, returns dihedral lists for the residues included in resran. If resran has 2 elements and the last is -1, RamaList will get all the dihedral for residues from resran[0] to the end of the chain. It only obtain dihedral lists for residues belonging to a chain included in chains.

func RamaResidueFilter Uses

func RamaResidueFilter(dihedrals []RamaSet, filterdata []string, shouldBePresent bool) ([]RamaSet, []int)

Filter the set of dihedral angles of a ramachandran plot by residue.(ex. only GLY, everything but GLY) The 3 letter code of the residues to be filtered in or out is in filterdata, whether they are filter in or out depends on shouldBePresent. It returns the filtered data and a slice containing the indexes in the new data of the residues in the old data, when they are included, or -1 when they are not included.

type Topology Uses

type Topology struct {
    Atoms []*Atom
    // contains filtered or unexported fields
}

Topology contains information about a molecule which is not expected to change in time (i.e. everything except for coordinates and b-factors)

func MergeAtomers Uses

func MergeAtomers(A, B Atomer) *Topology

Merges A and B in a single topology which is returned

func NewTopology Uses

func NewTopology(charge, multi int, ats ...[]*Atom) *Topology

NewTopology returns topology with ats atoms charge charge and multi multiplicity. It doesnt check for consitency across slices or correct charge or unpaired electrons.

func (*Topology) AppendAtom Uses

func (T *Topology) AppendAtom(at *Atom)

AppendAtom appends an atom at the end of the reference

func (*Topology) Atom Uses

func (T *Topology) Atom(i int) *Atom

Atom returns the Atom corresponding to the index i of the Atom slice in the Topology. Panics if out of range.

func (*Topology) Charge Uses

func (T *Topology) Charge() int

Charge returns the total charge of the topology

func (*Topology) CopyAtoms Uses

func (T *Topology) CopyAtoms(A Atomer)

Copy atoms into a topology. This is a deep copy, so T must have at least as many atoms as A.

func (*Topology) DelAtom Uses

func (T *Topology) DelAtom(i int)

DelAtom Deletes atom i by reslicing. This means that the copy still uses as much memory as the original T.

func (*Topology) Len Uses

func (T *Topology) Len() int

Len returns the length of the molecule.

func (*Topology) Masses Uses

func (T *Topology) Masses() ([]float64, error)

MassCol returns a slice of float64 with the masses of atoms, or nil and an error if they have not been calculated

func (*Topology) Multi Uses

func (T *Topology) Multi() int

Unpaired returns the multiplicity in the topology

func (*Topology) ResetIDs Uses

func (T *Topology) ResetIDs()

Sets the current order of atoms as Id and the order of molecules as MolID for all atoms

func (*Topology) SetAtom Uses

func (T *Topology) SetAtom(i int, at *Atom)

SetAtom sets the (i+1)th Atom of the topology to aT. Panics if out of range

func (*Topology) SetCharge Uses

func (T *Topology) SetCharge(i int)

SetCharge sets the total charge of the topology to i

func (*Topology) SetMulti Uses

func (T *Topology) SetMulti(i int)

SetUnpaired sets the multiplicity in the topology to i

func (*Topology) SomeAtoms Uses

func (R *Topology) SomeAtoms(T Atomer, atomlist []int)

SelectAtoms puts the atoms of T with indexes in atomlist into the receiver.

func (*Topology) SomeAtomsSafe Uses

func (R *Topology) SomeAtomsSafe(T Atomer, atomlist []int) error

SelectAtoms puts the atoms of T with indexes in atomlist into the receiver. Returns error if problem.

type Traj Uses

type Traj interface {

    //Is the trajectory ready to be read?
    Readable() bool

    //reads the next frame and returns it as DenseMatrix if keep==true, or discards it if false
    Next(output *v3.Matrix) error

    //Returns the number of atoms per frame
    Len() int
}

Traj is an interface for any trajectory object, including a Molecule Object

type TrajError Uses

type TrajError interface {
    Error
    Critical() bool
    FileName() string
    Format() string
}

TrajError is the nterface for errors in trajectories

Directories

PathSynopsis
chemjson
chemplot
dcd
qm
v3Package v3 implements a Matrix type representing a row-major 3D matrix (i.e.
xtc

Package chem imports 12 packages (graph) and is imported by 15 packages. Updated 2018-01-31. Refresh now. Tools for package owners.