RIF Documentation¶

Changelog¶

2017-08-22¶

added Changelog!
added Rif Package autodoc index
started on documenting the _rif extension module
- _rif._index._stripe_index.stripe_index_3d_Atom
added docs for rif.index.stripe_index_3d()

README¶

master

https://img.shields.io/travis/willsheffler/rif.svg

https://img.shields.io/codecov/c/github/willsheffler/rif.svg

devel

tell willsheffler@gmail.com to put something in this readme!¶

http://rif.readthedocs.io

INSTALL¶

universal¶

pip install -rrequirements.txt python tools/build_and_test.py

or

pip install . python3 -mpytest –pyargs rif

or

pip install -rrequirements.txt mkdir build && cd build CXXFLAGS=-Ipath-to-site-packages/numpy/core/include cmake .. PYTHONPATH=. python3 -mpytest –pyargs rif

crappy details¶

python requirements numpy, jinja2, pytest, pytest-cpp*, hypothesis, pandas

*use the pytest_cpp in external: “cd <rifdir>/external/hacked_packages/pytest-cpp-0.4 && pip install .”

for graphics requirements pymol + pip packages: PyOpenGL OpenGLContext PyVRML97 pydispatcher

SUBLIME SETUP¶

I use sublime plugins Anaconda, EasyClangComplete, Git, C++11, GitGutter, SublimeLinter

NOTES¶

iff you have pyrosetta4 on your path, functions using pyrosetta will be used/tested iff you have pymol on your path, functions using pymol will be used/tested

docs¶

install doxygen (–with-libclang probably best…) pip3 install sphinx_rtd_theme breathe

ubuntu 16.04¶

# better to use a virtualenv, no sudo # this ppa is to install gcc-6 sudo add-apt-repository ppa:ubuntu-toolchain-r/test # optional, for gcc6 sudo apt update sudo sudo apt install git git clone git@github.com:willsheffler/rif.git sudo apt-get install python-pip python3-pip cmake ninja-build clang g++-6 libboost-system-dev libboost-iostreams-dev sudo apt-get install clang g++-6 # optional

mac (tested on 10.10 and 10.12)¶

brew install python cmake ninja boost

CONTRIBUTING¶

conding conventions?¶

clang-format and pep8 everything please (hook it to your editor)
python style as laid out in Fluent Python
use libraries hdf5, numpy, pandas

Build System¶

This project should support

python setup.py
pip
cmake (needs help finding a numpy library)
tools/build_and_test.py
- should be re-factored to use click (or other)
- should dispatch on file pattern (use a library?) *.pybind.cpp *_test.py

Adding wrapped c++ code¶

This section needs an audit.

Create a file called path/to/my_module.pybind.cpp which defines a function

void RIFLIB_PYBIND_something_unique_path_is_good(py::module &m) {
    // pybind11 code here, adding to py::module m
}

http://pybind11.readthedocs.io/en/master/

This function will be auto-detected by the build system and called with a pybind module rif.path.to.my_module, i.e. module namespaces based on the file path. If your filename starts in caps: path/to/MyClass.pybind.cpp, the module m will instead be rif.path.to. That is, build system will assume this is a class that you will bind and should not be a module namespace.

contributors¶

will sheffler
alex ford (name and advice)
david baker (funding)

Introduction¶

RIF is a toolkit for protein structure modeling built around numpy. Name originally an acronyn for (R)otamer (I)nteraction (F)ield as the toolkit is growing from these ideas.

nothing to see here… yell at will. maybe you want the Legacy docs?

Concepts¶

Implementation¶

Examples¶

high level¶

fiber docking¶

This is prototype code to dock a monomer into helical protein fiber. Xform from \(B_0\) to \(B_N\) is sampled via hierarchical search. Possible placements a third body \(B_M\) are checked based on possible \(M/N\) roots of the transform \(X_N\):

\[X_M = X_N^{M/N}\]

"""
a prototype module for fiber docking. for illustration see:

https://docs.google.com/presentation/d/1qv6yZIy0b2dfmZAzA7bFNJD4pl9h7Nt1mOL4aqiVpOw/pub?start=false&loop=false&delayms=3000#slide=id.g1d6846b840_0_8

(see prev slide for simpler, less-general case)
"""

# from rif.io.pdbio import read_pdb
# from rif.nest import NestedXforms
# from rif.models import ImplicitEnsemble, ImplicitEnsembleScore
# from rif.search import HierBeamSearch

# import numpy as np


def rel_prime_flat_tri(Nmax):
    """return two same-size lists of relatively prime numbers

    mapped to tuples, would be: 1,2 1,3 2,3 1,4 3,4 1,5 2,5 3,5 4,5 ...
    replace with numpy one-liner, probably more clear
    """
    pass


def center_welzl(model):
    """center model based on welzl optimal bounding sphere, report radius"""
    pass


def fiber_dock(pdb_file):
    """turn monomer into helical fiber in "all" possible ways"""

    # setup
    resolutions = (16, 8, 4, 2, 1, 0.5)
    M, N = rel_prime_flat_tri(Nmax=20)
    powers = (M / N)[np.newaxis, :]  # shape (1, whatever)
    assert all(powers < 1)  # else Hsearch can't work
    prot, radius = center_welzl(read_pdb(pdb_file))
    bbframes = prot.bbframes(ss='HE', burial='exposed')
    samp_grid = NestedXforms(baseresl=resolutions[0], lever=radius)
    # ImplicitEnsemble represents a structural ensemble in which each member's
    # elements (atoms, rays, stubs, etc) are within well-defined translational
    # and rotational bounds of a reference structure. (rotation mag defined by
    # radius/lever)
    ensemble0 = ImplicitEnsemble(atoms=prot.atoms, rotslots=bbframes,
                                 radius=0.0, lever=radius, bvh_max_radius=1.0,
                                 voxel_convolv_radii=resolutions)
    hsearch = HierBeamSearch(resls=resolutions, grid=samp_grid, beam_size_m=10,
                             extra_fields=(('M', 'u1'), ('N', 'u1')), )
    score = ImplicitEnsembleScore(
        # todo
    )

    # assuming I can do what I want with numpy dtype/operators...
    def score_samples(stage, samples):
        resl = resolutions[stage]
        indices = samples['index']
        positions = samp_grid.position_of(stage, indices)
        ensembleN = positions * ensemble0.with_radius(resl)
        score0N = score(ensemble0, ensembleN)  # radius 0 + resl
        # better to extract sub-array of score0N < 0?
        if np.any(score0N < 0):
            # in addition to ensemble['position']**(M/N), operator
            # __pow__(ensemble,float) must scale ensemble['radius']
            # linearly by pow... is this well-defined algebra?
            # maybe just use an ensemble_power() func directly?
            ensembleM = np.power(ensembleN, powers)
            # radius 0 + resl * m / n (not constant over array)
            score0M = score(ensemble0, ensembleM)
            # - coupling='linear' means M ensemble motion tied exactly to
            #   N motion, so instead of radius = Arad + Brad,
            #   radius = min(Arad - Brad, Brad - Arad)
            # - this check may not be worth the cost here, skip or
            #   maybe filter on ensembleM good enough first
            scoreMN = score(ensembleM, ensembleN, coupling='linear')
            score0MN = score0N + score0M + scoreMN  # score0N bcast
            idxmin = score0MN.argmin(axis=2)  # best score forall M / N
            samples['M'] = M[idxmin]
            samples['N'] = N[idxmin]
            samples['score'] = score0MN[:, idxmin]
            # todo: eventually check all good M/N? will probably be rare
            #       to have more than one hit
            # todo: also score ensemble0 vs all 1..N) for thoroughness
        return samples

    # hsearch interface option 1
    # simple, but perhaps not as flexible as below
    hsearch.search(score_samples)
    print(hsearch.results())

    # hsearch interface option 2
    # not sure if using iteration as an interface to the search is a
    # good idea or not... note that samples['score'] gets communicated
    # back to the search, which is not a pattern I've seen used much
    # but it would be clearer for simpler cases... don't actually have to
    # define a function (ie score_samples)
    hsearch_accum = hsearch.accumulator()
    for iresl, samples in hsearch_accum:
        score_samples(iresl, samples)
    print(hsearch_accum.results())


if __name__ == '__main__':
    fiber_dock()

Rif Module Documentation¶

Subpackages¶

rif.apps package¶

Subpackages¶

rif.apps.examples package¶

Subpackages¶

rif.apps.examples.pyrosetta package¶

Submodules¶

rif.apps.examples.pyrosetta.170407_david_rbsample module¶

Module contents¶

Submodules¶

rif.apps.examples.170430_david_patch module¶

rif.apps.examples.db_helix_asp module¶

rif.apps.examples.db_stub_index module¶

rif.apps.examples.db_stub_index_test module¶

rif.apps.examples.geom_examples module¶

rif.apps.examples.glumpy_ex module¶

Module contents¶

docstring for rif.apps.examples

Submodules¶

rif.apps.cyclic_cyclic_dock module¶

rif.apps.fiber_dock module¶

a prototype module for fiber docking. for illustration see:

