NMR

pytraj.analysis.nmr.jcoupling(traj=None, mask='', top=None, kfile=None, dtype='dataset', frame_indices=None)

compute j-coupling

Parameters:

traj : any things that make frame_iter_master returning Frame object

command : str, default “”

cpptraj’s command/mask

kfile : str, default None, optional

Dir for Karplus file. If “None”, use $AMBERHOME dir

dtype : str, {‘dataset’, ...}, default ‘dataset’

Examples

>>> import pytraj as pt
>>> traj = pt.datafiles.load_tz2()
>>> data = pt.calc_jcoupling(traj, ':1-12', kfile='data/Karplus.txt')
pytraj.analysis.nmr.ired_vector_and_matrix(traj=None, mask='', frame_indices=None, order=2, dtype='tuple', top=None)

perform vector calculation and then calculate ired matrix

Parameters:

traj : Trajectory-like or iterable that produces pytraj.Frame

mask : str or a list of strings

frame_indices : array-like, optional, default None

only perform calculation for given frame indices

order : default 2

dtype : output’s dtype, {‘dataset’, ‘tuple’} default ‘dataset’

top : Topology, optional, default None

Returns:

out : if dtype is ‘dataset’, return pytraj.DatasetList with shape=(n_vectors+1,)

last index is a matrix, otherwise n_vectors. If dtype is ‘tuple’, return a a tuple of vectors and matrix

Examples

>>> import pytraj as pt
>>> traj = pt.datafiles.load_tz2_ortho()
>>> h = pt.select('@H', traj.top)
>>> n = h - 1
>>> nh = list(zip(n, h))
>>> vecs, mat = pt.ired_vector_and_matrix(traj, mask=nh)
>>> dslist = pt.ired_vector_and_matrix(traj, mask=nh, dtype='dataset')
pytraj.analysis.nmr.nh_order_parameters(traj, vector_pairs, order=2, tstep=1.0, tcorr=10000.0, n_cores=1, **kwargs)

compute NH order parameters

Parameters:

traj : Trajectory-like

vector_pairs : 2D array-like, shape (n_pairs, 2)

order : default 2

tstep : default 1.

tcorr : default 10000.

kwargs : additional keyword argument

Returns:

S2 : 1D array, order parameters

Examples

>>> import pytraj as pt
>>> traj = pt.datafiles.load_tz2_ortho()
>>> h_indices = pt.select('@H', traj.top)
>>> n_indices = h_indices - 1
>>> nh_pairs = list(zip(n_indices, h_indices))
>>> data = pt.nh_order_parameters(traj, nh_pairs)