In [ ]:

    

import h5py
    import condor
    import copy
    import numpy as np
    import condor.utils.linalg as linalg
    import matplotlib.pyplot as plt
    import matplotlib.colors as colors
    import matplotlib.cm as cm


    

In [ ]:

    

fullp_res_edge   = lambda det,w:  1./linalg.length(det.get_q_max(w, pos='edge')/(2*np.pi))
fullp_res_corner = lambda det,w:  1./linalg.length(det.get_q_max(w, pos='corner')/(2*np.pi))


    

In [ ]:

    

def plot_diffraction(data, mask=None):
    if mask is None:
        mask = np.ones_like(data).astype(np.bool)
    image = np.ma.array(np.copy(data).astype(np.float), mask=~mask)
    image.data[image.data<0.5] = 1e-5
    palette = copy.copy(cm.magma)
    palette.set_bad('w', 1.)
    palette.set_under('0.9',1.)
    plt.figure()
    plt.imshow(image, norm=colors.LogNorm(vmin=1), interpolation='none', cmap=palette)
    plt.axis('off')
    cb = plt.colorbar(pad=0)
    cb.ax.set_ylabel('Intensity [photons / px]')
    cb.outline.set_visible(False)
    plt.tight_layout()
    plt.show()


    

In [ ]:

    

# Detector (pnCCD)
pixelsize = 8*75e-6
nx,ny = (1024//8,1024//8)
detector_distance = 300e-3


    

In [ ]:

    

# Source
N = 10
photon_energy = np.random.normal(loc=3500, scale=10, size=(N,)) # [eV]
ph = condor.utils.photon.Photon(energy_eV=photon_energy)
wavelength = ph.get_wavelength()
fluence = 1e15 #[ph/um2]
focus_diameter = 0.2e-6
pulse_energy = fluence * ph.get_energy() * (np.pi*((1e6*focus_diameter/2.)**2)) # [J]


    

In [ ]:

    

# Sample
pdb_id = '1FFK'
sample_size = 18e-9


    

In [ ]:

    

angle_degrees = 72.5
angle = angle_degrees/360.*2*np.pi
rotation_axis = np.array([1.,1.,0.])/np.sqrt(2.)
quaternion = condor.utils.rotation.quat(angle,rotation_axis[0],rotation_axis[1], rotation_axis[2])
rotation_values = np.array([quaternion])
rotation_formalism = "quaternion"
rotation_mode = "extrinsic"


    

In [ ]:

    

incoherent_average = 0
det = condor.Detector(distance=detector_distance, pixel_size=pixelsize, nx=nx, ny=ny, noise="poisson")
par = condor.ParticleAtoms(pdb_id=pdb_id,
                           rotation_formalism=rotation_formalism, 
                           rotation_values=rotation_values,
                           rotation_mode=rotation_mode)
for i in range(N):
    src = condor.Source(wavelength=wavelength[i], pulse_energy=pulse_energy[i], focus_diameter=focus_diameter)
    E = condor.Experiment(source=src, particles={"particle_atoms":par}, detector=det)
    o = E.propagate()
    incoherent_average += o["entry_1"]["data_1"]["data"][:]
incoherent_average /= N
o["entry_1"]["data_1"]["data"] = incoherent_average


    

In [ ]:

    

print "Photon energy: %d eV" %photon_energy.mean()
print "Fluence: %g ph/um2" %fluence
print "Pulse energy: %.4f mJ" %(1e3*pulse_energy.mean())
print "Sample size: %d nm" %(1e9*sample_size)
print "Detector distance: %d mm" %(1e3*detector_distance)
print "Full period resolution (corner): %.2f nm" %(1e9*fullp_res_corner(det,wavelength.mean()))
print "Full period resolution (edge): %.2f nm" %(1e9*fullp_res_edge(det,wavelength.mean()))
print "Nr. of resolution elements: %.2f" %(sample_size / fullp_res_corner(det,wavelength.mean())*2)
plot_diffraction(o['entry_1']['data_1']['data'])
print "\n"


    

In [ ]:

    

W = condor.utils.cxiwriter.CXIWriter("../data/single_protein_polychromatic.h5")
W.write(o)
W.close()