https://docs.google.com/presentation/d/1qv6yZIy0b2dfmZAzA7bFNJD4pl9h7Nt1mOL4aqiVpOw/pub?start=false&loop=false&delayms=3000#slide=id.g1d6846b840_0_8

(see prev slide for simpler, less-general case)

rif.apps.fiber_dock.center_welzl(model)[source]¶: center model based on welzl optimal bounding sphere, report radius

rif.apps.fiber_dock.fiber_dock(pdb_file)[source]¶: turn monomer into helical fiber in “all” possible ways

rif.apps.fiber_dock.rel_prime_flat_tri(Nmax)[source]¶

return two same-size lists of relatively prime numbers

mapped to tuples, would be: 1,2 1,3 2,3 1,4 3,4 1,5 2,5 3,5 4,5 … replace with numpy one-liner, probably more clear

rif.apps.numpy_experiments module¶

Module contents¶

docstring for rif.apps

rif.geom package¶

Submodules¶

rif.geom.Ray module¶

rif.geom.Ray_test module¶

rif.geom.rand_geom module¶

rif.geom.rand_geom_test module¶

rif.geom.ray_hash module¶

rif.geom.ray_hash_error_test module¶

rif.geom.ray_hash_test module¶

Module contents¶

docstring for rif.geom

rif.hash package¶

Submodules¶

rif.hash.xform_hash_6d module¶

rif.hash.xform_hash_6d_test module¶

Module contents¶

docstring for rif.hash

rif.index package¶

Module contents¶

docstring for rif.index

rif.index.stripe_index_3d(maxdis, points, payloads=None)[source]¶

create a 3D stripe index of points[payloads] based on dtype of inputs

stripe index is a fast index of 3D coords (possible enhanced with metadata) which can quickly get all the neighbors within a certain maxdis, which must be specified at construction

Parameters:

maxdis (float) – maximum distance for points to be neighbors
points (array) – can map to array of C++ structs via dtype
payloads (array) – optional, additional data attached to each point

Returns:

Wrapper, for example, see the *_V3 or *_V3_V3 versions

Points only (*_V3) neighbors method returns points
With payloads (*_V3_V3) neighbors method returns payloads (only)

Return type:

stripe_index_3d_<points_t>_<payload_t=NULL>

Return value notes:: If payloads is a sequence of python objects, internally integers 1-N are used as the payloads, and the objects in the sequence are returned (not recommended, you should be using array style computing).

Examples

Example: simple points with python objects

>>> array = np.random.uniform(size=(1000,3))
>>> points = [(x[0], x[1], x[2]) for x in array]
>>> type(points[13])  # converted to tuples
<class 'tuple'>
>>> index = stripe_index_3d(maxdis=0.1, points=points)
>>> index
<_rif._index._stripe_index.stripe_index_3d_V3 object at ...>
>>> test = (0.1, 0.2, 0.3)
>>> neighbors_fast = index.neighbor_count(test)
>>> neighbors_slow = index.neighbor_count_brute(test)
>>> # speedup for this example would be nothing
>>> # because python operations dominate runtime

Example: simple points with numpy

>>> points = np.random.uniform(size=30000).astype('f').view(V3)
>>> points.shape
(10000,)
>>> index = stripe_index_3d(maxdis=0.05, points=points)
>>> index
<_rif._index._stripe_index.stripe_index_3d_V3 object at ...>
>>> testpts = np.random.uniform(low=-1, high=2, size=3000).astype('f').view(V3)
>>> tfast = clock()
>>> neighbors_fast = index.neighbor_count(testpts)
>>> tfast = clock() - tfast
>>> tslow = clock()
>>> neighbors_slow = index.neighbor_count_brute(testpts)
>>> tslow = clock() - tslow
>>> print('avg number of neighbors:', np.mean(neighbors_fast))  
avg number of neighbors: 0.881
>>> np.all(neighbors_fast == neighbors_slow)
True
>>> print('speedup is %fx vs brute force loop over points' % (tslow/tfast)) 
speedup is 81.904762x vs brute force loop over points
>>> # this example is faster because it uses a numpy vectorized query
>>> # this amortizes away the cost of python operations

Example: mapping points to arbitrary python payloads

>>> points = [(0,0,0), (1,1,1), (2,2,2), (3,3,3)]
>>> values = ['string payloads: ' + str(x) for x in points]
>>> index = stripe_index_3d(maxdis=1.0, points=points, payloads=values)
>>> index.neighbors((0,0,0))  # returns python objects
['string payloads: (0, 0, 0)']
>>> index.neighbors((0.5,0.5,0.5))
['string payloads: (0, 0, 0)', 'string payloads: (1, 1, 1)']
>>> # this example will also be slow... no reason to use this over brute force

Example: mapping point-like structs to structs via numpy dtypes

>>> N = 1000
>>> points = (np.arange(3*N).astype('f')/N).view(V3)
>>> atoms = np.empty(N, dtype=Atom)
>>> atoms['pos'] = points
>>> atoms['atype'] = np.arange(N) % 23
>>> atoms['rtype'] = np.arange(N) % 20
>>> atoms['anum'] = np.arange(N) % 7
>>> atoms[-3:]
array([(([ 2.99099994,  2.9920001 ,  2.99300003],),  8, 3, 17),
       (([ 2.99399996,  2.99499989,  2.99600005],),  9, 4, 18),
       (([ 2.99699998,  2.99799991,  2.99900007],), 10, 5, 19)],
      dtype=[('pos', [('raw', '<f4', (3,))]), ('atype', 'u1'), ('anum', 'u1'), ('rtype', '<u2')])
>>> index = stripe_index_3d(maxdis=0.05, points=points, payloads=atoms)
>>> index.neighbors((1, 1, 1))
array([(([ 0.972     ,  0.97299999,  0.97399998],),  2, 2,  4),
       (([ 0.97500002,  0.97600001,  0.977     ],),  3, 3,  5),
       (([ 0.97799999,  0.97899997,  0.98000002],),  4, 4,  6),
       (([ 0.98100001,  0.98199999,  0.98299998],),  5, 5,  7),
       (([ 0.98400003,  0.98500001,  0.986     ],),  6, 6,  8),
       (([ 0.98699999,  0.98799998,  0.98900002],),  7, 0,  9),
       (([ 0.99000001,  0.991     ,  0.99199998],),  8, 1, 10),
       (([ 0.99299997,  0.99400002,  0.995     ],),  9, 2, 11),
       (([ 0.99599999,  0.99699998,  0.99800003],), 10, 3, 12),
       (([ 0.99900001,  1.        ,  1.00100005],), 11, 4, 13),
       (([ 1.00199997,  1.00300002,  1.00399995],), 12, 5, 14),
       (([ 1.005     ,  1.00600004,  1.00699997],), 13, 6, 15),
       (([ 1.00800002,  1.00899994,  1.00999999],), 14, 0, 16),
       (([ 1.01100004,  1.01199996,  1.01300001],), 15, 1, 17),
       (([ 1.01400006,  1.01499999,  1.01600003],), 16, 2, 18),
       (([ 1.01699996,  1.01800001,  1.01900005],), 17, 3, 19),
       (([ 1.01999998,  1.02100003,  1.02199996],), 18, 4,  0),
       (([ 1.023     ,  1.02400005,  1.02499998],), 19, 5,  1),
       (([ 1.02600002,  1.02699995,  1.028     ],), 20, 6,  2)],
      dtype=[('pos', [('raw', '<f4', (3,))]), ('atype', 'u1'), ('anum', 'u1'), ('rtype', '<u2')])

C++ Module contents¶

_rif._index package¶

Submodules¶

rif.index.stripe_index module¶

_rif._index._stripe_index.make_huge_array(arg0: int) → array¶

class _rif._index._stripe_index.stripe_index_3d_Atom¶

Bases: pybind11_builtins.pybind11_object

Extension Class wrapping stripe_index template with Point typeYELL_AT_WILL

neighbor_count(*args, **kwargs)¶

Overloaded function.

neighbor_count(self: _rif._index._stripe_index.stripe_index_3d_Atom, query: numpy.ndarray[_rif._actor.Atom]) -> object

returns number of neighbors for input points

neighbor_count(self: _rif._index._stripe_index.stripe_index_3d_Atom, arg0: object) -> int

neighbor_count_brute(*args, **kwargs)¶

Overloaded function.

neighbor_count_brute(self: _rif._index._stripe_index.stripe_index_3d_Atom, query: numpy.ndarray[_rif._actor.Atom]) -> object

brute force version of neighbor_count

neighbor_count_brute(self: _rif._index._stripe_index.stripe_index_3d_Atom, arg0: object) -> int

neighbor_exists(*args, **kwargs)¶

Overloaded function.

neighbor_exists(self: _rif._index._stripe_index.stripe_index_3d_Atom, query: numpy.ndarray[_rif._actor.Atom]) -> object

return array true if a neighbor exists for corresponding input point

neighbor_exists(self: _rif._index._stripe_index.stripe_index_3d_Atom, arg0: object) -> bool

neighbor_exists_brute(*args, **kwargs)¶

Overloaded function.

neighbor_exists_brute(self: _rif._index._stripe_index.stripe_index_3d_Atom, query: numpy.ndarray[_rif._actor.Atom]) -> object

brute force version of neighbor_exists

neighbor_exists_brute(self: _rif._index._stripe_index.stripe_index_3d_Atom, arg0: object) -> int

neighbors(*args, **kwargs)¶

Overloaded function.

neighbors(self: _rif._index._stripe_index.stripe_index_3d_Atom, query: _rif._actor.Atom) -> array

return neighboring points for input point(s)? as custom JaggedArray

neighbors(self: _rif._index._stripe_index.stripe_index_3d_Atom, query: object) -> array

return neighboring points for input point

neighbors_brute(*args, **kwargs)¶

Overloaded function.

neighbors_brute(self: _rif._index._stripe_index.stripe_index_3d_Atom, arg0: _rif._actor.Atom) -> array
neighbors_brute(self: _rif._index._stripe_index.stripe_index_3d_Atom, arg0: object) -> array

translation¶: translation applied to internal storage (so that it is compact)

class _rif._index._stripe_index.stripe_index_3d_Atom_object¶

Bases: pybind11_builtins.pybind11_object

Extension Class wrapping stripe_index template mapping YELL_AT_WILLto payload YELL_AT_WILL

neighbor_count(*args, **kwargs)¶

Overloaded function.

neighbor_count(self: _rif._index._stripe_index.stripe_index_3d_Atom_object, query: numpy.ndarray[_rif._actor.Atom]) -> object

returns number of neighbors for input points

neighbor_count(self: _rif._index._stripe_index.stripe_index_3d_Atom_object, arg0: object) -> int

neighbor_count_brute(*args, **kwargs)¶

Overloaded function.

neighbor_count_brute(self: _rif._index._stripe_index.stripe_index_3d_Atom_object, query: numpy.ndarray[_rif._actor.Atom]) -> object

brute force version of neighbor_count

neighbor_count_brute(self: _rif._index._stripe_index.stripe_index_3d_Atom_object, arg0: object) -> int

neighbor_exists(*args, **kwargs)¶

Overloaded function.

neighbor_exists(self: _rif._index._stripe_index.stripe_index_3d_Atom_object, query: numpy.ndarray[_rif._actor.Atom]) -> object

return array true if a neighbor exists for corresponding input point

neighbor_exists(self: _rif._index._stripe_index.stripe_index_3d_Atom_object, arg0: object) -> bool

neighbor_exists_brute(*args, **kwargs)¶

Overloaded function.

neighbor_exists_brute(self: _rif._index._stripe_index.stripe_index_3d_Atom_object, query: numpy.ndarray[_rif._actor.Atom]) -> object

brute force version of neighbor_exists

neighbor_exists_brute(self: _rif._index._stripe_index.stripe_index_3d_Atom_object, arg0: object) -> int

neighbors(*args, **kwargs)¶

Overloaded function.

neighbors(self: _rif._index._stripe_index.stripe_index_3d_Atom_object, query: _rif._actor.Atom) -> array

return neighboring *payloads* for input point

neighbors(self: _rif._index._stripe_index.stripe_index_3d_Atom_object, arg0: object) -> array

neighbors_brute(*args, **kwargs)¶

Overloaded function.

neighbors_brute(self: _rif._index._stripe_index.stripe_index_3d_Atom_object, arg0: _rif._actor.Atom) -> array
neighbors_brute(self: _rif._index._stripe_index.stripe_index_3d_Atom_object, arg0: object) -> array

translation¶: translation applied to internal storage (so that it is compact)

class _rif._index._stripe_index.stripe_index_3d_V3¶

Bases: pybind11_builtins.pybind11_object

Extension Class wrapping stripe_index template with Point typeYELL_AT_WILL

neighbor_count(*args, **kwargs)¶

Overloaded function.

neighbor_count(self: _rif._index._stripe_index.stripe_index_3d_V3, query: numpy.ndarray[_rif._eigen_types.V3]) -> object

returns number of neighbors for input points

neighbor_count(self: _rif._index._stripe_index.stripe_index_3d_V3, arg0: object) -> int

neighbor_count_brute(*args, **kwargs)¶

Overloaded function.

neighbor_count_brute(self: _rif._index._stripe_index.stripe_index_3d_V3, query: numpy.ndarray[_rif._eigen_types.V3]) -> object

brute force version of neighbor_count

neighbor_count_brute(self: _rif._index._stripe_index.stripe_index_3d_V3, arg0: object) -> int

neighbor_exists(*args, **kwargs)¶

Overloaded function.

neighbor_exists(self: _rif._index._stripe_index.stripe_index_3d_V3, query: numpy.ndarray[_rif._eigen_types.V3]) -> object

return array true if a neighbor exists for corresponding input point

neighbor_exists(self: _rif._index._stripe_index.stripe_index_3d_V3, arg0: object) -> bool

neighbor_exists_brute(*args, **kwargs)¶

Overloaded function.

neighbor_exists_brute(self: _rif._index._stripe_index.stripe_index_3d_V3, query: numpy.ndarray[_rif._eigen_types.V3]) -> object

brute force version of neighbor_exists

neighbor_exists_brute(self: _rif._index._stripe_index.stripe_index_3d_V3, arg0: object) -> int

neighbors(*args, **kwargs)¶

Overloaded function.

neighbors(self: _rif._index._stripe_index.stripe_index_3d_V3, query: _rif._eigen_types.V3) -> array

return neighboring points for input point(s)? as custom JaggedArray

neighbors(self: _rif._index._stripe_index.stripe_index_3d_V3, query: object) -> array

return neighboring points for input point

neighbors_brute(*args, **kwargs)¶

Overloaded function.

neighbors_brute(self: _rif._index._stripe_index.stripe_index_3d_V3, arg0: _rif._eigen_types.V3) -> array
neighbors_brute(self: _rif._index._stripe_index.stripe_index_3d_V3, arg0: object) -> array

translation¶: translation applied to internal storage (so that it is compact)

class _rif._index._stripe_index.stripe_index_3d_V3_Atom¶

Bases: pybind11_builtins.pybind11_object

Extension Class wrapping stripe_index template mapping YELL_AT_WILLto payload YELL_AT_WILL

neighbor_count(*args, **kwargs)¶

Overloaded function.

neighbor_count(self: _rif._index._stripe_index.stripe_index_3d_V3_Atom, query: numpy.ndarray[_rif._eigen_types.V3]) -> object

returns number of neighbors for input points

neighbor_count(self: _rif._index._stripe_index.stripe_index_3d_V3_Atom, arg0: object) -> int

neighbor_count_brute(*args, **kwargs)¶

Overloaded function.

neighbor_count_brute(self: _rif._index._stripe_index.stripe_index_3d_V3_Atom, query: numpy.ndarray[_rif._eigen_types.V3]) -> object

brute force version of neighbor_count

neighbor_count_brute(self: _rif._index._stripe_index.stripe_index_3d_V3_Atom, arg0: object) -> int

neighbor_exists(*args, **kwargs)¶

Overloaded function.

neighbor_exists(self: _rif._index._stripe_index.stripe_index_3d_V3_Atom, query: numpy.ndarray[_rif._eigen_types.V3]) -> object

return array true if a neighbor exists for corresponding input point

neighbor_exists(self: _rif._index._stripe_index.stripe_index_3d_V3_Atom, arg0: object) -> bool

neighbor_exists_brute(*args, **kwargs)¶

Overloaded function.

neighbor_exists_brute(self: _rif._index._stripe_index.stripe_index_3d_V3_Atom, query: numpy.ndarray[_rif._eigen_types.V3]) -> object

brute force version of neighbor_exists

neighbor_exists_brute(self: _rif._index._stripe_index.stripe_index_3d_V3_Atom, arg0: object) -> int

neighbors(*args, **kwargs)¶

Overloaded function.

neighbors(self: _rif._index._stripe_index.stripe_index_3d_V3_Atom, query: _rif._eigen_types.V3) -> array

return neighboring *payloads* for input point

neighbors(self: _rif._index._stripe_index.stripe_index_3d_V3_Atom, arg0: object) -> array

neighbors_brute(*args, **kwargs)¶

Overloaded function.

neighbors_brute(self: _rif._index._stripe_index.stripe_index_3d_V3_Atom, arg0: _rif._eigen_types.V3) -> array
neighbors_brute(self: _rif._index._stripe_index.stripe_index_3d_V3_Atom, arg0: object) -> array

translation¶: translation applied to internal storage (so that it is compact)

class _rif._index._stripe_index.stripe_index_3d_V3_V3¶

Bases: pybind11_builtins.pybind11_object

Extension Class wrapping stripe_index template mapping YELL_AT_WILLto payload YELL_AT_WILL

neighbor_count(*args, **kwargs)¶

Overloaded function.

neighbor_count(self: _rif._index._stripe_index.stripe_index_3d_V3_V3, query: numpy.ndarray[_rif._eigen_types.V3]) -> object

returns number of neighbors for input points

neighbor_count(self: _rif._index._stripe_index.stripe_index_3d_V3_V3, arg0: object) -> int

neighbor_count_brute(*args, **kwargs)¶

Overloaded function.

neighbor_count_brute(self: _rif._index._stripe_index.stripe_index_3d_V3_V3, query: numpy.ndarray[_rif._eigen_types.V3]) -> object

brute force version of neighbor_count

neighbor_count_brute(self: _rif._index._stripe_index.stripe_index_3d_V3_V3, arg0: object) -> int

neighbor_exists(*args, **kwargs)¶

Overloaded function.

neighbor_exists(self: _rif._index._stripe_index.stripe_index_3d_V3_V3, query: numpy.ndarray[_rif._eigen_types.V3]) -> object

return array true if a neighbor exists for corresponding input point

neighbor_exists(self: _rif._index._stripe_index.stripe_index_3d_V3_V3, arg0: object) -> bool

neighbor_exists_brute(*args, **kwargs)¶

Overloaded function.

neighbor_exists_brute(self: _rif._index._stripe_index.stripe_index_3d_V3_V3, query: numpy.ndarray[_rif._eigen_types.V3]) -> object

brute force version of neighbor_exists

neighbor_exists_brute(self: _rif._index._stripe_index.stripe_index_3d_V3_V3, arg0: object) -> int

neighbors(*args, **kwargs)¶

Overloaded function.

neighbors(self: _rif._index._stripe_index.stripe_index_3d_V3_V3, query: _rif._eigen_types.V3) -> array

return neighboring *payloads* for input point

neighbors(self: _rif._index._stripe_index.stripe_index_3d_V3_V3, arg0: object) -> array

neighbors_brute(*args, **kwargs)¶

Overloaded function.

neighbors_brute(self: _rif._index._stripe_index.stripe_index_3d_V3_V3, arg0: _rif._eigen_types.V3) -> array
neighbors_brute(self: _rif._index._stripe_index.stripe_index_3d_V3_V3, arg0: object) -> array

translation¶: translation applied to internal storage (so that it is compact)

class _rif._index._stripe_index.stripe_index_3d_V3_object¶

Bases: pybind11_builtins.pybind11_object

Extension Class wrapping stripe_index template mapping YELL_AT_WILLto payload YELL_AT_WILL

neighbor_count(*args, **kwargs)¶

Overloaded function.

neighbor_count(self: _rif._index._stripe_index.stripe_index_3d_V3_object, query: numpy.ndarray[_rif._eigen_types.V3]) -> object

returns number of neighbors for input points

neighbor_count(self: _rif._index._stripe_index.stripe_index_3d_V3_object, arg0: object) -> int

neighbor_count_brute(*args, **kwargs)¶

Overloaded function.

neighbor_count_brute(self: _rif._index._stripe_index.stripe_index_3d_V3_object, query: numpy.ndarray[_rif._eigen_types.V3]) -> object

brute force version of neighbor_count

neighbor_count_brute(self: _rif._index._stripe_index.stripe_index_3d_V3_object, arg0: object) -> int

neighbor_exists(*args, **kwargs)¶

Overloaded function.

neighbor_exists(self: _rif._index._stripe_index.stripe_index_3d_V3_object, query: numpy.ndarray[_rif._eigen_types.V3]) -> object

return array true if a neighbor exists for corresponding input point

neighbor_exists(self: _rif._index._stripe_index.stripe_index_3d_V3_object, arg0: object) -> bool

neighbor_exists_brute(*args, **kwargs)¶

Overloaded function.

neighbor_exists_brute(self: _rif._index._stripe_index.stripe_index_3d_V3_object, query: numpy.ndarray[_rif._eigen_types.V3]) -> object

brute force version of neighbor_exists

neighbor_exists_brute(self: _rif._index._stripe_index.stripe_index_3d_V3_object, arg0: object) -> int

neighbors(*args, **kwargs)¶

Overloaded function.

neighbors(self: _rif._index._stripe_index.stripe_index_3d_V3_object, query: _rif._eigen_types.V3) -> array

return neighboring *payloads* for input point

neighbors(self: _rif._index._stripe_index.stripe_index_3d_V3_object, arg0: object) -> array

neighbors_brute(*args, **kwargs)¶

Overloaded function.

neighbors_brute(self: _rif._index._stripe_index.stripe_index_3d_V3_object, arg0: _rif._eigen_types.V3) -> array
neighbors_brute(self: _rif._index._stripe_index.stripe_index_3d_V3_object, arg0: object) -> array

translation¶: translation applied to internal storage (so that it is compact)

Module contents¶

rif.legacy package¶

Submodules¶

rif.legacy.xyzMath module¶

Easy 3D Linear Algebra, like xyz* in rosetta

class rif.legacy.xyzMath.Mat(xx=None, xy=None, xz=None, yx=None, yy=None, yz=None, zx=None, zy=None, zz=None)[source]¶

docstring for Mat

>>> m = Mat(2,0,0,0,1,0,0,0,1)
>>> print(m)
Mat[ (2.000000,0.000000,0.000000), (0.000000,1.000000,0.000000), (0.000000,0.000000,1.000000) ]
>>> print(m*m)
Mat[ (4.000000,0.000000,0.000000), (0.000000,1.000000,0.000000), (0.000000,0.000000,1.000000) ]
>>> print(Mat(*list(range(1,10))) * Mat(*list(range(10,19))))
Mat[ (84.000000,90.000000,96.000000), (201.000000,216.000000,231.000000), (318.000000,342.000000,366.000000) ]
>>> assert Mat(0.0,1.0,2.0,3,4,5,6,7,8) == Mat(-0,1,2,3,4,5.0,6.0,7.0,8.0)
>>> print(Mat(100,2,3,4,5,6,7,8,9).det())
-297.0
>>> m = Mat(100,2,3,4,5,6,7,8,9)
>>> assert m * ~m == Imat

rotation_axis(R)[source]¶

>>> axis ,ang  = randnorm(),uniform(-pi,pi)
>>> axis2,ang2 = rotation_matrix(axis,ang).rotation_axis()
>>> assert abs( abs(ang) - abs(ang2) ) < EPS
>>> assert axis == axis2 * copysign(1,ang*ang2)

class rif.legacy.xyzMath.Vec(x=0.0, y=None, z=None)[source]¶

a Vector like xyzVector<Real> in rosetta

>>> v = Vec(1,2,3)
>>> print(v, 10*v)
(1.000000,2.000000,3.000000) (10.000000,20.000000,30.000000)

multiplication is a dot prod at the moment >>> v*v 14.0 >>> assert Vec(1,0,-0) == Vec(1,-0,0)

class rif.legacy.xyzMath.Xform(R=None, t=None)[source]¶

Coordinate frame like rosetta Xform, behaves also as a rosetta Stub

>>> x = Xform(R=Imat,t=Uz)
>>> print(x)
Xform( Mat[ (1.000000,0.000000,0.000000), (0.000000,1.000000,0.000000), (0.000000,0.000000,1.000000) ], (0.000000,0.000000,1.000000) )
>>> assert (x*x) == Xform(R=Imat,t=2*Uz)
>>> x = Xform(R=rotation_matrix_degrees(Vec(1,0,0),90.0),t=Vec(0,0,0))
>>> print(x)
Xform( Mat[ (1.000000,0.000000,0.000000), (0.000000,0.000000,-1.000000), (0.000000,1.000000,0.000000) ], (0.000000,0.000000,0.000000) )
>>> assert x*x*x*x == Ixform
>>> x.t = Ux
>>> assert x*x*x*x == Xform(R=Imat,t=4*Ux)
>>> x.t = Uz
>>> print(x)
Xform( Mat[ (1.000000,0.000000,0.000000), (0.000000,0.000000,-1.000000), (0.000000,1.000000,0.000000) ], (0.000000,0.000000,1.000000) )
>>> assert x               == Xform(R=rotation_matrix_degrees(Ux, 90.0),t=Vec(0, 0,1))
>>> assert x*x             == Xform(R=rotation_matrix_degrees(Ux,180.0),t=Vec(0,-1,1))
>>> assert x*x*x           == Xform(R=rotation_matrix_degrees(Ux,270.0),t=Vec(0,-1,0))
>>> assert x*x*x*x         == Xform(R=rotation_matrix_degrees(Ux,  0.0),t=Vec(0, 0,0))
>>> assert x*x*x*x*x       == Xform(R=rotation_matrix_degrees(Ux, 90.0),t=Vec(0, 0,1))
>>> assert x*x*x*x*x*x     == Xform(R=rotation_matrix_degrees(Ux,180.0),t=Vec(0,-1,1))
>>> assert x*x*x*x*x*x*x   == Xform(R=rotation_matrix_degrees(Ux,270.0),t=Vec(0,-1,0))
>>> assert x*x*x*x*x*x*x*x == Xform(R=rotation_matrix_degrees(Ux,  0.0),t=Vec(0, 0,0))
>>> x = Xform(rotation_matrix_degrees(Vec(1,2,3),123),Vec(5,7,9))
>>> assert ~x *  x == Ixform
>>> assert  x * ~x == Ixform

Frames / RTs are interchangable:

