In [2]:
# In this tutorial we'll go over how to simulate
# x-ray scattering on a polar detector. This is
# a fast way to test out what correlations should
# look like.

from thor import xray
import mdtraj

In [3]:
# load in a PDB file
t = mdtraj.load('3LYZ.pdb') # your fav pdb

In [4]:
n_shots = 5                         # total number of shots to do
n_molecules = 5                     # the number of molecules to include per shot
q_values = np.arange(1.2,1.6,0.02)  # the |q| values of the rings to sim
n_phi = 360                         # number of pts around the rings

rings = xray.Rings.simulate(t, n_molecules, q_values, n_phi, n_shots)


18:17:46 - Finished polar shot 1/5 on device 0
18:17:47 - Finished polar shot 2/5 on device 0
18:17:49 - Finished polar shot 3/5 on device 0
18:17:50 - Finished polar shot 4/5 on device 0
18:17:52 - Finished polar shot 5/5 on device 0

In [5]:
# check the WAXS
figure()

waxs = rings.intensity_profile() # first col is |q|, second is intensity

plot(waxs[:,0], waxs[:,1], lw=2)
xlabel(r'$q / \AA^{-1}$')
ylabel('Intensity')
show()



In [7]:
# looks like we have a peak at |q| = 1.46, so let's correlate that
# I'll make a plot of each ring's correlation function independently

corr = rings.correlate_intra(1.46, 1.46, mean_only=False) # args are q_1, q_2

figure()

for i in range(rings.num_shots):
    plot(rings.phi_values, corr[i,:], lw=2)

xlabel(r'$\phi$ / Radians')
ylabel(r'$C(\phi)$')
show()


18:18:19 - Correlating rings at 1.460000 / 1.460000
18:18:19 - Correlating shots...
18:18:19 - Correlating shot 1/5
18:18:19 - ... complete

In [8]:
# let's check out the average of each of those shots
# also include a visualization of the std

figure()
plot(rings.phi_values, corr.mean(0), lw=2)
plot(rings.phi_values, corr.mean(0) + corr.std(0), lw=0.5, color='k')
plot(rings.phi_values, corr.mean(0) - corr.std(0), lw=0.5, color='k')
xlabel(r'$\phi$ / Radians')
ylabel(r'$C(\phi)$')
show()



In [ ]: