 codex/electron-density-phase-example.py

"""Plot an original one-dimensional phase example; write PNG to the cwd.

Uses the repository Python 3.14, NumPy and Matplotlib environment.
Profiles are stipulated functions, not measured crystallographic density.
"""
from pathlib import Path
import os
import matplotlib
matplotlib.use('Agg')
import matplotlib.pyplot as plt
import numpy as np


def main():
    x = np.linspace(0, 1, 501)
    y0 = 3 + 2 * np.cos(2 * np.pi * x)
    ypi = 3 - 2 * np.cos(2 * np.pi * x)
    samples = np.arange(64) / 64
    f0 = np.fft.rfft(3 + 2 * np.cos(2 * np.pi * samples)) / 64
    fpi = np.fft.rfft(3 - 2 * np.cos(2 * np.pi * samples)) / 64
    assert np.allclose(np.abs(f0), np.abs(fpi), atol=1e-12)
    fig, ax = plt.subplots(figsize=(7.6, 3.6), dpi=100, facecolor='white')
    fig.subplots_adjust(left=0.11, right=0.97, bottom=0.25, top=0.78)
    fig.suptitle('Same Fourier magnitudes, different density profiles',
                 fontsize=13, fontweight='bold', y=0.95)
    ax.plot(x, y0, color='#2166a5', linewidth=2.5, label=r'$3+2\cos(2\pi x)$: phase $0$')
    ax.plot(x, ypi, color='#ba562d', linewidth=2.5, linestyle='--',
            label=r'$3-2\cos(2\pi x)$: phase $\pi$')
    ax.set(xlim=(0, 1), ylim=(0.6, 5.4), xlabel='Fractional position x',
           ylabel='Density (model units)')
    ax.set_xticks([0, .25, .5, .75, 1])
    ax.set_yticks([1, 3, 5])
    ax.grid(alpha=.22)
    ax.spines[['top', 'right']].set_visible(False)
    ax.legend(loc='upper center', bbox_to_anchor=(.5, 1.28), ncol=2,
              frameon=False, fontsize=9)
    fig.text(.5, .065, 'Original one-dimensional model; not measured crystal data',
             ha='center', fontsize=9, color='#4d5962')
    out = Path('electron-density-phase-example.png')
    tmp = out.with_name('electron-density-phase-example.Turing-tmp.png')
    fig.savefig(tmp, facecolor='white')
    plt.close(fig)
    os.replace(tmp, out)
    print(f'{out}: 760 x 360; Fourier-magnitude equality checked')


if __name__ == '__main__':
    main()