>>> fr = Xform(rotation_matrix_degrees(Vec(1,2,3), 65.64),t=Vec(3,2,1))
>>> to = Xform(rotation_matrix_degrees(Vec(7,5,3),105.44),t=Vec(10,9,8))
>>> x = to/fr
>>> assert to/Ixform ==  to
>>> assert Ixform/fr == ~fr
>>> assert (to * ~fr) * fr == to
>>> assert x * fr == to

>>> a1 = randnorm()
>>> b1 = randnorm()
>>> ang = uniform(0,1)*360.0-180.0
>>> a2 = rotation_matrix_degrees(a1.cross(randnorm()),ang) * a1
>>> b2 = rotation_matrix_degrees(b1.cross(randnorm()),ang) * b1
>>> assert abs(angle(a1,a2) - angle(b1,b2)) < EPS
>>> xa = Xform().from_two_vecs(a1,a2)
>>> xb = Xform().from_two_vecs(b1,b2)
>>> assert xa.tolocal(a1) == xb.tolocal(b1)
>>> assert xa.tolocal(a2) == xb.tolocal(b2)
>>> assert ~xa*a1 == ~xb*b1
>>> assert ~xa*a2 == ~xb*b2
>>> assert xb/xa*a1 == b1
>>> assert xb/xa*a2 == b2

add/sub with Vecs:

>>> X = randxform()
>>> u,v = randvec(2)
>>> assert isxform(u+X) and isxform(X+u) and isxform(u-X) and isxform(X-u)
>>> assert X*v+u == (u+X)*v
>>> assert X*(v+u) == (X+u)*v
>>> assert Xform(u)*X*v == (u+X)*v
>>> assert X*Xform(u)*v == (X+u)*v
>>> assert X*v-u == (u-X)*v
>>> assert X*(v-u) == (X-u)*v

mul,div with Mats:

>>> R = randrot()
>>> assert isxform(R*X) and isxform(X*R)
>>> assert R*X*u == (R*X)*u == R*(X*u)
>>> assert X*R*u == (X*R)*u == X*(R*u)
>>> assert Xform(R)*X*u == Xform(R)*(X*u)
>>> assert X*Xform(R)*u == X*(Xform(R,V0)*u)
>>> assert X/X*v == v

mul/div Xforms:

>>> Y = randxform()
>>> assert isxform(X/Y) and isxform(X*Y)
>>> assert X/Y*v == X*~Y*v

these don’t work yet:

>>> axis,ang,cen = randnorm(),uniform(-pi,pi),randvec()      
>>> X = rotation_around(axis,ang,cen)                        
>>> axis2,ang2,cen2 = X.rotation_center()                    
>>> assert abs( abs(ang) - abs(ang2) ) < EPS                 
>>> assert axis == axis2 * copysign(1,ang*ang2)              
>>> print(cen)                                                
>>> print(cen2)                                               

>>> x =                Xform( Mat( Vec( 0.816587, -0.306018, 0.489427 ), Vec( 0.245040, 0.951487, 0.186086 ), Vec( -0.522629, -0.032026, 0.851959 ) ), Vec( 1.689794, 1.535762, -0.964428 ) )
>>> assert repr(x) == "Xform( Mat( Vec( 0.816587, -0.306018, 0.489427 ), Vec( 0.245040, 0.951487, 0.186086 ), Vec( -0.522629, -0.032026, 0.851959 ) ), Vec( 1.689794, 1.535762, -0.964428 ) )"

rif.legacy.xyzMath.align_skew_line_pairs(aa1, ac1, aa2, ac2, ba1, bc1, ba2, bc2)[source]¶

>>> aa1 = Ux
>>> ac1 = V0
>>> aa2 = Uy
>>> ac2 = V0
>>> ba1 = Ux
>>> bc1 = V0
>>> ba2 = Uy
>>> bc2 = V0
>>> print(align_skew_line_pairs(aa1,ac1,aa2,ac2,ba1,bc1,ba2,bc2))
Xform( Mat[ (1.000000,-0.000000,-0.000000), (-0.000000,1.000000,0.000000), (0.000000,-0.000000,1.000000) ], (0.000000,0.000000,0.000000) )

>>> aa1 = Uy
>>> ac1 = V0
>>> aa2 = Uz
>>> ac2 = V0
>>> ba1 = Ux
>>> bc1 = Uz
>>> ba2 = Uy
>>> bc2 = Uz
>>> print(align_skew_line_pairs(aa1,ac1,aa2,ac2,ba1,bc1,ba2,bc2))
Xform( Mat[ (0.000000,1.000000,0.000000), (0.000000,0.000000,1.000000), (1.000000,0.000000,0.000000) ], (0.000000,0.000000,1.000000) )

rif.legacy.xyzMath.alignaroundaxis(axis, u, v)[source]¶

>>> axis = randnorm()
>>> u = randvec()
>>> ang = uniform(-pi,pi)
>>> v = rotation_matrix(axis,ang)*u
>>> uprime = alignaroundaxis(axis,u,v)*u
>>> assert v.angle(uprime) < EPS
>>> v = randvec()
>>> uprime = alignaroundaxis(axis,u,v)*u
>>> assert coplanar(V0,axis,v,uprime)

rif.legacy.xyzMath.alignvector(a, b)[source]¶

>>> u = randvec()
>>> v = randvec()
>>> assert v.angle(alignvector(u,v)*u) < EPS

rif.legacy.xyzMath.alignvectors(a1, a2, b1, b2)[source]¶: same as alignvectors_minangle

rif.legacy.xyzMath.alignvectors_minangle(a1, a2, b1, b2)[source]¶

exact alignment:

>>> for i in range(10):
...         angdeg = uniform(-180,180)
...         a1 = randvec()
...         b1 = randnorm()*a1.length()
...         l2 = gauss(0,1)
...         a2 = rotation_matrix_degrees(a1.cross(randnorm()),angdeg) * a1 * l2
...         b2 = rotation_matrix_degrees(b1.cross(randnorm()),angdeg) * b1 * l2
...         assert abs(angle(a1,a2) - angle(b1,b2)) < EPS
...         Xa2b = alignvectors_minangle(a1,a2,b1,b2)
...         assert Xa2b.t.length() < EPS
...         assert (Xa2b*a1).distance(b1) < EPS
...         assert (Xa2b*a2).distance(b2) < EPS

if angle(a1,a2) != angle(b1,2b), minimize deviation

>>> a1,a2,b1,b2 = randvec(4)
>>> Xa2b = alignvectors_minangle(a1,a2,b1,b2)
>>> assert coplanar(b1,b2,Xa2b*a1,Xa2b*a2)
>>> assert (b1.angle(a1)+b2.angle(a2)) > (b1.angle(Xa2b*a1)+b2.angle(Xa2b*a2))

# >>> tgt1 = -Vec(0.816497,0.000000,0.577350) # >>> tgt2 = Vec(0.000000,0.000000,1.000000) # >>> orig1 = Vec(0.000000,0.000000,1.000000) # >>> orig2 = Vec(-0.723746,0.377967,-0.577350) # >>> print orig1.angle_degrees(orig2) # >>> print tgt1.angle_degrees(tgt2) # >>> x = alignvectors_minangle(orig1,orig2,tgt1,tgt2) # >>> print tgt1,x*orig1 # >>> print tgt2,x*orig2

rif.legacy.xyzMath.coplanar(x1, x2, x3, x4)[source]¶

>>> u,v,w = randvec(3)
>>> a,b,c = (gauss(0,10) for i in range(3))
>>> assert     coplanar(u, v, w, u + a*(u-v) + b*(v-w) + c*(w-u) )
>>> assert not coplanar(u, v, w, u + a*(u-v) + b*(v-w) + c*(w-u) + randvec().cross(u-v) )

rif.legacy.xyzMath.dihedral(p1, p2, p3, p4)[source]¶

>>> dihedral_degrees(Ux,Uy,V0,Uz)
90.0
>>> dihedral_degrees(Ux,V0,Uy,Uz)
-90.0

rif.legacy.xyzMath.expand_xforms(G, N=3, c=Vec( 1.000000, 3.000000, 10.000000 ), maxrad=9000000000.0)[source]¶

>>> G = get_test_generators1()
>>> for x in expand_xforms(G): print(x*Ux)
(-1.000000,0.000000,0.000000)
(1.000000,0.000000,0.000000)
(1.000000,-0.000000,0.000000)
(-1.000000,0.000000,0.000000)
(1.000000,-2.000000,0.000000)
(1.000000,0.000000,0.000000)
(-1.000000,2.000000,0.000000)
(-1.000000,0.000000,0.000000)
(1.000000,-2.000000,0.000000)
(1.000000,-0.000000,-2.000000)

rif.legacy.xyzMath.find_identities(G, n=6, c=Vec( 1.000000, 3.000000, 10.000000 ))[source]¶

>>> G = get_test_generators1()
>>> for I in find_identities(G): print(I.t)
(0.000000,0.000000,0.000000)
(-2.000000,2.000000,2.000000)
(2.000000,-2.000000,2.000000)

rif.legacy.xyzMath.get_cell_bounds_orthogonal_only(G, n=6, c=Vec( 1.000000, 3.000000, 10.000000 ))[source]¶

very slow… need to speed up

>>> G = get_test_generators1()                  
>>> get_cell_bounds_orthogonal_only(G[:2],12)   
(4.0, 4.0, 4.0)

