3D Molecules

Understanding molecular structure is a great challenge. In X-ray imaging, for instance, the experimentalist illuminates an object of interest, e.g. a molecule, and then collects the intensity of the diffracted rays, please see figure below for an illustrative setup.


The two figures below show the schematic representation and the corresponding electron density maps for the Caffeine and Nicotine molecules: the density map rho(x_1, x_2, x_3) is the 3D object we seek to infer. Here, we demonstrate the performance of the Wirtinger flow algorithm for recovering projections of 3D molecule density maps from simulated data.


Figure: Schematic representation and electron density map of the Caffeine molecule.


Figure: Schematic representation and electron density map of the Nicotine molecule.

Consider an experimental apparatus as in the first figure above. If we imagine that light propagates in the direction of the x_3-axis, an approximate model for the collected data reads

   I(f_1,f_2) = Big| intleft(int rho(x_1,x_2,x_3)       d x_3right)e^{-2ipi(f_1x_1+f_2x_2)}d x_1 d x_2Big|^2.

In other words, we collect the intensity of the diffraction pattern of the projection int rho(x_1,x_2,x_3) d x_3. The 2D image we wish to recover is thus the line integral of the density map along a given direction. As an example, the Caffeine molecule along with its projection on the x_1x_2-plane (line integral in the x_3 direction) is shown in the Figure below.


Figure: Electron density rho(x_1,x_2,x_3) of the Caffeine molecule along with its projection onto the x_1x_2-plane.

Now, if we let R be the Fourier transform of the density rho, one can re-express the identity above as


Therefore, by imputing the missing phase using phase retrieval algorithms, one can recover a slice of the 3D Fourier transform of the electron density map, i.e. R(f_1,f_2,0). Viewing the object from different angles or directions gives us different slices. One can then recover the 3D Fourier transform of the electron density map from all these slices and, in turn, the 3D electron density map.

We now generate 51 observation planes by rotating the x_1x_2-plane around the x_1-axis by equally spaced angles in the interval [0,2pi]. Each of these planes is associated with a 2D projection of size 1024times 1024, giving us 20 coded diffraction octanary patterns (we use the same patterns for all 51 projections). We run the Wirtinger flow algorithm. The figure below reports the average relative error over the 51 projections and the total computational time required for reconstructing all 51 images.

Caffeine molecule: Mean rel. error is 9.6times 10^{-6}; Total time is 5.4 hours.