Get error distribution#
Created with PtyRAD v1.0.0
Documentation: https://ptyrad.readthedocs.io/en/latest/
PtyRAD paper: https://doi.org/10.1093/mam/ozaf070
PtyRAD arXiv: https://arxiv.org/abs/2505.07814
Zenodo record: https://doi.org/10.5281/zenodo.15273176
Box folder: https://cornell.box.com/s/n5balzf88jixescp9l15ojx7di4xn1uo
Youtube channel: https://www.youtube.com/@ptyrad_official
Before running this notebook, you must first follow the instruction in README.md to:
Create the Python environment with all dependant Python packages like PyTorch
Activate that python environment
Install
ptyradpackage into your activated Python environement (only need to install once)
Note: This notebook is designed for visualization of the reconstruction error.
Author: Chia-Hao Lee, cl2696@cornell.edu
00. Setup working directory and imports#
import os
import numpy as np
import torch
import matplotlib.pyplot as plt
# Change this to the ABSOLUTE PATH to the ptyrad/ folder so you can correctly access data/ and params/
work_dir = "/home/cl2696/altas/workspace/ptyrad"
os.chdir(work_dir)
print("Current working dir: ", os.getcwd())
# The printed working dir should be ".../ptyrad/" to locate the example params files easily
# Note that the output/ directory will be automatically generated under your working directory
Current working dir: /home/cl2696/altas/workspace/ptyrad
from ptyrad.core.models import PtychoModel
from ptyrad.init import Initializer
from ptyrad.io.handlers import save_array
from ptyrad.io.load import load_ptyrad
from ptyrad.params import load_params
from ptyrad.plotting.model import plot_forward_pass
from ptyrad.runtime.device import set_gpu_device
from ptyrad.runtime.diagnostics import print_system_info
from ptyrad.runtime.logging import LoggingManager
from ptyrad.solver.reconstruction import make_batches, select_scan_indices
01. Setup file path and parse params#
model_path = "output/tBL_WSe2/20260606_full_N16384_dp128_flipT100_random32_p6_1obj_6slice_dz2_orblur0.4_ozblur1.0_oposc_sng1.0_spr0.1/model_iter0200.hdf5"
params_path = "output/tBL_WSe2/20260606_full_N16384_dp128_flipT100_random32_p6_1obj_6slice_dz2_orblur0.4_ozblur1.0_oposc_sng1.0_spr0.1/tBL_WSe2.yaml"
ckpt = load_ptyrad(model_path)
params = load_params(params_path, validate=True)
init_params = params.get('init_params')
model_params = params.get('model_params')
recon_params = params.get('recon_params')
# Modify the init_params so we're loading the reconstructed model
init_params['obj_source'] = 'PtyRAD'
init_params['probe_source'] = 'PtyRAD'
init_params['pos_source'] = 'PtyRAD'
init_params['obj_params'] = model_path
init_params['probe_params'] = model_path
init_params['pos_params'] = model_path
02. Initialize logger, device, initializer, and model#
LoggingManager(log_file='ptyrad_log.txt', log_dir='auto', prefix_time='datetime', show_timestamp=True)
print_system_info()
device = set_gpu_device(gpuid=0) # Pass in `gpuid = None` if you don't have access to a CUDA-compatible GPU. Note that running PtyRAD with CPU would be much slower than on GPU.
init = Initializer(init_params).init_all()
model = PtychoModel(init.init_variables, model_params, device=device)
2026-06-06 17:09:55,959 - ### PtyRAD LoggingManager configuration ###
2026-06-06 17:09:55,960 - log_file = 'ptyrad_log.txt'. If log_file = None, no log file will be created.
2026-06-06 17:09:55,960 - log_dir = 'auto'. If log_dir = 'auto', then log will be saved to `output_path` or 'logs/'.
2026-06-06 17:09:55,961 - flush_file = True. Automatically set to True if `log_file is not None`
2026-06-06 17:09:55,961 - prefix_time = datetime. If true, preset strings ('date', 'time', 'datetime'), or a string of time format, a datetime str is prefixed to the `log_file`.
2026-06-06 17:09:55,962 - prefix_jobid = '0'. If not 0, it'll be prefixed to the log file. This is used for hypertune mode with multiple GPUs.
2026-06-06 17:09:55,962 - append_to_file = True. If true, logs will be appended to the existing file. If false, the log file will be overwritten.
2026-06-06 17:09:55,962 - show_timestamp = True. If true, the printed information will contain a timestamp.
2026-06-06 17:09:55,962 -
2026-06-06 17:09:55,963 - ### System information ###
2026-06-06 17:09:55,963 - Platform: Linux-6.8.0-111-generic-x86_64-with-glibc2.35
2026-06-06 17:09:55,963 - Operating System: Linux 6.8.0-111-generic
2026-06-06 17:09:55,963 - OS Version: #111~22.04.1-Ubuntu SMP PREEMPT_DYNAMIC Tue Apr 14 17:13:45 UTC
2026-06-06 17:09:55,964 - Machine: x86_64
2026-06-06 17:09:55,964 - Processor: x86_64
2026-06-06 17:09:55,964 - Available CPU cores: 32
2026-06-06 17:09:55,964 - Total Memory: 62.05 GB
2026-06-06 17:09:55,965 - Available Memory: 48.67 GB
2026-06-06 17:09:55,965 -
2026-06-06 17:09:55,965 - ### GPU information ###
2026-06-06 17:09:55,996 - CUDA Available: True
2026-06-06 17:09:55,996 - CUDA Version: 12.8
2026-06-06 17:09:56,022 - Available CUDA GPUs: ['NVIDIA RTX 5000 Ada Generation', 'NVIDIA RTX 5000 Ada Generation']
2026-06-06 17:09:56,023 - CUDA Compute Capability: ['8.9', '8.9']
2026-06-06 17:09:56,024 - INFO: For torch.compile with Triton, you'll need CUDA GPU with Compute Capability >= 7.0.
2026-06-06 17:09:56,025 - In addition, Triton does not directly support Windows.
2026-06-06 17:09:56,025 - For Windows users, please follow the instruction and download `triton-windows` from https://github.com/woct0rdho/triton-windows.
2026-06-06 17:09:56,057 - MIG (Multi-Instance GPU) mode = False
2026-06-06 17:09:56,058 - INFO: MIG splits a physical GPU into multiple GPU slices, but multiGPU does not support these MIG slices.
2026-06-06 17:09:56,059 - In addition, multiGPU is currently only available on Linux due to the limited NCCL support.
2026-06-06 17:09:56,060 - -> If you're doing normal reconstruction/hypertune, you can safely ignore this.
2026-06-06 17:09:56,060 - -> If you want to do multiGPU, you must provide multiple 'full' GPUs that are not in MIG mode.
2026-06-06 17:09:56,061 -
2026-06-06 17:09:56,061 - ### Python information ###
2026-06-06 17:09:56,061 - Python Executable: /home/cl2696/miniforge3/envs/ptyrad/bin/python
2026-06-06 17:09:56,061 - Python Version: 3.12.12 | packaged by conda-forge | (main, Jan 26 2026, 23:51:32) [GCC 14.3.0]
2026-06-06 17:09:56,061 -
2026-06-06 17:09:56,061 - ### Packages information ###
2026-06-06 17:09:56,067 - Numpy Version (metadata): 2.3.5
2026-06-06 17:09:56,070 - PyTorch Version (metadata): 2.7.1+cu128
2026-06-06 17:09:56,074 - Optuna Version (metadata): 4.7.0
2026-06-06 17:09:56,077 - Accelerate Version (metadata): 1.12.0
2026-06-06 17:09:56,080 - PtyRAD Version (ptyrad/__init__.py): 0.1.0b13.post3
2026-06-06 17:09:56,080 - PtyRAD is located at: /home/cl2696/altas/workspace/ptyrad/src/ptyrad/__init__.py
2026-06-06 17:09:56,081 -
2026-06-06 17:09:56,081 - ### Setting GPU Device ###
2026-06-06 17:09:56,081 - Selected GPU device: cuda:0 (NVIDIA RTX 5000 Ada Generation)
2026-06-06 17:09:56,081 -
2026-06-06 17:09:56,082 - init_params are displayed below:
2026-06-06 17:09:56,082 - random_seed: None
2026-06-06 17:09:56,082 - probe_illum_type: electron
2026-06-06 17:09:56,082 - probe_kv: 80.0
2026-06-06 17:09:56,082 - probe_conv_angle: 24.9
2026-06-06 17:09:56,083 - probe_aberrations: {}
2026-06-06 17:09:56,083 - beam_kev: None
2026-06-06 17:09:56,083 - probe_dRn: None
2026-06-06 17:09:56,083 - probe_Rn: None
2026-06-06 17:09:56,083 - probe_D_H: None
2026-06-06 17:09:56,084 - probe_D_FZP: None
2026-06-06 17:09:56,084 - probe_Ls: None
2026-06-06 17:09:56,084 - meas_Npix: 128
2026-06-06 17:09:56,084 - pos_N_scans: 16384
2026-06-06 17:09:56,084 - pos_N_scan_slow: 128
2026-06-06 17:09:56,085 - pos_N_scan_fast: 128
2026-06-06 17:09:56,086 - pos_scan_step_size: 0.429
2026-06-06 17:09:56,086 - meas_calibration: {'mode': 'fitRBF', 'value': None}
2026-06-06 17:09:56,088 - probe_pmode_max: 6
2026-06-06 17:09:56,088 - probe_pmode_init_pows: [0.02]
2026-06-06 17:09:56,089 - obj_omode_max: 1
2026-06-06 17:09:56,089 - obj_omode_init_occu: {'occu_type': 'uniform', 'init_occu': None}
2026-06-06 17:09:56,089 - obj_Nlayer: 6
2026-06-06 17:09:56,089 - obj_slice_thickness: 2.0
2026-06-06 17:09:56,089 - simu_Npix: None
2026-06-06 17:09:56,089 - simu_match_mode: None
2026-06-06 17:09:56,090 - meas_permute: None
2026-06-06 17:09:56,090 - meas_reshape: None
2026-06-06 17:09:56,090 - meas_flipT: [1, 0, 0]
2026-06-06 17:09:56,090 - meas_crop: None
2026-06-06 17:09:56,091 - meas_pad: None
2026-06-06 17:09:56,091 - meas_resample: None
2026-06-06 17:09:56,091 - meas_add_source_size: None
2026-06-06 17:09:56,091 - meas_add_detector_blur: None
2026-06-06 17:09:56,091 - meas_remove_neg_values: {'mode': 'clip_neg', 'value': None, 'force': False}
2026-06-06 17:09:56,091 - meas_normalization: {'mode': 'max_at_one', 'value': None}
2026-06-06 17:09:56,092 - meas_add_poisson_noise: None
2026-06-06 17:09:56,092 - meas_export: None
2026-06-06 17:09:56,092 - probe_permute: None
2026-06-06 17:09:56,092 - probe_z_shift: None
2026-06-06 17:09:56,092 - probe_normalization: {'mode': 'mean_total_ints', 'value': None}
2026-06-06 17:09:56,093 - pos_scan_flipT: None
2026-06-06 17:09:56,093 - pos_scan_affine: None
2026-06-06 17:09:56,093 - pos_scan_rand_std: 0.15
2026-06-06 17:09:56,093 - obj_z_crop: None
2026-06-06 17:09:56,093 - obj_z_pad: None
2026-06-06 17:09:56,093 - obj_z_resample: None
2026-06-06 17:09:56,093 - meas_source: file
2026-06-06 17:09:56,094 - meas_params: {'path': 'data/tBL_WSe2/Panel_g-h_Themis/scan_x128_y128.raw', 'key': None, 'shape': None, 'offset': None, 'gap': None, 'selection': None, 'zarr_kwargs': None}
2026-06-06 17:09:56,094 - probe_source: PtyRAD
2026-06-06 17:09:56,094 - probe_params: output/tBL_WSe2/20260606_full_N16384_dp128_flipT100_random32_p6_1obj_6slice_dz2_orblur0.4_ozblur1.0_oposc_sng1.0_spr0.1/model_iter0200.hdf5
2026-06-06 17:09:56,094 - pos_source: PtyRAD
2026-06-06 17:09:56,094 - pos_params: output/tBL_WSe2/20260606_full_N16384_dp128_flipT100_random32_p6_1obj_6slice_dz2_orblur0.4_ozblur1.0_oposc_sng1.0_spr0.1/model_iter0200.hdf5
2026-06-06 17:09:56,095 - obj_source: PtyRAD
2026-06-06 17:09:56,095 - obj_params: output/tBL_WSe2/20260606_full_N16384_dp128_flipT100_random32_p6_1obj_6slice_dz2_orblur0.4_ozblur1.0_oposc_sng1.0_spr0.1/model_iter0200.hdf5
2026-06-06 17:09:56,095 - tilt_source: simu
2026-06-06 17:09:56,095 - tilt_params: {'tilt_type': 'all', 'init_tilts': [[0, 0]]}
2026-06-06 17:09:56,095 -
2026-06-06 17:09:56,095 - ### Initializing cache ###
2026-06-06 17:09:56,096 - Loading 'PtyRAD' file from output/tBL_WSe2/20260606_full_N16384_dp128_flipT100_random32_p6_1obj_6slice_dz2_orblur0.4_ozblur1.0_oposc_sng1.0_spr0.1/model_iter0200.hdf5 for caching
2026-06-06 17:09:56,307 - use_cached_obj = True
2026-06-06 17:09:56,309 - use_cached_probe = True
2026-06-06 17:09:56,309 - use_cached_pos = True
2026-06-06 17:09:56,309 -
2026-06-06 17:09:56,310 - ### Initializing measurements ###
2026-06-06 17:09:56,310 - Loading measurements from source = 'file'
2026-06-06 17:09:56,310 - Detected measurement file type = '.raw'
2026-06-06 17:09:56,311 - WARNING: Couldn't find the 'shape' in 'meas_params' with file type = '.raw'.
2026-06-06 17:09:56,311 - It is strongly recommended to provide an explicit shape to better load from .raw files
2026-06-06 17:09:56,312 - PtyRAD will still try to load the dataset based on the provided 'init_params', but you may consider setting 'shape': (N_scans, meas_Npix, meas_Npix) inside your 'meas_params' dict.
2026-06-06 17:10:05,584 - Original measurements dtype is float32, casting to float32 (single precision) for computational efficiency.
2026-06-06 17:10:05,586 - Imported meausrements shape / dtype = (16384, 128, 128), dtype = float32
2026-06-06 17:10:05,846 - Imported meausrements int. statistics (min, mean, max) = (-30.8785, 1814.7078, 99441.8047)
2026-06-06 17:10:05,848 - Flipping measurements with [flipud, fliplr, transpose] = [1, 0, 0]
2026-06-06 17:10:06,056 - Removing negative values in measurement with method = clip_neg and value = None
2026-06-06 17:10:06,136 - Minimum value = -30.8785, negative values are clipped to 0 due to the positive px value constraint of measurements
2026-06-06 17:10:07,226 - Normalizing measurements with mode = 'max_at_one' and value = 'None'
2026-06-06 17:10:07,317 - Normalizing by max of the 2D mean pattern intensity: 67097.391
2026-06-06 17:10:07,397 - meausrements shape / dtype = (16384, 128, 128), dtype = float32
2026-06-06 17:10:07,655 - meausrements int. statistics (min, mean, max) = (0.0000, 0.0271, 1.4821)
2026-06-06 17:10:07,788 - No negative values found in measurements. Skipping non-neg correction.
2026-06-06 17:10:07,982 - Pattern total int. statistics (min, mean, max) = (438.8489, 443.3112, 447.6965), with min/max = 98.0%
2026-06-06 17:10:08,239 - Global meausrements int. statistics (min, mean, max) = (0.0000, 0.0271, 1.4821)
2026-06-06 17:10:08,240 - measurements (N, Ky, Kx) = float32, (16384, 128, 128)
2026-06-06 17:10:08,241 -
2026-06-06 17:10:08,242 - ### Setting up calibration ###
2026-06-06 17:10:08,243 - meas_calibration mode = 'fitRBF', value = None
2026-06-06 17:10:08,244 - Using loaded raw averaged measurement (before crop/pad/resample) to fit RBF as a part of the meas calibration
2026-06-06 17:10:08,244 - Radius of fitted bright field disk (RBF) = 11.64 px with meas_Npix = 128
2026-06-06 17:10:08,244 - Suggested probe_mask_k radius (RBF*2/Npix) > 0.1819
2026-06-06 17:10:08,244 - Fitting raw averaged measurement with center, radius, and Gaussian blur std as a sanity check
2026-06-06 17:10:08,244 - Note that the fitted Gaussian blur std (detector blur) would be affected by overlapping Bragg disks
2026-06-06 17:10:08,325 - Initial guess: center=(63.78, 63.93), radius=11.64, Gaussian blur std=0.50
2026-06-06 17:10:08,334 - Final fit: center=(63.78, 63.93), radius=11.64, Gaussian blur std=0.54
2026-06-06 17:10:08,335 - The 'fitRBF' and 'RBF' calibration methods are highly dependent on the accuracy of user-provided experimental parameters and acquisition conditions,
2026-06-06 17:10:08,336 - including convergence angle, kV, dose, specimen thickness, and collection angle for the estimation of RBF.
2026-06-06 17:10:08,337 - For example, a 5-10% error in convergence angle is fairly common.
2026-06-06 17:10:08,338 - Users are strongly advised to perform proper microscope calibration to ensure accurate results.
2026-06-06 17:10:08,338 - These method should only be used as a rough estimate and not as a substitute for proper experimental calibration.
2026-06-06 17:10:08,338 - dx (real space pixel size of probe and object) set to 0.1526 Ang with Npix = 128
2026-06-06 17:10:08,338 -
2026-06-06 17:10:08,339 - ### Setting init_variables dict ###
2026-06-06 17:10:08,339 - Derived values given input init_params:
2026-06-06 17:10:08,339 - kv = 80.0 kV
2026-06-06 17:10:08,339 - wavelength = 0.0418 Ang
2026-06-06 17:10:08,340 - conv_angle = 24.9 mrad
2026-06-06 17:10:08,341 - Npix = 128 px
2026-06-06 17:10:08,342 - dk = 0.0512 Ang^-1
2026-06-06 17:10:08,342 - kMax = 3.2773 Ang^-1
2026-06-06 17:10:08,342 - da = 2.1383 mrad
2026-06-06 17:10:08,343 - angleMax = 136.8515 mrad
2026-06-06 17:10:08,343 - RBF = 11.6447 px (Inferred from the given calibration, NOT necessarily from the loaded measurement data)
2026-06-06 17:10:08,343 - n_alpha = 5.4960 (# conv_angle)
2026-06-06 17:10:08,343 - dx = 0.1526 Ang, Nyquist-limited dmin = 2*dx = 0.3051 Ang
2026-06-06 17:10:08,344 - Rayleigh-limited resolution = 1.0230 Ang (0.61*lambda/alpha for focused probe )
2026-06-06 17:10:08,344 - Real space probe extent = 19.5282 Ang
2026-06-06 17:10:08,344 -
2026-06-06 17:10:08,344 - ### Initializing probe ###
2026-06-06 17:10:08,345 - Loading probe from source = 'PtyRAD'
2026-06-06 17:10:08,345 - Loaded probe shape = (6, 128, 128), dtype = complex64
2026-06-06 17:10:08,345 - pmode_now: 6 and pmode_max: 6, leaving the pmode unchanged.
2026-06-06 17:10:08,346 - Orthogonalizing 6 pmodes
2026-06-06 17:10:08,349 - Sorting 6 pmodes by their intensities
2026-06-06 17:10:08,350 - Normalizing probe intensity with mode = 'mean_total_ints' and value = 'None'
2026-06-06 17:10:08,351 - sum(|probe_data|**2) = 443.31, while meas_total_ints (min, mean, max) = (438.8489, 443.3112, 447.6965)
2026-06-06 17:10:08,352 - probe (pmode, Ny, Nx) = complex64, (6, 128, 128)
2026-06-06 17:10:08,352 -
2026-06-06 17:10:08,352 - ### Initializing probe positions ###
2026-06-06 17:10:08,353 - Loading probe positions from source = 'PtyRAD'
2026-06-06 17:10:08,353 - Applying Gaussian distributed random displacement with std = 0.15 px to scan positions
2026-06-06 17:10:08,355 - crop_pos (N,2) = int16, (16384, 2)
2026-06-06 17:10:08,355 - crop_pos 1st and last px coords (y,x) = ([48, 49], [404, 404])
2026-06-06 17:10:08,356 - crop_pos extent (Ang) = [54.77040555 55.22809696]
2026-06-06 17:10:08,356 - probe_pos_shifts (N,2) = float32, (16384, 2)
2026-06-06 17:10:08,356 -
2026-06-06 17:10:08,357 - ### Initializing object ###
2026-06-06 17:10:08,357 - Loading object from source = 'PtyRAD'
2026-06-06 17:10:08,428 - omode_now: 1 and omode_max: 1, leaving the omode unchanged.
2026-06-06 17:10:08,434 - object (omode, Nz, Ny, Nx) = complex64, (1, 6, 583, 584)
2026-06-06 17:10:08,435 - object extent (Z, Y, X) (Ang) = [12. 88.9447 89.0973]
2026-06-06 17:10:08,435 -
2026-06-06 17:10:08,435 - ### Initializing omode_occu from 'uniform' ###
2026-06-06 17:10:08,436 - omode_occu (omode) = float32, (1,)
2026-06-06 17:10:08,436 -
2026-06-06 17:10:08,436 - ### Initializing H (Fresnel propagator) ###
2026-06-06 17:10:08,436 - Calculating H with probe_shape = (128, 128) px, dx = 0.1526 Ang, slice_thickness = 2.0000 Ang, lambd = 0.0418 Ang
2026-06-06 17:10:08,437 - H (Ky, Kx) = complex64, (128, 128)
2026-06-06 17:10:08,437 -
2026-06-06 17:10:08,437 - ### Initializing obj tilts from = 'simu' ###
2026-06-06 17:10:08,438 - Initialized obj_tilts with init_tilts = [[0, 0]] (theta_y, theta_x) mrad
2026-06-06 17:10:08,438 - obj_tilts (N, 2) = float32, (1, 2)
2026-06-06 17:10:08,438 -
2026-06-06 17:10:08,438 - ### Checking consistency between input params with the initialized variables ###
2026-06-06 17:10:08,439 - meas_Npix, simu_Npix, DP measurements, probe, and H shapes are consistent as '128'
2026-06-06 17:10:08,439 - N_scans, len(meas), N_scan_slow*N_scan_fast, len(crop_pos), and len(probe_pos_shifts) are consistent as '16384'
2026-06-06 17:10:08,439 - obj.shape[0] is consistent with len(omode_occu) as '1'
2026-06-06 17:10:08,439 - obj.shape[1] is consistent with Nlayer as '6'
2026-06-06 17:10:08,440 - crop positions (yx_min=[48 46], yx_max=[535 536]) are well contained inside object canvas (Ny,Nx) = (583, 584).
2026-06-06 17:10:08,441 - obj_tilts is consistent with either 1 or N_scans
2026-06-06 17:10:08,441 - Pass the consistency check of initialized variables, initialization is done!
2026-06-06 17:10:08,441 -
2026-06-06 17:10:08,441 - ### Collecting reconstruction provenance ###
2026-06-06 17:10:08,459 - Reconstruction provenance is collected and initialized.
2026-06-06 17:10:08,460 - ### Initializing PtychoModel model ###
2026-06-06 17:10:08,949 - ### PtychoModel optimizable variables ###
2026-06-06 17:10:08,950 - obja : torch.Size([1, 6, 583, 584]) , torch.float32 , device:cuda:0, grad:True , lr:5e-04
2026-06-06 17:10:08,951 - objp : torch.Size([1, 6, 583, 584]) , torch.float32 , device:cuda:0, grad:True , lr:5e-04
2026-06-06 17:10:08,951 - obj_tilts : torch.Size([1, 2]) , torch.float32 , device:cuda:0, grad:False, lr:0e+00
2026-06-06 17:10:08,951 - slice_thickness : torch.Size([]) , torch.float32 , device:cuda:0, grad:False, lr:0e+00
2026-06-06 17:10:08,952 - probe : torch.Size([6, 128, 128, 2]) , torch.float32 , device:cuda:0, grad:True , lr:1e-04
2026-06-06 17:10:08,952 - probe_pos_shifts: torch.Size([16384, 2]) , torch.float32 , device:cuda:0, grad:True , lr:5e-04
2026-06-06 17:10:08,952 -
2026-06-06 17:10:08,953 - ### Optimizable variables statitsics ###
2026-06-06 17:10:08,953 - Total measurement values : 268,435,456
2026-06-06 17:10:08,953 - Total optimizing variables: 4,315,040
2026-06-06 17:10:08,953 - Overdetermined ratio : 62.21
2026-06-06 17:10:08,954 -
2026-06-06 17:10:08,954 - ### Model behavior ###
2026-06-06 17:10:08,955 - Tilt propagator : False
2026-06-06 17:10:08,955 - Change slice thickness : False
2026-06-06 17:10:08,955 - Detector blur : False
2026-06-06 17:10:08,956 - Preload data : True
2026-06-06 17:10:08,956 - On-the-fly meas padding : False
2026-06-06 17:10:08,956 - On-the-fly meas resample : False
2026-06-06 17:10:08,957 - On-the-fly simu match mode: None
2026-06-06 17:10:08,957 -
2026-06-06 17:10:08,998 - ### Done initializing PtychoModel model ###
2026-06-06 17:10:09,000 -
# Run a forward pass to check if the reconstructed model (obj, probe, pos) was correctly loaded
indices = np.random.randint(0, init.init_variables["N_scans"], 2)
dp_power = 0.5
plot_forward_pass(model, indices, dp_power)
03. Run 1 full foward pass for the error distribution#
INDICES_MODE = recon_params.get('INDICES_MODE')
batch_config = recon_params.get('BATCH_SIZE', {})
grad_accumulation = batch_config.get("grad_accumulation", 1)
batch_size = batch_config.get('size') * grad_accumulation
GROUP_MODE = recon_params.get('GROUP_MODE')
pos = (model.crop_pos + model.opt_probe_pos_shifts).detach().cpu().numpy()
indices = select_scan_indices(
init.init_variables['N_scan_slow'],
init.init_variables['N_scan_fast'],
subscan_slow=INDICES_MODE.get('subscan_slow'),
subscan_fast=INDICES_MODE.get('subscan_fast'),
mode=INDICES_MODE.get('mode', 'random'),
)
batches = make_batches(indices, pos, batch_size, mode=GROUP_MODE)
2026-06-06 17:10:10,084 - Selecting indices with the 'full' mode
2026-06-06 17:10:10,481 - Generated 512 'random' groups of ~32 scan positions in 0.001 sec
# Initialize error_DP to hold the relative error patterns, note that this would take the same amount of GPU VRAM as the 4D-dataset
Npix = model.get_complex_probe_view().shape[-1]
error_DP = torch.zeros((len(indices), Npix, Npix), device=device, dtype=torch.float32)
dp_pow = 0.5
for batch in batches:
with torch.no_grad():
model_DP = model(batch)
measured_DP = model.get_measurements(batch)
data_mean = measured_DP.pow(dp_pow).mean(0).clamp(min=1e-8)
error_DP[batch] = (model_DP.pow(dp_pow) - measured_DP.pow(dp_pow)) / data_mean # error_DP corresponds to the relative error
error_DP_arr = error_DP.detach().cpu().numpy()
# Saving is optional but you can port it into py4DGUI for interactive visualization
save_array(error_DP_arr, file_dir='./', file_name='relative_error', file_format="hdf5", output_shape=None, append_shape=True,)
2026-06-06 17:10:11,711 - Saving array with shape = (16384, 128, 128) and dtype = float32
2026-06-06 17:10:11,714 - file path = './relative_error_16384_128_128.hdf5' already exists, the file will be overwritten.
N_scan_slow = init.init_variables['N_scan_slow']
N_scan_fast = init.init_variables['N_scan_fast']
fig, axs = plt.subplots(nrows=1,ncols=2, figsize = (12,6))
axs[0].set_title('Fourier Space Error Distribution')
im0 = axs[0].imshow(error_DP_arr.mean(0), cmap='seismic', vmin=-0.2, vmax=0.2)
fig.colorbar(im0, shrink=0.6)
axs[1].set_title('Real Space Error Distribution')
im1 = axs[1].imshow(error_DP_arr.mean((-1,-2)).reshape(N_scan_slow, N_scan_fast), cmap='seismic', vmin=-0.2, vmax=0.2)
fig.colorbar(im1, shrink=0.6)
plt.show()
Note#
Since it’s relative error, it naturally contains positive and negative values for each error_DP (N, ky, kx)
If we take the average along N, then inevitably some k pixels would have their errors cancel out.
Similarly, when we average along (ky, kx), errors for all k pixels would be averaged together and show an overall trend for each probe position.
Fourier Space Error Distribution#
Low k region (center) shows low relative error, which is more or less expected as it contains most of the signal and generally SNR ~ 1/sqrt(signal)
High k region has higher relative error, indicating model DP predicts stronger intensity than experiment pattern.
Real Space Error Distribution#
On atomic sites the relative error tend to be negative values, indicating model DP predicts overall lower intensity than experiment pattern.