rif.legacy.xyzMath.line_line_closest_points(A1, C1, A2, C2)[source]¶

>>> print(line_line_closest_points(Ux,Ux,Uy,Uy))
(Vec( 0.000000, 0.000000, 0.000000 ), Vec( 0.000000, 0.000000, 0.000000 ))

>>> print(line_line_closest_points(Ux,Uy,Uy,Uz))
(Vec( 0.000000, 1.000000, 0.000000 ), Vec( 0.000000, 1.000000, 1.000000 ))

rif.legacy.xyzMath.line_line_distance(a1, c1, a2, c2)[source]¶

>>> line_line_distance(Ux,V0,Uy,V0)
0.0
>>> round(line_line_distance(Ux,Vec(0,1,2),Ux,Vec(3,2,1)) , 8)
1.41421356
>>> line_line_distance(Ux,10*Uy,Uz,99.0*Ux)
10.0
>>> X = randxform()
>>> round(line_line_distance(X.R*Ux,X*Vec(0,1,2),X.R*Ux,X*Vec(3,2,1)) , 8)
1.41421356

rif.legacy.xyzMath.line_plane_intersection(l, l0, n, p0)[source]¶

>>> l  = Ux
>>> l0 = randvec()
>>> n  = Ux
>>> p0 = V0
>>> assert line_plane_intersection(l,l0,n,p0)[1] == Vec(0,l0.y,l0.z)
>>> n = randnorm()
>>> p0 = randvec().cross(n)
>>> l = randvec()
>>> l0 = p0+l*gauss(0,10)
>>> assert line_plane_intersection(l,l0,n,p0)[1] == p0

rif.legacy.xyzMath.point_line_distance(p, a, c)[source]¶

>>> point_line_distance(V0,Uy,V0)
0.0
>>> round(point_line_distance(V0,Uy,Ux+Uz),8)
1.41421356
>>> round(point_line_distance(Ux,Ux,Vec(3,2,1)) , 8)
2.23606798

rif.legacy.xyzMath.proj(u, v)[source]¶

>>> u = Vec(1,0,0); v = Vec(1,1,1)
>>> proj(u,v)
Vec( 1.000000, 0.000000, 0.000000 )

rif.legacy.xyzMath.projperp(u, v)[source]¶

>>> u = Vec(1,0,0); v = Vec(1,1,1)
>>> projperp(u,v)
Vec( 0.000000, 1.000000, 1.000000 )

rif.legacy.xyzMath.randnorm(n=None)[source]¶

>>> assert abs(randnorm().length()-1.0) < 0.0000001

rif.legacy.xyzMath.ray_sphere_intersection(lin, l0in, cen, r)[source]¶

>>> v = randnorm()
>>> v = Ux
>>> assert v.distance( ray_sphere_intersection(v,V0,V0,1) ) < 0.00001
>>> assert not ray_sphere_intersection(v,V0,-2*v,1.0)

rif.legacy.xyzMath.rmsd(l, m)[source]¶

>>> l,m = randvec(6),randvec(6)
>>> rmsd(l,l)
0.0

rif.legacy.xyzMath.rotation_around(axs, ang, cen=None)[source]¶

>>> x = rotation_around(Ux,1,Uy)
>>> x * Uy
Vec( 0.000000, 1.000000, 0.000000 )

rif.legacy.xyzMath.rotation_around_dof_to_set_vec_vec_angle(dofaxis, tgt0, v1, v2)[source]¶

>>> from random import random
>>> for i in range(10):
...         dof = randnorm()
...         tgt = random()*180
...         v1 = randnorm()
...         v2 = randnorm()
...         ANGs = rotation_around_dof_to_set_vec_vec_angle(dof,tgt,v1,v2)
...         for a in ANGs:
...                 R = rotation_around_degrees(dof,a,Vec(0,0,0))
...                 act = line_line_angle_degrees(v1,R*v2)
...                 if 0.000001 < min(abs(act-tgt),abs(act-180+tgt)):
...                         print(a,tgt,act)
>>> print("if no other output, correct")
if no other output, correct

rif.legacy.xyzMath.rotation_matrix_degrees(axis, angle)[source]¶

get a rotation matrix

>>> rx180 = rotation_matrix_degrees(Vec(1,0,0),180.0)
>>> rx90  = rotation_matrix_degrees(Vec(1,0,0),90.0)
>>> print(rx90*rx90 == rx180)
True
>>> r = rotation_matrix_degrees(Vec(1,0,0),45.0)
>>> print(r)
Mat[ (1.000000,0.000000,0.000000), (0.000000,0.707107,-0.707107), (0.000000,0.707107,0.707107) ]
>>> assert r*r == rx90
>>> assert r*r*r*r == rx180
>>> assert r*r*r*r*r*r*r*r == Imat
>>> assert ~r == r.transposed()

>>> ang = uniform(0,1)*360.0-180.0
>>> v = randvec()
>>> axs = randnorm()
>>> while(abs(v.dot(axs))>0.9): axs = randnorm()
>>> u = rotation_matrix_degrees(projperp(v,axs),ang)*v
>>> assert abs(u.angle_degrees(v)-abs(ang)) < 10*SQRTEPS

rif.legacy.xyzMath.skew_lines_center(A1, C1, A2, C2)[source]¶

>>> skew_lines_center(Ux,V0,Uy,V0)
Vec( 0.000000, 0.000000, 0.000000 )

>>> skew_lines_center(Ux,Uy,Uy,Ux)
Vec( 1.000000, 1.000000, 0.000000 )

>>> skew_lines_center(10*Ux,10*Uy,10*Uy,10*Uz)
Vec( 0.000000, 10.000000, 5.000000 )

# >>> skew_lines_center(Ux,Uy,Uz,Ux)

rif.legacy.xyzMath.slide_to_make_lines_intersect(dof, l, l0, m, m0)[source]¶

>>> v = randvec()
>>> assert abs(slide_to_make_lines_intersect(Ux,Uy,v,Uz,V0) + v.x ) < EPS
>>> dof,l,l0,m,m0 = randvec(5)
>>> d = slide_to_make_lines_intersect(dof,l,l0,m,m0)
>>> l0 = l0 + d*dof
>>> assert abs(line_line_distance(l,l0,m,m0)) < EPS

rif.legacy.xyzMath.symmetrize_xform(anchor, x, nf=None)[source]¶: # >>> x = rotation_around_degrees(Uz,180,V0) # >>> assert symmetrize_xform(Ux,x,2) == x #>>> x = rotation_around_degrees(Vec(0,0,1),121,Vec(1,1,1)) #>>> x,c = symmetrize_xform(Ux,x,3); print x.pretty(); print c # >>> print x.pretty()

rif.legacy.xyzMath_test module¶

Module contents¶

docstring for rif.legacy

rif.numeric package¶

Submodules¶

rif.numeric.jagged module¶

rif.numeric.jagged_test module¶

rif.numeric.lattice module¶

rif.numeric.lattice_test module¶

Module contents¶

docstring for rif.numeric

rif.rcl package¶

Submodules¶

rif.rcl.conversions module¶

rif.rcl.conversions.atoms(pose, sele='HEAVY', **kwargs)[source]¶: extract rif style atoms from rosetta pose

rif.rcl.conversions.bbstubs(pose, which_resi=None)[source]¶: extract rif style stubs from rosetta pose

rif.rcl.conversions.rays(pose, sele, shifts=None, protein_only=True, **kwargs)[source]¶: extract rif style rays from rosetta pose

rif.rcl.conversions_test module¶

rif.rcl.pyrosetta_util module¶

nRosetta Compatibility Layer

if pyrosetta is available, provides pyrosetta. check value of HAVE_PYROSETTA.

exception rif.rcl.pyrosetta_util.AlreadyInitError[source]¶: Rosetta already initialized but should not be

exception rif.rcl.pyrosetta_util.Error[source]¶: Rosetta Compatibility Exception

exception rif.rcl.pyrosetta_util.ReInitError[source]¶: Tried to reinitialize rosetta with different options

exception rif.rcl.pyrosetta_util.RosettaDBError[source]¶: Rosetta Database Error

rif.rcl.pyrosetta_util.add_files_to_rosettaDB(fname, contents)[source]¶: returns true if update was done

rif.rcl.pyrosetta_util.add_residue_types_to_rosettaDB(newfiles, comment='## RIF added restypes', typeset='fa_standard')[source]¶: return true if update was done

rif.rcl.pyrosetta_util.concatenate_pose(p, to_append, lb=1, ub=-1, out=None)[source]¶: put poses together aligning on end of p

rif.rcl.pyrosetta_util_test module¶

Module contents¶

rif.rotamer package¶

Submodules¶

rif.rotamer.richardson module¶

rif.rotamer.richardson.get_rotamer_index(rotspace)[source]¶: extract AA structural info via pyrosetta

rif.rotamer.richardson.get_rotamer_space[source]¶: B is disulfide, x4 -> x2other n num observations fabo total fraction/alpha/beta/other x#a chi # “act” x# chi # “com” x#r chi # range text x#w chi # peak width lb#/ub# chi # range

rif.rotamer.richardson.sample_rotamer_space(rotspace, resl=[10, 10, 10, 10])[source]¶: rotamer samples

rif.rotamer.richardson_test module¶

Module contents¶

docstring for rif.rotamers

rif.sampling package¶

Submodules¶

rif.sampling.orientations module¶

rif.sampling.orientations_test module¶

Module contents¶

docstring for rif.sampling

rif.sele package¶

Submodules¶

rif.sele.sele module¶

rif.sele.sele_test module¶

Module contents¶

rif.visual package¶

Subpackages¶

rif.visual.datavis package¶

Submodules¶

rif.visual.datavis.gridplot module¶

Module contents¶

Submodules¶

rif.visual.vispymol module¶

Module contents¶

Submodules¶

rif.dtypes module¶

overrides numpy operators to work with rif dpypes

override numpy operators using numpy.set_numeric_ops. keep track of overridden operators and call them if a custom operator is not found

class rif.dtypes.RifOperators[source]¶: contect manager for locally enabling rif ops

class rif.dtypes.RifOperatorsDisabled[source]¶: context manager to locally disable RifOperators

rif.dtypes.global_rif_operators_disable(quiet=False)[source]¶: disable rif operators via numpy.set_numeric_opts

rif.dtypes.global_rif_operators_enable(quiet=False)[source]¶: enable rif operators via numpy.set_numeric_opts

Module contents¶

docstring for rif

Legacy¶

This section documents the original RIF docking method and the legacy C++ implementation. Docs are here because (1) this is a convenient place for me to put them, and (2) the current library has grown out of the C++ implementation.

Overview of Prototype one-sided design implementation¶

Rotamer Interaction Field (RIF) Docking for de novo design of binders to small molecules and fixed protein targets works in two phases: RIF generation and scaffold docking.

RIF Generation¶

First, a RIF tailored to the small molecule target is generated. The RIF contains billions of rotamers that interact productively with the target in both polar and apolar fashion. Polar interactions with the ligand are generated based on hydrogen bond geometry and are tagged with the polar group(s) on the small molecule with which they hydrogen bond (Figure XXX, groups with numeric labels). All possible polar interactions are included in the RIF. Apolar interacting rotamers are generated with a docking process that generates all possible packing interactions above an configurable amino acid specific rosetta energy threshold (Figure XXX label APO). The polar and apolar RIF rotamers are stored ~0.5Å gridded representation (REF 1) of the six dimensional rigid body space, based on the backbone position of the rotamers generated. To facilitate search, described next, lower resolution RIF data structures are produced from the primary, gridded at 1.0Å, 2.0Å, 4.0Å, 8.0Å and 16.0Å resolutions.`

Implemented in binary rifgen

RIF Docking¶

Second, after RIF generation, a set of protein backbone scaffolds selected by the user is docked into the RIF using a hierarchical branch and bound style search (REF 2). Starting with the coarsest 16Å resolution RIF, we perform an enumerative search of possible 6D positions of each scaffold. Foreach scaffold placement, the designable scaffold backbone positions, as specified by the user, are checked against the RIF to determine a set of productively interacting rotamers which can be placed on the scaffold. If the rotamers found satisfy the all polar groups on the small molecule target (or a user-specified acceptable fraction), the scaffold position is assigned a score based on the quality of the polar and apolar interactions that can be formed, else it is rejected. All acceptable scaffold positions up to a configurable limit (typically 10 million) are sorted by score and promoted to the next search stage. Each promoted scaffold is split into 64 child positions by moving the scaffold in the 6D rigid body space, providing a finer sampling appropriate for the next RIF resolution. Search is done at 16Å, 8Å, 4Å, 2Å, 1Å, and 0.5Å resolutions. After the final 0.5Å search stage, a monte-carlo based combinatorial rotamer packing step is performed on a configurable fraction of the best scaffold placements (typically 1 million) to find internally consistent rotamer placements. The docking process is highly optimized and takes roughly one cpu hour per scaffold to generate thousands of possible binding modes.

Implemented in binary rif_dock_test

Slide Decks¶

IPD Presentation 2016/08¶

more “polish” than 2016/04, benchmark data here are accurate, to best of my knowledge

IPD Presentation 2016/04¶

more technical info vs 2016/08, protein docking DDG benchmark numbers inaccurate here!

Building¶

The current RIF docking protocol is not implemented in python and is not part of this repository.

Building Rosetta¶

There is only one dependency for the legacy rif implementation, but it’s a whopper: Rosetta. You will need a relatively current clone of rosetta main branch master from http://github.com/rosettacommons/main. Building with a rosetta release version should be possible, but has not been tested.

You must use rosetta’s cmake build system as follows:

cd <path_to_rosetta>/source/cmake
python make_project.py all
cd build_cxx11_omp
cmake .
make -j <num_processors>

The rosetta build will take a few cpu hours.

Obtain legacy scheme source¶

First, obtain the legacy rif source you want to build. There are a few choices. The ‘standard’ version most people in the bakerlab are using is from willsheffler@gmail.com. It has been largely unchanged for a year and should be relatively stable. There are also versions which have been enhanced with some additional features by other members of the baker lab, most notably one available from Brian Coventry https://github.com/bcov77/scheme/tree/brian.

Building RIF¶

The following should produce a build:

git clone git@github.com:willsheffler/scheme
cd scheme
git checkout sheffler/devel

The above will get you the ‘stable’ version from willsheffler@gmail.com. If you with to use another, substitute the repository and branch as appropriate. For example, for Brian Coventry’s version:

git clone git@github.com:bcov77/scheme
cd scheme
git checkout brian

Do the following to produce build config with cmake and build (ninja can be used instead of make). Note that the CMAKE_ROSETTA_PATH environment variable must be specified as this code links against rosetta.

export CMAKE_ROSETTA_PATH='<path_to_rosetta>'
mkdir build
cd build
cmake .. -DCMAKE_BUILD_TYPE=Release
make -j <num_processors> rifgen rif_dock_test

Usage¶

The build will product two binaries: rifgen for Rif generation, and rif_dock_test for RIF Docking. Both binaries are multi-threaded via OpenMP. You may control the number of threads used with the OMP_NUM_THREADS environment variable. The default is to use all available hardware threads (which is probably best). Running rifgen and rif_dock_test will in most require at least 64gb of ram per process, hence, it is best to use more threads per process instead of multiple processes. For some cases, particularly if you are not generating a full de novo RIF but are instead manually inputting hotspots, you may be able to get away with less ram.

You really really want to contact somebody who is using rif for a similar problem as you are and use their flags files as a starting point. Broadly speaking, people in the Baker lab are using rif for one-sided protein interface design (particularly denovo) and for small molecule binder design. If you contact willsheffler@gmail.com I can put you in touch with somebody.

rifgen¶

like most rosetta-based programs, usage is: rifgen @flags_file

rifgen makes use of a fair amount of NON-TARGET-SPECIFIC cached data. In most cases, this data will be generated if it is not available. This means the first time you run rifgen, it could take a long time. This data is stored and accessed in the location passed to ‘-rifgen:data_cache_dir’. As this data is not target specific, it can be shared between runs, projects, and users. Most people in the baker lab point to one such directory, for example.

please please contact somebody doing similar work and use their flags as a starting point. This documentation is likely incomplete and some may be out of date.

flags you will likely need to modify¶

-rifgen:data_cache_dir <location_for_cached_data>

-database <path to rosetta database, MUST be the same database as the rosetta you built>

-rifgen::rif_type RotScore  # use RotScoreSat for small molecules

-rifgen:target     <pdb file>
-rifgen:target_res <residue_numbers.txt> OPTIONAL
-rifgen:outdir     <output_dir_location>
-rifgen:outfile    <output_file_name.rif.gz>

flags you may want need to tweak¶

# apolar residues to generate in de novo RIF
-rifgen:apores VAL ILE LEU MET PHE TRP

# polar residues to generate in de novo RIF
-rifgen:donres SER THR TYR     GLN     ASN HIS HIS_D TRP
-rifgen:accres SER THR TYR GLU GLN ASP ASN HIS HIS_D

# multiplier for default hydrophobic residue score cut
# this will largely determine how long rif generation takes and
# how large the resulting rif is. lowering this number will cause
# more possible hydrophobic interactions to be found, but will make
# the rif bigger and make the search take longer. if you don't care
# much about hydrophobic interactions, or only want a few good ones
# raise this number to maybe 1.2 (fine to play with it) OR if you
# aren't getting good hydrophobic packing, lower it to maybe 0.8 or 0.7.
-rifgen:score_cut_adjust 1.0   # scale factor for hydrophobic score cutoffs

-rifgen::rosetta_field_resl 0.25          # grid spacing for vdw/lksol scoring
-rifgen::search_resolutions 3.0 1.5 0.75  # search resl for denovo hydrophobic interactions

-hash_cart_resl    0.7  # memory use will scale as (1/resl)**6 so be careful!
-hash_angle_resl  14.0  # should stay roughly 20x hash_cart_resl

-rif_accum_scratch_size_M 24000  # megabytes of memory to use as a scratch space

-rifgen:score_threshold 0          # score cut for rotamer inclusion in rif
-rifgen:hbond_weight 1.0           # max score per-hbond
-rifgen:upweight_multi_hbond 1.0   # extra score factor for bidentate hbonds

# for sanity-checking only, you can dump a (small!) fraction RIF rotamers:
-rifgen:rif_hbond_dump_fraction  0.000001
-rifgen:rif_apo_dump_fraction    0.000001

flags you must have but probably shouldn’t change¶

# same name as rif_dock_test flag, but not same!
# this roughly controls how many apolar rotamers are 'docked' for denovo apolar rif gen
# default 1 billion
-rifgen:beam_size_M 1000.0

-rifgen:extra_rotamers false
-rifgen:extra_rif_rotamers true

# rosetta stuff
-renumber_pdb
-add_orbitals

-hbond_cart_sample_hack_range 0.33
-hbond_cart_sample_hack_resl  0.33
-rifgen:tip_tol_deg        60.0 # for now, do either 60 or 36
-rifgen:rot_samp_resl       6.0

-rifgen:hash_preallocate_mult 0.125
-rifgen:max_rf_bounding_ratio 4.0

# geometry of bounding grids
-rifgen:hash_cart_resls   16.0   8.0   4.0   2.0   1.0
-rifgen:hash_cart_bounds   512   512   512   512   512
-rifgen:lever_bounds      16.0   8.0   4.0   2.0   1.0
-rifgen:hash_ang_resls     38.8  24.4  17.2  13.6  11.8 # yes worky worky
-rifgen:lever_radii        23.6 18.785501 13.324600  8.425850  4.855575

output¶

All products of rifgen will be in the directory specified by the ‘-output_dir’ flag. In addition, rifgen will dump output like the following at the end of the run. You can copy and paste this output into your rif_dock_test flags file to avoid some tedious and error-prone editing.

########################################### what you need for docking ###########################################
-rif_dock:target_pdb            rifgen_fz7_v4/rif_64_fz7_v4_sca0.7_noKR.rif_target.pdb.gz
-rif_dock:target_res            test_input/fz7/Fz7_model.res
-rif_dock:target_rf_resl        0.125
-rif_dock:target_rf_cache       rifgen_fz7_v4/__RF_Fz7_model.pdb_CEN_trhash12511471_resl0.125_osamp2_replonlybdry
-rif_dock:target_bounding_xmaps rifgen_fz7_v4/rif_64_fz7_v4_sca0.7_noKR.rif_BOUNDING_RIF_16.xmap.gz
-rif_dock:target_bounding_xmaps rifgen_fz7_v4/rif_64_fz7_v4_sca0.7_noKR.rif_BOUNDING_RIF_08.xmap.gz
-rif_dock:target_bounding_xmaps rifgen_fz7_v4/rif_64_fz7_v4_sca0.7_noKR.rif_BOUNDING_RIF_04.xmap.gz
-rif_dock:target_bounding_xmaps rifgen_fz7_v4/rif_64_fz7_v4_sca0.7_noKR.rif_BOUNDING_RIF_02.xmap.gz
-rif_dock:target_bounding_xmaps rifgen_fz7_v4/rif_64_fz7_v4_sca0.7_noKR.rif_BOUNDING_RIF_01.xmap.gz
-rif_dock:target_rif            rifgen_fz7_v4/rif_64_fz7_v4_sca0.7_noKR.rif.gz
-rif_dock:extra_rotamers        0
-rif_dock:extra_rif_rotamers    1
#################################################################################################################

rif_dock_test¶

like most rosetta-based programs, usage is: rif_dock_test @flags_file

like rifgen, rif_dock_test makes use of some NON-TARGET-SPECIFIC cached data. In most cases, this data will be generated if it is not available. This means the first time you run, it could take a long time. This data is stored and accessed in the location passed to ‘-rif_dock:data_cache_dir’. As this data is not target specific, it can be shared between runs, projects, and users. Most people in the baker lab point to one such directory, for example.

please please contact somebody doing similar work and use their flags as a starting point. This documentation is likely incomplete and some may be out of date.

data cache paths¶

# standard rosetta database flag, must match rosetta you compiled with
-database <path_to_rosetta_database>

# precomputed data for computing 2-body energy tables at hyperspeed
-rif_dock:rotrf_cache_dir <location>

# location to store per-scaffold non-target specific data
# mostly one and two body energy tables
-rif_dock:data_cache_dir  <scaffold data cache>
-rif_dock:cache_scaffold_data true  # set to false if you don't want to use/generate

flags you will likely need/want to modify¶

-rif_dock:scaffolds # list of scaffold pdb files

 # optional list of 'res' files for each scaffold
 # must be either one file which applies to all input scaffolds
 # or a list exactly matching the number of input scaffolds
-rif_dock:scaffold_res

# flags that should be copied from rifgen log output, just copy them
-rif_dock:target_pdb            <centered target pdb file>
-in:file:extra_res_fa           <ligand params, only if necessary>
-rif_dock:target_rf_resl        0.125
-rif_dock:target_rf_cache       <vdw/lksol target score grid files prefix>
-rif_dock:target_bounding_xmaps <bounding maps>
-rif_dock:target_rif            <primary rif file>
-rif_dock:extra_rotamers        0
-rif_dock:extra_rif_rotamers    1

# this is where the output will go, and how much
-rif_dock:outdir rifdock_v4_out/ads_test
-rif_dock:dokfile all.dok
-rif_dock:n_pdb_out 20 # max number of output pdbs

# optional flag to add extra output file tag
#-rif_dock:target_tag conf01


# set to true to align all output to scaffolds instead of target
# mostly useful for small molecule binder design
-rif_dock:align_output_to_scaffold false

# include some pikaa lines in output pdbs for rif residues(??)
-rif_dock:pdb_info_pikaa false # this is default I think

# *IFF* using a rif with satisfaction constraints (e.g. RotScoreSat)
# set this to the number of hbonds to the target which are required
# no way yet to explicity say which hbonds are required, but this gives
# some control. searches will be much faster if this number is higher
# of course, don't make it higher than the number of hbonds that can
# actually be made
-require_satisfaction 5

flags you may want/need to tweak¶

-beta
-score:weights beta_soft
-add_orbitals false # this breaks -beta scoring

# these flags control the overall time the search will take. a few alternte options are included
# setting the require_satisfaction flag above to a high vaule will make search faster across the
# board, so experiment with that also

# reasonable defaults:
-beam_size_M 5
-hsearch_scale_factor 1.2

# # very fast search, probably with low quality results
# -beam_size_M 1
# -hsearch_scale_factor 1.6

# # slower more thorough search
# -beam_size_M 30
# -hsearch_scale_factor 1.0


# make this number higher to have less redundant results or lower to have more similar results
# this is NOT a proper rmsd (yet), unfortunately, so if you want to tweak it you'll have to experiment
-rif_dock:redundancy_filter_mag 1.5


# hbond weight relative to vdw/lksol
-rif_dock:hbond_weight 2.0
# bonus term for bidentate hbonds
# value of 1.0 could up to double hbscore if bidentate, triple if tridentate...
# best in conjunction with low-ish starting hbweight
-rif_dock:upweight_multi_hbond 1.0

rosetta re-scoring rifdocks¶

Rosetta scoring is really slow compared to rif design, so usually only a small fraction of the rif designs are re-scored with the rosetta energy function.

-rif_dock:rosetta_score_fraction 0.015 # set to 0 if you don't want to re-score with rosetta
#-rif_dock:rosetta_score_then_min_below_thresh -20.0 # this is in "rif docking score units"
#-rif_dock:rosetta_score_at_least 3000
#-rif_dock:rosetta_score_at_most  30000
-rif_dock:rosetta_min_fraction 0.10 # fraction of re-scored designs to minimize
-rif_dock:rosetta_min_scaffoldbb false
-rif_dock:rosetta_min_targetbb   false
-rif_dock:rosetta_hard_min false
-rif_dock:global_score_cut -10.0

flags you must have but probably shouldn’t change¶

-rif_dock:replace_orig_scaffold_res false

-rif_dock::pack_iter_mult 1.0
-rif_dock:hack_pack_frac 0.025
-hack_pack true

-rif_dock::rf_resl 0.5
-rif_dock::rf_oversample 2

-rif_dock:rotrf_resl   0.25
-rif_dock:rotrf_spread 0.0
-rif_dock:rotrf_scale_atr 1.0

-rif_dock:use_scaffold_bounding_grids 0

-rif_dock:upweight_iface 1.3

-rif_dock:scaffold_to_ala true
-rif_dock:scaffold_to_ala_selonly false

# these are for encouraging
-add_native_scaffold_rots_when_packing 0 # 1
-bonus_to_native_scaffold_res          0 # -0.5

output¶

rif_dock_test will output up to -rif_dock:n_pdb_out pdb files per input scaffold with partial docks with interface-only design to the directory specified by -rif_dock:outdir. Additionally, a file with scores for each dock/design to the file specified by -rif_dock:dokfile, possibly with a suffix to prevent overwrite of a previous file.

test link: rif.legacy.xyzMath

test link: rif.legacy.xyzMath.Xform